QQ登录

只需要一步,快速开始

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

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

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

1

主题

2

听众

35

积分

升级  31.58%

该用户从未签到

自我介绍
20100103重要的日子
跳转到指定楼层
1#
发表于 2009-12-17 20:17 |只看该作者 |倒序浏览
|招呼Ta 关注Ta
%运用Jacobi迭代拟合出来的u关于x和y的矩阵  _$ {- R& q, v
function J=jac(A,b,u0,eps)
8 `" o6 w% X0 I: |+ J" Qif nargin==3
& ^1 u! H4 y. p. S+ q    eps=1.0e-8/ B# |( {5 o7 n2 `
elseif nargin<3
, J+ a5 i1 b7 ]4 o" }; `    'error'/ A7 w/ V, Y& J# k- X
    return
; S" }# [. X; Rend$ {9 y0 l# e' p/ c* k' F- a0 W5 M
" M: x& r/ J! J
%定义内部节点矩阵u0
. {6 }0 Q- J% o3 H4 z% ph=1;k=1;
# {9 j; e9 x4 B% @& d# _. px=0:h:17;y=0:k:10;
3 n' T3 K9 q4 e% e/ |) ie=length(x)-2;f=length(y)-2;7 Z' c& u) M$ b7 V) X
u=zeros(e,f);
9 {* ^9 x4 `5 F+ I& y4 n! `0 du0=u;
9 }: t' @* p* U
+ ?4 M# A& Q8 I6 M. y%定义外部节点p( l, D) Y3 D3 {+ u" S
p=zeros(e+2,e+2);
1 X# q* D" q- dp(1,1:f+2)=0;p(e+2,1:f+2)=0;       8 e- m: T# B: y( a$ V
p(1:e+2,1)=100;p(1:e+2,f+2)=100;
' i: \0 |% n6 M
8 h. C1 V/ V5 B%定义系数矩阵A6 y8 g$ V8 ~# _# M0 t% [! F, ?
A=zeros(e*f,e*f);
' o# D3 I4 ^/ g3 j3 A" dB=mat2cell(A,ones(e*f/e,1)*e,ones(e*f/e,1)*e);
( Z' Z# d6 z$ p' j5 qd1=ones(e,1);d2=ones(e-1,1);/ I4 I8 c( |" _! {" f6 @) I
M=4*diag(d1)-diag(d2,1)-diag(d2,-1);
( L, `* o5 f$ n) C; C4 rN=-eye(e);* E: v' M; j' i3 p2 n% `9 @
B{1,1}=M;B{e,e}=M
' l5 H. s" x' d  H4 A( M$ bfor i=2:e-1
7 X; G( E8 Z# R6 ^; h    B{i,i}=M;B{i-1,i}=N;B{i+1,i}=N
. ?, x' O1 q) }1 |* d5 K( Qend) X- p: U+ Z3 e! O
A=cell2mat(B);  
7 R7 z) O2 }; s; t  b7 Z* A+ ?这里总是显示
) b- I$ F8 Q* g) F' n??? function J=jac(A,b,u0,eps)$ R" e" \' z9 J
    |3 C. H& b( M3 J  t6 w' X
Error: Function definitions are not permitted at the prompt or in scripts.

2 w& N  H7 r6 b# ~* [( c& ^. m: o1 ^/ Y6 A5 b( C7 }' _) w9 t" `
%定义b. @8 w( ?9 s& c8 S0 ^
b=ones(e*f,1);  T) o5 T/ q5 t9 ^1 ~3 v2 c. j
for i=1:e' I- m% _* I' c) U) B" @9 M2 s
    for j=1:f
3 l! a4 i& ]' a7 d9 G$ j/ \* U% \3 {" @        b(i+j)=p(i,j+1)+p(i+1,j)+p(i+2,j+1)+p(i+1,j+2)
0 J2 v' p( \- T& j3 ~) v    end+ t# _. G* \% \& E/ N6 y* w3 T
end
% K2 }9 ^2 w' S6 x. a& t%运用Jacobi迭代法计算2 |) C# D( i1 {* V$ C
D=diag(diag(A));4 [: d; d$ B! l% h+ f, R
D=inv(D);
# w( S+ E/ C% }# i# z  TL=tril(A,-1);
3 G# J4 Z8 b+ u* X/ ^& QU=triu(A,1);6 h4 k4 ?) S7 i- U( h4 P6 s
B=-D*(L+U);
; [2 k: v- \5 o7 J1 @% Kf=D*b;' @" ^, y; d' M' F3 [2 J( C
J=B*u0+f;  Y" s0 Z: e3 P* v. V" I# [+ m) F
while norm(J-u0)>=eps# x: l" R# X% I' ?6 |! Y
x0=J;
. N% c2 S% A' l5 R4 wJ=B*u0+f;
- ~- S6 O) M: W8 A7 Cend
1 u$ I& V0 G; O* _0 C1 e* freturn
zan
转播转播0 分享淘帖0 分享分享0 收藏收藏0 支持支持0 反对反对0 微信微信

1

主题

2

听众

35

积分

升级  31.58%

该用户从未签到

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

使用道具 举报

BenCam 实名认证       

9

主题

6

听众

89

积分

该用户从未签到

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

使用道具 举报

madio        

3万

主题

1312

听众

5万

积分

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

    [LV.Master]伴坛终老

    自我介绍
    数学中国站长

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

    群组数学建模培训课堂1

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

    群组Matlab讨论组

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

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

    回复

    使用道具 举报

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

    qq
    收缩
    • 电话咨询

    • 04714969085
    fastpost

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

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

    蒙公网安备 15010502000194号

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

    GMT+8, 2026-4-14 11:40 , Processed in 0.438379 second(s), 73 queries .

    回顶部