Re: [問題] 解非線性聯立方程式(已爬文)已回收

看板MATLAB作者zupo (幫解MATLAB難題囧)時間17年前 (2008/11/10 17:34)推噓0(0推 0噓 0→)

留言0則, 0人參與討論串3/3 (看更多)

※ 引述《hagry (阿貴)》之銘言： : 我的問題是要解聯立方程式 : 在前面解的時候都不會有問題， : 但是用solve解時，卻出現此方程式非有效式子， : 請各位大大幫我看ㄧ下，給我ㄧ下意見， : 假如有數值解和逼近解怎麼寫?? : 我的程式 : K=0.1; : D=0.09; : S=0.2; : a=0.5+sqrt(0.25+2*D/S^2); : b=0.5-sqrt(0.25+2*D/S^2); : syms uP lP; : %uP = sym('uP','real'); : %uP = sym('lP','real'); : F ='0=((b*uP^a*(uP^(1-a)*lP^(b-a)-lP^(1-a)*uP^(b-a)))-(a*uP^b*(uP^(1-a)-lP^(1-a))))/(a*b*D*(lP^(b-a)-uP^(b-a))+(1-uP )/D+K'; : G ='0=((b*lP^a*(uP^(1-a)*lP^(b-a)-lP^(1-a)*uP^(b-a)))-(a*lP^b*(uP^(1-a)-lP^(1-a))))/(a*D*(lP^(b-a)-uP^(b-a))+(1-lP) /D-K '; : %以上在跑的時候都沒問題 : [uP,lP]=solve(F,G,'uP,lP') : %但是在解這步的時候卻出現is not a valid expression or equation : %我就換一個fsolve去解，結果他說FSOLVE only accepts inputs of data type double. : 我要探討的是在K不同變化下，uP和lP的變化， : 前提是要先解出來這個問題，我才能繼續下去，謝謝各位，感謝。 : 還是說我要去學其他程式?? 我用數值解的方式解了.畢竟這麼大串長的一團.符號解析解解出來大概也不是人類所能看的懂得. ------------------------------------------------------------------------- function pttex119 clc format long % uP == lP (positive) a1 = fsolve(@verylong,[2 2]) % uP ~= lP (positive) a2 = fsolve(@verylong,[90 10]) % uP == lP (negative) a3 = fsolve(@verylong,[-2 -2]) % uP ~= lP (negative) a4 = fsolve(@verylong,[-90 -10]) % uP ~= lP (final_1) a5 = fsolve(@verylong,[1 -10]) % uP ~= lP (final_2) a6 = fsolve(@verylong,[10 -10]) function f = verylong(x) K=0.1;D=0.09;S=0.2; a=0.5+sqrt(0.25+2*D/S^2); b=0.5-sqrt(0.25+2*D/S^2); f(1) = ((b*x(1)^a*(x(1)^(1-a)*x(2)^(b-a)-x(2)^(1-a)*x(1)^(b-a)))-... (a*x(1)^b*(x(1)^(1-a)-x(2)^(1-a))))/(a*b*D*(x(2)^(b-a)-x(1)^(b-a)))... +(1-x(1) )/D+K ; f(2) = ((b*x(2)^a*(x(1)^(1-a)*x(2)^(b-a)-x(2)^(1-a)*x(1)^(b-a)))-... (a*x(2)^b*(x(1)^(1-a)-x(2)^(1-a))))/(a*D*(x(2)^(b-a)-x(1)^(b-a)))... +(1-x(2)) /D-K ; -------------------------------------------------------------------------- 解大致上分為六種 1.兩個起始猜測值相等皆正數 uP = lP 2.兩個起始猜測值不相等皆正數唯一解 3.兩個起始猜值相等皆為負數 uP = lP 4.兩個起始猜值不相等皆為負數唯一複數解(不共厄) 5.兩個起始猜值異號絕對值不相等唯一複數解(不共厄) 且不跟4.同值 6.兩個起始猜測值異號絕對值相等 uP = lP 共厄複數以下是執行結果 Warning: Matrix is singular, close to singular or badly scaled. Results may be inaccurate. RCOND = NaN. > In optim\private\dogleg at 84 In optim\private\trustnleqn at 205 In fsolve at 295 In pttex119 at 5 Optimization terminated: no further progress can be made. Trust-region radius less than 2*eps. Problem may be ill-conditioned or Jacobian may be inaccurate. Try using exact Jacobian or check Jacobian for errors. a1 = 2 2 Optimization terminated: first-order optimality is less than options.TolFun. a2 = 1.015777627938991 1.002281634655273 Warning: Matrix is singular, close to singular or badly scaled. Results may be inaccurate. RCOND = NaN. > In optim\private\dogleg at 84 In optim\private\trustnleqn at 205 In fsolve at 295 In pttex119 at 9 Optimization terminated: no further progress can be made. Trust-region radius less than 2*eps. Problem may be ill-conditioned or Jacobian may be inaccurate. Try using exact Jacobian or check Jacobian for errors. a3 = -2 -2 Optimization terminated: first-order optimality is less than options.TolFun. a4 = 0.924516204796464 - 0.891463732532082i -0.415022055625012... + 0.806951489943288i Optimization terminated: first-order optimality is less than options.TolFun. a5 = 0.517994557698605 + 0.690508711744296i -0.591428251673179 ... + 0.336262180297366i Optimizer appears to be converging to a point which is not a root. Norm of relative change in X is less than max(options.TolX^2,eps) but sum-of-squares of function values is greater than or equal to sqrt... (options.TolFun) Try again with a new starting guess. a6 = 10.000000000000870 - 0.000000000002203i -9.999999999998328 ... - 0.000000000002198i -- 1.MATLAB programming 2.ASPEN process simulation package 3.FORTRN programming 4.Advance Engineering Mathematics 5.Process Control Theory 6.Chemical Engineering Basic Theory(Kinetic.thermodynamics.transport) 7.Numerical Method and Analysis 8.MATLAB Toolbox.Simulink system basic design -- ※ 發信站: 批踢踢實業坊(ptt.cc) ◆ From: 124.9.133.177