Question about the corrector step of symplectic scheme in the code

edited August 2014 in Old versions
Dear All

I am not quite clear about the corrector step of Symplectic step in the code JSphCpu.cpp.
The function shows ComputeSymplecticCorr(bool rhopbound, float dt) shows the the way to calculate the velocity and position in the second step which is different as described in the user guide v2.000 or the publish paper R.Vacondio(2013).

In the code:
Line 1394: Vel[p].x=VelPre[p].x + Ace[p].x * dt;
According to the paper, the velocity in the second step should be
Vel[p].x=Vel[p].x+Ace[p].x*dt05
(here, Vel[p].x and Ace[p].x is the velocity and acceleration calculated by the first step)

And, the position of the first half step(n+1/2) and the velocity in the step n+1 should be used to update the position in step n+1, however the code does not seem to calculate like that.
Line 1397: float dx=(VelPre[p].x+Vel[p].x) * dt05;
Line 1401: Pos[p].x=PosPre[p].x + dx;

I am not sure if I misunderstand the code or not, please help me.
Thank you.

Wang

Comments

  • In terms of the density in corrector step, it is updated by the function ComputeRhopEpsilon(), which is also different from the references. would you please tell me the paper referring to this method?

    Thank you.
  • Hi,

    In line-1394 that you refer above, the code is correct.
    I recommend you find a book with the predictor/corrector scheme (not paper, book an old trustworthy book).

    the scheme works as follows:

    Predictor,
    V(n+1/2)=V(n)+F(n)Dt/2
    (I have not checked that in the code yet, but I pretty sure it is correct, however you never know I might give it a look)

    Corrector:
    V(n+1/2)=V(n)+F(n+1/2)Dt/2

    Finally:
    V(n+1)=2V(n+1/2)-V(n)

    I hope we all agree to that :),
    Now substitute the corrector step to the 3rd equation (that has obviously been done to avoid extra arithmetic and memory) and you get

    V(n+1)=V(n)+DtF(n+1/2) as per line-1394

    I will check the other half of you question (lines 1397 and 1401) and get back as soon as I can.

    Regards,
    George

Sign In or Register to comment.