Re: Re: [Freepooma-devel] ODE System of arbitrary length

[Frank Gollas]

Richard Guenther wrote:
>Harris wrote:
>>Hi,
>>
>>I just stumbled on POOMA and I think it's pretty 
>>cool and better suited for my needs than blitz++.
>>
>>What I want to do is a numerical 4 step Runge-Kutta 
>>integration of an ODE system 
>>
>>dx1/dt = f1(x1,x2, x3 ..., Y )
>>dx2/dt = f2(x1,x2, x3 ..., Y )
>>dx3/dt = f3(x1,x2, x3 ..., Y )
>>...
>>
>>this should be rather straight forward.
>>But: the number of state equations is only 
>>known at _runtime_, the class that does 
>>the integration holds a dynamical list 
>>of arrays. Is there an efficient
>>way of doing this using stencils?
>
>I don't exactly understand what you mean by "dynamical list
>of arrays", but I can show you the code-snipped I'm using to
>integrate a gravitating n-Body system in time using a 4th order
>Runge-Kutta (attached).
>
>Basically, you need to use DynamicArray to store the state (and
>temporaries) for a variable number of bodies.  The integration
>part is data-parallel, while the force calculation (which is
>brute force n^2 algorithm) is using manual loops.
>
>If I understand you correctly, you want to use a stencil for
>calculating the integration in time.  This is possible, if
>your "force" can be calculated locally.  Maybe the code
>snipped helps answering more questions, otherwise just ask.
>
>Richard.


thank you very much for your fast answer and the code snipsets that were
pretty helpful. I'm also attaching some of my testing code to make clear my
intention. x is a container of Array's, each representing a state variable
field of the ODE system. f1 to fn are represented by stencils. Unfortunatly
I didn't manage to express f1(x1,x2,x3 ... ), f2(x1,x2,x3 ... ) in POOMA. In
my code only one equation f1(x[i]) with one array as parameter is applied.

The goal is to have as many Stencils as state equations in the ODE System,
each taking all arrays representing the state variables (and some constant
input ..) as argument, such that interaction between different state
variables can be implemented.
 
Is it possible to do this with layers? Does my code look ok up to this
point? BTW the Laplace Stencil example did compute nonsense until I replaced
"world(interior_2) = stencil(world);" with "temp(interior_2) =
stencil(world); world = temp;" what might be the reason for that?

Frank





-- 
GMX DSL = Maximale Leistung zum minimalen Preis!
2000 MB nur 2,99, Flatrate ab 4,99 Euro/Monat: http://www.gmx.net/de/go/dsl

pooma_test.cpp
Description: Text Data

MyStencils.h
Description: Text Data