Oasis3-MCT
 All Classes Files Functions Variables Macros Pages
mod_oasis_advance.F90
Go to the documentation of this file.
1 
2 !> Advances the OASIS coupling
3 
5 
10  USE mod_oasis_map
11  USE mod_oasis_part
12  USE mod_oasis_timer
13  USE mod_oasis_var
14  USE mod_oasis_sys
15  USE mod_oasis_io
16  USE mod_oasis_mpi
17  USE mct_mod
18 
19  IMPLICIT NONE
20 
21  private
22 
23  public oasis_advance_init
24  public oasis_advance_run
25 
26 
27 contains
28 
29 !---------------------------------------------------------------------
30 
31 !> Initializes the OASIS fields
32 
33  SUBROUTINE oasis_advance_init(kinfo)
34 
35 ! ----------------------------------------------------------------
36 ! This routine handles initial restart and communication
37 ! of data for fields with positive lags
38 ! ----------------------------------------------------------------
39 
40  IMPLICIT none
41 ! ----------------------------------------------------------------
42  INTEGER(kind=ip_i4_p), intent(inout) :: kinfo !< status, not used
43 ! ----------------------------------------------------------------
44  integer(kind=ip_i4_p) :: cplid,partid,varid
45  INTEGER(kind=ip_i4_p) :: nf,lsize,nflds,npc
46  integer(kind=ip_i4_p) :: dt,ltime,lag,getput
47  integer(kind=ip_i4_p) :: msec
48  real (kind=ip_r8_p), allocatable :: array(:) ! data
49  real (kind=ip_r8_p), allocatable :: array2(:) ! data
50  real (kind=ip_r8_p), allocatable :: array3(:) ! data
51  real (kind=ip_r8_p), allocatable :: array4(:) ! data
52  real (kind=ip_r8_p), allocatable :: array5(:) ! data
53  logical :: a2on,a3on,a4on,a5on ! data 2-5 logicals
54  integer(kind=ip_i4_p) :: mseclag ! model time + lag
55  character(len=ic_xl) :: rstfile ! restart filename
56  character(len=ic_med) :: vstring ! temporary string
57  type(prism_coupler_type),pointer :: pcpointer
58  character(len=*),parameter :: subname = '(oasis_advance_init)'
59 
60  call oasis_debug_enter(subname)
61 
62  !----------------------------------------------------------------
63  !> oasis_advance_init does the following
64  !> * Aborts if it's called from non-active tasks
65  !----------------------------------------------------------------
66 
67  if (mpi_comm_local == mpi_comm_null) then
68  write(nulprt,*) subname,estr,'called on non coupling task'
69  call oasis_abort()
70  endif
71 
72  kinfo = oasis_ok
73 
74  call oasis_timer_start('advance_init')
75 
76  if (oasis_debug >= 2) then
77  write(nulprt,*) ' subname at time time+lag act: field '
78  write(nulprt,*) ' diags : fldname min max sum '
79  endif
80 
81  !----------------------------------------------------------------
82  !> * Loop over all coupler connections,
83  !> Loop over get and put connections,
84  !> For valid connections
85  !----------------------------------------------------------------
86 
87  call oasis_debug_note(subname//' loop over cplid')
88  DO cplid = 1,prism_mcoupler
89  DO npc = 1,2
90  if (npc == 1) pcpointer => prism_coupler_put(cplid)
91  if (npc == 2) pcpointer => prism_coupler_get(cplid)
92  if (pcpointer%valid) then
93  dt = pcpointer%dt
94  lag = pcpointer%lag
95  ltime = pcpointer%ltime
96  getput= pcpointer%getput
97  rstfile=trim(pcpointer%rstfile)
98  partid= pcpointer%partID
99  msec = 0 ! reasonable default to start with
100  mseclag = msec
101 
102  IF (oasis_debug >= 2) THEN
103  WRITE(nulprt,*) '----------------------------------------------------------------'
104  WRITE(nulprt,*) subname,' Field cplid :',cplid
105  WRITE(nulprt,*) subname,' read restart:',npc,trim(pcpointer%fldlist)
106  WRITE(nulprt,*) subname,' lag prism_mcoupler :',lag,prism_mcoupler
107  WRITE(nulprt,*) '----------------------------------------------------------------'
108  CALL oasis_flush(nulprt)
109  ENDIF
110 
111  !------------------------------------------------
112  !> * Checks that lag is reasonable
113  !------------------------------------------------
114 
115  IF (lag > dt .OR. lag <= -dt) THEN
116  WRITE(nulprt,*) subname,estr,'lag out of dt range cplid/dt/lag=',cplid,dt,lag
117  CALL oasis_abort()
118  ENDIF
119 
120  !------------------------------------------------
121  !> * For put fields that need to read a restart file because of a lag,
122  !> call oasis_advance_run and read in the restart file.
123  !> Set readrest to true in oasis_advance_run to indicate it's called from init.
124  ! right now, do not know whether any hot map terms are there
125  ! assume they are if something is read, otherwise not
126  !------------------------------------------------
127  IF ( (getput == oasis3_put .AND. lag > 0) ) THEN
128  msec=0-lag
129  lsize = mct_avect_lsize(pcpointer%aVect1)
130  nflds = mct_avect_nrattr(pcpointer%aVect1)
131 
132  ALLOCATE(array(lsize))
133  ALLOCATE(array2(lsize))
134  ALLOCATE(array3(lsize))
135  ALLOCATE(array4(lsize))
136  ALLOCATE(array5(lsize))
137 
138  DO nf = 1,nflds
139  varid = pcpointer%varid(nf)
140  CALL oasis_advance_run(oasis_out,varid,msec,kinfo,nff=nf, &
141  namid=cplid,array1din=array,&
142  readrest=.true., array2=array2,array3=array3, &
143  array4=array4,array5=array5)
144  ENDDO
145  IF (oasis_debug >= 2) THEN
146  write(nulprt,*) subname,' advance_run ',cplid,trim(pcpointer%fldlist)
147  endif
148 
149  DEALLOCATE(array)
150  DEALLOCATE(array2)
151  DEALLOCATE(array3)
152  DEALLOCATE(array4)
153  DEALLOCATE(array5)
154  ENDIF ! put
155  ENDIF ! valid
156  enddo ! npc
157  ENDDO ! cplid
158 
159  !----------------------------------------------------------------
160  !> * Loop over all coupler connections,
161  !> Loop over get and put connections,
162  !> For valid connections
163  !----------------------------------------------------------------
164 
165  DO cplid=1, prism_mcoupler
166  DO npc = 1,2
167  if (npc == 1) pcpointer => prism_coupler_put(cplid)
168  if (npc == 2) pcpointer => prism_coupler_get(cplid)
169  IF (pcpointer%valid) then
170  dt = pcpointer%dt
171  lag = pcpointer%lag
172  ltime = pcpointer%ltime
173  getput= pcpointer%getput
174  rstfile=trim(pcpointer%rstfile)
175  partid= pcpointer%partID
176  msec = 0 ! reasonable default to start with
177  mseclag = msec
178 
179  IF (oasis_debug >= 2) THEN
180  WRITE(nulprt,*) '----------------------------------------------------------------'
181  WRITE(nulprt,*) subname,' Field cplid :',cplid
182  WRITE(nulprt,*) subname,' loctrans:',npc,trim(pcpointer%fldlist)
183  WRITE(nulprt,*) '----------------------------------------------------------------'
184  CALL oasis_flush(nulprt)
185  ENDIF
186 
187  !------------------------------------------------
188  !> * Read restart for LOCTRANS fields.
189  !> Do after restart and advance above because prism_advance_run
190  !> fills in the avect with the array info
191  !------------------------------------------------
192 
193  call oasis_debug_note(subname//' check for loctrans restart')
194  IF (getput == oasis3_put .AND. pcpointer%trans /= ip_instant) THEN
195  if (len_trim(rstfile) < 1) then
196  write(nulprt,*) subname,estr,'restart undefined'
197  call oasis_abort()
198  endif
199  if (oasis_debug >= 2) then
200  write(nulprt,*) subname,' at ',msec,mseclag,' RTRN: ',&
201  trim(pcpointer%fldlist),' ',trim(rstfile)
202  endif
203  lsize = mct_avect_lsize(pcpointer%aVect1)
204 
205  write(vstring,'(a,i6.6,a)') 'loc',pcpointer%namID,'_cnt'
206  call oasis_io_read_array(rstfile,prism_part(partid)%mpicom,iarray=pcpointer%avcnt,&
207  ivarname=trim(vstring),abort=.false.)
208 
209  write(vstring,'(a,i6.6,a)') 'loc',pcpointer%namID,'_'
210  call oasis_io_read_avfile(rstfile,pcpointer%avect1,&
211  prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
212  abort=.false.,nampre=trim(vstring))
213 
214  call mct_avect_init(pcpointer%aVect2,pcpointer%aVect1,lsize)
215  call mct_avect_zero(pcpointer%aVect2)
216  write(vstring,'(a,i6.6,a)') 'av2loc',pcpointer%namID,'_'
217  call oasis_io_read_avfile(rstfile,pcpointer%avect2,&
218  prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
219  abort=.false.,nampre=trim(vstring),&
220  didread=pcpointer%aVon(2))
221  if (.not. pcpointer%aVon(2)) then
222  call mct_avect_clean(pcpointer%avect2)
223  endif
224 
225  call mct_avect_init(pcpointer%aVect3,pcpointer%aVect1,lsize)
226  call mct_avect_zero(pcpointer%aVect3)
227  write(vstring,'(a,i6.6,a)') 'av3loc',pcpointer%namID,'_'
228  call oasis_io_read_avfile(rstfile,pcpointer%avect3,&
229  prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
230  abort=.false.,nampre=trim(vstring),&
231  didread=pcpointer%aVon(3))
232  if (.not. pcpointer%aVon(3)) then
233  call mct_avect_clean(pcpointer%avect3)
234  endif
235 
236  call mct_avect_init(pcpointer%aVect4,pcpointer%aVect1,lsize)
237  call mct_avect_zero(pcpointer%aVect4)
238  write(vstring,'(a,i6.6,a)') 'av4loc',pcpointer%namID,'_'
239  call oasis_io_read_avfile(rstfile,pcpointer%avect4,&
240  prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
241  abort=.false.,nampre=trim(vstring),&
242  didread=pcpointer%aVon(4))
243  if (.not. pcpointer%aVon(4)) then
244  call mct_avect_clean(pcpointer%avect4)
245  endif
246 
247  call mct_avect_init(pcpointer%aVect5,pcpointer%aVect1,lsize)
248  call mct_avect_zero(pcpointer%aVect5)
249  write(vstring,'(a,i6.6,a)') 'av5loc',pcpointer%namID,'_'
250  call oasis_io_read_avfile(rstfile,pcpointer%avect5,&
251  prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
252  abort=.false.,nampre=trim(vstring),&
253  didread=pcpointer%aVon(5))
254  if (.not. pcpointer%aVon(5)) then
255  call mct_avect_clean(pcpointer%avect5)
256  endif
257 
258  if (oasis_debug >= 20) then
259  write(nulprt,*) subname,' DEBUG read loctrans restart',&
260  cplid,pcpointer%avcnt
261  write(nulprt,*) subname,' DEBUG read loctrans restart',cplid,&
262  minval(pcpointer%avect1%rAttr),&
263  maxval(pcpointer%avect1%rAttr)
264  endif
265  endif
266  ENDIF ! valid
267  enddo ! npc
268  ENDDO ! cplid
269  call oasis_timer_stop('advance_init')
270 
271  call oasis_debug_exit(subname)
272 
273  end SUBROUTINE oasis_advance_init
274 !---------------------------------------------------------------------
275 
276 !> Advances the OASIS coupling
277 
278 !> Only one from array1din, array1dout, or array2dout can be passed in.
279 !> readrest is set to true when called by the oasis_advance_init method.
280 !> Arrays 2 to 5 are for the higher order terms (hot)
281 
282  SUBROUTINE oasis_advance_run(mop,varid,msec,kinfo,nff,namid,&
283  array1din,array1dout,array2dout,readrest,&
284  a2on,array2,a3on,array3,a4on,array4,a5on,array5)
285 
286  IMPLICIT none
287 ! ----------------------------------------------------------------
288  integer(kind=ip_i4_p), intent(in) :: mop !< OASIS_Out or OASIS_In
289  INTEGER(kind=ip_i4_p), intent(in) :: varid !< prism_var id
290  INTEGER(kind=ip_i4_p), intent(in) :: msec !< model time
291  INTEGER(kind=ip_i4_p), intent(inout) :: kinfo !< status
292 
293  INTEGER(kind=ip_i4_p), optional :: nff !< specify particular field for restart
294  INTEGER(kind=ip_i4_p), optional :: namid !< only do this namcouple
295  !< method for restart
296  REAL (kind=ip_r8_p), optional :: array1din(:) !< 1D put data
297  REAL (kind=ip_r8_p), optional :: array1dout(:) !< 1D get data
298  REAL (kind=ip_r8_p), optional :: array2dout(:,:) !< 2D get data
299  logical , optional :: readrest !< special flag to indicate this
300  !< is called from the advance_init
301  logical , optional :: a2on !< logical for array2
302  REAL (kind=ip_r8_p), optional :: array2(:) !< hot put data
303  logical , optional :: a3on !< logical for array3
304  REAL (kind=ip_r8_p), optional :: array3(:) !< hot put data
305  logical , optional :: a4on !< logical for array4
306  REAL (kind=ip_r8_p), optional :: array4(:) !< hot put data
307  logical , optional :: a5on !< logical for array5
308  REAL (kind=ip_r8_p), optional :: array5(:) !< hot put data
309 ! ----------------------------------------------------------------
310  character(len=ic_lvar):: vname
311  INTEGER(kind=ip_i4_p) :: cplid,rouid,mapid,partid
312  INTEGER(kind=ip_i4_p) :: nfav,nsav,nsa,n,nc,nf,npc
313  INTEGER(kind=ip_i4_p) :: lsize,nflds
314  integer(kind=ip_i4_p) :: tag,dt,ltime,lag,getput,maxtime,conserv
315  logical :: consbfb
316  logical :: sndrcv,output,input,unpack
317  logical :: snddiag,rcvdiag
318  logical :: arrayon(prism_coupler_avsmax)
319  LOGICAL :: a11on,a22on,a33on,a44on,a55on
320  real(kind=ip_double_p):: sndmult,sndadd,rcvmult,rcvadd
321  character(len=ic_xl) :: rstfile ! restart filename
322  character(len=ic_xl) :: inpfile ! input filename
323  integer(kind=ip_i4_p) :: nx,ny
324  integer(kind=ip_i4_p) :: mseclag ! model time + lag
325  real(kind=ip_r8_p) :: rcnt ! 1./cnt
326  character(len=ic_med) :: tstring ! timer label string
327  character(len=ic_med) :: fstring ! output file string
328  character(len=ic_med) :: cstring ! temporary string
329  character(len=ic_med) :: vstring ! temporary string
330  logical :: comm_now ! time to communicate
331  logical :: time_now ! coupling time
332  logical :: lreadrest ! local readrest
333  logical :: runit ! advance the variable
334  TYPE(mct_avect) :: avtest ! temporary
335  type(mct_avect) :: avtmp ! data read from restart
336  type(mct_avect) :: avtmp2 ! data read from restart
337  type(mct_avect) :: avtmp3 ! data read from restart
338  type(mct_avect) :: avtmp4 ! data read from restart
339  type(mct_avect) :: avtmp5 ! data read from restart
340  type(prism_coupler_type),pointer :: pcpointer
341  type(prism_coupler_type),pointer :: pcpointmp
342  logical, parameter :: allow_no_restart = .false.
343 !------------------------------------------------------------------
344 ! use this to set restart fields to zero if no restart file exists.
345 ! by default, code will abort if restart file is missing. this
346 ! logical avoids that error and allows the run to continue.
347 ! logical, parameter :: allow_no_restart = .true.
348 !------------------------------------------------------------------
349  character(len=*),parameter :: subname = '(oasis_advance_run)'
350  character(len=*),parameter :: F01 = '(a,i3.3)'
351 ! ----------------------------------------------------------------
352 
353  call oasis_debug_enter(subname)
354 
355  !----------------------------------------------------------------
356  !> oasis_advance_run does the following
357  !> * Aborts if it's called from non-active tasks
358  !----------------------------------------------------------------
359 
360  if (mpi_comm_local == mpi_comm_null) then
361  write(nulprt,*) subname,estr,'called on non coupling task'
362  call oasis_abort()
363  endif
364 
365  kinfo = oasis_ok
366  if (varid < 1 .or. varid > prism_nvar) then
367  write(nulprt,*) subname,estr,'invalid varid',varid
368  call oasis_abort()
369  endif
370  vname = prism_var(varid)%name
371 
372  lreadrest = .false.
373  if (present(readrest)) then
374  lreadrest = readrest
375  endif
376  if (lreadrest) kinfo = oasis_fromrest
377 
378  !------------------------------------------------
379  !> * Verify field (var) is either In or Out
380  !------------------------------------------------
381 
382  if (mop /= oasis_out .and. mop /= oasis_in) then
383  write(nulprt,*) subname,estr,'at ',msec,mseclag,' for var = ',trim(vname)
384  write(nulprt,*) subname,estr,'mop invalid expecting OASIS_Out or OASIS_In = ',mop
385  call oasis_abort()
386  endif
387 
388  call oasis_debug_note(subname//' loop over var ncpl')
389 
390  !------------------------------------------------
391  ! Check namid
392  !------------------------------------------------
393 
394  if (present(namid)) then
395  IF (oasis_debug >= 20) THEN
396  write(nulprt,*) subname,'namid',namid,prism_var(varid)%cpl(1:prism_var(varid)%ncpl)
397  endif
398  ! verify namid is associated with the variable or abort
399  runit = .false.
400  do nc = 1,prism_var(varid)%ncpl
401  cplid = prism_var(varid)%cpl(nc)
402  if (cplid == namid) runit = .true.
403  enddo
404  if (.not.runit) then
405  WRITE(nulprt,*) subname,estr,'namid not found for var = ',trim(vname)
406  WRITE(nulprt,*) subname,estr,'namid = ',namid
407  CALL oasis_abort()
408  endif
409  endif
410 
411  !------------------------------------------------
412  !> * Loop over all the couplers associated with this var
413  ! or if nam is present, just for nc == nam
414  !------------------------------------------------
415 
416  DO nc = 1,prism_var(varid)%ncpl
417  cplid = prism_var(varid)%cpl(nc)
418  runit = .true.
419  if (present(namid)) then
420  runit = .false.
421  if (cplid == namid) runit = .true.
422  endif
423  if (runit) then
424  cplid = prism_var(varid)%cpl(nc)
425  if (mop == oasis_out) pcpointer => prism_coupler_put(cplid)
426  if (mop == oasis_in ) pcpointer => prism_coupler_get(cplid)
427 
428  !------------------------------------------------
429  !> * check this prism_coupler is valid
430  !------------------------------------------------
431  if (.not.pcpointer%valid) then
432  WRITE(nulprt,*) subname,estr,'invalid prism_coupler for var = ',trim(vname)
433  CALL oasis_abort()
434  endif
435 
436  !------------------------------------------------
437  !> * check again that model op matches coupler op
438  !------------------------------------------------
439  getput = pcpointer%getput
440  if ((mop == oasis_out .and. getput == oasis3_put) .or. &
441  (mop == oasis_in .and. getput == oasis3_get)) then
442  !-continue
443  else
444  write(nulprt,*) subname,estr,'model def_var in-out does not match model get-put call for var = ',trim(vname)
445  call oasis_abort()
446  endif
447 
448  !------------------------------------------------
449  !> * set a bunch of local variables
450  !------------------------------------------------
451  rouid = pcpointer%routerid
452  mapid = pcpointer%mapperid
453  tag = pcpointer%tag
454  dt = pcpointer%dt
455  lag = pcpointer%lag
456  ltime = pcpointer%ltime
457  sndrcv = pcpointer%sndrcv
458  rstfile = trim(pcpointer%rstfile)
459  inpfile = trim(pcpointer%inpfile)
460  maxtime = pcpointer%maxtime
461  output = pcpointer%output
462  input = pcpointer%input
463  partid = pcpointer%partID
464  conserv = pcpointer%conserv
465  consbfb = .true.
466  IF (trim(pcpointer%consopt) == "opt") consbfb = .false.
467  snddiag = pcpointer%snddiag
468  rcvdiag = pcpointer%rcvdiag
469  sndadd = pcpointer%sndadd
470  sndmult = pcpointer%sndmult
471  rcvadd = pcpointer%rcvadd
472  rcvmult = pcpointer%rcvmult
473 
474  unpack = (sndrcv .OR. input)
475 
476  CALL oasis_debug_note(subname//' set nx and ny')
477  IF (prism_part(partid)%nx >= 1) THEN
478  nx = prism_part(partid)%nx
479  ny = prism_part(partid)%ny
480  ELSE
481  nx = prism_part(partid)%gsize
482  ny = 1
483  ENDIF
484 
485  IF (oasis_debug >= 20) THEN
486  WRITE(nulprt,*) subname,' DEBUG nx, ny = ',nx,ny
487  CALL oasis_flush(nulprt)
488  ENDIF
489 
490  !------------------------------------------------
491  !> * check that lag is reasonable
492  !------------------------------------------------
493 
494  IF (abs(lag) > dt) THEN
495  WRITE(nulprt,*) subname,estr,'lag setting greater than dt for var = ',trim(vname)
496  CALL oasis_abort()
497  ENDIF
498 
499  !------------------------------------------------
500  !> * read restart for call from init phase
501  ! right now, do not know whether any hot map terms are there
502  ! assume they are if something is read, otherwise not
503  !------------------------------------------------
504 
505  call oasis_debug_note(subname//' check for lag restart')
506  IF (getput == oasis3_put .AND. lag > 0 .AND. lreadrest) THEN
507  ! effective model time of restart : msec
508  mseclag = msec + lag
509  IF (.not.present(nff)) THEN
510  WRITE(nulprt,*) subname,estr,'nff optional argument not passed but expected for var = ',trim(vname)
511  CALL oasis_abort()
512  ENDIF
513  IF (len_trim(rstfile) < 1) THEN
514  WRITE(nulprt,*) subname,estr,'restart file undefined for var = ',trim(vname)
515  CALL oasis_abort()
516  ENDIF
517  lsize = mct_avect_lsize(pcpointer%aVect1)
518  IF (oasis_debug >= 2) THEN
519  WRITE(nulprt,*) subname,' at ',msec,mseclag,' RRST: ',&
520  trim(pcpointer%fldlist),' ',trim(rstfile)
521  ENDIF
522  CALL mct_avect_init(avtmp,rlist=pcpointer%fldlist,lsize=lsize)
523  if (allow_no_restart) then
524  avtmp%rAttr(nff,1:lsize) = 0.0
525  CALL oasis_io_read_avfile(trim(rstfile),avtmp,prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
526  abort=.false.,didread=a11on)
527  else
528  CALL oasis_io_read_avfile(trim(rstfile),avtmp,prism_part(partid)%pgsmap,prism_part(partid)%mpicom)
529  a11on = .true.
530  endif
531 
532  CALL mct_avect_init(avtmp2,rlist=pcpointer%fldlist,lsize=lsize)
533  CALL oasis_io_read_avfile(trim(rstfile),avtmp2,prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
534  abort=.false.,nampre='av2_',didread=a22on)
535 
536  CALL mct_avect_init(avtmp3,rlist=pcpointer%fldlist,lsize=lsize)
537  CALL oasis_io_read_avfile(trim(rstfile),avtmp3,prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
538  abort=.false.,nampre='av3_',didread=a33on)
539 
540  CALL mct_avect_init(avtmp4,rlist=pcpointer%fldlist,lsize=lsize)
541  CALL oasis_io_read_avfile(trim(rstfile),avtmp4,prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
542  abort=.false.,nampre='av4_',didread=a44on)
543 
544  CALL mct_avect_init(avtmp5,rlist=pcpointer%fldlist,lsize=lsize)
545  CALL oasis_io_read_avfile(trim(rstfile),avtmp5,prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
546  abort=.false.,nampre='av5_',didread=a55on)
547 
548  array1din(1:lsize) = avtmp%rAttr(nff,1:lsize)
549  IF (a22on) array2(1:lsize) = avtmp2%rAttr(nff,1:lsize)
550  IF (a33on) array3(1:lsize) = avtmp3%rAttr(nff,1:lsize)
551  IF (a44on) array4(1:lsize) = avtmp4%rAttr(nff,1:lsize)
552  IF (a55on) array5(1:lsize) = avtmp5%rAttr(nff,1:lsize)
553 
554  CALL mct_avect_clean(avtmp)
555  IF (a22on) CALL mct_avect_clean(avtmp2)
556  IF (a33on) CALL mct_avect_clean(avtmp3)
557  IF (a44on) CALL mct_avect_clean(avtmp4)
558  IF (a55on) CALL mct_avect_clean(avtmp5)
559  ENDIF
560 
561  !------------------------------------------------
562  !> * compute lag time, only on put side
563  !> * set time_now, is it a coupling period?
564  !------------------------------------------------
565 
566  call oasis_debug_note(subname//' set mseclag')
567  if (getput == oasis3_put) then
568  mseclag = msec + lag
569  elseif (getput == oasis3_get) then
570  mseclag = msec
571  endif
572 
573  if (oasis_debug >= 20) then
574  write(nulprt,*) subname,' DEBUG msec,mseclag = ',msec,mseclag
575  CALL oasis_flush(nulprt)
576  endif
577 
578  time_now = .false.
579  if (mod(mseclag,dt) == 0) time_now = .true.
580 
581  !------------------------------------------------
582  !> * check that model hasn't gone past maxtime
583  !------------------------------------------------
584 
585  if (msec >= maxtime) then
586  write(nulprt,*) subname,estr,'at ',msec,mseclag,' for var = ',trim(vname)
587  write(nulprt,*) subname,estr,'model time beyond namcouple maxtime = ',msec,maxtime
588  call oasis_abort()
589  endif
590 
591  !------------------------------------------------
592  !> * check that model isn't going backwards
593  ! msec >= 0 does the check only in run mode, not in initialization
594  !------------------------------------------------
595 
596  if (lcouplertime /= ispval .and. msec >= 0 .and. msec < lcouplertime) then
597  write(nulprt,*) subname,estr,'at ',msec,mseclag,' for var = ',trim(vname)
598  write(nulprt,*) subname,estr,'model seems to be running backwards = ',lcouplertime
599  call oasis_abort()
600  endif
601 
602  !------------------------------------------------
603  !> * check that variable didn't miss a coupling period
604  ! check only for send recv operations where deadlock
605  ! is possible
606  !------------------------------------------------
607 
608  do n = 1,prism_mcoupler
609  do npc = 1,2
610  if (npc == 1) pcpointmp => prism_coupler_put(n)
611  if (npc == 2) pcpointmp => prism_coupler_get(n)
612  if (pcpointmp%valid) then
613  if ((pcpointmp%ltime /= ispval) .and. &
614  (sndrcv .and. pcpointmp%sndrcv) .and. &
615  (msec > pcpointmp%ltime + pcpointmp%dt)) then
616  write(nulprt,*) subname,estr,'coupling skipped at earlier time, potential deadlock '
617  write(nulprt,*) subname,estr,'my coupler = ',cplid,' variable = ',&
618  trim(pcpointer%fldlist)
619  write(nulprt,*) subname,estr,'current time = ',msec,' mseclag = ',mseclag
620  write(nulprt,*) subname,estr,'skipped coupler = ',n,' variable = ',&
621  trim(pcpointmp%fldlist)
622  write(nulprt,*) subname,estr,'skipped coupler last time and dt = ',&
623  pcpointmp%ltime,pcpointmp%dt
624  WRITE(nulprt,*) subname,estr,'model timestep does not match coupling timestep'
625  call oasis_abort()
626  endif
627  endif ! valid
628  enddo ! npc
629  enddo ! prism_mcoupler
630 
631  !------------------------------------------------
632  !> * check that prior sequences weren't missed at this
633  !> step for get (recv) operation.
634  ! attempts to trap deadlocks before they happen
635  !------------------------------------------------
636 
637  if (sndrcv .and. getput == oasis3_get) then
638  if (lastseqtime /= ispval .and. msec == lastseqtime .and. pcpointer%seq < lastseq) then
639  write(nulprt,*) subname,estr,'coupling sequence out of order, potential deadlock '
640  write(nulprt,*) subname,estr,'my coupler = ',cplid,' variable = ',&
641  trim(pcpointer%fldlist)
642  write(nulprt,*) subname,' ERRRO: sequence number = ',pcpointer%seq
643  write(nulprt,*) subname,estr,'current time = ',msec,' mseclag = ',mseclag
644  write(nulprt,*) subname,estr,'last sequence and time = ',lastseq,lastseqtime
645  WRITE(nulprt,*) subname,estr,'model sequence does not match coupling sequence'
646  call oasis_abort()
647  else
648  lastseq = pcpointer%seq
649  lastseqtime = msec
650  endif
651  if (oasis_debug >= 20) then
652  write(nulprt,*) subname,' DEBUG sequence ',trim(vname),msec,lastseqtime,lastseq,pcpointer%seq
653  endif
654  endif
655 
656  !------------------------------------------------
657  !> * compute field index and check sizes
658  !------------------------------------------------
659 
660  call oasis_debug_note(subname//' compute field index and sizes')
661  nfav = mct_avect_indexra(pcpointer%avect1,trim(vname))
662  nsav = mct_avect_lsize(pcpointer%avect1)
663  if (lag > 0 .and. lreadrest) nsa=size(array1din)
664  if (present(array1din )) nsa = size(array1din )
665  if (present(array1dout)) nsa = size(array1dout)
666  if (present(array2dout)) nsa = size(array2dout)
667 
668  if (oasis_debug >= 20) then
669  write(nulprt,*) subname,' DEBUG nfav,nsav,nsa = ',nfav,nsav,nsa
670  endif
671 
672  if (nsav /= nsa) then
673  write(nulprt,*) subname,estr,'at ',msec,mseclag,' for var = ',trim(vname)
674  write(nulprt,*) subname,estr,'in field size passed into get/put compare to expected size ',nsav,nsa
675  call oasis_abort()
676  endif
677 
678  !------------------------------------------------
679  !> * check for higher order coupling fields
680  !> and get everything ready
681  ! arrayon is what's passed this time
682  ! optional args only on put side
683  !------------------------------------------------
684 
685  IF (.not.lreadrest) THEN
686  arrayon = .false.
687  arrayon(1) = .true.
688  if (present(a2on)) arrayon(2) = a2on
689  if (present(a3on)) arrayon(3) = a3on
690  if (present(a4on)) arrayon(4) = a4on
691  if (present(a5on)) arrayon(5) = a5on
692 
693  if ((getput == oasis3_get) .or. &
694  (getput == oasis3_put .and. trim(pcpointer%maploc) == "dst" )) then
695  if (arrayon(2) .or. arrayon(3) .or. &
696  arrayon(4) .or. arrayon(5)) then
697  write(nulprt,*) subname,estr,'at ',msec,mseclag,' for var = ',trim(vname)
698  write(nulprt,*) subname,estr,'higher order mapping not allowed on get side'
699  write(nulprt,*) subname,estr,'consider changing map location from dst to src'
700  call oasis_abort()
701  endif
702  endif
703 
704  if ((arrayon(2) .and. .not.present(array2)) .or. &
705  (arrayon(3) .and. .not.present(array3)) .or. &
706  (arrayon(4) .and. .not.present(array4)) .or. &
707  (arrayon(5) .and. .not.present(array5))) then
708  write(nulprt,*) subname,estr,'at ',msec,mseclag,' for var = ',trim(vname)
709  write(nulprt,*) subname,estr,'arrayon true but array not sent'
710  call oasis_abort()
711  ! With the current way of using oasis_advance_run, the above test is useless but we keep the test
712  ! as someone might be later adding an interface call that would violate the consistency
713  endif
714 
715  ! initialize aVect2-5 here if not already allocated
716 
717  if (arrayon(2) .and. .not. pcpointer%aVon(2)) then
718  call mct_avect_init(pcpointer%aVect2,pcpointer%aVect1,nsav)
719  call mct_avect_zero(pcpointer%aVect2)
720  pcpointer%aVon(2) = .true.
721  if (oasis_debug >= 2) then
722  write(nulprt,*) subname,' at ',msec,mseclag,' ALLO: ',&
723  trim(vname),' ','aVect2'
724  endif
725  endif
726 
727  if (arrayon(3) .and. .not. pcpointer%aVon(3)) then
728  call mct_avect_init(pcpointer%aVect3,pcpointer%aVect1,nsav)
729  call mct_avect_zero(pcpointer%aVect3)
730  pcpointer%aVon(3) = .true.
731  if (oasis_debug >= 2) then
732  write(nulprt,*) subname,' at ',msec,mseclag,' ALLO: ',&
733  trim(vname),' ','aVect3'
734  endif
735  endif
736 
737  if (arrayon(4) .and. .not. pcpointer%aVon(4)) then
738  call mct_avect_init(pcpointer%aVect4,pcpointer%aVect1,nsav)
739  call mct_avect_zero(pcpointer%aVect4)
740  pcpointer%aVon(4) = .true.
741  if (oasis_debug >= 2) then
742  write(nulprt,*) subname,' at ',msec,mseclag,' ALLO: ',&
743  trim(vname),' ','aVect4'
744  endif
745  endif
746 
747  if (arrayon(5) .and. .not. pcpointer%aVon(5)) then
748  call mct_avect_init(pcpointer%aVect5,pcpointer%aVect1,nsav)
749  call mct_avect_zero(pcpointer%aVect5)
750  pcpointer%aVon(5) = .true.
751  if (oasis_debug >= 2) then
752  write(nulprt,*) subname,' at ',msec,mseclag,' ALLO: ',&
753  trim(vname),' ','aVect5'
754  endif
755  endif
756  endif ! .not. lreadrest
757 
758  !------------------------------------------------
759  !> * update avect1-5 on put side and apply appropriate transform
760  !> * if its coupling time, set status of this var to ready
761  ! on restart, treat as instant value
762  !------------------------------------------------
763 
764  if (getput == oasis3_put) then
765 
766  call oasis_debug_note(subname//' loctrans operation')
767  write(tstring,f01) 'pcpy_',cplid
768  call oasis_timer_start(tstring)
769 
770  cstring = 'none'
771  if (lreadrest .or. pcpointer%trans == ip_instant) then
772  if (time_now) then
773  cstring = 'instant'
774  do n = 1,nsav
775  pcpointer%avect1%rAttr(nfav,n) = array1din(n)
776  if (pcpointer%aVon(2)) then
777  if (present(array2)) then
778  pcpointer%avect2%rAttr(nfav,n) = array2(n)
779  else
780  pcpointer%avect2%rAttr(nfav,n) = 0.0
781  endif
782  endif
783  if (pcpointer%aVon(3)) then
784  if (present(array3)) then
785  pcpointer%avect3%rAttr(nfav,n) = array3(n)
786  else
787  pcpointer%avect3%rAttr(nfav,n) = 0.0
788  endif
789  endif
790  if (pcpointer%aVon(4)) then
791  if (present(array4)) then
792  pcpointer%avect4%rAttr(nfav,n) = array4(n)
793  else
794  pcpointer%avect4%rAttr(nfav,n) = 0.0
795  endif
796  endif
797  if (pcpointer%aVon(5)) then
798  if (present(array5)) then
799  pcpointer%avect5%rAttr(nfav,n) = array5(n)
800  else
801  pcpointer%avect5%rAttr(nfav,n) = 0.0
802  endif
803  endif
804  enddo
805  pcpointer%avcnt(nfav) = 1
806  endif
807 
808  elseif (pcpointer%trans == ip_average) then
809  cstring = 'average'
810  if (kinfo == oasis_ok) kinfo = oasis_loctrans
811  do n = 1,nsav
812  pcpointer%avect1%rAttr(nfav,n) = &
813  pcpointer%avect1%rAttr(nfav,n) + array1din(n)
814  if (pcpointer%aVon(2)) then
815  if (present(array2)) then
816  pcpointer%avect2%rAttr(nfav,n) = &
817  pcpointer%avect2%rAttr(nfav,n) + array2(n)
818  endif
819  endif
820  if (pcpointer%aVon(3)) then
821  if (present(array3)) then
822  pcpointer%avect3%rAttr(nfav,n) = &
823  pcpointer%avect3%rAttr(nfav,n) + array3(n)
824  endif
825  endif
826  if (pcpointer%aVon(4)) then
827  if (present(array4)) then
828  pcpointer%avect4%rAttr(nfav,n) = &
829  pcpointer%avect4%rAttr(nfav,n) + array4(n)
830  endif
831  endif
832  if (pcpointer%aVon(5)) then
833  if (present(array5)) then
834  pcpointer%avect5%rAttr(nfav,n) = &
835  pcpointer%avect5%rAttr(nfav,n) + array5(n)
836  endif
837  endif
838  enddo
839  pcpointer%avcnt(nfav) = pcpointer%avcnt(nfav) + 1
840 
841  elseif (pcpointer%trans == ip_accumul) then
842  cstring = 'accumul'
843  if (kinfo == oasis_ok) kinfo = oasis_loctrans
844  do n = 1,nsav
845  pcpointer%avect1%rAttr(nfav,n) = &
846  pcpointer%avect1%rAttr(nfav,n) + array1din(n)
847  if (pcpointer%aVon(2)) then
848  if (present(array2)) then
849  pcpointer%avect2%rAttr(nfav,n) = &
850  pcpointer%avect2%rAttr(nfav,n) + array2(n)
851  endif
852  endif
853  if (pcpointer%aVon(3)) then
854  if (present(array3)) then
855  pcpointer%avect3%rAttr(nfav,n) = &
856  pcpointer%avect3%rAttr(nfav,n) + array3(n)
857  endif
858  endif
859  if (pcpointer%aVon(4)) then
860  if (present(array4)) then
861  pcpointer%avect4%rAttr(nfav,n) = &
862  pcpointer%avect4%rAttr(nfav,n) + array4(n)
863  endif
864  endif
865  if (pcpointer%aVon(5)) then
866  if (present(array5)) then
867  pcpointer%avect5%rAttr(nfav,n) = &
868  pcpointer%avect5%rAttr(nfav,n) + array5(n)
869  endif
870  endif
871  enddo
872  pcpointer%avcnt(nfav) = 1
873 
874  elseif (pcpointer%trans == ip_max) then
875  cstring = 'max'
876  if (kinfo == oasis_ok) kinfo = oasis_loctrans
877  if (pcpointer%aVon(2) .or. pcpointer%aVon(3) .or. &
878  pcpointer%aVon(4) .or. pcpointer%aVon(5)) then
879  write(nulprt,*) subname,estr,'at ',msec,mseclag,' for var = ',trim(vname)
880  write(nulprt,*) subname,estr,'higher order mapping with MAX transform not supported'
881  call oasis_abort()
882  endif
883  do n = 1,nsav
884  if (pcpointer%avcnt(nfav) == 0) then
885  pcpointer%avect1%rAttr(nfav,n) = array1din(n)
886  else
887  pcpointer%avect1%rAttr(nfav,n) = &
888  max(pcpointer%avect1%rAttr(nfav,n),array1din(n))
889  endif
890  enddo
891  pcpointer%avcnt(nfav) = 1
892 
893  elseif (pcpointer%trans == ip_min) then
894  cstring = 'min'
895  if (kinfo == oasis_ok) kinfo = oasis_loctrans
896  if (pcpointer%aVon(2) .or. pcpointer%aVon(3) .or. &
897  pcpointer%aVon(4) .or. pcpointer%aVon(5)) then
898  write(nulprt,*) subname,estr,'at ',msec,mseclag,' for var = ',trim(vname)
899  write(nulprt,*) subname,estr,'higher order mapping with MIN transform not supported'
900  call oasis_abort()
901  endif
902  do n = 1,nsav
903  if (pcpointer%avcnt(nfav) == 0) then
904  pcpointer%avect1%rAttr(nfav,n) = array1din(n)
905  else
906  pcpointer%avect1%rAttr(nfav,n) = &
907  min(pcpointer%avect1%rAttr(nfav,n),array1din(n))
908  endif
909  enddo
910  pcpointer%avcnt(nfav) = 1
911 
912  else
913  write(nulprt,*) subname,estr,'transform not known for var = ',trim(vname),pcpointer%trans
914  call oasis_abort()
915  endif
916  call oasis_timer_stop(tstring)
917 
918  if (oasis_debug >= 2 .and. trim(cstring) /= 'none') then
919  write(nulprt,*) subname,' at ',msec,mseclag,' PACK: ',&
920  trim(vname),' ',trim(cstring)
921  endif
922 
923  if (oasis_debug >= 20) then
924  write(nulprt,*) subname,' DEBUG loctrans update ',cplid,' ',&
925  trim(cstring),pcpointer%avcnt(nfav)
926  endif
927 
928  if (time_now) then
929  pcpointer%status(nfav) = oasis_comm_ready
930  endif
931  endif
932 
933  !------------------------------------------------
934  !> * decide if it's time to communicate based on time
935  ! also, on the put side, status of all vars
936  ! must be ready which means all vars have called put.
937  ! on get side, all ready means all vars have unpacked
938  ! from last get.
939  !------------------------------------------------
940 
941  call oasis_debug_note(subname//' comm_now compute')
942  comm_now = .false.
943  if (time_now) then
944  comm_now = .true.
945  do nf = 1,pcpointer%nflds
946  if (pcpointer%status(nf) /= oasis_comm_ready) then
947  comm_now = .false.
948  if (oasis_debug >= 15) then
949  write(nulprt,*) subname,' at ',msec,mseclag,' STAT: ',nf,' NOT READY'
950  endif
951  kinfo=oasis_waitgroup
952  else
953  if (oasis_debug >= 15) then
954  write(nulprt,*) subname,' at ',msec,mseclag,' STAT: ',nf,' READY'
955  endif
956  endif
957  enddo
958  endif
959 
960  !------------------------------------------------
961  !> * If it's time to communicate
962  !------------------------------------------------
963 
964  if (comm_now) then
965 
966  call oasis_debug_note(subname//' comm_now')
967 
968  !------------------------------------------------
969  !> * check again that time is correct
970  ! this is the time critical bit, we need to make sure the
971  ! model is truly advancing in time when comms are called.
972  ! must ignore the initial call, ltime = 0
973  !------------------------------------------------
974 
975  if (pcpointer%ltime /= ispval .and. msec <= pcpointer%ltime) then
976  write(nulprt,*) subname,estr,'model did not advance in time correctly for var = ',trim(vname)
977  write(nulprt,*) subname,estr,'msec, ltime = ',msec,pcpointer%ltime
978  call oasis_abort()
979  endif
980 
981  !------------------------------------------------
982  !> * average as needed for some transforms
983  ! (not cache friendly yet)
984  !------------------------------------------------
985 
986  IF (getput == oasis3_put) THEN
987  call oasis_debug_note(subname//' loctrans calc')
988  write(tstring,f01) 'pavg_',cplid
989  call oasis_timer_start(tstring)
990  do nf = 1,pcpointer%nflds
991  if (pcpointer%avcnt(nf) > 1) then
992  rcnt = 1.0/pcpointer%avcnt(nf)
993  do n = 1,nsav
994  pcpointer%avect1%rAttr(nf,n) = &
995  pcpointer%avect1%rAttr(nf,n) * rcnt
996  if (pcpointer%aVon(2)) then
997  pcpointer%avect2%rAttr(nf,n) = &
998  pcpointer%avect2%rAttr(nf,n) * rcnt
999  endif
1000  if (pcpointer%aVon(3)) then
1001  pcpointer%avect3%rAttr(nf,n) = &
1002  pcpointer%avect3%rAttr(nf,n) * rcnt
1003  endif
1004  if (pcpointer%aVon(4)) then
1005  pcpointer%avect4%rAttr(nf,n) = &
1006  pcpointer%avect4%rAttr(nf,n) * rcnt
1007  endif
1008  if (pcpointer%aVon(5)) then
1009  pcpointer%avect5%rAttr(nf,n) = &
1010  pcpointer%avect5%rAttr(nf,n) * rcnt
1011  endif
1012  enddo
1013  endif
1014  if (oasis_debug >= 20) then
1015  write(nulprt,*) subname,' DEBUG loctrans calc0 = ',cplid,nf,&
1016  pcpointer%avcnt(nf)
1017  write(nulprt,*) subname,' DEBUG loctrans calc1 = ',cplid,nf,&
1018  minval(pcpointer%avect1%rAttr(nf,:)),&
1019  maxval(pcpointer%avect1%rAttr(nf,:))
1020  call oasis_flush(nulprt)
1021  if (pcpointer%aVon(2)) &
1022  write(nulprt,*) subname,' DEBUG loctrans calc2 = ',cplid,nf,&
1023  minval(pcpointer%avect2%rAttr(nf,:)),&
1024  maxval(pcpointer%avect2%rAttr(nf,:))
1025  if (pcpointer%aVon(3)) &
1026  write(nulprt,*) subname,' DEBUG loctrans calc3 = ',cplid,nf,&
1027  minval(pcpointer%avect3%rAttr(nf,:)),&
1028  maxval(pcpointer%avect3%rAttr(nf,:))
1029  if (pcpointer%aVon(4)) &
1030  write(nulprt,*) subname,' DEBUG loctrans calc4 = ',cplid,nf,&
1031  minval(pcpointer%avect4%rAttr(nf,:)),&
1032  maxval(pcpointer%avect4%rAttr(nf,:))
1033  if (pcpointer%aVon(5)) &
1034  write(nulprt,*) subname,' DEBUG loctrans calc5 = ',cplid,nf,&
1035  minval(pcpointer%avect5%rAttr(nf,:)),&
1036  maxval(pcpointer%avect5%rAttr(nf,:))
1037  endif
1038  enddo
1039  call oasis_timer_stop(tstring)
1040  ENDIF
1041 
1042  !------------------------------------------------
1043  !> * write to restart file if put and at the end of the run,
1044  !> turn off communication
1045  ! past namcouple runtime (maxtime) no communication
1046  ! do restart if time+lag = maxtime, this assumes coupling
1047  ! period and lag and maxtime are all nicely consistent
1048  !------------------------------------------------
1049 
1050  if (mseclag >= maxtime) then
1051  sndrcv = .false. ! turn off communication
1052  unpack = .false. ! nothing to unpack
1053  if (getput == oasis3_put .and. lag > 0 .and. mseclag == maxtime) then
1054  kinfo = oasis_torest
1055  call oasis_debug_note(subname//' lag restart write')
1056  write(tstring,f01) 'wrst_',cplid
1057  call oasis_timer_start(tstring)
1058  call oasis_io_write_avfile(rstfile,pcpointer%avect1, &
1059  prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny)
1060  if (pcpointer%aVon(2)) &
1061  call oasis_io_write_avfile(rstfile,pcpointer%avect2, &
1062  prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre='av2_')
1063  if (pcpointer%aVon(3)) &
1064  call oasis_io_write_avfile(rstfile,pcpointer%avect3, &
1065  prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre='av3_')
1066  if (pcpointer%aVon(4)) &
1067  call oasis_io_write_avfile(rstfile,pcpointer%avect4, &
1068  prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre='av4_')
1069  if (pcpointer%aVon(5)) &
1070  call oasis_io_write_avfile(rstfile,pcpointer%avect5, &
1071  prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre='av5_')
1072  call oasis_timer_stop(tstring)
1073  if (oasis_debug >= 2) then
1074  write(nulprt,*) subname,' at ',msec,mseclag,' WRST: ', &
1075  trim(mct_avect_exportrlist2c(pcpointer%avect1)),' ',&
1076  trim(rstfile)
1077  call oasis_flush(nulprt)
1078  endif
1079  endif
1080  endif
1081 
1082  !------------------------------------------------
1083  !> * map and communicate operations
1084  !------------------------------------------------
1085 
1086  if (sndrcv) then
1087  if (getput == oasis3_put) then
1088  kinfo = oasis_sent
1089  call oasis_debug_note(subname//' put section')
1090  if (oasis_debug >= 2) then
1091  write(nulprt,*) subname,' at ',msec,mseclag,' SEND: ', &
1092  trim(mct_avect_exportrlist2c(pcpointer%avect1))
1093  call oasis_flush(nulprt)
1094  endif
1095  if (sndadd /= 0.0_ip_double_p .or. sndmult /= 1.0_ip_double_p) then
1096  call oasis_debug_note(subname//' apply sndmult sndadd')
1097  if (oasis_debug >= 20) then
1098  write(nulprt,*) subname,' DEBUG sndmult,add = ',sndmult,sndadd
1099  write(nulprt,*) subname,' DEBUG put b4 sndmult,add = ',cplid,&
1100  minval(pcpointer%avect1%rAttr),&
1101  maxval(pcpointer%avect1%rAttr)
1102  endif
1103  pcpointer%avect1%rAttr(:,:) = pcpointer%avect1%rAttr(:,:)*sndmult &
1104  + sndadd
1105  endif
1106  if (snddiag) call oasis_advance_avdiag(pcpointer%avect1,prism_part(partid)%mpicom)
1107  if (mapid > 0) then
1108  write(tstring,f01) 'pmap_',cplid
1109  call oasis_debug_note(subname//' put map')
1110  if (oasis_debug >= 20) then
1111  write(nulprt,*) subname,' DEBUG put av11 b4 map = ',cplid,&
1112  minval(pcpointer%avect1%rAttr),&
1113  maxval(pcpointer%avect1%rAttr)
1114  if (pcpointer%aVon(2)) &
1115  write(nulprt,*) subname,' DEBUG put av2 b4 map = ',cplid,&
1116  minval(pcpointer%avect2%rAttr),&
1117  maxval(pcpointer%avect2%rAttr)
1118  if (pcpointer%aVon(3)) &
1119  write(nulprt,*) subname,' DEBUG put av3 b4 map = ',cplid,&
1120  minval(pcpointer%avect3%rAttr),&
1121  maxval(pcpointer%avect3%rAttr)
1122  if (pcpointer%aVon(4)) &
1123  write(nulprt,*) subname,' DEBUG put av4 b4 map = ',cplid,&
1124  minval(pcpointer%avect4%rAttr),&
1125  maxval(pcpointer%avect4%rAttr)
1126  if (pcpointer%aVon(5)) &
1127  write(nulprt,*) subname,' DEBUG put av5 b4 map = ',cplid,&
1128  minval(pcpointer%avect5%rAttr),&
1129  maxval(pcpointer%avect5%rAttr)
1130  endif
1131  call oasis_timer_start(tstring)
1132  if (lucia_debug > 0) &
1133  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1134  'Balance: ',pcpointer%namID,' Before interpo ', mpi_wtime()
1135  call mct_avect_zero(pcpointer%avect1m)
1136  call oasis_advance_map(pcpointer%avect1, &
1137  pcpointer%avect1m,prism_mapper(mapid),conserv,consbfb, &
1138  pcpointer%aVon ,pcpointer%avect2, &
1139  pcpointer%avect3,pcpointer%avect4, &
1140  pcpointer%avect5)
1141  if (lucia_debug > 0) &
1142  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1143  'Balance: ',pcpointer%namID,' After interpo ', mpi_wtime()
1144  call oasis_timer_stop(tstring)
1145  write(tstring,f01) 'psnd_',cplid
1146  call oasis_debug_note(subname//' put send')
1147  if (oasis_debug >= 20) then
1148  write(nulprt,*) subname,' DEBUG put av1m b4 send = ',cplid,&
1149  minval(pcpointer%avect1m%rAttr),&
1150  maxval(pcpointer%avect1m%rAttr)
1151  endif
1152  call oasis_timer_start(tstring)
1153  if (lucia_debug > 0) &
1154  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1155  'Balance: ',pcpointer%namID,' Before MPI put ', mpi_wtime()
1156  call mct_waitsend(prism_router(rouid)%router)
1157  call mct_isend(pcpointer%avect1m,prism_router(rouid)%router,tag)
1158  if (lucia_debug > 0) &
1159  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1160  'Balance: ',pcpointer%namID,' After MPI put ', mpi_wtime()
1161  call oasis_timer_stop(tstring)
1162  ELSE
1163  write(tstring,f01) 'psnd_',cplid
1164  call oasis_debug_note(subname//' put send')
1165  if (oasis_debug >= 20) then
1166  write(nulprt,*) subname,' DEBUG put av1 b4 send = ',cplid,&
1167  minval(pcpointer%avect1%rAttr),&
1168  maxval(pcpointer%avect1%rAttr)
1169  endif
1170  call oasis_timer_start(tstring)
1171  if (lucia_debug > 0) &
1172  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1173  'Balance: ',pcpointer%namID,' Before MPI put ', mpi_wtime()
1174  call mct_waitsend(prism_router(rouid)%router)
1175  call mct_isend(pcpointer%avect1,prism_router(rouid)%router,tag)
1176  if (lucia_debug > 0) &
1177  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1178  'Balance: ',pcpointer%namID,' After MPI put ', mpi_wtime()
1179  call oasis_timer_stop(tstring)
1180  ENDIF
1181  elseif (getput == oasis3_get) then
1182  call oasis_debug_note(subname//' get section')
1183  if (oasis_debug >= 2 ) then
1184  write(nulprt,*) subname,' at ',msec,mseclag,' RECV: ', &
1185  trim(mct_avect_exportrlist2c(pcpointer%avect1))
1186  call oasis_flush(nulprt)
1187  endif
1188  if (mapid > 0) then
1189  call oasis_debug_note(subname//' get recv')
1190  write(tstring,f01) 'grcv_',cplid
1191  call oasis_timer_start(tstring)
1192  call mct_avect_zero(pcpointer%avect1m)
1193  if (lucia_debug > 0) &
1194  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1195  'Balance: ',pcpointer%namID,' Before MPI get ', mpi_wtime()
1196  call mct_recv(pcpointer%avect1m,prism_router(rouid)%router,tag)
1197  if (lucia_debug > 0) &
1198  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1199  'Balance: ',pcpointer%namID,' After MPI get ', mpi_wtime()
1200  call oasis_timer_stop(tstring)
1201  if (oasis_debug >= 20) then
1202  write(nulprt,*) subname,' DEBUG get af recv = ',cplid,&
1203  minval(pcpointer%avect1m%rAttr),&
1204  maxval(pcpointer%avect1m%rAttr)
1205  endif
1206  call oasis_debug_note(subname//' get map')
1207  write(tstring,f01) 'gmap_',cplid
1208  call oasis_timer_start(tstring)
1209  if (lucia_debug > 0) &
1210  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1211  'Balance: ',pcpointer%namID,' Before interpo ', mpi_wtime()
1212  call mct_avect_zero(pcpointer%avect1)
1213  call oasis_advance_map(pcpointer%avect1m, &
1214  pcpointer%avect1,prism_mapper(mapid),conserv,consbfb)
1215  if (lucia_debug > 0) &
1216  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1217  'Balance: ',pcpointer%namID,' After interpo ', mpi_wtime()
1218  call oasis_timer_stop(tstring)
1219  if (oasis_debug >= 20) then
1220  write(nulprt,*) subname,' DEBUG get af map = ',cplid,&
1221  minval(pcpointer%avect1%rAttr),&
1222  maxval(pcpointer%avect1%rAttr)
1223  endif
1224  else
1225  write(tstring,f01) 'grcv_',cplid
1226  call oasis_debug_note(subname//' get recv')
1227  call oasis_timer_start(tstring)
1228  if (lucia_debug > 0) &
1229  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1230  'Balance: ',pcpointer%namID,' Before MPI get ', mpi_wtime()
1231  call mct_recv(pcpointer%avect1,prism_router(rouid)%router,tag)
1232  if (lucia_debug > 0) &
1233  WRITE(nullucia, fmt='(A,I3.3,A,F16.5)') &
1234  'Balance: ',pcpointer%namID,' After MPI get ', mpi_wtime()
1235  call oasis_timer_stop(tstring)
1236  if (oasis_debug >= 20) then
1237  write(nulprt,*) subname,' DEBUG get af recv = ',cplid,&
1238  minval(pcpointer%avect1%rAttr),&
1239  maxval(pcpointer%avect1%rAttr)
1240  endif
1241  endif
1242  call oasis_debug_note(subname//' apply rcvmult rcvadd')
1243  if (rcvadd /= 0.0_ip_double_p .or. rcvmult /= 1.0_ip_double_p) then
1244  pcpointer%avect1%rAttr(:,:) = pcpointer%avect1%rAttr(:,:)*rcvmult &
1245  + rcvadd
1246  if (oasis_debug >= 20) then
1247  write(nulprt,*) subname,' DEBUG rcvmult,add = ',rcvmult,rcvadd
1248  write(nulprt,*) subname,' DEBUG get af rcvmult,add = ',cplid,&
1249  minval(pcpointer%avect1%rAttr),&
1250  maxval(pcpointer%avect1%rAttr)
1251  endif
1252  endif
1253  if (rcvdiag) call oasis_advance_avdiag(pcpointer%avect1,prism_part(partid)%mpicom)
1254  endif ! getput
1255  endif ! sndrcv
1256 
1257  !------------------------------------------------
1258  !> * write to output files if output is turned on
1259  !------------------------------------------------
1260 
1261  if (output) then
1262  if (kinfo == oasis_sent) then
1263  kinfo = oasis_sentout
1264  elseif (kinfo == oasis_torest) then
1265  kinfo = oasis_torestout
1266  else
1267  kinfo = oasis_output
1268  endif
1269  call oasis_debug_note(subname//' output')
1270  write(tstring,f01) 'wout_',cplid
1271  call oasis_timer_start(tstring)
1272  if (oasis_debug >= 2) then
1273  write(nulprt,*) subname,' at ',msec,mseclag,' WRIT: ', &
1274  trim(mct_avect_exportrlist2c(pcpointer%avect1))
1275  call oasis_flush(nulprt)
1276  endif
1277  write(fstring,'(A,I2.2)') '_'//trim(compnm)//'_',cplid
1278  call oasis_io_write_avfbf(pcpointer%avect1,prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
1279  nx,ny,msec,fstring)
1280  call oasis_timer_stop(tstring)
1281 
1282  if (oasis_debug >= 30) then
1283  call mct_avect_init(avtest,pcpointer%avect1,&
1284  mct_avect_lsize(pcpointer%avect1))
1285  write(tstring,f01) 'rinp_',cplid
1286  call oasis_timer_start(tstring)
1287  call oasis_io_read_avfbf(avtest,prism_part(partid)%pgsmap,prism_part(partid)%mpicom,msec,fstring)
1288  write(nulprt,*) subname,' DEBUG write/read test avfbf should be zero ',&
1289  sum(pcpointer%avect1%rAttr-avtest%rAttr)
1290  call mct_avect_clean(avtest)
1291  call oasis_timer_stop(tstring)
1292  endif
1293 
1294  endif
1295 
1296  !------------------------------------------------
1297  !> * set avcnt, avect1, ltime, and status
1298  !------------------------------------------------
1299 
1300  call oasis_debug_note(subname//' reset status')
1301  if (getput == oasis3_put) then
1302  if (.not.lreadrest) pcpointer%ltime = msec
1303  pcpointer%status(:) = oasis_comm_wait
1304  pcpointer%avcnt(:) = 0
1305  call mct_avect_zero(pcpointer%avect1)
1306  if (pcpointer%aVon(2)) &
1307  call mct_avect_zero(pcpointer%avect2)
1308  if (pcpointer%aVon(3)) &
1309  call mct_avect_zero(pcpointer%avect3)
1310  if (pcpointer%aVon(4)) &
1311  call mct_avect_zero(pcpointer%avect4)
1312  if (pcpointer%aVon(5)) &
1313  call mct_avect_zero(pcpointer%avect5)
1314  if (oasis_debug >= 20) then
1315  write(nulprt,*) subname,' DEBUG put reset status = '
1316  endif
1317  elseif (getput == oasis3_get) then
1318  if (.not.lreadrest) pcpointer%ltime = msec
1319  pcpointer%status(:) = oasis_comm_wait
1320  if (oasis_debug >= 20) then
1321  write(nulprt,*) subname,' DEBUG get reset status = '
1322  endif
1323  endif
1324 
1325  else
1326 
1327  !------------------------------------------------
1328  ! no action, document
1329  !------------------------------------------------
1330 
1331  if (oasis_debug >= 15) then
1332  if (getput == oasis3_put) then
1333  write(nulprt,*) subname,' at ',msec,mseclag,' SKIP: ', &
1334  trim(mct_avect_exportrlist2c(pcpointer%avect1))
1335  elseif (getput == oasis3_get) then
1336  write(nulprt,*) subname,' at ',msec,mseclag,' SKIP: ', &
1337  trim(mct_avect_exportrlist2c(pcpointer%avect1))
1338  endif
1339  call oasis_flush(nulprt)
1340  endif
1341 
1342  endif ! comm_now
1343 
1344  !------------------------------------------------
1345  !> * at the end of the run only, save fields associated
1346  !> with non-instant loctrans operations to restart files
1347  !------------------------------------------------
1348 
1349  IF (mseclag + dt >= maxtime .AND. &
1350  getput == oasis3_put .and. pcpointer%trans /= ip_instant) then
1351  call oasis_debug_note(subname//' loctrans restart write')
1352  write(tstring,f01) 'wtrn_',cplid
1353  call oasis_timer_start(tstring)
1354  WRITE(vstring,'(a,i6.6,a)') 'loc',pcpointer%namID,'_cnt'
1355  CALL oasis_io_write_array(rstfile,prism_part(partid)%mpicom,iarray=pcpointer%avcnt,&
1356  ivarname=trim(vstring))
1357  write(vstring,'(a,i6.6,a)') 'loc',pcpointer%namID,'_'
1358  CALL oasis_io_write_avfile(rstfile,pcpointer%avect1, &
1359  prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=trim(vstring))
1360  if (pcpointer%aVon(2)) then
1361  write(vstring,'(a,i6.6,a)') 'av2loc',pcpointer%namID,'_'
1362  CALL oasis_io_write_avfile(rstfile,pcpointer%avect2, &
1363  prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=trim(vstring))
1364  endif
1365  if (pcpointer%aVon(3)) then
1366  write(vstring,'(a,i6.6,a)') 'av3loc',pcpointer%namID,'_'
1367  CALL oasis_io_write_avfile(rstfile,pcpointer%avect3, &
1368  prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=trim(vstring))
1369  endif
1370  if (pcpointer%aVon(4)) then
1371  write(vstring,'(a,i6.6,a)') 'av4loc',pcpointer%namID,'_'
1372  CALL oasis_io_write_avfile(rstfile,pcpointer%avect4, &
1373  prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=trim(vstring))
1374  endif
1375  if (pcpointer%aVon(5)) then
1376  write(vstring,'(a,i6.6,a)') 'av5loc',pcpointer%namID,'_'
1377  CALL oasis_io_write_avfile(rstfile,pcpointer%avect5, &
1378  prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=trim(vstring))
1379  endif
1380  call oasis_timer_stop(tstring)
1381  if (oasis_debug >= 2) then
1382  write(nulprt,*) subname,' at ',msec,mseclag,' WTRN: ', &
1383  trim(mct_avect_exportrlist2c(pcpointer%avect1)),' ',trim(rstfile)
1384  call oasis_flush(nulprt)
1385  endif
1386  if (oasis_debug >= 20) then
1387  write(nulprt,*) subname,' DEBUG write loctrans restart',cplid,&
1388  pcpointer%avcnt
1389  write(nulprt,*) subname,' DEBUG write loctrans restart',cplid,&
1390  minval(pcpointer%avect1%rAttr),&
1391  maxval(pcpointer%avect1%rAttr)
1392  endif
1393  ENDIF
1394 
1395  !------------------------------------------------
1396  !> * GET only, unpack avect1 if it's newly received
1397  !------------------------------------------------
1398 
1399  if (getput == oasis3_get) then
1400  IF (time_now .AND. unpack) THEN
1401  if (kinfo == oasis_output) then
1402  kinfo = oasis_recvout
1403  elseif (kinfo == oasis_fromrest) then
1404  kinfo = oasis_fromrestout
1405  else
1406  kinfo = oasis_recvd
1407  endif
1408  if (input) then
1409  kinfo = oasis_input
1410  call oasis_debug_note(subname//' input')
1411  if (oasis_debug >= 2) then
1412  write(nulprt,*) subname,' at ',msec,mseclag,' READ: ', &
1413  trim(mct_avect_exportrlist2c(pcpointer%avect1))
1414  call oasis_flush(nulprt)
1415  endif
1416  write(tstring,f01) 'grin_',cplid
1417  call oasis_timer_start(tstring)
1418  if (trim(inpfile) /= trim(cspval)) then
1419  call oasis_io_read_avfbf(pcpointer%avect1,prism_part(partid)%pgsmap,prism_part(partid)%mpicom,&
1420  msec,filename=trim(inpfile))
1421  else
1422  fstring = '_'//trim(compnm)
1423  call oasis_io_read_avfbf(pcpointer%avect1,prism_part(partid)%pgsmap,prism_part(partid)%mpicom,&
1424  msec,f_string=fstring)
1425  endif
1426  call oasis_timer_stop(tstring)
1427  endif
1428  if (oasis_debug >= 2) then
1429  write(nulprt,*) subname,' at ',msec,mseclag,' UPCK: ',trim(vname)
1430  endif
1431  write(tstring,f01) 'gcpy_',cplid
1432  call oasis_debug_note(subname//' get copy to array')
1433  call oasis_timer_start(tstring)
1434  if (present(array1dout)) array1dout(:) = &
1435  pcpointer%avect1%rAttr(nfav,:)
1436  if (present(array2dout)) array2dout(:,:) = &
1437  reshape(pcpointer%avect1%rAttr(nfav,:),shape(array2dout))
1438  call oasis_timer_stop(tstring)
1439  if (oasis_debug >= 20) then
1440  if (present(array1dout)) write(nulprt,*) subname,' DEBUG array copy = ',&
1441  cplid,minval(array1dout),maxval(array1dout)
1442  if (present(array2dout)) write(nulprt,*) subname,' DEBUG array copy = ',&
1443  cplid,minval(array2dout),maxval(array2dout)
1444  endif
1445  ENDIF
1446  if (time_now) pcpointer%status(nfav) = oasis_comm_ready
1447  endif
1448 
1449  !------------------------------------------------
1450  !> * always remember last id and last coupler time
1451  !------------------------------------------------
1452 
1453  lcouplerid = cplid
1454  lcouplertime = msec
1455 
1456  if (oasis_debug >= 2) then
1457  write(nulprt,*) subname,' at ',msec,mseclag,' KINF: ',trim(vname),kinfo
1458  endif
1459 
1460  endif ! runit
1461  enddo ! nc = 1,var%ncpl
1462 
1463  call oasis_debug_exit(subname)
1464 
1465  END SUBROUTINE oasis_advance_run
1466 
1467 
1468 !-------------------------------------------------------------------
1469 
1470 !> Provides interpolation functionality
1471 
1472 !> Maps (regrids, interpolates) data from av1 to avd.
1473 !> av2-av5 are for higher order mapping (hot).
1474 
1475  SUBROUTINE oasis_advance_map(av1,avd,mapper,conserv,consbfb,&
1476  avon,av2,av3,av4,av5)
1477 
1478  ! NOTE: mask = 0 is active point according to oasis3 conserv.f
1479 
1480  implicit none
1481  type(mct_avect) ,intent(in) :: av1 !< source av
1482  type(mct_avect) ,intent(inout) :: avd !< dst av
1483  type(prism_mapper_type),intent(inout) :: mapper !< prism_mapper
1484  integer(kind=ip_i4_p) ,intent(in),optional :: conserv !< conserv flag
1485  logical ,intent(in),optional :: consbfb !< conserv bfb option
1486  logical ,intent(in),optional :: avon(:) !< which source hot are on
1487  type(mct_avect) ,intent(in),optional :: av2 !< source av2 hot
1488  type(mct_avect) ,intent(in),optional :: av3 !< source av3 hot
1489  type(mct_avect) ,intent(in),optional :: av4 !< source av4 hot
1490  type(mct_avect) ,intent(in),optional :: av5 !< source av5 hot
1491 
1492  integer(kind=ip_i4_p) :: fsize,lsizes,lsized,nf,ni,n,m
1493  real(kind=ip_r8_p) :: sumtmp, wts_sums, wts_sumd, zradi, zlagr
1494  integer(kind=ip_i4_p),allocatable :: imasks(:),imaskd(:)
1495  real(kind=ip_r8_p),allocatable :: areas(:),aread(:)
1496  real(kind=ip_r8_p),allocatable :: av_sums(:),av_sumd(:) ! local sums
1497  character(len=ic_med) :: tstring ! timer label string
1498  type(mct_avect) :: avdtmp ! for summing multiple mapping weights
1499  type(mct_avect) :: av2g ! for bfb sums
1500  logical :: lconsbfb
1501  integer(kind=ip_i4_p),parameter :: avsmax = prism_coupler_avsmax
1502  logical :: locavon(avsmax) ! local avon
1503  integer(kind=ip_i4_p) :: avonsize, nterm
1504  character(len=*),parameter :: subname = '(oasis_advance_map)'
1505 
1506  call oasis_debug_enter(subname)
1507 
1508  !> oasis_advance_map does the following
1509  !> * check for conservation flags
1510 
1511  lconsbfb = .true.
1512  if (present(consbfb)) then
1513  lconsbfb = consbfb
1514  endif
1515 
1516  !> * check for higher order terms
1517  !--- assume avon and av2-5 are not passed but av1 always is ---
1518  avonsize = 1
1519  nterm = 1
1520  locavon = .false.
1521  locavon(1) = .true.
1522 
1523  !--- but if avon is passed, use avon flags ---
1524  if (present(avon)) then
1525  avonsize = size(avon)
1526  nterm = min(avsmax,avonsize)
1527  locavon(1:nterm) = avon(1:nterm)
1528  else
1529  !--- if avon is not passed, av2-5 should not be ---
1530  if (present(av2) .or. present(av3) .or. present(av4) .or. present(av5)) then
1531  WRITE(nulprt,*) subname,estr,'av2-5 passed but avon not passed'
1532  CALL oasis_abort()
1533  endif
1534  endif
1535 
1536  !> * check consistency between weights and coupling terms
1537  do n = 1,nterm
1538  if (locavon(n) .and. n > mapper%nwgts) then
1539  WRITE(nulprt,*) subname,estr,'in nwgts and coupling terms',mapper%nwgts,n
1540  CALL oasis_abort()
1541  endif
1542  enddo
1543 
1544  !> * run mct sparse matrix mapper on data and separately on hot as needed
1545 
1546  if (locavon(1)) then
1547  if (mct_avect_nrattr(av1) /= mct_avect_nrattr(avd)) then
1548  WRITE(nulprt,*) subname,estr,'in av1 num of flds'
1549  CALL oasis_abort()
1550  endif
1551  call mct_smat_avmult(av1, mapper%sMatP(1), avd)
1552  endif
1553 
1554  if (locavon(2).or.locavon(3).or.locavon(4).or.locavon(5)) then
1555  lsized = mct_avect_lsize(avd)
1556  call mct_avect_init(avdtmp,avd,lsized)
1557 
1558  if (locavon(2)) then
1559  if (mct_avect_nrattr(av2) /= mct_avect_nrattr(avd)) then
1560  WRITE(nulprt,*) subname,estr,'in av2 num of flds'
1561  CALL oasis_abort()
1562  endif
1563  call mct_smat_avmult(av2, mapper%sMatP(2), avdtmp)
1564  avd%rAttr = avd%rAttr + avdtmp%rAttr
1565  endif
1566 
1567  if (locavon(3)) then
1568  if (mct_avect_nrattr(av3) /= mct_avect_nrattr(avd)) then
1569  WRITE(nulprt,*) subname,estr,'in av3 num of flds'
1570  CALL oasis_abort()
1571  endif
1572  call mct_smat_avmult(av3, mapper%sMatP(3), avdtmp)
1573  avd%rAttr = avd%rAttr + avdtmp%rAttr
1574  endif
1575 
1576  if (locavon(4)) then
1577  if (mct_avect_nrattr(av4) /= mct_avect_nrattr(avd)) then
1578  WRITE(nulprt,*) subname,estr,'in av4 num of flds'
1579  CALL oasis_abort()
1580  endif
1581  call mct_smat_avmult(av4, mapper%sMatP(4), avdtmp)
1582  avd%rAttr = avd%rAttr + avdtmp%rAttr
1583  endif
1584 
1585  if (locavon(5)) then
1586  if (mct_avect_nrattr(av5) /= mct_avect_nrattr(avd)) then
1587  WRITE(nulprt,*) subname,estr,'in av5 num of flds'
1588  CALL oasis_abort()
1589  endif
1590  call mct_smat_avmult(av5, mapper%sMatP(5), avdtmp)
1591  avd%rAttr = avd%rAttr + avdtmp%rAttr
1592  endif
1593 
1594  call mct_avect_clean(avdtmp)
1595  endif
1596 
1597 
1598  call oasis_debug_note(subname//' map')
1599 
1600  !> * enforce conservation
1601 
1602  IF (PRESENT(conserv)) THEN
1603  call oasis_debug_note(subname//' conserv')
1604  IF (conserv /= ip_cnone) THEN
1605  fsize = mct_avect_nrattr(av1)
1606  allocate(av_sums(fsize),av_sumd(fsize))
1607 
1608  zradi = 1./(eradius*eradius)
1609 
1610  !-------------------
1611  ! extract mask and area and compute sum of masked area for source
1612  !-------------------
1613  lsizes = mct_avect_lsize(mapper%av_ms)
1614  allocate(imasks(lsizes),areas(lsizes))
1615  nf = mct_avect_indexia(mapper%av_ms,'mask')
1616  imasks(:) = mapper%av_ms%iAttr(nf,:)
1617  nf = mct_avect_indexra(mapper%av_ms,'area')
1618  areas(:) = mapper%av_ms%rAttr(nf,:)*zradi
1619 
1620  if (lconsbfb) then
1621  call mct_avect_gather(mapper%av_ms,av2g,prism_part(mapper%spart)%gsmap,&
1622  0,mpi_comm_local)
1623  wts_sums = 0.0_ip_r8_p
1624  if (mpi_rank_local == 0) then
1625  ni = mct_avect_indexia(av2g,'mask')
1626  nf = mct_avect_indexra(av2g,'area')
1627  do n = 1,mct_avect_lsize(av2g)
1628  if (av2g%iAttr(ni,n) == 0) wts_sums = wts_sums + av2g%rAttr(nf,n)*zradi
1629  enddo
1630  endif
1631  call oasis_mpi_bcast(wts_sums,mpi_comm_local,subname//" bcast wts_sums")
1632  if (mpi_rank_local == 0) then
1633  call mct_avect_clean(av2g)
1634  endif
1635  else
1636  sumtmp = 0.0_ip_r8_p
1637  do n = 1,lsizes
1638  if (imasks(n) == 0) sumtmp = sumtmp + areas(n)
1639  enddo
1640  call oasis_mpi_sum(sumtmp,wts_sums,mpi_comm_local,string=subname//':wts_sums',&
1641  all=.true.)
1642  endif
1643 
1644  !-------------------
1645  ! extract mask and area and compute sum of masked area for destination
1646  !-------------------
1647  lsized = mct_avect_lsize(mapper%av_md)
1648  allocate(imaskd(lsized),aread(lsized))
1649  nf = mct_avect_indexia(mapper%av_md,'mask')
1650  imaskd(:) = mapper%av_md%iAttr(nf,:)
1651  nf = mct_avect_indexra(mapper%av_md,'area')
1652  aread(:) = mapper%av_md%rAttr(nf,:)*zradi
1653 
1654  if (lconsbfb) then
1655  call mct_avect_gather(mapper%av_md,av2g,prism_part(mapper%dpart)%gsmap,0,mpi_comm_local)
1656  wts_sumd = 0.0_ip_r8_p
1657  if (mpi_rank_local == 0) then
1658  ni = mct_avect_indexia(av2g,'mask')
1659  nf = mct_avect_indexra(av2g,'area')
1660  do n = 1,mct_avect_lsize(av2g)
1661  if (av2g%iAttr(ni,n) == 0) wts_sumd = wts_sumd + av2g%rAttr(nf,n)*zradi
1662  enddo
1663  endif
1664  call oasis_mpi_bcast(wts_sumd,mpi_comm_local,subname//" bcast wts_sumd")
1665  if (mpi_rank_local == 0) then
1666  call mct_avect_clean(av2g)
1667  endif
1668  else
1669  sumtmp = 0.0_ip_r8_p
1670  do n = 1,lsized
1671  if (imaskd(n) == 0) sumtmp = sumtmp + aread(n)
1672  enddo
1673  call oasis_mpi_sum(sumtmp,wts_sumd,mpi_comm_local,string=subname//':wts_sumd',all=.true.)
1674  endif
1675 
1676  if (oasis_debug >= 30) then
1677  write(nulprt,*) subname,' DEBUG conserve src mask ',minval(imasks),&
1678  maxval(imasks),sum(imasks)
1679  write(nulprt,*) subname,' DEBUG conserve dst mask ',minval(imaskd),&
1680  maxval(imaskd),sum(imaskd)
1681  write(nulprt,*) subname,' DEBUG conserve src area ',minval(areas),&
1682  maxval(areas),sum(areas)
1683  write(nulprt,*) subname,' DEBUG conserve dst area ',minval(aread),&
1684  maxval(aread),sum(aread)
1685  write(nulprt,*) subname,' DEBUG conserve wts_sum ',wts_sums,wts_sumd
1686  endif
1687 
1688  !-------------------
1689  ! compute global sums of av1
1690  ! assume av1 is the thing to be conserved
1691  !-------------------
1692  call oasis_advance_avsum(av1,av_sums,prism_part(mapper%spart)%gsmap,mpi_comm_local, &
1693  mask=imasks,wts=areas,consbfb=lconsbfb)
1694  call oasis_advance_avsum(avd,av_sumd,prism_part(mapper%dpart)%gsmap,mpi_comm_local, &
1695  mask=imaskd,wts=aread,consbfb=lconsbfb)
1696 
1697  if (oasis_debug >= 20) then
1698  write(nulprt,*) subname,' DEBUG src sum b4 conserve ',av_sums
1699  write(nulprt,*) subname,' DEBUG dst sum b4 conserve ',av_sumd
1700  endif
1701 
1702  if (conserv == ip_cglobal) then
1703  if (wts_sumd == 0.0_ip_r8_p) then
1704  WRITE(nulprt,*) subname,estr,'global conserve sums to zero '
1705  CALL oasis_abort()
1706  endif
1707  do m = 1,fsize
1708  zlagr = (av_sumd(m) - av_sums(m)) / wts_sumd
1709  do n = 1,lsized
1710  if (imaskd(n) == 0) avd%rAttr(m,n) = avd%rAttr(m,n) - zlagr
1711  enddo
1712  enddo
1713  elseif (conserv == ip_cglbpos) then
1714  do m = 1,fsize
1715  if (av_sumd(m) == 0.0_ip_r8_p .and. av_sums(m) /= 0.0_ip_r8_p) then
1716  WRITE(nulprt,*) subname,estr,'cglpos conserve one of the sums is zero'
1717  CALL oasis_abort()
1718  elseif (av_sumd(m) /= 0.0_ip_r8_p) then
1719  zlagr = av_sums(m) / av_sumd(m)
1720  do n = 1,lsized
1721  if (imaskd(n) == 0) avd%rAttr(m,n) = avd%rAttr(m,n) * zlagr
1722  enddo
1723  endif
1724  enddo
1725  elseif (conserv == ip_cbasbal) then
1726  if (wts_sumd == 0.0_ip_r8_p .or. wts_sums == 0.0_ip_r8_p) then
1727  WRITE(nulprt,*) subname,estr,'cbasbal conserve both sums are zero'
1728  CALL oasis_abort()
1729  endif
1730  do m = 1,fsize
1731  zlagr = (av_sumd(m) - (av_sums(m)*(wts_sumd/wts_sums))) / wts_sumd
1732  do n = 1,lsized
1733  if (imaskd(n) == 0) avd%rAttr(m,n) = avd%rAttr(m,n) - zlagr
1734  enddo
1735  enddo
1736  elseif (conserv == ip_cbaspos) then
1737  do m = 1,fsize
1738  if (av_sumd(m) == 0.0_ip_r8_p .and. av_sums(m) /= 0.0_ip_r8_p) then
1739  WRITE(nulprt,*) subname,estr,'cbaspos conserve one of the sums is zero'
1740  CALL oasis_abort()
1741  elseif (av_sumd(m) /= 0.0_ip_r8_p) then
1742  zlagr = (av_sums(m)/av_sumd(m)) * (wts_sumd/wts_sums)
1743  do n = 1,lsized
1744  if (imaskd(n) == 0) avd%rAttr(m,n) = avd%rAttr(m,n) * zlagr
1745  enddo
1746  endif
1747  enddo
1748  else
1749  WRITE(nulprt,*) subname,estr,'conserv option unknown = ',conserv
1750  CALL oasis_abort()
1751  endif
1752 
1753  if (oasis_debug >= 20) then
1754  call oasis_advance_avsum(av1,av_sums,prism_part(mapper%spart)%gsmap,mpi_comm_local, &
1755  mask=imasks,wts=areas,consbfb=lconsbfb)
1756  call oasis_advance_avsum(avd,av_sumd,prism_part(mapper%dpart)%gsmap,mpi_comm_local, &
1757  mask=imaskd,wts=aread,consbfb=lconsbfb)
1758  write(nulprt,*) subname,' DEBUG src sum af conserve ',av_sums
1759  write(nulprt,*) subname,' DEBUG dst sum af conserve ',av_sumd
1760  CALL oasis_flush(nulprt)
1761  endif
1762 
1763  deallocate(imasks,imaskd,areas,aread)
1764  deallocate(av_sums,av_sumd)
1765  ENDIF
1766  ENDIF ! present conserve
1767 
1768  call oasis_debug_exit(subname)
1769 
1770  END SUBROUTINE oasis_advance_map
1771 
1772 !-------------------------------------------------------------------
1773 
1774 !> A generic method for summing fields in an attribute vector
1775 
1776  SUBROUTINE oasis_advance_avsum(av,sum,gsmap,mpicom,mask,wts,consbfb)
1777 
1778  implicit none
1779  type(mct_avect) ,intent(in) :: av ! input av
1780  real(kind=ip_r8_p) ,intent(inout) :: sum(:) ! sum of av fields
1781  type(mct_gsmap) ,intent(in) :: gsmap ! gsmap associate with av
1782  integer(kind=ip_i4_p),intent(in) :: mpicom ! mpicom
1783  integer(kind=ip_i4_p),intent(in),optional :: mask(:) ! mask to apply to av
1784  real(kind=ip_r8_p) ,intent(in),optional :: wts(:) ! wts to apply to av
1785  logical ,intent(in),optional :: consbfb ! bfb conserve
1786 
1787  integer(kind=ip_i4_p) :: n,m,ierr,mytask
1788  integer(kind=ip_i4_p) :: lsize,fsize ! local size of av, number of flds in av
1789  real(kind=ip_r8_p),allocatable :: lsum(:) ! local sums
1790  real(kind=ip_r8_p),allocatable :: lwts(:) ! local wts taking into account mask and wts
1791  type(mct_avect) :: av1, av1g ! use av1,av1g for gather and bfb sum
1792  logical :: lconsbfb ! local conserve bfb
1793  character(len=*),parameter :: subname = '(oasis_advance_avsum)'
1794 
1795  call oasis_debug_enter(subname)
1796 
1797  lconsbfb = .true.
1798  if (present(consbfb)) then
1799  lconsbfb = consbfb
1800  endif
1801 
1802  fsize = mct_avect_nrattr(av)
1803  lsize = mct_avect_lsize(av)
1804 
1805  allocate(lsum(fsize))
1806  lsum = 0.0_ip_r8_p
1807  allocate(lwts(lsize))
1808  lwts = 1.0_ip_r8_p
1809 
1810  if (size(sum) /= fsize) then
1811  WRITE(nulprt,*) subname,estr,'size sum ne size av'
1812  CALL oasis_abort()
1813  endif
1814 
1815  if (present(mask)) then
1816  if (size(mask) /= lsize) then
1817  WRITE(nulprt,*) subname,estr,'size mask ne size av'
1818  CALL oasis_abort()
1819  endif
1820  do n = 1,lsize
1821  if (mask(n) /= 0) lwts(n) = 0.0_ip_r8_p
1822  enddo
1823  endif
1824 
1825  if (present(wts)) then
1826  if (size(wts) /= lsize) then
1827  WRITE(nulprt,*) subname,estr,'size wts ne size av'
1828  CALL oasis_abort()
1829  endif
1830  do n = 1,lsize
1831  lwts(n) = lwts(n) * wts(n)
1832  enddo
1833  endif
1834 
1835  if (lconsbfb) then
1836  call mct_avect_init(av1,av,lsize)
1837  do n = 1,lsize
1838  do m = 1,fsize
1839  av1%rAttr(m,n) = av%rAttr(m,n)*lwts(n)
1840  enddo
1841  enddo
1842  call mct_avect_gather(av1,av1g,gsmap,0,mpicom)
1843  call mpi_comm_rank(mpicom,mytask,ierr)
1844  sum = 0.0_ip_r8_p
1845  if (mytask == 0) then
1846  do n = 1,mct_avect_lsize(av1g)
1847  do m = 1,fsize
1848  sum(m) = sum(m) + av1g%rAttr(m,n)
1849  enddo
1850  enddo
1851  endif
1852  call oasis_mpi_bcast(sum,mpicom,subname//" bcast sum")
1853  call mct_avect_clean(av1)
1854  if (mytask == 0) then
1855  call mct_avect_clean(av1g)
1856  endif
1857  else
1858  lsum = 0.0_ip_r8_p
1859  do n = 1,lsize
1860  do m = 1,fsize
1861  lsum(m) = lsum(m) + av%rAttr(m,n)*lwts(n)
1862  enddo
1863  enddo
1864  call oasis_mpi_sum(lsum,sum,mpicom,string=trim(subname)//':sum',all=.true.)
1865  endif
1866 
1867  deallocate(lsum)
1868  deallocate(lwts)
1869 
1870  call oasis_debug_exit(subname)
1871 
1872  END SUBROUTINE oasis_advance_avsum
1873 
1874 !-------------------------------------------------------------------
1875 
1876 !> A generic method for writing the global sums of fields in an attribute vector
1877 
1878  SUBROUTINE oasis_advance_avdiag(av,mpicom,mask,wts)
1879 
1880  implicit none
1881  type(mct_avect) ,intent(in) :: av ! input av
1882  integer(kind=ip_i4_p),intent(in) :: mpicom ! mpi communicator
1883  integer(kind=ip_i4_p),intent(in),optional :: mask(:) ! mask to apply to av
1884  real(kind=ip_r8_p) ,intent(in),optional :: wts(:) ! wts to apply to av
1885 
1886  integer(kind=ip_i4_p) :: n,m,ierr,mype
1887  integer(kind=ip_i4_p) :: lsize,fsize ! local size of av, number of flds in av
1888  logical :: first_call
1889  real(kind=ip_r8_p) :: lval ! local temporary
1890  real(kind=ip_r8_p),allocatable :: lsum(:) ! local sum
1891  real(kind=ip_r8_p),allocatable :: lmin(:) ! local min
1892  real(kind=ip_r8_p),allocatable :: lmax(:) ! local max
1893  real(kind=ip_r8_p),allocatable :: gsum(:) ! global sum
1894  real(kind=ip_r8_p),allocatable :: gmin(:) ! global min
1895  real(kind=ip_r8_p),allocatable :: gmax(:) ! global max
1896  real(kind=ip_r8_p),allocatable :: lwts(:) ! local wts taking into account mask and wts
1897  type(mct_string) :: mstring ! mct char type
1898  character(len=64):: itemc ! string converted to char
1899  character(len=*),parameter :: subname = '(oasis_advance_avdiag)'
1900 
1901  call oasis_debug_enter(subname)
1902 
1903  fsize = mct_avect_nrattr(av)
1904  lsize = mct_avect_lsize(av)
1905 
1906  allocate(lsum(fsize))
1907  allocate(lmin(fsize))
1908  allocate(lmax(fsize))
1909  allocate(gsum(fsize))
1910  allocate(gmin(fsize))
1911  allocate(gmax(fsize))
1912 
1913  allocate(lwts(lsize))
1914  lwts = 1.0_ip_r8_p
1915 !!$ lmin=HUGE(lwts)
1916 !!$ lmax=-lmin
1917  if (present(mask)) then
1918  if (size(mask) /= lsize) then
1919  WRITE(nulprt,*) subname,estr,'size mask ne size av'
1920  CALL oasis_abort()
1921  endif
1922  do n = 1,lsize
1923  if (mask(n) /= 0) lwts(n) = 0.0_ip_r8_p
1924  enddo
1925  endif
1926 
1927  if (present(wts)) then
1928  if (size(wts) /= lsize) then
1929  WRITE(nulprt,*) subname,estr,'size wts ne size av'
1930  CALL oasis_abort()
1931  endif
1932  do n = 1,lsize
1933  lwts(n) = lwts(n) * wts(n)
1934  enddo
1935  endif
1936 
1937  lsum = 0.0_ip_r8_p
1938  lmin = 9.99e36 ! in case lsize is zero
1939  lmax = -9.99e36 ! in case lsize is zero
1940  do m = 1,fsize
1941  first_call = .true.
1942  do n = 1,lsize
1943  lval = av%rAttr(m,n)*lwts(n)
1944  lsum(m) = lsum(m) + lval
1945  if (lwts(n) /= 0.0_ip_r8_p) then
1946  if (first_call) then
1947  lmin(m) = lval
1948  lmax(m) = lval
1949  first_call = .false.
1950  else
1951  lmin(m) = min(lmin(m),lval)
1952  lmax(m) = max(lmax(m),lval)
1953  endif
1954  endif
1955  enddo
1956  enddo
1957 
1958  mype = -1
1959  if (mpicom /= mpi_comm_null) then
1960  call mpi_comm_rank(mpicom,mype,ierr)
1961  call oasis_mpi_sum(lsum,gsum,mpicom,string=trim(subname)//':sum',all=.false.)
1962  call oasis_mpi_min(lmin,gmin,mpicom,string=trim(subname)//':min',all=.false.)
1963  call oasis_mpi_max(lmax,gmax,mpicom,string=trim(subname)//':max',all=.false.)
1964  endif
1965  if (mype == 0) then
1966  do m = 1,fsize
1967  call mct_avect_getrlist(mstring,m,av)
1968  itemc = mct_string_tochar(mstring)
1969  call mct_string_clean(mstring)
1970  write(nulprt,'(a,a16,3g21.12)') ' diags: ',trim(itemc),gmin(m),gmax(m),gsum(m)
1971  enddo
1972  endif
1973 
1974  deallocate(lsum,lmin,lmax)
1975  deallocate(gsum,gmin,gmax)
1976  deallocate(lwts)
1977 
1978  call oasis_debug_exit(subname)
1979 
1980  END SUBROUTINE oasis_advance_avdiag
1981 
1982 !-------------------------------------------------------------------
1983 END MODULE mod_oasis_advance
1984 
Generic overloaded interface into MPI sum reduction.
subroutine, public oasis_debug_exit(string)
Used when a subroutine is exited, write info to log file at some debug level.
OASIS partition data and methods.
System type methods.
Generic overloaded interface into MPI max reduction.
subroutine, public oasis_io_write_array(rstfile, mpicom, iarray, ivarname, rarray, rvarname)
Writes a real or integer array to a file.
subroutine oasis_advance_map(av1, avd, mapper, conserv, consbfb, avon, av2, av3, av4, av5)
Provides interpolation functionality.
Provides a generic and simpler interface into MPI calls for OASIS.
Generic overloaded interface into MPI broadcast.
Advances the OASIS coupling.
subroutine, public oasis_timer_start(timer_label, barrier)
Start a timer.
subroutine, public oasis_io_read_avfile(rstfile, av, gsmap, mpicom, abort, nampre, didread)
Reads all fields for an attribute vector from a file.
subroutine, public oasis_timer_stop(timer_label)
Stop a timer.
subroutine, public oasis_advance_run(mop, varid, msec, kinfo, nff, namid, array1din, array1dout, array2dout, readrest, a2on, array2, a3on, array3, a4on, array4, a5on, array5)
Advances the OASIS coupling.
Provides a common location for several OASIS variables.
subroutine, public oasis_io_write_avfbf(av, gsmap, mpicom, nx, ny, msec, f_string, filename)
Write each field in an attribute vector to an individual files.
OASIS map (interpolation) data and methods.
subroutine, public oasis_debug_enter(string)
Used when a subroutine is entered, write info to log file at some debug level.
Defines kinds for OASIS.
subroutine, public oasis_debug_note(string)
Used to write information from a subroutine, write info to log file at some debug level...
subroutine, public oasis_flush(nu)
Flushes output to file.
Initialize the OASIS coupler infrastructure.
Coupler data for managing all aspects of coupling in OASIS.
Generic overloaded interface into MPI min reduction.
Performance timer methods.
subroutine, public oasis_io_write_avfile(rstfile, av, gsmap, mpicom, nx, ny, nampre)
Writes all fields from an attribute vector to a file.
subroutine oasis_advance_avdiag(av, mpicom, mask, wts)
A generic method for writing the global sums of fields in an attribute vector.
subroutine, public oasis_abort(id_compid, cd_routine, cd_message)
OASIS abort method, publically available to users.
subroutine oasis_advance_avsum(av, sum, gsmap, mpicom, mask, wts, consbfb)
A generic method for summing fields in an attribute vector.
subroutine, public oasis_io_read_array(rstfile, mpicom, iarray, ivarname, rarray, rvarname, abort)
Reads an integer or real field from a file into an array.
Provides reusable IO routines for OASIS.
Definition: mod_oasis_io.F90:4
subroutine, public oasis_advance_init(kinfo)
Initializes the OASIS fields.
Mapper data for interpolating data between grids.
OASIS variable data and methods.
Defines parameters for OASIS.
subroutine, public oasis_io_read_avfbf(av, gsmap, mpicom, msec, f_string, filename)
Read each field in an attribute vector from individual files.