getgrid.f

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

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

Variables on this page (8)

chi, cln, dchi, echi, es1, hpttol, json, ntot

Purpose

GETGRID builds the equally spaced grid on which the diffusion equations for rotational mixing and angular momentum transport are solved. The grid variable is the combination of $\log m$, $L$ and $\log P$ that HPOINT uses to place the model's points, with one difference: HPOINT scales $\Delta_L$ by the largest $L$ in the model, GETGRID by the surface $L$.

How it fits

Who calls it: MIXGRID and ROTGRID, for one region (IBEG..IEND).
Inputs: HL ($L$ in solar units), HP ($\log P$), HS ($\log m$).
What it calls: OSPLIN, a shape-preserving osculatory quadratic spline.
Outputs: in EGRID_MOD: CHI (the grid variable at the model points), ECHI (equally spaced values), ES1 (their masses, in g), DCHI, NTOT. EQ2MOD, SECULR, MIXCOM and CHECKJ read them.
Modules: NMLPAR_MOD gives the point-spacing tolerances HPTTOL (formerly in /CTOL/, set in PARMIN); EGRID_MOD (output); $\ln 10$ from CONSTANTS_MOD.

The method

With the model's spacing criteria $\Delta_m$ = HPTTOL(2), $\Delta_L$ = HPTTOL(9)$\cdot L_\mathrm{surf}$ (with $L_\mathrm{surf}$ = HL(M)), $\Delta_P$ = HPTTOL(11), $$\chi_i = \frac{\log m_i}{\Delta_m} + \frac{L_i}{\Delta_L} - \frac{\log P_i}{\Delta_P}.$$ The equally spaced grid runs from $\chi_1$ to $\chi_N$ with $N$ = NTOT points, $\Delta\chi = (\chi_N - \chi_1)/(N-1)$, and $\log m$ is interpolated onto it with OSPLIN, then unlogged.

Definitions

subroutineGETGRID

lines 2–94 · GitHub

Equally spaced grid in $\chi = \log m/\Delta_m + L/\Delta_L - \log P/\Delta_P$, with the mass of each grid point by spline interpolation.

Uses modules CONSTANTS_MOD, DIMENSIONS_MOD, EGRID_MOD, NMLPAR_MOD

Called by MIXGRID, ROTGRID

Flowchart of its calls, in order (1)

OSPLIN 88

Change log

The changes since rp_1

Changed since rp_1. Declarations only: COMMON blocks to modules, IMPLICIT NONE.

1C$$$$$$ 2 SUBROUTINE GETGRID(HL,HP,HS,IBEG,IEND,M) 3C 10/15/26 MHP via Claude 4C code changed to take JSON from DIMENSIONS_MOD instead of declaring 5C PARAMETER(JSON=5000) here. 6 USE DIMENSIONS_MOD 7C 10/14/26 MHP via Claude 8C code changed to convert COMMON/EGRID/ into EGRID_MOD. 9 USE EGRID_MOD 10C 10/14/26 MHP via Claude 11C code changed to convert COMMON/NMLPAR/ into NMLPAR_MOD. 12 USE NMLPAR_MOD 13C 10/14/26 MHP via Claude 14C code changed to convert COMMON/CONSTANTS/ into CONSTANTS_MOD. 15 USE CONSTANTS_MOD 16C 10/20/26 MHP via Claude 17C code changed to replace implicit typing with IMPLICIT NONE and 18C explicit declarations. 19 IMPLICIT NONE 20 REAL*8 DL,DM,DP,HL,HP,HS,YTAB 21 INTEGER I,IBEG,IEND,II,M,NTAB 22C 9/25/26 MHP via Claude 23C code changed to merge COMMON/CONST/, COMMON/CONST1/, COMMON/CONST2/, 24C and COMMON/CONST3/ into a single COMMON/CONSTANTS/ 25C 10/12/26 MHP via Claude 26C code changed to extend COMMON/NMLPAR/ to 234 members, adding the 27C namelist variables that the physics routines read from 27 COMMON 28C blocks which also held derived or working members. Those blocks keep 29C only their non-namelist members; five of them -- NEWMX, VMULT, 30C PMMWIND, CCOUT1 and ENVGEN -- had nothing left and are gone. 31C Derivation-only namelist variables from the same blocks went to 32C COMMON/NMLBYPASS/ instead, not here. 33C 34C Declared identically in all 121 files that hold any of these members. 35C Generated from one member list, so every copy is byte-identical; the 36C linker's 'size of symbol changed' warning is the check that they are. 37C 10/13/26 MHP via Claude: removed the unused COMMON/CTOL/ declaration 38C (no member referenced here; residue of the obsolete-option removals). 39C OUTPUT VARIABLES 40C INPUT VARIABLES 41 DIMENSION HL(JSON),HP(JSON),HS(JSON) 42C LOCAL VARIABLES 43 DIMENSION YTAB(JSON) 44 SAVE 45C THIS ROUTINE FORMS THE SUM OF CHI, DEFINED AS 46C CHI = LOG(M)/DM + L/(LTOT*DL) - LOG(P)/DP, WHERE DM, DP, AND DL ARE 47C THE POINT SPACINGS FOR LOG MASS, LUMINOSITY, AND PRESSURE RESPECTIVELY. 48C THIS ROUTINE TRANSFORMS THESE VARIABLES TO AN EQUALLY SPACED GRID FROM 49C THE FIRST POINT TO THE LAST POINT. IT RETURNS A SET OF EQUALLY SPACED 50C CHI VALUES AND THEIR LOCATION IN MASS. 51C 52C USE THE MODEL CRITERIA FOR ASSIGNING THE SPACING BETWEEN GRID POINTS. 53 DM = HPTTOL(2) 54 DL = HPTTOL(9)*HL(M) 55 DP = HPTTOL(11) 56 NTAB = IEND - IBEG + 1 57 DO II = 1, NTAB 58 I = IBEG + II - 1 59 CHI(II) = HS(I)/DM+HL(I)/DL-HP(I)/DP 60 END DO 61C DEFINE NTOT EQUAL TO NTAB 62 NTOT = NTAB 63C TOTAL NUMBER OF ZONES IS MODULUS OF FINAL CHI PLUS ONE. 64C NTOT = INT(CHI(NTAB)-CHI(1))+1 65C FOR ROTATION PURPOSES, DEFINE THE MINIMUM NUMBER OF 66C EQUALLY SPACED SHELLS AS 3. 67 NTOT = MAX(NTOT,3) 68C EQUALLY SPACED INCREMENT IN CHI 69C CHT = 0.5D0*(CHI(NTAB)+CHI(NTAB-1)) 70C CHB = 0.5D0*(CHI(2)+CHI(1)) 71C DCHI = (CHT-CHB)/FLOAT(NTOT-2) 72 DCHI = (CHI(NTAB)-CHI(1))/DFLOAT(NTOT-1) 73C ASSIGN VECTOR OF EQUALLY SPACED CHI 74C ECHI(1) = CHB - 0.5D0*DCHI 75C DO I = 2,NTOT 76C ECHI(I) = ECHI(I-1)+DCHI 77C END DO 78 ECHI(1) = CHI(1) 79 DO I = 2,NTOT 80 ECHI(I) = ECHI(I-1)+DCHI 81 END DO 82C NOW ASSIGN MASSES TO THE NEW EQUALLY SPACED GRID POINTS. 83C PERFORM INTERPOLATION IN LOG MASS - CONSISTENT WITH CHI DEFINED 84C IN LOG M, LOG P, ETC. 85 DO I = 1,NTAB 86 YTAB(I) = HS(IBEG+I-1) 87 END DO 88 CALL OSPLIN(ECHI,ES1,CHI,YTAB,NTAB,NTOT) 89C TRANSFORM TO PHYSICAL MASS (GM) 90 DO I = 1,NTOT 91 ES1(I) = EXP(CLN*ES1(I)) 92 END DO 93 RETURN 94 END

