关于 matlab 里的矩阵除法问题
这个问题简单来说就是求解不定或超定方程A*X=B,A为m*n的矩阵在matlab里直接用X=B\A就可以了
可我现在需要用C语言实现
所以请问大家matlab中是如何实现B\A的
或者在matlab 里的“ / ”除运算符调用的是哪个内部函数?
小弟不胜感激:@)
[ 本帖最后由 eight 于 2008-3-25 19:50 编辑 ] 解这样的方程,是不是要用到数值计算中的迭代方法阿? 在matlab中矩阵的左除"\"是一种自适应算法
对于有唯一解的线性方程组求出精确解
对于超定方程线性方程组求出最小二乘解
对于欠定线性方程组求其基本解
这个好像是个内部函数,看不了它的代码 原帖由 happy 于 2006-8-23 15:22 发表
在matlab中矩阵的左除"\"是一种自适应算法
对于有唯一解的线性方程组求出精确解
对于超定方程线性方程组求出最小二乘解
对于欠定线性方程组求其基本解
这个好像是个内部函数,看不了它的代码
首选感谢教授的热心指点:@)
我现在面临的问题是欠定线性方程组的求解,即A为M*(M+1)的长方阵,X为(M+1)*3的矩阵,B为M*3的矩阵
请问教授在Matlab 里是如何求此基本解的?
我最近也看了不少广义逆矩阵的资料,有A-A+A(1,3)A(1,4)等等,拭了几种都跟matlab 算的不一样,另外还有用QR分解的方法,现在正在做,不知道行不行
下面是matlab 里的mldivide文件,不知道是否是调用的这个函数来求"\"运行符,请指教
function X = mldivide(A, B)
%MLDIVIDE Symbolic matrix left division.
% MLDIVIDE(A,B) overloads symbolic A \ B.
% X = A\B solves the symbolic linear equations A*X = B.
% Warning messages are produced if X does not exist or is not unique.
% Rectangular matrices A are allowed, but the equations must be
% consistent; a least squares solution is not computed.
% Copyright 1993-2003 The MathWorks, Inc.
% $Revision: 1.18.4.2 $$Date: 2004/04/16 22:22:54 $
A = sym(A);
B = sym(B);
if all(size(A) == 1)
% Division by a scalar
X = ldivide(A,B);
elseif ndims(A) > 2 | ndims(B) > 2
error('symbolic:sym:mldivide:errmsg1','Input arguments must be 2-D.')
elseif size(A,1) ~= size(B,1)
error('symbolic:sym:mldivide:errmsg2','First dimensions must agree.')
else
% Matrix divided by matrix
X = maple('linsolve',char(A),char(B),'''_rank''');
% Solution does not exist.
if isempty(X)
warning('symbolic:sym:mldivide:warnmsg1','System is inconsistent. Solution does not exist.')
X = Inf;
X = sym(X(ones(size(A,2),size(B,2))));
maple('_rank := ''_rank'';');
return
end;
% Check rank and clear _rank in Maple workspace.
if str2double(maple('_rank')) < min(size(A))
warning('symbolic:sym:mldivide:warnmsg2','System is rank deficient. Solution is not unique.')
end
maple('_rank := ''_rank'';');
% Set any free parameters, _t, to zero.
t = findstr(X,'_t[');
s = findstr(X,']');
for k = fliplr(t)
r = s(s > k);
X(k:r(2)) = '0';
end
X = maple('',sym(X));
end I'm not sure...It seems like this???
b=round(10*rand(3));
A=round(10*rand(3,4));
X=pinv(A)*b 对于欠定方程组matlab一般采用两种方法求解
1. 采用除法求解,这种方法的求解结果是具有最多零元素的解
2. 基于伪逆pinv的求解方法,这种方法求解的具有最小长度或范数的解,楼上所说的属于这种情况
举一个简单的例子,比如
a=;b=;
则采用两种方法求解,x=a\b 和 x=pinv(a)b
第一种方法的结果是:
第二种方法的结果是: mldivide这个函数一直没有注意过,不过从说明上看好像这个函数就是左除 原帖由 happy 于 2006-8-23 19:10 发表
对于欠定方程组matlab一般采用两种方法求解
1. 采用除法求解,这种方法的求解结果是具有最多零元素的解
的确如此,比如取M=325,求解出来的X大多数行都是0,如以下截取的前13行,为表示方便,转置了一下
0 0-1.2593e+0090 0 0 0 0 0 0 0 0 -1.3489e+009
0 0-6.8914e+0090 0 0 0 0 0 0 0 0 -7.1698e+009
0 0-1.7833e+0100 0 0 0 0 0 0 0 0 -1.7629e+010
请教授大概讲一下计算原理,我得用C语言来实现
另外,今天我拭验了一下,似乎是采取如下的方法,下面是一段例程,可以直接运行,拷贝到matlab里看着清楚一些
clc
clear
m=15;
a=rand(m,m+1);
b=rand(m,3);
%相当于求解欠定线性方程组a*xx=b
=qr(a);%对a进行QR分解,即a=q*r,q为m*m的正交阵(重要性质:q*q'=I),r为m*(m+1)的上三角矩阵
rr=r;%保存r
for k=1:m+1
r=rr;
r(:,k)=[];%去掉r阵的第k列,r成为方阵,从而可以用inv求逆
X{1,k}=inv(r)*q'*b;
end
xx=a\b;%此为由matlab内部函数求出的方程解
%说明
%若xx矩阵中k行全为0,即表明去掉了r阵的k行,亦即 X{1,k}中的矩阵即为xx 请教授明查
%
%所以现在的问题就是不知道matlab是依据什么规则去掉r阵的某一行,
再次感谢教授的指导:@)
俺是研二学生,现在课题马上就要作完了,就卡在这一个问题上了,恳请教授不吝赐教
[ 本帖最后由 wwdjkp 于 2006-8-24 18:55 编辑 ] 今天又查了一下相关资料,可以确定的是mldivide就是左除,右除是mrdivide
其基本算法采用的QR因式分解,超定方程的整个求解方法你可以查看qr这个命令的帮助文件Example 2.
欠定的类似的,你看一下相关数值分析的书籍吧
受益匪浅
受益匪浅:handshake
页:
[1]