Scrambler
1
|
00001 !######################################################################### 00002 ! 00003 ! Copyright (C) 2003-2012 Department of Physics and Astronomy, 00004 ! University of Rochester, 00005 ! Rochester, NY 00006 ! 00007 ! time_step.f90 is part of AstroBEAR. 00008 ! 00009 ! AstroBEAR is free software: you can redistribute it and/or modify 00010 ! it under the terms of the GNU General Public License as published by 00011 ! the Free Software Foundation, either version 3 of the License, or 00012 ! (at your option) any later version. 00013 ! 00014 ! AstroBEAR is distributed in the hope that it will be useful, 00015 ! but WITHOUT ANY WARRANTY; without even the implied warranty of 00016 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00017 ! GNU General Public License for more details. 00018 ! 00019 ! You should have received a copy of the GNU General Public License 00020 ! along with AstroBEAR. If not, see <http://www.gnu.org/licenses/>. 00021 ! 00022 !######################################################################### 00025 00029 00030 00031 00032 !============================================================================================ 00033 ! Public Methods: RelaxTimeStep(), TestBadCFL(), MIN_TIMESTEP 00034 ! Created: 20100629 by Brandon D. Shroyer 00035 ! Notes: Not a strictly necessary module, but the functionality here is replicated 00036 ! across both MUSCL and sweep, so there's no reason to try and get it right 00037 ! twice. 00038 !============================================================================================ 00039 00042 MODULE TimeStep 00043 00044 USE DataDeclarations 00045 USE GlobalDeclarations 00046 USE EllipticDeclarations 00047 USE ExplicitDeclarations 00048 USE HyperbolicDeclarations 00049 USE ParticleDeclarations 00050 USE Timing 00051 IMPLICIT NONE 00052 SAVE 00053 PRIVATE 00054 00055 PUBLIC RelaxTimeStep, TestBadCFL, TimeStepInit, GetNextTimeStep, GetFirstTimeStep, PrintAdvance 00056 00057 REAL(KIND=qPREC), PUBLIC :: maxcfl 00058 REAL, PARAMETER :: BAD_CFL_THRESHOLD = 1d3 ! Kill the program if the CFL goes above this. 00059 REAL(KIND=qPrec), PARAMETER :: MIN_TIMESTEP = 1.d-32 00060 ! CFL-management variables and their index parameters. 00061 REAL(KIND=qPREC), PUBLIC :: overall_maxspeed, olddt 00062 INTEGER, PARAMETER :: MAX_CFL = 1 00063 INTEGER, PARAMETER :: TARGET_CFL = 2 00064 INTEGER, PARAMETER :: CFL_RELAXATION = 3 00065 ! INCLUDE 'mpif.h' 00066 CONTAINS 00067 00069 SUBROUTINE TimeStepInit() 00070 INTEGER :: iErr 00071 ! Open and read in data from the physics data file. 00072 ! READ(GLOBAL_DATA_HANDLE, NML=TimeStepData, IOSTAT=ierr) 00073 00074 ! IF (iErr /= 0) THEN 00075 ! PRINT *, "TimeStepInit() error: unable to read TimeStepData namelist." 00076 ! STOP 00077 ! END IF 00078 CLOSE(GLOBAL_DATA_HANDLE, IOStat=iErr) 00079 00080 ! IF (iErr /= 0) THEN 00081 ! PRINT *, "TimeStepInit() error: unable to close ", GLOBAL_DATA_FILE, "." 00082 ! END IF 00083 00084 00085 END SUBROUTINE TimeStepInit 00086 00088 SUBROUTINE GetFirstTimeStep 00089 overall_maxspeed=initial_maxspeed 00090 olddt=cfl_vars(TARGET_CFL) * (levels(ROOT_LEVEL)%dX/overall_maxspeed)!levels(ROOT_LEVEL)%dt 00091 maxcfl=cfl_vars(TARGET_CFL) 00092 END SUBROUTINE GetFirstTimeStep 00093 00094 00098 SUBROUTINE GetNextTimeStep(maxdt, lasttimestep) 00099 INTEGER :: level 00100 REAL(KIND=qPREC), OPTIONAL :: maxdt 00101 REAL(KIND=qPREC) :: newdt 00102 LOGICAL :: lasttimestep 00103 newdt=RelaxTimeStep() 00104 IF (PRESENT(maxdt)) THEN 00105 levels(ROOT_LEVEL)%dt = min(newdt, maxdt) 00106 if (newdt >= maxdt) THEN 00107 lasttimestep=.true. 00108 ELSEIF (newdt > .8*maxdt) THEN !should be able to guarantee two steps so shorten this one to avoid a tiny time step later. 00109 levels(ROOT_LEVEL)%dt = .5*maxdt 00110 ELSE 00111 olddt=newdt !create backup of full time step to use for future time step calculations 00112 END if 00113 ! IF (levels(ROOT_LEVEL)%dt <= maxdt) lasttimestep=.true. 00114 ELSE 00115 levels(ROOT_LEVEL)%dt = newdt 00116 olddt=newdt !create backup of full time step 00117 END IF 00118 DO level=ROOT_LEVEL+1,MaxLevel 00119 levels(level)%dt = levels(level-1)%dt / levels(level-1)%CoarsenRatio 00120 END DO 00121 RestartStep=.false. 00122 lRequestRestart=.false. 00123 END SUBROUTINE GetNextTimeStep 00124 00125 00126 00134 00135 FUNCTION RelaxTimeStep() 00136 REAL(KIND=qPrec) :: RelaxTimeStep 00137 REAL(KIND=qPrec) :: dt_wanted 00138 00139 00140 ! Establish a target timestep. 00141 dt_wanted = cfl_vars(TARGET_CFL) * levels(ROOT_LEVEL)%dX / overall_maxspeed 00142 00143 ! Create a timestep adjusted in the direction of dt_wanted. The relaxation parameter 00144 ! limits the degree to which the timestep can be shifted. 00145 IF (restartstep) THEN 00146 RelaxTimeStep=levels(ROOT_LEVEL)%dt*.5d0 00147 ELSE IF (cfl_vars(TARGET_CFL) > maxcfl) THEN 00148 RelaxTimeStep = (dt_wanted - olddt) * cfl_vars(CFL_RELAXATION) + olddt 00149 ELSE 00150 RelaxTimeStep = MAX((dt_wanted - olddt) / cfl_vars(CFL_RELAXATION) + olddt, levels(ROOT_LEVEL)%dt*cfl_vars(CFL_RELAXATION)) 00151 END IF 00152 00153 ! Make sure time step is above a minimum timestep 00154 RelaxTimeStep = MAX(RelaxTimeStep, MIN_TIMESTEP) 00155 maxspeed=0 00156 IF (lElliptic) elliptic_maxspeed=0 00157 IF (lExplicit) explicit_maxspeed=0 00158 maxcfl=0 00159 END FUNCTION RelaxTimeStep 00160 00164 SUBROUTINE TestBadCFL(lSuccess) 00165 LOGICAL :: lSuccess 00166 REAL(KIND=qPREC) :: LevelMaxSpeeds(4) 00167 INTEGER :: iErr 00168 CALL StartTimer(iTestBadCFL, -2) 00169 IF (lRequestRestart) THEN 00170 LevelMaxSpeeds(1)=1 00171 ELSE 00172 LevelMaxSpeeds(1)=0 00173 END IF 00174 00175 LevelMaxSpeeds(2)=maxval(maxspeed(0:MaxLevel)) 00176 LevelMaxSpeeds(3)=maxval(explicit_maxspeed(0:MaxLevel)) 00177 00178 # if defined HYPRE 00179 LevelMaxSpeeds(4)=maxval(elliptic_maxspeed(0:MaxLevel)) 00180 CALL StartTimer(iBarrier, 0) 00181 CALL MPI_ALLREDUCE(MPI_IN_PLACE, LevelMaxSpeeds, 4, MPI_DOUBLE_PRECISION, MPI_MAX, levels(0)%MPI_COMM, iErr) 00182 CALL StopTimer(iBarrier, 0) 00183 overall_maxspeed=maxval(LevelMaxSpeeds(2:4)) 00184 # else 00185 CALL StartTimer(iBarrier, 0) 00186 CALL MPI_ALLREDUCE(MPI_IN_PLACE, LevelMaxSpeeds, 3, MPI_DOUBLE_PRECISION, MPI_MAX, levels(0)%MPI_COMM, iErr) 00187 CALL StopTimer(iBarrier, 0) 00188 overall_maxspeed=maxval(LevelMaxSpeeds(2:3)) 00189 # endif 00190 00191 00192 IF (LevelMaxSpeeds(1) == 1) THEN !Someone requested a restart 00193 RestartStep=.true. 00194 END IF 00195 00196 IF (lParticles) overall_maxspeed=max(overall_maxspeed,ParticleMaxSpeed) 00197 00198 maxcfl=overall_maxspeed*levels(0)%dt/levels(0)%dx 00199 IF (maxcfl > cfl_vars(MAX_CFL) .OR. RestartStep) THEN 00200 IF (MPI_ID == 0) THEN 00201 IF (maxcfl > cfl_vars(MAX_CFL)) THEN 00202 write(*,'(A,E13.5,A,E13.5)') 'High CFL - restarting step. ', maxcfl, ' > ', cfl_vars(MAX_CFL) 00203 ELSE 00204 write(*,*) 'Restart requested' 00205 END IF 00206 END IF 00207 lSuccess=.false. 00208 olddt=levels(ROOT_LEVEL)%dt !Make sure we can correct the old time step 00209 ELSE 00210 lSuccess=.true. 00211 END IF 00212 CALL StopTimer(iTestBadCFL, -2) 00213 END SUBROUTINE TestBadCFL 00214 00216 SUBROUTINE PrintAdvance(n) 00217 INTEGER :: n 00218 INTEGER :: i 00219 ! IF (lRequestRestart) THEN 00220 ! write(*,*) 'processor', mpi_id, 'would like to restart', restartstep 00221 ! END IF 00222 IF (MPI_ID == 0 .AND. n <= FinestLevel) THEN 00223 WRITE(*,ADVANCE="NO",FMT='(8A2)')(' ',i=1,n) 00224 IF (lElliptic) THEN 00225 IF (lParticles) THEN 00226 WRITE(*,ADVANCE="YES",FMT=1003) n,& 00227 levels(n)%tnow, & 00228 levels(n)%dt, & 00229 max(maxspeed(n),elliptic_maxspeed(n))*levels(n)%dt/levels(n)%dx, & 00230 maxspeed(n), & 00231 elliptic_maxspeed(n), Particlemaxspeed 00232 00233 ELSE 00234 00235 WRITE(*,ADVANCE="YES",FMT=1001) n,& 00236 levels(n)%tnow, & 00237 levels(n)%dt, & 00238 max(maxspeed(n),elliptic_maxspeed(n))*levels(n)%dt/levels(n)%dx, & 00239 maxspeed(n), & 00240 elliptic_maxspeed(n) 00241 END IF 00242 ELSE 00243 IF (lParticles) THEN 00244 WRITE(*,ADVANCE="YES",FMT=1004) n,& 00245 levels(n)%tnow, & 00246 levels(n)%dt, & 00247 maxspeed(n)*levels(n)%dt/levels(n)%dx, & 00248 maxspeed(n), particlemaxspeed 00249 00250 ELSE 00251 WRITE(*,ADVANCE="YES",FMT=1002) n,& 00252 levels(n)%tnow, & 00253 levels(n)%dt, & 00254 maxspeed(n)*levels(n)%dt/levels(n)%dx, & 00255 maxspeed(n) 00256 END IF 00257 END IF 00258 1001 FORMAT(' Advanced level ',i2,' to tnext= ',e10.4,' with dt= ',e10.4,' CFL= ',e10.4,' max speed= ',e10.4,' elliptic_speed= ',e10.4) 00259 1002 FORMAT(' Advanced level ',i2,' to tnext= ',e10.4,' with dt= ',e10.4,' CFL= ',e10.4,' max speed= ',e10.4) 00260 1003 FORMAT(' Advanced level ',i2,' to tnext= ',e10.4,' with dt= ',e10.4,' CFL= ',e10.4,' max speed= ',e10.4,' elliptic_speed= ',e10.4, ' particle_speed= ',e10.4) 00261 1004 FORMAT(' Advanced level ',i2,' to tnext= ',e10.4,' with dt= ',e10.4,' CFL= ',e10.4,' max speed= ',e10.4, ' particle_speed= ',e10.4) 00262 END IF 00263 00264 END SUBROUTINE PrintAdvance 00265 00266 END MODULE TimeStep 00267