《matlab数学建模微分法与优化问题》由会员分享，可在线阅读，更多相关《matlab数学建模微分法与优化问题（45页珍藏版）》请在金锄头文库上搜索。
1、 M文件M文件建立方法：1. 在MATLAB中，点:FileNew M-file2. 在编辑窗口中输入程序内容3. 点File Save，存盘，M文件名必须与函数名一 致.M文件有两种类型：脚本M文件和函数M文件 l脚本M文件就是命令的简单叠加。 l有时为了研究某一个函数的各种性态，需要为MATLAB定义 新函数，为此必须编写函数文件。函数文件是文件名后缀为M 的文件，这类文件的第一行必须是一特殊字符function开始，格 式为：function 因变量名=函数名（自变量名）函数值的获得必须通过具体的运算实现，并赋给因变量。M文件名必须与函数名一致.例：定义函数 f(x1,x2)=funct
2、ion f=fun(x) f=100*(x(2)-x(1)2)2+(1-x(1)2；1.建立M文件：fun.m2. 可以直接使用函数fun.m例如：计算 f(1,2), 只需在MATLAB命令窗口键入命令：x=1 2fun(x)关系与逻辑运算1. 关系操作符2. 逻辑运算符for循环 允许一组命令以固定的和预定的次数重复for x=arraycommands；end在for和end语句之间的命令串commands按数组(array)中的每一列执行 一次. 在每一次迭代中，x被指定为数组的下一列，即在第n次循环中， x=array(：，n)控制流MATLAB提供三种决策或控制流结构：for循环、
3、while循环、if-else-end结构. 这些结构经常包含大量的MATLAB命令，故经常出现在MATLAB程序中 ，而不是直接加在MATLAB提示符下. 例 对n=1,2,10,求xn= 的值.while expressioncommands；end只要在表达式(expression)里的所有元素为真，就执行while 和end语句之间的命令串commands. While循环与for循环以固定次数求一组命令相反，while循环以不定的次 数求一组语句的值.例 设银行年利率为11.25%.将10000元钱存入银行，问 多长时间会连本带利翻一番？if-else-end结构 （1）一个选择的一
4、般形式是：if (expression)commands；end 如果在表达式(expression)里的所有元素为真，就执 行if和end语句之间的命令串commands. （2）多个选择的一般形式是：if （expression1）commands1；else if （expression2）commands2；else if （expression3）commands3；else if elsecommands； endendendendcontinue命令的作用是结束本次循环，即跳过循环体中 尚未执行的语句，接着进行下一次是否执行循环的判断。 (当条件满足时，continue语句才被调
5、用，跳出后面的语 句，循环继续，条件不满足则继续执行后面的程序。)break命令的作用是终止本次循环，跳出最内层循环，即 不必等到循环的结束而是根据条件退出循环，用法与 continue类似，常常和条件语句合用来强制终止循环。( 当条件满足时， break语句才被调用，系统将在循环尚未 结束时跳出当前循环，条件不满足则继续执行后面的程 序。)return命令可使正常运行的函数正常退出，并返回调用它 的函数继续运行，经常用于函数的末尾以正常结束函数的 运行，也可用在某条件满足时强行结束执行该函数。input命令用来提示用户从键盘输入数据、字符串或表达式 ，并接受输入值。 r=input(prom
6、pt)在屏幕上显示提示信息prompt，等待用户 输入，并将值赋给r，直接回车则输出空矩阵。 r=input(prompt,s)返回的字符串作为文本变量而不作为 变量名或者数值。pause命令用于暂时终止程序的运行，等待用户按任意键继 续进行。 pause暂停M文件的执行，按任意键继续 pause(n)终止程序n秒，n是任意实数 pause on允许后面的pause命令起作用 pause off使后面的pause命令不起作用代数方程的求解一般多项式方程的根可为实数，也可为复数。 可用MATLAB符号工具箱中的solve( )函数。 最简调用格式：S=solve(eqn1,eqn2,eqnn)
7、（返回一个结构体型变量S，如S.x表示方程的根。） 直接得出根： （变量返回到MATLAB工作空间）x,=solve(eqn1,eqn2,eqnn) 同上，并指定变量x,=solve(eqn1,eqn2,eqnn，x,)例：x,y,z=solve(x+3*y3+2*z2=1/2 , x2+3*y+z3 =2 ,x3+2*z+2*y2=2/4) ;函数 fzero 非线性方程的数值解格式 x = fzero (fun,x0) %用fun定义表达式f(x)，x0为初始 解。x = fzero (fun,x0,options) 说明 该函数采用数值解求方程f(x)=0的一个根。例: 求 的根fun=
8、x3-2*x-5; z=fzero(fun,2) %初始估计值为2非线性方程组数值解求解语句x=fsolve(Fun, x0)x=fsolve(Fun, x0,opt) opt 求解控制参数，结构体数据，有20多个选项。 获得默认的常用变量opt=optimset; 可以用这个命令显示出来。 用这两种方法修改参数（解误差控制量）opt.Tolx=1e-10; 或 set(opt.Tolx, 1e-10) 也可以调用optimset()函数来完成。例：OPT=optimset; OPT.Tolx=1e-10;x= fsolve(my2deq,1; 2,OPT)微分方程的解析解求微分方程（组）解析
9、解的命令:dsolve(方程1,方程2,方程n,初始条件,自变量)结 果：u = tg(t+c)y=dsolve(D2y+4*Dy+29*y=0,y(0)=0,Dy(0)=15,x)x,y,z=dsolve(Dx=2*x-3*y+3*z,Dy=4*x-5*y+3*z,Dz=4*x-4*y+2*z, t)；x=simple(x) % 将x化简y=simple(y)z=simple(z)微分方程的数值解例如：用差商代替导数 欧拉法若步长h较小，则有故有公式：使用泰勒公式，有龙格-库塔法、线性多步法等方法数值公式的精度当一个数值公式的截断误差可表示为O（hk+1）（其 中k为正整数，h为步长）时，称
10、它是一个k阶公式k越大，则数值公式的精度越高欧拉法是一阶公式，改进的欧拉法是二阶公式龙格-库塔法有二阶公式和四阶公式线性多步法有四阶亚当斯外插公式和内插公式t，x=solver(f,ts,x0,options)ode45 ode23 ode113 ode15s ode23s由待解 方程写 成的M 文件名ts=t0, tf，t0、tf为自 变量的初 值和终值函数 的初 值ode23：组合的2/3阶龙格库塔费尔贝格算法 ode45：运用组合的4/5阶龙格库塔费尔贝格算法自变 量值函数 值用于设定误差限(缺省时设定相对误差10-3, 绝对误差10-6), 命令为：options=odeset（rel
11、tol,rt,abstol,at）, rt，at：分别为设定的相对误差和绝对误差1在解含n个未知数的方程组时，x0和x均为n维向量，M文件中的待 解方程组应以x的分量形式写出2使用MATLAB软件求数值解时，高阶微分方程必须等价地变换 成一阶微分方程组注意:令 y1=x，y2=y11建立M文件vdp1000m如下：function dy=vdp1000(t,y)dy=zeros(2,1);dy(1)=y(2);dy(2)=1000*(1-y(1)2)*y(2)-y(1); 2取t0=0，tf=3000，输入命令：T,Y=ode45(vdp1000,0 3000,2 0); plot(T,Y(:
12、,1),-)3结果如图1建立M文件rigidm如下：function dy=rigid(t,y)dy=zeros(3,1);dy(1)=y(2)*y(3);dy(2)=-y(1)*y(3);dy(3)=-0.51*y(1)*y(2);2取t0=0，tf=12，输入命令：T,Y=ode45(rigid,0 12,0 1 1);plot(T,Y(:,1),-,T,Y(:,2),*,T,Y(:,3),+)3结果如图图中，y1的图形为实线，y2的图形为“*”线，y3的图形为“+”线导弹追踪问题设位于坐标原点的甲舰向位于x轴上点A(1, 0)处的乙舰发 射导弹，导弹头始终对准乙舰如果乙舰以最大的速度v0
13、( 常数)沿平行于y轴的直线行驶，导弹的速度是5v0，求导弹 运行的曲线方程乙舰行驶多远时，导弹将它击中？导弹追踪问题解法一（解析法）由(1),(2)消去t, 整理得模型:解法二(数值解法 )1建立M文件eq1mfunction dy=eq1(x,y)dy=zeros(2,1);dy(1)=y(2);dy(2)=1/5*sqrt(1+y(1)2)/(1-x); 2 x0=0;xf=0.9999;x,y=ode15s(eq1,x0 xf,0 0);plot(x,y(:,1),b.) hold ony=0:0.01:2;plot(1,y,b*)结论: 导弹大致在（1，02）处击中乙舰.令y1=y,
14、 y2=y1，将方程（3）化为一阶微分方程组解法三(建立参数方程求数值解)设时刻t乙舰的坐标为(X(t)，Y(t)，导弹的坐标为(x(t)，y(t)3因乙舰以速度v0沿直线x=1运动，设v0=1，则w=5， X=1，Y=t.4 解导弹运动轨迹的参数方程建立M文件eq2m如下：function dy=eq2(t,y)dy=zeros(2,1);dy(1)=5*(1-y(1)/sqrt(1-y(1)2+(t-y(2)2);dy(2)=5*(t-y(2)/sqrt(1-y(1)2+(t-y(2)2);取t0=0，tf=2：t,y=ode45(eq2,0 2,0 0);plot(y(:,1),y(:,
15、2),*)用MATLAB优化工具箱解线性规划min z=cX 1. 模型：命令：x=linprog（c, A, b） 2. 模型：min z=cX 命令：x=linprog（c，A，b，Aeq,beq）注意：若没有不等式： 存在，则令A= ，b= .3. 模型：min z=cX VLBXVUB4. 命令：x,fval=linprog()返回最优解及处的目标函数值fval.命令：x=linprog（c，A，b，Aeq,beq, VLB，VUB）解 编写M文件xxgh.m如下：c=-0.4 -0.28 -0.32 -0.72 -0.64 -0.6;A=0.01 0.01 0.01 0.03 0.03 0.03;0.02 0 0 0.05 0 0;0 0.02 0 0 0.05 0;0 0 0.03 0 0 0.08;b=850;700;100;900;Aeq=; beq=;vlb=0;0;0;0;0;0; vub=;x
