# Wronski matrix of a discrete evolution and condition

To compute the condition

$$\kappa_{\mathrm{abs}}^\infty := \max\{\| W(t,t_0, y_0) \|: t_0 \leq t \leq T \}$$

of an initial value problem on $[t_0, T]$, I need to compute the Wronski matrix $W$.

The approximation of the solution is denoted by the discrete evolution operator $\Psi^{t_0, t}$ and the exact solution is denoted $\Phi^{t_0, t}$ and the Wronski matrix is given as $W(t,t_0, y) := \frac{d}{dy} \Phi^{t_0, t}y \Big |_{y=z}$

My first question: is the condition really a function of the exact solution? According to my understanding numerical analysis deals with finding good approximation functions so the condition should be a measure of how good (stable?) a given approximation is, therefore the Wronski matrix should be a function of $\Psi$, not of $\Phi$.

My second question: how do I explicitly calculate the Wronski matrix for a given problem? Say, I'm using the explicit Euler method $y_{k+1} = y(t_k) + h*f(t_k, y(t_k))$ to solve the ODE $\dot{y} = (\alpha - \beta y)y$ or the system $\dot{u} = (\alpha - \beta v)u, \dot{v} = (\delta u - \gamma)v$.

-

First question: As you say, the condition number is a function of the exact solution. It measures how much the (exact) solution of the differential equation changes if you change the initial condition. The condition number depends only on the problem you are solving, not on the numerical method you are using.

So, why are we interested in this condition number in numerical analysis? Because the condition number gives a lower bound on the error committed by any numerical method. The idea is that the initial condition is inaccurate, either because it comes from some measurements which is not perfect or because of rounding error. This inaccuracy in the initial condition leads to an error in the solution as measured by the condition number, regardless of how good the numerical method is. On top of that comes the errors committed by the numerical method.

Second question: One idea for computing the condition number is to use that the numerical evolution operator $\Psi$ approximates the exact evolution operator $\Phi$. Thus, to compute the condition number, replace $\Phi$ in your definition by $\Psi$. You can compute the derivative $\frac{\partial}{\partial y} \Psi^{t_0,t}(y)$ by finite differences or by evaluating the derivative by hand.

-
Many thanks. I like your answer to my first question, it's very clear. As for my second question: I'm not entirely satisfied. For example: $\Psi^h y = y + hf(y)$ then trying to solve $\dot{y} = f(y) = (\alpha - \beta y)y$ yields $\frac{\partial}{\partial y} \Psi^h y = 1 + h \alpha - 2h\beta y$ which is a scalar so the maximal determinant is just the maximum value of the function $1 + h\alpha - 2h \beta y(t)$, for $t$ in the domain. Is this correct? –  Rudy the Reindeer Jan 26 '11 at 19:41
Second example: using $\Psi^h y = y + h f(y)$ again, this time to solve $\dot{u} = (\alpha - \beta v)u$, $\dot{v} = (\delta u - \gamma)v$, I get $$\Psi^h y = \left( \begin{array}{c} u + h(\alpha - \beta v) u \\\ v + h(\delta u - \gamma)v \end{array} \right)$$ but I don't know what to do now. Can you tell me how to proceed from there? –  Rudy the Reindeer Jan 26 '11 at 19:41