1"""Python version of getgrid.f (YREC 7be4133): the equally spaced grid 2for the rotational mixing and angular momentum diffusion.""" 3 4import math 5 6import numpy as np 7 8from pyrec import constants_mod 9from pyrec import dimensions_mod 10from pyrec import egrid_mod 11from pyrec import nmlpar_mod 12from pyrec.osplin import osplin 13 14# SAVE: the work array of log m on the region's model points. 15ytab = np.zeros(dimensions_mod.json + 1) 16 17 18def getgrid(hl: np.ndarray, hp: np.ndarray, hs: np.ndarray, ibeg: int, iend: int, m: int) -> None: 19 """Build the equally spaced grid in 20 chi = log m / D_m + L / D_L - log P / D_P over the region ibeg..iend, 21 with as many points as the region has, and the mass at each grid 22 point by OSPLIN's spline in chi. 23 24 Translates: getgrid.f#GETGRID (7be4133) 25 26 Arguments: hl (L, solar units), hp (log P), hs (log m), ibeg, iend 27 (the region), m (number of shells). 28 Changed in place: egrid_mod.chi, echi, es1; egrid_mod.dchi and 29 egrid_mod.ntot are set. 30 31 Returns nothing. 32 33 """ 34 35 e = egrid_mod 36 hpttol = nmlpar_mod.hpttol 37 38 # The point-spacing criteria, and chi at the model points (Fortran 39 # lines 53–60). 40 dm = hpttol[2] 41 dl = hpttol[9] * hl[m] 42 dp = hpttol[11] 43 ntab = iend - ibeg + 1 44 for ii in range(1, ntab + 1): 45 i = ibeg + ii - 1 46 e.chi[ii] = hs[i] / dm + hl[i] / dl - hp[i] / dp 47 48 # Question (see the page): as many grid points as model points 49 # (Fortran lines 61–67). 50 e.ntot = ntab 51 e.ntot = max(e.ntot, 3) 52 53 # The equally spaced chi (Fortran lines 68–81). 54 e.dchi = (e.chi[ntab] - e.chi[1]) / float(e.ntot - 1) 55 e.echi[1] = e.chi[1] 56 for i in range(2, e.ntot + 1): 57 e.echi[i] = e.echi[i - 1] + e.dchi 58 59 # log m at the grid points, then unlogged (Fortran lines 82–92). 60 for i in range(1, ntab + 1): 61 ytab[i] = hs[ibeg + i - 1] 62 osplin(e.echi, e.es1, e.chi, ytab, ntab, e.ntot) 63 for i in range(1, e.ntot + 1): 64 e.es1[i] = math.exp(constants_mod.cln * e.es1[i])