*/ Correction set for Version 2.3 of the LEM */ Latest updates made on 28 Nov 2002 */ Each set has comments from the sender */ Untested by me!! */--------------------------------------------------------- */--------------------------------------------------------- */ Mods by Marc to allow more portability */ */--------------------------------------------------------- *IDENT FORTE_CORR */--------------------------------------------------------- */ GC_IMAX is routine which is called from LEM but doesn't */ currently have a dummy fortran stub for the portable version */ as yet - marc 13/11/02 */--------------------------------------------------------- *I NO_GCOM.128 *IF DEF,PORT c--------------------------------------------------------- SUBROUTINE GC_IMAX(inum,nproc,istat,icommax) c--------------------------------------------------------- PRINT*,'ERROR: This model is meant for only one '// 1 'processor. GC_IMAX should never of have'// 2 'been accessed' STOP END *ENDIF */--------------------------------------------------------- */ Forte requires that the seed for random numbers be */ at least 4 integers long - marc 20/11/02 */--------------------------------------------------------- */APL*D PRAMETR.151 */APL! seed for forte compiler needs to be at least 4 integer */APL! long - marc 20/11/02 */APL! INTEGER, PARAMETER :: ISDSZP = 2 ! Maximum pos. size of seed array */APL INTEGER, PARAMETER :: ISDSZP = 4 ! Maximum pos. size of seed array */===================================================== */ Mods by A. Lock to stop problem with fort.15 being written */ */ Only a minor cosmetic problem */ with run_info (unit 15) being closed in BEGIN */ before the model's finished writing to it. *D BEGIN.1037,1038 *I BEGIN.1185 IF(MYPE.EQ.0)THEN ! Close unit 15 (opened in CHKPARMS) CLOSE(15) ENDIF */===================================================== */ Mods by A. Lock to call timestep checker every 10 */ steps rather than every NN which can cause things to */ blow up with a bigger NN. */ Note NNCLF gets added to the namelist but default will be as */ before *// *// Version 2: changes to make TESTCFL calls at the same steps as before, *// ie: reinstate calls every step for the first NNCFL steps. *// Will only give the same results as before if NN was constant (ie. the *// same in set-up and chain jobs). */ *IDENT NNCFL *// *// Mod to test CFL and change timestep every NNCFL steps, rather than NN *// NNCFL added to namelist CNTRL *// *I IPARMS.3 & ,NNCFL *I IPARMS.12 INTEGER NNCFL ! Number of steps between CFL calculations *D BEGIN.62 & ISTART,NN,NNDIAG,NNDUMP,NNCFL *I BEGIN.161 NNCFL = -1 *B BEGIN.166 IF (NNCFL .EQ. -1) NNCFL = NN *D NNSTEPS.411 IF (MOD(NSTEP+1,NNCFL) .EQ. 1 .OR. NSTEP+1 .LE. NNCFL) & CALL CVELGAL2(1) ! !complete calculation of UGAL,VGAL (NSTEP not yet updated) *D NNSTEPS.476 IF (MOD(NSTEP,NNCFL) .EQ. 1 .OR. NSTEP .LE. NNCFL) THEN *D TESTCFL.5,6 ! called from NNSTEPS if MOD(NSTEP,NNCFL) EQ 1 *D CVELGAL.4 ! called from STEPFLDS if MOD(NSTEP+1,NNCFL) EQ 1 *D CVELGAL2.5 ! called from NNSTEPS if MOD(NSTEP+1,NNCFL) EQ 1 *I STEPFLDS.30 *CALL IPARMS *D STEPFLDS.81 IF (MOD(NSTEP+1,NNCFL) .EQ. 1 .OR. NSTEP+1 .LE. NNCFL) & CALL CVELGAL(I,ZU,ZV,ZW) *I SMAG.20 *CALL IPARMS *D SMAG.213 IF (MOD(NSTEP+1,NNCFL) .EQ. 1 .OR. NSTEP+1 .LE. NNCFL) THEN *I SMAG2D.19 *CALL IPARMS *D SMAG2D.140 IF (MOD(NSTEP+1,NNCFL) .EQ. 1 .OR. NSTEP+1 .LE. NNCFL) THEN *I SMAG1D.19 *CALL IPARMS *D SMAG1D.121 IF (MOD(NSTEP+1,NNCFL) .EQ. 1 .OR. NSTEP+1 .LE. NNCFL) THEN */===================================================== */ Mods by Malcolm for varoius things */ */ *IDENT CORR_MKM */ previous code prevented a movie from starting at time=0.*** *D BEGIN.1103 IF( tminc.LE.0.0 .OR. time_movie.LT.0.0 )THEN */ */ allows the option of backscatter with SCT=0, which is useful in */ some circumstances *** *D BEGIN.804 IF(SCT.LT.0.0)THEN */ */*** avoids out-of-bounds indexing for array DIFF_DG, which is a */ potentially serious overwriting problem*** *I LOWERBC.197 IF(I.GE.0.and.I.le.iipep+1)THEN *I LOWERBC.200 ENDIF */ */*** counter incremented in the wrong place relative to an IPASQP.EQ.0 */ IF-block***** *D DIAG.748 *I DIAG.746 NPARMS=IHS+67 */ */*** sets the value of K, previously omitted (which made the values of */ the surface sensible and latent heat fluxes wrong if ISCBCP.EQ.2) *I DIAG.1147 K=1 */ */*** removes a potential spurious failure caused by previous */ STATUS='NEW' flag*** NOTE this is also protected by recent */ versions of lemsub *D CHKPARMS.26 OPEN(15,FORM='FORMATTED',FILE='run_info') */ */*** as far as I can tell, the error corrected by the following */ was only diagnostic*** *I BCKSCT.193 !------------------------------------------------------- ! Ensure that U,V boundary conditions remain satisfied !------------------------------------------------------- IF(INOVISP.EQ.0.AND.INOSURFP.EQ.0)THEN DO J=1,JJP SU(J,1)=-SU(J,2) SV(J,1)=-SV(J,2) ENDDO ELSE DO J=1,JJP SU(J,1)=SU(J,2) SV(J,1)=SV(J,2) ENDDO ENDIF */ */***corrects an omission in the OVRIDE1 namelist**** *I BEGIN.145 & FTSURF,FQSURF, */Here's another desirable correction to the LEM - stops you selecting */backscatter with the anelastic option , as Andy suspects the backscatter is */not coded properly in that case. */ *I CHKPARMS.321 IF(IANELP.EQ.1)THEN IF(MYPE.EQ.0)THEN WRITE(6,*)'****PARAMETER CONFLICT****' WRITE(6,*)' IANELP = ',IANELP,' IBSCATP = ',IBSCATP WRITE(6,*)' ANELASTIC OPTION PROBABLY NOT CODED CORRECTLY' WRITE(6,*)' IF BACKSCATTER IS SELECTED' ENDIF L_ERROR=.TRUE. ENDIF */! */!Glen Shutts adds: */! I mentioned an 'unfinalized MPI error' last week which actually */!prevented me from using the threads-enabled variant of the FORTRAN */!compiler. I fixed this by calling the routine GC_EXIT() at the end of */!each run using this update: */! *I EXXIT.379 IF(NPES.GT.1) THEN CALL GC_GSYNC(NPROC,INFO) CALL GC_EXIT() ! finalize MPI before quitting ENDIF */--------------------------------------------------------- */ GC_EXIT doesnt have a port routine so here is one */ Added by jon to compliment Glens mod 13/11/02 */--------------------------------------------------------- *B NO_GCOM.129 c--------------------------------------------------------- SUBROUTINE GC_EXIT(inum,nproc,istat,icommax) c--------------------------------------------------------- PRINT*,'ERROR: This model is meant for only one '// 1 'processor. GC_EXIT should never of have'// 2 'been accessed' STOP END */===================================================== */ These sets are pre November 2002 */ */ */ Taken from spec_cor */ */ Correction set for Version 2.3 of the LEM */ Correction to stop run if spectral file is */ not read correctly *I READLWSP.28 IF (ierr.ne.i_normal) THEN print*, cmessagein print*, 'Stopped in Read_lw_spec' call EXXIT(11) ! 11 is a radiation error ENDIF *I READSWSP.28 IF (ierr.ne.i_normal) THEN print*, cmessagein print*, 'Stopped in Read_sw_spec' call EXXIT(11) !11 is a radiation error ENDIF */ */------------------------------------------------------ *IDENT ALLSC230 *D ALLOCSP.26 CHARACTER*3 CTEMP ! String for spectra level *IDENT INDNC230 *D INDGNEW.45 DGSPNEW(JSPEC,N)=0.0 *IDENT INDAC230 *D INDGAV.34 DGSPAV(JSPEC,N)=0.0 *IDENT SET1C230 *D SET1D.658 IF (NMETP.GE.1.AND.IRAINP.GE.1.AND.IRAINP.LT.10) THEN *IDENT DIAGC230 */ request from A.Lock *D DIAG.527 APARMS(27)=REAL(NSERACT) */These updates correct the 2-D field diagnostics for MPP 2-D runs *I DIAG.151 REAL TEMPFLD(0:JJP+1,KKP) ! temp 2-D field for diagnostic fields *I DIAG.874 IF(IIP.GT.1)THEN *I DIAG.899 ELSE ! special case for 2-D constructed diagnostics L_INFIELD = .false. IF(IDGU.EQ.2)THEN DO K=1,KKP DO J=0,JJP+1 TEMPFLD(J,K) = U(J,K,1) + UGAL ENDDO ENDDO CALL WOUT2D(NUN, TEMPFLD, FIELD, L_INFIELD, & TEMPYZ2D, TEMPXZ2D, TEMPXY2D, & PACKXZ, PACKXY, BLK, LENBLK, ISTART ) ENDIF IF(IDGV.EQ.2)THEN DO K=1,KKP DO J=0,JJP+1 TEMPFLD(J,K) = V(J,K,1) + VGAL ENDDO ENDDO CALL WOUT2D(NUN, TEMPFLD, FIELD, L_INFIELD, & TEMPYZ2D, TEMPXZ2D, TEMPXY2D, & PACKXZ, PACKXY, BLK, LENBLK, ISTART ) ENDIF ENDIF ! IIP > 1 *I DIAG.1072 IF(IIP.GT.1)THEN *I DIAG.1083 ELSE DO K=1,KKP DO J=0,JJP+1 TEMPFLD(J,K) = TH(J,K,1) + THREF(K) ENDDO ENDDO L_INFIELD = .false. CALL WOUT2D(NUN, TEMPFLD, FIELD, L_INFIELD, & TEMPYZ2D, TEMPXZ2D, TEMPXY2D, & PACKXZ, PACKXY, BLK, LENBLK, ISTART ) ENDIF ! IIP >1 *I DIAG.1086 IF(IIP.GT.1)THEN *I DIAG.1109 ELSE DO K=1,KKP DO J=0,JJP+1 TEMPFLD(J,K) = TH(J,K,1) - OLZTHBAR(K) ENDDO ENDDO IF(IUSEQP.EQ.1)THEN DO IQ=1,NQP DO K=1,KKP DO J=0,JJP+1 TEMPFLD(J,K) = TEMPFLD(J,K) + & CQ(IQ)*THREF(K)*( Q(J,K,IQ,1) - OLZQBAR(K,IQ) ) ENDDO ENDDO ENDDO ENDIF L_INFIELD = .false. CALL WOUT2D(NUN, TEMPFLD, FIELD, L_INFIELD, & TEMPYZ2D, TEMPXZ2D, TEMPXY2D, & PACKXZ, PACKXY, BLK, LENBLK, ISTART ) ENDIF ! IIP > 1 *IDENT POISC230 */ For some reason this message passing no longer works under MPI. */ This change introduces same logic as in TIMSER which does work. *D POISSON.122,125 IF(RLOCDIVMAX.NE.DIVMAX)THEN ICOMAX(1)=-888 ICOMAX(2)=-888 ICOMAX(3)=-888 ENDIF ! This makes the assumption that the value of DIVMAX is ! unique to one PE - probably mostly true. CALL GC_IMAX(3,NPROC,INFO,ICOMAX) *IDENT DUMPC230 */ Certain diagnostics flags are not saved to the model dump *D DUMP.71 & ,IEXDG,IDGTHETA,IDGTHV *IDENT RETRC230 *D RETRIEVE.214 &,IEXDG,IDGTHETA,IDGTHV *IDENT NEWRC230 *D NEWRETR.228 &,IEXDG,IDGTHETA,IDGTHV *IDENT NEWDC230 *D NEWDUMP.78 & ,IEXDG,IDGTHETA,IDGTHV *IDENT MPHYC230 */ Add precipitation and column-classified diagnostics for warm rain *I MPHYS.64 ! Initialise the column classified and conditionally sampled ! diagnostics if IDGDEEP=1 IF ( idgdeep .NE. 0 ) CALL PRECDIAGSET( i, q ) *I MPHYS.87 ! Update the PRECIP diagnostic for warm rain IF ( idgdeep .NE. 0 ) CALL PRECDIAG( qdownflux, i, iqr, kkp ) *B MPHYS.94 !--------------------------------------- ! Reset RAINSINK and RAINSOURCE to zero !--------------------------------------- CALL SET0(RAINSINK,LENP) CALL SET0(RAINSOURCE,LENP) *IDENT MMOVC230 *D MMOVIES.123 slicemov(imov)=0.5*(v(j,k,i)+v(j-1,k,i))+vgal *IDENT FORIC230 *D FORINTER.97 & ( SLS_READ(KTMFOR,N)-SLS_READ(KTMFOR-1,N) ) */ This section is from Adrian Lock to correct an error in the */ calculation of Moist Ri no. *IDENT NEWMRI *// *// Mods to change LEM2.3 MOISTRI2 to work in terms of potential *// temperature throughout so that in cloud-free conditions neutrality *// (ie. RI=0) is now zero theta_v gradient. *// *// Simplest to strip it all out and start again! *// *D MOISTRI2.33,135 REAL EPS ! volume fraction which is "virtually exchanged" & ! (mixed) upper and lowere level & , REPSH ! 1/(2*eps) & , SCLTMP !Factor for Ri calculation & , THCONAP1 ! THCONA at z level K+1 & , QT_K, QT_KP1! QT (= QV+QL) at levels K and K+1 & , THLPR_K, THLPR_KP1 ! theta_l [=theta-(l/c_p)q_l/EXNER] at levels ! K and K+1 calculated at constant pressure & ! Before mixing: & , QLLI,QLUI !Liquid water at lower and upper levels & , THVLI,THVUI !Non Q_T part of "virt. pot. temperature" at & !lower/upper levels used to calculate buoyancy & !changes through mixing at lower and upper levels & !respectively [=Pot.Temp.-THREF*RATIO_MOL_WTS*QL] & ! AFTER mixing: & , QLLN,QLUN ! As QLLI,QLUI but after mixing a fraction EPS & , THVLN,THVUN ! AS THVLI,THVUI but after mixing a fraction EPS & , THLPR_UN,THLPR_LN ! Perturbation liq.water pot. temp. at lower/upper ! levels after mixing a fraction EPS & , QTUN,QTLN ! tot.water mix.ratio at lower/upper levels ! after mixing a fraction EPS EPS=0.01 REPSH=0.5/EPS *IF DEF,PORT ! ! FORTRAN EQUIVALENT OF ONEOVER_V CALL BELOW DO K=2,KKP-1 DO J=0,JJP+1 TMPINV(J,K) = 1./SSQ(J,K) ENDDO ENDDO *ELSE CALL ONEOVER_V((JJP+2)*(KKP-2),SSQ(0,2),TMPINV(0,2)) *ENDIF DO K=2,KKP-1 SCLTMP=RDZN(K+1)*BUOYCO(K) IF(IANELP.EQ.1)THEN THCONA=RATIO_MOL_WTS*THREF(K) THCONB=0.5*(RATIO_MOL_WTS-1.0)*(THREF(K)+THREF(K+1)) THCONAP1=RATIO_MOL_WTS*THREF(K+1) ELSE THCONAP1=THCONA ENDIF !DIR$ NOSPLIT DO J=JMINP,JMAXP QT_K = Q(J,K)+QL(J,K) QT_KP1 = Q(J,K+1)+QL(J,K+1) THLPR_K = TH(J,K) - RLVAP_ON_CP*QL(J,K)*PREFRCP(K) ! calculate theta_l at constant (level K) pressure THLPR_KP1 = TH(J,K+1) - RLVAP_ON_CP*QL(J,K+1)*PREFRCP(K) ! !.......CALCULATE QL AND BUOYANCY OF LEVELS K AND K+1 (LOWER AND UPPER) !.......QL CALC AT CONSTANT PRESSURE QLLI = MAX(0.0,( QT_K - ( QSAT(K) + & DQSATDT(K)*(RPREFRCP(K)*THLPR_K-TSTARPR(K)) ) & )*QSATFAC(K) ) QLUI = MAX(0.,( QT_KP1 - ( QSAT(K) + DQSATDT(K)* & (RPREFRCP(K)*(THLPR_KP1+THREF(K+1))-TSTARPR(K)-TREF(K)) ) & )*QSATFAC(K) ) THVLI = THLPR_K+THREF(K) + & ( RLVAP_ON_CP*PREFRCP(K)-THCONA ) * QLLI THVUI = THLPR_KP1+THREF(K+1) + & ( RLVAP_ON_CP*PREFRCP(K)-THCONAP1 ) * QLUI ! !.......MIX CONSERVATIVE QUANTITIES AT CONSTANT PRESSURE !.......NOTE THAT FOR COMPUTATIONAL CONVENIENCE THE RELEVANT LEVEL'S !.......THREF HAS BEEN OMITTED FROM ALL THV'S, THLPR_LN AND THLPR_UN AS !.......ONLY THE DIFFERENCE BETWEEN THV'S AT THE SAME LEVEL IS REQUIRED. THLPR_UN = THLPR_KP1 - EPS*( THLPR_KP1 - THLPR_K + DTHREF(K) ) THLPR_LN = THLPR_K + EPS*( THLPR_KP1 - THLPR_K + DTHREF(K) ) QTUN = QT_KP1 - EPS*( QT_KP1 - QT_K ) QTLN = QT_K + EPS*( QT_KP1 - QT_K ) !.......CALCULATE LIQUID WATER AND BUOYANCY OF MIXED AIR QLLN = MAX(0.0,( QTLN - ( QSAT(K) + & DQSATDT(K)*(RPREFRCP(K)*THLPR_LN-TSTARPR(K)) ) & )*QSATFAC(K) ) QLUN = MAX(0.,( QTUN - ( QSAT(K) + DQSATDT(K)* & (RPREFRCP(K)*(THLPR_UN+THREF(K+1))-TSTARPR(K)-TREF(K)) ) & )*QSATFAC(K) ) THVLN = THLPR_LN+THREF(K) + & ( RLVAP_ON_CP*PREFRCP(K)-THCONA ) * QLLN THVUN = THLPR_UN+THREF(K+1) + & ( RLVAP_ON_CP*PREFRCP(K)-THCONAP1 ) * QLUN ! !.....CHANGE IN BUOYANCY DUE TO FRAC MIXING IS USED TO DETERMINE RI RI(J,K) = SCLTMP*TMPINV(J,K)* ( THCONB*(QT_KP1-QT_K) & +((THVLN-THVLI)-(THVUN-THVUI))*REPSH ) ENDDO ENDDO */ Glen Richardson passed the following on to me. */ I havent tested them but they apparently work */ */ For the Kessler scheme (IRAINP=1) */ *I MPHYS.107 RAINSINK(J,K)=0. *I MPHYS.110 RAINSINK(J,K)=0. *I MPHYS.113 RAINSOURCE(J,K)=0. */ */For the Lee scheme (IRAINP=2) */ *I MPHYS.135 RAINSOURCE(J,K)=0. *I MPHYS.138 RAINSINK(J,K)=0. */ */ For the Swann scheme(IRAINP=3) */ *I MPHYS.165 RAINSINK(J,K)=0. *I MPHYS.171 RAINSOURCE(J,K)=0.