QQ登录

只需要一步,快速开始

 注册地址  找回密码
查看: 3041|回复: 3
打印 上一主题 下一主题

求助:关于分块矩阵的还原,为什么实行不了

[复制链接]
字体大小: 正常 放大

1

主题

2

听众

35

积分

升级  31.58%

该用户从未签到

自我介绍
20100103重要的日子
跳转到指定楼层
1#
发表于 2009-12-17 20:17 |只看该作者 |正序浏览
|招呼Ta 关注Ta
%运用Jacobi迭代拟合出来的u关于x和y的矩阵9 h+ o% G6 P4 X! U+ S, H
function J=jac(A,b,u0,eps)8 i$ \0 i/ j" E7 p* O2 l' ?2 j
if nargin==3
1 s9 z( W0 e2 X. f9 s' t    eps=1.0e-8
  P! _" P5 ]  }  x0 belseif nargin<3# D# a$ ^) [! A3 C; ?) ~+ [
    'error'
0 D, V8 i& H, W- W    return' ]  I# c& b3 @. L9 P- n
end
! M8 Q9 u! f5 G' k- t9 _! m# x. B5 }  c9 M
%定义内部节点矩阵u0. n$ r  A) E) M- B
h=1;k=1;
1 s; U) }0 B, l3 x3 jx=0:h:17;y=0:k:10;
- {& X  P4 e' |7 n) V# N2 ce=length(x)-2;f=length(y)-2;
+ n# S& e; V6 [! w8 ?0 m+ lu=zeros(e,f);# |, |: z+ T) h+ h
u0=u;& ~4 H# P. O2 m; h$ Q' S5 h5 ^' V
5 i. [. t' L3 _" |! O
%定义外部节点p
5 t6 P& K( S9 U3 R; xp=zeros(e+2,e+2);" f: ]( f- d# |
p(1,1:f+2)=0;p(e+2,1:f+2)=0;       ( \; w" N0 \+ C8 _  b9 X# d# l9 K: [
p(1:e+2,1)=100;p(1:e+2,f+2)=100;! ?0 p3 P7 g! Q( e* y
' T- e/ C8 ~+ V  v$ C/ Q
%定义系数矩阵A
' x/ i. n* a) y  n! `A=zeros(e*f,e*f);
% n4 t: m2 }+ b3 D8 o# r, WB=mat2cell(A,ones(e*f/e,1)*e,ones(e*f/e,1)*e);7 i; ?) S4 Z2 j1 `* @8 F0 s
d1=ones(e,1);d2=ones(e-1,1);  \3 P" i/ u5 I8 g4 _' q: c, [
M=4*diag(d1)-diag(d2,1)-diag(d2,-1);
" H) A, ~$ H; {% T( rN=-eye(e);
$ x- b% w4 X0 j% l& S# S" q: ]7 VB{1,1}=M;B{e,e}=M
* w) M: B; J: b( jfor i=2:e-17 a  W3 q3 `5 T8 R
    B{i,i}=M;B{i-1,i}=N;B{i+1,i}=N% N; n$ i& H/ M
end. {8 B& `9 w7 ]4 J, p" z
A=cell2mat(B);  / z- e1 u4 j; `9 |- v2 V! m4 _
这里总是显示
/ v/ }$ h0 Q0 |7 M* r* u* x+ ~9 e9 r??? function J=jac(A,b,u0,eps)
5 T0 ^1 t' J  Y) `  [    |
, t" {. y& x* F7 x9 oError: Function definitions are not permitted at the prompt or in scripts.

7 B! R2 q/ L- N! l  f9 L1 N
) ~$ }) @, q8 s%定义b
% P, l' r6 U9 F- Bb=ones(e*f,1);
% W1 T5 P  [  Y. r. Sfor i=1:e
( s; {# r, R+ B6 I8 w# \    for j=1:f
# i* [( d+ w. P        b(i+j)=p(i,j+1)+p(i+1,j)+p(i+2,j+1)+p(i+1,j+2)
; [7 a* H% U4 x7 s; R3 m( ]    end6 A$ r  G$ {# Z! J' m$ @9 J7 L
end$ `. X1 h  ^3 R% ?; [  y
%运用Jacobi迭代法计算
! h( Z$ S9 E: c. J9 R0 hD=diag(diag(A));
3 J3 c, y, {6 S$ g, tD=inv(D);
, l: p. }8 [2 ]# LL=tril(A,-1);8 m  G! X' O# [3 m0 H  {+ X% ]
U=triu(A,1);
$ S* Q  `$ j% |0 [+ ~+ \9 xB=-D*(L+U);" @" h% N# H( E  d; ^
f=D*b;
' T0 K6 V, ^: o1 n/ c, ~9 ^J=B*u0+f;
; j  T, U+ O1 y9 D* Zwhile norm(J-u0)>=eps
  l; E1 L8 F( v5 j1 S: H4 px0=J;
0 h/ {8 K2 _7 n. e/ Z" y2 M' GJ=B*u0+f;
" n+ v9 h8 Y. G  R8 Q: W6 _end6 `# a% e( {' h9 h
return
zan
转播转播0 分享淘帖0 分享分享0 收藏收藏0 支持支持0 反对反对0 微信微信
madio        

3万

主题

1312

听众

5万

积分

  • TA的每日心情
    奋斗
    2024-7-1 22:21
  • 签到天数: 2014 天

    [LV.Master]伴坛终老

    自我介绍
    数学中国站长

    社区QQ达人 邮箱绑定达人 优秀斑竹奖 发帖功臣 风雨历程奖 新人进步奖 最具活力勋章

    群组数学建模培训课堂1

    群组数学中国美赛辅助报名

    群组Matlab讨论组

    群组2013认证赛A题讨论群组

    群组2013认证赛C题讨论群组

    回复

    使用道具 举报

    BenCam 实名认证       

    9

    主题

    6

    听众

    89

    积分

    该用户从未签到

    自我介绍
    200 字节以内
    不支持自定义 Discuz! 代码
    回复

    使用道具 举报

    1

    主题

    2

    听众

    35

    积分

    升级  31.58%

    该用户从未签到

    自我介绍
    20100103重要的日子
    自己顶一下,拜托哪位高手指点一下,纠结这个矩阵的还原,想了好多方法还是不行~~实在想不出哪里出错了
    回复

    使用道具 举报

    您需要登录后才可以回帖 登录 | 注册地址

    qq
    收缩
    • 电话咨询

    • 04714969085
    fastpost

    关于我们| 联系我们| 诚征英才| 对外合作| 产品服务| QQ

    手机版|Archiver| |繁體中文 手机客户端  

    蒙公网安备 15010502000194号

    Powered by Discuz! X2.5   © 2001-2013 数学建模网-数学中国 ( 蒙ICP备14002410号-3 蒙BBS备-0002号 )     论坛法律顾问:王兆丰

    GMT+8, 2026-7-11 12:33 , Processed in 0.781396 second(s), 73 queries .

    回顶部