数学建模社区-数学中国
标题:
求助：关于分块矩阵的还原，为什么实行不了
[打印本页]
作者:
舒米牛牛
时间:
2009-12-17 20:17
标题:
求助：关于分块矩阵的还原，为什么实行不了
%运用Jacobi迭代拟合出来的u关于x和y的矩阵
9 g" v+ X2 F( f/ @4 D
function J=jac(A,b,u0,eps)
; u9 F0 S0 T2 H4 z/ @
if nargin==3
; w! S* ^& ~8 z# x' H- t3 s
eps=1.0e-8
8 L9 d1 j& S# m2 `& Y
elseif nargin<3
$ c1 ~/ O! D0 b. y$ d, e5 N
'error'
2 D$ t4 |5 C- |! \) v$ H" b u
return
) v+ j7 U& g) t- e. L
end
* D X7 X4 [4 v: |3 K$ o: v
{2 d$ c% e& `* z
%定义内部节点矩阵u0
) G; t% n& Y* Z$ [( U
h=1;k=1;
. A! e( T+ r! K B$ I S4 E& Z: N% C
x=0:h:17;y=0:k:10;
8 o" `+ z: G) t$ p
e=length(x)-2;f=length(y)-2;
8 g& D, N" v1 \% F
u=zeros(e,f);
/ ^% w) v$ }3 j! F" j
u0=u;
/ [2 \$ Q2 h! {
) `/ R: k& u, g" X8 N* S( p
%定义外部节点p
9 R# u T% r7 b! M& H
p=zeros(e+2,e+2);
; @6 m( [5 [% q/ M, ~% l" v
p(1,1:f+2)=0;p(e+2,1:f+2)=0;
( a2 U& F+ O' R
p(1:e+2,1)=100;p(1:e+2,f+2)=100;
& K7 f0 o$ p1 `4 B
. r2 d3 q3 h9 e& Y% j; k
%定义系数矩阵A
' T7 A; W5 @8 |8 l
A=zeros(e*f,e*f);
& \' Y* F6 l( u. i: F
B=mat2cell(A,ones(e*f/e,1)*e,ones(e*f/e,1)*e);
9 _ s. @: v: E9 v3 B
d1=ones(e,1);d2=ones(e-1,1);
6 m& I0 G0 X) @1 Y8 x( m
M=4*diag(d1)-diag(d2,1)-diag(d2,-1);
5 ~. H s3 y! C/ _/ c4 F
N=-eye(e);
o8 F- [; A( P4 p U
B{1,1}=M;B{e,e}=M
, s, [1 W' e7 H( {% O, _0 g
for i=2:e-1
) l2 ~; E7 }" D- y0 q9 o
B{i,i}=M;B{i-1,i}=N;B{i+1,i}=N
' N% t% G' o- m3 ~
end
) a/ w; X# _; q& F# e# F# h1 _) e* \
A=cell2mat(B);
" T" P" K/ o/ F) d5 N
这里总是显示
/ F* m+ {7 y7 t: S5 u
??? function J=jac(A,b,u0,eps)
4 o( w4 p, ^* s, S' c/ V
|
1 K. ?* U' f ^' B6 t4 [
Error: Function definitions are not permitted at the prompt or in scripts.
9 o- N2 k6 V1 E
/ }3 R5 N0 R' q
%定义b
9 a+ a: j) x/ d. f- `8 h! g3 ]+ k
b=ones(e*f,1);
. N- @6 a/ Z3 |* i
for i=1:e
: O* e. C' W) q" c) L, E
for j=1:f
5 |4 y. `" y- R* S" u6 l4 e
b(i+j)=p(i,j+1)+p(i+1,j)+p(i+2,j+1)+p(i+1,j+2)
; i) u- f2 |2 z; C0 z8 S& D) w
end
( V6 h5 m$ O9 T O6 U1 {
end
8 t8 b# W" z# j4 X4 C$ M
%运用Jacobi迭代法计算
+ B% e$ @- |% e0 Q$ Q& h
D=diag(diag(A));
# [6 ]; r# \& q2 ]+ @5 q. T
D=inv(D);
+ T/ t& Q+ s* Q) m' k/ y
L=tril(A,-1);
3 c1 m3 j. g- {$ O3 |
U=triu(A,1);
# g2 A1 _ _" m5 s. |
B=-D*(L+U);
, N( R+ l: ^! B. v9 w2 x# ?
f=D*b;
# e2 u+ ~0 D# o" ~
J=B*u0+f;
! x& O- Z5 f6 J( H' v
while norm(J-u0)>=eps
' k. {, ]' `' V$ ~+ D
x0=J;
& ^3 x$ `+ \6 i) T/ U" v( C
J=B*u0+f;
4 V- N @; z( C: B! T! r2 l
end
' Y9 g( B, T% p
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