17 Apr 2013 00:36

## Yampa integral function

Hello,
I'm trying to write some simulator of
physical systems using Yampa.
The fundamental idea is to integrate
the position and velocity, using
the following algorithm:

1. get the velocity using the formula,

v = v0 + integral ( 0.5 * force * recip mass)

where v0 is the previous velocity.

2. calculate the new position as follows,

x = x0 + integral (v)

and some other steps that are not relevant here.
Therefore I wrote the following piece of
code:

type Pos   = Double
type Vel   = Double
type Mass  = Double
type Force = Double

integralCoordSF :: Force -> Mass ->  SF (Pos,Vel) (Pos,Vel)
integralCoordSF  force mass = proc (x0,v0) -> do
v <- integralVelSF force mass -< v0
vdt2 <- integral -< v
let x = x0 + vdt2
returnA -< (x,v)

integralVelSF :: Force -> Mass ->  SF Double Double
integralVelSF force mass = proc v0 -> do
t1 <- integral -< 0.5 * (force/mass)
returnA -< v0 + t1

Now, let's the speed (velocity), force and the mass equal 1.0.  And
the initial potion equal zero. Then using the "embed" function we get that,

*Main>  embed (integralCoordSF  1.0 1.0) ((0.0,1.0),[(0.1,Nothing),(0.1,Nothing)])
[(0.0,1.0),(0.1,1.05),(0.20500000000000002,1.1)]

According to the algorithm, if  dt = 0.1 then

v = 1.0 + integral (0.5 * 1.0 * recip 1.0) = 1.05

in agreement with the algorithm.  But let's see about the position

x = 0 + integral (v)

Using the formula of step 2, the position should be equal to 0.105,
but the result given by embed is 0.1.

As a beginner certainly and missing a lot of concepts and details,
can someone please explain me what's going on here?

Thank you for your help!

Felipe Z

```_______________________________________________
Haskell-Cafe mailing list
Haskell-Cafe <at> haskell.org
http://www.haskell.org/mailman/listinfo/haskell-cafe
```
17 Apr 2013 01:48

### Re: Yampa integral function

Le 17/04/2013 00:36, felipe zapata a écrit :
Hello,
I'm trying to write some simulator of
physical systems using Yampa.
The fundamental idea is to integrate
the position and velocity, using
the following algorithm:
... ...
Nothing to do with Haskell.
Where do you get this factor 0.5 in your formula:  v = v0 + integral ( 0.5 * force * recip mass) ?
With constant acceleration v=v0+a*Dt =>  1.01, not 1.05

What is that:  let x = x0 + vdt2 ??
Again, a constant acceleration movement gives : x = x0 + v0*Dt + a*Dt^2/2.

Yampa or not Yampa, I suspect that your implementation of physics is simply wrong.
==

In general, even correcting all, you might have reasonable results in some trivial cases, but in general the extrapolating Euler schema is unstable, produces growing errors (e. g. in the oscillating case).

Jerzy Karczmarczuk

```_______________________________________________
Haskell-Cafe mailing list
Haskell-Cafe <at> haskell.org
http://www.haskell.org/mailman/listinfo/haskell-cafe
```
17 Apr 2013 01:55

### Re: Yampa integral function

```Le 17/04/2013 01:48, Jerzy Karczmarczuk a écrit :
> With constant acceleration v=v0+a*Dt =>  1.01, not 1.05
Gosh, trivial errors seem to be contagious. Of course I meant 1.1, not
1.01. Sorry. JK
```
17 Apr 2013 10:53

### Re: Yampa integral function

Warning: no knowledge of Yampa.

It seems that the integration for the position is always using the velocity from the previous step. Looking at the documentation in source code of Yampa this seems plausible as there is also a function called imIntegrate with the comment

-- "immediate" integration (using the function's value at the current time)

Which I assume would yield your answer. As I don't know Yampa I can only speculate about the reason for this difference, the two things I can think of now are better speed and to prevent trouble (in cases like x = integral y, y = integral x).

Lars

On Wed, Apr 17, 2013 at 1:55 AM, Jerzy Karczmarczuk wrote:
Le 17/04/2013 01:48, Jerzy Karczmarczuk a écrit :

With constant acceleration v=v0+a*Dt =>  1.01, not 1.05
Gosh, trivial errors seem to be contagious. Of course I meant 1.1, not 1.01. Sorry. JK

_______________________________________________
Haskell-Cafe mailing list
Haskell-Cafe <at> haskell.org
http://www.haskell.org/mailman/listinfo/haskell-cafe

```_______________________________________________
Haskell-Cafe mailing list
Haskell-Cafe <at> haskell.org
http://www.haskell.org/mailman/listinfo/haskell-cafe
```

Gmane