数学建模社区-数学中国
标题:
求助：关于分块矩阵的还原，为什么实行不了
[打印本页]
作者:
舒米牛牛
时间:
2009-12-17 20:17
标题:
求助：关于分块矩阵的还原，为什么实行不了
%运用Jacobi迭代拟合出来的u关于x和y的矩阵
! T4 {6 j x5 K8 `
function J=jac(A,b,u0,eps)
" D6 ^4 u1 q( Q) ~+ Z" l( N/ k3 z
if nargin==3
$ q( ?+ e/ d+ _
eps=1.0e-8
+ M4 G1 d6 J: L0 j5 f9 ^, D1 k% a
elseif nargin<3
. o6 d- o f) s' Y
'error'
: A/ v' j( K9 G
return
2 m! ]4 ?' R: }: g0 u
end
/ F( U' R( h3 ~+ z, I
, m" `7 G0 L. m" n# t5 b" i2 r5 Z
%定义内部节点矩阵u0
5 E5 s- Z1 D- p! |9 _7 J
h=1;k=1;
, h2 Y9 o5 P6 F1 ^
x=0:h:17;y=0:k:10;
0 k& W1 M- Q/ p4 |: U
e=length(x)-2;f=length(y)-2;
* U0 \: J2 r# m
u=zeros(e,f);
' R P- m% u; U
u0=u;
9 i& n' P/ `3 k0 Y5 q
3 V, G: X' P8 y/ ~
%定义外部节点p
, M X8 c& T$ |
p=zeros(e+2,e+2);
1 v- u% ?+ g( _' K
p(1,1:f+2)=0;p(e+2,1:f+2)=0;
7 l& K% i3 c" S p
p(1:e+2,1)=100;p(1:e+2,f+2)=100;
" Y- F9 z' W# n1 b
4 O7 j; Z8 s Q( E1 K* s1 t# j
%定义系数矩阵A
: Z; l+ I. a) q. o
A=zeros(e*f,e*f);
# r8 j ^5 Y7 E+ i8 K7 d. _
B=mat2cell(A,ones(e*f/e,1)*e,ones(e*f/e,1)*e);
' O7 U6 p6 C, z" v
d1=ones(e,1);d2=ones(e-1,1);
- r! E2 {. L% t8 v/ M" I% K
M=4*diag(d1)-diag(d2,1)-diag(d2,-1);
' g, ^3 S2 r3 V" O9 n N9 ?# l9 V4 T
N=-eye(e);
9 I' ], [: r/ ^4 z# ?" H0 I
B{1,1}=M;B{e,e}=M
7 z: ^; I W. n7 K
for i=2:e-1
/ p& J" u' P1 H+ X
B{i,i}=M;B{i-1,i}=N;B{i+1,i}=N
6 X7 c( b$ ` ^3 D4 S8 _" `, c$ W' X8 B
end
2 ~7 l5 t$ F% K* j' M
A=cell2mat(B);
2 p. h6 z$ I/ X; K
这里总是显示
1 Z1 r! V r, B; ~- ~; V! S4 A& j
??? function J=jac(A,b,u0,eps)
5 P( d/ O& ]# ]" [
|
$ r+ D+ q3 e& l0 n& V
Error: Function definitions are not permitted at the prompt or in scripts.
& S: Y) o, y ?+ `9 ]" l8 i
% v5 F& s+ C1 L. e1 W
%定义b
) U+ c. |1 i9 w3 ~. M) u$ d
b=ones(e*f,1);
( W7 v& s) L" i& T0 j1 m& w
for i=1:e
, _, A6 J- ~; J* M
for j=1:f
- G4 e/ I/ d* y# r) G; `9 ]1 d; J5 `
b(i+j)=p(i,j+1)+p(i+1,j)+p(i+2,j+1)+p(i+1,j+2)
/ n1 ?0 J# V: E
end
( a) B' ?9 f" D. j I0 p: E- i
end
' n) r9 t: s$ {$ ]% ~
%运用Jacobi迭代法计算
) E) q! |" _! N4 O3 m9 a+ h C9 @
D=diag(diag(A));
9 Y, B4 c: |) W# d
D=inv(D);
) E' S6 ]9 d4 @+ f
L=tril(A,-1);
) ]) t! |1 l+ K. G7 I7 Y( Q; O
U=triu(A,1);
% U( ^' f, e, O) e
B=-D*(L+U);
, _7 ]; W: }1 }& E* ^
f=D*b;
$ H- I3 L; B! Q i( G
J=B*u0+f;
: L- h, S* G: I! p- t/ i% Y
while norm(J-u0)>=eps
0 [5 g6 H# o7 T
x0=J;
, U* s c- v& I8 e
J=B*u0+f;
& x* v y1 |& C0 h' t* X
end
1 G& ^* W) F) W% g' Z
return
作者:
舒米牛牛
时间:
2009-12-17 20:18
自己顶一下，拜托哪位高手指点一下，纠结这个矩阵的还原，想了好多方法还是不行~~实在想不出哪里出错了
作者:
BenCam
时间:
2009-12-17 21:06
对不起,我也不知道,帮不了忙!
作者:
madio
时间:
2009-12-17 22:38
你是不是函数的定义没有放在M文件中？
欢迎光临 数学建模社区-数学中国 (http://www.madio.net/)
Powered by Discuz! X2.5