This posting will explore the “dynamic” parts of the Global Circulation Model GCM Model II.
We saw the overview of all subroutines, a look at the main loop at the core of it, plus how to download and compile it in earlier postings under the GCM category (see the “Category” list at the right).
While looking at the MAIN LOOP, we say it do a call to DYNAM, then do a list of other functions having to do with radiative and convective and similar activities:
>
So there is a ‘Dynamic” part “CALL DYNAM” to investigate, then we do the non-DYNAM parts. DYNAM is a set of subroutines in this program, but will be examined separately. I’m cutting the list of non-DYNAM steps down to just the most pertinent lines:C**** CONDENSTATION, SUPER SATURATION AND MOIST CONVECTION 140. CALL CONDSE 141. CALL PRECIP C**** RADIATION, SOLAR AND THERMAL 149. CALL RADIA 150. C**** SURFACE INTERACTION AND GROUND CALCULATION 157. 400 CALL SURFCE 158. CALL GROUND 159. CALL DRYCNV 160. C**** STRATOSPHERIC MOMENTUM DRAG 166. CALL SDRAG 167. C**** SEA LEVEL PRESSURE FILTER 177. CALL FILTER 181. c ** This code calls the ocean diffusion in the deep ocean mode CALL ODIFS C**** RESTRUCTURE THE OCEAN LAYERS AND ELIMINATE SMALL ICE BERGS 205. CALL OSTRUC 206.I’ve cut out the other bits from around this core just to make it easy to see what happens in what order.
This posting is looking only at that “CALL DYNAM” portion, but in more detail.
So what is in “DYNAM”? And why to comments say it is the new improved form using subroutines? To me, it looks like it iterates over the grid cells and drives the physical motion of the atmosphere. Advection.
First, it does the ‘setup the variables” part. It includes the BA94jalC9.COM COMMON block source, so we could look in there to see just what data items are being shared (so can be used to move data in and results out of the subroutine). Further down, it creates a set of arrays (DIMENSION statement) PA, PB, and PC. With three of them it implies 3-D geometry usage, but we’ll see when we get there.
It then loops over three different counters (L, J, and I) computing values of various arrays.
SUBROUTINE DYNAM
C****
C**** INTEGRATE THE DYNAMICS TERMS
C****
INCLUDE 'BA94jalC9.COM'
COMMON U,V,T,P,Q
COMMON/WORK6/WORKX(IM,JM,LM),UT(IM,JM,LM),VT(IM,JM,LM),
* TT(IM,JM,LM),PT(IM,JM),QT(IM,JM,LM)
COMMON/WORK2/UX(IM,JM,LM),VX(IM,JM,LM),TX(IM,JM,LM),PX(IM,JM)
DIMENSION PA(IM,JM),PB(IM,JM),PC(IM,JM)
DTFS=DT*2./3.
DTLF=2.*DT
NDYNO=MOD(NDYN,2)
DO 305 L=1,LM*3+1
DO 305 J=1,JM
DO 305 I=1,IM 91.
UX(I,J,L)=U(I,J,L) 92.
305 UT(I,J,L)=U(I,J,L) 93.
DO 310 L=1,LM 94.
DO 310 J=1,JM
DO 310 I=1,IM 95.
310 QT(I,J,L)=Q(I,J,L) 96.
DO 311 J=1,JM
DO 311 I=1,IM
PA(I,J)=P(I,J)
PB(I,J)=P(I,J)
311 PC(I,J)=P(I,J)
It is a very compact form, but sometimes looks a bit odd. The “DO 305” is one loop, then inside that loop is another doing “DO 305” (but with a different loop counter J and limit JM) and inside that another loop ending on 305 using yet another loop counter (index) and limit. Each of those passes through executes the text before 305 AND the text AT 305. Best Practice has that target be a CONTINUE statement to avoid ambiguity in the mind of the reader about the execution of that statement on each iteration, or at exit. So OK, they did the compact form not the most readable form.
It then does a set of calculations of air flow. We see that initially there was a “CALL DYNAM” that is now turned into a set of separate subroutines and that original line is commented out. The “flavor” of it is that the original may have been recursive (calling itself). There is a chance this commented line is just stating what the call to this program ought to look like, but it doesn’t say that, and the top comment says DYNAM was made into subroutines.
cat Mjal2cpdC9.f C DYNAM broken into subroutines
which strongly implies it was something else before…
They have a ‘forward step’ then odd and even ‘backward steps’. Overall, looks like moving the air from one cell to the next, then finding that effect and adjusting for that change into the neighbor cells.
OK, to me it doesn’t look all that efficient as a code block. LOTS of cost to set up subroutine call / return sequences and LOTS of passed parameters causing memory moves. It would be good to profile the code and see if this is where lots of time is spent. IF it isn’t that much time spent here, who cares about efficiency. IFF this is a long slow step, I think there’s room for improvement. (Though likely at the expense of readability).
I’m not going to break this up as it is pretty much the same in a few blocks one after the other. The exact advection happening changes in each pass, but the idea is the same. (P, PB, DTFS in one, P, PA, DT in another).
There is a “FLUX” component that needs defining, and a PGF subroutine that is a mystery, but the middle three are clearly advection (in three dimensions?).
C**** INITIAL FORWARD STEP, QX = Q + .667*DT*F(Q) 97.
NS=0 98.
MRCH=0 99.
C CALL DYNAM (UX,VX,TX,PX,Q,U,V,T,P,Q,DTFS) 100.
CALL AFLUX (U,V,P)
CALL ADVECM (P,PB,DTFS)
CALL ADVECV (P,UX,VX,PB,U,V,P,DTFS)
CALL ADVECT (P,TX,PB,T,DTFS)
CALL PGF (UX,VX,PB,U,V,T,P,DTFS)
IF(NDYNO.EQ.1) GO TO 320 101.
C**** INITIAL BACKWARD STEP IS ODD, QT = QT + DT*F(QX) 102.
MRCH=-1 103.
C CALL DYNAM (UT,VT,TT,PT,QT,UX,VX,TX,PX,Q,DT) 104.
CALL AFLUX (UX,VX,PB)
CALL ADVECM (P,PA,DT)
CALL ADVECV (P,UT,VT,PA,UX,VX,PB,DT)
CALL ADVECT (P,TT,PA,TX,DT)
CALL ADVECQ (P,QT,PA,Q,DT)
CALL PGF (UT,VT,PA,UX,VX,TX,PB,DT)
GO TO 360 105.
C**** INITIAL BACKWARD STEP IS EVEN, Q = Q + DT*F(QX) 106.
320 NS=1 107.
MODD5K=MOD(NSTEP+NS-NDYN+NDA5K,NDA5K) 108.
MRCH=1 109.
C CALL DYNAM (U,V,T,P,Q,UX,VX,TX,PX,QT,DT) 110.
CD DIAGA SHOULD BE CALLED HERE BUT THEN ARRAYS MUST BE CHANGED 111.
C**** ODD LEAP FROG STEP, QT = QT + 2*DT*F(Q) 112.
340 MRCH=-2 113.
C CALL DYNAM (UT,VT,TT,PT,QT,U,V,T,P,Q,DTLF) 114.
CALL AFLUX (U,V,P)
CALL ADVECM (PA,PB,DTLF)
CALL ADVECV (PA,UT,VT,PB,U,V,P,DTLF)
CALL ADVECT (PA,TT,PB,T,DTLF)
CALL ADVECQ (PA,QT,PB,Q,DTLF)
CALL PGF (UT,VT,PB,U,V,T,P,DTLF)
C**** LOAD PB TO PA
DO 341 J=1,JM
DO 341 I=1,IM
341 PA(I,J)=PB(I,J)
C**** EVEN LEAP FROG STEP, Q = Q + 2*DT*F(QT) 115.
360 NS=NS+2 116.
MODD5K=MOD(NSTEP+NS-NDYN+NDA5K,NDA5K) 117.
MRCH=2 118.
C CALL DYNAM (U,V,T,P,Q,UT,VT,TT,PT,QT,DTLF) 119.
CALL AFLUX (UT,VT,PA)
CALL ADVECM (PC,P,DTLF)
CALL ADVECV (PC,U,V,P,UT,VT,PA,DTLF)
CALL ADVECT (PC,T,P,TT,DTLF)
CALL ADVECQ (PC,Q,P,QT,DTLF)
CALL PGF (U,V,P,UT,VT,TT,PA,DTLF)
C**** LOAD P TO PC
DO 371 J=1,JM
DO 371 I=1,IM
371 PC(I,J)=P(I,J)
IF(MOD(NSTEP+NS-NDYN+NDAA,NDAA).LT.MRCH) THEN 120.
CALL DIAGA (UT,VT,TT,PB,QT)
CALL DIAGB (UT,VT,TT,PB,QT)
END IF 121.6
IF(NS.LT.NDYN) GO TO 340 122.
RETURN
END
So the main loop starts with the wind, does the rain and clouds, then moves to the radiative et. al. effects, fiddles with the stratosphere and surface pressures, then repeats. Got it.
We have an outer loop over time (years, months, days, hours) and inner looping over grid cells (Lat, Long, Alt layer) with some very inner looping over the sides of cells as ‘stuff moves’.
That’s about it.
What is AFLUX?
It calculates “Air Mass Fluxes”… but I thought we had advection routines… we’ll need to see what makes them different…
SUBROUTINE AFLUX (U,V,P) 1001.
C**** 1002.
C**** THIS SUBROUTINE CALCULATES THE HORIZONTAL AIR MASS FLUXES 1003.
C**** AND VERTICAL AIR MASS FLUXES AS DETERMINED BY U, V AND P. 1004.
C**** 1005.
Gee, it would be nice if they said what U, V, and P were. I presume P is pressure. V is likely volume. U? Need to work that one out. Thermal property? Specific heat? Who knows… Or I could be all wrong and it might be Z X and Y axis names… That’s the trouble with not having a Data Dictionary or decent comments with variables defined.
FWIW, in my very first programming class in FORTRAN IV, in the introductory comment block, we were required to list all variables used and define them…
This chunk of code does some nice things. Has ‘convergence’ calculations and handles the pole special cases.
THE biggest thing I see wrong with it is just that it is LINEAR oriented. Linear and BOXES. Watch earth.nullschool and you rapidly see that the air is dominated by circular patterns, not linear, and by bands not boxes. So the Polar Vortex, the near polar storm vortexes (sometimes making a nice hexagon). Equatorial band rise, rollover, and descend at the desert band. I can’t help the feeling that trying to map vortex circles and bands of a sphere onto boxes and lines is not going to be quite right…
With that, here’s the rest of this program with little hand holding (as I suspect it is done well and doesn’t really need a detail exploration).
INCLUDE 'BA94jalC9.COM'
COMMON/WORK1/PIT(IM,JM),SD(IM,JM,LM-1),PU(IM,JM,LM) 1007.
* ,PV(IM,JM,LM) 1008.
COMMON/WORK3/PHI(IM,JM,LM),SPA(IM,JM,LM) 1009.
COMMON/WORK4/FD(IM,0:JM+1),FLUXQ(IM),DUMMYS(IM),DUMMYN(IM) 1010.
COMMON/WORK5/DUT(IM,JM,LM),DVT(IM,JM,LM) 1011.
DIMENSION UT(IM,JM,LM),VT(IM,JM,LM),TT(IM,JM,LM), 1012.
* PA(IM,JM),PB(IM,JM),QT(IM,JM,LM),CONV(IM,JM,LM) 1013.
EQUIVALENCE (CONV,PIT) 1014.
C**** 1015.
C**** 1017.
C**** BEGINNING OF LAYER LOOP 1018.
C**** 1019.
L=LM 1020.
C**** 1021.
C**** COMPUTATION OF MASS FLUXES P,T PU PRIMARY GRID ROW 1022.
C**** ARAKAWA'S SCHEME B PV U,V SECONDARY GRID ROW 1023.
C**** 1024.
C**** COMPUTE PU, THE WEST-EAST MASS FLUX, AT NON-POLAR POINTS 1025.
2150 DO 2154 J=2,JM-1 1026.
DO 2154 I=1,IM 1027.
2154 SPA(I,J,L)=U(I,J,L)+U(I,J+1,L) 1028.
CALL AVRX (SPA(1,1,L)) 1029.
I=IM 1030.
DO 2166 IP1=1,IM 1031.
DO 2165 J=2,JM-1 1032.
2165 PU(I,J,L)=.25*DYP(J)*SPA(I,J,L)*(P(I,J)+P(IP1,J)) 1033.
2166 I=IP1 1034.
C**** COMPUTE PV, THE SOUTH-NORTH MASS FLUX 1035.
IM1=IM 1036.
DO 2172 I=1,IM 1037.
DO 2170 J=2,JM 1038.
2170 PV(I,J,L)=.25*DXV(J)*(V(I,J,L)+V(IM1,J,L))*(P(I,J)+P(I,J-1)) 1039.
2172 IM1=I 1040.
C**** COMPUTE PU*3 AT THE POLES 1041.
PUS=0. 1042.
PUN=0. 1043.
PVS=0. 1044.
PVN=0. 1045.
DO 1110 I=1,IM 1046.
PUS=PUS+U(I,2,L) 1047.
PUN=PUN+U(I,JM,L) 1048.
PVS=PVS+PV(I,2,L) 1049.
1110 PVN=PVN+PV(I,JM,L) 1050.
PUS=.25*DYP(2)*PUS*P(1,1)/FIM 1051.
PUN=.25*DYP(JM-1)*PUN*P(1,JM)/FIM 1052.
PVS=PVS/FIM 1053.
PVN=PVN/FIM 1054.
DUMMYS(1)=0. 1055.
DUMMYN(1)=0. 1056.
DO 1120 I=2,IM 1057.
DUMMYS(I)=DUMMYS(I-1)+(PV(I,2,L)-PVS) 1058.
1120 DUMMYN(I)=DUMMYN(I-1)+(PV(I,JM,L)-PVN) 1059.
PBS=0. 1060.
PBN=0. 1061.
DO 1130 I=1,IM 1062.
PBS=PBS+DUMMYS(I) 1063.
1130 PBN=PBN+DUMMYN(I) 1064.
PBS=PBS/FIM 1065.
PBN=PBN/FIM 1066.
DO 1140 I=1,IM 1067.
SPA(I,1,L)=4.*(PBS-DUMMYS(I)+PUS)/(DYP(2)*P(1,1)) 1068.
SPA(I,JM,L)=4.*(DUMMYN(I)-PBN+PUN)/(DYP(JM-1)*P(1,JM)) 1069.
PU(I,1,L)=3.*(PBS-DUMMYS(I)+PUS) 1070.
1140 PU(I,JM,L)=3.*(DUMMYN(I)-PBN+PUN) 1071.
C**** 1072.
C**** CONTINUITY EQUATION 1073.
C**** 1074.
C**** COMPUTE CONV, THE HORIZONTAL MASS CONVERGENCE 1075.
DO 1510 J=2,JM-1 1076.
IM1=IM 1077.
DO 1510 I=1,IM 1078.
CONV(I,J,L)=(PU(IM1,J,L)-PU(I,J,L)+PV(I,J,L)-PV(I,J+1,L))*DSIG(L) 1079.
1510 IM1=I 1080.
CONV(1,1,L)=-PVS*DSIG(L) 1081.
CONV(1,JM,L)=PVN*DSIG(L) 1082.
L=L-1 1083.
IF(L.GE.1) GO TO 2150 1084.
C**** 1085.
C**** END OF HORIZONTAL ADVECTION LAYER LOOP 1086.
C**** 1087.
C**** COMPUTE PIT, THE PRESSURE TENDENCY 1088.
C PIT(I,J)=CONV(I,J,1) 1089.
DO 2420 LX=2,LM 1090.
L=2+LM-LX 1091.
PIT(1,1)=PIT(1,1)+CONV(1,1,L) 1092.
PIT(1,JM)=PIT(1,JM)+CONV(1,JM,L) 1093.
DO 2420 J=2,JM-1 1094.
DO 2420 I=1,IM 1095.
2420 PIT(I,J)=PIT(I,J)+CONV(I,J,L) 1096.
C**** COMPUTE SD, SIGMA DOT 1097.
SD(1, 1,LM-1)=CONV(1, 1,LM)-DSIG(LM)*PIT(1, 1) 1098.
SD(1,JM,LM-1)=CONV(1,JM,LM)-DSIG(LM)*PIT(1,JM) 1099.
DO 2430 J=2,JM-1 1100.
DO 2430 I=1,IM 1101.
2430 SD(I,J,LM-1)=CONV(I,J,LM)-DSIG(LM)*PIT(I,J) 1102.
DO 2440 LX=2,LM-1 1103.
L=LM-LX 1104.
SD(1, 1,L)=SD(1, 1,L+1)+CONV(1, 1,L+1)-DSIG(L+1)*PIT(1, 1) 1105.
SD(1,JM,L)=SD(1,JM,L+1)+CONV(1,JM,L+1)-DSIG(L+1)*PIT(1,JM) 1106.
DO 2440 J=2,JM-1 1107.
DO 2440 I=1,IM 1108.
2440 SD(I,J,L)=SD(I,J,L+1)+CONV(I,J,L+1)-DSIG(L+1)*PIT(I,J) 1109.
DO 2450 L=1,LM-1 1110.
SDSP=SD(1,1,L) 1111.
SDNP=SD(1,JM,L) 1112.
DO 2450 I=1,IM 1113.
SD(I,1,L)=SDSP 1114.
2450 SD(I,JM,L)=SDNP 1115.
C**** 1116.
C**** 1117.
RETURN 1118.
END 1119.
PGF
But what is that PGF subroutine?
Momentum.
SUBROUTINE PGF (UT,VT,PB,U,V,T,P,DT1) 2501.
C**** 2502.
C**** THIS SUBROUTINE ADDS TO MOMENTUM THE TENDENCIES DETERMINED BY 2503.
C**** THE PRESSURE GRADIENT FORCE 2504.
I’m going to include the rest of that subroutine here, but not comment on it much. It already said what it does.
C**** 2505.
INCLUDE 'BA94jalC9.COM'
COMMON/WORK1/PIT(IM,JM),SD(IM,JM,LM-1),PU(IM,JM,LM) 2507.
1 ,PV(IM,JM,LM) 2508.
COMMON/WORK3/PHI(IM,JM,LM),SPA(IM,JM,LM) 2509.
COMMON/WORK4/FD(IM,0:JM+1),FLUXQ(IM),DUMMYS(IM),DUMMYN(IM) 2510.
COMMON/WORK5/DUT(IM,JM,LM),DVT(IM,JM,LM) 2511.
DIMENSION UT(IM,JM,LM),VT(IM,JM,LM),TT(IM,JM,LM), 2512.
* PA(IM,JM),PB(IM,JM),QT(IM,JM,LM),CONV(IM,JM,LM) 2513.
REAL*8 KAPAP1 2514.
C**** 2515.
SHA=RGAS/KAPA 2516.
KAPAP1=KAPA+1. 2517.
JMM2=JM-2 2518.
DT4=DT1/4. 2520.
C**** 2521.
C**** VERTICAL DIFFERENCING 2522.
C**** 2523.
IMAX=1 2524.
DO 3071 J=1,JM 2525.
IF(J.EQ.JM) IMAX=1 2526.
DO 3070 I=1,IMAX 2527.
SUM1=0. 2528.
SUM2=0. 2529.
SP=P(I,J) 2530.
PDN=SIG(1)*SP+PTOP 2531.
PKDN=EXPBYK(PDN) 2532.
DO 3040 L=1,LM-1 2533.
LP1=L+1 2534.
C**** CALCULATE SPA 2535.
SPA(I,J,L)=SIG(L)*SP*RGAS*T(I,J,L)*PKDN/PDN 2536.
SUM1=SUM1+SPA(I,J,L)*DSIG(L) 2537.
PUP=SIG(LP1)*SP+PTOP 2538.
PKUP=EXPBYK(PUP) 2539.
THETA=THBAR(T(I,J,LP1),T(I,J,L)) 2540.
C**** CALCULATE THE DIFFERENCE IN PHI BETWEEN ODD LEVELS LP1 AND L 2541.
PHI(I,J,LP1)=SHA*THETA*(PKDN-PKUP) 2542.
SUM2=SUM2+SIGE(LP1)*PHI(I,J,LP1) 2543.
PDN=PUP 2544.
3040 PKDN=PKUP 2545.
SPA(I,J,LM)=SIG(LM)*SP*RGAS*T(I,J,LM)*PKDN/PDN 2546.
SUM1=SUM1+SPA(I,J,LM)*DSIG(LM) 2547.
3050 PHI(I,J,1)=FDATA(I,J,1)+SUM1-SUM2 2548.
DO 3070 L=2,LM 2549.
3070 PHI(I,J,L)=PHI(I,J,L)+PHI(I,J,L-1) 2550.
3071 IMAX=IM 2551.
DO 3080 L=1,LM 2552.
DO 3080 I=2,IM 2553.
SPA(I,1,L)=SPA(1,1,L) 2554.
SPA(I,JM,L)=SPA(1,JM,L) 2555.
PHI(I,1,L)=PHI(1,1,L) 2556.
3080 PHI(I,JM,L)=PHI(1,JM,L) 2557.
C**** 2558.
C**** PRESSURE GRADIENT FORCE 2559.
C**** 2560.
DO 3340 L=1,LM 2561.
C**** NORTH-SOUTH DERIVATIVE AFFECTS THE V-COMPONENT OF MOMENTUM 2562.
IM1=IM 2563.
DO 3236 I=1,IM 2564.
DO 3234 J=2,JM 2565.
FLUX=DT4*((P(I,J)+P(I,J-1))*(PHI(I,J,L)-PHI(I,J-1,L))+ 2566.
* (SPA(I,J,L)+SPA(I,J-1,L))*(P(I,J)-P(I,J-1)))*DXV(J) 2567.
DVT(I,J,L)=DVT(I,J,L)-FLUX 2568.
3234 DVT(IM1,J,L)=DVT(IM1,J,L)-FLUX 2569.
3236 IM1=I 2570.
C**** SMOOTHED EAST-WEST DERIVATIVE AFFECTS THE U-COMPONENT 2571.
I=IM 2572.
DO 3293 IP1=1,IM 2573.
PU(I,1,L)=0. 2574.
PU(I,JM,L)=0. 2575.
DO 3280 J=2,JM-1 2576.
3280 PU(I,J,L)=(P(IP1,J)+P(I,J))*(PHI(IP1,J,L)-PHI(I,J,L))+ 2577.
* (SPA(IP1,J,L)+SPA(I,J,L))*(P(IP1,J)-P(I,J)) 2578.
3293 I=IP1 2579.
CALL AVRX (PU(1,1,L)) 2580.
DO 3294 J=2,JM 2581.
DO 3294 I=1,IM 2582.
3294 DUT(I,J,L)=DUT(I,J,L)-DT4*DYV(J)*(PU(I,J,L)+PU(I,J-1,L)) 2583.
3340 CONTINUE 2584.
C**** CALL DIAGNOSTICS 2585.
IF(MRCH.LE.0) GO TO 500 2586.
IF(MODD5K.LT.MRCH) CALL DIAG5A (6,MRCH) 2587.
IF(MODD5K.LT.MRCH) CALL DIAG9D (3,DT1,U,V) 2588.
C**** 2589.
C**** 2590.
C**** UNDO SCALING PERFORMED AT BEGINNING OF DYNAM 2591.
C**** 2592.
500 CONTINUE 2593.
DO 3410 J=2,JM-1 2594.
DO 3410 I=1,IM 2595.
3410 FD(I,J)=PB(I,J)*DXYP(J) 2596.
FDSP=PB(1, 1)*DXYP( 1) 2597.
FDNP=PB(1,JM)*DXYP(JM) 2598.
FDSP=FDSP+FDSP 2599.
FDNP=FDNP+FDNP 2600.
DO 3520 I=1,IM 2601.
FD(I, 1)=FDSP 2602.
3520 FD(I,JM)=FDNP 2603.
I=IM 2604.
DO 3540 IP1=1,IM 2605.
DO 3530 J=2,JM 2606.
RFDU=4./(FD(I,J)+FD(IP1,J)+FD(I,J-1)+FD(IP1,J-1)) 2607.
DO 3530 L=1,LM+LM 2608.
3530 UT(I,J,L)=UT(I,J,L)+DUT(I,J,L)*RFDU 2609.
3540 I=IP1 2610.
RETURN 2611.
END
ADVECTION- ADVECM, ADVECV, ADVECT, ADVECQ
I’m not going to quote them all here. Just enough to know what they do.
SUBROUTINE ADVECM (P,PA,DT1) 1501.
C**** 1502.
C**** THIS SUBROUTINE CALCULATES UPDATED COLUMN PRESSURES AS 1503.
C**** DETERMINED BY DT1 AND THE CURRENT AIR MASS FLUXES. 1504.
C****
SUBROUTINE ADVECV (PA,UT,VT,PB,U,V,P,DT1) 2001.
C**** 2002.
C**** THIS SUBROUTINE ADVECTS MOMENTUM (INCLUDING THE CORIOLIS FORCE) 2003.
C**** AS DETERMINED BY DT1 AND THE CURRENT AIR MASS FLUXES. 2004.
C**** 2005.
SUBROUTINE ADVECT (PA,TT,PB,T,DT1) 3001.
C**** 3002.
C**** THIS SUBROUTINE ADVECTS POTENTIAL TEMPERATURE 3003.
C**** 3004.
C**** 3005.
SUBROUTINE ADVECQ (PA,QT,PB,Q,DT1) 3501.
C**** 3502.
C**** THIS SUBROUTINE ADVECTS HUMIDITY AS DETERMINED BY DT1 AND THE 3503.
C**** CURRENT AIR MASS FLUXES 3504.
C**** 3505.
Air mass, momentum, temperature, and humidity. OK. Wonder why they can’t be advected in one batch, but ‘later’.
“Someday” it might be fun to go back and look at the ADVEC* set in more detail, but for now, just knowing what they move is enough. That humidity one, though, might well be where it doesn’t move enough humidity out of the upper atmosphere. Worth a ‘go fish’ at some point.
In Conclusion
With that, I’m done for the day ;-0
Anyone wants to dig around in it, now you have some pointers …
Musings from the Chiefio 
On the variable U if I had to make a wild guess I think it would be potential energy as noted in this index of common physics variables
https://en.wikipedia.org/wiki/List_of_common_physics_notations
If so its units should resolve to Joules or one of its related values (kilojoules etc.).
In heating calculations the variable U is also used to measure heat loss through a wall of known insulation value as in:
http://www.engineeringtoolbox.com/heat-loss-transmission-d_748.html
Am I missing a lot? The advection part seems to be but some airmass motion from grid location to adjacent grid location (finite element analysis), rather than compressible fluid dynamics, with strict continuum mechanics imposed (streamlines). How can such possibly be claimed as some model of Earth’s atmosphere? The intersect voltage and branch currents for a large array of electrical components, including reactive components can be modeled in this way as the various current paths are indeed finite. But Earth’s atmosphere, with 1000MPH eastern tangential momentum, near the surface at the equator?, give me a break!! This is the difference between a circuit designer and an aerodynamicist like Lockheed’s Kelly Johnson. As his boss once remarked ‘that damn Swede can see the air’! Do academic meteorologists actually have any comprehension of this atmosphere?
Larry Ledwick says: 3 January 2017 at 12:41 am
“On the variable U if I had to make a wild guess I think it would be potential energy as noted in this index of common physics variables.”
Perhaps lower case (u) energy density of the atmosphere. But what of latent heat within that same volume?
E.M. commented on the main loop
I would think that somewhere in there is an underlying assumption of the inter-dependencies of the variables. For example, The sun rising and setting affects the wind, both in intensity and direction. Would one get a different outcome if it started with the clouds, fiddled with the stratosphere, then calculated the wind, rain, and radiative effects?
Just askin’ out of ignorance.
@Larry:
Thanks for that idea. I’ll run with it and see if it helps clarity.
@Will:
Nope, not missing it, near as I can tell. I think you have grasped it clearly…
@H.R.:
I’m pretty sure that there are ordering dependencies… one hopes they got it right.
One example: Say you do rain first, then wind. That ought to get more evaporative cooling and higher moisture in air, so lower air density and more rise. The way they do it, the air gets moving, then you add water. That ought to have less upward convection from lower density and then greater wind strength from higher temperature / pressure in cells (that otherwise would be wetter and lower pressure).
Would a better cycle be found? Is it a significant difference? Who knows… One hopes the cycle time step is short enough to have it not matter much, but the time steps in this model are pretty big.
@All:
I put a comment on the prior Model E thread to this same effect:
I’ve done the compile and assemble of the parts and gotten the Makefile working. I’ve made a run job directory and populated it with a simple command and a place to put input data and results and I’ve run the program (with broken inputs so not expecting much). I’ve given it a couple of the data input files it expected (and demanded with error messages) and then found I needed to write a “rundeck”.
So if anyone has an example “Model II Rundeck” I’d love to have a copy.
I found a couple of descriptions on line, but they look to be from different eras (or are just different) so I’m in the land of reading the input READ FORMAT line, decoding the variables, and writing a rundeck. It is mostly just setting a bunch of switches and starting parameters to the model (so you can create the fictional world as you like it…) so it won’t be too hard to create one, but it takes time.
I expect to get back to this in a week or so after Valentines Day and some literal fence mending (the “drought” rained so much with wind that I’ve got between 30 and 60 feet of fence to repair…)
Early indications are that I have a working model, I just to to get all the input data sets and parameters and input rundeck worked out.