42 INTEGER(kind=ip_i4_p),
intent(inout) :: kinfo
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(:)
49 real (kind=ip_r8_p),
allocatable :: array2(:)
50 real (kind=ip_r8_p),
allocatable :: array3(:)
51 real (kind=ip_r8_p),
allocatable :: array4(:)
52 real (kind=ip_r8_p),
allocatable :: array5(:)
53 logical :: a2on,a3on,a4on,a5on
54 integer(kind=ip_i4_p) :: mseclag
55 character(len=ic_xl) :: rstfile
56 character(len=ic_med) :: vstring
58 character(len=*),
parameter :: subname =
'(oasis_advance_init)'
67 if (mpi_comm_local == mpi_comm_null)
then
68 write(nulprt,*) subname,estr,
'called on non coupling task'
76 if (oasis_debug >= 2)
then
77 write(nulprt,*)
' subname at time time+lag act: field '
78 write(nulprt,*)
' diags : fldname min max sum '
88 DO cplid = 1,prism_mcoupler
90 if (npc == 1) pcpointer => prism_coupler_put(cplid)
91 if (npc == 2) pcpointer => prism_coupler_get(cplid)
92 if (pcpointer%valid)
then
95 ltime = pcpointer%ltime
96 getput= pcpointer%getput
97 rstfile=trim(pcpointer%rstfile)
98 partid= pcpointer%partID
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,*)
'----------------------------------------------------------------'
115 IF (lag > dt .OR. lag <= -dt)
THEN
116 WRITE(nulprt,*) subname,estr,
'lag out of dt range cplid/dt/lag=',cplid,dt,lag
127 IF ( (getput == oasis3_put .AND. lag > 0) )
THEN
129 lsize = mct_avect_lsize(pcpointer%aVect1)
130 nflds = mct_avect_nrattr(pcpointer%aVect1)
132 ALLOCATE(array(lsize))
133 ALLOCATE(array2(lsize))
134 ALLOCATE(array3(lsize))
135 ALLOCATE(array4(lsize))
136 ALLOCATE(array5(lsize))
139 varid = pcpointer%varid(nf)
141 namid=cplid,array1din=array,&
142 readrest=.true., array2=array2,array3=array3, &
143 array4=array4,array5=array5)
145 IF (oasis_debug >= 2)
THEN
146 write(nulprt,*) subname,
' advance_run ',cplid,trim(pcpointer%fldlist)
165 DO cplid=1, prism_mcoupler
167 if (npc == 1) pcpointer => prism_coupler_put(cplid)
168 if (npc == 2) pcpointer => prism_coupler_get(cplid)
169 IF (pcpointer%valid)
then
172 ltime = pcpointer%ltime
173 getput= pcpointer%getput
174 rstfile=trim(pcpointer%rstfile)
175 partid= pcpointer%partID
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,*)
'----------------------------------------------------------------'
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'
199 if (oasis_debug >= 2)
then
200 write(nulprt,*) subname,
' at ',msec,mseclag,
' RTRN: ',&
201 trim(pcpointer%fldlist),
' ',trim(rstfile)
203 lsize = mct_avect_lsize(pcpointer%aVect1)
205 write(vstring,
'(a,i6.6,a)')
'loc',pcpointer%namID,
'_cnt'
207 ivarname=trim(vstring),abort=.false.)
209 write(vstring,
'(a,i6.6,a)')
'loc',pcpointer%namID,
'_'
211 prism_part(partid)%pgsmap,prism_part(partid)%mpicom, &
212 abort=.false.,nampre=trim(vstring))
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,
'_'
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)
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,
'_'
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)
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,
'_'
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)
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,
'_'
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)
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)
283 array1din,array1dout,array2dout,readrest,&
284 a2on,array2,a3on,array3,a4on,array4,a5on,array5)
288 integer(kind=ip_i4_p),
intent(in) :: mop
289 INTEGER(kind=ip_i4_p),
intent(in) :: varid
290 INTEGER(kind=ip_i4_p),
intent(in) :: msec
291 INTEGER(kind=ip_i4_p),
intent(inout) :: kinfo
293 INTEGER(kind=ip_i4_p),
optional :: nff
294 INTEGER(kind=ip_i4_p),
optional :: namid
296 REAL (kind=ip_r8_p),
optional :: array1din(:)
297 REAL (kind=ip_r8_p),
optional :: array1dout(:)
298 REAL (kind=ip_r8_p),
optional :: array2dout(:,:)
299 logical ,
optional :: readrest
301 logical ,
optional :: a2on
302 REAL (kind=ip_r8_p),
optional :: array2(:)
303 logical ,
optional :: a3on
304 REAL (kind=ip_r8_p),
optional :: array3(:)
305 logical ,
optional :: a4on
306 REAL (kind=ip_r8_p),
optional :: array4(:)
307 logical ,
optional :: a5on
308 REAL (kind=ip_r8_p),
optional :: array5(:)
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
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
322 character(len=ic_xl) :: inpfile
323 integer(kind=ip_i4_p) :: nx,ny
324 integer(kind=ip_i4_p) :: mseclag
325 real(kind=ip_r8_p) :: rcnt
326 character(len=ic_med) :: tstring
327 character(len=ic_med) :: fstring
328 character(len=ic_med) :: cstring
329 character(len=ic_med) :: vstring
334 TYPE(mct_avect) :: avtest
335 type(mct_avect) :: avtmp
336 type(mct_avect) :: avtmp2
337 type(mct_avect) :: avtmp3
338 type(mct_avect) :: avtmp4
339 type(mct_avect) :: avtmp5
342 logical,
parameter :: allow_no_restart = .false.
349 character(len=*),
parameter :: subname =
'(oasis_advance_run)'
350 character(len=*),
parameter :: F01 =
'(a,i3.3)'
360 if (mpi_comm_local == mpi_comm_null)
then
361 write(nulprt,*) subname,estr,
'called on non coupling task'
366 if (varid < 1 .or. varid > prism_nvar)
then
367 write(nulprt,*) subname,estr,
'invalid varid',varid
370 vname = prism_var(varid)%name
373 if (
present(readrest))
then
376 if (lreadrest) kinfo = oasis_fromrest
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
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)
400 do nc = 1,prism_var(varid)%ncpl
401 cplid = prism_var(varid)%cpl(nc)
402 if (cplid == namid) runit = .true.
405 WRITE(nulprt,*) subname,estr,
'namid not found for var = ',trim(vname)
406 WRITE(nulprt,*) subname,estr,
'namid = ',namid
416 DO nc = 1,prism_var(varid)%ncpl
417 cplid = prism_var(varid)%cpl(nc)
419 if (
present(namid))
then
421 if (cplid == namid) runit = .true.
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)
431 if (.not.pcpointer%valid)
then
432 WRITE(nulprt,*) subname,estr,
'invalid prism_coupler for var = ',trim(vname)
439 getput = pcpointer%getput
440 if ((mop == oasis_out .and. getput == oasis3_put) .or. &
441 (mop == oasis_in .and. getput == oasis3_get))
then
444 write(nulprt,*) subname,estr,
'model def_var in-out does not match model get-put call for var = ',trim(vname)
451 rouid = pcpointer%routerid
452 mapid = pcpointer%mapperid
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
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
474 unpack = (sndrcv .OR. input)
477 IF (prism_part(partid)%nx >= 1)
THEN
478 nx = prism_part(partid)%nx
479 ny = prism_part(partid)%ny
481 nx = prism_part(partid)%gsize
485 IF (oasis_debug >= 20)
THEN
486 WRITE(nulprt,*) subname,
' DEBUG nx, ny = ',nx,ny
494 IF (abs(lag) > dt)
THEN
495 WRITE(nulprt,*) subname,estr,
'lag setting greater than dt for var = ',trim(vname)
506 IF (getput == oasis3_put .AND. lag > 0 .AND. lreadrest)
THEN
509 IF (.not.
present(nff))
THEN
510 WRITE(nulprt,*) subname,estr,
'nff optional argument not passed but expected for var = ',trim(vname)
513 IF (len_trim(rstfile) < 1)
THEN
514 WRITE(nulprt,*) subname,estr,
'restart file undefined for var = ',trim(vname)
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)
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)
528 CALL
oasis_io_read_avfile(trim(rstfile),avtmp,prism_part(partid)%pgsmap,prism_part(partid)%mpicom)
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)
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)
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)
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)
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)
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)
567 if (getput == oasis3_put)
then
569 elseif (getput == oasis3_get)
then
573 if (oasis_debug >= 20)
then
574 write(nulprt,*) subname,
' DEBUG msec,mseclag = ',msec,mseclag
579 if (mod(mseclag,dt) == 0) time_now = .true.
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
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
608 do n = 1,prism_mcoupler
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'
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'
648 lastseq = pcpointer%seq
651 if (oasis_debug >= 20)
then
652 write(nulprt,*) subname,
' DEBUG sequence ',trim(vname),msec,lastseqtime,lastseq,pcpointer%seq
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)
668 if (oasis_debug >= 20)
then
669 write(nulprt,*) subname,
' DEBUG nfav,nsav,nsa = ',nfav,nsav,nsa
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
685 IF (.not.lreadrest)
THEN
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
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'
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'
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'
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'
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'
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'
764 if (getput == oasis3_put)
then
767 write(tstring,f01)
'pcpy_',cplid
771 if (lreadrest .or. pcpointer%trans == ip_instant)
then
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)
780 pcpointer%avect2%rAttr(nfav,n) = 0.0
783 if (pcpointer%aVon(3))
then
784 if (
present(array3))
then
785 pcpointer%avect3%rAttr(nfav,n) = array3(n)
787 pcpointer%avect3%rAttr(nfav,n) = 0.0
790 if (pcpointer%aVon(4))
then
791 if (
present(array4))
then
792 pcpointer%avect4%rAttr(nfav,n) = array4(n)
794 pcpointer%avect4%rAttr(nfav,n) = 0.0
797 if (pcpointer%aVon(5))
then
798 if (
present(array5))
then
799 pcpointer%avect5%rAttr(nfav,n) = array5(n)
801 pcpointer%avect5%rAttr(nfav,n) = 0.0
805 pcpointer%avcnt(nfav) = 1
808 elseif (pcpointer%trans == ip_average)
then
810 if (kinfo == oasis_ok) kinfo = oasis_loctrans
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)
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)
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)
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)
839 pcpointer%avcnt(nfav) = pcpointer%avcnt(nfav) + 1
841 elseif (pcpointer%trans == ip_accumul)
then
843 if (kinfo == oasis_ok) kinfo = oasis_loctrans
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)
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)
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)
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)
872 pcpointer%avcnt(nfav) = 1
874 elseif (pcpointer%trans == ip_max)
then
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'
884 if (pcpointer%avcnt(nfav) == 0)
then
885 pcpointer%avect1%rAttr(nfav,n) = array1din(n)
887 pcpointer%avect1%rAttr(nfav,n) = &
888 max(pcpointer%avect1%rAttr(nfav,n),array1din(n))
891 pcpointer%avcnt(nfav) = 1
893 elseif (pcpointer%trans == ip_min)
then
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'
903 if (pcpointer%avcnt(nfav) == 0)
then
904 pcpointer%avect1%rAttr(nfav,n) = array1din(n)
906 pcpointer%avect1%rAttr(nfav,n) = &
907 min(pcpointer%avect1%rAttr(nfav,n),array1din(n))
910 pcpointer%avcnt(nfav) = 1
913 write(nulprt,*) subname,estr,
'transform not known for var = ',trim(vname),pcpointer%trans
918 if (oasis_debug >= 2 .and. trim(cstring) /=
'none')
then
919 write(nulprt,*) subname,
' at ',msec,mseclag,
' PACK: ',&
920 trim(vname),
' ',trim(cstring)
923 if (oasis_debug >= 20)
then
924 write(nulprt,*) subname,
' DEBUG loctrans update ',cplid,
' ',&
925 trim(cstring),pcpointer%avcnt(nfav)
929 pcpointer%status(nfav) = oasis_comm_ready
945 do nf = 1,pcpointer%nflds
946 if (pcpointer%status(nf) /= oasis_comm_ready)
then
948 if (oasis_debug >= 15)
then
949 write(nulprt,*) subname,
' at ',msec,mseclag,
' STAT: ',nf,
' NOT READY'
951 kinfo=oasis_waitgroup
953 if (oasis_debug >= 15)
then
954 write(nulprt,*) subname,
' at ',msec,mseclag,
' STAT: ',nf,
' READY'
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
986 IF (getput == oasis3_put)
THEN
988 write(tstring,f01)
'pavg_',cplid
990 do nf = 1,pcpointer%nflds
991 if (pcpointer%avcnt(nf) > 1)
then
992 rcnt = 1.0/pcpointer%avcnt(nf)
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
1000 if (pcpointer%aVon(3))
then
1001 pcpointer%avect3%rAttr(nf,n) = &
1002 pcpointer%avect3%rAttr(nf,n) * rcnt
1004 if (pcpointer%aVon(4))
then
1005 pcpointer%avect4%rAttr(nf,n) = &
1006 pcpointer%avect4%rAttr(nf,n) * rcnt
1008 if (pcpointer%aVon(5))
then
1009 pcpointer%avect5%rAttr(nf,n) = &
1010 pcpointer%avect5%rAttr(nf,n) * rcnt
1014 if (oasis_debug >= 20)
then
1015 write(nulprt,*) subname,
' DEBUG loctrans calc0 = ',cplid,nf,&
1017 write(nulprt,*) subname,
' DEBUG loctrans calc1 = ',cplid,nf,&
1018 minval(pcpointer%avect1%rAttr(nf,:)),&
1019 maxval(pcpointer%avect1%rAttr(nf,:))
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,:))
1050 if (mseclag >= maxtime)
then
1053 if (getput == oasis3_put .and. lag > 0 .and. mseclag == maxtime)
then
1054 kinfo = oasis_torest
1056 write(tstring,f01)
'wrst_',cplid
1059 prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny)
1060 if (pcpointer%aVon(2)) &
1062 prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=
'av2_')
1063 if (pcpointer%aVon(3)) &
1065 prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=
'av3_')
1066 if (pcpointer%aVon(4)) &
1068 prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=
'av4_')
1069 if (pcpointer%aVon(5)) &
1071 prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=
'av5_')
1073 if (oasis_debug >= 2)
then
1074 write(nulprt,*) subname,
' at ',msec,mseclag,
' WRST: ', &
1075 trim(mct_avect_exportrlist2c(pcpointer%avect1)),
' ',&
1087 if (getput == oasis3_put)
then
1090 if (oasis_debug >= 2)
then
1091 write(nulprt,*) subname,
' at ',msec,mseclag,
' SEND: ', &
1092 trim(mct_avect_exportrlist2c(pcpointer%avect1))
1095 if (sndadd /= 0.0_ip_double_p .or. sndmult /= 1.0_ip_double_p)
then
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)
1103 pcpointer%avect1%rAttr(:,:) = pcpointer%avect1%rAttr(:,:)*sndmult &
1108 write(tstring,f01)
'pmap_',cplid
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)
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)
1137 pcpointer%avect1m,prism_mapper(mapid),conserv,consbfb, &
1138 pcpointer%aVon ,pcpointer%avect2, &
1139 pcpointer%avect3,pcpointer%avect4, &
1141 if (lucia_debug > 0) &
1142 WRITE(nullucia, fmt=
'(A,I3.3,A,F16.5)') &
1143 'Balance: ',pcpointer%namID,
' After interpo ', mpi_wtime()
1145 write(tstring,f01)
'psnd_',cplid
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)
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()
1163 write(tstring,f01)
'psnd_',cplid
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)
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()
1181 elseif (getput == oasis3_get)
then
1183 if (oasis_debug >= 2 )
then
1184 write(nulprt,*) subname,
' at ',msec,mseclag,
' RECV: ', &
1185 trim(mct_avect_exportrlist2c(pcpointer%avect1))
1190 write(tstring,f01)
'grcv_',cplid
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()
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)
1207 write(tstring,f01)
'gmap_',cplid
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)
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()
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)
1225 write(tstring,f01)
'grcv_',cplid
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()
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)
1243 if (rcvadd /= 0.0_ip_double_p .or. rcvmult /= 1.0_ip_double_p)
then
1244 pcpointer%avect1%rAttr(:,:) = pcpointer%avect1%rAttr(:,:)*rcvmult &
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)
1262 if (kinfo == oasis_sent)
then
1263 kinfo = oasis_sentout
1264 elseif (kinfo == oasis_torest)
then
1265 kinfo = oasis_torestout
1267 kinfo = oasis_output
1270 write(tstring,f01)
'wout_',cplid
1272 if (oasis_debug >= 2)
then
1273 write(nulprt,*) subname,
' at ',msec,mseclag,
' WRIT: ', &
1274 trim(mct_avect_exportrlist2c(pcpointer%avect1))
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, &
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
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)
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 = '
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 = '
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))
1349 IF (mseclag + dt >= maxtime .AND. &
1350 getput == oasis3_put .and. pcpointer%trans /= ip_instant)
then
1352 write(tstring,f01)
'wtrn_',cplid
1354 WRITE(vstring,
'(a,i6.6,a)')
'loc',pcpointer%namID,
'_cnt'
1356 ivarname=trim(vstring))
1357 write(vstring,
'(a,i6.6,a)')
'loc',pcpointer%namID,
'_'
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,
'_'
1363 prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=trim(vstring))
1365 if (pcpointer%aVon(3))
then
1366 write(vstring,
'(a,i6.6,a)')
'av3loc',pcpointer%namID,
'_'
1368 prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=trim(vstring))
1370 if (pcpointer%aVon(4))
then
1371 write(vstring,
'(a,i6.6,a)')
'av4loc',pcpointer%namID,
'_'
1373 prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=trim(vstring))
1375 if (pcpointer%aVon(5))
then
1376 write(vstring,
'(a,i6.6,a)')
'av5loc',pcpointer%namID,
'_'
1378 prism_part(partid)%pgsmap,prism_part(partid)%mpicom,nx,ny,nampre=trim(vstring))
1381 if (oasis_debug >= 2)
then
1382 write(nulprt,*) subname,
' at ',msec,mseclag,
' WTRN: ', &
1383 trim(mct_avect_exportrlist2c(pcpointer%avect1)),
' ',trim(rstfile)
1386 if (oasis_debug >= 20)
then
1387 write(nulprt,*) subname,
' DEBUG write loctrans restart',cplid,&
1389 write(nulprt,*) subname,
' DEBUG write loctrans restart',cplid,&
1390 minval(pcpointer%avect1%rAttr),&
1391 maxval(pcpointer%avect1%rAttr)
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
1411 if (oasis_debug >= 2)
then
1412 write(nulprt,*) subname,
' at ',msec,mseclag,
' READ: ', &
1413 trim(mct_avect_exportrlist2c(pcpointer%avect1))
1416 write(tstring,f01)
'grin_',cplid
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))
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)
1428 if (oasis_debug >= 2)
then
1429 write(nulprt,*) subname,
' at ',msec,mseclag,
' UPCK: ',trim(vname)
1431 write(tstring,f01)
'gcpy_',cplid
1434 if (
present(array1dout)) array1dout(:) = &
1435 pcpointer%avect1%rAttr(nfav,:)
1436 if (
present(array2dout)) array2dout(:,:) = &
1437 reshape(pcpointer%avect1%rAttr(nfav,:),shape(array2dout))
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)
1446 if (time_now) pcpointer%status(nfav) = oasis_comm_ready
1456 if (oasis_debug >= 2)
then
1457 write(nulprt,*) subname,
' at ',msec,mseclag,
' KINF: ',trim(vname),kinfo
1476 avon,av2,av3,av4,av5)
1481 type(mct_avect) ,
intent(in) :: av1
1482 type(mct_avect) ,
intent(inout) :: avd
1484 integer(kind=ip_i4_p) ,
intent(in),
optional :: conserv
1485 logical ,
intent(in),
optional :: consbfb
1486 logical ,
intent(in),
optional :: avon(:)
1487 type(mct_avect) ,
intent(in),
optional :: av2
1488 type(mct_avect) ,
intent(in),
optional :: av3
1489 type(mct_avect) ,
intent(in),
optional :: av4
1490 type(mct_avect) ,
intent(in),
optional :: av5
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(:)
1497 character(len=ic_med) :: tstring
1498 type(mct_avect) :: avdtmp
1499 type(mct_avect) :: av2g
1501 integer(kind=ip_i4_p),
parameter :: avsmax = prism_coupler_avsmax
1502 logical :: locavon(avsmax)
1503 integer(kind=ip_i4_p) :: avonsize, nterm
1504 character(len=*),
parameter :: subname =
'(oasis_advance_map)'
1512 if (
present(consbfb))
then
1524 if (
present(avon))
then
1525 avonsize =
size(avon)
1526 nterm = min(avsmax,avonsize)
1527 locavon(1:nterm) = avon(1:nterm)
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'
1538 if (locavon(n) .and. n > mapper%nwgts)
then
1539 WRITE(nulprt,*) subname,estr,
'in nwgts and coupling terms',mapper%nwgts,n
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'
1551 call mct_smat_avmult(av1, mapper%sMatP(1), avd)
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)
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'
1563 call mct_smat_avmult(av2, mapper%sMatP(2), avdtmp)
1564 avd%rAttr = avd%rAttr + avdtmp%rAttr
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'
1572 call mct_smat_avmult(av3, mapper%sMatP(3), avdtmp)
1573 avd%rAttr = avd%rAttr + avdtmp%rAttr
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'
1581 call mct_smat_avmult(av4, mapper%sMatP(4), avdtmp)
1582 avd%rAttr = avd%rAttr + avdtmp%rAttr
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'
1590 call mct_smat_avmult(av5, mapper%sMatP(5), avdtmp)
1591 avd%rAttr = avd%rAttr + avdtmp%rAttr
1594 call mct_avect_clean(avdtmp)
1602 IF (
PRESENT(conserv))
THEN
1604 IF (conserv /= ip_cnone)
THEN
1605 fsize = mct_avect_nrattr(av1)
1606 allocate(av_sums(fsize),av_sumd(fsize))
1608 zradi = 1./(eradius*eradius)
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
1621 call mct_avect_gather(mapper%av_ms,av2g,prism_part(mapper%spart)%gsmap,&
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
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)
1636 sumtmp = 0.0_ip_r8_p
1638 if (imasks(n) == 0) sumtmp = sumtmp + areas(n)
1640 call
oasis_mpi_sum(sumtmp,wts_sums,mpi_comm_local,string=subname//
':wts_sums',&
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
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
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)
1669 sumtmp = 0.0_ip_r8_p
1671 if (imaskd(n) == 0) sumtmp = sumtmp + aread(n)
1673 call
oasis_mpi_sum(sumtmp,wts_sumd,mpi_comm_local,string=subname//
':wts_sumd',all=.true.)
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
1693 mask=imasks,wts=areas,consbfb=lconsbfb)
1695 mask=imaskd,wts=aread,consbfb=lconsbfb)
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
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 '
1708 zlagr = (av_sumd(m) - av_sums(m)) / wts_sumd
1710 if (imaskd(n) == 0) avd%rAttr(m,n) = avd%rAttr(m,n) - zlagr
1713 elseif (conserv == ip_cglbpos)
then
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'
1718 elseif (av_sumd(m) /= 0.0_ip_r8_p)
then
1719 zlagr = av_sums(m) / av_sumd(m)
1721 if (imaskd(n) == 0) avd%rAttr(m,n) = avd%rAttr(m,n) * zlagr
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'
1731 zlagr = (av_sumd(m) - (av_sums(m)*(wts_sumd/wts_sums))) / wts_sumd
1733 if (imaskd(n) == 0) avd%rAttr(m,n) = avd%rAttr(m,n) - zlagr
1736 elseif (conserv == ip_cbaspos)
then
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'
1741 elseif (av_sumd(m) /= 0.0_ip_r8_p)
then
1742 zlagr = (av_sums(m)/av_sumd(m)) * (wts_sumd/wts_sums)
1744 if (imaskd(n) == 0) avd%rAttr(m,n) = avd%rAttr(m,n) * zlagr
1749 WRITE(nulprt,*) subname,estr,
'conserv option unknown = ',conserv
1753 if (oasis_debug >= 20)
then
1755 mask=imasks,wts=areas,consbfb=lconsbfb)
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
1763 deallocate(imasks,imaskd,areas,aread)
1764 deallocate(av_sums,av_sumd)
1779 type(mct_avect) ,
intent(in) :: av
1780 real(kind=ip_r8_p) ,
intent(inout) :: sum(:)
1781 type(mct_gsmap) ,
intent(in) :: gsmap
1782 integer(kind=ip_i4_p),
intent(in) :: mpicom
1783 integer(kind=ip_i4_p),
intent(in),
optional :: mask(:)
1784 real(kind=ip_r8_p) ,
intent(in),
optional :: wts(:)
1785 logical ,
intent(in),
optional :: consbfb
1787 integer(kind=ip_i4_p) :: n,m,ierr,mytask
1788 integer(kind=ip_i4_p) :: lsize,fsize
1789 real(kind=ip_r8_p),
allocatable :: lsum(:)
1790 real(kind=ip_r8_p),
allocatable :: lwts(:)
1791 type(mct_avect) :: av1, av1g
1793 character(len=*),
parameter :: subname =
'(oasis_advance_avsum)'
1798 if (
present(consbfb))
then
1802 fsize = mct_avect_nrattr(av)
1803 lsize = mct_avect_lsize(av)
1805 allocate(lsum(fsize))
1807 allocate(lwts(lsize))
1810 if (
size(sum) /= fsize)
then
1811 WRITE(nulprt,*) subname,estr,
'size sum ne size av'
1815 if (
present(mask))
then
1816 if (
size(mask) /= lsize)
then
1817 WRITE(nulprt,*) subname,estr,
'size mask ne size av'
1821 if (mask(n) /= 0) lwts(n) = 0.0_ip_r8_p
1825 if (
present(wts))
then
1826 if (
size(wts) /= lsize)
then
1827 WRITE(nulprt,*) subname,estr,
'size wts ne size av'
1831 lwts(n) = lwts(n) * wts(n)
1836 call mct_avect_init(av1,av,lsize)
1839 av1%rAttr(m,n) = av%rAttr(m,n)*lwts(n)
1842 call mct_avect_gather(av1,av1g,gsmap,0,mpicom)
1843 call mpi_comm_rank(mpicom,mytask,ierr)
1845 if (mytask == 0)
then
1846 do n = 1,mct_avect_lsize(av1g)
1848 sum(m) = sum(m) + av1g%rAttr(m,n)
1853 call mct_avect_clean(av1)
1854 if (mytask == 0)
then
1855 call mct_avect_clean(av1g)
1861 lsum(m) = lsum(m) + av%rAttr(m,n)*lwts(n)
1864 call
oasis_mpi_sum(lsum,sum,mpicom,string=trim(subname)//
':sum',all=.true.)
1881 type(mct_avect) ,
intent(in) :: av
1882 integer(kind=ip_i4_p),
intent(in) :: mpicom
1883 integer(kind=ip_i4_p),
intent(in),
optional :: mask(:)
1884 real(kind=ip_r8_p) ,
intent(in),
optional :: wts(:)
1886 integer(kind=ip_i4_p) :: n,m,ierr,mype
1887 integer(kind=ip_i4_p) :: lsize,fsize
1888 logical :: first_call
1889 real(kind=ip_r8_p) :: lval
1890 real(kind=ip_r8_p),
allocatable :: lsum(:)
1891 real(kind=ip_r8_p),
allocatable :: lmin(:)
1892 real(kind=ip_r8_p),
allocatable :: lmax(:)
1893 real(kind=ip_r8_p),
allocatable :: gsum(:)
1894 real(kind=ip_r8_p),
allocatable :: gmin(:)
1895 real(kind=ip_r8_p),
allocatable :: gmax(:)
1896 real(kind=ip_r8_p),
allocatable :: lwts(:)
1897 type(mct_string) :: mstring
1898 character(len=64):: itemc
1899 character(len=*),
parameter :: subname =
'(oasis_advance_avdiag)'
1903 fsize = mct_avect_nrattr(av)
1904 lsize = mct_avect_lsize(av)
1906 allocate(lsum(fsize))
1907 allocate(lmin(fsize))
1908 allocate(lmax(fsize))
1909 allocate(gsum(fsize))
1910 allocate(gmin(fsize))
1911 allocate(gmax(fsize))
1913 allocate(lwts(lsize))
1917 if (
present(mask))
then
1918 if (
size(mask) /= lsize)
then
1919 WRITE(nulprt,*) subname,estr,
'size mask ne size av'
1923 if (mask(n) /= 0) lwts(n) = 0.0_ip_r8_p
1927 if (
present(wts))
then
1928 if (
size(wts) /= lsize)
then
1929 WRITE(nulprt,*) subname,estr,
'size wts ne size av'
1933 lwts(n) = lwts(n) * wts(n)
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
1949 first_call = .false.
1951 lmin(m) = min(lmin(m),lval)
1952 lmax(m) = max(lmax(m),lval)
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.)
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)
1974 deallocate(lsum,lmin,lmax)
1975 deallocate(gsum,gmin,gmax)
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.
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.
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.
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.