function Cha2demo3 P = 9.33; % atm T = 300.2; % K n = 2; % mol a = 4.17; b = 0.0371; R = 0.08206; V0 = n*R*T/P % 5.2807 [V,fval] = fzero(@PVTeq,V0,[],P,T,n,a,b,R) % ------------------------------------------function f = PVTeq(V,P,T,n,a,b,R) f = (P + a*n^2/V^2) * (V-n*b) - n*R*T
例题：


在945.36kPa（9.33atm）、300.2K时，容器中充以 2mol氮气，试求容器体积。已知此状态下氮气的P-V-T 关系符合范德华方程，其范德华常数为a＝ 4.17atm•L/mol2，b＝0.0371L/mol。 数学模型：范德华方程为：
an2 f (V ) ( p 2 )(V nb) nRT 0 V
脚本文件
非线性方程求解函数fzero

fzero


对于一般的含单个未知数的超越方程，可以采用 fzero函数求解 fzero函数结合使用二分法、割线法和可逆二次内插 法
函数fzero

[x,fval,exitflag,output] = fzero(fun,x0,options, p1, p2, ...) 此函数的作用求函数fun在x0附件的零值点，x0是 标量
P [a0 , a1 ,, an1 , an ]
这样就把多项式问题转化为向量问题
函数roots


r = roots(P)，用于求解多项式的根 其中，行向量P的元素是多项式的系数，按多 项式次数降序排列 如果P中含有n+1个元素，则多项式为n次 roots可以获得多项式的所有根
例题：
1. Fun函数如何编写？ 2. X0如何选取？
说明
1) 采用句柄函数定义函数
function Cha2demo1 x=fzero(@fun,1) function y=fun(x) y=x^3-2*sin(x);
2) 初值的选择对于解有影响，不同的初值可能 获得不同的解

可以根据感兴趣的解的区间确定初值范围 可以作出函数在一定范围内的曲线，直观的确定解的大致 范围
课堂练习


某蒸馏釜的操作压力位106.7kpa，其中溶液含苯摩尔分 数为0.2，甲苯为0.8，求此溶液的泡点和平衡的气相组 成。 苯甲苯溶液可以看做理想溶液，组分蒸汽压为：a为苯， b为甲苯
1211 log Pa 6.301 t 220.8 1345 log Pb 6.080 t 219.5
function y=fun(x) y(1)=sin(x(1))+x(2)^2+log(x(3))-7; y(2)=3*x(1)+2^x(2)-x(3)^3+1; y(3)=x(1)+x(2)+x(3)-5;
x= 0.5991
2.3959
2.0050
fsolve函数的应用
在铜管内在1 atm下将异丙醇加热到400℃。已知铜是生产丙酮和丙 醛的催化剂，或许还有某些异丙醇异构化为正丙醇。这三种产物的 生成可用如下三个独立反应表示： iC3H7OH(IP) = n C3H7OH(NP) K1 = 0.064 iC3H7OH(IP) = (CH3)CO(AC)+H2 K2 = 0.076 iC3H7OH(IP) = C２H５ＣHO(PR) +H2 K3 = 0.00012 问，反应达到平衡时，产物各自的mol百分含量是多少。 数学模型：各反应的化学平衡方程如下
function Cha2demo7 x0 = [0.05 0.2 0.01]; x = fsolve(@EquiC3,x0); function f = EquiC3(x) f1 = x(1)-0.064*(1-x(1)-x(2)-x(3)); f2 = x(2)*(x(2)+x(3))-0.076*(1-x(1)-x(2)-x(3))*(1+x(2)+x(3)); f3 = x(3)*(x(2)+x(3))-0.00012*(1-x(1)-x(2)-x(3))*(1+x(2)+x(3)); f = [f1 f2 f3];
roots，求多项式的值，polyval；多项式微分，polyder； 多项式拟合，polyfit；多项式乘法，conv；多项式除法， deconv；
例题：


在945.36kPa（9.33atm）、300.2K时，容器中充以 2mol氮气，试求容器体积。已知此状态下氮气的P-V-T 关系符合范德华方程，其范德华常数为a＝ 4.17atm•L/mol2，b＝0.0371L/mol。 数学模型：范德华方程变形可得

