| 查看: 1553 | 回复: 8 | ||
[求助]
matlab实现牛顿迭代法求非线性方程组程序运行出现问题 已有3人参与
|
|
已知非线性方程组如下： log(rhomcA*x1*(1-yA*rhob/rhomcA-yB*rhob/rhomcB)/(yA*rhob*(1-x1-x2)))+EAA/T*(7*x1-8*yA*rhob/rhomcA)+(EAA*EBB)^0.5/T*(7*x2-8*yB*rhob/rhomcB)+EAs/T=0 f2=log(rhomcB*x2*(1-yA*rhob/rhomcA-yB*rhob/rhomcB)/(yB*rhob*(1-x1-x2)))+EBB/T*(7*x2-8*yB*rhob/rhomcB)+(EAA*EBB)^0.5/T*(7*x1-8*yA*rhob/rhomcA)+EBs/T=0 首先建立函数fun function f=fun(x); %定义非线性方程组如下 %变量x1 x2 %函数f1 f2 T=318.2; %(K) rhomcA=0.34; %(g/cm3) CA=3.26; %(mmol/g) EAA=64; %(K) EAs=-1385; %(K) rhomcB=0.98; %(g/cm3) CB=4.53; %(mmol/g) EBB=82; %(K) EBs=-1690; %(K) rhob=[0.0047 0.0090 0.0177 0.0264 0.0367 0.0466 0.0566 0.0674 0.0787 0.0890]'; yA=[0.89 0.886 0.884 0.883 0.877 0.871 0.866 0.861 0.856 0.852]'; yB=[0.11 0.114 0.116 0.117 0.123 0.129 0.134 0.139 0.144 0.148]'; syms x1 x2 x3 x4 f1=log(rhomcA*x1.*(1-yA.*rhob/rhomcA-yB.*rhob/rhomcB)./(yA.*rhob.*(1-x1-x2)))+EAA/T*(7*x1-8*yA.*rhob/rhomcA)+(EAA*EBB)^0.5/T*(7*x2-8*yB.*rhob/rhomcB)+EAs/T; f2=log(rhomcB*x2.*(1-yA.*rhob/rhomcA-yB.*rhob/rhomcB)./(yB.*rhob.*(1-x1-x2)))+EBB/T*(7*x2-8*yB.*rhob/rhomcB)+(EAA*EBB)^0.5/T*(7*x1-8*yA.*rhob/rhomcA)+EBs/T; %f3=nexA-2*CA*(x1-x3*rhob/rhomcA); %f4=nexB-2*CB*(x2-x4*rhob/rhomcB); f=[f1 f2]; 建立函数dfun function df=dfun(x); %用来求解方程组的雅克比矩阵储存在dfun中 f=fun(x); df=[diff(f,'x1');diff(f,'x2')]; df=conj(df'); 编程牛顿法求解非线性方程组将newton.m保存到工作路径中: function x=newton(x0,eps,N); con=0; %其中x0为迭代初值,eps为精度要求N为最大迭代步数,con用来记录结果是否收敛 for i=1:N; f=subs(fun(x0),{'x1' 'x2'},{x0(1) x0(2)}); df=subs(dfun(x0),{'x1' 'x2'},{x0(1) x0(2)}); x=x0-f./df; for j=1:length(x0); il(i,j)=x(j); end if norm(x-x0)<eps con=1; break; end x0=x; end %以下是将迭代过程写入txt文档文件名为iteration.txt fid=fopen('iteration.txt','w'); fprintf(fid,'iteration'); for j=1:length(x0) fprintf(fid,' x%d',j); end for j=1:i fprintf(fid,'\n%6d ',j); for k=1:length(x0) fprintf(fid,' %10.6f',il(j,k)); end end if con==1 fprintf(fid,'\n计算结果收敛！'); end if con==0 fprintf(fid,'\n迭代步数过多可能不收敛！'); end fclose(fid); 在matlab中命令窗口中输入： newton([0.8 0.8],0.00001,20) 结果出现： ??? Error using ==> rdivide Matrix dimensions must agree. Error in ==> newton at 7 x=x0-f./df; 检查了几遍程序没发现问题出现在哪，请高手指点  |
回复此楼» 猜你喜欢
国自科面上基金字体
已经有6人回复
-
304求调剂
已经有3人回复
-
304求调剂
已经有6人回复
-
材料与化工一志愿南昌大学327求调剂推荐
已经有8人回复
-
0703 物理化学调剂
已经有5人回复
-
化学工程321分求调剂
已经有7人回复
-
0703化学调剂 290分有科研经历，论文在投
已经有7人回复
-
308求调剂
已经有3人回复
-
本人考085602 化学工程 专硕
已经有10人回复
-
焦虑
已经有9人回复
» 本主题相关价值贴推荐，对您同样有帮助:
- 求助如何用matlab求非线性方程组的解
已经有8人回复
- 求解决matlab解非线性方程组出现的问题
已经有6人回复
- MATLAB解非线性方程组
已经有9人回复
- 解线性方程组的时候如何取限定条件使解是正解？
已经有6人回复
- matlab 解非线性偏微分方程组
已经有6人回复
- 病态方程组求解问题，望指导啊
已经有7人回复
- 求用牛顿迭代法解非线性方程组的matlab程序
已经有4人回复
- 非线性二阶微分方程组求解，matlab
已经有12人回复
- 有限差分法化之后化简为线性方程组
已经有50人回复
- matlab 非线性微分方程求解
已经有3人回复
- 请教一个matlab求解非线性方程组的问题
已经有9人回复
- 用matlab求解非线性方程组说无解，一定是方程组本身无解，还是有可能程序有问题呢？
已经有11人回复
- matlab如何求解一个非线性微分方程组
已经有8人回复
- 用fortran程序遗传算法解非线性方程组
已经有7人回复
- 牛顿迭代法解含对数的非线性方程组
已经有11人回复
- MATLAB求解非线性方程组
已经有5人回复
- 求Matlab高手解决线性方程组的迭代求解问题
已经有17人回复
- Matlab：Newton-Raphson 解线非线性性方程组的问题
已经有12人回复
- matlab求解非线性方程组，错误提示怎么解决
已经有5人回复
- 求助matlab---fsolve解非线性方程组
已经有6人回复
- matlab求解非线性方程组
已经有16人回复
- 求一个mathmetica的程序
已经有10人回复
- 【求助】求教matlab解非线性方程组
已经有9人回复
感谢参与，应助指数 +1
|
本帖内容被屏蔽 |
2楼2015-04-20 15:18:07
3楼2015-04-20 15:33:42
|
本帖内容被屏蔽 |
4楼2015-04-20 15:42:05
|
按着你说的可能出现的问题，我对程序又改了改，发现还是没结果，然后捋了捋思路，问题可能出现在如果我按着上述程序计算，即使能得出结果，对应的x1和x2也应该是数，而不是数组。 我把我要解决的问题给你描述一下： 以下参数已知： T=318.2; rhomcA=0.34; CA=3.26; EAA=64; EAs=-1385; rhomcB=0.98; CB=4.53; EBB=82; EBs=-1690; rhob=[0.0047 0.0090 0.0177 0.0264 0.0367 0.0466 0.0566 0.0674 0.0787 0.0890]'; yA=[0.89 0.886 0.884 0.883 0.877 0.871 0.866 0.861 0.856 0.852]'; yB=[0.11 0.114 0.116 0.117 0.123 0.129 0.134 0.139 0.144 0.148]'; 已知公式如下： log(rhomcA*x1*(1-yA*rhob/rhomcA-yB*rhob/rhomcB)/(yA*rhob*(1-x1-x2)))+EAA/T*(7*x1-8*yA*rhob/rhomcA)+(EAA*EBB)^0.5/T*(7*x2-8*yB*rhob/rhomcB)+EAs/T=0 log(rhomcB*x2*(1-yA*rhob/rhomcA-yB*rhob/rhomcB)/(yB*rhob*(1-x1-x2)))+EBB/T*(7*x2-8*yB*rhob/rhomcB)+(EAA*EBB)^0.5/T*(7*x1-8*yA*rhob/rhomcA)+EBs/T=0 想使用newton迭代法求x1和x2，要求的x1和x2应该分别为10×1数据 请大神帮帮忙呗，这个问题困扰我一周了，matlab初学者，一遇到循环就蒙  |
5楼2015-04-20 16:33:16
dingd
铁杆木虫 (职业作家)
- 应助: 1641 (讲师)
- 金币: 15037.3
- 散金: 101
- 红花: 234
- 帖子: 3410
- 在线: 1223.7小时
- 虫号: 291104
- 注册: 2006-10-28
【答案】应助回帖
感谢参与，应助指数 +1
|
Matlab感觉好复杂的样子。1stOpt求解看看： Loop constant rhob Loop constant ya Loop constant yb x1 x2 0.0047 0.89 0.11 0.349788429554512 0.0364739534111791 0.009 0.886 0.114 0.460684483692858 0.0489828447137933 0.0177 0.884 0.116 0.580741252753925 0.0617304175836032 0.0264 0.883 0.117 0.649820821824047 0.0690991828120868 0.0367 0.877 0.123 0.699996476771262 0.0783646455649231 0.0466 0.871 0.129 0.732280709438354 0.0863910296157365 0.0566 0.866 0.134 0.756187405093584 0.0931836288289366 0.0674 0.861 0.139 0.775229538588392 0.0997840413222308 0.0787 0.856 0.144 0.789855870165985 0.106177522100049 0.089 0.852 0.148 0.800097262910873 0.111365586073744 |
6楼2015-04-20 17:32:40
cooooldog
铁杆木虫 (著名写手)
ส็็็
- 应助: 237 (大学生)
- 金币: 6101.9
- 散金: 1114
- 红花: 39
- 帖子: 1380
- 在线: 553.8小时
- 虫号: 506699
- 注册: 2008-02-18
- 专业: 模式识别
7楼2015-04-20 19:55:59
8楼2015-04-20 21:02:48
9楼2015-04-20 21:11:40
