Jose A. Rodriguez, Project 3, AMS 530

Jose A. Rodriguez, Project 3, AMS 530

Problem 3.1

I.  Algorithm

The program MOLECDYN calculates the forces, velocities and positions of a large

number of particles. The forces are calculated from simple Coulomb interactions

between the particles. The positions are determined from the Newton’s second

law of movement. The program calculates each step using P processors in parallel.

This is the sequence of steps:

·  Assign charges to the N particles using P processors. Done in the

subroutine charges. Random charges are generated in the [-1,1]

interval. Half of them positive, half negative.

·  Assign initial coordinates to the N particles using P processors. Done in the

subroutines coordx, coordy and coordz. The particles are randomly located

in a square of 10x10x10 (10-8 meter units). This in the central or working

box. In boxes around, the particles have similar relative positions which

can be obtained by applying the corresponding translation vectors

{Tx,Ty,Tz}.

·  Assign initial velocities to the N particles using P processors. Done in the

subroutine vel.The velocities are assigned assuming a Gaussian

distribution centered in (0,0,0) with a variance of 4.

·  Compute forces for each of the N particles by P processors. The charges

and coordinates of each particle are broadcasted to P processors from the

root. Each processor calculates the total force on (N/P) particles in a

central or working box. Done in subroutine forces assuming periodic

boundary conditions. The particles in the central box interact with

particles in 26 nearest-neighbor boxes located around. Interactions with

particles at longer distances are negligible.

* Compute new positions and velocities for N particles using P processors.

Each processor calculates the position and velocities for (N/P) particles

in the central box. It is assumed that particles in the surrounding boxes

move in an identical way (i.e. number of particles in each box constant).

Newton’s equation is solved numerically and the new X,Y,Z positions

are corrected if necessary assuming periodic boundary conditions. Done

in subroutine velnp. The differential equations are solved numerically,

using finite differences. To find the new coordinates a second-order

Taylor (or Runge-Kutta) approximation is used. For the velocities a

first-order Taylor (or Euler) approximation is used.

·  The last three steps are repeated 100 times for a given time step or

increment, in this case 10-6 seconds.

II.  The program

C23456789123456789123456789

PROGRAM MOLECDYN

PARAMETER (LENGTH = 10000)

INCLUDE 'mpif.h'

DIMENSION Q(LENGTH), COX(LENGTH), COY(LENGTH),

&COZ(LENGTH),FTX(LENGTH),FTY(LENGTH),FTZ(LENGTH),

&FTSX(LENGTH),FTSY(LENGTH),FTSZ(LENGTH),QW(LENGTH),

&COXW(LENGTH), COYW(LENGTH), COZW(LENGTH),

&FTTX(LENGTH),FTTY(LENGTH),FTTZ(LENGTH),VX(LENGTH),

&VY(LENGTH),VZ(LENGTH),NCX(LENGTH),NCY(LENGTH),

&NCZ(LENGTH)

INTEGER SIZE, MY_RANK, ROOT

INTEGER N_SIZE

DOUBLE PRECISION STARTTIME, ENDTIME, EXETIME,

&FTX,FTY,FTZ,FTSX,FTSY,FTSZ,FTTX,FTTY,FTTZ,

&COX,COY,COZ,COXW,COYW,COZW,VX,VY,VZ

N= LENGTH

ROOT= 0

C

C

CALL MPI_Init(ierr)

STARTTIME= MPI_Wtime()

CALL MPI_Comm_rank(MPI_comm_world,MY_RANK,ierr)

CALL MPI_Comm_size(MPI_comm_world,SIZE,ierr)

N_SIZE= N/SIZE

C

IF(MY_RANK.EQ.ROOT)THEN

C

C OBTAIN CHARGES, N_SIZE PARTICLES PER PROCESSOR

C

CALL CHARGES (N_SIZE,Q)

C

C GATHER CHARGES

C

