Choose the best method of integration of equations of movement about sibernetic HOT 9 CLOSED

Fortunately, I've also tried another one symplectic method (1st order) - semi-implicit Euler method.
Surprisingly, it does not suffer from finite precision problem - the solution is very close to analytical one at timeStep values up to 1e-04 s (previously we were able to use timeStep up to 5e-6 s).

In this simple test it really works at this timeStep, but in a simulation including liquid I see stable behavior of the system only at timeStep values up to 1e-05 s, but even this is twice better than 5e-06. Our recent worm body model and liquid around it works fine with this 1st order semi-implicit Euler method. So, I consider it as a good solution which singificanly improved adequacy of Sibernetic and I plan to update the master branch of the Sibernetic code soon with this new method.
Also there is one more variant - we should keep it in mind for the case if even better precision will be necessary. It is 2nd order symplectic LeapFrog algorithm, which works only if we use double values for accelerations, velocities and positions, but right now, in my opinion, it is not necessary.

from sibernetic.

What is the computation cost per step of semi-implicit Euler vs the existing method? The total computation time to advance the simulation by time delta-t is the important metric.

Double precision calculations do take longer than 32 bit, but if we are able to use a simpler numerical method, or to significantly increase the time step, then the trade-off may be worth it.

from sibernetic.

@a-palyanov believes we have the best option for floats but we may need a different option for doubles

from sibernetic.

@a-palyanov thanks for the update on this today in the hangout -- feel free to write a bit more to illustrate where this is so other folks can get the update too.

from sibernetic.

Float_max = 3.402823e+038
Float_min = 1.175494e-038
Double_max = 1.797693e+308
Double_min = 2.225074e-308

Double vs float:

Float:

As we can see at this example, -118.625 is composed of [-1] * [118625] * [1e-03].
So, for fraction part we can work with values from 0 to some MAX_VALUE.

from sibernetic.

And here is an illustration of how this influences on precision which is available for us in case of calculations (numerical integration).

from sibernetic.

My calculations showed that for 32bit float coordinates (in all three dimensions) inside the interval ~ (-100000_r0...100000_r0) (r0 = 'radius' of the SPH particle) the maximal integral error connected with this effect will be limited by ~0.002_r0 per 100 integration steps (or ~0.00002_r0 per 1 integration step), which will not be noticeable. This numbers are for integration time step = 3.00e-05 s.
We currently use in simulation 5e-06 s. I hope I will soon recalculate for this value too, but anyway results shouldn't be worse.

from sibernetic.

Remember to ensure that every stage of calculation must stay within the
bounds 1E38 1E-38. There are some places in the original kernel that
violate this restriction.

Also when summing you need to ensure that the operands are of a similar
magnitude.
On Jun 19, 2014 12:54 AM, "a-palyanov" [email protected] wrote:

My calculations showed that for 32bit float coordinates (in all three
dimensions) inside the interval ~ -100000_r0...100000_r0
http://r0%20=%20'radius'%20of%20the%20SPH%20particle the maximal
integral error connected with this effect will be limited by ~0.002_r0 per
100 integration steps (or ~0.00002_r0 per 1 integration step), which will
not be noticeable.

—
Reply to this email directly or view it on GitHub
#19 (comment)
.

from sibernetic.

More precisely, we must stay within the bounds (3.40e+38...1.17e-38, 0.00e+00 , -1.17e-38...-3.40e+38)
And I hope in this new code (https://github.com/openworm/Smoothed-Particle-Hydrodynamics/tree/development) which I committed today this restrictions are satisfied.

from sibernetic.