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
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
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 < jerzy.karczmarczuk@unicaen.fr> 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@haskell.org http://www.haskell.org/**mailman/listinfo/haskell-cafehttp://www.haskell.org/mailman/listinfo/haskell-cafe
participants (3)
-
felipe zapata -
Jerzy Karczmarczuk -
L Corbijn