help-octave
[Top][All Lists]

## Re: solving odes

 From: Thomas Shores Subject: Re: solving odes Date: Thu, 16 Jun 2011 14:01:43 -0500 User-agent: Mozilla/5.0 (X11; U; Linux x86_64; en-US; rv:1.9.2.17) Gecko/20110428 Fedora/3.1.10-1.fc14 Thunderbird/3.1.10

```On 06/15/2011 06:44 AM, Piotr wrote:
```
```Hello. I want to solve numerically the system of ODEs:

\$y_0(t)=some\_given\_function,\$

\$y_n'(t)=y_n(t)a(t)+y_{n-1}(t)b(t)+y_n(t)x(t),  y_n(0)=y_n, n>0,\$

```
where \$a, b\$ are given functions and \$x\$ is the solution of the auxiliary problem:
```
\$x'(t)=f(t,x(t)), x(0)=x_0,\$

and \$n>0\$ is some natural number, which can be sometimes huge.

```
I approximated solution of the auxiliary problem by lsode. Then I interpolated it by a piecewise linear function and solve the main equation for \$n=1\$ by lsode. For \$n=2\$ I tried to apply the same procedure with interpolation of \$x\$ and \$y_1,\$ but it failed. Could somebody suggest me a better (simpler) approach to this problem? I'm a newbie in Octave, so my idea is not sophisticated. I would be grateful for any help. Piotr.
```
_______________________________________________
Help-octave mailing list
https://mailman.cae.wisc.edu/listinfo/help-octave
```
Looks like you're trying to do some sort of method of lines. There are any number of reasons for which lsode might fail, one of the leading being that your system might not have a solution over the interval under consideration (consider the innocuous looking y'(t) = 1/(1+y^2), y(0) = 0, over the interval [0,2]). Without further information it's hard to tell what the problem is. One puzzling thing, though, is why you go to the trouble to uncouple this system. You could express it as a single vector IVP, of the form dY/dt = F(t,Y(t)), Y(0)=Y0, where, in your notation (treating x(0), y_n(0), and y_0(t) as known),
```
```
Y(t) = [x(t), y_1(t), ... , y_N(t)]', (here ' means transpose, not derivative)
```Y0 = [x(0), y_1(0), ... , y_N(0)]', and
F(t,Y(t)) = [f(t,x(t)), (a(t)+x(t))*y_1(t)+b(t)*y_0(t), ... ,
(a(t)+x(t))y_N(t)+b(t)*y_{N-1}(t)]'

Hope this helps,

Thomas Shores

```