SUBROUTINE SFLX ( I ICE,DT,Z,NSOIL,SLDPTH, I LWDN,SOLDN,SFCPRS,PRCP,SFCTMP,TH2,Q2,SFCSPD,Q2SAT,DQSDT2, I VEGTYP,SOILTYP,SLOPETYP, I SHDFAC,PTU,TBOT,ALB,SNOALB, 2 CMC,T1,STC,SMC,SH2O,SNOWH,SNEQV,ALBEDO,CH,CM, O ETP,ETA,H,S,RUNOFF1,RUNOFF2,Q1,SNMAX, O SOILW,SOILM, SMCWLT,SMCDRY,SMCREF,SMCMAX ) C IMPLICIT NONE CC C ---------------------------------------------------------------------- CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC CC PURPOSE: SUB-DRIVER FOR "NOAH/OSU LSM" FAMILY OF PHYSICS SUBROUTINES CC FOR A SOIL/VEG/SNOWPACK LAND-SURFACE MODEL TO UPDATE SOIL CC MOISTURE, SOIL ICE, SOIL TEMPERATURE, SKIN TEMPERATURE, CC SNOWPACK WATER CONTENT, SNOWDEPTH, AND ALL TERMS CC OF THE SURFACE ENERGY BALANCE AND SURFACE WATER CC BALANCE (EXCLUDING INPUT ATMOSPHERIC FORCINGS OF CC DOWNWARD RADIATION AND PRECIP) CC CC VERSION 2.2 07 FEBRUARY 2001 CC CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC CC C ---------------------------------------------------------------------- C ------------ FROZEN GROUND VERSION ---------------------------- C ADDED STATES: SH2O(NSOIL) - UNFROZEN SOIL MOISTURE C SNOWH - SNOW DEPTH C C ---------------------------------------------------------------------- C C NOTE ON SNOW STATE VARIABLES: C SNOWH = actual physical snow depth in m C SNEQV = liquid water-equivalent snow depth in m C (time-dependent snow density is obtained from SNEQV/SNOWH) C C NOTE ON ALBEDO FRACTIONS: C Input: C ALB = BASELINE SNOW-FREE ALBEDO, FOR JULIAN DAY OF YEAR C (USUALLY FROM TEMPORAL INTERPOLATION OF MONTHLY MEAN VALUES) C (CALLING PROG MAY OR MAY NOT INCLUDE DIURNAL SUN ANGLE EFFECT) C SNOALB = UPPER BOUND ON MAXIMUM ALBEDO OVER DEEP SNOW C (E.G. FROM ROBINSON AND KUKLA, 1985, J. CLIM. & APPL. METEOR.) C Output: C ALBEDO = COMPUTED ALBEDO WITH SNOWCOVER EFFECTS C (COMPUTED USING ALB, SNOALB, SNEQV, AND SHDFAC->green veg frac) C C ARGUMENT LIST IN THE CALL TO SFLX C C ---------------------------------------------------------------------- C 1. CALLING STATEMENT C C SUBROUTINE SFLX C I (ICE,DT,Z,NSOIL,SLDPTH, C I LWDN,SOLDN,SFCPRS,PRCP,SFCTMP,TH2,Q2,Q2SAT,DQSDT2, C I VEGTYP,SOILTYP,SLOPETYP, C I SHDFAC,PTU,TBOT,ALB,SNOALB, C I SFCSPD, C 2 CMC,T1,STC,SMC,SH2O,SNOWH,SNEQV,CH,CM, C O ETP,ETA,H,S,RUNOFF1,RUNOFF2,Q1,SNMAX,ALBEDO, C O SOILW,SOILM,SMCWLT,SMCDRY,SMCREF,SMCMAX) C C 2. INPUT (denoted by "I" in column six of argument list at top of routine) C ### GENERAL PARAMETERS ### C C ICE: SEA-ICE FLAG (=1: SEA-ICE, =0: LAND) C DT: TIMESTEP (SEC) C (DT SHOULD NOT EXCEED 3600 SECS, RECOMMEND 1800 SECS OR LESS) C Z: HEIGHT (M) ABOVE GROUND OF ATMOSPHERIC FORCING VARIABLES C NSOIL: NUMBER OF SOIL LAYERS C (at least 2, and not greater than parameter NSOLD set below) C SLDPTH: THE THICKNESS OF EACH SOIL LAYER (M) C C ### ATMOSPHERIC VARIABLES ### C C LWDN: LW DOWNWARD RADIATION (W M-2; POSITIVE, not net longwave) C SOLDN: SOLAR DOWNWARD RADIATION (W M-2; POSITIVE, not net shortwave) C SFCPRS: PRESSURE AT HEIGHT Z ABOVE GROUND (PASCALS) C PRCP: PRECIP RATE (KG M-2 S-1) (note, this is a rate) C SFCTMP: AIR TEMPERATURE (K) AT HEIGHT Z ABOVE GROUND C TH2: AIR POTENTIAL TEMPERATURE (K) AT HEIGHT Z ABOVE GROUND C Q2: MIXING RATIO AT HEIGHT Z ABOVE GROUND (KG KG-1) C SFCSPD: WIND SPEED (M S-1) AT HEIGHT Z ABOVE GROUND C Q2SAT: SAT MIXING RATIO AT HEIGHT Z ABOVE GROUND (KG KG-1) C DQSDT2: SLOPE OF SAT SPECIFIC HUMIDITY CURVE AT T=SFCTMP (KG KG-1 K-1) C C ### CANOPY/SOIL CHARACTERISTICS ### C C VEGTYP: VEGETATION TYPE (INTEGER INDEX) C SOILTYP: SOIL TYPE (INTEGER INDEX) C SLOPETYP: CLASS OF SFC SLOPE (INTEGER INDEX) C SHDFAC: AREAL FRACTIONAL COVERAGE OF GREEN VEGETATION (range 0.0-1.0) C PTU: PHOTO THERMAL UNIT (PLANT PHENOLOGY FOR ANNUALS/CROPS) C (not yet used, but passed to REDPRM for future use in veg parms) C TBOT: BOTTOM SOIL TEMPERATURE (LOCAL YEARLY-MEAN SFC AIR TEMPERATURE) C ALB: BACKROUND SNOW-FREE SURFACE ALBEDO (FRACTION) C SNOALB: ALBEDO UPPER BOUND OVER DEEP SNOW (FRACTION) C C 3. STATE VARIABLES: BOTH INPUT AND OUTPUT C (NOTE: OUTPUT USUALLY MODIFIED FROM INPUT BY PHYSICS) C C (denoted by "2" in column six of argument list at top of routine) C C !!! ########### STATE VARIABLES ############## !!! C C CMC: CANOPY MOISTURE CONTENT (M) C T1: GROUND/CANOPY/SNOWPACK) EFFECTIVE SKIN TEMPERATURE (K) C C STC(NSOIL): SOIL TEMP (K) C SMC(NSOIL): TOTAL SOIL MOISTURE CONTENT (VOLUMETRIC FRACTION) C SH2O(NSOIL): UNFROZEN SOIL MOISTURE CONTENT (VOLUMETRIC FRACTION) C NOTE: FROZEN SOIL MOISTURE = SMC - SH2O C C SNOWH: SNOW DEPTH (M) C SNEQV: WATER-EQUIVALENT SNOW DEPTH (M) C NOTE: SNOW DENSITY = SNEQV/SNOWH C ALBEDO: SURFACE ALBEDO INCLUDING SNOW EFFECT (UNITLESS FRACTION) C CH: SFC EXCH COEF FOR HEAT AND MOISTURE (M S-1) C CM: SFC EXCH COEF FOR MOMENTUM (M S-1) C NOTE: CH AND CM ARE TECHNICALLY CONDUCTANCES SINCE THEY C HAVE BEEN MULTIPLIED BY THE WIND SPEED. C C 4. OUTPUT (denoted by "O" in column six of argument list at top of routine) C C NOTE-- SIGN CONVENTION OF SFC ENERGY FLUXES BELOW IS: NEGATIVE IF C SINK OF ENERGY TO SURFACE C C ETP: POTENTIAL EVAPORATION (W M-2) C ETA: ACTUAL LATENT HEAT FLUX (W M-2: NEGATIVE, IF UP FROM SURFACE) C H: SENSIBLE HEAT FLUX (W M-2: NEGATIVE, IF UPWARD FROM SURFACE) C S: SOIL HEAT FLUX (W M-2: NEGATIVE, IF DOWNWARD FROM SURFACE) C RUNOFF1: SURFACE RUNOFF (M S-1), NOT INFILTRATING THE SURFACE C RUNOFF2: SUBSURFACE RUNOFF (M S-1), DRAINAGE OUT BOTTOM OF LAST SOIL LYR C Q1: EFFECTIVE MIXING RATIO AT GRND SFC (KG KG-1) C (NOTE: Q1 IS NUMERICAL EXPENDIENCY FOR EXPRESSING ETA C EQUIVALENTLY IN A BULK AERODYNAMIC FORM) C SNMAX: SNOW MELT (M) (WATER EQUIVALENT) C SOILW: AVAILABLE SOIL MOISTURE IN ROOT ZONE (UNITLESS FRACTION BETWEEN C SOIL SATURATION AND WILTING POINT) C SOILM: TOTAL SOIL COLUMN MOISTURE CONTENT (M) (FROZEN + UNFROZEN) C C FOR DIAGNOSTIC PURPOSES, RETURN SOME PRIMARY PARAMETERS NEXT C (SET IN ROUTINE REDPRM) C C SMCWLT: WILTING POINT (VOLUMETRIC) C SMCDRY: DRY SOIL MOISTURE THRESHOLD WHERE DIRECT EVAP FRM TOP LYR ENDS C SMCREF: SOIL MOISTURE THRESHOLD WHERE TRANSPIRATION BEGINS TO STRESS C SMCMAX: POROSITY, I.E. SATURATED VALUE OF SOIL MOISTURE INTEGER NSOLD PARAMETER (NSOLD = 20) C LOGICAL SNOWNG LOGICAL FRZGRA LOGICAL SATURATED C INTEGER K INTEGER KZ INTEGER ICE INTEGER NSOIL,VEGTYP,SOILTYP,NROOT INTEGER SLOPETYP C REAL ALBEDO REAL ALB REAL B REAL BETA REAL CFACTR C..................CH IS SFC EXCHANGE COEF FOR HEAT/MOIST C..................CM IS SFC MOMENTUM DRAG (NOT NEEDED IN SFLX) REAL CH REAL CM C REAL CMC REAL CMCMAX REAL CP REAL CSNOW REAL CSOIL REAL CZIL REAL DEW REAL DF1 REAL DKSAT REAL DT REAL DWSAT REAL DQSDT2 REAL DSOIL REAL DTOT REAL DRIP REAL EC REAL EDIR REAL ETT REAL EXPSNO REAL EXPSOI REAL EPSCA REAL ETA REAL ETP REAL EDIR1 REAL EC1 REAL ETT1 REAL F REAL F1 REAL FLX1 REAL FLX2 REAL FLX3 REAL FXEXP REAL FRZX REAL H REAL HS REAL KDT REAL LWDN REAL LVH2O REAL PC REAL PRCP REAL PTU REAL PRCP1 REAL PSISAT REAL Q1 REAL Q2 REAL Q2SAT REAL QUARTZ REAL R REAL RCH REAL REFKDT REAL RR REAL RTDIS (NSOLD) REAL RUNOFF1 REAL RUNOFF2 REAL RGL REAL RUNOF REAL RIB REAL RUNOFF3 REAL RSMAX REAL RC REAL RCMIN REAL RSNOW REAL SNDENS REAL SNCOND REAL S REAL SBETA REAL SFCPRS REAL SFCSPD REAL SFCTMP REAL SHDFAC REAL SH2O(NSOIL) REAL SLDPTH(NSOIL) REAL SMCDRY REAL SMCMAX REAL SMCREF REAL SMCWLT REAL SMC(NSOIL) REAL SNEQV REAL SNOWH REAL SNOFAC REAL SN_NEW REAL SLOPE REAL SNUP REAL SALP REAL SNOALB REAL STC(NSOIL) REAL SOLDN REAL SNMAX REAL SOILM REAL SOILW REAL SOILWM REAL SOILWW REAL T1 REAL T1V REAL T24 REAL T2V REAL TBOT REAL TH2 REAL TH2V REAL TOPT REAL TFREEZ REAL XLAI REAL Z REAL ZBOT REAL Z0 REAL ZSOIL(NSOLD) C PARAMETER ( TFREEZ = 273.15 ) PARAMETER ( LVH2O = 2.501000E+6 ) PARAMETER ( R = 287.04 ) PARAMETER ( CP = 1004.5 ) C C COMMON BLK "RITE" CARRIES DIAGNOSTIC QUANTITIES FOR PRINTOUT, C BUT IS NOT INVOLVED IN MODEL PHYSICS AND IS NOT PRESENT IN C PARENT MODEL THAT CALLS SFLX C COMMON/RITE/ BETA,DRIP,EC,EDIR,ETT,FLX1,FLX2,FLX3,RUNOF, & DEW,RIB,RUNOFF3 C INITIALIZATION RUNOFF1 = 0.0 RUNOFF2 = 0.0 RUNOFF3 = 0.0 SNMAX = 0.0 C C THE VARIABLE "ICE" IS A FLAG DENOTING SEA-ICE CASE IF(ICE .EQ. 1) THEN C SEA-ICE LAYERS ARE EQUAL THICKNESS AND SUM TO 3 METERS DO KZ = 1, NSOIL ZSOIL(KZ)=-3.*FLOAT(KZ)/FLOAT(NSOIL) END DO ELSE C CALCULATE DEPTH (NEGATIVE) BELOW GROUND FROM TOP SKIN SFC TO C BOTTOM OF EACH SOIL LAYER. C NOTE:!!! SIGN OF ZSOIL IS NEGATIVE (DENOTING BELOW GROUND) ZSOIL(1)=-SLDPTH(1) DO KZ = 2, NSOIL ZSOIL(KZ)=-SLDPTH(KZ)+ZSOIL(KZ-1) END DO ENDIF C ---------------------------------------------------------------------- CC CC NEXT IS CRUCIAL CALL TO SET THE LAND-SURFACE PARAMETERS, CC INCLUDING SOIL-TYPE AND VEG-TYPE DEPENDENT PARAMETERS. CC CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC CC CALL REDPRM(VEGTYP,SOILTYP, SLOPETYP, + CFACTR, CMCMAX, RSMAX, TOPT, REFKDT, KDT, SBETA, O SHDFAC, RCMIN, RGL, HS, ZBOT, FRZX, PSISAT, SLOPE, + SNUP, SALP, B, DKSAT, DWSAT, SMCMAX, SMCWLT, SMCREF, O SMCDRY, F1, QUARTZ, FXEXP, RTDIS, SLDPTH, ZSOIL, + NROOT, NSOIL, Z0, CZIL, XLAI, CSOIL, PTU) C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC CC CC NEXT CALL ROUTINE SFCDIF TO CALCULATE CC THE SFC EXCHANGE COEF (CH) FOR HEAT AND MOISTURE CC CC NOTE NOTE !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! CC CC COMMENT OUT CALL SFCDIF, IF SFCDIF ALREADY CALLED CC IN CALLING PROGRAM (SUCH AS IN COUPLED ATMOSPHERIC MODEL) CC CC NOTE !! DO NOT CALL SFCDIF UNTIL AFTER ABOVE CALL TO REDPRM, CC IN CASE ALTERNATIVE VALUES OF ROUGHNESS LENGTH (Z0) AND CC ZILINTINKEVICH COEF (CZIL) ARE SET THERE VIA NAMELIST I/O CC CC NOTE !! ROUTINE SFCDIF RETURNS A CH THAT REPRESENTS THE WIND SPD CC TIMES THE "ORIGINAL" NONDIMENSIONAL "Ch" TYPICAL IN LITERATURE. CC HENCE THE CH RETURNED FROM SFCDIF HAS UNITS OF M/S. CC THE IMPORTANT COMPANION COEFFICIENT OF CH, CARRIED HERE AS "RCH", CC IS THE CH FROM SFCDIF TIMES AIR DENSITY AND PARAMETER "CP". CC "RCH" IS COMPUTED IN "CALL PENMAN". RCH RATHER THAN CH IS THE C COEFF USUALLY INVOKED LATER IN EQNS. CC CC NOTE !! SFCDIF ALSO RETURNS THE SURFACE EXCHANGE COEFFICIENT FOR C MOMENTUM, CM, ALSO KNOWN AS THE SURFACE DRAGE COEFFICIENT, C BUT CM IS NOT USED HERE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C ---------------------------------------------------------------------- CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C CALC VIRTUAL TEMPS AND VIRTUAL POTENTIAL TEMPS NEEDED BY C SUBROUTINES SFCDIF AND PENMAN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC T2V = SFCTMP * (1.0 + 0.61 * Q2 ) c comment out below 2 lines if CALL SFCDIF is commented out, i.e. in c the coupled model c T1V = T1 * (1.0 + 0.61 * Q2 ) c TH2V = TH2 * (1.0 + 0.61 * Q2 ) C C CALL SFCDIF ( Z, Z0, T1V, TH2V, SFCSPD, CZIL, CM, CH ) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C INITIALIZE MISC VARIABLES. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC SNOWNG = .FALSE. FRZGRA = .FALSE. C IF SEA-ICE CASE, ASSIGN DEFAULT WATER-EQUIV SNOW ON TOP IF(ICE .EQ. 1) THEN SNEQV = 0.01 SNOWH = 0.05 ENDIF C C IF INPUT SNOWPACK IS NONZERO, THEN COMPUTE SNOW DENSITY "SNDENS" C AND SNOW THERMAL CONDUCTIVITY "SNCOND" C (NOTE THAT CSNOW IS A FUNCTION SUBROUTINE) C IF(SNEQV .EQ. 0.0) THEN SNDENS = 0.0 SNOWH = 0.0 SNCOND = 1.0 ELSE SNDENS=SNEQV/SNOWH SNCOND = CSNOW (SNDENS) ENDIF C ---------------------------------------------------------------------- CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C DETERMINE IF IT'S PRECIPITATING AND WHAT KIND OF PRECIP IT IS. C IF IT'S PRCPING AND THE AIR TEMP IS COLDER THAN 0 C, IT'S SNOWING! C IF IT'S PRCPING AND THE AIR TEMP IS WARMER THAN 0 C, BUT THE GRND C TEMP IS COLDER THAN 0 C, FREEZING RAIN IS PRESUMED TO BE FALLING. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IF ( PRCP .GT. 0.0 ) THEN IF ( SFCTMP .LE. TFREEZ ) THEN SNOWNG = .TRUE. ELSE IF ( T1 .LE. TFREEZ ) FRZGRA = .TRUE. ENDIF ENDIF C ---------------------------------------------------------------------- C If either prcp flag is set, determine new snowfall (converting prcp C rate from kg m-2 s-1 to a liquid equiv snow depth in meters) and add C it to the existing snowpack. C Note that since all precip is added to snowpack, no precip infiltrates C into the soil so that PRCP1 is set to zero. IF ( ( SNOWNG ) .OR. ( FRZGRA ) ) THEN SN_NEW = PRCP * DT * 0.001 SNEQV = SNEQV + SN_NEW PRCP1 = 0.0 C ---------------------------------------------------------------------- C Update snow density based on new snowfall, using old and new snow. CALL SNOW_NEW (SFCTMP,SN_NEW,SNOWH,SNDENS) C --- debug ------------------------------------------------------------ c SNDENS = 0.2 c SNOWH = SNEQV/SNDENS C --- debug ------------------------------------------------------------ C ---------------------------------------------------------------------- C Update snow thermal conductivity SNCOND = CSNOW (SNDENS) C ---------------------------------------------------------------------- ELSE C C PRECIP IS LIQUID (RAIN), HENCE SAVE IN THE PRECIP VARIABLE THAT C LATER CAN WHOLELY OR PARTIALLY INFILTRATE THE SOIL (ALONG WITH C ANY CANOPY "DRIP" ADDED TO THIS LATER) C PRCP1 = PRCP ENDIF C ---------------------------------------------------------------------- C Thermal conductivity for sea-ice case IF (ICE .EQ. 1) THEN DF1=2.2 ELSE C C NEXT CALCULATE THE SUBSURFACE HEAT FLUX, WHICH FIRST REQUIRES C CALCULATION OF THE THERMAL DIFFUSIVITY. TREATMENT OF THE C LATTER FOLLOWS THAT ON PAGES 148-149 FROM "HEAT TRANSFER IN C COLD CLIMATES", BY V. J. LUNARDINI (PUBLISHED IN 1981 C BY VAN NOSTRAND REINHOLD CO.) I.E. TREATMENT OF TWO CONTIGUOUS C "PLANE PARALLEL" MEDIUMS (NAMELY HERE THE FIRST SOIL LAYER C AND THE SNOWPACK LAYER, IF ANY). THIS DIFFUSIVITY TREATMENT C BEHAVES WELL FOR BOTH ZERO AND NONZERO SNOWPACK, INCLUDING THE C LIMIT OF VERY THIN SNOWPACK. THIS TREATMENT ALSO ELIMINATES C THE NEED TO IMPOSE AN ARBITRARY UPPER BOUND ON SUBSURFACE C HEAT FLUX WHEN THE SNOWPACK BECOMES EXTREMELY THIN. C C ---------------------------------------------------------------------- C FIRST CALCULATE THERMAL DIFFUSIVITY OF TOP SOIL LAYER, USING C BOTH THE FROZEN AND LIQUID SOIL MOISTURE, FOLLOWING THE C SOIL THERMAL DIFFUSIVITY FUNCTION OF PETERS-LIDARD ET AL. C (1998,JAS, VOL 55, 1209-1224), WHICH REQUIRES THE SPECIFYING C THE QUARTZ CONTENT OF THE GIVEN SOIL CLASS (SEE ROUTINE REDPRM) C CALL TDFCND ( DF1, SMC(1),QUARTZ,SMCMAX,SH2O(1) ) C ---------------------------------------------------------------------- C NEXT ADD SUBSURFACE HEAT FLUX REDUCTION EFFECT FROM THE C OVERLYING GREEN CANOPY, ADAPTED FROM SECTION 2.1.2 OF C PETERS-LIDARD ET AL. (1997, JGR, VOL 102(D4)) C DF1 = DF1 * EXP(SBETA*SHDFAC) ENDIF C ---------------------------------------------------------------------- C FINALLY "PLANE PARALLEL" SNOWPACK EFFECT FOLLOWING C V.J. LINARDINI REFERENCE CITED ABOVE. NOTE THAT DTOT IS C COMBINED DEPTH OF SNOWDEPTH AND THICKNESS OF FIRST SOIL LAYER C DSOIL = -(0.5 * ZSOIL(1)) IF (SNEQV .EQ. 0.) THEN S = DF1 * (T1 - STC(1) ) / DSOIL ELSE DTOT = SNOWH + DSOIL EXPSNO = SNOWH/DTOT EXPSOI = DSOIL/DTOT c 1. harmonic mean (series flow) c DF1 = (SNCOND*DF1)/(EXPSOI*SNCOND+EXPSNO*DF1) c 2. arithmetic mean (parallel flow) c DF1 = EXPSNO*SNCOND + EXPSOI*DF1 c 3. geometric mean (intermediate between c harmonic and arithmetic mean) DF1 = (SNCOND**EXPSNO)*(DF1**EXPSOI) C ---------------------------------------------------------------------- C CALCULATE SUBSURFACE HEAT FLUX, S, FROM FINAL THERMAL DIFFUSIVITY C OF SURFACE MEDIUMS, DF1 ABOVE, AND SKIN TEMPERATURE AND TOP C MID-LAYER SOIL TEMPERATURE S = DF1 * (T1 - STC(1) ) / DTOT ENDIF C ---------------------------------------------------------------------- C Update albedo, except over sea-ice IF (ICE .EQ. 0) THEN C ---------------------------------------------------------------------- C NEXT IS TIME-DEPENDENT SURFACE ALBEDO MODIFICATION DUE TO C TIME-DEPENDENT SNOWDEPTH STATE AND TIME-DEPENDENT CANOPY GREENNESS c IF ( (SNEQV .EQ. 0.0) .OR. (ALB .GE. SNOALB) ) THEN IF (SNEQV .EQ. 0.0) THEN ALBEDO = ALB ELSE C ---------------------------------------------------------------------- C SNUP IS VEG-CLASS DEPENDENT SNOWDEPTH THRESHHOLD (SET IN ROUTINE C REDPRM)WHERE MAX SNOW ALBEDO EFFECT IS FIRST ATTAINED IF (SNEQV .LT. SNUP) THEN RSNOW = SNEQV/SNUP SNOFAC = 1. - ( EXP(-SALP*RSNOW) - RSNOW*EXP(-SALP)) ELSE SNOFAC = 1.0 ENDIF C ---------------------------------------------------------------------- C SNOALB IS ARGUMENT REPRESENTING MAXIMUM ALBEDO OVER DEEP SNOW, C AS PASSED INTO SFLX, AND ADAPTED FROM THE SATELLITE-BASED MAXIMUM C SNOW ALBEDO FIELDS PROVIDED BY D. ROBINSON AND G. KUKLA C (1985, JCAM, VOL 24, 402-411) ALBEDO = ALB + (1.0-SHDFAC)*SNOFAC*(SNOALB-ALB) IF (ALBEDO .GT. SNOALB) ALBEDO=SNOALB ENDIF ELSE C ---------------------------------------------------------------------- C albedo over sea-ice ALBEDO = 0.60 SNOFAC = 1.0 ENDIF C ---------------------------------------------------------------------- C CALCULATE TOTAL DOWNWARD RADIATION (SOLAR PLUS LONGWAVE) C NEEDED IN PENMAN EP SUBROUTINE THAT FOLLOWS F = SOLDN*(1.0-ALBEDO) + LWDN C ---------------------------------------------------------------------- CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C CALL PENMAN SUBROUTINE TO CALCULATE POTENTIAL EVAPORATION (ETP) C (AND OTHER PARTIAL PRODUCTS AND SUMS SAVE IN COMMON/RITE FOR C LATER CALCULATIONS) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC CALL PENMAN ( SFCTMP,SFCPRS,CH,T2V,TH2,PRCP,F,T24,S,Q2, & Q2SAT,ETP,RCH,EPSCA,RR,SNOWNG,FRZGRA,DQSDT2) C C following old constraint is disabled C.....IF(SATURATED) ETP = 0.0 C ---------------------------------------------------------------------- CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C CALL CANRES TO CALCULATE THE CANOPY RESISTANCE AND CONVERT IT C INTO PC IF MORE THAN TRACE AMOUNT OF CANOPY GREENNESS FRACTION CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC c IF(SHDFAC .GT. 1.E-6) THEN c make this threshold consistent with the one in SMFLX for TRANSP c and EC(anopy) IF(SHDFAC .GT. 0.) THEN C FROZEN GROUND EXTENSION: TOTAL SOIL WATER "SMC" WAS REPLACED C BY UNFROZEN SOIL WATER "SH2O" IN CALL TO CANRES BELOW C CALL CANRES(SOLDN,CH,SFCTMP,Q2,SFCPRS,SH2O,ZSOIL,NSOIL, & SMCWLT,SMCREF,RCMIN,RC,PC,NROOT,Q2SAT,DQSDT2, & TOPT,RSMAX,RGL,HS,XLAI) ENDIF C ---------------------------------------------------------------------- CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C NOW DECIDE MAJOR PATHWAY BRANCH TO TAKE DEPENDING ON WHETHER C SNOWPACK EXISTS OR NOT CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IF ( SNEQV .EQ. 0.0 ) THEN CALL NOPAC ( ETP, ETA, PRCP, SMC, SMCMAX, SMCWLT, & SMCREF,SMCDRY, CMC, CMCMAX, NSOIL, DT, SHDFAC, & SBETA,Q1,Q2,T1,SFCTMP,T24,TH2,F,F1,S,STC, & EPSCA, B, PC, RCH, RR, CFACTR, + SH2O, SLOPE, KDT, FRZX, PSISAT, ZSOIL, & DKSAT, DWSAT, TBOT, ZBOT, RUNOFF1,RUNOFF2, & RUNOFF3, EDIR1, EC1, ETT1,NROOT,ICE,RTDIS, & QUARTZ, FXEXP,CSOIL) ELSE CALL SNOPAC ( ETP,ETA,PRCP,PRCP1,SNOWNG,SMC,SMCMAX,SMCWLT, & SMCREF, SMCDRY, CMC, CMCMAX, NSOIL, DT, & SBETA,Q1,DF1, & Q2,T1,SFCTMP,T24,TH2,F,F1,S,STC,EPSCA,SFCPRS, c & B, PC, RCH, RR, CFACTR, SALP, SNEQV, & B, PC, RCH, RR, CFACTR, SNOFAC, SNEQV,SNDENS, + SNOWH, SH2O, SLOPE, KDT, FRZX, PSISAT, SNUP, & ZSOIL, DWSAT, DKSAT, TBOT, ZBOT, SHDFAC,RUNOFF1, & RUNOFF2,RUNOFF3,EDIR1,EC1,ETT1,NROOT,SNMAX,ICE, & RTDIS,QUARTZ, FXEXP,CSOIL) ENDIF C ---------------------------------------------------------------------- C PREPARE SENSIBLE HEAT (H) FOR RETURN TO PARENT MODEL CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC H = -(CH * CP * SFCPRS)/(R * T2V) * ( TH2 - T1 ) C ---------------------------------------------------------------------- C CONVERT UNITS AND/OR SIGN OF TOTAL EVAP (ETA), POTENTIAL EVAP (ETP), C SUBSURFACE HEAT FLUX (S), AND RUNOFFS FOR WHAT PARENT MODEL EXPECTS CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C CONVERT ETA FROM KG M-2 S-1 TO W M-2 C ETA = ETA*LVH2O ETP = ETP*LVH2O C CONVERT THE SIGN OF SOIL HEAT FLUX SO THAT: C S>0: WARM THE SURFACE (NIGHT TIME) C S<0: COOL THE SURFACE (DAY TIME) S=-1.0*S C C CONVERT RUNOFF3 (INTERNAL LAYER RUNOFF FROM SUPERSAT) FROM M TO M S-1 C AND ADD TO SUBSURFACE RUNOFF/DRAINAGE/BASEFLOW C RUNOFF3 = RUNOFF3/DT RUNOFF2 = RUNOFF2+RUNOFF3 C C TOTAL COLUMN SOIL MOISTURE IN METERS (SOILM) AND ROOT-ZONE C SOIL MOISTURE AVAILABILITY (FRACTION) RELATIVE TO POROSITY/SATURATION SOILM=-1.0*SMC(1)*ZSOIL(1) DO K = 2, NSOIL SOILM=SOILM+SMC(K)*(ZSOIL(K-1)-ZSOIL(K)) END DO SOILWM=-1.0*(SMCMAX-SMCWLT)*ZSOIL(1) SOILWW=-1.0*(SMC(1)-SMCWLT)*ZSOIL(1) DO K = 2, NROOT SOILWM=SOILWM+(SMCMAX-SMCWLT)*(ZSOIL(K-1)-ZSOIL(K)) SOILWW=SOILWW+(SMC(K)-SMCWLT)*(ZSOIL(K-1)-ZSOIL(K)) END DO SOILW=SOILWW/SOILWM C RETURN END