如何定义函数文件 数学软件与数学实验 教学课件.ppt

如何定义函数文件 如何定义函数文件 举例说明 举例说明（单个方程） 举例说明（方程组） 思考 说明 例 高阶常微分方程 参数传递 * * [T,Y] = ode45(odefun,tspan,y0) [T,Y] = ode23(odefun,tspan,y0) odefun？ ode45、ode23 能解什么样的 ODE？ Oh，I see 微分方程组 当 是函数向量时呢？ 那么odefun就是 odefun [T,Y] = ode45(odefun,tspan,y0) [T,Y] = ode23(odefun,tspan,y0) ode45、ode23 等函数可用于求解显式常微分方程 当 是向量函数时，所对应的方程即为微分方程组 odefun fun=inline(-2*y+2*x^2+2*x,x,y); [x,y]=ode23(fun,[0,0.5],1); 例：求初值问题 的数值解。 解法一：使用 inline 定义微分方程 odefun odefun 为方程右端项 f(t,y) 可以用 inline 定义（只适合于单个方程的情形） 通过函数文件定义，然后用函数句柄调用（适合所有情形） 注：自变量必须在前面，因变量在后面！ function dy = myfun1(x,y) dy = -2*y+2*x^2+2*x; 解法二：通过函数文件定义微分方程 odefun 1、先编写函数文件 myfun1.m clear; [x,y]=ode23(@myfun1,[0,0.5],1); 2、编写主文件 main1.m 或直接在 Matlab 命令窗口输入上面的语句。 解：此时只能通过函数文件定义微分方程 odefun 例：求 , , 的数值解。 function dy = myfun2(t,y) dy = zeros(3,1); % dy must be a column vector! dy(1) = y(2) * y(3); dy(2) = -y(1) * y(3); dy(3) = -0.51 * y(1) * y(2); 1、先编写函数文件 myfun2.m clear; [T,Y]=ode45(@myfun2,[0,12],[0,1,1]); 2、编写主文件 main2.m dy = [ y(2)*y(3); -y(1)*y(3); ... -0.51*y(1)*y(2)]; function dy = myfun2(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); 1、函数文件 myfun2.m 能不能写成下面形式？ function dy = myfun2(t,x,y,z) dy = zeros(3,1); dy(1) = y * z; dy(2) = -x * z; dy(3) = -0.51 * x * y; X odefun 变量属性必须一一对应！ function dy = myfun2(t,y) 如果是常微分方程组，y 就是列向量！ dy 必须是列向量，长度为方程组的个数，通常与y的长度相同！ 函数中的输入参数和输出参数是形参，名字可以任意取，但必须满足上述条件。即输入参数有两个，第一个表示自变量，第二个是由因变量组成的列向量，输出参数必须是列向量。 function dy = myfun3(t,y) dy = zeros(2,1); dy(1) = y(2) + t; dy(2) = t - 2; 例：解初值问题： ， ， function out = myfun3(t,y) out = [y(2) + t; t – 2]; function yprime = myfun3(x,y) yprime = [y(2) + x; x – 2]; clear; [T,Y]=ode45(@myfun3,[0,10],[1,1]); 2、主文件 main3.m 1、函数文件 myfun3.m 高阶常微分方程 例： Van der Pol 初值问题 令 ，则原方程可化为 一阶常微分方程组 变量替换 化为 参数怎么处理？ 用全局变量传递 1、函数文件 verderpol.m function xprime = verderpol(t,x) global mu;