24小時熱門版塊排行榜

>論壇更新日志 (5558)
>考研 (1102)
>導師招生 (1010)
>文獻求助 (410)
>蟲友互識 (264)
>基金申請 (202)
>休閑灌水 (180)
>論文投稿 (140)
>碩博家園 (123)
>考博 (120)
>招聘信息布告欄 (96)
>博后之家 (67)
>教師之家 (65)
>綠色求助(高懸賞) (40)
>公派出國 (37)
>留學生活 (32)

返回列表

小木蟲zb

木蟲 (正式寫手)

應助: 11 (小學生)
金幣: 1560.9
紅花: 5
帖子: 345
在線: 142.1小時
蟲號: 1227944
注冊: 2011-03-10
性別: GG
專業(yè): 分離過程

[求助] 不知道m(xù)atlab程序為什么錯了（錯誤應該存在非線性方程組求解中）

求大神幫忙，我碰到個問題已經苦苦折騰了一個多月了，不知道怎么解決，matlab用線上法求偏微分方程組（對床層高度z離散，將偏微分方程組轉換為常微分方程組，然后用ode15s求解），在這個方程組中，有兩個參數(shù)ce1、ce2需要需要用兩個隱函數(shù)方程求解帶入，下面是我的程序，運行不出，不知道具體錯在哪，只知道可能在隱函數(shù)方程的編程中有些錯誤，
function DBT_MNA
clear all;clc
c10=1e-10*ones(1,20);q10=zeros(1,20);q20=zeros(1,20);c20=zeros(1,20);
y0=[q10 c10 q20 c20];c0=7.03;
options=odeset('relTol',1e-6);tspan=[0:0.2:20];
[t,y]=ode15s(@fangcheng,tspan,y0,options);
plot(t*7.03,y(:,40)/c0,'ro',t*7.03,y(:,80)/c0,'bo','LineWidth',3);xlabel('Amount of treated MDF(g-MDF/g-AC)'),ylabel('c/c0'),axis([-0.05 120 -0.05 1.05]),
legend('DBT','MNA','Location','best') ,title('DBT+MNA-RS')
grid,hold off
%------------odefun---------------------------------------------
function dydt=fangcheng(t,y)
kf=200;a=0.363;u=0.361032/a;c0=7.03;L=8.31e-2;N=20;dz=L/N;p=436.8/(1-0.363);
q1=y(1:N);c1=y(N+1:2*N);q2=y(2*N+1:3*N);c2=y(3*N+1:4*N);
%-----------定義dq/dt-----------------------------
for j=1:N
dq1dt(j)=kf/p*(c1(j)-myfun1(q1(j),q2(j)));%兩個未知參數(shù)調用其他函數(shù)求解ce1(j)=myfun1(q1(j),q2(j))
dq2dt(j)=kf/p*(c2(j)-myfun2(q1(j),q2(j)));%兩個未知參數(shù)調用其他函數(shù)求解ce2(j)=myfun2(q1(j),q2(j))
end
%---------定義dC/dt-----------------------------------
dc1dt(1)=-u*(c1(1)-c0)/dz-(1-a)/a*p*dq1dt(1);
dc2dt(1)=-u*(c2(1)-c0)/dz-(1-a)/a*p*dq2dt(1);
for j=2:N
dc1dt(j)=-u*(c1(j)-c1(j-1))/dz-(1-a)/a*p*dq1dt(j);
dc2dt(j)=-u*(c2(j)-c2(j-1))/dz-(1-a)/a*p*dq2dt(j);
end
dydt=[dq1dt dc1dt dq2dt dc2dt]';
%隱函數(shù)方程數(shù)值求解,隱函數(shù)方程為q1=0.45967*ce1^0.3839/(ce1^0.0457+0.3*ce2^0.0003206);
%q2=0.17356*ce2^(-1.55201)/(ce2^(-2)+0.1*ce1^(-1)),下面用x(1)表示ce1,x(2)表示ce2，上面的程序可以確定沒錯誤，下面是可能出現(xiàn)錯誤的地方
function y1=myfun1(q1,q2)
options=optimset('MaxFunEvals',1e6);
y=fsolve(@(x)fun(x,q1,q2),[0.1,0.01],options);y1=y(1);
function y2=myfun2(q1,q2)
options=optimset('MaxFunEvals',1e6);
y=fsolve(@(x)fun(x,q1,q2),[0.1,0.01],options);y2=y(2);
function y=fun(x,q1,q2)
y(1)=q1-0.45967*x(1)^0.3839/(x(1)^0.0457+0.3*x(2)^0.0003206);
y(2)=q2-0.17356*x(2)^(-1.55201)/(x(2)^(-2)+0.1*x(1)^(-1));
上面是我編的程序，但是運行不出來，不知道錯誤在哪，希望好心人士幫我看下，這是我之前就該問題在小木蟲上的提問http://www.gaoyang168.com/bbs/viewthread.php?tid=5856359

[ Last edited by 小木蟲zb on 2013-5-28 at 22:41 ]

回復此樓

» 猜你喜歡

一志愿中國海洋大學，生物學，301分，求調劑已經有4人回復
一志愿武理材料305分求調劑已經有6人回復
材料專碩英一數(shù)二306 已經有3人回復
085600材料與化工調劑 324分已經有7人回復
311求調劑已經有3人回復
26調劑/材料/英一數(shù)二/總分289/已過A區(qū)線已經有7人回復
295求調劑已經有5人回復
354求調劑已經有5人回復
材料專業(yè)求調劑已經有3人回復
求材料調劑已經有7人回復

1樓 2013-05-28 22:17:16

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

智能機器人

Robot (super robot)

我們都愛小木蟲

找到一些相關的精華帖子，希望有用哦~

求助matlab寫循環(huán) 已經有6人回復
求助MATLAB做三維圖問題已經有4人回復
求大神幫我改一下下面的程序！matlab迭代法的不知道錯在哪里！已經有7人回復
求助��！求一程序，用matlab程序做，用蒙特卡洛方法模擬已經有10人回復
求助matlab編程已經有13人回復
【求助】用matlab編程，需要將數(shù)據(jù)自檢分類（高手進來指點下）已經有5人回復
【求助】請教一個matlab程序，謝謝幫忙！！已經有15人回復
【求助】求助MATLAB程序已經有6人回復
【求助】初學matlab 已經有11人回復
【求助】使用MATLAB怎么實現(xiàn)擬合力場參數(shù)的程序？【已解決】已經有15人回復
【求助】求MAtlab求解程序�。�！已經有20人回復

點擊這里搜索更多相關資源

科研從小木蟲開始，人人為我，我為人人

相關版塊跳轉我要訂閱樓主小木蟲zb 的主題更新