decic

为 ode15i 计算一致的初始条件

语法

[y0_new,yp0_new] = decic(odefun,t0,y0,fixed_y0,yp0,fixed_yp0)

[y0_new,yp0_new] = decic(odefun,t0,y0,fixed_y0,yp0,fixed_yp0,options)

[y0_new,yp0_new,resnrm] = decic(___)

说明

[y0_new,yp0_new] = decic(odefun,t0,y0,fixed_y0,yp0,fixed_yp0) 使用 y0 和 yp0 作为完全隐函数 odefun 的初始条件的估计值，保留由 fixed_y0 和 fixed_yp0 指定的分量作为固定分量，然后计算非固定分量的值。结果将获得一套完整的一致初始条件。新值 yo_new 和 yp0_new 满足 odefun(t0,y0_new,yp0_new) = 0，适合用作 ode15i 的初始条件。

示例

[y0_new,yp0_new] = decic(odefun,t0,y0,fixed_y0,yp0,fixed_yp0,options) 还使用 options 结构体 options 为 AbsTol 和 RelTol 指定值。可以使用 odeset 创建 options 结构体。

[y0_new,yp0_new,resnrm] = decic(___) 返回 odefun(t0,y0_new,yp0_new) 的范数作为 resnrm。如果范数看起来异常大，则使用 options 减小相对误差容限 RelTol，其默认值为 1e-3。

示例

全部折叠

为隐式方程计算一致初始条件

打开实时脚本

考虑隐式方程组

