ptime.f

138 lines · 1 definitions · on GitHub (7be4133)

Python version: ptime.py, 84 lines, on the right (switch to Python above the code).

Variables on this page (6)

atime, hlo, hpo, hro, hto, json

Purpose

Limits the time step by how much the structure changed from the previous model to this one. Its comment: the changes "in T, P, R, and L [must] be less than the user parameters ATIME(8)–ATIME(11) at any model point"; this criterion "dominates in the pre-main sequence, and this SR must be used (LPTIME=T) to perform pre-main sequence evolution."

How it fits

Inputs: the last step DELTS; the current model (HL, HP, HR, HT, M) and the previous one (HLO, HPO, HRO, HTO in OLDMOD_MOD); the limits ATIME in NMLPAR_MOD (set by PARMIN, possibly through REMAP).
Who calls it: HTIMER, when LPTIME; otherwise HTIMER uses $10^{20}$.
Output: PDELT, the proposed step; HTIMER takes the smallest of its candidates.

The method

Over all shells, the largest changes $\Delta\log P$, $\Delta\log T$, $\Delta\log r$, and the relative luminosity change $\Delta L / \bar L$ (with $\bar L$ the mean of the shell and its inner neighbour). Then $$f = \max\left(\frac{\Delta\log P}{a_9},\ \frac{\Delta\log T}{a_8},\ \frac{\Delta\log r}{a_{10}},\ \frac{\Delta L/\bar L}{a_{11}\ln 10}\right),\qquad \frac{1}{a_{13}} \le f \le a_{13},\qquad \texttt{PDELT} = \frac{\Delta t}{f},$$ with $a_k$ = ATIME(k), and $f = 1$ if nothing changed. The shells where the maxima occur are recorded in NHMAX, which nothing reads (unused?). The header still calls ATIME(12) the "maximum percentage change in the timestep"; since MHP 10/14 the code uses ATIME(13) as the factor, and REMAP now sets ATIME(12) from TIME_DY_SHELL, the helium-shell limit that YTIME reads, and ATIME(13) from TIME_MAX_DT_FRAC: The header's line is older wording (settled from the code, Nov 9, 2026: line 68 ↗; lines 128–135 ↗; remap.f, lines 274–275 ↗; remap.f, lines 291–294 ↗).

Definitions

subroutinePTIME

lines 6–138 · GitHub

Uses modules DIMENSIONS_MOD, NMLPAR_MOD, OLDMOD_MOD

Called by HTIMER

Compiler warnings (1)

gfortran -Wall -Wextra, compiling only; nothing was run.

line 127: Equality comparison for REAL(8) [-Wcompare-reals]

Change log

The changes since rp_1

Changed since rp_1. Declarations only: COMMON blocks to modules (/OLDMOD/ to OLDMOD_MOD, ATIME of /CTLIM/ to NMLPAR_MOD; the unused /CONST/–/CONST3/ dropped), JSON from DIMENSIONS_MOD, IMPLICIT NONE.

