-
Notifications
You must be signed in to change notification settings - Fork 166
/
Copy pathmodule_cplscalars.F90
203 lines (165 loc) · 8.65 KB
/
module_cplscalars.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
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
!> @file
!> @brief Manage cpl_scalars
!> @author [email protected]
!> @author modified for FV3atm by [email protected] @date 03-03-2024
!> Manage scalars in import and export states. Called by realizeConnectedCplFields
!> to set the required scalar data into a state. The scalar_value will be set into
!> a field with name flds_scalar_name. The scalar_id identifies which dimension in
!> the scalar field is given by the scalar_value. The number of scalars is used to
!> ensure that the scalar_id is within the bounds of the scalar field
module module_cplscalars
use ESMF, only : ESMF_Field, ESMF_Distgrid, ESMF_Grid, ESMF_State
use ESMF, only : ESMF_VM, ESMF_VMGetCurrent, ESMF_VMGet, ESMF_VMBroadCast
use ESMF, only : ESMF_LogFoundError, ESMF_LOGERR_PASSTHRU, ESMF_FAILURE, ESMF_SUCCESS
use ESMF, only : ESMF_LOGMSG_INFO, ESMF_LOGWRITE, ESMF_TYPEKIND_R8, ESMF_KIND_R8
use ESMF, only : ESMF_GridCreate, ESMF_FieldCreate, ESMF_StateGet, ESMF_DistGridCreate
use ESMF, only : ESMF_FieldGet
implicit none
private
public SetScalarField
public State_SetScalar
public State_GetScalar
! set from config
integer, public :: flds_scalar_num, flds_scalar_index_nx
integer, public :: flds_scalar_index_ny, flds_scalar_index_ntile
character(len=80), public :: flds_scalar_name
contains
!================================================================================
!> Create a scalar field
!>
!> @param[inout] field an ESMF_Field
!> @param[in] flds_scalar_name the name of the scalar
!> @param[in] flds_scalar_num the number of scalars
!> @param[inout] rc a return code
!>
!> @author [email protected], [email protected]
!> @date 03-03-2024
subroutine SetScalarField(field, flds_scalar_name, flds_scalar_num, rc)
type(ESMF_Field) , intent(inout) :: field
character(len=*) , intent(in) :: flds_scalar_name
integer , intent(in) :: flds_scalar_num
integer , intent(inout) :: rc
! local variables
type(ESMF_Distgrid) :: distgrid
type(ESMF_Grid) :: grid
character(len=*), parameter :: subname='(module_cplscalars:SetScalarField)'
! ----------------------------------------------
rc = ESMF_SUCCESS
! create a DistGrid with a single index space element, which gets mapped onto DE 0.
distgrid = ESMF_DistGridCreate(minIndex=(/1/), maxIndex=(/1/), rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
grid = ESMF_GridCreate(distgrid, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
field = ESMF_FieldCreate(name=trim(flds_scalar_name), grid=grid, typekind=ESMF_TYPEKIND_R8, &
ungriddedLBound=(/1/), ungriddedUBound=(/flds_scalar_num/), gridToFieldMap=(/2/), rc=rc) ! num of scalar values
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
end subroutine SetScalarField
!================================================================================
!> Set scalar data into a state
!>
!> @param[inout] State an ESMF_State
!> @param[in] scalar_value the value of the scalar
!> @param[in] scalar_id the identity of the scalar
!> @param[in] flds_scalar_name the name of the scalar
!> @param[in] flds_scalar_num the number of scalars
!> @param[inout] rc a return code
!>
!> @author [email protected], [email protected]
!> @date 03-02-2024
subroutine State_SetScalar(scalar_value, scalar_id, State, flds_scalar_name, flds_scalar_num, rc)
! input/output arguments
real(ESMF_KIND_R8), intent(in) :: scalar_value
integer, intent(in) :: scalar_id
type(ESMF_State), intent(inout) :: State
character(len=*), intent(in) :: flds_scalar_name
integer, intent(in) :: flds_scalar_num
integer, intent(inout) :: rc
! local variables
integer :: mytask
type(ESMF_Field) :: lfield
type(ESMF_VM) :: vm
real(ESMF_KIND_R8), pointer :: farrayptr(:,:)
character(len=*), parameter :: subname = ' (module_cplscalars:state_setscalar) '
! ----------------------------------------------
rc = ESMF_SUCCESS
call ESMF_VMGetCurrent(vm, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
call ESMF_VMGet(vm, localPet=mytask, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
call ESMF_StateGet(State, itemName=trim(flds_scalar_name), field=lfield, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
if (mytask == 0) then
call ESMF_FieldGet(lfield, farrayPtr = farrayptr, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
if (scalar_id <= 0 .or. scalar_id > flds_scalar_num) then
call ESMF_LogWrite(trim(subname)//": ERROR in scalar_id", ESMF_LOGMSG_INFO)
rc = ESMF_FAILURE
return
endif
farrayptr(scalar_id,1) = scalar_value
endif
end subroutine State_SetScalar
!===============================================================================
!> Get scalar data from a state
!>
!> @details Obtain the field flds_scalar_name from a State and broadcast and
!> it to all PEs
!>
!> @param[in] State an ESMF_State
!> @param[in] scalar_value the value of the scalar
!> @param[in] scalar_id the identity of the scalar
!> @param[in] flds_scalar_name the name of the scalar
!> @param[in] flds_scalar_num the number of scalars
!> @param[out] rc a return code
!>
!> @author [email protected], [email protected]
!> @date 03-02-2024
subroutine State_GetScalar(state, scalar_id, scalar_value, flds_scalar_name, flds_scalar_num, rc)
! ----------------------------------------------
! Get scalar data from State for a particular name and broadcast it to all other pets
! ----------------------------------------------
! input/output variables
type(ESMF_State), intent(in) :: state
integer, intent(in) :: scalar_id
real(ESMF_KIND_R8), intent(out) :: scalar_value
character(len=*), intent(in) :: flds_scalar_name
integer, intent(in) :: flds_scalar_num
integer, intent(inout) :: rc
! local variables
integer :: mytask, ierr, icount
type(ESMF_VM) :: vm
type(ESMF_Field) :: field
real(ESMF_KIND_R8), pointer :: farrayptr(:,:)
real(ESMF_KIND_R8) :: tmp(1)
character(len=*), parameter :: subname = ' (module_cplscalars:state_getscalar) '
! ----------------------------------------------
rc = ESMF_SUCCESS
call ESMF_VMGetCurrent(vm, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
call ESMF_VMGet(vm, localPet=mytask, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
! check item exist or not?
call ESMF_StateGet(State, itemSearch=trim(flds_scalar_name), itemCount=icount, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
if (icount > 0) then
call ESMF_StateGet(State, itemName=trim(flds_scalar_name), field=field, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
if (mytask == 0) then
call ESMF_FieldGet(field, farrayPtr = farrayptr, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
if (scalar_id < 0 .or. scalar_id > flds_scalar_num) then
call ESMF_LogWrite(trim(subname)//": ERROR in scalar_id", ESMF_LOGMSG_INFO, line=__LINE__, file=__FILE__)
rc = ESMF_FAILURE
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
endif
tmp(:) = farrayptr(scalar_id,:)
endif
call ESMF_VMBroadCast(vm, tmp, 1, 0, rc=rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
scalar_value = tmp(1)
else
scalar_value = 0.0_ESMF_KIND_R8
call ESMF_LogWrite(trim(subname)//": no ESMF_Field found named: "//trim(flds_scalar_name), ESMF_LOGMSG_INFO)
end if
end subroutine State_GetScalar
end module module_cplscalars