Simulation for the n-body problem

A small java application is shown below. The equivalent MoPL-3 code is appended.

#define G 9.81
#define x0 SW1[1,1]
#define y0 SW1[2,1]
#define z0 SW1[3,1]
#define vx0 SW1[1,2]
#define vy0 SW1[2,2]
#define vz0 SW1[2,3]
#define m0 SW1[1,3]
#define pos0 SW1[2,3]
#define ax0 SW1[1,4]
#define ay0 SW1[2,4]
#define az0 SW1[3,4]
#define aax0 SW1[1,5]
#define aay0 SW1[2,5]
#define aaz0 SW1[3,5] #define x1 SW2[1,1]
#define y1 SW2[2,1]
#define z1 SW2[3,1]
#define vx1 SW2[1,2]
#define vy1 SW2[2,2]
#define vz1 SW2[2,3]
#define m1 SW2[1,3]
#define pos1 SW2[2,3]
#define ax1 SW2[1,4]
#define ay1 SW2[2,4]
#define az1 SW2[3,4]
#define aax1 SW2[1,5]
#define aay1 SW2[2,5]
#define aaz1 SW2[3,5]
#define r sqrt ((x0-x1)*(x0-x1) + (y0-y1)*(y0-y1) + (z0-z1)*(z0-z1))
#define r3 r * r * r
#define r5 r * r * r * r * r
array		Data		[1:3, 1:500] of float;	/* data 100x(3x5) */
scanPattern	InnerScan		is 100 steps [5,0];
OuterScan is 100 steps [5,0];
window nBodyW is
SW1, SW2 [1:5, 1:3] handle [1,1] of float;

rALUsubnet nBody of nBodyW is begin
if (pos0 != pos1) then begin
ax1 = ax1 + G * (m1 / r3) * (x0 - x1);
ay1 = ay1 + G * (m1 / r3) * (y0 - y1);
az1 = az1 + G * (m1 / r3) * (z0 - z1);
aax1 = aax1 + G * m1 *(vx0 - vx1) / (r3))
- (3 * (vx0 - vx1) * (x0 - x1) * (x0 - x1) / (r5);
aay1 = aay1 + G * m1 *(vy0 - vy1) / (r3))
- (3 * (vy0 - vy1) * (y0 - y1) * (y0 - y1) / (r5);
aaz1 = aaz1 + G * m1 *(vz0 - vz1) / (r3))
- (3 * (vz0 - vz1) * (z0 - z1) * (z0 - z1) / (r5);
end;
end;
begin
with BodyW do begin
apply nBody;
move SW1 to Data [1,1];
move SW2 to Data [1,1];
OuterScan (InnerScan [SW2]) [SW1];
end;
end;