-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathw3agcmmd.ftn
339 lines (337 loc) · 13.2 KB
/
w3agcmmd.ftn
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
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
#include "w3macros.h"
!/ ------------------------------------------------------------------- /
MODULE W3AGCMMD
!/
!/ +-----------------------------------+
!/ | WAVEWATCH III NOAA/NCEP |
!/ | J. Pianezze |
!/ | FORTRAN 90 |
!/ | Last update : April-2016 |
!/ +-----------------------------------+
!/
!/ Mar-2014 : Origination. ( version 4.18 )
!/ For upgrades see subroutines.
!/ Apr-2016 : Add comments (J. Pianezze) ( version 5.07 )
!/
!/ Copyright 2009-2012 National Weather Service (NWS),
!/ National Oceanic and Atmospheric Administration. All rights
!/ reserved. WAVEWATCH III is a trademark of the NWS.
!/ No unauthorized use without permission.
!/
! 1. Purpose :
!
! Module used for coupling applications between atmospheric model and WW3 with OASIS3-MCT
!
! 2. Variables and types :
!
! 3. Subroutines and functions :
!
! Name Type Scope Description
! ----------------------------------------------------------------
! SND_FIELDS_TO_ATMOS Subr. Public Send fields to atmos model
! RCV_FIELDS_FROM_ATMOS Subr. Public Receive fields from atmos model
! ----------------------------------------------------------------
!
! 4. Subroutines and functions used :
!
! Name Type Module Description
! ----------------------------------------------------------------
! CPL_OASIS_SEND Subr. W3OACPMD Send fields
! CPL_OASIS_RECV Subr. W3OACPMD Receive fields
! ----------------------------------------------------------------
!
! 5. Remarks
! 6. Switches :
! 7. Source code :
!
!/ ------------------------------------------------------------------- /
!
IMPLICIT NONE
!
INCLUDE "mpif.h"
!
PRIVATE
!
! * Accessibility
PUBLIC SND_FIELDS_TO_ATMOS
PUBLIC RCV_FIELDS_FROM_ATMOS
!
CONTAINS
!/ ------------------------------------------------------------------- /
SUBROUTINE SND_FIELDS_TO_ATMOS()
!/
!/ +-----------------------------------+
!/ | WAVEWATCH III NOAA/NCEP |
!/ | J. Pianezze |
!/ | FORTRAN 90 |
!/ | Last update : Apr-2016 |
!/ +-----------------------------------+
!/
!/ Mar-2014 : Origination. ( version 4.18 )
!/ Apr-2016 : Add comments (J. Pianezze) ( version 5.07 )
!/
! 1. Purpose :
!
! Send coupling fields to atmospheric model
!
! 2. Method :
! 3. Parameters :
! 4. Subroutines used :
!
! Name Type Module Description
! -------------------------------------------------------------------
! CPL_OASIS_SND Subr. W3OACPMD Send field to atmos/ocean model
! -------------------------------------------------------------------
!
! 5. Called by :
!
! Name Type Module Description
! ------------------------------------------------------------------
! W3WAVE Subr. W3WAVEMD Wave model
! ------------------------------------------------------------------
!
! 6. Error messages :
! 7. Remarks :
! 8. Structure :
! 9. Switches :
! 10. Source code :
!
!/ ------------------------------------------------------------------- /
!
USE W3OACPMD, ONLY: ID_OASIS_TIME, IL_NB_SND, SND_FLD, CPL_OASIS_SND
USE W3GDATMD, ONLY: NSEAL, NSEA
USE W3ADATMD, ONLY: CX, CY, CHARN, HS, FP0, TWS
USE W3ODATMD, ONLY: UNDEF, NAPROC, IAPROC
!
!/ ------------------------------------------------------------------- /
!/ Local parameters
!/
REAL(kind=8), DIMENSION(NSEAL,1) :: RLA_OASIS_SND
INTEGER :: IB_DO
LOGICAL :: LL_ACTION
REAL(kind=8), DIMENSION(NSEAL) :: TMP
INTEGER :: JSEA, ISEA
!
!----------------------------------------------------------------------
! * Executable part
!
DO IB_DO = 1, IL_NB_SND
!
! Ocean sea surface current (m.s-1) (u-component)
! ---------------------------------------------------------------------
IF (SND_FLD(IB_DO)%CL_FIELD_NAME == 'WW3_WSSU') THEN
TMP(1:NSEAL) = 0.0
DO JSEA=1, NSEAL
ISEA=IAPROC+(JSEA-1)*NAPROC
IF(CX(ISEA) /= UNDEF) TMP(JSEA)=CX(ISEA)
END DO
RLA_OASIS_SND(:,1) = DBLE(TMP(1:NSEAL))
CALL CPL_OASIS_SND(IB_DO, ID_OASIS_TIME, RLA_OASIS_SND, LL_ACTION)
ENDIF
!
! Ocean sea surface current (m.s-1) (v-component)
! ---------------------------------------------------------------------
IF (SND_FLD(IB_DO)%CL_FIELD_NAME == 'WW3_WSSV') THEN
TMP(1:NSEAL) = 0.0
DO JSEA=1, NSEAL
ISEA=IAPROC+(JSEA-1)*NAPROC
IF(CY(ISEA) /= UNDEF) TMP(JSEA)=CY(ISEA)
TMP(JSEA)=CY(ISEA)
END DO
RLA_OASIS_SND(:,1) = DBLE(TMP(1:NSEAL))
CALL CPL_OASIS_SND(IB_DO, ID_OASIS_TIME, RLA_OASIS_SND, LL_ACTION)
ENDIF
!
! Charnock Coefficient (-)
! ---------------------------------------------------------------------
IF (SND_FLD(IB_DO)%CL_FIELD_NAME == 'WW3_ACHA') THEN
TMP(1:NSEAL) = 0.0
WHERE(CHARN /= UNDEF) TMP=CHARN
RLA_OASIS_SND(:,1) = DBLE(TMP(1:NSEAL))
CALL CPL_OASIS_SND(IB_DO, ID_OASIS_TIME, RLA_OASIS_SND, LL_ACTION)
ENDIF
!
! Significant wave height (m)
! ---------------------------------------------------------------------
IF (SND_FLD(IB_DO)%CL_FIELD_NAME == 'WW3__AHS') THEN
TMP(1:NSEAL) = 0.0
WHERE(HS /= UNDEF) TMP=HS
RLA_OASIS_SND(:,1) = DBLE(TMP(1:NSEAL))
CALL CPL_OASIS_SND(IB_DO, ID_OASIS_TIME, RLA_OASIS_SND, LL_ACTION)
ENDIF
!
! Peak frequency (s-1)
! ---------------------------------------------------------------------
IF (SND_FLD(IB_DO)%CL_FIELD_NAME == 'WW3___FP') THEN
TMP(1:NSEAL) = 0.0
WHERE(FP0 /= UNDEF) TMP=FP0
RLA_OASIS_SND(:,1) = DBLE(TMP(1:NSEAL))
CALL CPL_OASIS_SND(IB_DO, ID_OASIS_TIME, RLA_OASIS_SND, LL_ACTION)
ENDIF
!
! Peak period (s)
! ---------------------------------------------------------------------
IF (SND_FLD(IB_DO)%CL_FIELD_NAME == 'WW3___TP') THEN
TMP(1:NSEAL) = 0.0
WHERE(FP0 /= UNDEF) TMP=1./FP0
RLA_OASIS_SND(:,1) = DBLE(TMP(1:NSEAL))
CALL CPL_OASIS_SND(IB_DO, ID_OASIS_TIME, RLA_OASIS_SND, LL_ACTION)
ENDIF
!
! Wind sea Mean period (s)
! ---------------------------------------------------------------------
IF (SND_FLD(IB_DO)%CL_FIELD_NAME == 'WW3__FWS') THEN
TMP(1:NSEAL) = 0.0
WHERE(TWS /= UNDEF) TMP=TWS
RLA_OASIS_SND(:,1) = DBLE(TMP(1:NSEAL))
CALL CPL_OASIS_SND(IB_DO, ID_OASIS_TIME, RLA_OASIS_SND, LL_ACTION)
ENDIF
!
ENDDO
!
!/ ------------------------------------------------------------------- /
END SUBROUTINE SND_FIELDS_TO_ATMOS
!/ ------------------------------------------------------------------- /
SUBROUTINE RCV_FIELDS_FROM_ATMOS(ID_LCOMM, IDFLD, FXN, FYN, FAN)
!/
!/ +-----------------------------------+
!/ | WAVEWATCH III NOAA/NCEP |
!/ | J. Pianezze |
!/ | FORTRAN 90 |
!/ | Last update : April-2016 |
!/ +-----------------------------------+
!/
!/ Mar-2014 : Origination. ( version 4.18 )
!/ Apr-2015 : Modification (M. Accensi) ( version 5.07 )
!/ Apr-2016 : Add comments (J. Pianezze) ( version 5.07 )
!/
! 1. Purpose :
!
! Receive coupling fields from atmospheric model
!
! 2. Method :
! 3. Parameters :
!
! Parameter list
! ----------------------------------------------------------------
! ID_LCOMM Char. I MPI communicator
! IDFLD Int. I Name of the exchange fields
! FXN Int. I/O First exchange field
! FYN Int. I/O Second exchange field
! FAN Int. I/O Third exchange field
! ----------------------------------------------------------------
!
! 4. Subroutines used :
!
! Name Type Module Description
! -------------------------------------------------------------------
! CPL_OASIS_RCV Subr. W3OACPMD Receive fields from atmos/ocean model
! W3S2XY Subr. W3SERVMD Convert from storage (NSEA) to spatial grid (NX, NY)
! -------------------------------------------------------------------
!
! 5. Called by :
!
! Name Type Module Description
! ------------------------------------------------------------------
! W3FLDG Subr. W3FLDSMD Manage input fields of depth,
! current, wind and ice concentration
! ------------------------------------------------------------------
!
! 6. Error messages :
! 7. Remarks :
! 8. Structure :
! 9. Switches :
! 10. Source code :
!
!/ ------------------------------------------------------------------- /
!
USE W3OACPMD, ONLY: ID_OASIS_TIME, IL_NB_RCV, RCV_FLD, CPL_OASIS_RCV
USE W3GDATMD, ONLY: NX, NY, NSEAL, NSEA, MAPSF
USE W3ODATMD, ONLY: NAPROC, IAPROC
USE W3SERVMD, ONLY: W3S2XY
!
!/ ------------------------------------------------------------------- /
!/ Parameter list
!/
INTEGER, INTENT(IN) :: ID_LCOMM
CHARACTER(LEN=3), INTENT(IN) :: IDFLD
REAL, INTENT(INOUT) :: FXN(:,:), FYN(:,:), FAN(:,:)
!
!/ ------------------------------------------------------------------- /
!/ Local parameters
!/
LOGICAL :: LL_ACTION
INTEGER :: IB_DO, IB_I, IB_J, IL_ERR
REAL(kind=8), DIMENSION(NSEAL,1) :: RLA_OASIS_RCV
REAL(kind=8), DIMENSION(NSEAL) :: TMP
REAL, DIMENSION(1:NSEA) :: SND_BUFF,RCV_BUFF
!
!----------------------------------------------------------------------
! * Executable part
!
RLA_OASIS_RCV(:,:) = 0.0
!
DO IB_DO = 1, IL_NB_RCV
IF (IDFLD == 'WND') THEN
!
! Wind speed at 10m (m.s-1) (u-component)
! ----------------------------------------------------------------------
IF (RCV_FLD(IB_DO)%CL_FIELD_NAME == 'WW3__U10') THEN
CALL CPL_OASIS_RCV(IB_DO, ID_OASIS_TIME, RLA_OASIS_RCV, LL_ACTION)
IF (LL_ACTION) THEN
TMP(1:NSEAL) = RLA_OASIS_RCV(1:NSEAL,1)
SND_BUFF(1:NSEA) = 0.0
DO IB_I = 1, NSEAL
IB_J = IAPROC + (IB_I-1)*NAPROC
SND_BUFF(IB_J) = TMP(IB_I)
ENDDO
!
CALL MPI_ALLREDUCE(SND_BUFF(1:NSEA), &
RCV_BUFF(1:NSEA), &
NSEA, &
MPI_REAL, &
MPI_SUM, &
ID_LCOMM, &
IL_ERR)
!
! Convert from storage (NSEA) to spatial grid (NX, NY)
CALL W3S2XY(NSEA,NSEA,NX,NY,RCV_BUFF(1:NSEA),MAPSF,FXN)
!
ENDIF
ENDIF
!
! Wind speed at 10m (m.s-1) (v-component)
! ----------------------------------------------------------------------
IF (RCV_FLD(IB_DO)%CL_FIELD_NAME == 'WW3__V10') THEN
CALL CPL_OASIS_RCV(IB_DO, ID_OASIS_TIME, RLA_OASIS_RCV, LL_ACTION)
IF (LL_ACTION) THEN
TMP(1:NSEAL) = RLA_OASIS_RCV(1:NSEAL,1)
SND_BUFF(1:NSEA) = 0.0
DO IB_I = 1, NSEAL
IB_J = IAPROC + (IB_I-1)*NAPROC
SND_BUFF(IB_J) = TMP(IB_I)
END DO
!
CALL MPI_ALLREDUCE(SND_BUFF(1:NSEA), &
RCV_BUFF(1:NSEA), &
NSEA, &
MPI_REAL, &
MPI_SUM, &
ID_LCOMM, &
IL_ERR)
!
! Convert from storage (NSEA) to spatial grid (NX, NY)
CALL W3S2XY(NSEA,NSEA,NX,NY,RCV_BUFF(1:NSEA),MAPSF,FYN)
!
ENDIF
ENDIF
!
ENDIF
ENDDO
!/ ------------------------------------------------------------------- /
END SUBROUTINE RCV_FIELDS_FROM_ATMOS
!/ ------------------------------------------------------------------- /
!/
END MODULE W3AGCMMD
!/
!/ ------------------------------------------------------------------- /