24小时热门版块排行榜    

查看: 2119  |  回复: 5
当前只显示满足指定条件的回帖,点击这里查看本话题的所有回帖

eudemonia333

新虫 (初入文坛)

[求助] 病态矩阵求解

病态线性方程组阵Ax=B,A为19*19的矩阵,条件数E11,行列式2.3E90。
(1)是不是这类病态矩阵直接采用A\B或者inv(A)*B结果偏差都比较大?试过用符号解再转换为数值和直接inv(A)*B,发现二者结果误差很小。对此问题比较疑惑。
(2)矩阵较大,不好贴出,A和B均在附件中,如果有其它高精度的求解方法(如tikhonov正则化方法,在网上下载了tikhonov源代码,采用后结果和A\B相比增大了1000倍,不知道是我自己写的调用代码有问题还是其它问题)求大神赐教。

[ Last edited by eudemonia333 on 2013-4-22 at 21:03 ]
回复此楼

» 本帖附件资源列表

  • 欢迎监督和反馈:小木虫仅提供交流平台,不对该内容负责。
    本内容由用户自主发布,如果其内容涉及到知识产权问题,其责任在于用户本人,如对版权有异议,请联系邮箱:xiaomuchong@tal.com
  • 附件 1 : A.xls
  • 2013-04-22 20:58:46, 16.5 K
  • 附件 2 : B.xls
  • 2013-04-22 20:58:49, 13.5 K

» 收录本帖的淘帖专辑推荐

精彩

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

chunbendan

铜虫 (初入文坛)

【答案】应助回帖

寒,怎么变成表情了哦。 就是一个:加一个)
6楼2013-04-27 11:55:12
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 6 个回答

chunbendan

铜虫 (初入文坛)

【答案】应助回帖


lovibond: 金币+1, 鼓励新虫应助 2013-04-27 12:12:00
病态矩阵求解需要特别注意。matlab里直接用inv(A)*B肯定是不行的,误差超级大。用A\B不知道有没有内部进行处理,也不建议用。建议自己写算法,奇异值修正或者是迭代法都可以。我这里写两个简单的迭代法程序供楼主参考:

Jacobi迭代的函数文件:XJGJacobi.m
function [x1,n]=XJGJacobi(A,b,x0,eps)
%   [x1,n]=XJGJacobi(A,b,x0,eps)使用JACOBI迭代法进行迭代求解方程A*x=b
%   A,b     输入量,方程的矩阵和向量
%   x0      输入量,迭代初始值
%   eps     精度,比如0.0001
%   x1      输出,结果
%   n       输出,迭代次数

D=diag(diag(A));
B=D\(D-A);
n=1;
% 初始迭代
x1=B*x0+D\b;
while norm(x1-x0,inf)>eps   %终止条件
n=n+1;
x0=x1;
x1=B*x0+D\b;
end

GAUSS_SEIDIAL迭代的函数文件:XJGGauss.m
function [x1,n]=XJGGauss(A,b,x0,eps)
%   [x1,n]=XJGGauss(A,b,x0,eps)
%   本函数使用GAUSS-SEDIAL迭代法进行迭代求解线性方程A*x=b
%   A,b     方程中的矩阵和向量
%   x0      初始迭代值
%   eps     精度

% 2011年12月5日 18:47:56

LD=tril(A,0);
U=-triu(A,1);
B=LD\U;
f=LD\b;
x1=B*x0+f;
n=1;    %迭代次数

while norm(x1-x0,inf)>eps   %终止条件
n=n+1;
x0=x1;
x1=B*x0+f;
end
return

%我随便写了A和B,你的矩阵太大了
%用于求解的脚本文件:
A=[6,2,1,-2
2,5,0,-2
-2,0,8,5
1,3,2,7];
b=[4;7;-1;0];
x0=zeros(4,1);  %初始值
eps=1e-5;       %精度
% Jacobi迭代法求解并显示结果
[x1,n]=XJGJacobi(A,b,x0,eps);
disp('Jacobi迭代的结果为:')
disp(x1)
disp('Jacobi迭代的次数为:')
disp(n)
% GAUSS-SEDIAL迭代法求解并显示结果
[x1,n]=XJGGauss(A,b,x0,eps);
disp('Gauss-Sedial迭代的结果为:')
disp(x1)
disp('Guass-Sedial迭代次数为:')
disp(n)
2楼2013-04-26 09:27:39
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

eudemonia333

新虫 (初入文坛)

引用回帖:
2楼: Originally posted by chunbendan at 2013-04-26 09:27:39
病态矩阵求解需要特别注意。matlab里直接用inv(A)*B肯定是不行的,误差超级大。用A\B不知道有没有内部进行处理,也不建议用。建议自己写算法,奇异值修正或者是迭代法都可以。我这里写两个简单的迭代法程序供楼主参 ...

谢谢!
上面两个程序不太适应我的矩阵,计算完出现了NaN,Warning: Matrix is singular to working precision. 直接采用inv没有上面两个错误。
3楼2013-04-26 16:49:33
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

chunbendan

铜虫 (初入文坛)

【答案】应助回帖

引用回帖:
3楼: Originally posted by eudemonia333 at 2013-04-26 16:49:33
谢谢!
上面两个程序不太适应我的矩阵,计算完出现了NaN,Warning: Matrix is singular to working precision. 直接采用inv没有上面两个错误。...

我没用你的数据测试,应该修改一下能用,或者用列主元消去法也可以。我明天帮你弄一下
用inv处理奇异矩阵不行的

[ 发自手机版 http://muchong.com/3g ]
4楼2013-04-26 23:24:09
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
信息提示
请填处理意见