$\begin{array}{cc} 0 & = 2 y_{1}^{'} - y_{2} \\ 0 & = y_{1} + y_{2} \end{array}$

这些方程很简单，很容易读取变量的一致初始条件。例如，如果固定 $y_{1} = 1$ ，则根据第二个方程 $y_{2} = - 1$ ，根据第一个方程 $y_{1}^{'} = - 1 / 2$ 。由于 $y_{1}$ 、 $y_{1}^{'}$ 和 $y_{2}$ 的这些值满足方程，因此它们是一致的。

通过使用 decic 计算方程的一致初始条件并固定值 $y_{1} = 1$ ，可以确认这些值。使用 y0 = [1 0] 和 yp0 = [0 0] 的估计值，它们不满足方程，因此不一致。

odefun = @(t,y,yp) [2*yp(1)-y(2); y(1)+y(2)];
t0 = 0;
y0 = [1 0];
yp0 = [0 0];
[y0,yp0] = decic(odefun,t0,y0,[1 0],yp0,[])

求解 Weissinger 隐式 ODE

打开实时脚本

计算一致的初始条件，并用 ode15i 求解隐式 ODE。

Weissinger 方程是

${t y}^{2} {(y^{'})}^{3} - y^{3} {(y^{'})}^{2} + t (t^{2} + 1) y^{'} - t^{2} y = 0$ .

由于该方程采用一般形式 $f (t, y, y^{'}) = 0$ ，您可以使用 ode15i 函数来求解隐式微分方程。

编写方程代码

要以适合 ode15i 的形式对方程进行编码，您需要编写一个函数，该函数使用 $t$ 、 $y$ 和 $y^{'}$ 的输入并返回方程的残差值。可以使用函数 @weissinger 对此方程进行编码。查看函数文件。

type weissinger

function res = weissinger(t,y,yp)
%WEISSINGER  Evaluate the residual of the Weissinger implicit ODE
%
%   See also ODE15I.

%   Jacek Kierzenka and Lawrence F. Shampine
%   Copyright 1984-2014 The MathWorks, Inc.

res = t*y^2 * yp^3 - y^3 * yp^2 + t*(t^2 + 1)*yp - t^2 * y;

计算一致的初始条件

ode15i 求解器需要一致的初始条件，即提供给求解器的初始条件必须满足

$f (t_{0}, y, y^{'}) = 0$ .

由于可能会提供不一致的初始条件，并且 ode15i 不会检查一致性，因此建议您使用辅助函数 decic 来计算这些条件。decic 会固定一些指定的变量，并为不固定的变量计算一致的初始值。

在本例中，固定初始值 $y (t_{0}) = \sqrt{\frac{3}{2}}$ ，并让 decic 从 $y^{'} (t_{0}) = 0$ 的初始估计值开始，计算导数 $y^{'} (t_{0})$ 的一致初始值。

t0 = 1;
y0 = sqrt(3/2);
yp0 = 0;
[y0,yp0] = decic(@weissinger,t0,y0,1,yp0,0)

y0 = 
1.2247

yp0 = 
0.8165

求解方程

使用 decic 加上 ode15i 返回的一致初始条件，在时间区间 $[1 10]$ 上求解 ODE。

[t,y] = ode15i(@weissinger,[1 10],y0,yp0);

绘制结果

此 ODE 的精确解为

$y (t) = \sqrt{t^{2} + \frac{1}{2}}$ .

绘制 ode15i 计算的数值解 y 对解析解 ytrue 的图。

ytrue = sqrt(t.^2 + 0.5);
plot(t,y,'*',t,ytrue,'-o')
legend('ode15i', 'exact')

Figure contains an axes object. The axes object contains 2 objects of type line. One or more of the lines displays its values using only markers These objects represent ode15i, exact.

输入参数

全部折叠

`odefun` — 要求解的函数
函数句柄

要求解的函数，指定为函数句柄，用于定义要求积分的函数。odefun 代表您要使用 ode15i 求解的隐式微分方程组。

对于标量 t 以及列向量 y 和 yp 来说，函数 f = odefun(t,y,yp) 必须返回数据类型为 single 或 double 的列向量 f，这些类型与 $f (t, y, y')$ 相对应。odefun 必须同时接受 t、y 和 yp 这三个输入参量，即使其中某个参量未在函数中使用也是如此。

例如，要求解 $y' - y = 0$ ，请使用此函数。

function f = odefun(t,y,yp)
f = yp - y;

对于方程组，odefun 的输出为向量。每个方程变成解向量中的一个元素。例如，要求解

$\begin{array}{l} y'_{1} - y_{2} = 0 \\ y'_{2} + 1 = 0, \end{array}$

，请使用以下函数。

function dy = odefun(t,y,yp)
dy = zeros(2,1);
dy(1) = yp(1)-y(2);
dy(2) = yp(2)+1;

有关如何为函数 odefun 提供其他参数的信息，请参阅参数化函数。

示例: @myFcn

数据类型: function_handle

`t0` — 初始时间
标量

初始时间，指定为标量。decic 使用初始时间计算满足 odefun(t0,y0_new,yp0_new) = 0 的一致初始条件。

数据类型: single | double

`y0` — `y` 分量的初始估计值
向量

y 分量的初始估计值，指定为向量。y0 中的每个元素为 odefun 定义的方程组中的一个因变量 $y_{n}$ 指定一个初始条件。

数据类型: single | double

`fixed_y0` — 要保持固定的 `y` 分量
向量 1 和 0 | `[]`

要保持固定的 y 分量，指定为向量 1 和 0，或者指定为 []。

如果 y0(i) 的估计值中不允许进行任何更改，则设置 fixed_y0(i) = 1。
如果可以更改任何条目，则设置 fixed_y0 = []。

您最多只能固定 length(yp0) 个分量。根据具体问题的不同，可能无法始终固定 y0 或 yp0 的某些分量。最好不要固定不必要的多余分量。

`yp0` — `y'` 分量的初始估计值
向量

y' 分量的初始估计值，指定为向量。yp0 中的每个元素为 odefun 定义的方程组中的一个微分因变量 $y'_{n}$ 指定一个初始条件。

数据类型: single | double

`fixed_yp0` — 要保持固定的 `y'` 分量
向量 1 和 0 | `[]`

要保持固定的 y' 分量，指定为向量 1 和 0，或者指定为 []。

如果 yp0(i) 的估计值中不允许进行任何更改，则设置 fixed_yp0(i) = 1。
如果可以更改任何条目，则设置 fixed_yp0 = []。

您最多只能固定 length(yp0) 个分量。根据具体问题的不同，可能无法始终固定 y0 或 yp0 的某些分量。最好不要固定不必要的多余分量。

`options` — Options 结构体
结构体数组

Options 结构体，指定为结构体数组。使用 odeset 函数创建或修改 options 结构体。可与 decic 函数一起使用的相关选项包括 RelTol 和 AbsTol，它们控制用于计算初始条件的误差阈值。

示例: options = odeset('RelTol',1e-5)

数据类型: struct

输出参量

全部折叠

`y0_new` — `y0` 的一致初始条件
向量

y0 的一致初始条件，以向量形式返回。如果 resnrm 的值较小，则 yo_new 和 yp0_new 满足 odefun(t0,y0_new,yp0_new) = 0，适合用作 ode15i 的初始条件。

`yp0_new` — `yp0` 的一致初始条件
向量

yp0 的一致初始条件，以向量形式返回。如果 resnrm 的值较小，则 yo_new 和 yp0_new 满足 odefun(t0,y0_new,yp0_new) = 0，适合用作 ode15i 的初始条件。

`resnrm` — 残差范数
向量

残差范数，以向量形式返回。resnrm 是 odefun(t0,y0_new,yp0_new) 的范数。

如果 resnrm 的值较小，则表示 decic 已成功计算满足 odefun(t0,y0_new,yp0_new) = 0 的一致初始条件。
如果 resnrm 的值很大，请使用 options 输入尝试调整误差阈值 RelTol 和 AbsTol。

提示

在使用 ode15i 求解之前，ihb1dae 和 iburgersode 示例文件使用 decic 计算一致初始条件。可以键入 edit ihb1dae 或 edit iburgersode 查看代码。
还可以使用 decic 计算通过 ode15s 或 ode23t 求解的 DAE 的一致初始条件。为此，请按以下步骤操作。
1. 采用完全隐式的形式 f(t,y,y') = 0 重写方程组。
2. 调用 decic 计算方程的一致初始条件。
3. 指定 y0_new 作为调用求解器的初始条件，并指定 yp_new 作为 odeset 的 InitialSlope 选项的值。

版本历史记录

在 R2006a 之前推出

另请参阅

ode15i | odeget | odeset

decic

语法

说明

示例

为隐式方程计算一致初始条件

求解 Weissinger 隐式 ODE

输入参数

odefun — 要求解的函数 函数句柄

t0 — 初始时间 标量

y0 — y 分量的初始估计值 向量

fixed_y0 — 要保持固定的 y 分量 向量 1 和 0 | []

yp0 — y' 分量的初始估计值 向量

fixed_yp0 — 要保持固定的 y' 分量 向量 1 和 0 | []

options — Options 结构体 结构体数组

输出参量

y0_new — y0 的一致初始条件 向量

yp0_new — yp0 的一致初始条件 向量

resnrm — 残差范数 向量

提示

版本历史记录

另请参阅

`odefun` — 要求解的函数
函数句柄

`t0` — 初始时间
标量

`y0` — `y` 分量的初始估计值
向量

`fixed_y0` — 要保持固定的 `y` 分量
向量 1 和 0 | `[]`

`yp0` — `y'` 分量的初始估计值
向量

`fixed_yp0` — 要保持固定的 `y'` 分量
向量 1 和 0 | `[]`

`options` — Options 结构体
结构体数组

`y0_new` — `y0` 的一致初始条件
向量

`yp0_new` — `yp0` 的一致初始条件
向量

`resnrm` — 残差范数
向量