-
Notifications
You must be signed in to change notification settings - Fork 312
/
lnd_import_export_utils.F90
177 lines (135 loc) · 7.17 KB
/
lnd_import_export_utils.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
module lnd_import_export_utils
use shr_kind_mod , only : r8 => shr_kind_r8
use shr_infnan_mod , only : isnan => shr_infnan_isnan
use shr_sys_mod , only : shr_sys_abort
use clm_varctl , only : iulog
use decompmod , only : bounds_type
use atm2lndType , only : atm2lnd_type
use Wateratm2lndBulkType , only : wateratm2lndbulk_type
implicit none
private ! except
public :: derive_quantities
public :: check_for_errors
public :: check_for_nans
!=============================================================================
contains
!=============================================================================
!===========================================================================
subroutine derive_quantities( bounds, atm2lnd_inst, wateratm2lndbulk_inst, &
forc_rainc, forc_rainl, forc_snowc, forc_snowl )
!-------------------------------------------------------------------------
! Convert the input data from the mediator to the land model
!-------------------------------------------------------------------------
use clm_varcon, only: rair, o2_molar_const
use QSatMod, only: QSat
! input/output variabes
type(bounds_type), intent(in) :: bounds ! bounds
type(atm2lnd_type), intent(inout) :: atm2lnd_inst ! clm internal input data type
type(wateratm2lndbulk_type), intent(inout) :: wateratm2lndbulk_inst
real(r8), intent(in) :: forc_rainc(bounds%begg:bounds%endg) ! convective rain (mm/s)
real(r8), intent(in) :: forc_rainl(bounds%begg:bounds%endg) ! large scale rain (mm/s)
real(r8), intent(in) :: forc_snowc(bounds%begg:bounds%endg) ! convective snow (mm/s)
real(r8), intent(in) :: forc_snowl(bounds%begg:bounds%endg) ! large scale snow (mm/s)
! local variables
integer :: g ! indices
integer :: begg, endg ! bounds
real(r8) :: qsat_kg_kg ! saturation specific humidity (kg/kg)
real(r8) :: forc_t ! atmospheric temperature (Kelvin)
real(r8) :: forc_q ! atmospheric specific humidity (kg/kg)
real(r8) :: forc_pbot ! atmospheric pressure (Pa)
character(len=*), parameter :: subname='(cpl:utils:derive_quantities)'
!-------------------------------------------------------------------------
! Set bounds
begg = bounds%begg; endg=bounds%endg
!--------------------------
! Derived quantities
!--------------------------
do g = begg, endg
forc_t = atm2lnd_inst%forc_t_not_downscaled_grc(g)
forc_q = wateratm2lndbulk_inst%forc_q_not_downscaled_grc(g)
forc_pbot = atm2lnd_inst%forc_pbot_not_downscaled_grc(g)
atm2lnd_inst%forc_hgt_u_grc(g) = atm2lnd_inst%forc_hgt_grc(g) !observational height of wind [m]
atm2lnd_inst%forc_hgt_t_grc(g) = atm2lnd_inst%forc_hgt_grc(g) !observational height of temperature [m]
atm2lnd_inst%forc_hgt_q_grc(g) = atm2lnd_inst%forc_hgt_grc(g) !observational height of humidity [m]
atm2lnd_inst%forc_vp_grc(g) = forc_q * forc_pbot / (0.622_r8 + 0.378_r8 * forc_q)
atm2lnd_inst%forc_rho_not_downscaled_grc(g) = &
(forc_pbot - 0.378_r8 * atm2lnd_inst%forc_vp_grc(g)) / (rair * forc_t)
atm2lnd_inst%forc_po2_grc(g) = o2_molar_const * forc_pbot
atm2lnd_inst%forc_wind_grc(g) = sqrt(atm2lnd_inst%forc_u_grc(g)**2 + atm2lnd_inst%forc_v_grc(g)**2)
atm2lnd_inst%forc_solar_grc(g) = atm2lnd_inst%forc_solad_grc(g,1) + atm2lnd_inst%forc_solai_grc(g,1) + &
atm2lnd_inst%forc_solad_grc(g,2) + atm2lnd_inst%forc_solai_grc(g,2)
wateratm2lndbulk_inst%forc_rain_not_downscaled_grc(g) = forc_rainc(g) + forc_rainl(g)
wateratm2lndbulk_inst%forc_snow_not_downscaled_grc(g) = forc_snowc(g) + forc_snowl(g)
call QSat(forc_t, forc_pbot, qsat_kg_kg)
! modify specific humidity if precip occurs
if (1==2) then
if ((forc_rainc(g) + forc_rainl(g)) > 0._r8) then
forc_q = 0.95_r8 * qsat_kg_kg
!forc_q = qsat_kg_kg
wateratm2lndbulk_inst%forc_q_not_downscaled_grc(g) = forc_q
endif
endif
wateratm2lndbulk_inst%forc_rh_grc(g) = 100.0_r8*(forc_q / qsat_kg_kg)
end do
end subroutine derive_quantities
!===========================================================================
subroutine check_for_errors( bounds, atm2lnd_inst, wateratm2lndbulk_inst )
! input/output variabes
type(bounds_type), intent(in) :: bounds ! bounds
type(atm2lnd_type), intent(inout) :: atm2lnd_inst ! clm internal input data type
type(wateratm2lndbulk_type), intent(inout) :: wateratm2lndbulk_inst
! local variables
integer :: g ! indices
integer :: begg, endg ! bounds
character(len=*), parameter :: subname='(cpl:utils:check_for_errors)'
!-------------------------------------------------------------------------
! Set bounds
begg = bounds%begg; endg=bounds%endg
!--------------------------
! Error checks
!--------------------------
! Check that solar, specific-humidity, and LW downward aren't negative
do g = begg, endg
if ( atm2lnd_inst%forc_lwrad_not_downscaled_grc(g) <= 0.0_r8 ) then
call shr_sys_abort( subname//&
' ERROR: Longwave down sent from the atmosphere model is negative or zero' )
end if
if ( (atm2lnd_inst%forc_solad_grc(g,1) < 0.0_r8) .or. &
(atm2lnd_inst%forc_solad_grc(g,2) < 0.0_r8) .or. &
(atm2lnd_inst%forc_solai_grc(g,1) < 0.0_r8) .or. &
(atm2lnd_inst%forc_solai_grc(g,2) < 0.0_r8) ) then
call shr_sys_abort( subname//&
' ERROR: One of the solar fields (indirect/diffuse, vis or near-IR)'// &
' from the atmosphere model is negative or zero' )
end if
if ( wateratm2lndbulk_inst%forc_q_not_downscaled_grc(g) < 0.0_r8 )then
call shr_sys_abort( subname//&
' ERROR: Bottom layer specific humidty sent from the atmosphere model is less than zero' )
end if
end do
! Make sure relative humidity is properly bounded
! atm2lnd_inst%forc_rh_grc(g) = min( 100.0_r8, atm2lnd_inst%forc_rh_grc(g) )
! atm2lnd_inst%forc_rh_grc(g) = max( 0.0_r8, atm2lnd_inst%forc_rh_grc(g) )
end subroutine check_for_errors
!=============================================================================
subroutine check_for_nans(array, fname, begg)
! input/output variables
real(r8) , intent(in) :: array(:)
character(len=*) , intent(in) :: fname
integer , intent(in) :: begg
! local variables
integer :: i
!---------------------------------------------------------------------------
! Check if any input from mediator or output to mediator is NaN
if (any(isnan(array))) then
write(iulog,*) '# of NaNs = ', count(isnan(array))
write(iulog,*) 'Which are NaNs = ', isnan(array)
do i = 1, size(array)
if (isnan(array(i))) then
write(iulog,*) "NaN found in field ", trim(fname), ' at gridcell index ',begg+i-1
end if
end do
call shr_sys_abort(' ERROR: One or more of the output from CLM to the coupler are NaN ' )
end if
end subroutine check_for_nans
end module lnd_import_export_utils