solve_nmm.F
References to this file elsewhere.
1 !-----------------------------------------------------------------------
2 !
3 !NCEP_MESO:MEDIATION_LAYER:SOLVER
4 !
5 !-----------------------------------------------------------------------
6 #include "nmm_loop_basemacros.h"
7 #include "nmm_loop_macros.h"
8 !-----------------------------------------------------------------------
9 !
10 SUBROUTINE SOLVE_NMM(GRID,CONFIG_FLAGS &
11 !
12 #include "nmm_dummy_args.inc"
13 !
14 & )
15 !-----------------------------------------------------------------------
16 USE MODULE_DOMAIN
17 USE MODULE_CONFIGURE
18 USE MODULE_MODEL_CONSTANTS
19 USE MODULE_STATE_DESCRIPTION
20 USE MODULE_CTLBLK
21 USE MODULE_DM
22 USE MODULE_IGWAVE_ADJUST, ONLY: PDTE,PFDHT,DDAMP,VTOA
23 USE MODULE_ADVECTION, ONLY: ADVE,VAD2,HAD2,VAD2_SCAL,HAD2_SCAL
24 USE MODULE_NONHY_DYNAM, ONLY: EPS,VADZ,HADZ
25 USE MODULE_DIFFUSION_NMM, ONLY: HDIFF
26 USE MODULE_BNDRY_COND, ONLY: BOCOH,BOCOV
27 USE MODULE_PHYSICS_CALLS
28 USE MODULE_EXT_INTERNAL
29 USE MODULE_PRECIP_ADJUST
30 USE MODULE_NEST_UTIL
31 #ifdef WRF_CHEM
32 USE MODULE_INPUT_CHEM_DATA, ONLY: GET_LAST_GAS
33 #endif
34 !-----------------------------------------------------------------------
35 !
36 IMPLICIT NONE
37 !
38 !-----------------------------------------------------------------------
39 INCLUDE "mpif.h"
40 !-----------------------------------------------------------------------
41 !
42 !*** INPUT DATA
43 !
44 !-----------------------------------------------------------------------
45 !
46 TYPE(DOMAIN),TARGET :: GRID
47 !
48 !*** DEFINITIONS OF DUMMY ARGUMENTS TO THIS ROUTINE (GENERATED FROM REGISTRY)
49 !
50 ! NOTE, REGISTRY NO LONGER GENERATES DUMMY ARGUMENTS OR DUMMY ARGUMENT
51 ! DECLARATIONS FOR RCONFIG ENTRIES. THEY ARE STILL PART OF STATE. ACCESS
52 ! TO THESE VARIABLES IS NOW THROUGH GRID STRUCTURE, AS MODIFIED BELOW.
53 ! AFFECTED VARIABLES: SIGMA, DT, NPHS, IDTAD, NRADS, NRADL, JULDAY,
54 ! JULYR, NUM_SOIL_LAYERS, NCNVC, ENSDIM, DY, AND SPEC_BDY_WIDTH.
55 ! JM, 20050819
56 !
57 !----------------------------
58 #include <nmm_dummy_decl.inc>
59 !----------------------------
60 !
61 !*** STRUCTURE THAT CONTAINS RUN-TIME CONFIGURATION (NAMELIST) DATA FOR DOMAIN
62 !
63 TYPE(GRID_CONFIG_REC_TYPE),INTENT(IN) :: CONFIG_FLAGS
64 #ifdef WRF_CHEM
65 INTEGER :: NUMGAS
66 #endif
67 !
68 !-----------------------------------------------------------------------
69 !
70 !*** LOCAL VARIABLES
71 !
72 !-----------------------------------------------------------------------
73 INTEGER :: IDS,IDE,JDS,JDE,KDS,KDE &
74 & ,IMS,IME,JMS,JME,KMS,KME &
75 & ,IPS,IPE,JPS,JPE,KPS,KPE &
76 & ,ITS,ITE,JTS,JTE,KTS,KTE
77 !
78 INTEGER :: I,ICLTEND,IDF,IJDE,IJDS,IRTN,J,JC,JDF,K,KDF,LB,N_MOIST &
79 & ,NTSD_current
80 INTEGER,SAVE :: NTSD_restart
81 !dusan INTEGER :: MPI_COMM_COMP,MYPE,MYPROC,NPES
82 INTEGER :: MYPROC
83 INTEGER :: KVH,NTSD_rad,RC
84 REAL :: WC, QI, QR, QW, FICE, FRAIN
85 INTEGER :: NUM_OZMIXM,NUM_AEROSOLC
86 !
87 CHARACTER*80 :: MESSAGE
88 !
89 ! For precip assimilation:
90 INTEGER :: ISTAT
91 REAL,ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: PPTDAT
92
93 ! For physics compatibility with other packages
94 REAL,ALLOCATABLE,DIMENSION(:,:,:) :: TTEN,QTEN
95 REAL,ALLOCATABLE,DIMENSION(:,:,:) :: RTHRATEN,RTHBLTEN,RQVBLTEN
96 !
97 REAL :: DT_INV,GPS
98 !
99 LOGICAL :: LAST_TIME
100 !
101 LOGICAL wrf_dm_on_monitor
102 EXTERNAL wrf_dm_on_monitor
103 !
104 !-----------------------------------------------------------------------
105 !*** TIMING VARIABLES
106 !-----------------------------------------------------------------------
107 real,save :: solve_tim,exch_tim,pdte_tim,adve_tim,vtoa_tim &
108 &, vadz_tim,hadz_tim,eps_tim,vad2_tim,had2_tim &
109 &, radiation_tim,rdtemp_tim,turbl_tim,cltend_tim &
110 &, cucnvc_tim,gsmdrive_tim,hdiff_tim,bocoh_tim &
111 &, pfdht_tim,ddamp_tim,bocov_tim,uv_htov_tim,sum_tim &
112 &, adjppt_tim
113 real,save :: exch_tim_max
114 real :: btim,btimx
115 real :: et_max,this_tim
116 integer :: n_print_time
117 !
118 #ifdef RSL
119 integer rsl_internal_milliclock
120 external rsl_internal_milliclock
121 # define timef rsl_internal_milliclock
122 #else
123 ! real*8,dimension(1) :: timef
124 real*8 :: timef
125 #endif
126 !-----------------------------------------------------------------------
127 !
128 #ifdef DEREF_KLUDGE
129 ! SEE http://www.mmm.ucar.edu/wrf/WG2/topics/deref_kludge.htm
130 INTEGER :: SM31,EM31,SM32,EM32,SM33,EM33
131 INTEGER :: SM31X,EM31X,SM32X,EM32X,SM33X,EM33X
132 INTEGER :: SM31Y,EM31Y,SM32Y,EM32Y,SM33Y,EM33Y
133 #endif
134 !
135 !-----------------------------------------------------------------------
136 !
137 ! LIMIT THE NUMBER OF ARGUMENTS IF COMPILED WITH -DLIMIT_ARGS BY COPYING
138 ! SCALAR (NON-ARRAY) ARGUMENTS OUT OF THE GRID DATA STRUCTURE INTO LOCALLY
139 ! DEFINED COPIES (DEFINED IN EM_DUMMY_DECL.INC, ABOVE, AS THEY ARE IF THEY
140 ! ARE ARGUMENTS). AN EQUIVALENT INCLUDE OF EM_SCALAR_DEREFS.INC APPEARS
141 ! AT THE END OF THE ROUTINE TO COPY BACK ANY CHNAGED NON-ARRAY VALUES.
142 ! THE DEFINITION OF COPY_IN OR COPY_OUT BEFORE THE INCLUDE DEFINES THE
143 ! DIRECTION OF THE COPY. NMM_SCALAR_DEREFS IS GENERATED FROM REGISTRY.
144 !
145 !-----------------------------------------------------------------------
146 #define COPY_IN
147 #include <nmm_scalar_derefs.inc>
148 !-----------------------------------------------------------------------
149 !
150 ! TRICK PROBLEMATIC COMPILERS INTO NOT PERFORMING COPY-IN/COPY-OUT BY ADDING
151 ! INDICES TO ARRAY ARGUMENTS IN THE CALL STATEMENTS IN THIS ROUTINE.
152 ! IT HAS THE EFFECT OF PASSING ONLY THE FIRST ELEMENT OF THE ARRAY, RATHER
153 ! THAN THE ENTIRE ARRAY. SEE:
154 ! http://www.mmm.ucar.edu/wrf/WG2/topics/deref_kludge.htm
155 !
156 !-----------------------------------------------------------------------
157 #include "deref_kludge.h"
158 !-----------------------------------------------------------------------
159 !
160 ! NEEDED BY SOME COMM LAYERS, E.G. RSL. IF NEEDED, nmm_data_calls.inc IS
161 ! GENERATED FROM THE REGISTRY. THE DEFINITION OF REGISTER_I1 ALLOWS
162 ! I1 DATA TO BE COMMUNICATED IN THIS ROUTINE IF NECESSARY.
163 !
164 !-----------------------------------------------------------------------
165 #ifdef DM_PARALLEL
166 # define REGISTER_I1
167 # include <nmm_data_calls.inc>
168 #endif
169 !-----------------------------------------------------------------------
170 !***********************************************************************
171 !-----------------------------------------------------------------------
172 CALL WRF_GET_MYPROC(MYPROC)
173 MYPE=MYPROC
174 !-----------------------------------------------------------------------
175 !
176 !*** OBTAIN DIMENSION INFORMATION STORED IN THE GRID DATA STRUCTURE.
177 !
178 CALL GET_IJK_FROM_GRID(GRID &
179 & ,IDS,IDE,JDS,JDE,KDS,KDE &
180 & ,IMS,IME,JMS,JME,KMS,KME &
181 & ,IPS,IPE,JPS,JPE,KPS,KPE )
182 !-----------------------------------------------------------------------
183 !
184 !*** COMPUTE THESE STARTING AND STOPPING LOCATIONS FOR EACH TILE AND
185 !*** NUMBER OF TILES.
186 !*** SEE: http://www.mmm.ucar.edu/wrf/WG2/topics/settiles
187 !
188 CALL SET_TILES(GRID,IDS,IDE,JDS,JDE,IPS,IPE,JPS,JPE)
189 !-----------------------------------------------------------------------
190 !
191 !*** TTEN, QTEN are used by GD convection scheme
192 !
193 ALLOCATE(TTEN(IMS:IME,KMS:KME,JMS:JME),STAT=ISTAT)
194 ALLOCATE(QTEN(IMS:IME,KMS:KME,JMS:JME),STAT=ISTAT)
195 ALLOCATE(RTHBLTEN(IMS:IME,KMS:KME,JMS:JME),STAT=ISTAT)
196 ALLOCATE(RQVBLTEN(IMS:IME,KMS:KME,JMS:JME),STAT=ISTAT)
197 ALLOCATE(RTHRATEN(IMS:IME,KMS:KME,JMS:JME),STAT=ISTAT)
198 !
199 #ifdef WRF_CHEM
200 NUMGAS = GET_LAST_GAS(CONFIG_FLAGS%CHEM_OPT)
201 #endif
202 DO J=JMS,JME
203 DO K=KMS,KME
204 DO I=IMS,IME
205 TTEN(I,K,J)=T(I,K,J)
206 QTEN(I,K,J)=Q(I,K,J)
207 ENDDO
208 ENDDO
209 ENDDO
210 !
211 GRID%SIGMA=1
212 HYDRO=.FALSE.
213 !
214 IJDS=MIN(IDS,JDS)
215 IJDE=MAX(IDE,JDE)
216 !
217 IDF=IDE-1
218 JDF=JDE-1
219 KDF=KDE-1
220 !
221 !-----------------------------------------------------------------------
222 !
223 !*** FOR NOW SET CONTROLS FOR TILES TO PATCHES
224 !
225 !-----------------------------------------------------------------------
226 ITS=IPS
227 ITE=MIN(IPE,IDF)
228 JTS=JPS
229 JTE=MIN(JPE,JDF)
230 KTS=KPS
231 KTE=MIN(KPE,KDF)
232 if(ntsd==0)then
233 write(0,*)' its=',its,' ite=',ite
234 write(0,*)' jts=',jts,' jte=',jte
235 write(0,*)' kts=',kts,' kte=',kte
236 endif
237 !-----------------------------------------------------------------------
238 !*** SET TIMING VARIABLES TO ZERO AT START OF FORECAST.
239 !-----------------------------------------------------------------------
240 if(ntsd==0)then
241 solve_tim=0.
242 exch_tim=0.
243 pdte_tim=0.
244 adve_tim=0.
245 vtoa_tim=0.
246 vadz_tim=0.
247 hadz_tim=0.
248 eps_tim=0.
249 vad2_tim=0.
250 had2_tim=0.
251 radiation_tim=0.
252 rdtemp_tim=0.
253 turbl_tim=0.
254 cltend_tim=0.
255 cucnvc_tim=0.
256 gsmdrive_tim=0.
257 hdiff_tim=0.
258 bocoh_tim=0.
259 pfdht_tim=0.
260 ddamp_tim=0.
261 bocov_tim=0.
262 uv_htov_tim=0.
263 exch_tim_max=0.
264 adjppt_tim=0.
265 endif
266 !-----------------------------------------------------------------------
267 N_MOIST=NUM_MOIST
268 !
269 DO J=MYJS_P4,MYJE_P4
270 IHEG(J)=MOD(J+1,2)
271 IHWG(J)=IHEG(J)-1
272 IVEG(J)=MOD(J,2)
273 IVWG(J)=IVEG(J)-1
274 ENDDO
275
276 DO J=MYJS_P4,MYJE_P4
277 IVW(J)=IVWG(J)
278 IVE(J)=IVEG(J)
279 IHE(J)=IHEG(J)
280 IHW(J)=IHWG(J)
281 ENDDO
282 !
283 !*** LATERAL POINTS IN THE BOUNDARY ARRAYS
284 !
285 LB=2*(IDF-IDS+1)+(JDF-JDS+1)-3
286 !
287 !*** APPROXIMATE GRIDPOINT SPACING (METERS)
288 !
289 JC=JMS+(JME-JMS)/2
290 GPS=SQRT(DX_NMM(IMS,JC)**2+DY_NMM**2)
291 !
292 !*** TIMESTEPS PER HOUR
293 !
294 TSPH=3600./GRID%DT
295 !
296 n_print_time=nint(3600./grid%dt) ! Print stats once per hour
297 !-----------------------------------------------------------------------
298 !
299 NBOCO=0
300 !-----------------------------------------------------------------------
301 !***
302 !*** THE MAIN TIME INTEGRATION LOOP
303 !***
304 !-----------------------------------------------------------------------
305 !
306 !*** NTSD IS THE TIMESTEP COUNTER (Number of Time Steps Done)
307 !
308 !-----------------------------------------------------------------------
309 !***
310 !*** ADVANCE_count STARTS AT ZERO FOR ALL RUNS (REGULAR AND RESTART).
311 !***
312 !-----------------------------------------------------------------------
313 !
314 CALL DOMAIN_CLOCK_GET(GRID,ADVANCEcOUNT=NTSD_current)
315 !
316 IF(NTSD_current==0)THEN
317 IF(GRID%RESTART.AND.GRID%TSTART>0.)THEN
318 IHRST=NSTART_HOUR
319 NTSD_restart=NTSD+1
320 ELSE
321 IHRST=GRID%GMT
322 NSTART_HOUR=IHRST
323 NTSD_restart=0
324 ENDIF
325 ENDIF
326 !
327 NTSD=NTSD_restart+NTSD_current
328 LAST_TIME=domain_last_time_step(GRID)
329 !
330 !-----------------------------------------------------------------------
331 !
332 IF(WRF_DM_ON_MONITOR() )THEN
333 WRITE(MESSAGE,125)NTSD,NTSD*GRID%DT/3600.
334 125 FORMAT(' SOLVE_NMM: TIMESTEP IS ',I5,' TIME IS ',F7.3,' HOURS')
335 CALL WRF_MESSAGE(TRIM(MESSAGE))
336 ENDIF
337 !
338 !-----------------------------------------------------------------------
339 !
340 CALL WRF_GET_DM_COMMUNICATOR(MPI_COMM_COMP)
341 CALL WRF_GET_NPROC(NPES)
342 !
343 #if (NMM_NEST == 1)
344 !-----------------------------------------------------------------------------
345 !*** PATCHING NESTED BOUNDARIES.
346 !-----------------------------------------------------------------------------
347 !
348 CALL wrf_debug ( 100 , 'nmm: in patch' )
349
350 btimx=timef()
351 #ifdef DM_PARALLEL
352 # include "HALO_NMM_ZZ.inc"
353 #endif
354
355 IF(GRID%ID/=1)THEN
356 !
357 CALL NESTBC_PATCH (PD_B,T_B,Q_B,U_B,V_B,Q2_B,CWM_B &
358 ,PD_BT,T_BT,Q_BT,U_BT,V_BT,Q2_BT,CWM_BT &
359 ,PDNEST_B,TNEST_B,QNEST_B,UNEST_B,VNEST_B,Q2NEST_B,CWMNEST_B &
360 ,PDNEST_BT,TNEST_BT,QNEST_BT,UNEST_BT,VNEST_BT,Q2NEST_BT,CWMNEST_BT &
361 ,IJDS,IJDE,GRID%SPEC_BDY_WIDTH &
362 ,IDS,IDF,JDS,JDF,KDS,KDE &
363 ,IMS,IME,JMS,JME,KMS,KME &
364 ,ITS,ITE,JTS,JTE,KTS,KTE )
365 CALL wrf_debug ( 100 , 'nmm: out of patch' )
366 !
367 #ifdef MOVE_NESTS
368
369 IF(GRID%ID/=1.AND.MOD(NTSD,1)==0.AND.GRID%NUM_MOVES==-99)THEN
370 XLOC_1=(IDE-1)/2 ! This maneuvers the storm to the center of the nest quickly
371 YLOC_1=(JDE-1)/2 ! This maneuvers the storm to the center of the nest quickly
372 ENDIF
373 #endif
374
375 ENDIF
376 #endif
377 !
378 !-----------------------------------------------------------------------
379 !*** ALLOCATE PPTDAT ARRAY (PRECIP ASSIM):
380 !-----------------------------------------------------------------------
381 !
382 IF(GRID%PCPFLG.AND..NOT.ALLOCATED(PPTDAT))THEN
383 ALLOCATE(PPTDAT(IMS:IME,JMS:JME,3),STAT=ISTAT)
384 ENDIF
385 !
386 !-----------------------------------------------------------------------
387 !***
388 !*** Call READPCP to
389 !*** 1) READ IN PRECIPITATION FOR HOURS 1, 2 and 3;
390 !*** 2) Initialize DDATA to 999. (this is the amount
391 !*** of input precip allocated to each physics time step
392 !*** in ADJPPT; TURBL/SURFCE, which uses DDATA, is called
393 !*** before ADJPPT)
394 !*** 3) Initialize LSPA to zero
395 !***
396 !-----------------------------------------------------------------------
397 IF (NTSD==0) THEN
398 IF (GRID%PCPFLG) THEN
399 CALL READPCP(PPTDAT,DDATA,LSPA &
400 & ,IDS,IDE,JDS,JDE,KDS,KDE &
401 & ,IMS,IME,JMS,JME,KMS,KME &
402 & ,ITS,ITE,JTS,JTE,KTS,KTE)
403 ENDIF
404 ENDIF
405 !-----------------------------------------------------------------------
406 !
407 btim=timef()
408 !
409 !-----------------------------------------------------------------------
410 !*** ZERO OUT ACCUMULATED QUANTITIES WHEN NEEDED.
411 !-----------------------------------------------------------------------
412 !
413 CALL BUCKETS(NTSD,NPREC,NSRFC,NRDSW,NRDLW &
414 & ,GRID%RESTART,GRID%TSTART &
415 & ,NCLOD,NHEAT,GRID%NPHS,TSPH &
416 & ,ACPREC,CUPREC,ACSNOW,ACSNOM,SSROFF,BGROFF &
417 & ,SFCEVP,POTEVP,SFCSHX,SFCLHX,SUBSHX,SNOPCX &
418 & ,SFCUVX,POTFLX &
419 & ,ARDSW,ASWIN,ASWOUT,ASWTOA &
420 & ,ARDLW,ALWIN,ALWOUT,ALWTOA &
421 & ,ACFRST,NCFRST,ACFRCV,NCFRCV &
422 & ,AVCNVC,AVRAIN,TCUCN,TRAIN &
423 & ,ASRFC &
424 & ,T,TLMAX,TLMIN &
425 & ,IDS,IDE,JDS,JDE,KDS,KDE &
426 & ,IMS,IME,JMS,JME,KMS,KME &
427 & ,ITS,ITE,JTS,JTE,KTS,KTE)
428 !-----------------------------------------------------------------------
429 !
430 IF(NTSD==0)THEN
431 FIRST=.TRUE.
432 ! call hpm_init()
433 btimx=timef()
434 !
435 !-----------------------------------------------------------------------
436 #ifdef DM_PARALLEL
437 # include "HALO_NMM_A.inc"
438 #endif
439 !
440 #ifdef DM_PARALLEL
441 IF(CONFIG_FLAGS%MP_PHYSICS/=ETAMPNEW)THEN
442 # include "HALO_NMM_A_3.inc"
443 ENDIF
444 #endif
445 !
446 !*** Only for chemistry:
447 !
448 #ifdef WRF_CHEM
449 #ifdef DM_PARALLEL
450 # include "HALO_NMM_A_2.inc"
451 #endif
452 #endif
453 !
454 !-----------------------------------------------------------------------
455 !*** USE THE FOLLOWING VARIABLES TO KEEP TRACK OF EXCHANGE TIMES.
456 !-----------------------------------------------------------------------
457 exch_tim=exch_tim+timef()-btimx
458 ! this_tim=timef()-btimx
459 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
460 ! & ,mpi_comm_comp,irtn)
461 ! exch_tim_max=exch_tim_max+et_max
462 !-----------------------------------------------------------------------
463 !
464 GO TO 2003
465 ENDIF
466 !
467 !-----------------------------------------------------------------------
468 !-----------------------------------------------------------------------
469 2000 CONTINUE
470 !-----------------------------------------------------------------------
471 !-----------------------------------------------------------------------
472 !
473 !-----------------------------------------------------------------------
474 !*** PRESSURE TENDENCY, SIGMA DOT, VERTICAL PART OF OMEGA-ALPHA
475 !-----------------------------------------------------------------------
476 !
477 btimx=timef()
478 !-----------------
479 #ifdef DM_PARALLEL
480 # include "HALO_NMM_D.inc"
481 #endif
482 !-----------------
483 exch_tim=exch_tim+timef()-btimx
484 ! this_tim=timef()-btimx
485 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
486 ! & ,mpi_comm_comp,irtn)
487 ! exch_tim_max=exch_tim_max+et_max
488 !
489 btimx=timef()
490 !
491 CALL PDTE( &
492 #ifdef DM_PARALLEL
493 & GRID, &
494 #endif
495 & NTSD,GRID%DT,PT,ETA2,RES,HYDRO &
496 & ,HTM,HBM2 &
497 & ,PD,PDSL,PDSLO &
498 & ,PETDT,DIV,PSDT &
499 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
500 & ,IDS,IDF,JDS,JDF,KDS,KDE &
501 & ,IMS,IME,JMS,JME,KMS,KME &
502 & ,ITS,ITE,JTS,JTE,KTS,KTE)
503
504 pdte_tim=pdte_tim+timef()-btimx
505 !
506 !-----------------------------------------------------------------------
507 !*** ADVECTION OF T, U, AND V
508 !-----------------------------------------------------------------------
509 !
510 btimx=timef()
511 !-----------------
512 #ifdef DM_PARALLEL
513 # include "HALO_NMM_F.inc"
514 # include "HALO_NMM_F1.inc"
515 #endif
516 !-----------------
517 exch_tim=exch_tim+timef()-btimx
518 ! this_tim=timef()-btimx
519 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
520 ! & ,mpi_comm_comp,irtn)
521 ! exch_tim_max=exch_tim_max+et_max
522 !
523 btimx=timef()
524 !
525 CALL ADVE(NTSD,GRID%DT,DETA1,DETA2,PDTOP &
526 & ,CURV,F,FAD,F4D,EM_LOC,EMT_LOC,EN,ENT,DX_NMM,DY_NMM &
527 & ,HTM,HBM2,VTM,VBM2,LMH,LMV &
528 & ,T,U,V,PDSLO,TOLD,UOLD,VOLD &
529 & ,PETDT,UPSTRM &
530 & ,FEW,FNS,FNE,FSE &
531 & ,ADT,ADU,ADV &
532 & ,N_IUP_H,N_IUP_V &
533 & ,N_IUP_ADH,N_IUP_ADV &
534 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
535 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
536 & ,IDS,IDF,JDS,JDF,KDS,KDE &
537 & ,IMS,IME,JMS,JME,KMS,KME &
538 & ,ITS,ITE,JTS,JTE,KTS,KTE)
539 !
540 adve_tim=adve_tim+timef()-btimx
541 !
542 !-----------------------------------------------------------------------
543 !*** PRESSURE TENDENCY, ETA/SIGMADOT, VERTICAL PART OF OMEGA-ALPHA TERM
544 !-----------------------------------------------------------------------
545 !
546 btimx=timef()
547 !
548 CALL VTOA( &
549 #ifdef DM_PARALLEL
550 & GRID, &
551 #endif
552 & NTSD,GRID%DT,PT,ETA2 &
553 & ,HTM,HBM2,EF4T &
554 & ,T,DWDT,RTOP,OMGALF &
555 & ,PINT,DIV,PSDT,RES &
556 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
557 & ,IDS,IDF,JDS,JDF,KDS,KDE &
558 & ,IMS,IME,JMS,JME,KMS,KME &
559 & ,ITS,ITE,JTS,JTE,KTS,KTE)
560 !
561 vtoa_tim=vtoa_tim+timef()-btimx
562 !
563 !-----------------------------------------------------------------------
564 !*** VERTICAL ADVECTION OF HEIGHT
565 !-----------------------------------------------------------------------
566 !
567 btimx=timef()
568 !
569 CALL VADZ(NTSD,GRID%DT,FIS,GRID%SIGMA,DFL,HTM,HBM2 &
570 & ,DETA1,DETA2,PDTOP &
571 & ,PINT,PDSL,PDSLO,PETDT &
572 & ,RTOP,T,Q,CWM,Z,W,DWDT,PDWDT &
573 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
574 & ,IDS,IDF,JDS,JDF,KDS,KDE &
575 & ,IMS,IME,JMS,JME,KMS,KME &
576 & ,ITS,ITE,JTS,JTE,KTS,KTE)
577
578 vadz_tim=vadz_tim+timef()-btimx
579 !
580 !-----------------------------------------------------------------------
581 !*** HORIZONTAL ADVECTION OF HEIGHT
582 !-----------------------------------------------------------------------
583 !
584 btimx=timef()
585 !-----------------
586 #ifdef DM_PARALLEL
587 # include "HALO_NMM_G.inc"
588 #endif
589 !-----------------
590 exch_tim=exch_tim+timef()-btimx
591 ! this_tim=timef()-btimx
592 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
593 ! & ,mpi_comm_comp,irtn)
594 ! exch_tim_max=exch_tim_max+et_max
595 !
596 btimx=timef()
597 !
598 CALL HADZ(NTSD,GRID%DT,HYDRO,HTM,HBM2,DETA1,DETA2,PDTOP &
599 & ,DX_NMM,DY_NMM,FAD &
600 & ,FEW,FNS,FNE,FSE &
601 & ,PDSL,U,V,W,Z &
602 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
603 & ,IDS,IDF,JDS,JDF,KDS,KDE &
604 & ,IMS,IME,JMS,JME,KMS,KME &
605 & ,ITS,ITE,JTS,JTE,KTS,KTE)
606 !
607 hadz_tim=hadz_tim+timef()-btimx
608 !
609 !-----------------------------------------------------------------------
610 !*** ADVECTION OF W
611 !-----------------------------------------------------------------------
612 !
613 btimx=timef()
614 !-----------------
615 #ifdef DM_PARALLEL
616 # include "HALO_NMM_H.inc"
617 #endif
618 !-----------------
619 exch_tim=exch_tim+timef()-btimx
620 ! this_tim=timef()-btimx
621 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
622 ! & ,mpi_comm_comp,irtn)
623 ! exch_tim_max=exch_tim_max+et_max
624 !
625 btimx=timef()
626 !
627 CALL EPS(NTSD,GRID%DT,HYDRO,DX_NMM,DY_NMM,FAD &
628 & ,DETA1,DETA2,PDTOP,PT &
629 & ,HTM,HBM2,HBM3,LMH &
630 & ,PDSL,PDSLO,PINT,RTOP,PETDT,PDWDT &
631 & ,DWDT,DWDTMN,DWDTMX &
632 & ,FNS,FEW,FNE,FSE &
633 & ,T,U,V,W,Q,CWM &
634 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
635 & ,IDS,IDF,JDS,JDF,KDS,KDE &
636 & ,IMS,IME,JMS,JME,KMS,KME &
637 & ,ITS,ITE,JTS,JTE,KTS,KTE)
638 !
639 eps_tim=eps_tim+timef()-btimx
640 !
641 !-----------------------------------------------------------------------
642 !*** VERTICAL ADVECTION OF Q, TKE, AND CLOUD WATER
643 !-----------------------------------------------------------------------
644 !
645 IF(MOD(NTSD,GRID%IDTAD)==0)THEN
646 btimx=timef()
647 !
648 vad2_micro_check: IF(CONFIG_FLAGS%MP_PHYSICS==ETAMPNEW)THEN
649 CALL VAD2(NTSD,GRID%DT,GRID%IDTAD,DX_NMM,DY_NMM &
650 & ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
651 & ,HBM2,LMH &
652 & ,Q,Q2,CWM,PETDT &
653 & ,N_IUP_H,N_IUP_V &
654 & ,N_IUP_ADH,N_IUP_ADV &
655 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
656 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
657 & ,IDS,IDF,JDS,JDF,KDS,KDE &
658 & ,IMS,IME,JMS,JME,KMS,KME &
659 & ,ITS,ITE,JTS,JTE,KTS,KTE)
660 !
661 ELSE vad2_micro_check
662 CALL VAD2_SCAL(NTSD,GRID%DT,GRID%IDTAD,DX_NMM,DY_NMM &
663 & ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
664 & ,HBM2,LMH &
665 & ,Q2,PETDT &
666 & ,N_IUP_H,N_IUP_V &
667 & ,N_IUP_ADH,N_IUP_ADV &
668 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
669 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
670 & ,1,1 &
671 & ,IDS,IDF,JDS,JDF,KDS,KDE &
672 & ,IMS,IME,JMS,JME,KMS,KME &
673 & ,ITS,ITE,JTS,JTE,KTS,KTE)
674
675 CALL VAD2_SCAL(NTSD,GRID%DT,GRID%IDTAD,DX_NMM,DY_NMM &
676 & ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
677 & ,HBM2,LMH &
678 & ,MOIST,PETDT &
679 & ,N_IUP_H,N_IUP_V &
680 & ,N_IUP_ADH,N_IUP_ADV &
681 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
682 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
683 & ,NUM_MOIST,2 &
684 & ,IDS,IDF,JDS,JDF,KDS,KDE &
685 & ,IMS,IME,JMS,JME,KMS,KME &
686 & ,ITS,ITE,JTS,JTE,KTS,KTE)
687 !
688 CALL VAD2_SCAL(NTSD,GRID%DT,GRID%IDTAD,DX_NMM,DY_NMM &
689 & ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
690 & ,HBM2,LMH &
691 & ,SCALAR,PETDT &
692 & ,N_IUP_H,N_IUP_V &
693 & ,N_IUP_ADH,N_IUP_ADV &
694 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
695 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
696 & ,NUM_SCALAR,2 &
697 & ,IDS,IDF,JDS,JDF,KDS,KDE &
698 & ,IMS,IME,JMS,JME,KMS,KME &
699 & ,ITS,ITE,JTS,JTE,KTS,KTE)
700 !
701
702 DO J=MYJS,MYJE
703 DO K=KTS,KTE
704 DO i=MYIS,MYIE
705 Q(I,K,J)=MOIST(I,K,J,P_QV)/(1.+MOIST(I,K,J,P_QV))
706 ENDDO
707 ENDDO
708 ENDDO
709 !
710 ENDIF vad2_micro_check
711 !
712 vad2_tim=vad2_tim+timef()-btimx
713 !
714 ENDIF
715 !
716 !-----------------------------------------------------------------------
717 !*** VERTICAL ADVECTION OF CHEMISTRY
718 !-----------------------------------------------------------------------
719 !
720 #ifdef WRF_CHEM
721 IF(MOD(NTSD,GRID%IDTAD)==0)THEN
722 #ifdef IBM
723 btimx=timef()
724 #endif
725 !
726 CALL VAD2_SCAL(NTSD,GRID%DT,GRID%IDTAD,DX_NMM,DY_NMM &
727 & ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
728 & ,HBM2,LMH &
729 & ,CHEM,PETDT &
730 & ,N_IUP_H,N_IUP_V &
731 & ,N_IUP_ADH,N_IUP_ADV &
732 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
733 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
734 & ,NUM_CHEM,1 &
735 & ,IDS,IDF,JDS,JDF,KDS,KDE &
736 & ,IMS,IME,JMS,JME,KMS,KME &
737 & ,ITS,ITE,JTS,JTE,KTS,KTE)
738 !
739 ENDIF
740 #endif
741 !
742 !-----------------------------------------------------------------------
743 !*** HORIZONTAL ADVECTION OF Q, TKE, AND CLOUD WATER
744 !-----------------------------------------------------------------------
745 !
746 IF(MOD(NTSD,GRID%IDTAD)==0)THEN
747 btimx=timef()
748 !-----------------
749 #ifdef DM_PARALLEL
750 # include "HALO_NMM_I.inc"
751 #endif
752 !
753 #ifdef DM_PARALLEL
754 IF(CONFIG_FLAGS%MP_PHYSICS/=ETAMPNEW)THEN
755 # include "HALO_NMM_I_3.inc"
756 ENDIF
757 #endif
758 !
759 !-----------------
760 exch_tim=exch_tim+timef()-btimx
761 ! this_tim=timef()-btimx
762 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
763 ! & ,mpi_comm_comp,irtn)
764 ! exch_tim_max=exch_tim_max+et_max
765 !
766 btimx=timef()
767 !
768 !-----------------------------------------------------------------------
769 had2_micro_check: IF(CONFIG_FLAGS%MP_PHYSICS==ETAMPNEW)THEN
770 !-----------------------------------------------------------------------
771 !
772 CALL HAD2( &
773 #if defined(DM_PARALLEL)
774 & GRID%DOMDESC, &
775 #endif
776 & NTSD,GRID%DT,GRID%IDTAD,DX_NMM,DY_NMM &
777 & ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
778 & ,HTM,HBM2,HBM3,LMH &
779 & ,Q,Q2,CWM,U,V,Z,HYDRO &
780 & ,N_IUP_H,N_IUP_V &
781 & ,N_IUP_ADH,N_IUP_ADV &
782 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
783 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
784 & ,IDS,IDF,JDS,JDF,KDS,KDE &
785 & ,IMS,IME,JMS,JME,KMS,KME &
786 & ,ITS,ITE,JTS,JTE,KTS,KTE)
787 !
788 !*** UPDATE MOIST ARRAY
789 !
790 DO J=MYJS,MYJE
791 DO K=KTS,KTE
792 DO I=MYIS,MYIE
793 MOIST(I,K,J,P_QV)=Q(I,K,J)/(1.-Q(I,K,J))
794 WC = CWM(I,K,J)
795 QI = 0.
796 QR = 0.
797 QW = 0.
798 FICE=F_ICE(I,K,J)
799 FRAIN=F_RAIN(I,K,J)
800 !
801 IF(FICE>=1.)THEN
802 QI=WC
803 ELSEIF(FICE<=0.)THEN
804 QW=WC
805 ELSE
806 QI=FICE*WC
807 QW=WC-QI
808 ENDIF
809 !
810 IF(QW>0..AND.FRAIN>0.)THEN
811 IF(FRAIN>=1.)THEN
812 QR=QW
813 QW=0.
814 ELSE
815 QR=FRAIN*QW
816 QW=QW-QR
817 ENDIF
818 ENDIF
819 !
820 MOIST(I,K,J,P_QC)=QW
821 MOIST(I,K,J,P_QR)=QR
822 MOIST(I,K,J,P_QI)=0.
823 MOIST(I,K,J,P_QS)=QI
824 MOIST(I,K,J,P_QG)=0.
825 ENDDO
826 ENDDO
827 ENDDO
828 !
829 !-----------------------------------------------------------------------
830 ELSE had2_micro_check
831 !-----------------------------------------------------------------------
832 !
833 CALL HAD2_SCAL( &
834 #if defined(DM_PARALLEL)
835 & GRID%DOMDESC, &
836 #endif
837 & NTSD,GRID%DT,GRID%IDTAD,DX_NMM,DY_NMM &
838 & ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
839 & ,HTM,HBM2,HBM3,LMH &
840 & ,Q2,U,V,Z,HYDRO &
841 & ,N_IUP_H,N_IUP_V &
842 & ,N_IUP_ADH,N_IUP_ADV &
843 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
844 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
845 & ,1,1 &
846 & ,IDS,IDF,JDS,JDF,KDS,KDE &
847 & ,IMS,IME,JMS,JME,KMS,KME &
848 & ,ITS,ITE,JTS,JTE,KTS,KTE)
849 !
850 CALL HAD2_SCAL( &
851 #if defined(DM_PARALLEL)
852 & GRID%DOMDESC, &
853 #endif
854 & NTSD,GRID%DT,GRID%IDTAD,DX_NMM,DY_NMM &
855 & ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
856 & ,HTM,HBM2,HBM3,LMH &
857 & ,MOIST,U,V,Z,HYDRO &
858 & ,N_IUP_H,N_IUP_V &
859 & ,N_IUP_ADH,N_IUP_ADV &
860 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
861 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
862 & ,NUM_MOIST,2 &
863 & ,IDS,IDF,JDS,JDF,KDS,KDE &
864 & ,IMS,IME,JMS,JME,KMS,KME &
865 & ,ITS,ITE,JTS,JTE,KTS,KTE)
866 !
867 CALL HAD2_SCAL( &
868 #if defined(DM_PARALLEL)
869 & GRID%DOMDESC, &
870 #endif
871 & NTSD,GRID%DT,GRID%IDTAD,DX_NMM,DY_NMM &
872 & ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
873 & ,HTM,HBM2,HBM3,LMH &
874 & ,SCALAR,U,V,Z,HYDRO &
875 & ,N_IUP_H,N_IUP_V &
876 & ,N_IUP_ADH,N_IUP_ADV &
877 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
878 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
879 & ,NUM_SCALAR,2 &
880 & ,IDS,IDF,JDS,JDF,KDS,KDE &
881 & ,IMS,IME,JMS,JME,KMS,KME &
882 & ,ITS,ITE,JTS,JTE,KTS,KTE)
883 !
884 DO J=MYJS,MYJE
885 DO K=KTS,KTE
886 DO I=MYIS,MYIE
887 Q(I,K,J)=MOIST(I,K,J,P_QV)/(1.+MOIST(I,K,J,P_QV))
888 ENDDO
889 ENDDO
890 ENDDO
891 !
892 !-----------------------------------------------------------------------
893 ENDIF had2_micro_check
894 !-----------------------------------------------------------------------
895 !
896 had2_tim=had2_tim+timef()-btimx
897 ENDIF
898 !
899 !-----------------------------------------------------------------------
900 !*** HORIZONTAL ADVECTION OF CHEMISTRY
901 !-----------------------------------------------------------------------
902 !
903 #ifdef WRF_CHEM
904 IF(MOD(NTSD,GRID%IDTAD)==0)THEN
905 btimx=timef()
906 #ifdef DM_PARALLEL
907 # include "HALO_NMM_I_2.inc"
908 #endif
909 exch_tim=exch_tim+timef()-btimx
910 ! this_tim=timef()-btimx
911 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
912 ! & ,mpi_comm_comp,irtn)
913 ! exch_tim_max=exch_tim_max+et_max
914 !
915 btimx=timef()
916 !
917 CALL HAD2_SCAL( &
918 #if defined(DM_PARALLEL)
919 & GRID%DOMDESC, &
920 #endif
921 & NTSD,GRID%DT,GRID%IDTAD,DX_NMM,DY_NMM &
922 & ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
923 & ,HTM,HBM2,HBM3,LMH &
924 & ,CHEM,U,V,Z,HYDRO &
925 & ,N_IUP_H,N_IUP_V &
926 & ,N_IUP_ADH,N_IUP_ADV &
927 & ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
928 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
929 & ,NUM_CHEM,1 &
930 & ,IDS,IDF,JDS,JDF,KDS,KDE &
931 & ,IMS,IME,JMS,JME,KMS,KME &
932 & ,ITS,ITE,JTS,JTE,KTS,KTE)
933
934 ENDIF
935 #endif
936 !
937 !----------------------------------------------------------------------
938 !*** RADIATION
939 !----------------------------------------------------------------------
940 !
941 !*** When allocating CAM radiation 4d arrays (ozmixm, aerosolc),
942 !*** the following two scalars are not needed.
943 !
944 NUM_OZMIXM=1
945 NUM_AEROSOLC=1
946 !
947 IF(NTSD<=0)THEN
948 NTSD_rad=NTSD
949 ELSE
950 !
951 !*** Call radiation just BEFORE the top of the hour
952 !*** so that updated fields are written to history files.
953 !
954 NTSD_rad=NTSD+1
955 ENDIF
956 !
957 IF(MOD(NTSD_rad,GRID%NRADS)==0.OR. &
958 & MOD(NTSD_rad,GRID%NRADL)==0)THEN
959 btimx=timef()
960 !
961 CALL RADIATION(NTSD_rad,GRID%DT,GRID%JULDAY,GRID%JULYR &
962 & ,GRID%XTIME,GRID%JULIAN &
963 & ,IHRST,GRID%NPHS &
964 & ,GLAT,GLON,GRID%NRADS,GRID%NRADL &
965 & ,DETA1,DETA2,AETA1,AETA2,ETA1,ETA2,PDTOP,PT &
966 & ,PD,RES,PINT,T,Q,MOIST,THS,ALBEDO,EPSR &
967 & ,F_ICE,F_RAIN &
968 #ifdef WRF_CHEM
969 & ,GD_CLOUD,GD_CLOUD2 &
970 #endif
971 & ,SM,HBM2,LMH,CLDFRA,N_MOIST,RESTRT &
972 & ,RLWTT,RSWTT,RLWIN,RSWIN,RSWINC,RSWOUT &
973 & ,RLWTOA,RSWTOA,CZMEAN &
974 & ,CFRACL,CFRACM,CFRACH,SIGT4 &
975 & ,ACFRST,NCFRST,ACFRCV,NCFRCV &
976 & ,CUPPT,VEGFRC,SNO,HTOP,HBOT &
977 & ,Z,SICE,NUM_AEROSOLC,NUM_OZMIXM &
978 & ,GRID,CONFIG_FLAGS &
979 & ,RTHRATEN &
980 #ifdef WRF_CHEM
981 & ,PM2_5_DRY, PM2_5_WATER, PM2_5_DRY_EC &
982 & ,TAUAER1, TAUAER2, TAUAER3, TAUAER4 &
983 & ,GAER1, GAER2, GAER3, GAER4 &
984 & ,WAER1, WAER2, WAER3, WAER4 &
985 #endif
986 & ,IDS,IDF,JDS,JDF,KDS,KDE &
987 & ,IMS,IME,JMS,JME,KMS,KME &
988 & ,ITS,ITE,JTS,JTE,KTS,KTE)
989 !
990 DO J=JMS,JME
991 DO I=IMS,IME
992 GSW(I,J)=RSWIN(I,J)-RSWOUT(I,J)
993 ENDDO
994 ENDDO
995 !
996 ! *** NOTE ***
997 ! RLWIN/RSWIN - downward longwave/shortwave at the surface (formerly TOTLWDN/TOTSWDN)
998 ! RSWINC - CLEAR-SKY downward shortwave at the surface (new for AQ)
999 ! *** NOTE ***
1000 !
1001 radiation_tim=radiation_tim+timef()-btimx
1002 ENDIF
1003 !
1004 !----------------------------------------------------------------------
1005 !*** APPLY TEMPERATURE TENDENCY DUE TO RADIATION
1006 !----------------------------------------------------------------------
1007 !
1008 btimx=timef()
1009 !
1010 CALL RDTEMP(NTSD,GRID%DT,GRID%JULDAY,GRID%JULYR,IHRST,GLAT,GLON &
1011 & ,CZEN,CZMEAN,T,RSWTT,RLWTT,HTM,HBM2 &
1012 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1013 & ,IMS,IME,JMS,JME,KMS,KME &
1014 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1015 !
1016 rdtemp_tim=rdtemp_tim+timef()-btimx
1017 !
1018 !----------------------------------------------------------------------
1019 !*** TURBULENT PROCESSES
1020 !----------------------------------------------------------------------
1021 !
1022 IF(MOD(NTSD,GRID%NPHS)==0)THEN
1023 !
1024 btimx=timef()
1025 !
1026 CALL TURBL(NTSD,GRID%DT,GRID%NPHS,RESTRT &
1027 & ,N_MOIST,GRID%NUM_SOIL_LAYERS,SLDPTH,DZSOIL &
1028 & ,DETA1,DETA2,AETA1,AETA2,ETA1,ETA2,PDTOP,PT &
1029 ! & ,SM,LMH,HTM,VTM,HBM2,VBM2,DX_NMM,DFL &
1030 & ,SM,LMH,HTM,VTM,HBM2,VBM2,DX_NMM,DFRLG &
1031 & ,CZEN,CZMEAN,SIGT4,RLWIN,RSWIN,RADOT &
1032 & ,PD,RES,PINT,T,Q,CWM,F_ICE,F_RAIN,SR &
1033 & ,Q2,U,V,THS,NMM_TSK,SST,PREC,SNO,ZERO_3D &
1034 & ,FIS,Z0,Z0BASE,USTAR,PBLH,LPBL,EL_MYJ &
1035 & ,MOIST,RMOL &
1036 & ,EXCH_H,AKHS,AKMS,AKHS_OUT,AKMS_OUT &
1037 & ,THZ0,QZ0,UZ0,VZ0,QSH,MAVAIL &
1038 & ,STC,SMC,CMC,SMSTAV,SMSTOT,SSROFF,BGROFF &
1039 & ,IVGTYP,ISLTYP,VEGFRC,SHDMIN,SHDMAX,GRNFLX &
1040 & ,SFCEXC,ACSNOW,ACSNOM,SNOPCX,SICE,TG,SOILTB &
1041 & ,ALBASE,MXSNAL,ALBEDO,SH2O,SI,EPSR &
1042 & ,U10,V10,TH10,Q10,TSHLTR,QSHLTR,PSHLTR &
1043 & ,T2,QSG,QVG,QCG,SOILT1,TSNAV,SMFR3D,KEEPFR3DFLAG &
1044 & ,TWBS,QWBS,SFCSHX,SFCLHX,SFCEVP &
1045 & ,POTEVP,POTFLX,SUBSHX &
1046 & ,APHTIM,ARDSW,ARDLW,ASRFC &
1047 & ,RSWOUT,RSWTOA,RLWTOA &
1048 & ,ASWIN,ASWOUT,ASWTOA,ALWIN,ALWOUT,ALWTOA &
1049 & ,UZ0H,VZ0H,DUDT,DVDT &
1050 & ,RTHBLTEN,RQVBLTEN &
1051 & ,GRID%PCPFLG,DDATA &
1052 & ,GRID,CONFIG_FLAGS &
1053 & ,IHE,IHW,IVE,IVW &
1054 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1055 & ,IMS,IME,JMS,JME,KMS,KME &
1056 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1057 !
1058 ! *** NOTE ***
1059 ! RLWIN/RSWIN - downward longwave/shortwave at the surface
1060 ! *** NOTE ***
1061 !
1062 turbl_tim=turbl_tim+timef()-btimx
1063 !
1064 btimx=timef()
1065 !-----------------
1066 #ifdef DM_PARALLEL
1067 # include "HALO_NMM_TURBL_A.inc"
1068 #endif
1069 !-----------------
1070 #ifdef DM_PARALLEL
1071 # include "HALO_NMM_TURBL_B.inc"
1072 #endif
1073 !-----------------
1074 exch_tim=exch_tim+timef()-btimx
1075 ! this_tim=timef()-btimx
1076 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
1077 ! & ,mpi_comm_comp,irtn)
1078 ! exch_tim_max=exch_tim_max+et_max
1079 !
1080 !*** INTERPOLATE WINDS FROM H POINTS BACK TO V POINTS.
1081 !
1082 btimx=timef()
1083 CALL UV_H_TO_V(NTSD,GRID%DT,GRID%NPHS,UZ0H,VZ0H,UZ0,VZ0 &
1084 & ,DUDT,DVDT,U,V,HBM2,VTM,IVE,IVW &
1085 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1086 & ,IMS,IME,JMS,JME,KMS,KME &
1087 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1088 uv_htov_tim=uv_htov_tim+timef()-btimx
1089 !
1090 !----------------------------------------------------------------------
1091 !*** STORE ORIGINAL TEMPERATURE ARRAY
1092 !----------------------------------------------------------------------
1093 !
1094 btimx=timef()
1095 !-----------------
1096 #ifdef DM_PARALLEL
1097 # include "HALO_NMM_J.inc"
1098 #endif
1099 !
1100 #ifdef DM_PARALLEL
1101 IF(CONFIG_FLAGS%MP_PHYSICS/=ETAMPNEW)THEN
1102 # include "HALO_NMM_J_3.inc"
1103 ENDIF
1104 #endif
1105 !
1106 #ifdef WRF_CHEM
1107 #ifdef DM_PARALLEL
1108 # include "HALO_NMM_J_2.inc"
1109 #endif
1110 #endif
1111 !-----------------
1112 exch_tim=exch_tim+timef()-btimx
1113 ! this_tim=timef()-btimx
1114 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
1115 ! & ,mpi_comm_comp,irtn)
1116 ! exch_tim_max=exch_tim_max+et_max
1117 !
1118 ICLTEND=-1
1119 btimx=timef()
1120 !
1121 CALL CLTEND(ICLTEND,GRID%NPHS,T,T_OLD,T_ADJ &
1122 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1123 & ,IMS,IME,JMS,JME,KMS,KME &
1124 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1125 !
1126 cltend_tim=cltend_tim+timef()-btimx
1127 ENDIF
1128 !
1129 !----------------------------------------------------------------------
1130 !*** CONVECTIVE PRECIPITATION
1131 !----------------------------------------------------------------------
1132 !
1133 IF(MOD(NTSD,GRID%NCNVC)==0.AND. &
1134 & CONFIG_FLAGS%CU_PHYSICS==KFETASCHEME)THEN
1135 !
1136 btimx=timef()
1137 !-----------------
1138 #ifdef DM_PARALLEL
1139 # include "HALO_NMM_C.inc"
1140 #endif
1141 !-----------------
1142 exch_tim=exch_tim+timef()-btimx
1143 ! this_tim=timef()-btimx
1144 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
1145 ! & ,mpi_comm_comp,irtn)
1146 ! exch_tim_max=exch_tim_max+et_max
1147 ENDIF
1148 !
1149 convection: IF(CONFIG_FLAGS%CU_PHYSICS/=0)THEN
1150 btimx=timef()
1151 !
1152 !*** GET TENDENCIES FOR GD SCHEME.
1153 !
1154 IF(CONFIG_FLAGS%CU_PHYSICS==GDSCHEME)THEN
1155 DT_INV=1./GRID%DT
1156 DO J=JMS,JME
1157 DO K=KMS,KME
1158 DO I=IMS,IME
1159 TTEN(I,K,J)=(T(I,K,J)-TTEN(I,K,J))*DT_INV
1160 QTEN(I,K,J)=(Q(I,K,J)-QTEN(I,K,J))*DT_INV
1161 ENDDO
1162 ENDDO
1163 ENDDO
1164 ENDIF
1165
1166 CALL CUCNVC(NTSD,GRID%DT,GRID%NCNVC,GRID%NRADS,GRID%NRADL &
1167 & ,GPS,RESTRT,HYDRO,CLDEFI,LMH,N_MOIST,GRID%ENSDIM &
1168 & ,MOIST &
1169 & ,DETA1,DETA2,AETA1,AETA2,ETA1,ETA2 &
1170 & ,F_ICE,F_RAIN &
1171 !*** Changes for other cu schemes, most for GD scheme
1172 & ,APR_GR,APR_W,APR_MC,TTEN,QTEN &
1173 & ,APR_ST,APR_AS,APR_CAPMA &
1174 & ,APR_CAPME,APR_CAPMI &
1175 & ,MASS_FLUX,XF_ENS &
1176 & ,PR_ENS,GSW &
1177 #ifdef WRF_CHEM
1178 & ,GD_CLOUD,GD_CLOUD2,RAINCV &
1179 #endif
1180 !
1181 & ,PDTOP,PT,PD,RES,PINT,T,Q,CWM,TCUCN &
1182 & ,OMGALF,U,V,VTM,W,Z,FIS,W0AVG &
1183 & ,PREC,ACPREC,CUPREC,CUPPT,CPRATE &
1184 & ,SM,HBM2,LPBL,CNVBOT,CNVTOP &
1185 & ,HTOP,HBOT,HTOPD,HBOTD,HTOPS,HBOTS &
1186 & ,RTHBLTEN,RQVBLTEN,RTHRATEN &
1187 & ,AVCNVC,ACUTIM,ZERO_3D,IHE,IHW &
1188 & ,GRID,CONFIG_FLAGS &
1189 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1190 & ,IMS,IME,JMS,JME,KMS,KME &
1191 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1192 !
1193 cucnvc_tim=cucnvc_tim+timef()-btimx
1194 !
1195 ENDIF convection
1196 !
1197 !----------------------------------------------------------------------
1198 !*** GRIDSCALE MICROPHYSICS (CONDENSATION & PRECIPITATION)
1199 !----------------------------------------------------------------------
1200 !
1201 IF(MOD(NTSD,GRID%NPHS)==0)THEN
1202 btimx=timef()
1203 !
1204 CALL GSMDRIVE(NTSD,GRID%DT,GRID%NPHS,N_MOIST &
1205 & ,DX_NMM(ITS,JC),GRID%DY,LMH,SM,HBM2,FIS &
1206 & ,DETA1,DETA2,AETA1,AETA2,ETA1,ETA2 &
1207 & ,PDTOP,PT,PD,RES,PINT,T,Q,CWM,TRAIN &
1208 & ,MOIST,SCALAR,NUM_SCALAR &
1209 & ,F_ICE,F_RAIN,F_RIMEF,SR &
1210 & ,PREC,ACPREC,AVRAIN,ZERO_3D &
1211 & ,MP_RESTART_STATE &
1212 & ,TBPVS_STATE &
1213 & ,TBPVS0_STATE &
1214 & ,GRID,CONFIG_FLAGS &
1215 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1216 & ,IMS,IME,JMS,JME,KMS,KME &
1217 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1218 !
1219 gsmdrive_tim=gsmdrive_tim+timef()-btimx
1220 !
1221 !-----------------------------------------------------------------------
1222 !---------PRECIPITATION ASSIMILATION------------------------------------
1223 !-----------------------------------------------------------------------
1224 !
1225 IF (GRID%PCPFLG) THEN
1226 btimx=timef()
1227 !
1228 CALL CHKSNOW(NTSD,GRID%DT,GRID%NPHS,SR,PPTDAT &
1229 & ,IDS,IDE,JDS,JDE,KDS,KDE &
1230 & ,IMS,IME,JMS,JME,KMS,KME &
1231 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1232 CALL ADJPPT(NTSD,GRID%DT,GRID%NPHS,PREC,LSPA,PPTDAT,DDATA &
1233 & ,IDS,IDE,JDS,JDE,KDS,KDE &
1234 & ,IMS,IME,JMS,JME,KMS,KME &
1235 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1236 !
1237 adjppt_tim=adjppt_tim+timef()-btimx
1238 ENDIF
1239 !
1240 !----------------------------------------------------------------------
1241 !*** CALCULATE TEMP TENDENCIES AND RESTORE ORIGINAL TEMPS
1242 !----------------------------------------------------------------------
1243 !
1244 ICLTEND=0
1245 btimx=timef()
1246 !
1247 CALL CLTEND(ICLTEND,GRID%NPHS,T,T_OLD,T_ADJ &
1248 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1249 & ,IMS,IME,JMS,JME,KMS,KME &
1250 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1251 !
1252 cltend_tim=cltend_tim+timef()-btimx
1253 ENDIF
1254 !
1255 !----------------------------------------------------------------------
1256 !*** UPDATE TEMP TENDENCIES FROM CLOUD PROCESSES EVERY TIME STEP
1257 !----------------------------------------------------------------------
1258 !
1259 ICLTEND=1
1260 btimx=timef()
1261 !
1262 CALL CLTEND(ICLTEND,GRID%NPHS,T,T_OLD,T_ADJ &
1263 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1264 & ,IMS,IME,JMS,JME,KMS,KME &
1265 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1266 !
1267 cltend_tim=cltend_tim+timef()-btimx
1268 !
1269 !----------------------------------------------------------------------
1270 !*** LATERAL DIFFUSION
1271 !----------------------------------------------------------------------
1272 !
1273 btimx=timef()
1274 !-----------------
1275 #ifdef DM_PARALLEL
1276 # include "HALO_NMM_K.inc"
1277 #endif
1278 !-----------------
1279 exch_tim=exch_tim+timef()-btimx
1280 ! this_tim=timef()-btimx
1281 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
1282 ! & ,mpi_comm_comp,irtn)
1283 ! exch_tim_max=exch_tim_max+et_max
1284 !
1285 btimx=timef()
1286 !
1287 CALL HDIFF(NTSD,GRID%DT,FIS,DY_NMM,HDAC,HDACV &
1288 & ,HTM,HBM2,VTM,DETA1,GRID%SIGMA &
1289 & ,T,Q,U,V,Q2,Z,W,SM,SICE &
1290 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
1291 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1292 & ,IMS,IME,JMS,JME,KMS,KME &
1293 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1294 !
1295 DO J=MYJS,MYJE
1296 DO K=KTS,KTE
1297 DO I=MYIS,MYIE
1298 !!! MOIST(I,K,J,P_QV)=MAX(0.,Q(I,K,J)/(1.-Q(I,K,J)))
1299 MOIST(I,K,J,P_QV)=Q(I,K,J)/(1.-Q(I,K,J))
1300 ENDDO
1301 ENDDO
1302 ENDDO
1303 !
1304 hdiff_tim=hdiff_tim+timef()-btimx
1305 !
1306 !----------------------------------------------------------------------
1307 !*** UPDATING BOUNDARY VALUES AT HEIGHT POINTS
1308 !----------------------------------------------------------------------
1309 !
1310 btimx=timef()
1311 !-----------------
1312 #ifdef DM_PARALLEL
1313 # include "HALO_NMM_L.inc"
1314 #endif
1315 !
1316 #ifdef DM_PARALLEL
1317 # include "HALO_NMM_L_3.inc"
1318 #endif
1319 !
1320 #ifdef WRF_CHEM
1321 #ifdef DM_PARALLEL
1322 # include "HALO_NMM_L_2.inc"
1323 #endif
1324 #endif
1325 !-----------------
1326 exch_tim=exch_tim+timef()-btimx
1327 ! this_tim=timef()-btimx
1328 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
1329 ! & ,mpi_comm_comp,irtn)
1330 ! exch_tim_max=exch_tim_max+et_max
1331 !
1332 btimx=timef()
1333 !
1334 CALL BOCOH(GRID%ID,NTSD,GRID%DT,NEST,NUNIT_NBC,NBOCO,LAST_TIME,TSPH &
1335 & ,LB,ETA1,ETA2,PDTOP,PT,RES,HTM &
1336 & ,PD_B,T_B,Q_B,U_B,V_B,Q2_B,CWM_B &
1337 & ,PD_BT,T_BT,Q_BT,U_BT,V_BT,Q2_BT,CWM_BT &
1338 & ,PD,T,Q,Q2,CWM,PINT,MOIST,N_MOIST,SCALAR,NUM_SCALAR &
1339 #ifdef WRF_CHEM
1340 & ,CHEM,NUMGAS,CONFIG_FLAGS &
1341 #endif
1342 & ,IJDS,IJDE,GRID%SPEC_BDY_WIDTH,Z &
1343 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
1344 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1345 & ,IMS,IME,JMS,JME,KMS,KME &
1346 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1347 !
1348 bocoh_tim=bocoh_tim+timef()-btimx
1349 ! if(mod(ntsd,n_print_time)==0)then
1350 ! call twr(t,0,'t',ntsd,mype,npes,mpi_comm_comp &
1351 ! & ,ids,ide,jds,jde,kds,kde &
1352 ! & ,ims,ime,jms,jme,kms,kme &
1353 ! & ,its,ite,jts,jte,kts,kte)
1354 ! endif
1355 !
1356 !----------------------------------------------------------------------
1357 !*** IS IT TIME FOR A CHECK POINT ON THE MODEL HISTORY FILE?
1358 !----------------------------------------------------------------------
1359 !
1360 2003 CONTINUE
1361 !
1362 !----------------------------------------------------------------------
1363 !*** PRESSURE GRD, CORIOLIS, DIVERGENCE, AND HORIZ PART OF OMEGA-ALPHA
1364 !----------------------------------------------------------------------
1365 !
1366 btimx=timef()
1367 !-----------------
1368 #ifdef DM_PARALLEL
1369 # include "HALO_NMM_A.inc"
1370 #endif
1371 !
1372 #ifdef DM_PARALLEL
1373 IF(CONFIG_FLAGS%MP_PHYSICS/=ETAMPNEW)THEN
1374 # include "HALO_NMM_A_3.inc"
1375 ENDIF
1376 #endif
1377 !
1378 #ifdef WRF_CHEM
1379 #ifdef DM_PARALLEL
1380 # include "HALO_NMM_A_2.inc"
1381 #endif
1382 #endif
1383 !-----------------
1384 exch_tim=exch_tim+timef()-btimx
1385 ! this_tim=timef()-btimx
1386 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
1387 ! & ,mpi_comm_comp,irtn)
1388 ! exch_tim_max=exch_tim_max+et_max
1389 !
1390 btimx=timef()
1391 !
1392 CALL PFDHT(NTSD,LAST_TIME,PT,DETA1,DETA2,PDTOP,RES,FIS &
1393 & ,HYDRO,GRID%SIGMA,FIRST,DX_NMM,DY_NMM &
1394 & ,HTM,HBM2,VTM,VBM2,VBM3 &
1395 & ,FDIV,FCP,WPDAR,DFL,CPGFU,CPGFV &
1396 & ,PD,PDSL,T,Q,U,V,CWM,OMGALF,PINT,DWDT &
1397 & ,RTOP,DIV,FEW,FNS,FNE,FSE &
1398 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
1399 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1400 & ,IMS,IME,JMS,JME,KMS,KME &
1401 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1402 !
1403 pfdht_tim=pfdht_tim+timef()-btimx
1404 !
1405 !----------------------------------------------------------------------
1406 !*** DIVERGENCE DAMPING
1407 !----------------------------------------------------------------------
1408 !
1409 btimx=timef()
1410 !-----------------
1411 #ifdef DM_PARALLEL
1412 # include "HALO_NMM_B.inc"
1413 #endif
1414 !-----------------
1415 exch_tim=exch_tim+timef()-btimx
1416 ! this_tim=timef()-btimx
1417 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
1418 ! & ,mpi_comm_comp,irtn)
1419 ! exch_tim_max=exch_tim_max+et_max
1420 !
1421 btimx=timef()
1422 !
1423 CALL DDAMP(NTSD,GRID%DT,DETA1,DETA2,PDSL,PDTOP,DIV,HBM2,VTM &
1424 & ,T,U,V,DDMPU,DDMPV &
1425 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
1426 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1427 & ,IMS,IME,JMS,JME,KMS,KME &
1428 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1429 !
1430 ddamp_tim=ddamp_tim+timef()-btimx
1431 !----------------------------------------------------------------------
1432 !----------------------------------------------------------------------
1433 !
1434 IF(FIRST.AND.NTSD==0)THEN
1435 FIRST=.FALSE.
1436 btimx=timef()
1437 !-----------------
1438 #ifdef DM_PARALLEL
1439 # include "HALO_NMM_A.inc"
1440 #endif
1441 #ifdef WRF_CHEM
1442 #ifdef DM_PARALLEL
1443 # include "HALO_NMM_A_2.inc"
1444 #endif
1445 #endif
1446 !-----------------
1447 exch_tim=exch_tim+timef()-btimx
1448 ! this_tim=timef()-btimx
1449 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
1450 ! & ,mpi_comm_comp,irtn)
1451 ! exch_tim_max=exch_tim_max+et_max
1452 !
1453 GO TO 2000
1454 ENDIF
1455 !
1456 !----------------------------------------------------------------------
1457 !*** UPDATING BOUNDARY VALUES AT VELOCITY POINTS
1458 !----------------------------------------------------------------------
1459 !
1460 btimx=timef()
1461 !-----------------
1462 #ifdef DM_PARALLEL
1463 # include "HALO_NMM_C.inc"
1464 #endif
1465 !-----------------
1466 exch_tim=exch_tim+timef()-btimx
1467 ! this_tim=timef()-btimx
1468 ! call mpi_allreduce(this_tim,et_max,1,mpi_real,mpi_max &
1469 ! & ,mpi_comm_comp,irtn)
1470 ! exch_tim_max=exch_tim_max+et_max
1471 !
1472 btimx=timef()
1473 !
1474 CALL BOCOV(GRID%ID,NTSD,GRID%DT,LB,VTM,U_B,V_B,U_BT,V_BT &
1475 & ,U,V &
1476 & ,IJDS,IJDE,GRID%SPEC_BDY_WIDTH &
1477 & ,IHE,IHW,IVE,IVW,INDX3_WRK &
1478 & ,IDS,IDF,JDS,JDF,KDS,KDE &
1479 & ,IMS,IME,JMS,JME,KMS,KME &
1480 & ,ITS,ITE,JTS,JTE,KTS,KTE )
1481 !
1482 bocov_tim=bocov_tim+timef()-btimx
1483 !
1484 !----------------------------------------------------------------------
1485 !*** COPY THE NMM VARIABLE Q2 TO THE WRF VARIABLE TKE_MYJ
1486 !----------------------------------------------------------------------
1487 !
1488 DO J=JTS,JTE
1489 DO K=KTS,KTE
1490 DO I=ITS,ITE
1491 TKE_MYJ(I,K,J)=0.5*Q2(I,K,J) !TKE is q squared over 2
1492 ENDDO
1493 ENDDO
1494 ENDDO
1495 !
1496 !----------------------------------------------------------------------
1497 !
1498 IF(LAST_TIME.AND.ALLOCATED(PPTDAT))THEN
1499 DEALLOCATE(PPTDAT,STAT=ISTAT)
1500 ENDIF
1501 !
1502 !----------------------------------------------------------------------
1503 !
1504 solve_tim=solve_tim+timef()-btim
1505 !
1506 !----------------------------------------------------------------------
1507 !*** PRINT TIMING VARIABLES WHEN DESIRED.
1508 !----------------------------------------------------------------------
1509 !
1510 sum_tim=pdte_tim+adve_tim+vtoa_tim+vadz_tim+hadz_tim+eps_tim &
1511 & +vad2_tim+had2_tim+radiation_tim+rdtemp_tim+turbl_tim &
1512 & +cltend_tim+cucnvc_tim+gsmdrive_tim+hdiff_tim &
1513 & +bocoh_tim+pfdht_tim+ddamp_tim+bocov_tim+uv_htov_tim &
1514 & +exch_tim+adjppt_tim
1515 !
1516 if(mod(ntsd,n_print_time)==0)then
1517 write(0,*)' ntsd=',ntsd,' solve_tim=',solve_tim*1.e-3 &
1518 & ,' sum_tim=',sum_tim*1.e-3
1519 write(0,*)' pdte_tim=',pdte_tim*1.e-3,' pct=',pdte_tim/sum_tim*100.
1520 write(0,*)' adve_tim=',adve_tim*1.e-3,' pct=',adve_tim/sum_tim*100.
1521 write(0,*)' vtoa_tim=',vtoa_tim*1.e-3,' pct=',vtoa_tim/sum_tim*100.
1522 write(0,*)' vadz_tim=',vadz_tim*1.e-3,' pct=',vadz_tim/sum_tim*100.
1523 write(0,*)' hadz_tim=',hadz_tim*1.e-3,' pct=',hadz_tim/sum_tim*100.
1524 write(0,*)' eps_tim=',eps_tim*1.e-3,' pct=',eps_tim/sum_tim*100.
1525 write(0,*)' vad2_tim=',vad2_tim*1.e-3,' pct=',vad2_tim/sum_tim*100.
1526 write(0,*)' had2_tim=',had2_tim*1.e-3,' pct=',had2_tim/sum_tim*100.
1527 write(0,*)' radiation_tim=',radiation_tim*1.e-3,' pct=',radiation_tim/sum_tim*100.
1528 write(0,*)' rdtemp_tim=',rdtemp_tim*1.e-3,' pct=',rdtemp_tim/sum_tim*100.
1529 write(0,*)' turbl_tim=',turbl_tim*1.e-3,' pct=',turbl_tim/sum_tim*100.
1530 write(0,*)' cltend_tim=',cltend_tim*1.e-3,' pct=',cltend_tim/sum_tim*100.
1531 write(0,*)' cucnvc_tim=',cucnvc_tim*1.e-3,' pct=',cucnvc_tim/sum_tim*100.
1532 write(0,*)' gsmdrive_tim=',gsmdrive_tim*1.e-3,' pct=',gsmdrive_tim/sum_tim*100.
1533 write(0,*)' adjppt_tim=',adjppt_tim*1.e-3,' pct=',adjppt_tim/sum_tim*100.
1534 write(0,*)' hdiff_tim=',hdiff_tim*1.e-3,' pct=',hdiff_tim/sum_tim*100.
1535 write(0,*)' bocoh_tim=',bocoh_tim*1.e-3,' pct=',bocoh_tim/sum_tim*100.
1536 write(0,*)' pfdht_tim=',pfdht_tim*1.e-3,' pct=',pfdht_tim/sum_tim*100.
1537 write(0,*)' ddamp_tim=',ddamp_tim*1.e-3,' pct=',ddamp_tim/sum_tim*100.
1538 write(0,*)' bocov_tim=',bocov_tim*1.e-3,' pct=',bocov_tim/sum_tim*100.
1539 write(0,*)' uv_h_to_v_tim=',uv_htov_tim*1.e-3,' pct=',uv_htov_tim/sum_tim*100.
1540 write(0,*)' exch_tim=',exch_tim*1.e-3,' pct=',exch_tim/sum_tim*100.
1541 ! call time_stats(exch_tim,'exchange',ntsd,mype,npes,mpi_comm_comp)
1542 ! write(0,*)' exch_tim_max=',exch_tim_max*1.e-3
1543 !
1544 call field_stats(t,mype,mpi_comm_comp &
1545 & ,ids,ide,jds,jde,kds,kde &
1546 & ,ims,ime,jms,jme,kms,kme &
1547 & ,its,ite,jts,jte,kts,kte)
1548 endif
1549 !
1550 ! if(last_time)then
1551 DEALLOCATE(TTEN,STAT=ISTAT)
1552 DEALLOCATE(QTEN,STAT=ISTAT)
1553 DEALLOCATE(RTHRATEN,STAT=ISTAT)
1554 DEALLOCATE(RTHBLTEN,STAT=ISTAT)
1555 DEALLOCATE(RQVBLTEN,STAT=ISTAT)
1556 !
1557 #define COPY_OUT
1558 #include <nmm_scalar_derefs.inc>
1559 Return
1560 !----------------------------------------------------------------------
1561 !**********************************************************************
1562 !**********************************************************************
1563 !************* EXIT FROM THE TIME LOOP **************************
1564 !**********************************************************************
1565 !**********************************************************************
1566 !----------------------------------------------------------------------
1567 END SUBROUTINE SOLVE_NMM
1568 !----------------------------------------------------------------------
1569 !**********************************************************************
1570 !----------------------------------------------------------------------
1571 SUBROUTINE TWR(ARRAY,KK,FIELD,NTSD,MYPE,NPES,MPI_COMM_COMP &
1572 & ,IDS,IDE,JDS,JDE,KDS,KDE &
1573 & ,IMS,IME,JMS,JME,KMS,KME &
1574 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1575 !----------------------------------------------------------------------
1576 !**********************************************************************
1577 USE MODULE_EXT_INTERNAL
1578 !
1579 IMPLICIT NONE
1580 INCLUDE "mpif.h"
1581 !----------------------------------------------------------------------
1582 !----------------------------------------------------------------------
1583 INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
1584 & ,IMS,IME,JMS,JME,KMS,KME &
1585 & ,ITS,ITE,JTS,JTE,KTS,KTE &
1586 & ,KK,MPI_COMM_COMP,MYPE,NPES,NTSD
1587 !
1588 REAL,DIMENSION(IMS:IME,KMS:KME+KK,JMS:JME),INTENT(IN) :: ARRAY
1589 !
1590 CHARACTER(*),INTENT(IN) :: FIELD
1591 !
1592 !*** LOCAL VARIABLES
1593 !
1594 INTEGER,DIMENSION(MPI_STATUS_SIZE) :: JSTAT
1595 INTEGER,DIMENSION(MPI_STATUS_SIZE,4) :: STATUS_ARRAY
1596 INTEGER,DIMENSION(2) :: IM_REM,JM_REM,IT_REM,JT_REM
1597 !
1598 INTEGER :: I,IENDX,IER,IPE,IRECV,IRTN,ISEND,IUNIT &
1599 & ,J,K,N,NLEN,NSIZE
1600 INTEGER :: ITS_REM,ITE_REM,JTS_REM,JTE_REM
1601 !
1602 REAL,DIMENSION(IDS:IDE,JDS:JDE) :: TWRITE
1603 REAL,ALLOCATABLE,DIMENSION(:) :: VALUES
1604 CHARACTER(5) :: TIMESTEP
1605 CHARACTER(6) :: FMT
1606 CHARACTER(12) :: FILENAME
1607 !----------------------------------------------------------------------
1608 !**********************************************************************
1609 !----------------------------------------------------------------------
1610 !
1611 IF(NTSD<=9)THEN
1612 FMT='(I1.1)'
1613 NLEN=1
1614 ELSEIF(NTSD<=99)THEN
1615 FMT='(I2.2)'
1616 NLEN=2
1617 ELSEIF(NTSD<=999)THEN
1618 FMT='(I3.3)'
1619 NLEN=3
1620 ELSEIF(NTSD<=9999)THEN
1621 FMT='(I4.4)'
1622 NLEN=4
1623 ELSEIF(NTSD<=99999)THEN
1624 FMT='(I5.5)'
1625 NLEN=5
1626 ENDIF
1627 WRITE(TIMESTEP,FMT)NTSD
1628 FILENAME=FIELD//'_'//TIMESTEP(1:NLEN)
1629 !
1630 IF(MYPE==0)THEN
1631 CALL INT_GET_FRESH_HANDLE(IUNIT)
1632 CLOSE(IUNIT)
1633 OPEN(UNIT=IUNIT,FILE=FILENAME,FORM='UNFORMATTED',IOSTAT=IER)
1634 ENDIF
1635 !
1636 !----------------------------------------------------------------------
1637 !!!! DO 500 K=KTS,KTE+KK !Unflipped
1638 !!!! DO 500 K=KTE+KK,KTS,-1
1639 DO 500 K=KDE-1,KDS,-1 !Write LM layers top down for checking
1640 !----------------------------------------------------------------------
1641 !
1642 IF(MYPE==0)THEN
1643 DO J=JTS,JTE
1644 DO I=ITS,ITE
1645 TWRITE(I,J)=ARRAY(I,K,J)
1646 ENDDO
1647 ENDDO
1648 !
1649 DO IPE=1,NPES-1
1650 CALL MPI_RECV(IT_REM,2,MPI_INTEGER,IPE,IPE &
1651 & ,MPI_COMM_COMP,JSTAT,IRECV)
1652 CALL MPI_RECV(JT_REM,2,MPI_INTEGER,IPE,IPE &
1653 & ,MPI_COMM_COMP,JSTAT,IRECV)
1654 !
1655 ITS_REM=IT_REM(1)
1656 ITE_REM=IT_REM(2)
1657 JTS_REM=JT_REM(1)
1658 JTE_REM=JT_REM(2)
1659 !
1660 NSIZE=(ITE_REM-ITS_REM+1)*(JTE_REM-JTS_REM+1)
1661 ALLOCATE(VALUES(1:NSIZE))
1662 !
1663 CALL MPI_RECV(VALUES,NSIZE,MPI_REAL,IPE,IPE &
1664 & ,MPI_COMM_COMP,JSTAT,IRECV)
1665 N=0
1666 DO J=JTS_REM,JTE_REM
1667 DO I=ITS_REM,ITE_REM
1668 N=N+1
1669 TWRITE(I,J)=VALUES(N)
1670 ENDDO
1671 ENDDO
1672 !
1673 DEALLOCATE(VALUES)
1674 !
1675 ENDDO
1676 !
1677 !----------------------------------------------------------------------
1678 ELSE
1679 NSIZE=(ITE-ITS+1)*(JTE-JTS+1)
1680 ALLOCATE(VALUES(1:NSIZE))
1681 !
1682 N=0
1683 DO J=JTS,JTE
1684 DO I=ITS,ITE
1685 N=N+1
1686 VALUES(N)=ARRAY(I,K,J)
1687 ENDDO
1688 ENDDO
1689 !
1690 IT_REM(1)=ITS
1691 IT_REM(2)=ITE
1692 JT_REM(1)=JTS
1693 JT_REM(2)=JTE
1694 !
1695 CALL MPI_SEND(IT_REM,2,MPI_INTEGER,0,MYPE &
1696 & ,MPI_COMM_COMP,ISEND)
1697 CALL MPI_SEND(JT_REM,2,MPI_INTEGER,0,MYPE &
1698 & ,MPI_COMM_COMP,ISEND)
1699 !
1700 CALL MPI_SEND(VALUES,NSIZE,MPI_REAL,0,MYPE &
1701 & ,MPI_COMM_COMP,ISEND)
1702 !
1703 DEALLOCATE(VALUES)
1704 !
1705 ENDIF
1706 !----------------------------------------------------------------------
1707 !
1708 CALL MPI_BARRIER(MPI_COMM_COMP,IRTN)
1709 !
1710 IF(MYPE==0)THEN
1711 !
1712 DO J=JDS,JDE-1
1713 IENDX=IDE-1
1714 IF(MOD(J,2)==0)IENDX=IENDX-1
1715 WRITE(IUNIT)(TWRITE(I,J),I=1,IENDX)
1716 ENDDO
1717 !
1718 ENDIF
1719 !
1720 !----------------------------------------------------------------------
1721 500 CONTINUE
1722 !
1723 IF(MYPE==0)CLOSE(IUNIT)
1724 !----------------------------------------------------------------------
1725 !
1726 END SUBROUTINE TWR
1727 !----------------------------------------------------------------------
1728 !**********************************************************************
1729 !----------------------------------------------------------------------
1730 SUBROUTINE EXIT(NAME,T,Q,U,V,Q2,NTSD,MYPE,MPI_COMM_COMP &
1731 & ,IDS,IDE,JDS,JDE,KDS,KDE &
1732 & ,IMS,IME,JMS,JME,KMS,KME &
1733 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1734 !----------------------------------------------------------------------
1735 !**********************************************************************
1736 USE MODULE_EXT_INTERNAL
1737 !
1738 !----------------------------------------------------------------------
1739 IMPLICIT NONE
1740 !----------------------------------------------------------------------
1741 INCLUDE "mpif.h"
1742 !----------------------------------------------------------------------
1743 INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
1744 & ,IMS,IME,JMS,JME,KMS,KME &
1745 & ,ITS,ITE,JTS,JTE,KTS,KTE &
1746 & ,MYPE,MPI_COMM_COMP,NTSD
1747 !
1748 REAL,DIMENSION(IMS:IME,KMS:KME,JMS:JME),INTENT(IN) :: T,Q,U,V,Q2
1749 CHARACTER(*),INTENT(IN) :: NAME
1750 !
1751 INTEGER :: I,J,K,IEND,IERR,IRET
1752 CHARACTER(256) :: ERRMESS
1753 LOGICAL :: E_BDY,S_BDY
1754 !----------------------------------------------------------------------
1755 IRET=0
1756 100 FORMAT(' EXIT ',A,' AT NTSD=',I5)
1757 IEND=ITE
1758 S_BDY=(JTS==JDS)
1759 E_BDY=(ITE==IDE-1)
1760 !
1761 DO J=JTS,JTE
1762 IEND=ITE
1763 DO K=KTS,KTE
1764 IF(E_BDY.AND.MOD(J,2)==0)IEND=ITE-1
1765 !
1766 DO I=ITS,IEND
1767 IF(T(I,K,J)>330..OR.T(I,K,J)<180..OR.T(I,K,J)/=T(I,K,J))THEN
1768 WRITE(0,100)NAME,NTSD
1769 WRITE(0,200)I,J,K,T(I,K,J),MYPE,NTSD
1770 200 FORMAT(' BAD VALUE I=',I3,' J=',I3,' K=',I2,' T=',E12.5 &
1771 &, ' MYPE=',I3,' NTSD=',I5)
1772 IRET=666
1773 return
1774 ! WRITE(ERRMESS,205)NAME,T(I,K,J),I,K,J,MYPE
1775 205 FORMAT(' EXIT ',A,' TEMPERATURE=',E12.5 &
1776 &, ' AT (',I3,',',I2,',',I3,')',' MYPE=',I3)
1777 ! CALL WRF_ERROR_FATAL(ERRMESS)
1778 ! CALL MPI_ABORT(MPI_COMM_WORLD,1,IERR)
1779 ELSEIF(Q(I,K,J)<-1.E-4.OR.Q(I,K,J)>30.E-3 &
1780 & .OR.Q(I,K,J)/=Q(I,K,J))THEN
1781 WRITE(0,100)NAME,NTSD
1782 WRITE(0,300)I,J,K,Q(I,K,J),MYPE,NTSD
1783 300 FORMAT(' BAD VALUE I=',I3,' J=',I3,' K=',I2,' Q=',E12.5 &
1784 &, ' MYPE=',I3,' NTSD=',I5)
1785 IRET=666
1786 return
1787 ! WRITE(ERRMESS,305)NAME,Q(I,K,J),I,K,J,MYPE
1788 305 FORMAT(' EXIT ',A,' SPEC HUMIDITY=',E12.5 &
1789 &, ' AT (',I3,',',I2,',',I3,')',' MYPE=',I3)
1790 ! CALL WRF_ERROR_FATAL(ERRMESS)
1791 ! CALL MPI_ABORT(MPI_COMM_WORLD,1,IERR)
1792 ENDIF
1793 ENDDO
1794 ENDDO
1795 ENDDO
1796 !
1797 DO J=JTS,JTE
1798 IEND=ITE
1799 DO K=KTS,KTE
1800 IF(E_BDY.AND.MOD(J,2)==1)IEND=ITE-1
1801 DO I=ITS,IEND
1802 IF(ABS(U(I,K,J))>125..OR.ABS(V(I,K,J))>125. &
1803 & .OR.U(I,K,J)/=U(I,K,J).OR.V(I,K,J)/=V(I,K,J))THEN
1804 WRITE(0,100)NAME,NTSD
1805 WRITE(0,400)I,J,K,U(I,K,J),V(I,K,J),MYPE,NTSD
1806 400 FORMAT(' BAD VALUE I=',I3,' J=',I3,' K=',I2,' U=',E12.5 &
1807 &, ' V=',E12.5,' MYPE=',I3,' NTSD=',I5)
1808 IRET=666
1809 return
1810 ! WRITE(ERRMESS,405)NAME,U(I,K,J),V(I,K,J),I,K,J,MYPE
1811 405 FORMAT(' EXIT ',A,' U=',E12.5,' V=',E12.5 &
1812 &, ' AT (',I3,',',I2,',',I3,')',' MYPE=',I3)
1813 ! CALL WRF_ERROR_FATAL(ERRMESS)
1814 ! CALL MPI_ABORT(MPI_COMM_WORLD,1,IERR)
1815 ENDIF
1816 ENDDO
1817 ENDDO
1818 ENDDO
1819 !----------------------------------------------------------------------
1820 END SUBROUTINE EXIT
1821 !----------------------------------------------------------------------
1822 !**********************************************************************
1823 !----------------------------------------------------------------------
1824 SUBROUTINE TIME_STATS(TIME_LCL,NAME,NTSD,MYPE,NPES,MPI_COMM_COMP)
1825 !----------------------------------------------------------------------
1826 !**********************************************************************
1827 USE MODULE_EXT_INTERNAL
1828 !
1829 !----------------------------------------------------------------------
1830 IMPLICIT NONE
1831 !----------------------------------------------------------------------
1832 INCLUDE "mpif.h"
1833 !----------------------------------------------------------------------
1834 INTEGER,INTENT(IN) :: MPI_COMM_COMP,MYPE,NPES,NTSD
1835 REAL,INTENT(IN) :: TIME_LCL
1836 !
1837 CHARACTER(*),INTENT(IN) :: NAME
1838 !
1839 !*** LOCAL VARIABLES
1840 !
1841 INTEGER,DIMENSION(MPI_STATUS_SIZE) :: JSTAT
1842 INTEGER,DIMENSION(MPI_STATUS_SIZE,4) :: STATUS_ARRAY
1843 INTEGER,ALLOCATABLE,DIMENSION(:) :: ID_PE,IPE_SORT
1844 !
1845 INTEGER :: IPE,IPE_MAX,IPE_MEDIAN,IPE_MIN,IRECV,IRTN,ISEND &
1846 & ,N,N_MEDIAN,NLEN
1847 !
1848 REAL,ALLOCATABLE,DIMENSION(:) :: TIME,SORT_TIME
1849 REAL,DIMENSION(2) :: REMOTE
1850 REAL :: TIME_MAX,TIME_MEAN,TIME_MEDIAN,TIME_MIN
1851 !
1852 CHARACTER(5) :: TIMESTEP
1853 CHARACTER(6) :: FMT
1854 CHARACTER(25) :: TITLE
1855 !----------------------------------------------------------------------
1856 !**********************************************************************
1857 !----------------------------------------------------------------------
1858 !
1859 IF(NTSD<=9)THEN
1860 FMT='(I1.1)'
1861 NLEN=1
1862 ELSEIF(NTSD<=99)THEN
1863 FMT='(I2.2)'
1864 NLEN=2
1865 ELSEIF(NTSD<=999)THEN
1866 FMT='(I3.3)'
1867 NLEN=3
1868 ELSEIF(NTSD<=9999)THEN
1869 FMT='(I4.4)'
1870 NLEN=4
1871 ELSEIF(NTSD<=99999)THEN
1872 FMT='(I5.5)'
1873 NLEN=5
1874 ENDIF
1875 WRITE(TIMESTEP,FMT)NTSD
1876 TITLE=NAME//'_'//TIMESTEP(1:NLEN)
1877 !
1878 !----------------------------------------------------------------------
1879 !
1880 IF(MYPE==0)THEN
1881 ALLOCATE(TIME(1:NPES))
1882 ALLOCATE(SORT_TIME(1:NPES))
1883 ALLOCATE(ID_PE(1:NPES))
1884 ALLOCATE(IPE_SORT(1:NPES))
1885 !
1886 TIME(1)=TIME_LCL
1887 ID_PE(1)=MYPE
1888 !
1889 !*** COLLECT TIMES AND PE VALUES FROM OTHER PEs
1890 !
1891 DO IPE=1,NPES-1
1892 CALL MPI_RECV(REMOTE,2,MPI_REAL,IPE,IPE &
1893 & ,MPI_COMM_COMP,JSTAT,IRECV)
1894 !
1895 TIME(IPE+1)=REMOTE(1)
1896 ID_PE(IPE+1)=NINT(REMOTE(2))
1897 ENDDO
1898 !
1899 !*** NOW GET STATS.
1900 !*** FIRST THE MAX, MIN, AND MEAN TIMES.
1901 !
1902 TIME_MEAN=0.
1903 TIME_MAX=-1.
1904 TIME_MIN=1.E10
1905 IPE_MAX=-1
1906 IPE_MIN=-1
1907 !
1908 DO N=1,NPES
1909 TIME_MEAN=TIME_MEAN+TIME(N)
1910 !
1911 IF(TIME(N)>TIME_MAX)THEN
1912 TIME_MAX=TIME(N)
1913 IPE_MAX=ID_PE(N)
1914 ENDIF
1915 !
1916 IF(TIME(N)<TIME_MIN)THEN
1917 TIME_MIN=TIME(N)
1918 IPE_MIN=ID_PE(N)
1919 ENDIF
1920 !
1921 ENDDO
1922 !
1923 TIME_MAX=TIME_MAX*1.E-3
1924 TIME_MIN=TIME_MIN*1.E-3
1925 TIME_MEAN=TIME_MEAN*1.E-3/REAL(NPES)
1926 !
1927 !*** THEN THE MEDIAN TIME.
1928 !
1929 CALL SORT(TIME,NPES,SORT_TIME,IPE_SORT)
1930 N_MEDIAN=(NPES+1)/2
1931 TIME_MEDIAN=SORT_TIME(N_MEDIAN)*1.E-3
1932 IPE_MEDIAN=IPE_SORT(N_MEDIAN)
1933 !
1934 !----------------------------------------------------------------------
1935 ELSE
1936 !
1937 !*** SEND TIME AND PE VALUE TO PE0.
1938 !
1939 REMOTE(1)=TIME_LCL
1940 REMOTE(2)=REAL(MYPE)
1941 !
1942 CALL MPI_SEND(REMOTE,2,MPI_REAL,0,MYPE,MPI_COMM_COMP,ISEND)
1943 !
1944 ENDIF
1945 !----------------------------------------------------------------------
1946 !
1947 CALL MPI_BARRIER(MPI_COMM_COMP,IRTN)
1948 !
1949 !*** WRITE RESULTS
1950 !
1951 IF(MYPE==0)THEN
1952 WRITE(0,100)TITLE
1953 WRITE(0,105)TIME_MAX,IPE_MAX
1954 WRITE(0,110)TIME_MIN,IPE_MIN
1955 WRITE(0,115)TIME_MEDIAN,IPE_MEDIAN
1956 WRITE(0,120)TIME_MEAN
1957 100 FORMAT(' Time for ',A)
1958 105 FORMAT(' Maximum=',G11.5,' for PE ',I2.2)
1959 110 FORMAT(' Minimum=',G11.5,' for PE ',I2.2)
1960 115 FORMAT(' Median =',G11.5,' for PE ',I2.2)
1961 120 FORMAT(' Mean =',G11.5)
1962 ENDIF
1963 !----------------------------------------------------------------------
1964 !
1965 END SUBROUTINE TIME_STATS
1966 !
1967 !----------------------------------------------------------------------
1968 !**********************************************************************
1969 !----------------------------------------------------------------------
1970 SUBROUTINE SORT(DATA,NPES,DATA_SORTED,IPE_SORTED)
1971 !----------------------------------------------------------------------
1972 !***
1973 !*** SORT DATA FROM MULTIPLE PEs. SEND BACK THE SORTED DATA ITEMS
1974 !*** ALONG WITH THE ASSOCIATED TASK IDs.
1975 !***
1976 !----------------------------------------------------------------------
1977 IMPLICIT NONE
1978 !----------------------------------------------------------------------
1979 INTEGER,INTENT(IN) :: NPES
1980 REAL,DIMENSION(NPES),INTENT(IN) :: DATA
1981 !
1982 INTEGER,DIMENSION(NPES),INTENT(OUT) :: IPE_SORTED
1983 REAL,DIMENSION(NPES),INTENT(OUT) :: DATA_SORTED
1984 !----------------------------------------------------------------------
1985 TYPE :: DATA_LINK
1986 REAL :: VALUE
1987 INTEGER :: IPE
1988 TYPE(DATA_LINK),POINTER :: NEXT_VALUE
1989 END TYPE
1990 !----------------------------------------------------------------------
1991 !
1992 !*** LOCAL VARIABLES
1993 !
1994 !----------------------------------------------------------------------
1995 INTEGER :: ISTAT,N
1996 !
1997 TYPE(DATA_LINK),POINTER :: HEAD,TAIL ! Smallest, largest
1998 TYPE(DATA_LINK),POINTER :: PTR_NEW ! Each new value
1999 TYPE(DATA_LINK),POINTER :: PTR1,PTR2 ! Working pointers
2000 !----------------------------------------------------------------------
2001 !**********************************************************************
2002 !----------------------------------------------------------------------
2003 pe_loop: DO N=1,NPES
2004 ALLOCATE(PTR_NEW,STAT=ISTAT) ! Location for next data items
2005 PTR_NEW%VALUE=DATA(N)
2006 PTR_NEW%IPE=N-1
2007 !
2008 !----------------------------------------------------------------------
2009 !*** DETERMINE WHERE IN LIST TO INSERT VALUE.
2010 !*** FIRST THE INITIAL DATA VALUE.
2011 !----------------------------------------------------------------------
2012 !
2013 ! main: IF(.NOT.ASSOCIATED(HEAD))THEN
2014 main: IF(N==1)THEN
2015 HEAD=>PTR_NEW
2016 TAIL=>HEAD
2017 NULLIFY(PTR_NEW%NEXT_VALUE)
2018 !
2019 !----------------------------------------------------------------------
2020 !*** THE NEW VALUE IS LESS THAN THE SMALLEST VALUE ALREADY SORTED.
2021 !----------------------------------------------------------------------
2022 !
2023 ELSE
2024 check: IF(PTR_NEW%VALUE<HEAD%VALUE)THEN
2025 PTR_NEW%NEXT_VALUE=>HEAD
2026 HEAD=>PTR_NEW
2027 !
2028 !----------------------------------------------------------------------
2029 !*** THE NEW VALUE IS GREATER THAN THE LARGEST VALUE ALREADY SORTED.
2030 !----------------------------------------------------------------------
2031 !
2032 ELSEIF(PTR_NEW%VALUE>=TAIL%VALUE)THEN
2033 TAIL%NEXT_VALUE=>PTR_NEW ! This is what connects the former
2034 ! final value in the list to
2035 ! the new value being appended.
2036 TAIL=>PTR_NEW
2037 NULLIFY(TAIL%NEXT_VALUE)
2038 !
2039 !----------------------------------------------------------------------
2040 !*** THE NEW VALUE IS IN BETWEEN VALUES ALREADY SORTED.
2041 !----------------------------------------------------------------------
2042 !
2043 ELSE
2044 PTR1=>HEAD
2045 PTR2=>PTR1%NEXT_VALUE
2046 !
2047 search: DO
2048 IF((PTR_NEW%VALUE>=PTR1%VALUE).AND. &
2049 & (PTR_NEW%VALUE<PTR2%VALUE))THEN
2050 PTR_NEW%NEXT_VALUE=>PTR2
2051 PTR1%NEXT_VALUE=>PTR_NEW
2052 EXIT search
2053 ENDIF
2054 !
2055 PTR1=>PTR2
2056 PTR2=>PTR2%NEXT_VALUE
2057 ENDDO search
2058 !
2059 ENDIF check
2060 !
2061 ENDIF main
2062 !
2063 ENDDO pe_loop
2064 !
2065 !----------------------------------------------------------------------
2066 !*** COLLECT THE SORTED NUMBERS FROM THE LINKED LIST.
2067 !----------------------------------------------------------------------
2068 !
2069 PTR1=>HEAD
2070 !
2071 DO N=1,NPES
2072 ! IF(.NOT.ASSOCIATED(PTR_NEW))EXIT
2073 DATA_SORTED(N)=PTR1%VALUE
2074 IPE_SORTED(N)=PTR1%IPE
2075 PTR1=>PTR1%NEXT_VALUE
2076 ENDDO
2077 !
2078 DEALLOCATE(PTR_NEW)
2079 NULLIFY (HEAD,TAIL,PTR1,PTR2)
2080 !----------------------------------------------------------------------
2081 END SUBROUTINE SORT
2082 !----------------------------------------------------------------------
2083 !**********************************************************************
2084 !----------------------------------------------------------------------
2085 SUBROUTINE FIELD_STATS(FIELD,MYPE,MPI_COMM_COMP &
2086 & ,IDS,IDE,JDS,JDE,KDS,KDE &
2087 & ,IMS,IME,JMS,JME,KMS,KME &
2088 & ,ITS,ITE,JTS,JTE,KTS,KTE)
2089 !----------------------------------------------------------------------
2090 !***
2091 !*** GENERATE STANDARD LAYER STATISTICS FOR THE DESIRED FIELD.
2092 !***
2093 !----------------------------------------------------------------------
2094 IMPLICIT NONE
2095 !----------------------------------------------------------------------
2096 INCLUDE "mpif.h"
2097 !----------------------------------------------------------------------
2098 !
2099 INTEGER,INTENT(IN) :: MPI_COMM_COMP,MYPE
2100 INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
2101 & ,IMS,IME,JMS,JME,KMS,KME &
2102 & ,ITS,ITE,JTS,JTE,KTS,KTE
2103 !
2104 REAL,DIMENSION(IMS:IME,KMS:KME,JMS:JME),INTENT(IN) :: FIELD
2105 !
2106 !----------------------------------------------------------------------
2107 !*** LOCAL
2108 !----------------------------------------------------------------------
2109 !
2110 INTEGER,PARAMETER :: DOUBLE=SELECTED_REAL_KIND(15,300)
2111 !
2112 INTEGER :: I,IEND,IRTN,I_BY_J,J,K,KFLIP
2113 !
2114 REAL :: FIKJ,FMAXK,FMINK
2115 REAL(KIND=DOUBLE) :: F_MEAN,POINTS,RMS,ST_DEV,SUMFK,SUMF2K
2116 REAL,DIMENSION(KTS:KTE) :: FMAX,FMAX_0,FMIN,FMIN_0
2117 REAL(KIND=DOUBLE),DIMENSION(KTS:KTE) :: SUMF,SUMF_0,SUMF2,SUMF2_0
2118 !----------------------------------------------------------------------
2119 !
2120 I_BY_J=(IDE-IDS)*(JDE-JDS)-(JDE-JDS-1)/2 ! This assumes that
2121 ! IDE and JDE are each
2122 ! one greater than the
2123 ! true grid size.
2124 !
2125 layer_loop: DO K=KTS,KTE
2126 !
2127 FMAXK=-1.E10
2128 FMINK=1.E10
2129 SUMFK=0.
2130 SUMF2K=0.
2131 !
2132 DO J=JTS,JTE
2133 IEND=MIN(ITE,IDE-1)
2134 IF(MOD(J,2)==0.AND.ITE==IDE-1)IEND=IEND-1
2135 DO I=ITS,IEND
2136 FIKJ=FIELD(I,K,J)
2137 FMAXK=MAX(FMAXK,FIKJ)
2138 FMINK=MIN(FMINK,FIKJ)
2139 SUMFK=SUMFK+FIKJ
2140 SUMF2K=SUMF2K+FIKJ*FIKJ
2141 ENDDO
2142 ENDDO
2143 !
2144 FMAX(K)=FMAXK
2145 FMIN(K)=FMINK
2146 SUMF(K)=SUMFK
2147 SUMF2(K)=SUMF2K
2148 !
2149 ENDDO layer_loop
2150 !
2151 !----------------------------------------------------------------------
2152 !*** GLOBAL STATS
2153 !----------------------------------------------------------------------
2154 !
2155 CALL MPI_REDUCE(SUMF,SUMF_0,KTE,MPI_REAL8,MPI_SUM,0 &
2156 & ,MPI_COMM_COMP,IRTN)
2157 CALL MPI_REDUCE(SUMF2,SUMF2_0,KTE,MPI_REAL8,MPI_SUM,0 &
2158 & ,MPI_COMM_COMP,IRTN)
2159 CALL MPI_REDUCE(FMAX,FMAX_0,KTE,MPI_REAL,MPI_MAX,0 &
2160 & ,MPI_COMM_COMP,IRTN)
2161 CALL MPI_REDUCE(FMIN,FMIN_0,KTE,MPI_REAL,MPI_MIN,0 &
2162 & ,MPI_COMM_COMP,IRTN)
2163 !
2164 IF(MYPE==0)THEN
2165 POINTS=I_BY_J
2166 DO K=KTE,KTS,-1
2167 F_MEAN=SUMF_0(K)/POINTS
2168 ST_DEV=SQRT((POINTS*SUMF2_0(K)-SUMF_0(K)*SUMF_0(K))/ &
2169 & (POINTS*(POINTS-1)))
2170 RMS=SQRT(SUMF2_0(K)/POINTS)
2171 KFLIP=KTE-K+1
2172 WRITE(0,101)KFLIP,FMAX_0(K),FMIN_0(K)
2173 WRITE(0,102)F_MEAN,ST_DEV,RMS
2174 101 FORMAT(' LAYER=',I2,' MAX=',E13.6,' MIN=',E13.6)
2175 102 FORMAT(9X,' MEAN=',E13.6,' STDEV=',E13.6,' RMS=',E13.6)
2176 ENDDO
2177 ENDIF
2178 !----------------------------------------------------------------------
2179 END SUBROUTINE FIELD_STATS
2180 !----------------------------------------------------------------------
2181 FUNCTION TIMEF()
2182 REAL*8 TIMEF
2183 INTEGER IC, IR
2184 CALL SYSTEM_CLOCK(COUNT=IC, COUNT_RATE=IR)
2185 TIMEF=REAL(IC)/REAL(IR) * 1000.0
2186 END