88
89
98
99# ifdef FLOAT_BIOLOGY
101# endif
102# ifdef DISTRIBUTE
104# endif
106# if defined SOLVE3D && defined FLOAT_VWALK
108# endif
109
110
111
112 integer, intent(in) :: ng, Lstr, Lend
113 integer, intent(in) :: knew, nnew, nfm3, nfm2, nfm1, nf, nfp1
114
115# ifdef ASSUMED_SHAPE
116 integer, intent(in) :: Ftype(:)
117 real(r8), intent(in) :: Tinfo(0:,:)
118 real(r8), intent(in) :: Fz0(:)
119
120 logical, intent(inout) :: bounded(:)
121# if defined SOLVE3D && defined FLOAT_STICKY
122 logical, intent(inout) :: stuck(:)
123# endif
124 real(r8), intent(inout) :: track(:,0:,:)
125# else
126 integer, intent(in) :: Ftype(Nfloats(ng))
127 real(r8), intent(in) :: Tinfo(0:izrhs,Nfloats(ng))
128 real(r8), intent(in) :: Fz0(Nfloats(ng))
129
130 logical, intent(inout) :: bounded(Nfloats(ng))
131# if defined SOLVE3D && defined FLOAT_STICKY
132 logical, intent(inout) :: stuck(Nfloats(ng))
133# endif
134 real(r8), intent(inout) :: track(NFV(ng),0:NFT,Nfloats(ng))
135# endif
136
137
138
139 logical, parameter :: Gmask = .false.
140# ifdef MASKING
141 logical, parameter :: Lmask = .true.
142# else
143 logical, parameter :: Lmask = .false.
144# endif
145 logical, dimension(Lstr:Lend) :: my_thread
146
147 integer :: LBi, UBi, LBj, UBj
148 integer :: Ir, Jr, Npts, i, i1, i2, j, j1, j2, itrc, l, k
149
150 real(r8), parameter :: Fspv = 0.0_r8
151
152 real(r8) :: cff1, cff2, cff3, cff4, cff5, cff6, cff7, cff8, cff9
153 real(r8) :: oHz, p1, p2, q1, q2, xrhs, yrhs, zrhs, zfloat
154 real(r8) :: HalfDT
155
156 real(r8), dimension(Lstr:Lend) :: nudg
157
158# ifdef DISTRIBUTE
159 real(r8) :: Xstr, Xend, Ystr, Yend
160 real(r8), dimension(Nfloats(ng)*NFV(ng)*(NFT+1)) :: Fwrk
161# endif
162
163
164
165 lbi=lbound(
grid(ng)%h,dim=1)
166 ubi=ubound(
grid(ng)%h,dim=1)
167 lbj=lbound(
grid(ng)%h,dim=2)
168 ubj=ubound(
grid(ng)%h,dim=2)
169
170# ifdef DISTRIBUTE
171
172
173
174
175
176
177
178
179
180
181
182
183
184
186
191 DO l=lstr,lend
192 my_thread(l)=.false.
193 IF ((xstr.le.track(
ixgrd,nf,l)).and. &
194 & (track(
ixgrd,nf,l).lt.xend).and. &
195 & (ystr.le.track(
iygrd,nf,l)).and. &
196 & (track(
iygrd,nf,l).lt.yend))
THEN
197 my_thread(l)=.true.
198 ELSE IF (
master.and.(.not.bounded(l)))
THEN
199 my_thread(l)=.true.
200 ELSE
203 track(i,j,l)=fspv
204 END DO
205 END DO
206 END IF
207 END DO
208# else
209
210
211
212
213
214 DO l=lstr,lend
215 my_thread(l)=.true.
216 END DO
217# endif
218# if !(defined SOLVE3D && defined FLOAT_VWALK)
219 DO l=lstr,lend
220 nudg(l)=0.0_r8
221 END DO
222# endif
223# if defined SOLVE3D && defined FLOAT_VWALK
224
225
226
227
228
229
230 CALL vwalk_floats (ng, lstr, lend, .true., my_thread, nudg)
231# endif
232
233
234
235
236
237
238 cff1=8.0_r8/3.0_r8
239 cff2=4.0_r8/3.0_r8
240 DO l=lstr,lend
241 IF (my_thread(l).and.bounded(l)) THEN
243 &
dt(ng)*(cff1*track(
ixrhs,nf ,l)- &
244 & cff2*track(
ixrhs,nfm1,l)+ &
245 & cff1*track(
ixrhs,nfm2,l))
247 &
dt(ng)*(cff1*track(
iyrhs,nf ,l)- &
248 & cff2*track(
iyrhs,nfm1,l)+ &
249 & cff1*track(
iyrhs,nfm2,l))
250
251# if defined SOLVE3D && !defined FLOAT_VWALK
252
253
254
256# ifdef FLOAT_BIOLOGY
258 &
dt(ng)*(cff1*track(
izrhs,nf ,l)- &
259 & cff2*track(
izrhs,nfm1,l)+ &
260 & cff1*track(
izrhs,nfm2,l)+ &
261 & cff1*track(
iwbio,nf ,l)* &
262 & track(
i1ohz,nf ,l)- &
263 & cff2*track(
iwbio,nfm1,l)* &
264 & track(
i1ohz,nfm1,l)+ &
265 & cff1*track(
iwbio,nfm2,l)* &
266 & track(
i1ohz,nfm2,l))
267# else
269 &
dt(ng)*(cff1*track(
izrhs,nf ,l)- &
270 & cff2*track(
izrhs,nfm1,l)+ &
271 & cff1*track(
izrhs,nfm2,l))
272# endif
273
274
275
276
277
280 ir=int(track(
ixgrd,nfp1,l))
281 jr=int(track(
iygrd,nfp1,l))
282
283 i1=min(max(ir ,0),
lm(ng)+1)
284 i2=min(max(ir+1,1),
lm(ng)+1)
285 j1=min(max(jr ,0),
mm(ng)+1)
286 j2=min(max(jr+1,0),
mm(ng)+1)
287
288 p2=real(i2-i1,r8)*(track(
ixgrd,nfp1,l)-real(i1,r8))
289 q2=real(j2-j1,r8)*(track(
iygrd,nfp1,l)-real(j1,r8))
290 p1=1.0_r8-p2
291 q1=1.0_r8-q2
292# ifdef MASKING
293 cff7=p1*q1*
grid(ng)%z_w(i1,j1,
n(ng))*
grid(ng)%rmask(i1,j1)+ &
294 & p2*q1*
grid(ng)%z_w(i2,j1,
n(ng))*
grid(ng)%rmask(i2,j1)+ &
295 & p1*q2*
grid(ng)%z_w(i1,j2,
n(ng))*
grid(ng)%rmask(i1,j2)+ &
296 & p2*q2*
grid(ng)%z_w(i2,j2,
n(ng))*
grid(ng)%rmask(i2,j2)
297 cff8=p1*q1*
grid(ng)%rmask(i1,j1)+ &
298 & p2*q1*
grid(ng)%rmask(i2,j1)+ &
299 & p1*q2*
grid(ng)%rmask(i1,j2)+ &
300 & p2*q2*
grid(ng)%rmask(i2,j2)
301 cff9=0.0_r8
302 IF (cff8.gt.0.0_r8) cff9=cff7/cff8
303# else
304 cff9=p1*q1*
grid(ng)%z_w(i1,j1,
n(ng))+ &
305 & p2*q1*
grid(ng)%z_w(i2,j1,
n(ng))+ &
306 & p1*q2*
grid(ng)%z_w(i1,j2,
n(ng))+ &
307 & p2*q2*
grid(ng)%z_w(i2,j2,
n(ng))
308# endif
309 cff6=cff9
310
312 zfloat=fz0(l)
314 zfloat=fz0(l)+cff9
315 END IF
316
318# ifdef MASKING
319 cff7=p1*q1*
grid(ng)%z_w(i1,j1,k)*
grid(ng)%rmask(i1,j1)+ &
320 & p2*q1*
grid(ng)%z_w(i2,j1,k)*
grid(ng)%rmask(i2,j1)+ &
321 & p1*q2*
grid(ng)%z_w(i1,j2,k)*
grid(ng)%rmask(i1,j2)+ &
322 & p2*q2*
grid(ng)%z_w(i2,j2,k)*
grid(ng)%rmask(i2,j2)
323 cff8=p1*q1*
grid(ng)%rmask(i1,j1)+ &
324 & p2*q1*
grid(ng)%rmask(i2,j1)+ &
325 & p1*q2*
grid(ng)%rmask(i1,j2)+ &
326 & p2*q2*
grid(ng)%rmask(i2,j2)
327 IF (cff8.gt.0.0_r8) THEN
328 cff5=cff7/cff8
329 ELSE
330 cff5=0.0_r8
331 END IF
332# else
333 cff5=p1*q1*
grid(ng)%z_w(i1,j1,k)+ &
334 & p2*q1*
grid(ng)%z_w(i2,j1,k)+ &
335 & p1*q2*
grid(ng)%z_w(i1,j2,k)+ &
336 & p2*q2*
grid(ng)%z_w(i2,j2,k)
337# endif
338 IF ((zfloat-cff5)*(cff6-zfloat).ge.0.0_r8) THEN
339 track(
izgrd,nfp1,l)=real(k,r8)+(zfloat-cff5)/(cff6-cff5)
340 END IF
341 cff6=cff5
342 END DO
343 END IF
344# endif
345 END IF
346 END DO
347
348
349
350
351
352# ifdef SOLVE3D
359# ifdef MASKING
360 &
grid(ng) % rmask, &
361# endif
362 &
ocean(ng) % u(:,:,:,nnew), &
363 & my_thread, bounded, track)
364
371# ifdef MASKING
372 &
grid(ng) % rmask, &
373# endif
374 &
ocean(ng) % v(:,:,:,nnew), &
375 & my_thread, bounded, track)
376
377# if !defined FLOAT_VWALK
384# ifdef MASKING
385 &
grid(ng) % rmask, &
386# endif
388 & my_thread, bounded, track)
389# endif
390# if defined FLOAT_BIOLOGY
393 &
r3dvar, gmask, fspv, nudg, &
397# ifdef MASKING
398 &
grid(ng) % rmask, &
399# endif
401 & my_thread, bounded, track)
402 DO l=lstr,lend
403 IF (my_thread(l).and.bounded(l)) THEN
404 track(
i1ohz,nfp1,l)=1.0_r8/track(
i1ohz,nfp1,l)
405 END IF
406 END DO
407# endif
408# else
414# ifdef MASKING
415 &
grid(ng) % rmask, &
416# endif
417 &
ocean(ng) % ubar(:,:,knew), &
418 & my_thread, bounded, track)
419
425# ifdef MASKING
426 &
grid(ng) % rmask, &
427# endif
428 &
ocean(ng) % vbar(:,:,knew), &
429 & my_thread, bounded, track)
430# endif
431
432# ifdef FLOAT_BIOLOGY
433
434
435
436
437
438
439
440
443 & lstr, lend, nfp1,
iftvar(itrc), &
448# ifdef MASKING
449 &
grid(ng) % rmask, &
450# endif
451 &
ocean(ng) % t(:,:,:,nnew,itrc), &
452 & my_thread, bounded, track)
453 END DO
454
455
456
458# endif
459
460
461
462
463
464
465 cff1=9.0_r8/8.0_r8
466 cff2=1.0_r8/8.0_r8
467 cff3=3.0_r8/8.0_r8
468 cff4=6.0_r8/8.0_r8
469 DO l=lstr,lend
470 IF (my_thread(l).and.bounded(l)) THEN
471 track(
ixgrd,nfp1,l)=cff1*track(
ixgrd,nf ,l)- &
472 & cff2*track(
ixgrd,nfm2,l)+ &
473 &
dt(ng)*(cff3*track(
ixrhs,nfp1,l)+ &
474 & cff4*track(
ixrhs,nf ,l)- &
475 & cff3*track(
ixrhs,nfm1,l))
476 track(
iygrd,nfp1,l)=cff1*track(
iygrd,nf ,l)- &
477 & cff2*track(
iygrd,nfm2,l)+ &
478 &
dt(ng)*(cff3*track(
iyrhs,nfp1,l)+ &
479 & cff4*track(
iyrhs,nf ,l)- &
480 & cff3*track(
iyrhs,nfm1,l))
481
482# if defined SOLVE3D && !defined FLOAT_VWALK
483
484
485
487# if defined FLOAT_BIOLOGY
488 track(
izgrd,nfp1,l)=cff1*track(
izgrd,nf ,l)- &
489 & cff2*track(
izgrd,nfm2,l)+ &
490 &
dt(ng)*(cff3*track(
izrhs,nfp1,l)+ &
491 & cff4*track(
izrhs,nf ,l)- &
492 & cff3*track(
izrhs,nfm1,l)+ &
493 & cff3*track(
iwbio,nfp1,l)* &
494 & track(
i1ohz,nfp1,l)+ &
495 & cff4*track(
iwbio,nf ,l)* &
496 & track(
i1ohz,nf ,l)- &
497 & cff3*track(
iwbio,nfm1,l)* &
498 & track(
i1ohz,nf ,l))
499# else
500 track(
izgrd,nfp1,l)=cff1*track(
izgrd,nf ,l)- &
501 & cff2*track(
izgrd,nfm2,l)+ &
502 &
dt(ng)*(cff3*track(
izrhs,nfp1,l)+ &
503 & cff4*track(
izrhs,nf ,l)- &
504 & cff3*track(
izrhs,nfm1,l))
505# endif
506
507
508
509
510
513 ir=int(track(
ixgrd,nfp1,l))
514 jr=int(track(
iygrd,nfp1,l))
515
516 i1=min(max(ir ,0),
lm(ng)+1)
517 i2=min(max(ir+1,1),
lm(ng)+1)
518 j1=min(max(jr ,0),
mm(ng)+1)
519 j2=min(max(jr+1,0),
mm(ng)+1)
520
521 p2=real(i2-i1,r8)*(track(
ixgrd,nfp1,l)-real(i1,r8))
522 q2=real(j2-j1,r8)*(track(
iygrd,nfp1,l)-real(j1,r8))
523 p1=1.0_r8-p2
524 q1=1.0_r8-q2
525# ifdef MASKING
526 cff7=p1*q1*
grid(ng)%z_w(i1,j1,
n(ng))*
grid(ng)%rmask(i1,j1)+ &
527 & p2*q1*
grid(ng)%z_w(i2,j1,
n(ng))*
grid(ng)%rmask(i2,j1)+ &
528 & p1*q2*
grid(ng)%z_w(i1,j2,
n(ng))*
grid(ng)%rmask(i1,j2)+ &
529 & p2*q2*
grid(ng)%z_w(i2,j2,
n(ng))*
grid(ng)%rmask(i2,j2)
530 cff8=p1*q1*
grid(ng)%rmask(i1,j1)+ &
531 & p2*q1*
grid(ng)%rmask(i2,j1)+ &
532 & p1*q2*
grid(ng)%rmask(i1,j2)+ &
533 & p2*q2*
grid(ng)%rmask(i2,j2)
534 IF (cff8.gt.0.0_r8) THEN
535 cff9=cff7/cff8
536 ELSE
537 cff9=0.0_r8
538 END IF
539# else
540 cff9=p1*q1*
grid(ng)%z_w(i1,j1,
n(ng))+ &
541 & p2*q1*
grid(ng)%z_w(i2,j1,
n(ng))+ &
542 & p1*q2*
grid(ng)%z_w(i1,j2,
n(ng))+ &
543 & p2*q2*
grid(ng)%z_w(i2,j2,
n(ng))
544# endif
545 cff6=cff9
546
548 zfloat=fz0(l)
550 zfloat=fz0(l)+cff9
551 END IF
552
554# ifdef MASKING
555 cff7=p1*q1*
grid(ng)%z_w(i1,j1,k)*
grid(ng)%rmask(i1,j1)+ &
556 & p2*q1*
grid(ng)%z_w(i2,j1,k)*
grid(ng)%rmask(i2,j1)+ &
557 & p1*q2*
grid(ng)%z_w(i1,j2,k)*
grid(ng)%rmask(i1,j2)+ &
558 & p2*q2*
grid(ng)%z_w(i2,j2,k)*
grid(ng)%rmask(i2,j2)
559 cff8=p1*q1*
grid(ng)%rmask(i1,j1)+ &
560 & p2*q1*
grid(ng)%rmask(i2,j1)+ &
561 & p1*q2*
grid(ng)%rmask(i1,j2)+ &
562 & p2*q2*
grid(ng)%rmask(i2,j2)
563 cff5=0.0_r8
564 IF (cff8.gt.0.0_r8) cff5=cff7/cff8
565# else
566 cff5=p1*q1*
grid(ng)%z_w(i1,j1,k)+ &
567 & p2*q1*
grid(ng)%z_w(i2,j1,k)+ &
568 & p1*q2*
grid(ng)%z_w(i1,j2,k)+ &
569 & p2*q2*
grid(ng)%z_w(i2,j2,k)
570# endif
571 IF ((zfloat-cff5)*(cff6-zfloat).ge.0.0_r8) THEN
572 track(
izgrd,nfp1,l)=real(k,r8)+(zfloat-cff5)/(cff6-cff5)
573 END IF
574 cff6=cff5
575 END DO
576 END IF
577# endif
578 END IF
579 END DO
580
581
582
583
584
587 DO l=lstr,lend
588 IF (my_thread(l).and.bounded(l)) THEN
589 IF (track(
ixgrd,nfp1,l).ge.real(
lm(ng)+1,r8)-0.5_r8)
THEN
595 ELSE IF (track(
ixgrd,nfp1,l).lt.0.5_r8)
THEN
601 END IF
602 END IF
603 END DO
604# ifdef DISTRIBUTE
606 fwrk=reshape(track,(/npts/))
609 DO l=lstr,lend
610 IF ((xstr.le.track(
ixgrd,nfp1,l)).and. &
611 & (track(
ixgrd,nfp1,l).lt.xend).and. &
612 & (ystr.le.track(
iygrd,nfp1,l)).and. &
613 & (track(
iygrd,nfp1,l).lt.yend))
THEN
614 my_thread(l)=.true.
615 ELSE IF (
master.and.(.not.bounded(l)))
THEN
616 my_thread(l)=.true.
617 ELSE
618 my_thread(l)=.false.
621 track(i,j,l)=fspv
622 END DO
623 END DO
624 END IF
625 END DO
626 END IF
627# endif
628 ELSE
629 DO l=lstr,lend
630 IF (my_thread(l).and.bounded(l)) THEN
631 IF ((track(
ixgrd,nfp1,l).ge.real(
lm(ng)+1,r8)-0.5_r8).or. &
632 & (track(
ixgrd,nfp1,l).lt.0.5_r8))
THEN
633 bounded(l)=.false.
634 END IF
635 END IF
636 END DO
637 END IF
638
641 DO l=lstr,lend
642 IF (my_thread(l).and.bounded(l)) THEN
643 IF (track(
iygrd,nfp1,l).ge.real(
mm(ng)+1,r8)-0.5_r8)
THEN
649 ELSE IF (track(
iygrd,nfp1,l).lt.0.5_r8)
THEN
655 END IF
656 END IF
657 END DO
658# ifdef DISTRIBUTE
660 fwrk=reshape(track,(/npts/))
663 DO l=lstr,lend
664 IF ((xstr.le.track(
ixgrd,nfp1,l)).and. &
665 & (track(
ixgrd,nfp1,l).lt.xend).and. &
666 & (ystr.le.track(
iygrd,nfp1,l)).and. &
667 & (track(
iygrd,nfp1,l).lt.yend))
THEN
668 my_thread(l)=.true.
669 ELSE IF (
master.and.(.not.bounded(l)))
THEN
670 my_thread(l)=.true.
671 ELSE
672 my_thread(l)=.false.
675 track(i,j,l)=fspv
676 END DO
677 END DO
678 END IF
679 END DO
680 END IF
681# endif
682 ELSE
683 DO l=lstr,lend
684 IF (my_thread(l).and.bounded(l)) THEN
685 IF ((track(
iygrd,nfp1,l).ge.real(
mm(ng)+1,r8)-0.5_r8).or. &
686 & (track(
iygrd,nfp1,l).lt.0.5_r8))
THEN
687 bounded(l)=.false.
688 END IF
689 END IF
690 END DO
691 END IF
692
693
694
695
696
697
699
700 DO l=lstr,lend
701 IF (.not.bounded(l).and. &
702 & (
time(ng)-halfdt.le.tinfo(
itstr,l).and. &
703 &
time(ng)+halfdt.gt.tinfo(
itstr,l)))
THEN
704 bounded(l)=.true.
705 IF ((tinfo(
ixgrd,l).lt.0.5_r8).or. &
706 & (tinfo(
iygrd,l).lt.0.5_r8).or. &
707 & (tinfo(
ixgrd,l).gt.real(
lm(ng),r8)+0.5_r8).or. &
708 & (tinfo(
iygrd,l).gt.real(
mm(ng),r8)+0.5_r8))
THEN
709 bounded(l)=.false.
710 END IF
711# if defined SOLVE3D && defined FLOAT_STICKY
712 stuck(l)=.false.
713# endif
714# ifdef DISTRIBUTE
715 IF ((xstr.le.tinfo(
ixgrd,l)).and. &
716 & (tinfo(
ixgrd,l).lt.xend).and. &
717 & (ystr.le.tinfo(
iygrd,l)).and. &
718 & (tinfo(
iygrd,l).lt.yend).and.bounded(l))
THEN
722# ifdef SOLVE3D
724# endif
725 END DO
726 my_thread(l)=.true.
727 ELSE
728 my_thread(l)=.false.
731 track(i,j,l)=fspv
732 END DO
733 END DO
734 END IF
735# else
736 IF (bounded(l)) THEN
737 my_thread(l)=.true.
741# ifdef SOLVE3D
743# endif
744 END DO
745 END IF
746# endif
747 END IF
748 END DO
749
750
751
752
753
754# ifdef SOLVE3D
761# ifdef MASKING
762 &
grid(ng) % rmask, &
763# endif
764 &
ocean(ng) % u(:,:,:,nnew), &
765 & my_thread, bounded, track)
766
773# ifdef MASKING
774 &
grid(ng) % rmask, &
775# endif
776 &
ocean(ng) % v(:,:,:,nnew), &
777 & my_thread, bounded, track)
778
779# if !defined FLOAT_VWALK
786# ifdef MASKING
787 &
grid(ng) % rmask, &
788# endif
790 & my_thread, bounded, track)
791# endif
792# ifdef FLOAT_BIOLOGY
795 &
r3dvar, gmask, fspv, nudg, &
799# ifdef MASKING
800 &
grid(ng) % rmask, &
801# endif
803 & my_thread, bounded, track)
804 DO l=lstr,lend
805 IF (my_thread(l).and.bounded(l)) THEN
806 track(
i1ohz,nfp1,l)=1.0_r8/track(
i1ohz,nfp1,l)
807 END IF
808 END DO
809# endif
810# else
816# ifdef MASKING
817 &
grid(ng) % rmask, &
818# endif
819 &
ocean(ng) % ubar(:,:,knew), &
820 & my_thread, bounded, track)
821
827# ifdef MASKING
828 &
grid(ng) % rmask, &
829# endif
830 &
ocean(ng) % vbar(:,:,knew), &
831 & my_thread, bounded, track)
832# endif
833
834
835
836 DO l=lstr,lend
837 IF (my_thread(l).and.bounded(l).and. &
838 & (
time(ng)-halfdt.le.tinfo(
itstr,l).and. &
839 &
time(ng)+halfdt.gt.tinfo(
itstr,l)))
THEN
840 xrhs=track(
ixrhs,nfp1,l)
841 yrhs=track(
iyrhs,nfp1,l)
842# ifdef SOLVE3D
843 zrhs=track(
izrhs,nfp1,l)
844# endif
845# ifdef FLOAT_BIOLOGY
846 ohz=track(
i1ohz,nfp1,l)
847# endif
849 track(
ixrhs,i,l)=xrhs
850 track(
iyrhs,i,l)=yrhs
851# ifdef SOLVE3D
852 track(
izrhs,i,l)=zrhs
853# endif
854# ifdef FLOAT_BIOLOGY
856# endif
857 END DO
858 END IF
859 END DO
860
861
862
863
864
871# ifdef SOLVE3D
873# endif
874# ifdef MASKING
875 &
grid(ng) % rmask, &
876# endif
878 & my_thread, bounded, track)
879
885# ifdef SOLVE3D
887# endif
888# ifdef MASKING
889 &
grid(ng) % rmask, &
890# endif
892 & my_thread, bounded, track)
893 ELSE
899# ifdef SOLVE3D
901# endif
902# ifdef MASKING
903 &
grid(ng) % rmask, &
904# endif
906 & my_thread, bounded, track)
907
913# ifdef SOLVE3D
915# endif
916# ifdef MASKING
917 &
grid(ng) % rmask, &
918# endif
920 & my_thread, bounded, track)
921 END IF
922# ifdef SOLVE3D
929# ifdef MASKING
930 &
grid(ng) % rmask, &
931# endif
933 & my_thread, bounded, track)
934
941# ifdef MASKING
942 &
grid(ng) % rmask, &
943# endif
945 & my_thread, bounded, track)
946
949 & lstr, lend, nfp1,
iftvar(itrc), &
954# ifdef MASKING
955 &
grid(ng) % rmask, &
956# endif
957 &
ocean(ng) % t(:,:,:,nnew,itrc), &
958 & my_thread, bounded, track)
959 END DO
960# endif
961# ifdef FLOAT_BIOLOGY
962
963
964
965
966
968# endif
969# if defined SOLVE3D && defined FLOAT_VWALK && !defined VWALK_FORWARD
970
971
972
973
974
975
976 CALL vwalk_floats (ng, lstr, lend, .false., my_thread, nudg)
977# endif
978# ifdef SOLVE3D
979# ifdef FLOAT_STICKY
980
981
982
983
984 DO l=lstr,lend
985 IF (my_thread(l).and.bounded(l)) THEN
986 IF (stuck(l)) THEN
990 ELSE
991 IF (track(
izgrd,nfp1,l).gt.real(
n(ng),r8))
THEN
993 track(
izgrd,j,l)=2.0_r8*real(
n(ng),r8)- &
995 END DO
996 ELSE IF (track(
izgrd,nfp1,l).lt.0.0_r8)
THEN
998 track(
izgrd,j,l)=0.0_r8
999 END DO
1000 stuck(l)=.true.
1001 END IF
1002 END IF
1003 END IF
1004 END DO
1005# else
1006
1007
1008
1009 DO l=lstr,lend
1010 IF (my_thread(l).and.bounded(l)) THEN
1011 IF (track(
izgrd,nfp1,l).gt.real(
n(ng),r8))
THEN
1013 track(
izgrd,j,l)=2.0_r8*real(
n(ng),r8)-track(
izgrd,j,l)
1014 END DO
1015 ELSE IF (track(
izgrd,nfp1,l).lt.0.0_r8)
THEN
1018 END DO
1019 END IF
1020 END IF
1021 END DO
1022# endif
1023# endif
1024# ifdef DISTRIBUTE
1025
1026
1027
1028
1029
1030 fwrk=reshape(track,(/npts/))
1033
1034
1035
1036 fwrk=fspv
1038 IF (bounded(l)) THEN
1039 fwrk(l)=1.0_r8
1040 END IF
1041 END DO
1044 IF (fwrk(l).ne.fspv) THEN
1045 bounded(l)=.true.
1046 ELSE
1047 bounded(l)=.false.
1048 END IF
1049 END DO
1050# endif
1051
1052 RETURN
subroutine, public biology_floats(ng, lstr, lend, predictor, my_thread)
subroutine, public interp_floats(ng, lbi, ubi, lbj, ubj, lbk, ubk, lstr, lend, itime, ifield, isbval, gtype, maskit, fspval, nudg, pm, pn, hz, amask, a, my_thread, bounded, track)
integer, parameter flt_isobar
integer, parameter flt_geopot
integer, dimension(:), allocatable iftvar
integer, parameter flt_lagran
type(t_grid), dimension(:), allocatable grid
integer, dimension(:), allocatable istvar
type(t_ocean), dimension(:), allocatable ocean
integer, dimension(:), allocatable nfv
integer, dimension(:), allocatable nfloats
integer, dimension(:), allocatable n
type(t_bounds), dimension(:), allocatable bounds
integer, parameter r3dvar
integer, parameter u3dvar
integer, dimension(:), allocatable lm
integer, parameter u2dvar
integer, dimension(:), allocatable ntilei
integer, parameter w3dvar
integer, dimension(:), allocatable nt
integer, dimension(:), allocatable mm
integer, parameter r2dvar
integer, parameter v2dvar
integer, parameter v3dvar
integer, dimension(:), allocatable ntilej
real(dp), parameter spval
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(dp), dimension(:), allocatable time
subroutine, public vwalk_floats(ng, lstr, lend, predictor, my_thread, nudg)