SUBROUTINE HEAT_CAP (NPNTS,SOIL_PTS,SOIL_INDEX,B,DZ,HCAP, 2,2HEATCP5A.23
& SATHH,SMCL,STHF,TSOIL,V_SAT,HCAPS HEATCP5A.24
! LOGICAL LTIMER HEATCP5A.25
+,LTIMER HEATCP5A.26
+) HEATCP5A.27
HEATCP5A.28
IMPLICIT NONE HEATCP5A.29
! HEATCP5A.30
! Description: HEATCP5A.31
! Calculates the effective heat capacity of a soil layer HEATCP5A.32
! including the effects of ice, water and phase changes. HEATCP5A.33
! (Cox, 6/95) HEATCP5A.34
! HEATCP5A.35
! Documentation : UM Documentation Paper 25 HEATCP5A.36
! HEATCP5A.37
! Current Code Owner : David Gregory HEATCP5A.38
! HEATCP5A.39
! History: HEATCP5A.40
! Version Date Comment HEATCP5A.41
! ------- ---- ------- HEATCP5A.42
! 4.1 New deck. Peter Cox HEATCP5A.43
! 4.5 26/5/98 Correction to stop overflow when SMCL very small. ACB2F405.1
! C.Bunton ACB2F405.2
! HEATCP5A.44
! Adds the vector shared, private directives. HEATCP5A.45
! Code Description: HEATCP5A.46
! Language: FORTRAN 77 + common extensions. HEATCP5A.47
! HEATCP5A.48
! System component covered: P25 HEATCP5A.49
! System Task: P25 HEATCP5A.50
! HEATCP5A.51
HEATCP5A.52
! Global variables: HEATCP5A.53
*CALL C_SOILH 
HEATCP5A.54
HEATCP5A.55
! Subroutine arguments HEATCP5A.56
! Scalar arguments with intent(IN) : HEATCP5A.57
INTEGER HEATCP5A.58
& NPNTS ! IN Number of gridpoints. HEATCP5A.59
&,SOIL_PTS ! IN Number of soil points. HEATCP5A.60
HEATCP5A.61
! Array arguments with intent(IN) : HEATCP5A.62
INTEGER HEATCP5A.63
& SOIL_INDEX(NPNTS) ! IN Array of soil points. HEATCP5A.64
HEATCP5A.65
REAL HEATCP5A.66
& B(NPNTS) ! IN Clapp-Hornberger exponent. HEATCP5A.67
&,DZ ! IN Thickness of the soil layer (m). HEATCP5A.68
&,HCAP(NPNTS) ! IN Soil heat capacity (J/K/m3). HEATCP5A.69
&,SATHH(NPNTS) ! IN Saturated soil water pressure (m). HEATCP5A.70
&,SMCL(NPNTS) ! IN Soil moisture content of the HEATCP5A.71
! ! layer (kg/m2). HEATCP5A.72
&,STHF(NPNTS) ! IN Frozen soil moisture content of HEATCP5A.73
! ! the layer as a fraction of HEATCP5A.74
! ! saturation. HEATCP5A.75
&,TSOIL(NPNTS) ! IN Layer temperature (K). HEATCP5A.76
&,V_SAT(NPNTS) ! IN Volumetric soil moisture HEATCP5A.77
! ! concentration at saturation HEATCP5A.78
! ! (m3 H2O/m3 soil). HEATCP5A.79
! HEATCP5A.80
LOGICAL LTIMER ! Logical switch for TIMER diags HEATCP5A.81
HEATCP5A.82
! Array arguments with intent(OUT) : HEATCP5A.83
REAL HEATCP5A.84
& HCAPS(NPNTS) ! OUT The total volumetric heat capacity HEATCP5A.85
! ! (soil+water) of the layer (J/m3/K). HEATCP5A.86
HEATCP5A.87
! Local scalars: HEATCP5A.88
INTEGER HEATCP5A.89
& I,J ! WORK Loop counter. HEATCP5A.90
HEATCP5A.91
! Local arrays: HEATCP5A.92
REAL HEATCP5A.93
& DTHU(NPNTS) ! WORK Rate of change of volumetric unfrozen HEATCP5A.94
! ! soil moisture concentration with HEATCP5A.95
! ! temperature (m3 liquid H2O/m3 soil/K) HEATCP5A.96
&,SMCLU(NPNTS) ! WORK Unfrozen moisture content of each HEATCP5A.97
! ! soil layer (kg/m2). HEATCP5A.98
&,SMCLSAT(NPNTS) ! WORK The saturation moisture content of HEATCP5A.99
! ! each layer (kg/m2). HEATCP5A.100
&,SMCLF(NPNTS) ! WORK Frozen moisture content of each HEATCP5A.101
! ! soil layer (kg/m2). HEATCP5A.102
&,TMAX(NPNTS) ! WORK Temperature above which all water is HEATCP5A.103
! ! unfrozen (Celsius) HEATCP5A.104
&,TSL(NPNTS) ! WORK Soil layer temperature (Celsius) HEATCP5A.105
HEATCP5A.106
*CALL C_DENSTY HEATCP5A.107
*CALL C_LHEAT HEATCP5A.108
*CALL C_PERMA HEATCP5A.109
*CALL C_0_DG_C HEATCP5A.110
HEATCP5A.111
IF (LTIMER) THEN HEATCP5A.112
CALL TIMER('HEATCAP ',3) HEATCP5A.113
ENDIF HEATCP5A.114
!---------------------------------------------------------------------- HEATCP5A.115
! Initialise all points HEATCP5A.116
!---------------------------------------------------------------------- HEATCP5A.117
DO I=1,NPNTS HEATCP5A.118
IF (V_SAT(I).GT.0.0) THEN ! Soil points HEATCP5A.119
HCAPS(I)=HCAP(I) HEATCP5A.120
ELSE ! Ice points HEATCP5A.121
HCAPS(I)=SNOW_HCAP HEATCP5A.122
ENDIF HEATCP5A.123
ENDDO HEATCP5A.124
HEATCP5A.125
!----------------------------------------------------------------------- HEATCP5A.126
! Calculate the temperature in Celsius and diagnose the saturation, HEATCP5A.127
! frozen and unfrozen soil moisture contents of the layer HEATCP5A.128
!----------------------------------------------------------------------- HEATCP5A.129
CDIR$ IVDEP HEATCP5A.130
! Fujitsu vectorization directive GRB0F405.341
!OCL NOVREC GRB0F405.342
DO J=1,SOIL_PTS HEATCP5A.131
I=SOIL_INDEX(J) HEATCP5A.132
TSL(I)=TSOIL(I)-ZERODEGC HEATCP5A.133
SMCLSAT(I)=RHO_WATER*DZ*V_SAT(I) HEATCP5A.134
SMCLF(I)=SMCLSAT(I)*STHF(I) HEATCP5A.135
SMCLU(I)=SMCL(I)-SMCLF(I) HEATCP5A.136
ENDDO HEATCP5A.137
HEATCP5A.138
HEATCP5A.139
CDIR$ IVDEP HEATCP5A.140
! Fujitsu vectorization directive GRB0F405.343
!OCL NOVREC GRB0F405.344
DO J=1,SOIL_PTS HEATCP5A.141
I=SOIL_INDEX(J) HEATCP5A.142
!----------------------------------------------------------------------- HEATCP5A.143
! Calculate TMAX, the temperature above which all soil water is HEATCP5A.144
! unfrozen. Check that (SMCLSAT/SMCL)**B will not overflow when SMCL is ACB2F405.3
! very small. The function EPSILON gives the number of type (real) of ACB2F405.4
! SMCL that is negligeable compared to 1. ACB2F405.5
!----------------------------------------------------------------------- ACB2F405.6
*IF DEF,SCMA,AND,-DEF,T3E ACB2F405.7
IF (SMCL(I) .GT. SMCLSAT(I)*1E-20) THEN ACB2F405.8
*ELSE ACB2F405.9
IF (SMCL(I) .GT. SMCLSAT(I)*(EPSILON(SMCL(I)) ) )THEN ACB2F405.10
*ENDIF ACB2F405.11
ACB2F405.12
TMAX(I)=-SATHH(I)/DPSIDT HEATCP5A.148
& *(SMCLSAT(I)/SMCL(I))**(B(I)) HEATCP5A.149
ELSE HEATCP5A.150
TMAX(I)=-ZERODEGC HEATCP5A.151
ENDIF HEATCP5A.152
HEATCP5A.153
!----------------------------------------------------------------------- HEATCP5A.154
! Calculate the rate of change of volumetric unfrozen soil moisture HEATCP5A.155
! concentration with temperature HEATCP5A.156
!----------------------------------------------------------------------- HEATCP5A.157
IF (TSL(I).GE.TMAX(I)) THEN HEATCP5A.158
DTHU(I)=0.0 HEATCP5A.159
ELSE HEATCP5A.160
DTHU(I)=DPSIDT*SMCLSAT(I) HEATCP5A.161
& /(B(I)*SATHH(I)*RHO_WATER*DZ) HEATCP5A.162
& *(-DPSIDT*TSL(I)/SATHH(I))**(-1.0/B(I)-1.0) HEATCP5A.163
ENDIF HEATCP5A.164
HEATCP5A.165
!----------------------------------------------------------------------- HEATCP5A.166
! Calculate the effective heat capacity HEATCP5A.167
!----------------------------------------------------------------------- HEATCP5A.168
HCAPS(I)=HCAP(I)+(HCAPW-HCAPI)*SMCLU(I)/DZ HEATCP5A.169
& +HCAPI*SMCL(I)/DZ HEATCP5A.170
& +RHO_WATER*DTHU(I)*((HCAPW-HCAPI)*TSL(I)+LF) HEATCP5A.171
ENDDO HEATCP5A.172
HEATCP5A.173
IF (LTIMER) THEN HEATCP5A.174
CALL TIMER('HEATCAP ',4) HEATCP5A.175
ENDIF HEATCP5A.176
HEATCP5A.177
RETURN HEATCP5A.178
END HEATCP5A.179
*ENDIF HEATCP5A.180