\ ob.fpgravity
\
\ hp_fpgravity
\ gravity class that uses floating point calculations.
\ part of the collection of software componenets: "henri poincare."
\
\ based on ideas presented by Joel Ryan.
\
\ code: Han-earl Park
\ copyright 2004 buster & friends' C-ALTO Labs
\ (Den haag, December 1997 -
\ (Valencia, October 1998 -
\ (Southampton, February 2002 -
\
\ MOD: HeP 03-19-04 First working floating point version.
\ Fixed the kludgy stack behavior in ABS<.
\ Remove redundant calcualtions when parallel updating of
\ particles.
\ MOD: HeP 03-20-04 Cleanup code: FPGRAV.CALC.DISP calls GRAV.CALC.DIM.DISP in
\ its loop.
\ Split off floating point version into new class.
\
\ ToDo: store vectored displacement between particles internally
\ ToDo: NO internal scaling...?
include? task-hp_fputil myt:hp_fputil
include? task-hp_force myt:hp_gravity
anew task-hp_fpgravity
:class OB.FPGRAVITY abs
LOOP
2drop
fvec->abs \ x y z... -f- |r|
;
: FPGRAV.CALC.FORCE ( disp m1 m2 -f- force , gravitational pull )
f* \ -f- disp m1*m2
iv-force-constant float f* \ -f- disp m1*m2*k
fswap \ -f- m1*m2*k disp
\
#dimensions@ 2
DO
fdup f*
LOOP
f/ \ -f- force
;
fvariable fpgrav-force \ temp stoarge to avoid stack manipulations
fvariable fpgrav-disp
: FPGRAV.UPDATE.FORCE { p1 p2 -- , disp force -f- }
fpgrav-force f!
fpgrav-disp f!
\
#dimensions@ 0
DO
p1 p2 i GRAV.CALC.DIM.DISP float \ -f- dim_disp
\
fpgrav-force f@ f* \ -f- dim_disp*force
fpgrav-disp f@ f/ \ -f- dim_force
\
FIX
hp_force i ed.at: p1 + \ sum value to current
hp_force i ed.to: p1 \ update
LOOP
;
:m (ELM.CALC.FORCE): ( -- )
iv-force-#pucks 0
DO
i iv-force-cur-puck# = not \ ignore if current particle
IF
iv-force-cur-puck \ -- p1
i iv-force-cur-space get: [] \ -- p1 p2
\
2dup FPGRAV.CALC.DISP \ -- p1 p2 , -f- disp
2dup fdup \ -- p1 p2 p1 p2 , -f- disp disp
\
get.mass: [] i>f \ -- p1 p2 p1 , -f- disp disp m2
get.mass: [] i>f \ -- p1 p2 , -f- disp disp m2 m1
FPGRAV.CALC.FORCE \ -- p1 p2 , -f- disp force
\
FPGRAV.UPDATE.FORCE \ -- , -f-
\
\ FIX FIX swap
\ GRAV.UPDATE.FORCE
THEN
LOOP
;m
\ update particles in parallel
: PARALLEL.FPGRAV.UPDATE.FORCE { p1 p2 -- , disp force -f- }
fpgrav-force f!
fpgrav-disp f!
\
#dimensions@ 0
DO
p1 p2 i GRAV.CALC.DIM.DISP float \ -f- dim_disp
\
fpgrav-force f@ f* \ -f- dim_disp*force
fpgrav-disp f@ f/ \ -f- dim_force
\
FIX
dup
hp_force i ed.at: p1 + \ sum value to current
hp_force i ed.to: p1 \ update
\
negate
hp_force i ed.at: p2 + \ sum value to current
hp_force i ed.to: p2 \ update
LOOP
;
: PARALLEL.FPGRAV.CALC.FORCE ( -- )
iv-force-#pucks iv-force-cur-puck# 1+ \ update only those needed
DO
iv-force-cur-puck \ -- p1
i iv-force-cur-space get: [] \ -- p1 p2
\
2dup FPGRAV.CALC.DISP \ -- p1 p2 , -f- disp
2dup fdup \ -- p1 p2 p1 p2 , -f- disp disp
\
get.mass: [] i>f \ -- p1 p2 p1 , -f- disp disp m2
get.mass: [] i>f \ -- p1 p2 , -f- disp disp m2 m1
FPGRAV.CALC.FORCE \ -- p1 p2 , -f- disp force
\
PARALLEL.FPGRAV.UPDATE.FORCE \ -- , -f-
LOOP
;
:m (CALC.FORCE): ( -- , calculate forces for all particles )
iv-force-#pucks 0
DO
i iv-force-cur-space get: [] \ -- p
\
dup get.type: [] self get.type: [] \ -- p t1 t2
particle.type=
IF
iv=> iv-force-cur-puck
i iv=> iv-force-cur-puck#
\
PARALLEL.FPGRAV.CALC.FORCE \ --
ELSE
drop
THEN
LOOP
;m
;class