wrtlst.f

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

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

Variables on this page (24)

ilast, iowr, ishort, json, kttau, ldh, ldify, ldifz, ldisk, lexcom, linstb, ljdot0, lopal95, lopale06, lovstc, lovste, lovstm, lpurez, lrot, lscv, lsemic, pdisk, tdisk, wmax

Purpose

WRTLST writes a converged model to a model file, in the MODEL2 format, which a later run can read back as its starting model (STARIN through GETMODEL2). Its comment: it writes the last converged model, and the stored models every NPUNCH models. Before writing, it turns the run's physics switches into the four header labels.

How it fits

Inputs: the model and its envelope triangle (arguments), the unit IWRITE, the mixing length CMIXL and COMPMIX (the initial composition's code, carried from input to output), both now arguments; and the physics switches in NMLPAR_MOD: KTTAU, LDH, LSCV, LOPALE06, LOPAL95 for the labels, and LEXCOM, LROT, LDIFY, LDIFZ, LDISK, LINSTB, LJDOT0, LOVSTC, LOVSTE, LOVSTM, LPUREZ, LSEMIC, PDISK, TDISK, WMAX, whose values are written as they are. The units ILAST, IOWR, ISHORT come from IOUNITS_MOD.
Who calls it: the main program main.f (unit ILAST) and WRTOUT (ILAST, and IJLAST).
What it calls: PUTMODEL2, which writes the file.
Outputs: the model file on IWRITE; a line on IOWR and one on ISHORT, which omits the age when IWRITE is 11, the unit number PARMIN gives ILAST (parmin.f, line 658 ↗) (settled from the code, Nov 9, 2026: lines 150–156 ↗). The unit ILAST is rewound, so it holds only the latest model.

The labels.

Atmosphere: KTTAU 0–5 → EDD (Eddington), KS, HRA, KUR (Kurucz), ALL (Allard), K/C (Kurucz–Castelli); any other value leaves ATM unset (bug 2).
Equation of state: SAHA, plus DH (Debye–Hückel), SCV, and the OPAL 2006 tables (O6) in six-character combinations: SAHA, SAH+DH, SCV, SCV+DH, SCDHO6, OPAL06, OP6+DH; the seventh, SCV+O6, is never written (bug 1).
Low-temperature opacity ALOK: always NONE now; high-temperature HIK: OP95 or NONE.

For users (from YREC's Quick Start). Through the unit ILAST, this is the .last file (FLAST in the control namelist): the last converged model, in the same format as a starting model, so that a later job can start from it (FFIRST).

Suspected bugs, checked by running the code

Each bug is marked by what running the code showed: (changes results, large) or (small), the fix moving a run clearly beyond a null change; (stops the run); (wrong output), the fix changing only what is written; (harmless in the runs), within the noise; (never occurred), the code ran but the fault did not; (not reached); (cannot be reached), a guard prevents it; (withdrawn), the suspicion was wrong.

1. (never occurred) A lost label.

(lines 110–118 ↗) In the SCV branch, IF (LOPALE06) EOS='SCV+O6' sits inside IF (LDH) and is overwritten on the next line by 'SCV+DH', so SCV with OPAL 2006 and without Debye–Hückel is labelled SCV. What follows: a wrong equation-of-state label in the model file's flags line; GETMODEL2 reads the label back, but nothing acts on it. Only with LSCV and LOPALE06 on and LDH off. PUTSTORE has the same lines (putstore.f, suspected bug 1).

A possible fix:

wrtlst.f

@@ lines 107–118 @@
 C Determine equation of state flag, EOS
       EOS='SAHA  '
       IF (LDH) EOS='SAH+DH'
+C 11/09/26 KC via Claude
+C code changed to label SCV with OPAL 2006 and without Debye-Huckel
+C SCV+O6; the test sat inside IF (LDH), where it was overwritten.
       IF (LSCV) THEN
          EOS='SCV   '
 C 9/24/26 MHP via Claude
 C code changed to remove the OPAL 1995/2001 equations of state
-         IF (LDH) THEN
          IF (LOPALE06) EOS='SCV+O6'
+         IF (LDH) THEN
             EOS='SCV+DH'
             IF (LOPALE06) EOS='SCDHO6'
          ENDIF

Checked by running YREC (Nov 9, 2026; a line added at each suspected path counted how often it was taken, the outputs otherwise unchanged): Of 34 runs (26 of YREC's cases, from a 0.04 $M_\odot$ brown dwarf and the solar models to 9 $M_\odot$, rotation, and helium burning, and 8 made-up ones), the SCV branch ran at every model in 32, but never with LOPALE06 on and LDH off: YREC's cases use SCV only with Debye–Hückel, and the documentation recommends LDH on (docs/source/microphysics.md, line 17). Not run: LSCV and LOPALE06 on with LDH off. Intent: none found.

2. (never occurred) Unset atmosphere label.

(lines 93–106 ↗) If KTTAU is outside 0–5, ATM keeps its value from the previous call (it is SAVEd) or, at the first call, is undefined, and that is written to the model file. Only with an invalid KTTAU (the default is 0); nothing checks its range.

A possible fix:

wrtlst.f

@@ lines 103–108 @@
 C JNT 06/14
       ELSEIF (KTTAU .EQ. 5) THEN
          ATM='K/C '
+C 11/09/26 KC via Claude
+C code changed to label an atmosphere option outside 0-5.
+      ELSE
+         ATM='UNK '
       ENDIF
 C Determine equation of state flag, EOS
       EOS='SAHA  '

Checked by running YREC (Nov 9, 2026; a line added at each suspected path counted how often it was taken, the outputs otherwise unchanged): WRTLST ran at every model in 33 runs, with KTTAU from 0 to 5 (5 in a made-up run: YREC's testsuite/Test_solarGS98_kurucz with KTTAU = 5 and the Castelli–Kurucz table, as input/atmos/CastelliKurucz/README asks), never outside 0–5. Not run: an invalid KTTAU, which nothing in YREC offers. Intent: none found.

Definitions

subroutineWRTLST

lines 7–159 · GitHub

Sets the header labels for atmosphere, EOS and opacities, and writes the model to IWRITE in MODEL2 format through PUTMODEL2.

Uses modules DIMENSIONS_MOD, IOUNITS_MOD, NMLPAR_MOD

Called by (file level), WRTOUT ×2

Flowchart of its calls, in order (1)

PUTMODEL2 137

Change log

The changes since rp_1

Changed since rp_1. COMMON blocks to modules, IMPLICIT NONE; CMIXL and COMPMIX (formerly COMMON CONST3 and I2O) are now arguments; the labels of the removed options are gone (OPAL 1995 and 2001 EOS, Alexander 1995 and Kurucz 1990 low-temperature opacities, OPAL92 and LAOL89 high-temperature opacities); WRTMONTE, a former caller, is deleted.

1C 2C 3C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ 4C 06 5C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ 6C MHP 10/02 NKK removed from declared variable list 7 SUBROUTINE WRTLST(IWRITE,HCOMP,HD,HL,HP,HR,HS,HT,LC,TRIT,TRIL,PS, 8 * *TS,RS,CFENV,FTRI,TLUMX,JCORE,JENV,MODEL,M,SMASS,TEFFL,BL,HSTOT, 9 * *DAGE,DDAGE,OMEGA,CMIXL,COMPMIX) 10C 10/15/26 MHP via Claude 11C code changed to take JSON from DIMENSIONS_MOD instead of declaring 12C PARAMETER(JSON=5000) here. 13 USE DIMENSIONS_MOD 14C 10/14/26 MHP via Claude 15C code changed to convert COMMON/NMLPAR/ into NMLPAR_MOD. 16 USE NMLPAR_MOD 17C 10/14/26 MHP via Claude 18C code changed to convert COMMON/IOUNITS/ into IOUNITS_MOD. 19 USE IOUNITS_MOD 20C 10/20/26 MHP via Claude 21C code changed to replace implicit typing with IMPLICIT NONE and 22C explicit declarations. 23 IMPLICIT NONE 24 REAL*8 BL,CFENV,CMIXL,DAGE,DDAGE,FTRI,HCOMP,HD,HL,HP,HR,HS,HSTOT, 25 * * HT,OMEGA,PS,RS,SMASS,TEFFL,TLUMX,TRIL,TRIT,TS 26 INTEGER IWRITE,JCORE,JENV,M,MODEL 27 LOGICAL*4 LC 28 29C WRTLST WRITES THE CONVERGED MODEL TO LAST MODEL A(STORES LAST 30C CONVERGED MODEL) AND STORE MODELS D(EVERY NPUNCH MODELS) 31 32C WRITE MODEL OUT IN ASCII FORMAT 33 CHARACTER*6 EOS 34c CHARACTER*4 ATM, LOK, HIK, COMPMIX 35c MHP 4/25 changed LOK name to make it unique, used elsewhere 36 CHARACTER*4 ATM, ALOK, HIK, COMPMIX 37C MHP 8/25 Removed unused variables 38C CHARACTER*256 FLAOL, FPUREZ 39C CHARACTER*256 FOPALE,FOPALE01,FOPALE06 ! FcondOpacP 40 41C 9/25/26 MHP via Claude 42C code changed to consolidate all I/O logical unit numbers, previously 43C scattered across COMMON/LUOUT/, COMMON/LUNUM/, COMMON/ALEX06/, 44C COMMON/ACDPTH/, COMMON/NEWOPAC/, COMMON/LOPAL95/, COMMON/ATMOS2/, 45C COMMON/ALATM03/, COMMON/OPALEOS/, and COMMON/SCV2/, into a single 46C COMMON/IOUNITS/ 47C 10/12/26 MHP via Claude 48C code changed to extend COMMON/NMLPAR/ to 234 members, adding the 49C namelist variables that the physics routines read from 27 COMMON 50C blocks which also held derived or working members. Those blocks keep 51C only their non-namelist members; five of them -- NEWMX, VMULT, 52C PMMWIND, CCOUT1 and ENVGEN -- had nothing left and are gone. 53C Derivation-only namelist variables from the same blocks went to 54C COMMON/NMLBYPASS/ instead, not here. 55C 56C Declared identically in all 121 files that hold any of these members. 57C Generated from one member list, so every copy is byte-identical; the 58C linker's 'size of symbol changed' warning is the check that they are. 59 DIMENSION HCOMP(15,JSON),HD(JSON),HL(JSON),HP(JSON),HR(JSON), 60 * * HS(JSON),HT(JSON),LC(JSON),TRIT(3),TRIL(3),PS(3),TS(3),RS(3), 61 * * CFENV(9),TLUMX(8),OMEGA(JSON) 62C llp 3/19/03 Add COMMON block /I2O/ for info directly transferred from 63C input to output model - starting with a code for th initial model 64C compostion (COMPMIX) 65C 10/13/26 MHP via Claude 66C code changed to replace COMMON/I2O/ with an explicit argument. COMPMIX 67C is CHARACTER, which does not belong in a COMMON block. GETMODEL2, 68C GETYREC7 and PUTMODEL2 already took it as an argument; only these three 69C files shared it through the block. 70 71C llp 3/19/03 Add required COMMON blocks such that header flags 72C ATM, EOS, HIK and LOK can be determined. 73C 10/13/26 MHP via Claude: removed the unused COMMON/ATMOS/ declaration 74C (no member referenced here; residue of the obsolete-option removals). 75C 10/13/26 MHP via Claude: removed the unused COMMON/CWIND/ declaration 76C (no member referenced here; residue of the obsolete-option removals). 77C 10/13/26 MHP via Claude: removed the unused COMMON/DEBHU/ declaration 78C (no member referenced here; residue of the obsolete-option removals). 79C 10/13/26 MHP via Claude: removed the unused COMMON/DISK/ declaration 80C (no member referenced here; residue of the obsolete-option removals). 81C 10/13/26 MHP via Claude: removed the unused COMMON/DPMIX/ declaration 82C (no member referenced here; residue of the obsolete-option removals). 83C 10/13/26 MHP via Claude: removed the unused COMMON/NEWOPAC/ declaration 84C (no member referenced here; residue of the obsolete-option removals). 85C 10/13/26 MHP via Claude: removed the unused COMMON/SCVEOS/ declaration 86C (no member referenced here; residue of the obsolete-option removals). 87 88 89 SAVE 90 91C physics flags: 92C Determine atmosphere flag, ATM 93 IF (KTTAU .EQ. 0) THEN 94 ATM='EDD ' 95 ELSEIF (KTTAU .EQ. 1) THEN 96 ATM='KS ' 97 ELSEIF (KTTAU .EQ. 2) THEN 98 ATM='HRA ' 99 ELSEIF (KTTAU .EQ. 3) THEN 100 ATM='KUR ' 101 ELSEIF (KTTAU .EQ. 4) THEN 102 ATM='ALL ' 103C JNT 06/14 104 ELSEIF (KTTAU .EQ. 5) THEN 105 ATM='K/C ' 106 ENDIF 107C Determine equation of state flag, EOS 108 EOS='SAHA ' 109 IF (LDH) EOS='SAH+DH' 110 IF (LSCV) THEN 111 EOS='SCV ' 112C 9/24/26 MHP via Claude 113C code changed to remove the OPAL 1995/2001 equations of state 114 IF (LDH) THEN 115 IF (LOPALE06) EOS='SCV+O6' 116 EOS='SCV+DH' 117 IF (LOPALE06) EOS='SCDHO6' 118 ENDIF 119 ELSE 120 IF (LOPALE06) THEN 121 EOS='OPAL06' 122 IF (LDH) EOS='OP6+DH' 123 ENDIF 124 ENDIF 125C Determine low temperature opacities flag, ALOK 126 ALOK='NONE' 127C 9/23/26 MHP via Claude 128C code changed to remove the Kurucz 1990 opacity tables 129C 9/24/26 MHP via Claude 130C code changed to remove the Alexander 1995 opacity tables 131C Determine high temperature opacities flag, HIK 132 HIK='NONE' 133 IF (LOPAL95) HIK='OP95' 134C 9/23/26 MHP via Claude 135C code changed to remove the LAOL89 opacity tables 136 137 CALL PUTMODEL2(BL,CFENV,CMIXL,DAGE,DDAGE,FTRI,HCOMP,HD,HL, 138 * * HP,HR,HS,HSTOT,HT,IWRITE,ISHORT,JCORE,JENV,LC,LEXCOM,LROT,M, 139 * * MODEL,OMEGA,PS,RS,SMASS,TEFFL,TLUMX,TRIL,TRIT,TS, 140 & & ATM,EOS,HIK,LDIFY,LDIFZ,LDISK,LINSTB,LJDOT0,ALOK, 141 & & LOVSTC,LOVSTE,LOVSTM,LPUREZ,LSEMIC,COMPMIX,PDISK,TDISK,WMAX) 142C First three lines above are YREC7 inputs 143C Last two lines are MODEL2 add-ons 144 145 146 147 WRITE(IOWR,360) MODEL,M,TEFFL,BL,DAGE 148 360 360 FORMAT(I6,' #SHELLS=', I4, ' LogTeff=',F8.5, 149 * * ' Log(L/Lsun)=',F8.5,' Age=',F12.5) 150 IF(IWRITE.EQ.11) THEN 151 WRITE(ISHORT,330) MODEL,IWRITE 152 ELSE 153 WRITE(ISHORT,340) MODEL,DAGE,IWRITE 154 ENDIF 155 330 330 FORMAT(' DUMPED MODEL',I5,' FILE',I3) 156 340 340 FORMAT(' DUMPED MODEL',I5,' AGE',F13.9,' FILE',I3) 157 IF(IWRITE.EQ.ILAST) REWIND ILAST 158 RETURN 159 END

1"""Python version of wrtlst.f (YREC 7be4133): write a converged model in 2MODEL2 format.""" 3 4import numpy as np 5 6from pyrec import iounits_mod 7from pyrec import nmlpar_mod 8from pyrec.putmodel2 import putmodel2 9 10# ATM is SAVEd: it keeps its last value when KTTAU is outside 0-5. 11_saved = {'atm': None} 12 13_ATM_LABELS = {0: 'EDD ', 1: 'KS ', 2: 'HRA ', 3: 'KUR ', 4: 'ALL ', 14 5: 'K/C '} # 5: Kurucz/Castelli (JNT 06/14) 15 16 17def wrtlst(iwrite: int, hcomp: np.ndarray, hd: np.ndarray, hl: np.ndarray, hp: np.ndarray, hr: np.ndarray, 18 hs: np.ndarray, ht: np.ndarray, lc: np.ndarray, trit: np.ndarray, tril: np.ndarray, ps: np.ndarray, 19 ts: np.ndarray, rs: np.ndarray, cfenv: np.ndarray, ftri: float, tlumx: np.ndarray, jcore: int, jenv: int, 20 model: int, m: int, smass: float, teffl: float, bl: float, hstot: float, dage: float, ddage: float, 21 omega: np.ndarray, cmixl: float, compmix: str) -> None: 22 """Set the header labels for the atmosphere, EOS and opacities, and 23 write the model to unit iwrite in MODEL2 format through putmodel2. 24 25 Translates: wrtlst.f#WRTLST (7be4133) 26 27 Arguments: the unit iwrite; the model (arrays hcomp ... omega, the 28 envelope triangle trit, tril, ps, ts, rs, cfenv, ftri, tlumx) and its 29 scalars; the mixing length cmixl; compmix, the initial composition's 30 code (a string). Nothing is returned. 31 32 Units are Fortran unit numbers; writes to the unit numbers iowr and 33 ishort are shown as print(). 34 35 """ 36 37 n = nmlpar_mod 38 39 # Physics flags. Atmosphere flag ATM (Fortran lines 91–106). 40 # Suspected bug (see the page): KTTAU outside 0-5 leaves ATM unchanged. 41 if n.kttau in _ATM_LABELS: 42 _saved['atm'] = _ATM_LABELS[n.kttau] 43 atm = _saved['atm'] 44 45 # Equation of state flag EOS (Fortran lines 107–124). 46 eos = 'SAHA ' 47 if n.ldh: 48 eos = 'SAH+DH' 49 if n.lscv: 50 eos = 'SCV ' 51 if n.ldh: 52 # Suspected bug (see the page): 'SCV+O6' is overwritten at once. 53 if n.lopale06: 54 eos = 'SCV+O6' 55 eos = 'SCV+DH' 56 if n.lopale06: 57 eos = 'SCDHO6' 58 else: 59 if n.lopale06: 60 eos = 'OPAL06' 61 if n.ldh: 62 eos = 'OP6+DH' 63 64 # Low-temperature opacity flag ALOK (no options left) (Fortran lines 65 # 125–130). 66 alok = 'NONE' 67 68 # High-temperature opacity flag HIK (Fortran lines 131–135). 69 hik = 'NONE' 70 if n.lopal95: 71 hik = 'OP95' 72 73 # The first three lines are YREC7 inputs, the last two MODEL2 add-ons 74 # (Fortran lines 137–143). 75 putmodel2(bl, cfenv, cmixl, dage, ddage, ftri, hcomp, hd, hl, 76 hp, hr, hs, hstot, ht, iwrite, iounits_mod.ishort, jcore, jenv, 77 lc, n.lexcom, n.lrot, m, 78 model, omega, ps, rs, smass, teffl, tlumx, tril, trit, ts, 79 atm, eos, hik, n.ldify, n.ldifz, n.ldisk, n.linstb, n.ljdot0, 80 alok, n.lovstc, n.lovste, n.lovstm, n.lpurez, n.lsemic, 81 compmix, n.pdisk, n.tdisk, n.wmax) 82 83 # To unit IOWR (format 360) (Fortran lines 147–149). 84 print(f"{model:6d} #SHELLS={m:4d} LogTeff={teffl:8.5f}" 85 f" Log(L/Lsun)={bl:8.5f} Age={dage:12.5f}") 86 87 # To unit ISHORT; unit 11 is ILAST (Fortran lines 150–156). 88 if iwrite == 11: 89 print(f" DUMPED MODEL{model:5d} FILE{iwrite:3d}") # 330 90 else: 91 print(f" DUMPED MODEL{model:5d} AGE{dage:13.9f} FILE{iwrite:3d}") 92 93 # Fortran line 157. 94 if iwrite == iounits_mod.ilast: 95 pass # REWIND ILAST: the unit holds only the latest model.