# decic

Compute consistent initial conditions for `ode15i`

## Syntax

``````[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(___)``````

## Description

example

``````[y0_new,yp0_new] = decic(odefun,t0,y0,fixed_y0,yp0,fixed_yp0)``` uses `y0` and `yp0` as guesses for the initial conditions of the fully implicit function `odefun`, holds the components specified by `fixed_y0` and `fixed_yp0` as fixed, then computes values for the nonfixed components. The result is a complete set of consistent initial conditions. The new values `yo_new` and `yp0_new` satisfy ```odefun(t0,y0_new,yp0_new) = 0``` and are suitable to be used as initial conditions with `ode15i`.```
``````[y0_new,yp0_new] = decic(odefun,t0,y0,fixed_y0,yp0,fixed_yp0,options)``` also uses the options structure `options` to specify values for `AbsTol` and `RelTol`. Create the options structure using `odeset`.```
``````[y0_new,yp0_new,resnrm] = decic(___)``` returns the norm of `odefun(t0,y0_new,yp0_new)` as `resnrm`. If the norm seems unduly large, then use `options` to decrease the relative error tolerance `RelTol`, which has a default value of `1e-3`.```

## Examples

collapse all

Consider the implicit system of equations

`$\begin{array}{cc}0& =2{y}_{1}^{\prime }-{y}_{2}\\ 0& ={y}_{1}+{y}_{2}\end{array}$`

These equations are straightforward enough that it is simple to read off consistent initial conditions for the variables. For example, if you fix ${y}_{1}=1$, then ${y}_{2}=-1$ according to the second equation and ${y}_{1}^{\prime }=-1/2$ according to the first equation. Since these values of ${y}_{1}$, ${y}_{1}^{\prime }$, and ${y}_{2}$ satisfy the equations, they are consistent.

Confirm these values by using `decic` to compute consistent initial conditions for the equations, fixing the value ${y}_{1}=1$. Use guesses of `y0 = [1 0]` and `yp0 = [0 0]`, which do not satisfy the equations and are thus inconsistent.

```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,[])```
```y0 = 2×1 1 -1 ```
```yp0 = 2×1 -0.5000 0 ```

Calculate consistent initial conditions and solve an implicit ODE with `ode15i`.

Weissinger's equation is

${\mathrm{ty}}^{2}{\left({\mathit{y}}^{\prime }\right)}^{3}-{\mathit{y}}^{3}{\left({\mathit{y}}^{\prime }\right)}^{2}+\mathit{t}\left({\mathit{t}}^{2}+1\right){\mathit{y}}^{\prime }-{\mathit{t}}^{2}\mathit{y}=0$.

Since the equation is in the generic form $\mathit{f}\left(\mathit{t},\mathit{y},{\mathit{y}}^{\prime }\right)=0$, you can use the `ode15i` function to solve the implicit differential equation.

Code Equation

To code the equation in a form suitable for `ode15i`, you need to write a function with inputs for $\mathit{t}$, $\mathit{y}$, and ${\mathit{y}}^{\prime }$ that returns the residual value of the equation. The function `@weissinger` encodes this equation. View the function file.

`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; ```

Calculate Consistent Initial Conditions

The `ode15i` solver requires consistent initial conditions, that is, the initial conditions supplied to the solver must satisfy

$\mathit{f}\left({\mathit{t}}_{0},\mathit{y},{\mathit{y}}^{\prime }\right)=0$.

Since it is possible to supply inconsistent initial conditions, and `ode15i` does not check for consistency, it is recommended that you use the helper function `decic` to compute such conditions. `decic` holds some specified variables fixed and computes consistent initial values for the unfixed variables.

In this case, fix the initial value $\mathit{y}\left({\mathit{t}}_{0}\right)=\sqrt{\frac{3}{2}}$ and let `decic` compute a consistent initial value for the derivative ${\mathit{y}}^{\prime }\left({\mathit{t}}_{0}\right)$, starting from an initial guess of ${\mathit{y}}^{\prime }\left({\mathit{t}}_{0}\right)=0$.

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

Solve Equation

Use the consistent initial conditions returned by `decic` with `ode15i` to solve the ODE over the time interval $\left[1\text{\hspace{0.17em}}10\right]$.

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

Plot Results

The exact solution of this ODE is

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

Plot the numerical solution `y` computed by `ode15i` against the analytical solution `ytrue`.

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

## Input Arguments

collapse all

Functions to solve, specified as a function handle that defines the functions to be integrated. `odefun` represents the system of implicit differential equations that you want to solve using `ode15i`.

The function `f = odefun(t,y,yp)`, for a scalar `t` and column vectors `y` and `yp`, must return a column vector `f` of data type `single` or `double` that corresponds to $f\left(t,y,y\text{'}\right)$. `odefun` must accept all three input arguments, `t`, `y`, and `yp` even if one of the arguments is not used in the function.

For example, to solve $y\text{'}-y=0$, use this function.

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

For a system of equations, the output of `odefun` is a vector. Each equation becomes an element in the solution vector. For example, to solve

`$\begin{array}{l}y{\text{'}}_{1}-{y}_{2}=0\\ y{\text{'}}_{2}+1=0\text{\hspace{0.17em}},\end{array}$`

use this function.

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

For information on how to provide additional parameters to the function `odefun`, see Parameterizing Functions.

Example: `@myFcn`

Data Types: `function_handle`

Initial time, specified as a scalar. `decic` uses the initial time to compute consistent initial conditions that satisfy ```odefun(t0,y0_new,yp0_new) = 0```.

Data Types: `single` | `double`

Initial guesses for `y`-components, specified as a vector. Each element in `y0` specifies an initial condition for one dependent variable ${y}_{n}$ in the system of equations defined by `odefun`.

Data Types: `single` | `double`

`y`-components to hold fixed, specified as a vector of 1s and 0s, or as `[]`.

• Set `fixed_y0(i) = 1` if no change is permitted in the guess for `y0(i)`.

• Set `fixed_y0 = []` if any entry can be changed.

You cannot fix more than `length(yp0)` components. Depending on the specific problem, it is not always possible to fix certain components of `y0` or `yp0`. It is a best practice not to fix more components than is necessary.

Initial guesses for `y'`-components, specified as a vector. Each element in `yp0` specifies an initial condition for one differentiated dependent variable $y{\text{'}}_{n}$ in the system of equations defined by `odefun`.

Data Types: `single` | `double`

`y'`-components to hold fixed, specified as a vector of 1s and 0s, or as `[]`.

• Set `fixed_yp0(i) = 1` if no change is permitted in the guess for `yp0(i)`.

• Set `fixed_yp0 = []` if any entry can be changed.

You cannot fix more than `length(yp0)` components. Depending on the specific problem, it is not always possible to fix certain components of `y0` or `yp0`. It is a best practice not to fix more components than is necessary.

Options structure, specified as a structure array. Use the `odeset` function to create or modify the options structure. The relevant options for use with the `decic` function are `RelTol` and `AbsTol`, which control the error thresholds used to compute the initial conditions.

Example: `options = odeset('RelTol',1e-5)`

Data Types: `struct`

## Output Arguments

collapse all

Consistent initial conditions for `y0`, returned as a vector. If the value of `resnrm` is small, then `yo_new` and `yp0_new` satisfy ```odefun(t0,y0_new,yp0_new) = 0``` and are suitable to be used as initial conditions with `ode15i`.

Consistent initial conditions for `yp0`, returned as a vector. If the value of `resnrm` is small, then `yo_new` and `yp0_new` satisfy ```odefun(t0,y0_new,yp0_new) = 0``` and are suitable to be used as initial conditions with `ode15i`.

Norm of residual, returned as a vector. `resnrm` is the norm of `odefun(t0,y0_new,yp0_new)`.

• A small value of `resnrm` indicates that `decic` successfully computed consistent initial conditions that satisfy `odefun(t0,y0_new,yp0_new) = 0`.

• If the value of `resnrm` is large, try adjusting the error thresholds `RelTol` and `AbsTol` using the `options` input.

## Tips

• The `ihb1dae` and `iburgersode` example files use `decic` to compute consistent initial conditions before solving with `ode15i`. Type ```edit ihb1dae``` or `edit iburgersode` to view the code.

• You can additionally use `decic` to compute consistent initial conditions for DAEs solved by `ode15s` or `ode23t`. To do this, follow these steps.

1. Rewrite the system of equations in fully implicit form `f(t,y,y') = 0`.

2. Call `decic` to compute consistent initial conditions for the equations.

3. Specify `y0_new` as the initial condition in the call to the solver, and specify `yp_new` as the value of the `InitialSlope` option of `odeset`.