1C 2C 3C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ 4C PTIME 5C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ 6 SUBROUTINE PTIME(DELTS,HL,HP,HR,HT,M,PDELT) 7C 10/15/26 MHP via Claude 8C code changed to take JSON from DIMENSIONS_MOD instead of declaring 9C PARAMETER(JSON=5000) here. 10 USE DIMENSIONS_MOD 11C 10/14/26 MHP via Claude 12C code changed to convert COMMON/OLDMOD/ into OLDMOD_MOD. 13 USE OLDMOD_MOD 14C 10/14/26 MHP via Claude 15C code changed to convert COMMON/NMLPAR/ into NMLPAR_MOD. 16 USE NMLPAR_MOD 17C 10/20/26 MHP via Claude 18C code changed to replace implicit typing with IMPLICIT NONE and 19C explicit declarations. 20 IMPLICIT NONE 21 REAL*8 DELTS,FAC,HL,HMAX,HP,HR,HT,PDELT,TEMP,TESTL,TESTP,TESTR, 22 * * TESTT 23 INTEGER I,M,NHMAX 24 25 26C 10/13/26 MHP via Claude: removed the unused COMMON/CONSTANTS/ declaration 27C (no member referenced here; residue of the obsolete-option removals). 28C 9/24/26 MHP via Claude 29C code changed to remove dead CMIXL from COMMON/CONST3/; it is now 30C passed as an explicit argument where actually used 31C 10/12/26 MHP via Claude 32C code changed to extend COMMON/NMLPAR/ to 234 members, adding the 33C namelist variables that the physics routines read from 27 COMMON 34C blocks which also held derived or working members. Those blocks keep 35C only their non-namelist members; five of them -- NEWMX, VMULT, 36C PMMWIND, CCOUT1 and ENVGEN -- had nothing left and are gone. 37C Derivation-only namelist variables from the same blocks went to 38C COMMON/NMLBYPASS/ instead, not here. 39C 40C Declared identically in all 121 files that hold any of these members. 41C Generated from one member list, so every copy is byte-identical; the 42C linker's 'size of symbol changed' warning is the check that they are. 43C 10/13/26 MHP via Claude 44C code changed to drop COMMON/CTLIM/, declared here but never used. 45 DIMENSION HT(JSON),HL(JSON),HP(JSON),HR(JSON) 46 DIMENSION NHMAX(4),HMAX(4) 47 SAVE 48 49C THIS SR LIMITS THE TIMESTEP BASED ON THE REQUIREMENT THAT THE CHANGES FROM 50C ONE MODEL TO THE NEXT IN T,P,R,AND L BE LESS THAN THE USER PARAMETERS 51C ATIME(8)-ATIME(11) AT ANY MODEL POINT. 52C THIS CRITERION DOMINATES IN THE PRE-MAIN SEQUENCE, AND THIS SR MUST BE 53C USED (LPTIME=T) TO PERFORM PRE-MAIN SEQUENCE EVOLUTION. 54C 55C INPUT VARIABLES 56C DECLARED VARIABLES : 57C DELTS : PREVIOUS MODEL TIMESTEP. 58C HL : RUN OF LUMINOSITY(AS A FUNCTION OF MASS) IN THE CURRENT MODEL. 59C HP : RUN OF PRESSURE IN THE CURRENT MODEL. 60C HR : RUN OF RADIUS IN THE CURRENT MODEL. 61C HT : RUN OF TEMPERATURE IN THE CURRENT MODEL. 62C M : NUMBER OF POINTS IN THE MODEL. 63C VARIABLES IN COMMON BLOCKS : 64C ATIME(8) : MAX DELTA LOG T (= 0.02 DEFAULT) 65C ATIME(9) : " " " P (= 0.04 " ) 66C ATIME(10): " " " R (= 0.02 " ) 67C ATIME(11): " DELTA L/L (= 0.02 " ) 68C ATIME(12): maximum percentage change in the timestep from the previous one 69C HLO : RUN OF LUMINOSITY(AS A FUNCTION OF MASS) IN THE PREVIOUS MODEL. 70C HPO : RUN OF PRESSURE IN THE PREVIOUS MODEL. 71C HRO : RUN OF RADIUS IN THE PREVIOUS MODEL. 72C HTO : RUN OF TEMPERATURE IN THE PREVIOUS MODEL. 73C 74C OUTPUT VARIABLES : 75C PDELT : TIME STEP FOR REQUESTED CHANGE IN P,T,R,L. 76C 77C FIND MAXIMUM ABSOLUTE TIME DIFFERENCES FOR EACH QUANTITY 78C 79C PRESSURE 80 HMAX(1)=ABS(HPO(1)-HP(1)) 81C TEMPERATURE 82 HMAX(2)=ABS(HTO(1)-HT(1)) 83C RADIUS 84 HMAX(3)=ABS(HRO(1)-HR(1)) 85C LUMINSOITY 86 IF(HL(1)+HL(2).GT.0.0D0) THEN 87 HMAX(4)=ABS((HLO(1)-HL(1))*2.D0/(HL(2)+HL(1))) 88 ELSE 89 HMAX(4) = 0.0D0 90 ENDIF 91 DO 40 I = 2,M 92 TESTP = ABS(HPO(I)-HP(I)) 93 IF(HMAX(1).LE.TESTP) THEN 94 HMAX(1) = TESTP 95 NHMAX(1) = I 96 ENDIF 97 TESTT = ABS(HTO(I)-HT(I)) 98 IF(HMAX(2).LE.TESTT) THEN 99 HMAX(2) = TESTT 100 NHMAX(2) = I 101 ENDIF 102 TESTR = ABS(HRO(I)-HR(I)) 103 IF(HMAX(3).LE.TESTR) THEN 104 HMAX(3) = TESTR 105 NHMAX(3) = I 106 ENDIF 107 IF(HL(I)+HL(I-1).GT.0.0D0) THEN 108 TESTL = ABS((HLO(I)-HL(I))*2.0D0/(HL(I)+HL(I-1))) 109 ELSE 110 TESTL = 0.0D0 111 ENDIF 112 IF(HMAX(4).LE.TESTL) THEN 113 HMAX(4) = TESTL 114 NHMAX(4) = I 115 ENDIF 116 40 40 CONTINUE 117C NOW ACTUALLY LIMIT THE TIMESTEP BY A FACTOR THAT REDUCES THE 118C TIME CHANGES IN ALL QUANTITIES TO THE PS VALUES OR LESS 119 FAC = HMAX(1)/ATIME(9) 120 IF (HMAX(2)/ATIME(8).GT.FAC) FAC=HMAX(2)/ATIME(8) 121 IF (HMAX(3)/ATIME(10).GT.FAC) FAC=HMAX(3)/ATIME(10) 122 IF (HMAX(4)/(ATIME(11)*2.3026D0).GT.FAC) THEN 123 FAC=HMAX(4)/(ATIME(11)*2.3026D0) 124 ENDIF 125C IF NO CHANGE FROM PREVIOUS MODEL,SET PDELT TO TIMESTEP 126C STORED IN THE PREVIOUS MODEL. 127 IF (FAC.EQ.0.D0) FAC=1.0D0 128C RESTRICT CHANGE IN TIMESTEP TO NO MORE THAN A FACTOR OF ATIME(12)% 129C UP OR DOWN. 130C TEMP=1.0D0+ATIME(12) 131C MHP 10/14 USE ATIME(13) AS THE GLOBAL FACTOR FOR LIMITING TIMESTEP CHANGES 132 TEMP = ATIME(13) 133 IF (FAC.GT.TEMP) FAC=TEMP 134 IF (FAC.LT.1.0D0/TEMP) FAC=1.0D0/TEMP 135 PDELT = DELTS/FAC 136 137 RETURN 138 END

