> y <- c(-44.394, -45.3049, -44.9594)
> a <- c(50.0, 43.0, 52.0)
> b <- c(37.0, 39.0, 41.0)
> regression = lm(y ~ a + b)
> regression
Call:
lm(formula = y ~ a + b)
Coefficients:
(Intercept) a b
-41.63759 0.07852 -0.18061
function x = LinSolve(A,y)
%
% Recursive Solution of Linear System Ax=y
% matlab equivalent: x = A\y
% x = n x 1
% A = n x n
% y = n x 1
% Uses stack space extensively. Not efficient.
% C allows recursion, so convert it into C.
% ----------------------------------------------
n=length(y);
x=zeros(n,1);
if(n>1)
x(1:n-1,1) = LinSolve( A(1:n-1,1:n-1) - (A(1:n-1,n)*A(n,1:n-1))./A(n,n) , ...
y(1:n-1,1) - A(1:n-1,n).*(y(n,1)/A(n,n)));
x(n,1) = (y(n,1) - A(n,1:n-1)*x(1:n-1,1))./A(n,n);
else
x = y(1,1) / A(1,1);
end
10 回答
Cramer's Rule和Gaussian Elimination是两个很好的通用算法(另见Simultaneous Linear Equations) . 如果您正在寻找代码,请查看GiNaC,Maxima和SymbolicC++(当然,取决于您的许可要求) .
编辑:我知道你在C领域工作,但我也必须为SymPy(Python中的计算机代数系统)提供一个好词 . 你可以从它的算法中学到很多东西(如果你能读懂一些python) . 此外,它是在新的BSD许可下,而大多数免费数学包是GPL .
您可以使用与手动解决方法完全相同的程序来解决此问题(使用乘法和减法,然后将结果反馈到方程式中) . 这是非常标准的中学数学 .
所以你最终得到:
如果将这些值插回A,B和C,您会发现它们是正确的 .
诀窍是使用一个简单的4x3矩阵,它依次减少到3x2矩阵,然后是2x1,即“a = n”,n是实际数字 . 一旦你有了它,你将它提供到下一个矩阵中以获得另一个值,然后将这两个值放入下一个矩阵,直到你解决了所有变量 .
如果你有N个不同的方程式,你总是可以求解N个变量 . 我说的很明显,因为这两个不是:
它们是相同的等式乘以2,所以你无法从它们得到一个解决方案 - 将第一个乘以2然后减去离开你的真实但无用的声明:
举个例子,这里's some C code that works out the simultaneous equations that you'重新放在你的问题中 . 首先是一些必要的类型,变量,用于打印方程的支持函数,以及
main
的开头:接下来,将具有三个未知数的三个方程式减少到具有两个未知数的两个方程:
接下来,将具有两个未知数的两个方程式减少到一个未知的方程式:
现在我们有一个
number1 = unknown * number2
类型的公式,我们可以使用unknown <- number1 / number2
简单地计算出未知值 . 然后,一旦你计算出该值,将其替换为具有两个未知数的方程之一并计算出第二个值 . 然后将这些(现在已知的)未知数替换为原始方程之一,现在您拥有所有三个未知数的值:该代码的输出与此答案中的早期计算相匹配:
对于3x3线性方程组,我想可以推出自己的算法 .
但是,您可能不得不担心准确性,除以零或非常小的数字以及如何处理无限多的解决方案 . 我的建议是使用标准的数字线性代数包,如LAPACK .
看看Microsoft Solver Foundation .
有了它,你可以写这样的代码:
Here is the output:
=== Solver Foundation Service Report ===
日期时间:2009年4月20日23:29:55
型号名称:默认
要求的能力:LP
求解时间(ms):1027
总时间(毫秒):1414
解决完成状态:最佳
选择的解算器:Microsoft.SolverFoundation.Solvers.SimplexSolver
指令:
Microsoft.SolverFoundation.Services.Directive
算法:原始
算术:混合
定价(准确):默认
定价(双倍):SteepestEdge
基础:松弛
枢轴数:3
===解决方案详情===
目标:
决定:
a:0.0785250000000004
b:-0.180612500000001
c:-41.6375875
您是在寻找能够完成工作还是实际执行矩阵操作的软件包,并且每个步骤都执行?
第一个,我的一个同事刚刚使用Ocaml GLPK . 它只是GLPK的包装器,但它消除了很多设置步骤 . 看起来你不得不坚持使用C语言中的GLPK . 对于后者,感谢美味的保存旧文章,我曾经学过LP一段时间,PDF . 如果您需要进一步设置特定帮助,请告诉我们,我可以从这里直接进行操作 . 祝好运!
来自NIST的Template Numerical Toolkit有这样做的工具 .
一种更可靠的方法是使用QR Decomposition .
这是一个包装器的例子,所以我可以在我的代码中调用“GetInverse(A,InvA)”,它会将反转INVA .
Array2D在库中定义 .
根据你的问题的措辞,你似乎有更多的方程而不是未知数,你想要最小化不一致性 . 这通常通过线性回归来完成,其最小化不一致的平方和 . 根据数据的大小,您可以在电子表格或统计包中执行此操作 . R是一个高质量的免费软件包,可以在很多其他方面进行线性回归 . 线性回归(以及很多问题)有很多,但是对于简单的情况来说这很简单 . 这是使用您的数据的R示例 . 请注意,“tx”是模型的截距 .
在运行时效率方面,其他人的回答都比I好 . 如果你总是会有与变量相同数量的方程式,我喜欢Cramer's rule,因为它已经写好了,我相信你能找到一个),并划分两个矩阵的决定因素 .
就个人而言,我偏爱Numerical Recipes的算法 . (我很喜欢C版 . )
本书将教您算法的工作原理,并向您展示这些算法的一些非常好的调试实现 .
当然,你可以盲目地使用CLAPACK(我已经非常成功地使用了它),但我首先手动输入高斯消除算法,至少对于制作这些算法所做的工作有一个微弱的想法 . 稳定 .
之后,如果你正在做更有趣的线性代数,那么查看Octave的源代码将回答很多问题 .