File: C:\NOAA\NEMS_11731\src\ENS_Cpl\ENS_Sto_Per_Scheme_Step1.f

1      SUBROUTINE ENS_Sto_Per_Scheme_Step1(cst)
2     
3     ! March 2009 Weiyu Yang  Modified for the NEMS model.
4     !----------------------------------------------------
5     
6     ! This subroutine is used to compute the first step of the 
7     ! stochastic perturbation scheme, in which X_i_dot = T_i + S_i 
8     ! and S_i ~ SUM(W_i,j P_j.
9     !-------------------------------------------------------------
10     
11     
12     ! In the current code, assume each ensemble member uses the same 
13     ! number of the processors.
14     !---------------------------------------------------------------
15      USE ENS_Cpl_InternalState_ESMFMod
16     
17      IMPLICIT none
18     
19      TYPE(ENS_Cpl_InternalState),                         INTENT(inout) :: cst
20     
21     ! Working arrays.
22     !----------------
23      INTEGER                                             :: tb
24      INTEGER                                             :: i, j, k, l
25      REAL(KIND=KIND_EVOD), DIMENSION(cst%ARRAY_ONE_SIZ4) :: pstm, psmtm
26      REAL(KIND=KIND_EVOD), DIMENSION(cst%ARRAY_TOT_SIZ4) :: ttm, tmtm, utm, umtm, vtm, vmtm
27     
28      REAL(KIND=KIND_EVOD), DIMENSION(cst%ARRAY_TOT_SIZ4, cst%ntrac) :: tracertm, tracermtm
29     
30     ! Sto_Coef(k,j) is the weight of the k'th perturabtion in member j's formulation. 
31     !--------------------------------------------------------------------------------
32      tb = cst%Total_member
33     
34      DO i = 1, cst%ARRAY_ONE_SIZ4
35          pstm (i) = cst%ps  (i, tb) - cst%ps6  (i, tb)
36          psmtm(i) = cst%psm (i, tb) - cst%ps6m (i, tb)
37      END DO
38     
39      DO i = 1, cst%ARRAY_TOT_SIZ4
40          ttm (i) = cst%t (i, tb) - cst%t6 (i, tb)
41          tmtm(i) = cst%tm(i, tb) - cst%t6m(i, tb)
42          utm (i) = cst%u (i, tb) - cst%u6 (i, tb)
43          umtm(i) = cst%um(i, tb) - cst%u6m(i, tb)
44          vtm (i) = cst%v (i, tb) - cst%v6 (i, tb)
45          vmtm(i) = cst%vm(i, tb) - cst%v6m(i, tb)
46      END DO
47     
48      DO j = 1, cst%ntrac
49          DO i = 1, cst%ARRAY_TOT_SIZ4
50              tracertm (i, j) = cst%tracer (i, tb, j) - cst%tracer6 (i, tb, j)
51              tracermtm(i, j) = cst%tracerm(i, tb, j) - cst%tracer6m(i, tb, j)
52          END DO
53      END DO
54     
55      DO i = 1, cst%ARRAY_ONE_SIZ4
56          cst%ps_step1 (i, tb) = 0.0
57          cst%psm_step1(i, tb) = 0.0
58      END DO
59      DO i = 1, cst%ARRAY_TOT_SIZ4
60          cst%t_step1   (i, tb) = 0.0
61          cst%u_step1   (i, tb) = 0.0
62          cst%v_step1   (i, tb) = 0.0
63          cst%tm_step1  (i, tb) = 0.0
64          cst%um_step1  (i, tb) = 0.0
65          cst%vm_step1  (i, tb) = 0.0
66      END DO
67     
68      DO j = 1, cst%ntrac
69          DO i = 1, cst%ARRAY_TOT_SIZ4
70              cst%tracer_step1 (i, tb, j) = 0.0
71              cst%tracerm_step1(i, tb, j) = 0.0
72          END DO
73      END DO
74     
75      DO j = 1, tb - 1     !for the forcing of each member j
76          DO i = 1, cst%ARRAY_ONE_SIZ4
77              cst%ps_step1 (i, j) = 0.0
78              cst%psm_step1(i, j) = 0.0
79          END DO
80          DO k = 1, tb - 1
81              DO i = 1, cst%ARRAY_ONE_SIZ4
82                  cst%ps_step1  (i, j) = cst%ps_step1  (i, j) + cst%Sto_Coef(k, j) *     &
83                      (cst%ps (i, k) - cst%ps6 (i, k) - pstm(i))
84                  cst%psm_step1 (i, j) = cst%psm_step1 (i, j) + cst%Sto_Coef(k, j) *     &
85                      (cst%psm(i, k) - cst%ps6m(i, k) - psmtm(i))
86              END DO
87          END DO
88     
89          DO i = 1, cst%ARRAY_TOT_SIZ4
90              cst%t_step1   (i, j) = 0.0
91              cst%u_step1   (i, j) = 0.0
92              cst%v_step1   (i, j) = 0.0
93              cst%tm_step1  (i, j) = 0.0
94              cst%um_step1  (i, j) = 0.0
95              cst%vm_step1  (i, j) = 0.0
96          END DO
97          DO k = 1, tb - 1
98              DO i = 1, cst%ARRAY_TOT_SIZ4
99                  cst%t_step1   (i, j) = cst%t_step1   (i, j) + cst%Sto_Coef(k, j) *     &
100                      (cst%t    (i, k) - cst%t6   (i, k) - ttm(i))
101                  cst%tm_step1  (i, j) = cst%tm_step1  (i, j) + cst%Sto_Coef(k, j) *     &
102                      (cst%tm   (i, k) - cst%t6m  (i, k) - tmtm(i))
103     
104                  cst%u_step1   (i, j) = cst%u_step1   (i, j) + cst%Sto_Coef(k, j) *     &
105                      (cst%u    (i, k) - cst%u6   (i, k) - utm(i))
106                  cst%um_step1  (i, j) = cst%um_step1  (i, j) + cst%Sto_Coef(k, j) *     &
107                      (cst%um   (i, k) - cst%u6m  (i, k) -umtm(i))
108     
109                  cst%v_step1   (i, j) = cst%v_step1   (i, j) + cst%Sto_Coef(k, j) *     &
110                      (cst%v    (i, k) - cst%v6  (i, k) - vtm(i))
111                  cst%vm_step1  (i, j) = cst%vm_step1  (i, j) + cst%Sto_Coef(k, j) *     &
112                      (cst%vm   (i, k) - cst%v6m  (i, k) - vmtm(i))
113              END DO
114          END DO
115      END DO
116     
117      DO l = 1, cst%ntrac
118          DO j = 1, tb - 1     !for the forcing of each member j
119              DO i = 1, cst%ARRAY_TOT_SIZ4
120                  cst%tracer_step1 (i, j, l) = 0.0
121                  cst%tracerm_step1(i, j, l) = 0.0
122              END DO
123          
124              DO k = 1, tb - 1
125                  DO i = 1, cst%ARRAY_TOT_SIZ4
126                      cst%tracer_step1 (i, j, l) = cst%tracer_step1 (i, j, l) + cst%Sto_Coef(k, j) *     &
127                          (cst%tracer (i, k, l) - cst%tracer6 (i, k, l) - tracertm(i, l))
128                      cst%tracerm_step1(i, j, l) = cst%tracerm_step1(i, j, l) + cst%Sto_Coef(k, j) *     &
129                          (cst%tracerm(i, k, l) - cst%tracer6m(i, k, l) - tracermtm(i, l))
130                  END DO
131              END DO
132          END DO
133      END DO
134     
135     ! New feature in B series  DHOU 03-28-2007
136     ! Centralize the perturbations w_step1
137     ! Find the average of the perturbations.
138     !-----------------------------------------
139      IF(cst%CENTER == 1) THEN
140          DO i = 1, cst%ARRAY_ONE_SIZ4
141              cst%ps_mean(i)  = 0.0
142              cst%psm_mean(i) = 0.0
143          END DO
144          DO j = 1, tb - 1     !for the forcing of each member j
145              DO i = 1, cst%ARRAY_ONE_SIZ4
146                  cst%ps_mean(i)  = cst%ps_mean(i)  + cst%ps_step1 (i, j)
147                  cst%psm_mean(i) = cst%psm_mean(i) + cst%psm_step1(i, j)
148              END DO
149          END  DO
150          DO i = 1, cst%ARRAY_ONE_SIZ4
151              cst%ps_mean (i) = cst%ps_mean (i) / FLOAT(tb - 1)
152              cst%psm_mean(i) = cst%psm_mean(i) / FLOAT(tb - 1)
153          END DO
154     
155     ! subtract the PS average from the perturbations.
156     !------------------------------------------------
157          DO j = 1, tb - 1  
158              DO i = 1, cst%ARRAY_ONE_SIZ4
159                  cst%ps_step1 (i, j) = cst%ps_step1 (i, j) - cst%ps_mean(i)
160                  cst%psm_step1(i, j) = cst%psm_step1(i, j) - cst%psm_mean(i)
161              END DO
162          END  DO
163     
164          DO i = 1, cst%ARRAY_TOT_SIZ4
165              cst%t_mean (i) = 0.0
166              cst%tm_mean(i) = 0.0
167              cst%u_mean (i) = 0.0
168              cst%um_mean(i) = 0.0
169              cst%v_mean (i) = 0.0
170              cst%vm_mean(i) = 0.0
171          END DO
172          DO j = 1, cst%ntrac
173              DO i = 1, cst%ARRAY_TOT_SIZ4
174                  cst%tracer_mean (i, j) = 0.0
175                  cst%tracerm_mean(i, j) = 0.0
176              END DO
177          END DO
178     
179     
180          DO j = 1, tb - 1     !for the forcing of each member j
181              DO i = 1, cst%ARRAY_TOT_SIZ4
182                  cst%t_mean (i) = cst%t_mean (i) + cst%t_step1 (i, j)
183                  cst%tm_mean(i) = cst%tm_mean(i) + cst%tm_step1(i, j)
184                  cst%u_mean (i) = cst%u_mean (i) + cst%u_step1 (i, j)
185                  cst%um_mean(i) = cst%um_mean(i) + cst%um_step1(i, j)
186                  cst%v_mean (i) = cst%v_mean (i) + cst%v_step1 (i, j)
187                  cst%vm_mean(i) = cst%vm_mean(i) + cst%vm_step1(i, j)
188              END DO
189          END DO
190     
191          DO l = 1, cst%ntrac
192              DO j = 1, tb - 1 
193                  DO i = 1, cst%ARRAY_TOT_SIZ4
194                      cst%tracer_mean (i, l) = cst%tracer_mean (i, l) + cst%tracer_step1 (i, j, l)
195                      cst%tracerm_mean(i, l) = cst%tracerm_mean(i, l) + cst%tracerm_step1(i, j, l)
196                  END DO
197              END DO
198          END DO
199          
200          DO i = 1, cst%ARRAY_TOT_SIZ4
201              cst%t_mean (i) = cst%t_mean (i) / FLOAT(tb - 1)
202              cst%tm_mean(i) = cst%tm_mean(i) / FLOAT(tb - 1)
203              cst%u_mean (i) = cst%u_mean (i) / FLOAT(tb - 1)
204              cst%um_mean(i) = cst%um_mean(i) / FLOAT(tb - 1)
205              cst%v_mean (i) = cst%v_mean (i) / FLOAT(tb - 1)
206              cst%vm_mean(i) = cst%vm_mean(i) / FLOAT(tb - 1)
207          END DO
208     
209          DO l = 1, cst%ntrac
210              DO i = 1, cst%ARRAY_TOT_SIZ4
211                  cst%tracer_mean (i, l) = cst%tracer_mean (i, l) / FLOAT(tb - 1)
212                  cst%tracerm_mean(i, l) = cst%tracerm_mean(i, l) / FLOAT(tb - 1)
213              END DO
214          END DO
215          
216     ! subtract the average from the perturbations.
217     !---------------------------------------------
218          DO j = 1, tb - 1  
219              DO i = 1, cst%ARRAY_TOT_SIZ4
220                  cst%t_step1 (i, j) = cst%t_step1 (i, j) - cst%t_mean (i)
221                  cst%tm_step1(i, j) = cst%tm_step1(i, j) - cst%tm_mean(i)
222                  cst%u_step1 (i, j) = cst%u_step1 (i, j) - cst%u_mean (i)
223                  cst%um_step1(i, j) = cst%um_step1(i, j) - cst%um_mean(i)
224                  cst%v_step1 (i, j) = cst%v_step1 (i, j) - cst%v_mean (i)
225                  cst%vm_step1(i, j) = cst%vm_step1(i, j) - cst%vm_mean(i)
226              END DO
227          END  DO
228     
229          DO l = 1, cst%ntrac
230              DO j = 1, tb - 1 
231                  DO i = 1, cst%ARRAY_TOT_SIZ4
232                      cst%tracer_step1 (i, j, l) = cst%tracer_step1 (i, j, l) - cst%tracer_mean (i, l)
233                      cst%tracerm_step1(i, j, l) = cst%tracerm_step1(i, j, l) - cst%tracerm_mean(i, l)
234                  END DO
235              END DO
236          END DO
237      END IF
238     
239      END SUBROUTINE ENS_Sto_Per_Scheme_Step1
240     
241     
242     
243     
244     
245      SUBROUTINE ENS_Sto_Per_Scheme_Step1_2(impENS, cst, rc)
246     
247     ! This subroutine is used to compute the first step of the 
248     ! stochastic perturbation scheme, in which X_i_dot = T_i + S_i 
249     ! and S_i ~ SUM(W_i,j P_j.
250     !-------------------------------------------------------------
251     
252     ! In the current code, assume each ensemble member uses the same 
253     ! number of the processors.
254     !---------------------------------------------------------------
255     
256      USE ESMF_Mod
257      USE ENS_Cpl_InternalState_ESMFMod
258      USE Lib_ESMFStateAddGetMod
259     
260      IMPLICIT none
261     
262      TYPE(ESMF_State),             INTENT(inout) :: impENS
263      TYPE(ENS_Cpl_InternalState),  INTENT(inout) :: cst
264      INTEGER,                      INTENT(out)   :: rc
265     
266     ! !WORKING ARRAYS AND LOCAL PARAMETERS.
267     !--------------------------------------
268      TYPE(ESMF_Field)                            :: Field
269      TYPE(ESMF_FieldBundle)                      :: Bundle
270      INTEGER                                     :: i, j, k, l, i1, j1
271      INTEGER                                     :: rc1
272     
273      rc1 = ESMF_SUCCESS
274      rc  = ESMF_SUCCESS
275     
276      CALL GetF90ArrayFromState(impENS, 'pps', cst%work1, 0, rc = rc1)
277     
278          IF(ESMF_LogMsgFoundAllocError(rc1, "Get From impENS to the Cpl Internal State -- PS.")) THEN
279              rc = ESMF_FAILURE
280              PRINT*, 'Error Happened When Getting From impENS to the Cpl Internal State -- PS, rc = ', rc1
281              rc1 = ESMF_SUCCESS
282          END IF
283     
284      DO j = 1, cst%arraysize_2
285          DO i = 1, cst%arraysize_1
286              cst%psw(i, j) = cst%work1(i, j) + cst%RS_GLOBAL * cst%psw(i, j)
287          END DO
288      END DO
289     
290      CALL GetF90ArrayFromState(impENS, 'psm', cst%work1, 0, rc = rc1)
291     
292          IF(ESMF_LogMsgFoundAllocError(rc1, "Get From impENS to the Cpl Internal State -- PS(-DT).")) THEN
293              rc = ESMF_FAILURE
294              PRINT*, 'Error Happened When Getting From impENS to the Cpl Internal State -- PS(-DT), rc = ', rc1
295              rc1 = ESMF_SUCCESS
296          END IF
297     
298      DO j = 1, cst%arraysize_2
299          DO i = 1, cst%arraysize_1
300              cst%pswm(i, j) = cst%work1(i, j) + cst%RS_GLOBAL * cst%pswm(i, j)
301          END DO
302      END DO
303     
304      CALL GetF90ArrayFromState(impENS, 'tt', cst%work2, 0, rc = rc1)
305     
306          IF(ESMF_LogMsgFoundAllocError(rc1, "Get From impENS to the Cpl Internal State -- T.")) THEN
307              rc = ESMF_FAILURE
308              PRINT*, 'Error Happened When Getting From impENS to the Cpl Internal State -- T, rc = ', rc1
309              rc1 = ESMF_SUCCESS
310          END IF
311     
312      DO k = 1, cst%arraysize_3
313          DO j = 1, cst%arraysize_2
314              DO i = 1, cst%arraysize_1
315                  cst%tw(i, j, k) = cst%work2(i, j, k) + cst%RS_GLOBAL * cst%tw(i, j, k)
316              END DO
317          END DO
318      END DO
319     
320      CALL GetF90ArrayFromState(impENS, 'tm', cst%work2, 0, rc = rc1)
321     
322          IF(ESMF_LogMsgFoundAllocError(rc1, "Get From impENS to the Cpl Internal State -- T(-DT).")) THEN
323              rc = ESMF_FAILURE
324              PRINT*, 'Error Happened When Getting From impENS to the Cpl Internal State -- T(-DT), rc = ', rc1
325              rc1 = ESMF_SUCCESS
326          END IF
327     
328      DO k = 1, cst%arraysize_3
329          DO j = 1, cst%arraysize_2
330              DO i = 1, cst%arraysize_1
331                  cst%twm(i, j, k) = cst%work2(i, j, k) + cst%RS_GLOBAL * cst%twm(i, j, k)
332              END DO
333          END DO
334      END DO
335     
336      CALL GetF90ArrayFromState(impENS, 'uu', cst%work2, 0, rc = rc1)
337     
338          IF(ESMF_LogMsgFoundAllocError(rc1, "Get From impENS to the Cpl Internal State -- U.")) THEN
339              rc = ESMF_FAILURE
340              PRINT*, 'Error Happened When Getting From impENS to the Cpl Internal State -- U, rc = ', rc1
341              rc1 = ESMF_SUCCESS
342          END IF
343     
344      DO k = 1, cst%arraysize_3
345          DO j = 1, cst%arraysize_2
346              DO i = 1, cst%arraysize_1
347                  cst%uw(i, j, k) = cst%work2(i, j, k) + cst%RS_GLOBAL * cst%uw(i, j, k)
348              END DO
349          END DO
350      END DO
351     
352      CALL GetF90ArrayFromState(impENS, 'um', cst%work2, 0, rc = rc1)
353     
354          IF(ESMF_LogMsgFoundAllocError(rc1, "Get From impENS to the Cpl Internal State -- U(-DT).")) THEN
355              rc = ESMF_FAILURE
356              PRINT*, 'Error Happened When Getting From impENS to the Cpl Internal State -- U(-DT), rc = ', rc1
357              rc1 = ESMF_SUCCESS
358          END IF
359     
360      DO k = 1, cst%arraysize_3
361          DO j = 1, cst%arraysize_2
362              DO i = 1, cst%arraysize_1
363                  cst%uwm(i, j, k) = cst%work2(i, j, k) + cst%RS_GLOBAL * cst%uwm(i, j, k)
364              END DO
365          END DO
366      END DO
367     
368      CALL GetF90ArrayFromState(impENS, 'vv', cst%work2, 0, rc = rc1)
369     
370          IF(ESMF_LogMsgFoundAllocError(rc1, "Get From impENS to the Cpl Internal State -- V.")) THEN
371              rc = ESMF_FAILURE
372              PRINT*, 'Error Happened When Getting From impENS to the Cpl Internal State -- V, rc = ', rc1
373              rc1 = ESMF_SUCCESS
374          END IF
375     
376      DO k = 1, cst%arraysize_3
377          DO j = 1, cst%arraysize_2
378              DO i = 1, cst%arraysize_1
379                  cst%vw(i, j, k) = cst%work2(i, j, k) + cst%RS_GLOBAL * cst%vw(i, j, k)
380              END DO
381          END DO
382      END DO
383     
384      CALL GetF90ArrayFromState(impENS, 'vm', cst%work2, 0, rc = rc1)
385     
386          IF(ESMF_LogMsgFoundAllocError(rc1, "Get From impENS to the Cpl Internal State -- V(-DT).")) THEN
387              rc = ESMF_FAILURE
388              PRINT*, 'Error Happened When Getting From impENS to the Cpl Internal State -- V(-DT), rc = ', rc1
389              rc1 = ESMF_SUCCESS
390          END IF
391     
392      DO k = 1, cst%arraysize_3
393          DO j = 1, cst%arraysize_2
394              DO i = 1, cst%arraysize_1
395                  cst%vwm(i, j, k) = cst%work2(i, j, k) + cst%RS_GLOBAL * cst%vwm(i, j, k)
396              END DO
397          END DO
398      END DO
399     
400      CALL ESMF_StateGet(impENS, 'tracers', Bundle, rc = rc1)
401     
402          IF(ESMF_LogMsgFoundAllocError(rc1, "Retrieve Bundle from impENS in ENS Cpl.")) THEN
403              rc = ESMF_FAILURE
404              PRINT*, 'Error Happened When Retrieving Bundle from impENS in ENS Cpl, rc = ', rc1
405              rc1 = ESMF_SUCCESS
406          END IF
407     
408      DO j1 = 2, 3
409          DO i1 = 1, cst%ntrac
410              CALL ESMF_FieldBundleGet(Bundle,                    &
411                                       cst%vname(i1, j1),         &
412                                       Field,                     &
413                                       rc = rc1)
414     
415              IF(ESMF_LogMsgFoundAllocError(rc1, "Retrieve Field from Bundle in ENS Cpl.")) THEN
416                  rc = ESMF_FAILURE
417                  PRINT*, 'Error Happened When Retrieving Field from Bundle in ENS Cpl, rc = ', rc1
418                  rc1 = ESMF_SUCCESS
419              END IF
420     
421              NULLIFY(cst%work2)
422              CALL ESMF_FieldGet(Field, FArray = cst%work2, localDE = 0, rc = rc1)
423     
424              IF(ESMF_LogMsgFoundAllocError(rc1, "Retrieve FArray from Field in ENS Cpl.")) THEN
425                  rc = ESMF_FAILURE
426                  PRINT*, 'Error Happened When Retrieving FArray from Field in ENS Cpl, rc = ', rc1
427                  rc1 = ESMF_SUCCESS
428              END IF
429     
430              SELECT CASE(j1)
431                  CASE(2)
432                      DO k = 1, cst%arraysize_3
433                          DO j = 1, cst%arraysize_2
434                              DO i = 1, cst%arraysize_1
435                                  cst%tracerw (i, j, k, i1) = cst%work2(i, j, k) + cst%RS_GLOBAL * cst%tracerw (i, j, k, i1)
436                              END DO
437                          END DO
438                      END DO
439                  CASE(3)
440                      DO k = 1, cst%arraysize_3
441                          DO j = 1, cst%arraysize_2
442                              DO i = 1, cst%arraysize_1
443                                  cst%tracerwm(i, j, k, i1) = cst%work2(i, j, k) + cst%RS_GLOBAL * cst%tracerwm(i, j, k, i1)
444                              END DO
445                          END DO
446                      END DO
447              END SELECT
448     
449              IF(ESMF_LogMsgFoundAllocError(rc1, "Distribute for the Step1 Scheme -- tracer.")) THEN
450                  rc = ESMF_FAILURE
451                  PRINT*, 'Error Happened When Distributing for the Step1 Scheme -- tracer, rc = ', rc1
452                  rc1 = ESMF_SUCCESS
453              END IF
454          END DO
455      END DO
456     
457      IF(rc == ESMF_SUCCESS) THEN
458          PRINT*, "PASS: ENS_Sto_Per_Scheme_Step1_2"
459      ELSE
460          PRINT*, "FAIL: ENS_Sto_Per_Scheme_Step1_2"
461      END IF
462     
463      END SUBROUTINE ENS_Sto_Per_Scheme_Step1_2
464