求方程 x x 1 的根
3 2
>>c = [1 -1 0 -1]; >>r = roots(c) r= 1.4656 -0.2328 + 0.7926i -0.2328 - 0.7926i
>>polyval(c, r(1)) ans = -2.5535e-015
Matlab提供了多种多项式计算函数，如多项式求根函数
第二讲 误差与非线性方程求解
数值计算的误差
模型误差 误 差 的 分 类 观测误差 截断误差 舍入误差 计算时只截 取有限项 计算机对所 储存的数据 位数有限制
ln( x 1) x 1 2 1 3 1 4 1 x x x (1) n 1 x n 2 3 4 n
ε
( x) x * x
定义待求解方程时，必须首先将方程组变换成 F(X)=0的形式！
例题10：
sin x y 2 ln z 7 y 3 3 x 2 z 1 0 x y z 5
在命令窗口输入： x0=[1 1 1]; x=fsolve(@fun,x0)
y (1) sin x(1) x(2) 2 ln x(3) 7 0 x (2) 3 y (2) 3x(1) 2 x(3) 1 0 y (3) x(1) x(2) x(3) 5 0
函数fsolve
与fzero函数只能求解单个方程的根不同，fsolve 函数可求解非线性方程组的解。其算法采用的是 最小二乘法。 调用格式： [x,fval,exitflag,output,jacobian] = fsolve(fun,x0,options, p1, p2, ...) 输入输出变量的意义同fzero函数 输出变量中的jacobian为函数fun在x处的 Jacobian矩阵
3) fzero不能获得多项式的多重根，尤其是复数 根。而roots函数求解，则可获得所有根
例题7：

计算以下方程的根 1) 求sinx在3附近的零点; 2) 求cosx在[1,2]范围内的零点； 3 3) x 2 x 5 0 3 4) x 2 sin x 0
本例较简单，可直接在命令窗口输入命令求解： 1) fzero(@sin,3) 2) fzero(@cos,[1,2]) 3) fzero(@(x) x^3-2*x-5,1); roots([1 0 -2 -5]) 4) fzero(@(x) x^3-2*sin(x),1)
绝对误差 误差的分类 相对误差
( x)
r ( x) ( x * x) / x r ( x ) ( x * x) / x *
r ( x)
有效数字
定义：设X*是数X的近似值，如果X*的绝对误差限是他 的某一位数位的半个单位，并且从X*左起第一个非零数 字到该数位共有N位，则称这个N个数字为X*的有效数字， 也称X*近似X时有N位有效数字
an2 f (V ) ( p 2 )(V nb) nRT 0 V
pV 3 ( pnb nRT )V 2 an2V an3b 0

这是关于V的三次方程，可以由roots
P = 9.33; % atm T = 300.2; % K n = 2; % mol a = 4.17; b = 0.0371; R = 0.08206; Eq=[P,-(P*n*b+n*R*T),a*n^2,a*n^3*b]; roots(Eq)
数学模型：范德华方程变形可得关于V的非线性 方程
an2 f (V ) ( p 2 )(V nb) nRT 0 V
非线性方程

非线性方程包括：高次代数方程、超越方程(具有未知量 的对数函数、指数函数、三角函数、反三角函数等的方程)及其它 们的组合 与线性方程相比，非线性方程求解问题无论从理论上 还是从计算公式上都要复杂得多 对于高次代数方程，当次数>4时，则没有通解公式可 用，对于超越方程既不知有几个根，也没有同样的求 解方式。实际上，对于n≥3代数方程以及超越方程都 采用数值方法求近似根。
x1 0.064 1 x1 x 2 x 3
x2 ( x2 x3 ) 0.076 (1 x1 x 2 x 3 )(1 x 2 x 3 )
x3 ( x2 x3 ) 0.00012 (1 x1 x 2 x 3 )(1 x 2 x 3 )


x 所求解 fval 函数在解x处的值 exitflag 程序结束情况: >0，程序收敛于解；<0，程序没有收敛； ＝0，计算达到了最大次数 output 是一个结构体，提供程序运行的信息； output.iterations，迭代次数；output.functions，函数fun的计算 次数；output.algorithm，使用的算法 options 选项，可用optimset函数设定选项的新值 fun可以是函数句柄或匿名函数。 P1,P2是传递的参数