版塊導(dǎo)航: 正在加載中...

應(yīng)《網(wǎng)絡(luò)安全法》要求，自2017年10月1日起，未進(jìn)行實(shí)名認(rèn)證將不得使用互聯(lián)網(wǎng)跟帖服務(wù)。為保障您的帳號(hào)能夠正常使用，請(qǐng)盡快對(duì)帳號(hào)進(jìn)行手機(jī)號(hào)驗(yàn)證，感謝您的理解與支持！

24小時(shí)熱門版塊排行榜

返回列表

holmescn

金蟲(chóng) (正式寫手)

程序強(qiáng)帖: 37
應(yīng)助: 1 (幼兒園)
金幣: 1918.8
帖子: 699
在線: 102.6小時(shí)
蟲(chóng)號(hào): 913482

[交流] 【求助】Matlab 程序結(jié)果不相同，求解

程序的一部分，目標(biāo)是求解一個(gè)類似泊松方程的場(chǎng)方程。原作者不是我，程序?qū)懙煤苤苯�，但很緩慢。�?jīng)過(guò)我修改后，速度雖然提上去了，但貌似兩個(gè)程序的結(jié)果不是很對(duì)，按說(shuō)我都使用的相同的算法，區(qū)別只在寫法上。

方程： $C_1\frac{\partial^2 \Phi}{\partial x^2}+C_2\frac{\partial^2 \Phi}{\partial z^2} =C_3\frac{\partial P}{\partial z}$
周期性邊界條件。

原始代碼：

CODE:

C1 = (2+2*D44*dx^2/(D11*dz^2));

C2 = (2+2*D11*dz^2/(D44*dx^2));

C3 = D44*abs(P0)*sqrt(D11)*dx^2*dz/ ...

    (2*D11*af*sqrt(abs(A))*dz^2+2*D44* ...

    af*sqrt(abs(A))*dx^2);

Nx=101;Nz=101;

P1 = zeros(Nx,Nz);

FI = zeros(Nx,Nz);Fnew=zeros(Nx,Nz);Fold=zeros(Nx,Nz)+1;

while max(max(abs(FI-Fold)))>=0.001

    for k=1:Nz

        for i=1:Nx

            if i~=1&&i~=Nx&&k~=1&&k~=Nz

                Fnew(i,k)=(FI(i+1,k)+FI(i-1,k))/C1+(FI(i,k+1)+FI(i,k-1))/C2- ...

                          (P1(i,k)-P1(i,k-1))*C3;

            end

            if k==1

                if i~=1&&i~=Nx

                    Fnew(i,k)=(FI(i+1,k)+FI(i-1,k))/C1+ ...

                              (FI(i,k+1)+FI(i,Nz))/C2- ...

                              (P1(i,k)-P1(i,Nz))*C3;

                end

                if i==1

                    Fnew(i,k)=(FI(i+1,k)+FI(Nx,k))/C1+ ...

                              (FI(i,k+1)+FI(i,Nz))/C2- ...

                              (P1(i,k)-P1(i,Nz))*C3;

                end

                if i==Nx

                    Fnew(i,k)=(FI(1,k)+FI(i-1,k))/C1+ ...

                              (FI(i,k+1)+FI(i,Nz))/C2- ...

                              (P1(i,k)-P1(i,Nz))*C3;

                end

            end

            if k==Nz

                if i~=1&&i~=Nx

                    Fnew(i,k)=(FI(i+1,k)+FI(i-1,k))/C1+ ...

                              (FI(i,1)+FI(i,k-1))/C2- ...

                              (P1(i,k)-P1(i,k-1))*C3;

                end

                if i==1

                    Fnew(i,k)=(FI(i+1,k)+FI(Nx,k))/C1+ ...

                              (FI(i,1)+FI(i,k-1))/C2- ...

                              (P1(i,k)-P1(i,k-1))*C3;

                end

                if i==Nx

                    Fnew(i,k)=(FI(1,k)+FI(i-1,k))/C1+ ...

                              (FI(i,1)+FI(i,k-1))/C2- ...

                              (P1(i,k)-P1(i,k-1))*C3;

                end

            end

            if i==1&&k~=1&&k~=Nz

                Fnew(i,k)=(FI(i+1,k)+FI(Nx,k))/C1+ ...

                          (FI(i,k+1)+FI(i,k-1))/C2- ...

                          (P1(i,k)-P1(i,k-1))*C3;

            end

            if i==Nx&&k~=1&&k~=Nz

                Fnew(i,k)=(FI(1,k)+FI(i-1,k))/C1+ ...

                          (FI(i,k+1)+FI(i,k-1))/C2- ...

                          (P1(i,k)-P1(i,k-1))*C3;

            end

        end

    end

    Fold=FI;FI=Fnew;