CALL MPI_Gather(Q,N_SIZE,MPI_REAL,Q,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

C

C OBTAIN INITIAL COORDINATES, N_SIZE PARTICLES PER

C PROCESSOR

C

CALL COORDX(N_SIZE,COX)

CALL COORDY(N_SIZE,COY)

CALL COORDZ(N_SIZE,COZ)

C

C GATHER COORDINATES

C

CALL MPI_Gather(COX,N_SIZE,MPI_REAL,COX,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(COY,N_SIZE,MPI_REAL,COY,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(COZ,N_SIZE,MPI_REAL,COZ,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

C

DO 23 M=1,N

QW(M)= Q(M)

COXW(M)= COX(M)

COYW(M)= COY(M)

COZW(M)= COZ(M)

23 CONTINUE

C

C OBTAIN INITIAL VELOCITIES, N_SIZE PARTICLES

C PER PROCESSOR

C

CALL VEL(N_SIZE,VX)

CALL VEL(N_SIZE,VY)

CALL VEL(N_SIZE,VZ)

C

C GATHER VELOCITIES

C

CALL MPI_Gather(VX,N_SIZE,MPI_REAL,VX,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(VY,N_SIZE,MPI_REAL,VY,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(VZ,N_SIZE,MPI_REAL,VZ,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

C

C INTEGRATION STEPS

C

NSTEPS= 100

DO 555 ML= 1,NSTEPS

C

C

C CALCULATE THE FORCES

C

C First broacast charges, velocities and coordinates

C to the processors

C

CALL MPI_Bcast(Q,N,MPI_REAL,root,MPI_Comm_World,ierr)

CALL MPI_Bcast(COX,N,MPI_REAL,root,MPI_Comm_World,ierr)

CALL MPI_Bcast(COY,N,MPI_REAL,root,MPI_Comm_World,ierr)

CALL MPI_Bcast(COZ,N,MPI_REAL,root,MPI_Comm_World,ierr)

CALL MPI_Scatter(QW,N_SIZE,MPI_REAL,QW,N_SIZE,MPI_REAL,

&ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Scatter(COXW,N_SIZE,MPI_REAL,COXW,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Scatter(COYW,N_SIZE,MPI_REAL,COYW,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Scatter(COZW,N_SIZE,MPI_REAL,COZW,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Scatter(FTTX,N_SIZE,MPI_REAL,FTTX,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Scatter(FTTY,N_SIZE,MPI_REAL,FTTY,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Scatter(FTTZ,N_SIZE,MPI_REAL,FTTZ,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Scatter(VX,N_SIZE,MPI_REAL,VX,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Scatter(VY,N_SIZE,MPI_REAL,VY,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Scatter(VZ,N_SIZE,MPI_REAL,VZ,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

C

C Calculate forces for N_SIZE particles in each processor

C

CALL FORCES(N,N_SIZE,COX,COY,COZ,Q,COXW,COYW,COZW,QW,

&FTTX,FTTY,FTTZ)

C

C Gather the results

C

CALL MPI_Gather(FTTX,N_SIZE,MPI_REAL,FTX,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(FTTY,N_SIZE,MPI_REAL,FTY,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(FTTZ,N_SIZE,MPI_REAL,FTZ,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

C

C Calculate velocities and new X,Y,Z positions for

C N_SIZE particles in each processor

C

CALL VELNP(N_SIZE,COXW,COYW,COZW,

&FTTX,FTTY,FTTZ,VX,VY,VZ,NCX,NCZ,NCY)

C

C Gather results

C

CALL MPI_Gather(COXW,N_SIZE,MPI_REAL,COXW,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(COXW,N_SIZE,MPI_REAL,COX,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(COYW,N_SIZE,MPI_REAL,COYW,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(COYW,N_SIZE,MPI_REAL,COY,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(COZW,N_SIZE,MPI_REAL,COZW,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(COZW,N_SIZE,MPI_REAL,COZ,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(VX,N_SIZE,MPI_REAL,VX,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(VY,N_SIZE,MPI_REAL,VY,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

CALL MPI_Gather(VZ,N_SIZE,MPI_REAL,VZ,N_SIZE,

&MPI_REAL,ROOT,MPI_COMM_WORLD,ierr)

555 CONTINUE

ENDTIME= MPI_Wtime()

IF(MY_RANK.EQ.ROOT) CAll FWRITE(STARTTIME,ENDTIME,N,

&SIZE,COX,COY,COZ,Q,FTX,FTY,FTZ)

END IF

call MPI_Finalize(ierr)

END

SUBROUTINE VELNP(N_SIZE,COXW,COYW,COZW,

&FTTX,FTTY,FTTZ,VX,VY,VZ,NCX,NCZ,NCY)

DIMENSION COXW(N_SIZE),COYW(N_SIZE),COZW(N_SIZE),

&FTTX(N_SIZE),FTTY(N_SIZE),FTTZ(N_SIZE),VX(N_SIZE),

&VY(N_SIZE),VZ(N_SIZE),NCX(N_SIZE),NCY(N_SIZE),

&NCZ(N_SIZE)

DOUBLE PRECISION FTTX,FTTY,FTTZ,DH,COXW,COYW,

&COZW,CX,ACX,CY,ACY,CZ,ACZ,VX,VY,VZ

C

C The time interval is 10-6 sec

C

DH=1.0E-06

C

C Calculate new positions for each particle. V=(dr/dt).

C According to Taylor's approximation r(t+1)= r(t)+

C h*r'(t) + (h**2)*r''(t)/2. V(t)= r'(t), F(t)=mr''(t),

C m= 1. Thus, r(t+1)= r(t) + h*V(t) + (h**2)*F(t)/2

C

DO 25 K=1,N_SIZE

CX= COXW(K) + DH*VX(K)+ (DH**2)*FTTX(K)/2

C

C Correct position according to periodic BC

C

COXW(K)= CX

ACX= ABS(CX)

IF(ACX.LE.10) GO TO 11

DIV= CX/10

IDIV= DIV

IF(CX.LT.0) GO TO 12

COXW(K)= CX - (IDIV)*10

GO TO 14

12 COXW(K)= (IDIV-1)*(-10) + CX

GO TO 14

11 IF(CX.GE.0) GO TO 14

COXW(K)= 10 + CX

14 CY= COYW(K) + DH*VY(K) + (DH**2)*FTTY(K)/2

C

C Correct position according to periodic BC

C

COYW(K)= CY

ACY= ABS(CY)

IF(ACY.LE.10) GO TO 15

DIV= CY/10

IDIV= DIV

IF(CY.LT.0) GO TO 16

COYW(K)= CY - (IDIV)*10

GO TO 17

16 COYW(K)= (IDIV-1)*(-10) + CY

GO TO 17

15 IF(CY.GE.0) GO TO 17

COYW(K)= 10 + CY

17 CZ= COZW(K) + DH*VZ(K) + (DH**2)*FTTZ(K)/2

C

C Correct position according to periodic BC

C

COZW(K)= CZ

ACZ= ABS(CZ)

IF(ACZ.LE.10) GO TO 18

DIV= CZ/10

IDIV= DIV

IF(CZ.LT.0) GO TO 19

COZW(K)= CZ - (IDIV)*10

GO TO 20

19 COZW(K)= (IDIV-1)*(-10) + CZ

GO TO 20

18 IF(CZ.GE.0) GO TO 20

COZW(K)= 10 + CZ

20 XX= CZ

25 CONTINUE

C

C Calculate velocities for each particle. F=m(dV/dt).

C According to Euler's approximation V(t+1)= V(t) +

C h*F(t)

C

DO 50 I=1,N_SIZE

VX(I)= VX(I) + DH*FTTX(I)

VY(I)= VY(I) + DH*FTTY(I)

VZ(I)= VZ(I) + DH*FTTZ(I)

50 CONTINUE

RETURN

END

SUBROUTINE FORCES(N,N_SIZE,COX,COY,COZ,Q,COXW,COYW,

&COZW,QW,FTX,FTY,FTZ)

DIMENSION COX(N),COY(N),Q(N),X(3),QW(N_SIZE),

&COZ(N),COXW(N_SIZE),COYW(N_SIZE),COZW(N_SIZE),

&FTX(N_SIZE),FTY(N_SIZE),FTZ(N_SIZE)

DOUBLE PRECISION FIX,FIY,FIZ,FTX,FTY,FTZ,CC,

&COX,COY,COZ,COXW,COYW,COZW

C

C Force= CQ*Q/(r**2), where force is calculated in

C in newton units. C= 8.9874 E+9 newton*(meter**2)/

C (coul**2)

C

CC= 2.3 E-12

DO 30 K= 1,N_SIZE

FTX(K)= 0.00

FTY(K)= 0.00

FTZ(K)= 0.00

C

C We will assume periodic boundary conditions.

C The particles in a box interact with particles

C inside the box and particles in 26 boxes

C surrounding the central box. We will assume

C that Coulomb interactions with particles in

C boxes located at larger distances are

C negligible. NA, NB and NC are translation

C vectors.

C

DO 32 IA=1,3

NA= (IA - 2)*10

DO 34 IB=1,3

NB= (IB - 2)*10

DO 36 IC=1,3

NC= (IC - 2)*10

DO 40 L= 1,N

IF(L.EQ.K) GO TO 40

X(1)= COXW(K) - COX(L) - NA

X(2)= COYW(K) - COY(L) - NB

X(3)= COZW(K) - COZ(L) - NC

DIS= X(1)*X(1) + X(2)*X(2) + X(3)*X(3)

IF(DIS.EQ.0) GO TO 40

FM= CC*QW(K)*Q(L)/DIS

RNORM= SQRT(DIS)

FIX = FM*X(1)/RNORM

FIY = FM*X(2)/RNORM

FIZ = FM*X(3)/RNORM

FTX(K)= FTX(K) + FIX

FTY(K)= FTY(K) + FIY

FTZ(K)= FTZ(K) + FIZ

40 CONTINUE

36 CONTINUE

34 CONTINUE

32 CONTINUE

30 CONTINUE

RETURN

END

SUBROUTINE CHARGES(N,Q)

DIMENSION Q(N)

C

C charges are calculated in units of "e",

C where 1e= 1.60219 E-19 coulomb

C

DO 10 I=1,N

SIGN=(-1)**I

Q(I)=SIGN*Rand(I+1)

10 CONTINUE

RETURN

END

SUBROUTINE COORDX(N,CO)

DIMENSION CO(N)

DOUBLE PRECISION CO

C

C coordinates in units of "A", where

C 1A= E-8 meters

C

DO 10 I=1,N

CO(I)=10*Rand(I)

10 CONTINUE

RETURN

END

SUBROUTINE COORDY(N,CO)

DIMENSION CO(N)

DOUBLE PRECISION CO

C

C coordinates in units of "A", where

C 1A= E-8 meters

C

DO 10 I=1,N

CO(I)=10*Rand(I+2)

10 CONTINUE

RETURN

END

SUBROUTINE COORDZ(N,CO)

DIMENSION CO(N)

DOUBLE PRECISION CO

C

C coordinates in units of "A", where

C 1A= E-8 meters

C

DO 10 I=1,N

CO(I)= 10*Rand(I+4)

10 CONTINUE

RETURN

END

SUBROUTINE VEL(N,V)

DIMENSION V(N)

DOUBLE PRECISION V

C

C The initial velocities are assumed to follow a

C gaussian distribution, f(v)=A*exp(-v**2/d**2),

C centered in 0 and with variance (d) equal to 4.

C The f(v) values were calculated for velocities

C of 0, 0.5, -0.5, 1.0, -1.0, 1.5, -1.5 ......

C 4.5, -4.5, 5.0, -5.0. Groups of particles are

C assigned these velocities following a Gaussian

C distribution.

C

C v=0

PV1= 0.18*N

NV1= PV1

DO 10 I= 1, NV1

V(I)= 0.00

10 CONTINUE

C v=0.5 or -0.5

PV2= 0.16*N

NV2= NV1 + PV2

DO 20 J= NV1+1,NV2

A=(-1)**J

V(J)= A*0.5

20 CONTINUE

C v=1.0 or -1.0

PV3= 0.16*N

NV3= NV2 + PV3

DO 30 IJ= NV2+1,NV3

A=(-1)**IJ

V(IJ)= A*1.0

30 CONTINUE

C v=1.5 or -1.5

PV4= 0.14*N

NV4= NV3 + PV4

DO 40 JI= NV3+1,NV4

A=(-1)**JI

V(JI)= A*1.5

40 CONTINUE

C v=2.0 or -2.0

PV5= 0.12*N

NV5= NV4 + PV5

DO 50 K= NV4+1,NV5

A=(-1)**K

V(K)= A*2.0

50 CONTINUE

C v=2.5 or -2.5

PV6= 0.08*N

NV6= NV5 + PV6

DO 60 KI= NV5+1, NV6

A=(-1)**KI

V(KI)= A*2.5

60 CONTINUE

C v=3.0 or 3.0

PV7= 0.06*N

NV7= NV6 + PV7

DO 70 IK= NV6+1,NV7

A=(-1)**IK

V(IK)= A*3.0

70 CONTINUE

C v=3.5 or -3.5

PV8= 0.04*N

NV8= NV7 + PV8

DO 80 L= NV7+1,NV8

A=(-1)**L

V(L)= A*3.5

80 CONTINUE

C v=4.0 or -4.0

PV9= 0.02*N

NV9= NV8 + PV9

DO 90 LI= NV8+1,NV9

A=(-1)**LI

V(LI)= A*4.0

90 CONTINUE

C v=4.5 or -4.5

PV10= 0.02*N

NV10= NV9 + PV10

DO 100 IL= NV9+1,NV10

A=(-1)**IL

V(IL)= A*4.5

100 CONTINUE

C v=5.0 or -5.0

PV11= 0.02*N

NV11= NV10 + PV11

DO 110 M= NV10+1,NV11

A=(-1)**M

V(M)= A*5.0

110 CONTINUE

RETURN

END

SUBROUTINE FWRITE(STARTTIME,ENDTIME,N,SIZE,

&COX,COY,COZ,Q,FTX,FTY,FTZ)

DIMENSION COX(N),COY(N),COZ(N),Q(N),FTX(N),

&FTY(N),FTZ(N)

DOUBLE PRECISION STARTTIME,ENDTIME,EXETIME,

&COX,COY,COZ,Q,FTX,FTY,FTZ

INTEGER SIZE

EXETIME= ENDTIME - STARTTIME

WRITE(*,100) STARTTIME, EXETIME, N, SIZE

100 FORMAT('To=', E10.5, ' Texe=', E10.5, ' N=',

& I5, ' NPro=', I4)

C DO 10 I= 1,10

C WRITE(*,200) I, COX(I), FTX(I)

C 200 FORMAT(I4, ' X=', F10.5,

C &' FX=', E15.5)

C 10 CONTINUE

RETURN

END

III.  Results

Next are shown (as an example) the calculated x coordinated and Fx force for one

of the particles in a system with a total of 21600 particles. A negative force is

acting over the particle and as a result one sees a monotonous decrease of the

x coordinate. This force→reaction response is typical of the calculations.

N= 21600 particles, P=16 processors, DT= 10-6 sec

STEP= 1

x= 5.29476 FX= -0.71793E+01

STEP= 2

x= 5.29476 FX= -0.71793E+01

STEP= 3

x= 5.29476 FX= -0.71792E+01

STEP= 4

x= 5.29475 FX= -0.71792E+01

STEP= 5

x= 5.29475 FX= -0.71792E+01

STEP= 6

x= 5.29475 FX= -0.71791E+01

STEP= 7

x= 5.29475 FX= -0.71791E+01

STEP= 8

x= 5.29475 FX= -0.71791E+01

STEP= 9

x= 5.29475 FX= -0.71791E+01

STEP= 10

x= 5.29475 FX= -0.71790E+01

STEP= 11

x= 5.29475 FX= -0.71790E+01

STEP= 12

x= 5.29474 FX= -0.71790E+01

STEP= 13

x= 5.29474 FX= -0.71790E+01

STEP= 14

x= 5.29474 FX= -0.71789E+01

STEP= 15

x= 5.29474 FX= -0.71789E+01

STEP= 16

x= 5.29474 FX= -0.71789E+01

STEP= 17

x= 5.29474 FX= -0.71789E+01

STEP= 18

x= 5.29474 FX= -0.71788E+01

STEP= 19

x= 5.29474 FX= -0.71788E+01

STEP= 20

x= 5.29474 FX= -0.71788E+01

STEP= 21

x= 5.29474 FX= -0.71788E+01

STEP= 22

x= 5.29474 FX= -0.71787E+01

STEP= 23

x= 5.29474 FX= -0.71787E+01

STEP= 24

x= 5.29474 FX= -0.71787E+01

STEP= 25

x= 5.29474 FX= -0.71787E+01

STEP= 26

x= 5.29474 FX= -0.71786E+01

STEP= 27

x= 5.29474 FX= -0.71786E+01

STEP= 28

x= 5.29474 FX= -0.71786E+01

STEP= 29

x= 5.29474 FX= -0.71786E+01

STEP= 30

x= 5.29474 FX= -0.71785E+01

STEP= 31

x= 5.29474 FX= -0.71785E+01

STEP= 32

x= 5.29474 FX= -0.71785E+01

STEP= 33

x= 5.29474 FX= -0.71785E+01

STEP= 34

x= 5.29473 FX= -0.71784E+01

STEP= 35

x= 5.29473 FX= -0.71784E+01

STEP= 36

x= 5.29473 FX= -0.71784E+01

STEP= 37

x= 5.29473 FX= -0.71784E+01

STEP= 38

x= 5.29473 FX= -0.71783E+01

STEP= 39

x= 5.29473 FX= -0.71783E+01

STEP= 40

x= 5.29473 FX= -0.71783E+01

STEP= 41

x= 5.29473 FX= -0.71783E+01

STEP= 42

x= 5.29473 FX= -0.71782E+01

STEP= 43

x= 5.29473 FX= -0.71782E+01

STEP= 44

x= 5.29473 FX= -0.71782E+01

STEP= 45

x= 5.29473 FX= -0.71782E+01

STEP= 46

x= 5.29473 FX= -0.71781E+01

STEP= 47

x= 5.29473 FX= -0.71781E+01

STEP= 48

x= 5.29473 FX= -0.71781E+01

STEP= 49

x= 5.29473 FX= -0.71781E+01

STEP= 50

x= 5.29473 FX= -0.71780E+01

STEP= 51

x= 5.29472 FX= -0.71780E+01

STEP= 52

x= 5.29472 FX= -0.71780E+01

STEP= 53

x= 5.29472 FX= -0.71779E+01

STEP= 54

x= 5.29472 FX= -0.71779E+01

STEP= 55

x= 5.29472 FX= -0.71779E+01

STEP= 56

x= 5.29472 FX= -0.71779E+01

STEP= 57

x= 5.29472 FX= -0.71778E+01

STEP= 58

x= 5.29472 FX= -0.71778E+01

STEP= 59

x= 5.29472 FX= -0.71778E+01

STEP= 60

x= 5.29472 FX= -0.71778E+01

STEP= 61

x= 5.29472 FX= -0.71777E+01

STEP= 62

x= 5.29472 FX= -0.71777E+01

STEP= 63

x= 5.29472 FX= -0.71777E+01

STEP= 64

x= 5.29472 FX= -0.71777E+01

STEP= 65

x= 5.29472 FX= -0.71776E+01

STEP= 66

x= 5.29472 FX= -0.71776E+01

STEP= 67

x= 5.29472 FX= -0.71776E+01

STEP= 68

x= 5.29472 FX= -0.71776E+01

STEP= 69

x= 5.29471 FX= -0.71775E+01

STEP= 70

x= 5.29471 FX= -0.71775E+01

STEP= 71

x= 5.29471 FX= -0.71775E+01

STEP= 72

x= 5.29471 FX= -0.71775E+01

STEP= 73

x= 5.29471 FX= -0.71774E+01

STEP= 74

x= 5.29471 FX= -0.71774E+01

STEP= 75

x= 5.29471 FX= -0.71774E+01

STEP= 76

x= 5.29471 FX= -0.71774E+01

STEP= 77

x= 5.29471 FX= -0.71773E+01

STEP= 78

x= 5.29471 FX= -0.71773E+01

STEP= 79

x= 5.29471 FX= -0.71773E+01

STEP= 80

x= 5.29471 FX= -0.71773E+01

STEP= 81

x= 5.29471 FX= -0.71772E+01

STEP= 82

x= 5.29471 FX= -0.71772E+01

STEP= 83

x= 5.29471 FX= -0.71772E+01

STEP= 84

x= 5.29471 FX= -0.71772E+01

STEP= 85

x= 5.29471 FX= -0.71771E+01

STEP= 86

x= 5.29470 FX= -0.71771E+01

STEP= 87

x= 5.29470 FX= -0.71771E+01

STEP= 88

x= 5.29470 FX= -0.71771E+01

STEP= 89

x= 5.29470 FX= -0.71770E+01

STEP= 90

x= 5.29470 FX= -0.71770E+01

STEP= 91

x= 5.29470 FX= -0.71770E+01

STEP= 92

x= 5.29470 FX= -0.71770E+01

STEP= 93

x= 5.29470 FX= -0.71769E+01

STEP= 94

x= 5.29470 FX= -0.71769E+01

STEP= 95

x= 5.29470 FX= -0.71769E+01

STEP= 96

x= 5.29470 FX= -0.71769E+01

STEP= 97

x= 5.29469 FX= -0.71768E+01

STEP= 98

x= 5.29469 FX= -0.71768E+01

STEP= 99

x= 5.29469 FX= -0.71768E+01

STEP= 100

x= 5.29469 FX= -0.71767E+01

Next are shown the computational times necessary in Galaxy to do 100