home *** CD-ROM | disk | FTP | other *** search
/ Aminet 10 / aminetcdnumber101996.iso / Aminet / misc / sci / StarCollapse.lha / pull.s < prev    next >
Text File  |  1995-09-11  |  5KB  |  234 lines
  1. *****************************************************************************
  2. *
  3. * CalcPull                                            (C) L. Vanhelsuwé
  4. * --------
  5. *
  6. * Assembler module for GR gravity simulation.
  7. * Calculate Gravitational Pull on mass Cloud[m].
  8. *
  9. * HISTORY
  10. * -------
  11. * 19-AUG-95: started this file out of frustration from Lattice's inability
  12. *            to generate correct 68881 code itself.
  13. *
  14. *
  15. *****************************************************************************
  16.  
  17.  
  18.         XDEF    _asm_calc_pull                ;what we export
  19.  
  20.         XREF    _CloudSize,_Cloud            ;what we import
  21.  
  22.     ; Mass STRUCTURE field offsets                ;**!! track C's idea of structure!
  23.  
  24. mass        equ        0*8
  25. x            equ        1*8
  26. y            equ        2*8
  27. z            equ        3*8
  28. Vx            equ        4*8
  29. Vy            equ        5*8
  30. Vz            equ        6*8
  31. mass_sizeof    equ        7*8
  32.  
  33.  
  34. saved        equ        7*4                        ;number of regs saved
  35.  
  36. star        equ        saved+4
  37. F_ptr        equ        saved+8
  38. Fx_ptr        equ        saved+12
  39. Fy_ptr        equ        saved+16
  40. Fz_ptr        equ        saved+20
  41.  
  42. *--------------------------------------------------------------------------
  43. * void    calc_pull (int star, num *F_ptr, num *Fx_ptr, num *Fy_ptr, num *Fz_ptr);
  44. *
  45. *--------------------------------------------------------------------------
  46.  
  47.     ; Receives C arguments on stack as follows:
  48.     ;  4(SP) = mass index (0.. CloudSize -1)
  49.     ;  8(SP) = F_ptr
  50.     ; 12(SP) = Fx_ptr
  51.     ; 16(SP) = Fy_ptr
  52.     ; 20(SP) = Fz_ptr
  53.  
  54.  
  55. _asm_calc_pull:    movem.l    d0/d7/a0-a4,-(SP)    ;we're called from C, so save all !
  56.  
  57.                 move.l    _CloudSize,d7        ;scan entire Cloud
  58.                 subq.l    #1,d7
  59.  
  60.                 move.l    _Cloud,a0            ;-> array of Mass structures (the Cloud)
  61.  
  62.                 move.l    star(SP),d0
  63.                 mulu    #mass_sizeof,d0
  64.                 lea        0(a0,d0.L),a1        ;-> mass under consideration
  65.  
  66.                 move.l    Fx_ptr(SP),a2        ; = Fx ptr
  67.                 move.l    Fy_ptr(SP),a3        ; = Fy ptr
  68.                 move.l    Fz_ptr(SP),a4        ; = Fz ptr
  69.  
  70.                 fmovecr    #$0F,fp0            ; constant 0.0
  71.                 fmove.d    fp0,(a2)            ; Fx = 0.0    (in double format)
  72.                 fmove.d    fp0,(a3)            ; Fy = 0.0    (in double format)
  73.                 fmove.d    fp0,(a4)            ; Fz = 0.0    (in double format)
  74.  
  75.                 fmove.d    x(a1),fp0            ; Cloud[star].x        ! CACHED !
  76.                 fmove.d    y(a1),fp1            ; Cloud[star].y
  77.                 fmove.d    z(a1),fp2            ; Cloud[star].z
  78.  
  79. sum_vectors        cmp.l    a0,a1
  80.                 beq        skip_self
  81.  
  82.                 fmove.d    x(a0),fp3
  83.                 fsub    fp0,fp3                ; dx = Cloud[m].x - Cloud[star].x;
  84.  
  85.                 fmove.d    y(a0),fp4
  86.                 fsub    fp1,fp4                ; dy = Cloud[m].y - Cloud[star].y;
  87.  
  88.                 fmove.d    z(a0),fp5
  89.                 fsub    fp2,fp5                ; dz = Cloud[m].z - Cloud[star].z;
  90.  
  91.                 fmove    fp3,fp7
  92.                 fmul    fp7,fp7                ; dx^2
  93.  
  94.                 fmove    fp4,fp6
  95.                 fmul    fp6,fp6                ; dy^2
  96.                 fadd    fp6,fp7                ; dx^2 + dy^2
  97.  
  98.                 fmove    fp5,fp6
  99.                 fmul    fp6,fp6                ; dz^2
  100.                 fadd    fp6,fp7                ; dx^2 + dy^2 + dz^2     (= R^2)
  101.  
  102.                 ftst    fp7                    ; if (R^2 == 0) continue
  103.                 fbeq    skip_self
  104.  
  105.                 fmove.d    mass(a0),fp6        ; mass 1
  106.                 fmul.d    mass(a1),fp6        ; mass 1 * mass 2
  107.                 fmul.l    #100,fp6            ; GRAV * m1 * m2
  108.                 fdiv    fp7,fp6                ; F = G*m1*m2/R^2
  109.  
  110.                 fsqrt    fp7                    ; R = SQRT( R^2 )
  111.  
  112. ;                ftst    fp7                    ; if (R == 0) continue
  113. ;                fbeq    skip_self
  114.  
  115.                 fmul    fp6,fp3                ; F*dx
  116.                 fdiv    fp7,fp3                ; Fx = F*dx/R
  117.  
  118.                 fmul    fp6,fp4                ; F*dy
  119.                 fdiv    fp7,fp4                ; Fy = F*dy/R
  120.  
  121.                 fmul    fp6,fp5                ; F*dz
  122.                 fdiv    fp7,fp5                ; Fz = F*dz/R
  123.  
  124.                 fmove.d    (a2),fp6            ; accumulate into Fx,Fy,Fz
  125.                 fadd    fp3,fp6                ; Fx += Fxcouple
  126.                 fmove.d    fp6,(a2)
  127.  
  128.                 fmove.d    (a3),fp6
  129.                 fadd    fp4,fp6                ; Fy += Fycouple
  130.                 fmove.d    fp6,(a3)
  131.  
  132.                 fmove.d    (a4),fp6
  133.                 fadd    fp5,fp6                ; Fz += Fzcouple
  134.                 fmove.d    fp6,(a4)
  135.  
  136. skip_self        add.w    #mass_sizeof,a0        ; next Cloud element
  137.                 dbra    d7,sum_vectors
  138.  
  139.                 fmove.d    (a2),fp0            ;Fx
  140.                 fmul    fp0,fp0                ;Fx^2
  141.  
  142.                 fmove.d    (a3),fp1            ;Fy
  143.                 fmul    fp1,fp1                ;Fy^2
  144.  
  145.                 fmove.d    (a4),fp2            ;Fz
  146.                 fmul    fp2,fp2                ;Fz^2
  147.  
  148.                 fadd    fp2,fp1
  149.                 fadd    fp1,fp0                ; Fx^2 + Fy^2 + Fz^2
  150.                 fsqrt    fp0                    ;SQRT ( Fx^2 + Fy^2 + Fz^2 )
  151.  
  152.                 move.l    F_ptr(SP),a0
  153.                 fmove.d    fp0,(a0)            ; return F
  154.  
  155.                 movem.l    (SP)+,d0/d7/a0-a4
  156.  
  157.                 rts
  158.  
  159.         END
  160.  
  161. *--------------------------------------------------------------------------
  162. * Below is the original C source which this assembler routine replaces.
  163. *--------------------------------------------------------------------------
  164.  
  165.  
  166.     IFD     THIS_IS_C_SOURCE
  167.  
  168. typedef struct Mass {
  169.     num    mass;
  170.     num    x,y,z;
  171.     num    Vx,Vy,Vz;
  172.  
  173. } Mass;
  174.  
  175. //--------------------------------------------------------------------------
  176. //--------------------------------------------------------------------------
  177.  
  178. void    calc_pull (int star, num *F_ptr, num *Fx_ptr, num *Fy_ptr, num *Fz_ptr) {
  179.  
  180. int m;
  181. num    F,Fx,Fy,Fz;
  182. num Fcouple;                    // force magnitude between 2 masses
  183. num    Fxcouple,Fycouple,Fzcouple;    // force vector
  184. num    dx,dy,dz;
  185. num    R;
  186.  
  187.         // reset force accumulators
  188.  
  189.     Fx = Fy = Fz = 0;
  190.  
  191.         // forall bodies (except body under consideration)
  192.  
  193.     for (m=0; m < CloudSize; m++) {
  194.  
  195.         if (m == star) continue;            // skip oneself...
  196.  
  197.             // calc distance between bodies
  198.  
  199.         dx = Cloud[m].x - Cloud[star].x;
  200.         dy = Cloud[m].y - Cloud[star].y;
  201.         dz = Cloud[m].z - Cloud[star].z;
  202.     
  203.         R = dx*dx + dy*dy + dz*dz;
  204.  
  205.             // If body collides with another, no force results
  206.     
  207.         if (R == 0) continue;
  208.     
  209.         Fcouple = GRAVITY_CONST * Cloud[star].mass * Cloud[m].mass / R;
  210.  
  211.         R = sqrt(R);
  212.     
  213.         if (R == 0) continue;
  214.     
  215.         Fxcouple = Fcouple * dx/R;
  216.         Fycouple = Fcouple * dy/R;
  217.         Fzcouple = Fcouple * dz/R;
  218.  
  219.         Fx += Fxcouple;
  220.         Fy += Fycouple;
  221.         Fz += Fzcouple;
  222.     }
  223.  
  224.     F = sqrt( Fx*Fx + Fy*Fy + Fz*Fz );
  225.  
  226.     *F_ptr    = F;
  227.     *Fx_ptr    = Fx;
  228.     *Fy_ptr    = Fy;
  229.     *Fz_ptr    = Fz;
  230. }
  231. //--------------------------------------------------------------------------
  232.     ENDIF
  233.  
  234.