版塊導航: 正在加載中...

登錄注冊

應《網絡安全法》要求，自2017年10月1日起，未進行實名認證將不得使用互聯網跟帖服務。為保障您的帳號能夠正常使用，請盡快對帳號進行手機號驗證，感謝您的理解與支持！

24小時熱門版塊排行榜

返回列表

當前只顯示滿足指定條件的回帖，點擊這里查看本話題的所有回帖

iMitsuha

新蟲 (初入文壇)

應助: 0 (幼兒園)
帖子: 7
在線: 13.5小時
蟲號: 8077957
注冊: 2018-02-24
專業(yè): 傳遞過程

[求助] Matlab求解非線性偏微分方程組已有1人參與

CODE:

在課題研究中建立了一個非線性偏微分方程組，用Matlab中pdepe函數求解時，無法求出結果，Matlab一直正忙，我不是數學專業(yè)的，一直找不出問題，懇請各位指導一下。謝謝！

程序如下
function pdexx11
m = 0;
x = linspace(0,0.02,3);
t = linspace(0,1,3);
sol = pdepe(m,@pdefun,@icfun,@bcfun,x,t);
u1 = sol(:,:,1);
u2 = sol(:,:,2);
u3 = sol(:,:,3);

figure;

surf(x,t,u1)

title('u1(x,t)')

xlabel('Distance x')

ylabel('Time t')

figure;

surf(x,t,u2)

title('u2(x,t)')

xlabel('Distance x')

ylabel('Time t')

figure;

surf(x,t,u3)

title('u3(x,t)')

xlabel('Distance x')

ylabel('Time t')

function [ c,f,s ] = pdefun( x,t,u,DuDx )
c = [1;1;1];
f = [1.25e-9;
   1.828e-9;
   1.856e-9].*DuDx+[0.59*u(1);0.59*u(2);0.59*u(3)];
s = [-1.03e-4/60*1.9*((u(1)*u(2)*u(3)/(5.51e-14))^(1/3)-1);
   -1.03e-4/60*1.9*((u(1)*u(2)*u(3)/(5.51e-14))^(1/3)-1);
   -1.03e-4/60*1.9*((u(1)*u(2)*u(3)/(5.51e-14))^(1/3)-1)];

function u0 = icfun(x)
u0 = [0;0;0];

function [ pl,ql,pr,qr ] = bcfun( xl,ul,xr,ur,t )
pl = [ul(1)-0.0493;
   ul(2)-0.07142;
   ul(3)-0.0387];
ql = [0;0;0];
pr = [ur(1)-0.0106;ur(2)-0.03272;ur(3)];
qr = [0;0;0];

回復此樓

» 猜你喜歡

伙伴們，祝我生日快樂吧已經有5人回復
274求調劑0856材料化工已經有10人回復
337一志愿華南理工0805材料求調劑已經有6人回復
（081700）化學工程與技術-298分求調劑已經有7人回復
材料與化工304求B區(qū)調劑已經有5人回復
NSFC申報書里申請人簡歷中代表性論著還需要在申報書最后的附件里面再上傳一遍嗎已經有9人回復
一志愿北京化工大學材料與化工296分求調劑已經有17人回復
0703化學一志愿211 總分320求調劑已經有3人回復
0703求調劑已經有6人回復
一志愿天津大學材料與化工275求調劑已經有22人回復

1樓 2019-02-25 10:09:10

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

ofwhy

至尊木蟲 (著名寫手)

應助: 3 (幼兒園)
金幣: 13429.6
散金: 93
紅花: 7
帖子: 2309
在線: 177.9小時
蟲號: 2777064
注冊: 2013-11-04
性別: GG
專業(yè): 工程熱力學

【答案】應助回帖

可能公式推導有錯誤，如果按照你導出的公式計算，會計算出虛值，然后不停的死循環(huán)。如將pu3/pt的最后k*Ac前面的負號改為正號則有解。從你的邊界條件來看，似乎這個分量也應該跟其他兩個應該是不一樣的。
function [ output_args ] = PDEs1( input_args )
m=0;
x=linspace(0,0.02,3);
t=linspace(0,1,3);
sol=pdepe(m,@pdes1pde,@pdes1ic,@pdes1bc,x,t);
u1=sol(:,:,1);
u2=sol(:,:,2);
u3=sol(:,:,3);

figure
surf(x,t,u1)
title('u1(x,t)')
xlabel('Distance x')
ylabel('Time t')

figure
surf(x,t,u2)
title('u2(x,t)')
xlabel('Distance x')
ylabel('Time t')

figure
surf(x,t,u3)
title('u3(x,t)')
xlabel('Distance x')
ylabel('Time t')

function[c,f,s]=pdes1pde(x,t,u,DuDx)
      c=[1;1;1];
      D1=1.25*10^(-9);
      D2=1.838*10^(-9);
      D3=1.856*10^(-9);
      V=0.59;
      k=1.03*10^(-4)/60;
      Ac=1.9;
      Ksp=5.51*10^(-14);
      f=[D1/(1+V);D2/(1+V);D3/(1+V)].*DuDx;
      s=[-k*Ac/(1+V)*((u(1)*u(2)*u(3)/Ksp)^(1/3)-1);-k*Ac/(1+V)*((u(1)*u(2)*u(3)/Ksp)^(1/3)-1);k*Ac/(1+V)*((u(1)*u(2)*u(3)/Ksp)^(1/3)-1);];
end
function u0 = pdes1ic(x)
      u0 = [0;0;0];
end
function [pl,ql,pr,qr] = pdes1bc(xl,ul,xr,ur,t)
      pl = [ul(1)-0.0493; ul(2)-0.07142;ul(3)-0.0387];
      ql = [0;0;0];
      pr = [ur(1)-0.0106;ur(2)-0.03272;ur(3)];
      qr = [0;0;0];
end

end

贊一下

回復此樓

4樓2019-03-31 12:08:01

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

相關版塊跳轉我要訂閱樓主 iMitsuha 的主題更新