!*********************************************************************** ! ! PROGRAM FLY ! ! ! FLY Ver. 5.0: December , 2012 ! ! Multi Platform one-side comm. based code ! for TREE cosmological LSS evolution ! ! (U. BECCIANI et al.) ! ! At variance with TREECODE we adopt the standard leapfrog integrator ! of Efstathiou and Eastwood (Hockney and Eastwood, eqs. 11.59a,b). ! ! Thanks to Lars Hernquist, Institute for Advanced Study, for the ! original source code of the serial version ! !*********************************************************************** ! !======================================================================= ! USE fly_h implicit none INCLUDE 'mpif.h' ! Declaration of local variables. ! ------------------------------- INTEGER(KIND=8) :: n,i,j REAL(KIND=8) :: init_rt, init_tot,init_tot1,init_tot2,step_rt REAL(KIND=8) :: psec,tot_sec,prev_sec REAL(KIND=8) :: min_rmt_time,max_rmt_time INTEGER :: istatus(MPI_STATUS_SIZE),ind_pe,tag0=0,tag1=1,tag2=2,tag3=3,lname INTEGER :: index_min, index_max character*(MPI_MAX_PROCESSOR_NAME) hostname_me #ifdef FLASH character*(MPI_MAX_PORT_NAME) portName #endif ! ! !======================================================================= ! !----------------------------------------------------------------------------- ! Initialize state of the system. !----------------------------------------------------------------------------- ! tot_sec=0.0 prev_sec=0.0 CALL MPI_Init(ierror) CALL MPI_COMM_RANK(MPI_COMM_WORLD, me, ierror) CALL MPI_COMM_SIZE(MPI_COMM_WORLD, npes, ierror) CALL MPI_GET_PROCESSOR_NAME(hostname_me, lname, ierror) WRITE(uterm, *) "FLY - RUN. PE=",me," HOSTNAME:",hostname_me(1:lname),' npes=',npes CALL flush(uterm) init_tot1=MPI_WTIME() ! start time counter of GLOBAL code !----------------------------------------------------------------------------- ! Verify the system consistence !----------------------------------------------------------------------------- IF(N_PES.NE.npes) THEN write(uterm,*)'Error: N_PES=',N_PES,' and number of PEs= ',npes, & ' DOES NOT match',' me=',me STOP ENDIF IF(nb_loc*N_PES .lt. nbodsmax) THEN IF(me.eq.0) write(uterm,*)'Error: nbodsmax=',nbodsmax,' and number of PEs= ',npes, & ' DOES NOT match' STOP ENDIF IF(nc_loc*N_PES .lt. nbodsmax) THEN IF(me.eq.0) write(uterm,*)'Error: ncells=',ncells,' and number of PEs= ',npes, & ' DOES NOT match' STOP ENDIF IF(nb_loc .gt. 2147483646/3) THEN IF(me.eq.0) write(uterm,*)'Error: nbodsmax=',nbodsmax,' too big nb_loc=',nb_loc STOP ENDIF DO i=1,15 IF(npes .eq. 2**i) EXIT ENDDO IF(i.eq.16) THEN POW2=.FALSE. IF(me.eq.0) write(6,*)"Power of two not enabled" ELSE POW2=.TRUE. IF(me.eq.0) write(6,*)"Power of two enabled" ENDIF ! !----------------------------------------------------------------------------- ! Start FLY inizialization phase !----------------------------------------------------------------------------- init_rt = MPI_WTIME() ! start time counter of the initialization + initial load balance phase CALL null ! routine to initialize variables CALL sys_init ! system initialization: read bodies, parameter, and initilize variables !----------------------------------------------------------------------------- ! WRITE OUT INFORMATION !----------------------------------------------------------------------------- ! IF(me.EQ.0) THEN write(uterm,*) 'Running with ',N_PES,' processors and tol=',tol write(uterm,*) 'GROUPING LEVEL=',ncrit,' Body in cell group=',nbodcrit call flush(uterm) ENDIF init_rt = MPI_WTIME()-init_rt CALL MPI_BARRIER(MPI_COMM_WORLD,ierror) IF(me.eq.0) THEN write(uterm,*)'INIT: executed =',init_rt,' sec.' call flush(uterm) ENDIF !----------------------------------------------------------------------------- ! START THE STEP SYSTEM EVOLUTION FOR nsteps !----------------------------------------------------------------------------- DO 100 n=1,nsteps step_rt = MPI_WTIME() !----------------------------------------------------------------------- ! Dynamic input parameter - START !----------------------------------------------------------------------- !----------------------------------------------------------------------- ! Open data files, read input parameters and initial system state, ! and reset system parameters. !----------------------------------------------------------------------- CALL inpar_dyn !----------------------------------------------------------------------- ! Open out_32.tab file: List of programmed output (in redshift value) ! and reset system parameters. !----------------------------------------------------------------------- CALL read_redsh ! reading out32 !----------------------------------------------------------------------- ! Dynamic input parameter - STOP !----------------------------------------------------------------------- IF(halt_sim.eq.1 .or. halt_sim.eq.2) EXIT !Stop Simulation !----------------------------------------------------------------------- ! System Evolution !----------------------------------------------------------------------- CALL step(n) CALL MPI_BARRIER(MPI_COMM_WORLD,ierror) step_rt=MPI_WTIME()-step_rt psec=step_rt IF(me.eq.0) THEN write(uterm,*) '--------------- step ',tstep, & ' executed sec =',step_rt call flush(uterm) ENDIF CALL MPI_BARRIER(MPI_COMM_WORLD,ierror) !------------------------------------------------- ! WRITE STATISTICAL TIMING !------------------------------------------------- IF(me.eq.0) THEN write(uterm,*) '--------------- Local_computation: ----------------- ' write(uterm,1000) 'PE',me,': ACCGRAV ', ctota call flush(uterm) DO ind_pe=1,npes-1 CALL MPI_RECV(ctota, 1, MPI_REAL8, ind_pe, tag0 ,MPI_COMM_WORLD , istatus,ierror) write(uterm,1000) 'PE', ind_pe,': ACCGRAV ', ctota call flush(uterm) ENDDO ELSE CALL MPI_SEND(ctota,1,MPI_REAL8,PE0,tag0,MPI_COMM_WORLD, ierror) ENDIF CALL MPI_BARRIER(MPI_COMM_WORLD,ierror) CALL MPI_BARRIER(MPI_COMM_WORLD,ierror) IF(me.eq.0) THEN write(uterm,*) '--------------- Remote_computation: ----------------- ' write(uterm,1000) 'PE',me,': ACC_EX ',ctotc call flush(uterm) DO ind_pe=1,npes-1 CALL MPI_RECV(ctotc, 1, MPI_REAL8, ind_pe, tag0 ,MPI_COMM_WORLD , istatus,ierror) write(uterm,1000) 'PE', ind_pe,': ACC_EX ',ctotc call flush(uterm) ENDDO write(uterm,*) '------------------------------------------------------ ' write(uterm,*)'PARTICLE/seconds =',nbodies/psec call flush(uterm) ELSE CALL MPI_SEND(ctotc,1,MPI_REAL8,PE0,tag0,MPI_COMM_WORLD, ierror) ENDIF 1000 format(x,a,i3,2x,3(a,g18.8,2X)) 1100 format(x,3(a,g18.8,2X)) !----------------------------------------------------------------------------- ! Check for stop simulation !----------------------------------------------------------------------------- IF(halt_sim.ne.0) EXIT !Stop Simulation tot_sec=tot_sec+psec prev_sec=tot_sec+2*psec*1.20 IF(max_time .gt. 0) THEN IF(prev_sec .ge. max_time) THEN halt_sim=3 IF(me.eq.0) THEN write(uterm,*)' ' write(uterm,*)'--------------------------------------------' write(uterm,*)' CPU TIME LIMIT REACHED NEXT STEP ' write(uterm,*)'--------------------------------------------' write(uterm,*)'Actual cpu time = ',tot_sec write(uterm,*)'Forecast cpu time = ',tot_sec+2*psec*1.20 write(uterm,*)'Maximum cpu time = ',max_time write(uterm,*)'--------------------------------------------' write(uterm,*)' ' ENDIF ELSE IF(me.eq.0) THEN write(uterm,*)' ' write(uterm,*)'--------------------------------------------' write(uterm,*)' CPU TIME LIMIT COUNTING ' write(uterm,*)'--------------------------------------------' write(uterm,*)'Actual cpu time = ',tot_sec write(uterm,*)'Forecast cpu time = ',tot_sec+2*psec*1.20 write(uterm,*)'Maximum cpu time = ',max_time write(uterm,*)'--------------------------------------------' write(uterm,*)' ' ENDIF ENDIF ENDIF 100 CONTINUE ! continue to the next step system evolution !----------------------------------------------------------------------------- ! END OF STEP SYSTEM EVOLUTION and Terminate the run !----------------------------------------------------------------------------- CALL MPI_BARRIER(MPI_COMM_WORLD,ierror) init_tot2 = MPI_WTIME() init_tot = init_tot2 - init_tot1 IF((tst_max.le.tstep).and.(halt_sim.eq.0) ) halt_sim=3 IF(me.eq.0) THEN write(uterm,*)'TOTAL: executed =',init_tot,' sec.' IF(tst_max.le.tstep) THEN write(uterm,*)'SIMULATION STOPED: reached MAX STEP value' ENDIF IF(halt_sim.eq.1 .or. halt_sim.eq.2) THEN write(uterm,*)'SIMULATION STOPED: reached final readshift value' ENDIF IF(halt_sim.eq.3) THEN write(uterm,*)'SIMULATION STOPED: reached maximum cpu limit' ENDIF ENDIF #ifdef FLASH if(me == 0) then call MPI_UNPUBLISH_NAME('particleInterface', MPI_INFO_NULL, & portName, ierror) call MPI_CLOSE_PORT(portName, ierror) end if call MPI_COMM_DISCONNECT(flashComm, ierror) #endif call MPI_Finalize(ierror) STOP END