end

注意這里的P1只是給了一個(gè)形式，為的是說(shuō)明他的結(jié)構(gòu)，計(jì)算中的P1是有一個(gè)分布的，不是零。

修改后的代碼：

CODE:

C1 = (2+2*D44*dx^2/(D11*dz^2));

C2 = (2+2*D11*dz^2/(D44*dx^2));

C3 = D44*abs(P0)*sqrt(D11)*dx^2*dz/ ...

    (2*D11*af*sqrt(abs(A))*dz^2+2*D44* ...

    af*sqrt(abs(A))*dx^2);

Nx=101;Nz=101;

P1=zeros(Nx+2,Nz+2);

FI=zeros(Nx+2,Nz+2);Fnew=zeros(Nx+2,Nz+2);Fold=zeros(Nx+2,Nz+2)+1;

i=2:Nx+1;k=2:Nz+1;

while max(max(abs(FI(i,k)-Fold(i,k))))>=0.001

    P1(1,:)=P1(Nx+1,:);P1(Nx+2,:)=P1(2,:);

    P1(:,1)=P1(:,Nz+1);P1(:,Nx+2)=P1(:,;2);

    FI(1,:)=FI(Nx+1,:);FI(Nx+2,:)=FI(2,:);

    FI(:,1)=FI(:,Nz+1);FI(:,Nx+2)=P1(:,2);

    Fnew(i,k)=(FI(i+1,k)+FI(i-1,k))/C1+ ...

              (FI(i,k+1)+FI(i,k-1))/C2- ...

              (P1(i,k)-P1(i,k-1))*C3;

    Fold=FI;FI=Fnew;

end

原始代碼迭代次數(shù)為498次，修改后要多迭代100次，而且結(jié)果還對(duì)不上。注意程序里的C1，C2，C3和方程里的不是一回事。

問(wèn)題：
1、第二段程序什么地方不正確，導(dǎo)致結(jié)果不同。
2、這樣的方程是不是有什么方法使用庫(kù)函數(shù)求解。
3、這種迭代格式是否正確。
4、還沒(méi)想到……

回復(fù)此樓

» 猜你喜歡

求調(diào)劑一志愿海大，0703化學(xué)學(xué)碩304分，有大創(chuàng)項(xiàng)目，四級(jí)已過(guò) 已經(jīng)有9人回復(fù)
289材料與化工（085600）B區(qū)求調(diào)劑已經(jīng)有4人回復(fù)
一志愿武理材料工程348求調(diào)劑已經(jīng)有6人回復(fù)
08工科 320總分求調(diào)劑已經(jīng)有11人回復(fù)
一志愿華中農(nóng)業(yè)071010，總分320求調(diào)劑已經(jīng)有6人回復(fù)
323求調(diào)劑已經(jīng)有5人回復(fù)
284求調(diào)劑已經(jīng)有5人回復(fù)
317求調(diào)劑已經(jīng)有16人回復(fù)
287求調(diào)劑已經(jīng)有9人回復(fù)
求調(diào)劑已經(jīng)有6人回復(fù)

» 搶金幣啦！回帖就可以得到:

查看全部散金貼

1樓 2010-12-09 18:31:36

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

holmescn

金蟲(chóng) (正式寫手)

程序強(qiáng)帖: 37
應(yīng)助: 1 (幼兒園)
金幣: 1918.8
帖子: 699
在線: 102.6小時(shí)
蟲(chóng)號(hào): 913482

★
余澤成(金幣+1):謝謝分享，金幣由版面來(lái)發(fā)！ 2010-12-12 20:46:36

看來(lái)程序太長(zhǎng)的話，大家都很頭痛。這個(gè)問(wèn)題已經(jīng)被我?guī)煹芙鉀Q了，是我的錯(cuò)誤。
見(jiàn)代碼段2倒數(shù)第6行末尾的P1應(yīng)該是FI。又是一次復(fù)制粘貼的錯(cuò)誤。

不過(guò)，我的問(wèn)題雖然解決了，本帖也很好地展現(xiàn)了快速matlab程序的寫法。代碼段2的運(yùn)行效率大約是代碼段1的10倍。此貼可供初學(xué)者圍觀。

金幣將發(fā)給有價(jià)值的回帖討論。

贊一下(1人)

