Among other things, Greg M. Johnson wrote:

> Could someone pseudo-code it out for me.
> 
> First, I'll share my non-energy conservation system that I used in my
> flocking algorithm ( http://www.geocities.com/pterandon/boids.html ), and
> in planetary orbiting simulations.
> 
> dt=1  (we're stepping along povray frame by povray frame)
> p= vector for current position
> v=velocity
> a=acceleration, or sum of forces at current position (based on gravity of
> sun, "repulsion" of neighbors in flock),
> 
> ----- for each frame----
> p2=p1+v1 *dt
> a2= constant * (sum of new forces based on new location).
> v2=v1+a2
> --- repeat ad nauseum ---
> 

Hmm... I thought I had replied to this... Was my post lost somehow? I can't 
find it anywhere :/

Well, I'll write it again.

First, shouldn't it be:

  v2=v1+a2*dt ?

Second, as I mentioned, there are several algorithms used for integrating 
the equations of motion (that's what you're trying). An algorithm that 
minimizes numerical errors is the so-called "velocity Verlet" [Swope et al. 
J. Chem. Phys. 76 (1982) 637]:

-------------
  Calculate all forces on all particles, for each particle:
    a1 = (total force at p1)/(real or fictitious mass)
1 Calculate new positions for all particles:
    p2 = p1 + v1*dt + 1/2*a1*(dt^2)
  Calculate intermediate velocities (at t1+1/2*dt) for all particles:
    v' = v1 + 1/2*a1*dt
  Calculate all forces on all particles at their new position:
    a2 = (total force at p2)/(real or fictitious mass)
  Caluclate the new velocities for all particles:
    v2 = v' + 1/2*a2*dt
  Goto 1 (p2 -> p1, v2 -> v1, a2 -> a1)
-------------

This should conserve energy. What is not conserved is kinetic energy 
(temperature). If a simulation at constant temperature is desired, there 
are some "tricks" available, the simplest one is just re-scaling all 
velocities so that kinetic energy remains constant. After calculating the 
new velocities, do:

  f = sqrt(2*Ek/sum(mass*v2^2))
  v2(scaled) = f*v2

where Ek is the desired constant kinetic energy.

-- 
light_source{9+9*x,1}camera{orthographic look_at(1-y)/4angle 30location
9/4-z*4}light_source{-9*z,1}union{box{.9-z.1+x clipped_by{plane{2+y-4*x
0}}}box{z-y-.1.1+z}box{-.1.1+x}box{.1z-.1}pigment{rgb<.8.2,1>}}//Jellby

Post a reply to this message