|
|
[求助]
小弟求教几个程序-关于迭代法的 已有2人参与
小弟编写了Table4.1SOR中的r=0.01及n=8的程序,现有如下几个疑问:1、我运行出it为8,cpu=0.042,w(exp)=1.02,为什么与作者不一样?2、作者终止条件是F范数,我的是2范数,应该怎么修改?3.当n=64时出现‘Out of memory. Type HELP MEMORY for your options.’请问是我程序太复杂了吗?
注:test为确定w(exp)程序,SOR为调用程序,test2为根据w(exp)计算it程序,各位大神不吝赐教。
test.m:
n=8;
r=0.01;
a=ones(n,1);
b=ones(n-1,1);
I=eye(n);
MM=2*diag(a)-diag(b,-1)-diag(b,1);
N=0.5*diag(b,-1)-0.5*diag(b,1);
A=MM+2*r*N+100*I/(n+1);
B=A;
X=rand(n,n);
F=A*X+X*B;
E=kron(I,A)+kron(B',I);
b=F( ;
%x=E\b;
x0=b;
ww=0.01:0.01:1.99;
D=diag(diag(E)); %求E的对角矩阵
L=-tril(E,-1); %求E的下三角矩阵
U=-triu(E,1); %求E的上三角矩阵
for i=1:length(ww)
w=ww(i);
B=inv(D-L*w)*((1-w)*D+w*U);
p(i)=max(max(abs(eig(B))));
%plot(ww,p);
[x,n1]=SOR(E,b,x0,w);
mm(i)=n1;
end
plot(ww,mm);
Xe=reshape(x,n,n);
errror=norm((Xe-X),2)
SOR.m:
function [x,n1]=SOR(E,b,x0,w,eps,M)
if nargin==4
eps=1.0e-6;
M=2000;
elseif nargin<4
error
return
elseif nargin==5
M=200;
end
if(w<=0 || w>=2) %收敛条件要求
error;
return;
end
D=diag(diag(E)); %求E的对角矩阵
L=-tril(E,-1); %求E的下三角矩阵
U=-triu(E,1); %求E的上三角矩阵
B=inv(D-L*w)*((1-w)*D+w*U);
f=w*inv((D-L*w))*b;
x=B*x0+f;
n1=1; %迭代次数
%迭代过程
while norm(x-x0)>=eps
x0=x;
x=B*x0+f;
n1=n1+1;
if(n1>=M)
disp('Warining:迭代次数太多,可能不收敛!');
return;
end
end
test2.m
w=1.02
n=8;
r=0.01;
a=ones(n,1);
b=ones(n-1,1);
I=eye(n);
M=2*diag(a)-diag(b,-1)-diag(b,1);
N=0.5*diag(b,-1)-0.5*diag(b,1);
A=M+2*r*N+100*I/(n+1);
B=A;
X=rand(n,n);
F=A*X+X*B;
E=kron(I,A)+kron(B',I);
b=F( ;
%x=E\b;
x0=b;
[x,n1]=SOR(E,b,x0,w)
Xe=reshape(x,n,n);
errror=norm((Xe-X),2)
![小弟求教几个程序-关于迭代法的]()
1.png
![小弟求教几个程序-关于迭代法的-1]()
2.png |
» 猜你喜欢
十年后又回来了,论文投稿求助
已经有5人回复
售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急
已经有6人回复
售SCI一区T0P文章,我:8.O.55.1.O.54,科目齐全,可+急
已经有4人回复
fileCode有新解读?
已经有17人回复
这样的filecode谁见过
已经有8人回复
售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急
已经有6人回复
售SCI一区T0P文章,我:8O.55.1.O.54,科目全,可伽急
已经有6人回复
售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急
已经有6人回复
售SCI一区T0P文章,我:8.O.55.1.O54,科目全,可伽急
已经有5人回复
售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急
已经有4人回复
» 本主题相关价值贴推荐,对您同样有帮助:
|