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