-
Notifications
You must be signed in to change notification settings - Fork 2
/
cpl_land.f90
executable file
·123 lines (94 loc) · 3.31 KB
/
cpl_land.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
subroutine ini_land(istart)
! subroutine ini_land (istart)
!
! Input : istart = restart flag ( 0 = no, 1 = yes)
use mod_atparam
use mod_var_land, only: stlcl_ob, stl_lm
implicit none
integer, intent(in) :: istart
! 1. Compute climatological fields for initial date
call atm2land(0)
! 2. Initialize prognostic variables of land model
! in case of no restart or no coupling
if (istart.le.0 .or. istart == 2) then
stl_lm(:) = stlcl_ob(:) ! land sfc. temperature
end if
! 3. Compute additional land variables
call land2atm(0)
end
subroutine atm2land(jday)
use mod_cpl_flags, only: icland
use mod_atparam
use mod_cpl_land_model, only: vland_input
use mod_flx_land, only: hflux_l
use mod_cli_land, only: stl12, snowd12, soilw12
use mod_date, only: imont1, tmonth
use mod_var_land, only: stlcl_ob, snowdcl_ob, soilwcl_ob, stl_lm
implicit none
integer, intent(in) :: jday
integer, parameter :: nlon=ix, nlat=il, ngp=nlon*nlat
! 1. Interpolate climatological fields to actual date
! Climatological land sfc. temperature
call forin5(ngp,imont1,tmonth,stl12,stlcl_ob)
! Climatological snow depth
call forint(ngp,imont1,tmonth,snowd12,snowdcl_ob)
! Climatological soil water availability
call forint(ngp,imont1,tmonth,soilw12,soilwcl_ob)
if (jday.le.0) return
! 2. Set input variables for mixed-layer/ocean model
if (icland.gt.0) then
vland_input(:,1) = stl_lm(:)
vland_input(:,2) = hflux_l(:)
vland_input(:,3) = stlcl_ob(:)
end if
! 3. Call message-passing routines to send data (if needed)
end
subroutine land2atm(jday)
use mod_cpl_flags, only: icland
use mod_atparam
use mod_cpl_land_model, only: land_model, vland_output
use mod_var_land
implicit none
integer, intent(in) :: jday
if (jday.gt.0.and.icland.gt.0) then
! 1. Run ocean mixed layer or
! call message-passing routines to receive data from ocean model
call land_model
! 2. Get updated variables for mixed-layer/ocean model
stl_lm(:) = vland_output(:,1) ! land sfc. temperature
end if
! 3. Compute land-sfc. fields for atm. model
! 3.1 Land sfc. temperature
if (icland.le.0) then
! Use observed climatological field
stl_am(:) = stlcl_ob(:)
else
! Use land model sfc. temperature
stl_am(:) = stl_lm(:)
end if
! 3.2 Snow depth and soil water availability
snowd_am(:) = snowdcl_ob(:)
soilw_am(:) = soilwcl_ob(:)
end
subroutine rest_land(imode)
! subroutine rest_land (imode)
! Purpose : read/write land variables from/to a restart file
! Input : IMODE = 0 : read model variables from a restart file
! = 1 : write model variables to a restart file
use mod_cpl_flags, only: icland
use mod_atparam
use mod_var_land, only: stl_am, stl_lm
implicit none
integer, intent(in) :: imode
if (imode.eq.0) then
read (3) stl_lm(:) ! Land sfc. temperature
else
! Write land model variables from coupled runs,
! otherwise write fields used by atmospheric model
if (icland.gt.0) then
write (10) stl_lm(:)
else
write (10) stl_am(:)
end if
end if
end