回復(fù)此樓

2樓2010-12-11 10:51:12

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

相關(guān)版塊跳轉(zhuǎn) 我要訂閱樓主 holmescn 的主題更新

返回列表

普通表情龍兔虎貓高級(jí)回復(fù) (可上傳附件)

最具人氣熱帖推薦 [查看全部]		作者	回/看	最后發(fā)表

[考研] 廣西大學(xué)材料導(dǎo)師推薦 +3	夏夏夏小正 2026-03-17	5/250	2026-03-21 22:20 by 金昊ML
[考研] 材料工程專碩 348分求調(diào)劑 +3	冬辭. 2026-03-17	5/250	2026-03-21 18:47 by 學(xué)員8dgXkO
[考研] 278求調(diào)劑 +9	煙火先于春 2026-03-17	9/450	2026-03-21 17:47 by 學(xué)員8dgXkO
[考研] 307求調(diào)劑 +3	余意卿 2026-03-18	3/150	2026-03-21 17:31 by ColorlessPI
[考研] 299求調(diào)劑 +5	shxchem 2026-03-20	7/350	2026-03-21 17:09 by ColorlessPI
[考研] 求調(diào)劑 +3	.m.. 2026-03-21	4/200	2026-03-21 16:25 by barlinike
[考研] 299求調(diào)劑 +6	△小透明* 2026-03-17	6/300	2026-03-21 02:42 by JourneyLucky
[考研] 一志愿中國(guó)石油大學(xué)（華東）本科齊魯工業(yè)大學(xué) +3	石能偉 2026-03-17	3/150	2026-03-21 02:22 by JourneyLucky
[考研] 一志愿華南師大 070300（化學(xué)）304分求調(diào)劑 +3	0703武芊慧雪304 2026-03-18	3/150	2026-03-21 00:48 by JourneyLucky
[考研] 材料專碩英一數(shù)二306 +7	z1z2z3879 2026-03-18	7/350	2026-03-20 23:48 by JourneyLucky
[考研] 317求調(diào)劑 +5	申子申申 2026-03-19	9/450	2026-03-20 22:26 by JourneyLucky
[考研] 材料與化工 322求調(diào)劑 +4	然11 2026-03-19	4/200	2026-03-20 22:12 by luoyongfeng
[考研] 295復(fù)試調(diào)劑 +8	簡(jiǎn)木ChuFront 2026-03-19	8/400	2026-03-20 20:44 by zhukairuo
[考研] 一志愿南理工085701環(huán)境302求調(diào)劑院校 +3	葵梓衛(wèi)隊(duì) 2026-03-20	3/150	2026-03-20 19:28 by zhukairuo
[考研] 320求調(diào)劑0856 +3	不想起名字112 2026-03-19	3/150	2026-03-19 22:53 by 學(xué)員8dgXkO
[考博] 申博26年 +3	八6八68 2026-03-19	3/150	2026-03-19 19:43 by nxgogo
[考研] 一志愿福大288有機(jī)化學(xué)，求調(diào)劑 +3	小木蟲(chóng)200408204 2026-03-18	3/150	2026-03-19 13:31 by houyaoxu
[考研] 085601求調(diào)劑 +4	Du.11 2026-03-16	4/200	2026-03-17 17:08 by ruiyingmiao
[考研] 考研調(diào)劑 +3	淇ya_~ 2026-03-17	5/250	2026-03-17 09:25 by Winj1e
[考研] 070303 總分349求調(diào)劑 +3	LJY9966 2026-03-15	5/250	2026-03-16 14:24 by xwxstudy

亭亭五月天在线观看,亭亭五月天在线观看,国产最新av一区二区,国产 高清 中文字幕,99re热久久亚洲综合精品成人,熟妇 一区二区三区,一级做a爰片性色毛片武则天,美女的骚穴视频播放,国产美女午夜免费视频

24小時(shí)熱門版塊排行榜

holmescn

[交流] 【求助】Matlab 程序結(jié)果不相同，求解

» 猜你喜歡

» 搶金幣啦！回帖就可以得到: 查看全部散金貼

holmescn

亭亭五月天在线观看,亭亭五月天在线观看,国产最新av一区二区,国产高清中文字幕,99re热久久亚洲综合精品成人,熟妇一区二区三区,一级做a爰片性色毛片武则天,美女的骚穴视频播放,国产美女午夜免费视频

[交流] 【求助】Matlab 程序結(jié)果不相同，求解

» 搶金幣啦！回帖就可以得到:

查看全部散金貼