1"""Python version of ptime.f (YREC 7be4133): the time step limited by the 2structure's change since the previous model.""" 3 4import numpy as np 5 6from pyrec import nmlpar_mod 7from pyrec import oldmod_mod 8 9# Locals kept between calls (SAVE) (Fortran lines 46–47). 10hmax = np.zeros(5) 11nhmax = np.zeros(5, dtype=int) 12 13 14def ptime(delts: float, hl: np.ndarray, hp: np.ndarray, hr: np.ndarray, ht: np.ndarray, m: int, pdelt: float) -> float: 15 """Propose a time step from the largest changes, over all shells, of 16 log P, log T, log r and the relative luminosity since the previous 17 model (OLDMOD_MOD), against the limits ATIME(8)-ATIME(11), with the 18 factor bounded by ATIME(13). 19 20 Translates: ptime.f#PTIME (7be4133) 21 22 Arguments: delts, the last step; hl, hp, hr, ht, the current model; 23 m, the number of shells; pdelt (set). Returns pdelt, the proposed step. 24 25 """ 26 27 o = oldmod_mod 28 atime = nmlpar_mod.atime 29 30 # The changes at the first shell (Fortran lines 80–90). 31 hmax[1] = abs(o.hpo[1] - hp[1]) 32 hmax[2] = abs(o.hto[1] - ht[1]) 33 hmax[3] = abs(o.hro[1] - hr[1]) 34 if hl[1] + hl[2] > 0.0: 35 hmax[4] = abs((o.hlo[1] - hl[1]) * 2.0 / (hl[2] + hl[1])) 36 else: 37 hmax[4] = 0.0 38 39 # The largest changes over all shells; NHMAX, where each occurs, is 40 # read nowhere (Fortran lines 91–116). 41 for i in range(2, m + 1): # loop 40 42 testp = abs(o.hpo[i] - hp[i]) 43 if hmax[1] <= testp: 44 hmax[1] = testp 45 nhmax[1] = i 46 testt = abs(o.hto[i] - ht[i]) 47 if hmax[2] <= testt: 48 hmax[2] = testt 49 nhmax[2] = i 50 testr = abs(o.hro[i] - hr[i]) 51 if hmax[3] <= testr: 52 hmax[3] = testr 53 nhmax[3] = i 54 if hl[i] + hl[i - 1] > 0.0: 55 testl = abs((o.hlo[i] - hl[i]) * 2.0 / (hl[i] + hl[i - 1])) 56 else: 57 testl = 0.0 58 if hmax[4] <= testl: 59 hmax[4] = testl 60 nhmax[4] = i 61 62 # The factor that brings every change within its limit (Fortran lines 63 # 119–127). 64 fac = hmax[1] / atime[9] 65 if hmax[2] / atime[8] > fac: 66 fac = hmax[2] / atime[8] 67 if hmax[3] / atime[10] > fac: 68 fac = hmax[3] / atime[10] 69 # Question (see the page): ATIME(11) times ln 10, a limit on d log L. 70 if hmax[4] / (atime[11] * 2.3026) > fac: 71 fac = hmax[4] / (atime[11] * 2.3026) 72 if fac == 0.0: 73 fac = 1.0 74 75 # The factor kept within [1/atime[13], atime[13]] (Fortran lines 76 # 132–135). 77 temp = atime[13] 78 if fac > temp: 79 fac = temp 80 if fac < 1.0 / temp: 81 fac = 1.0 / temp 82 pdelt = delts / fac 83 84 return pdelt