257
258
259
260
261 integer, intent(in ) :: ng, tile
262 integer, intent(in ) :: LBi, UBi, LBj, UBj, UBk
263 integer, intent(in ) :: IminS, ImaxS, JminS, JmaxS
264 integer, intent(in ) :: krhs, kstp, knew
265#ifdef SOLVE3D
266 integer, intent(in ) :: nstp, nnew
267#endif
268
269#ifdef ASSUMED_SHAPE
270# ifdef MASKING
271 real(r8), intent(in ) :: pmask(LBi:,LBj:)
272 real(r8), intent(in ) :: rmask(LBi:,LBj:)
273 real(r8), intent(in ) :: umask(LBi:,LBj:)
274 real(r8), intent(in ) :: vmask(LBi:,LBj:)
275# endif
276 real(r8), intent(in ) :: fomn(LBi:,LBj:)
277 real(r8), intent(in ) :: h(LBi:,LBj:)
278 real(r8), intent(in ) :: om_u(LBi:,LBj:)
279 real(r8), intent(in ) :: om_v(LBi:,LBj:)
280 real(r8), intent(in ) :: on_u(LBi:,LBj:)
281 real(r8), intent(in ) :: on_v(LBi:,LBj:)
282 real(r8), intent(in ) :: omn(LBi:,LBj:)
283 real(r8), intent(in ) :: pm(LBi:,LBj:)
284 real(r8), intent(in ) :: pn(LBi:,LBj:)
285# if defined CURVGRID && defined UV_ADV
286 real(r8), intent(in ) :: dndx(LBi:,LBj:)
287 real(r8), intent(in ) :: dmde(LBi:,LBj:)
288# endif
289# if defined UV_VIS2 || defined UV_VIS4
290 real(r8), intent(in ) :: pmon_r(LBi:,LBj:)
291 real(r8), intent(in ) :: pnom_r(LBi:,LBj:)
292 real(r8), intent(in ) :: pmon_p(LBi:,LBj:)
293 real(r8), intent(in ) :: pnom_p(LBi:,LBj:)
294 real(r8), intent(in ) :: om_r(LBi:,LBj:)
295 real(r8), intent(in ) :: on_r(LBi:,LBj:)
296 real(r8), intent(in ) :: om_p(LBi:,LBj:)
297 real(r8), intent(in ) :: on_p(LBi:,LBj:)
298# ifdef UV_VIS2
299 real(r8), intent(in ) :: visc2_p(LBi:,LBj:)
300 real(r8), intent(in ) :: visc2_r(LBi:,LBj:)
301# endif
302# ifdef UV_VIS4
303 real(r8), intent(in ) :: visc4_p(LBi:,LBj:)
304 real(r8), intent(in ) :: visc4_r(LBi:,LBj:)
305# endif
306# endif
307# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
308 real(r8), intent(in ) :: tl_bed_thick(LBi:,LBj:,:)
309# endif
310# ifdef WEC_MELLOR
311 real(r8), intent(in ) :: ubar_stokes(LBi:,LBj:)
312 real(r8), intent(in ) :: vbar_stokes(LBi:,LBj:)
313# endif
314# if defined TIDE_GENERATING_FORCES && !defined SOLVE3D
315 real(r8), intent(in ) :: eq_tide(LBi:,LBj:)
316 real(r8), intent(in ) :: tl_eq_tide(LBi:,LBj:)
317# endif
318 real(r8), intent(in ) :: rubar(LBi:,LBj:,:)
319 real(r8), intent(in ) :: rvbar(LBi:,LBj:,:)
320 real(r8), intent(in ) :: rzeta(LBi:,LBj:,:)
321 real(r8), intent(in ) :: ubar(LBi:,LBj:,:)
322 real(r8), intent(in ) :: vbar(LBi:,LBj:,:)
323 real(r8), intent(in ) :: zeta(LBi:,LBj:,:)
324# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
325 real(r8), intent(inout) :: tl_h(LBi:,LBj:)
326# else
327 real(r8), intent(in ) :: tl_h(LBi:,LBj:)
328# endif
329# ifndef SOLVE3D
330 real(r8), intent(in ) :: tl_bustr(LBi:,LBj:)
331 real(r8), intent(in ) :: tl_bvstr(LBi:,LBj:)
332# ifdef ATM_PRESS
333 real(r8), intent(in ) :: Pair(LBi:,LBj:)
334# endif
335# else
336# ifdef VAR_RHO_2D_NOT_YET
337 real(r8), intent(in ) :: rhoA(LBi:,LBj:)
338 real(r8), intent(in ) :: rhoS(LBi:,LBj:)
339 real(r8), intent(in ) :: tl_rhoA(LBi:,LBj:)
340 real(r8), intent(in ) :: tl_rhoS(LBi:,LBj:)
341# endif
342 real(r8), intent(in ) :: Zt_avg1(LBi:,LBj:)
343
344 real(r8), intent(inout) :: tl_DU_avg1(LBi:,LBj:)
345 real(r8), intent(inout) :: tl_DU_avg2(LBi:,LBj:)
346 real(r8), intent(inout) :: tl_DV_avg1(LBi:,LBj:)
347 real(r8), intent(inout) :: tl_DV_avg2(LBi:,LBj:)
348 real(r8), intent(inout) :: tl_Zt_avg1(LBi:,LBj:)
349 real(r8), intent(inout) :: tl_rufrc(LBi:,LBj:)
350 real(r8), intent(inout) :: tl_rvfrc(LBi:,LBj:)
351 real(r8), intent(inout) :: tl_ru(LBi:,LBj:,0:,:)
352 real(r8), intent(inout) :: tl_rv(LBi:,LBj:,0:,:)
353# endif
354# ifdef WEC_MELLOR
355 real(r8), intent(inout) :: tl_rustr2d(LBi:,LBj:)
356 real(r8), intent(inout) :: tl_rvstr2d(LBi:,LBj:)
357 real(r8), intent(inout) :: tl_rulag2d(LBi:,LBj:)
358 real(r8), intent(inout) :: tl_rvlag2d(LBi:,LBj:)
359 real(r8), intent(inout) :: tl_ubar_stokes(LBi:,LBj:)
360 real(r8), intent(inout) :: tl_vbar_stokes(LBi:,LBj:)
361# endif
362# if defined TIDE_GENERATING_FORCES && !defined SOLVE3D
363 real(r8), intent(in ) :: eq_tide(LBi:UBi,LBj:UBj)
364 real(r8), intent(in ) :: tl_eq_tide(LBi:UBi,LBj:UBj)
365# endif
366# ifdef WET_DRY_NOT_YET
367 real(r8), intent(inout) :: pmask_full(LBi:,LBj:)
368 real(r8), intent(inout) :: rmask_full(LBi:,LBj:)
369 real(r8), intent(inout) :: umask_full(LBi:,LBj:)
370 real(r8), intent(inout) :: vmask_full(LBi:,LBj:)
371
372 real(r8), intent(inout) :: pmask_wet(LBi:,LBj:)
373 real(r8), intent(inout) :: rmask_wet(LBi:,LBj:)
374 real(r8), intent(inout) :: umask_wet(LBi:,LBj:)
375 real(r8), intent(inout) :: vmask_wet(LBi:,LBj:)
376# ifdef SOLVE3D
377 real(r8), intent(inout) :: rmask_wet_avg(LBi:,LBj:)
378# endif
379# endif
380# ifdef DIAGNOSTICS_UV
381
382
383
384
385# ifdef SOLVE3D
386
387
388
389
390# endif
391# endif
392 real(r8), intent(inout) :: tl_rubar(LBi:,LBj:,:)
393 real(r8), intent(inout) :: tl_rvbar(LBi:,LBj:,:)
394 real(r8), intent(inout) :: tl_rzeta(LBi:,LBj:,:)
395 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
396 real(r8), intent(inout) :: tl_vbar(LBi:,LBj:,:)
397 real(r8), intent(inout) :: tl_zeta(LBi:,LBj:,:)
398
399#else
400
401# ifdef MASKING
402 real(r8), intent(in ) :: pmask(LBi:UBi,LBj:UBj)
403 real(r8), intent(in ) :: rmask(LBi:UBi,LBj:UBj)
404 real(r8), intent(in ) :: umask(LBi:UBi,LBj:UBj)
405 real(r8), intent(in ) :: vmask(LBi:UBi,LBj:UBj)
406# endif
407 real(r8), intent(in ) :: fomn(LBi:UBi,LBj:UBj)
408 real(r8), intent(in ) :: h(LBi:UBi,LBj:UBj)
409 real(r8), intent(in ) :: om_u(LBi:UBi,LBj:UBj)
410 real(r8), intent(in ) :: om_v(LBi:UBi,LBj:UBj)
411 real(r8), intent(in ) :: on_u(LBi:UBi,LBj:UBj)
412 real(r8), intent(in ) :: on_v(LBi:UBi,LBj:UBj)
413 real(r8), intent(in ) :: omn(LBi:UBi,LBj:UBj)
414 real(r8), intent(in ) :: pm(LBi:UBi,LBj:UBj)
415 real(r8), intent(in ) :: pn(LBi:UBi,LBj:UBj)
416# if defined CURVGRID && defined UV_ADV
417 real(r8), intent(in ) :: dndx(LBi:UBi,LBj:UBj)
418 real(r8), intent(in ) :: dmde(LBi:UBi,LBj:UBj)
419# endif
420# if defined UV_VIS2 || defined UV_VIS4
421 real(r8), intent(in ) :: pmon_r(LBi:UBi,LBj:UBj)
422 real(r8), intent(in ) :: pnom_r(LBi:UBi,LBj:UBj)
423 real(r8), intent(in ) :: pmon_p(LBi:UBi,LBj:UBj)
424 real(r8), intent(in ) :: pnom_p(LBi:UBi,LBj:UBj)
425 real(r8), intent(in ) :: om_r(LBi:UBi,LBj:UBj)
426 real(r8), intent(in ) :: on_r(LBi:UBi,LBj:UBj)
427 real(r8), intent(in ) :: om_p(LBi:UBi,LBj:UBj)
428 real(r8), intent(in ) :: on_p(LBi:UBi,LBj:UBj)
429# ifdef UV_VIS2
430 real(r8), intent(in ) :: visc2_p(LBi:UBi,LBj:UBj)
431 real(r8), intent(in ) :: visc2_r(LBi:UBi,LBj:UBj)
432# endif
433# ifdef UV_VIS4
434 real(r8), intent(in ) :: visc4_p(LBi:UBi,LBj:UBj)
435 real(r8), intent(in ) :: visc4_r(LBi:UBi,LBj:UBj)
436# endif
437# endif
438# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
439 real(r8), intent(in ) :: tl_bed_thick(LBi:UBi,LBj:UBj,3)
440# endif
441# ifdef WEC_MELLOR
442 real(r8), intent(in ) :: ubar_stokes(LBi:UBi,LBj:UBj)
443 real(r8), intent(in ) :: vbar_stokes(LBi:UBi,LBj:UBj)
444# endif
445# if defined TIDE_GENERATING_FORCES && !defined SOLVE3D
446 real(r8), intent(in ) :: eq_tide(LBi:UBi,LBj:UBj)
447 real(r8), intent(in ) :: tl_eq_tide(LBi:UBi,LBj:UBj)
448# endif
449 real(r8), intent(in ) :: rubar(LBi:UBi,LBj:UBj,2)
450 real(r8), intent(in ) :: rvbar(LBi:UBi,LBj:UBj,2)
451 real(r8), intent(in ) :: rzeta(LBi:UBi,LBj:UBj,2)
452 real(r8), intent(in ) :: ubar(LBi:UBi,LBj:UBj,:)
453 real(r8), intent(in ) :: vbar(LBi:UBi,LBj:UBj,:)
454 real(r8), intent(in ) :: zeta(LBi:UBi,LBj:UBj,:)
455# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
456 real(r8), intent(inout) :: tl_h(LBi:UBi,LBj:UBj)
457# else
458 real(r8), intent(in ) :: tl_h(LBi:UBi,LBj:UBj)
459# endif
460# ifndef SOLVE3D
461 real(r8), intent(in ) :: tl_bustr(LBi:UBi,LBj:UBj)
462 real(r8), intent(in ) :: tl_bvstr(LBi:UBi,LBj:UBj)
463# ifdef ATM_PRESS
464 real(r8), intent(in ) :: Pair(LBi:UBi,LBj:UBj)
465# endif
466# else
467# ifdef VAR_RHO_2D_NOT_YET
468 real(r8), intent(in ) :: rhoA(LBi:UBi,LBj:UBj)
469 real(r8), intent(in ) :: rhoS(LBi:UBi,LBj:UBj)
470 real(r8), intent(in ) :: tl_rhoA(LBi:UBi,LBj:UBj)
471 real(r8), intent(in ) :: tl_rhoS(LBi:UBi,LBj:UBj)
472# endif
473 real(r8), intent(in ) :: Zt_avg1(LBi:UBi,LBj:UBj)
474
475 real(r8), intent(inout) :: tl_DU_avg1(LBi:UBi,LBj:UBj)
476 real(r8), intent(inout) :: tl_DU_avg2(LBi:UBi,LBj:UBj)
477 real(r8), intent(inout) :: tl_DV_avg1(LBi:UBi,LBj:UBj)
478 real(r8), intent(inout) :: tl_DV_avg2(LBi:UBi,LBj:UBj)
479 real(r8), intent(inout) :: tl_Zt_avg1(LBi:UBi,LBj:UBj)
480 real(r8), intent(inout) :: tl_rufrc(LBi:UBi,LBj:UBj)
481 real(r8), intent(inout) :: tl_rvfrc(LBi:UBi,LBj:UBj)
482 real(r8), intent(inout) :: tl_ru(LBi:UBi,LBj:UBj,0:UBk,2)
483 real(r8), intent(inout) :: tl_rv(LBi:UBi,LBj:UBj,0:UBk,2)
484# endif
485# ifdef WEC_MELLOR
486 real(r8), intent(inout) :: tl_rustr2d(LBi:UBi,LBj:UBj)
487 real(r8), intent(inout) :: tl_rvstr2d(LBi:UBi,LBj:UBj)
488 real(r8), intent(inout) :: tl_rulag2d(LBi:UBi,LBj:UBj)
489 real(r8), intent(inout) :: tl_rvlag2d(LBi:UBi,LBj:UBj)
490 real(r8), intent(inout) :: tl_ubar_stokes(LBi:UBi,LBj:UBj)
491 real(r8), intent(inout) :: tl_vbar_stokes(LBi:UBi,LBj:UBj)
492# endif
493# ifdef WET_DRY_NOT_YET
494 real(r8), intent(inout) :: pmask_full(LBi:UBi,LBj:UBj)
495 real(r8), intent(inout) :: rmask_full(LBi:UBi,LBj:UBj)
496 real(r8), intent(inout) :: umask_full(LBi:UBi,LBj:UBj)
497 real(r8), intent(inout) :: vmask_full(LBi:UBi,LBj:UBj)
498
499 real(r8), intent(inout) :: pmask_wet(LBi:UBi,LBj:UBj)
500 real(r8), intent(inout) :: rmask_wet(LBi:UBi,LBj:UBj)
501 real(r8), intent(inout) :: umask_wet(LBi:UBi,LBj:UBj)
502 real(r8), intent(inout) :: vmask_wet(LBi:UBi,LBj:UBj)
503# ifdef SOLVE3D
504 real(r8), intent(inout) :: rmask_wet_avg(LBi:UBi,LBj:UBj)
505# endif
506# endif
507# ifdef DIAGNOSTICS_UV
508
509
510
511
512# ifdef SOLVE3D
513
514
515
516
517# endif
518# endif
519 real(r8), intent(inout) :: tl_rubar(LBi:UBi,LBj:UBj,2)
520 real(r8), intent(inout) :: tl_rvbar(LBi:UBi,LBj:UBj,2)
521 real(r8), intent(inout) :: tl_rzeta(LBi:UBi,LBj:UBj,2)
522 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
523 real(r8), intent(inout) :: tl_vbar(LBi:UBi,LBj:UBj,:)
524 real(r8), intent(inout) :: tl_zeta(LBi:UBi,LBj:UBj,:)
525#endif
526
527
528
529 logical :: CORRECTOR_2D_STEP
530
531 integer :: i, is, j, ptsk
532#ifdef DIAGNOSTICS_UV
533
534#endif
535
536 real(r8) :: cff, cff1, cff2, cff3, cff4, cff5, cff6, cff7
537 real(r8) :: fac, fac1, fac2, fac3
538 real(r8) :: tl_cff, tl_cff1, tl_cff2, tl_cff3, tl_cff4
539 real(r8) :: tl_fac, tl_fac1
540
541 real(r8), parameter :: IniVal = 0.0_r8
542
543 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Dgrad
544 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Dnew
545 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Drhs
546 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Drhs_p
547 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Dstp
548 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: DUon
549 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: DVom
550#ifdef WEC_MELLOR
551 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: DUSon
552 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: DVSom
553#endif
554#ifdef UV_VIS4
555 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: LapU
556 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: LapV
557#endif
558 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFe
559 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFx
560 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFe
561 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFx
562 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: grad
563 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: gzeta
564 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: gzeta2
565#if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
566 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: gzetaSA
567#endif
568 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: rhs_ubar
569 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: rhs_vbar
570 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: rhs_zeta
571 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: zeta_new
572 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: zwrk
573#ifdef WET_DRY_NOT_YET
574 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: wetdry
575#endif
576#ifdef DIAGNOSTICS_UV
577
578
579
580
581#endif
582
583 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Dgrad
584 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Dnew
585 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Drhs
586 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Drhs_p
587 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Dstp
588 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_DUon
589 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_DVom
590#ifdef WEC_MELLOR
591 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_DUSon
592 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_DVSom
593#endif
594#ifdef UV_VIS4
595 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_LapU
596 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_LapV
597#endif
598 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_UFe
599 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_UFx
600 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_VFe
601 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_VFx
602 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_grad
603 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_gzeta
604 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_gzeta2
605#if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
606 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_gzetaSA
607#endif
608 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_rhs_ubar
609 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_rhs_vbar
610 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_rhs_zeta
611 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_zeta_new
612 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_zwrk
613#ifdef WET_DRY_NOT_YET
614 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_wetdry
615#endif
616
617#include "set_bounds.h"
618
619 ptsk=3-kstp
620 corrector_2d_step=.not.predictor_2d_step(ng)
621#ifdef DEBUG
622 IF (master) THEN
623 WRITE (20,20) iic(ng), predictor_2d_step(ng), &
624 & kstp, krhs, knew, ptsk
625 20 FORMAT (' iic = ',i5.5,' predictor = ',l1,' kstp = ',i1, &
626 & ' krhs = ',i1,' knew = ',i1,' ptsk = ',i1)
627 END IF
628#endif
629#ifdef INITIALIZE_AUTOMATIC
630
631
632
633
634
635
636 DO j=jmins,jmaxs
637 DO i=imins,imaxs
638 dgrad(i,j)=inival
639 dnew(i,j)=inival
640 drhs(i,j)=inival
641 drhs_p(i,j)=inival
642 dstp(i,j)=inival
643 duon(i,j)=inival
644 dvom(i,j)=inival
645# ifdef WEC_MELLOR
646 duson(i,j)=inival
647 dvsom(i,j)=inival
648# endif
649# ifdef UV_VIS4
650 lapu(i,j)=inival
651 lapv(i,j)=inival
652# endif
653 ufe(i,j)=inival
654 ufx(i,j)=inival
655 vfe(i,j)=inival
656 vfx(i,j)=inival
657 grad(i,j)=inival
658 gzeta(i,j)=inival
659 gzeta2(i,j)=inival
660# if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
661 gzetasa(i,j)=inival
662# endif
663 rhs_ubar(i,j)=inival
664 rhs_vbar(i,j)=inival
665 rhs_zeta(i,j)=inival
666 zeta_new(i,j)=inival
667 zwrk(i,j)=inival
668
669 tl_dgrad(i,j)=inival
670 tl_dnew(i,j)=inival
671 tl_drhs(i,j)=inival
672 tl_drhs_p(i,j)=inival
673 tl_dstp(i,j)=inival
674 tl_duon(i,j)=inival
675 tl_dvom(i,j)=inival
676# ifdef WEC_MELLOR
677 tl_duson(i,j)=inival
678 tl_dvsom(i,j)=inival
679# endif
680# ifdef UV_VIS4
681 tl_lapu(i,j)=inival
682 tl_lapv(i,j)=inival
683# endif
684 tl_ufe(i,j)=inival
685 tl_ufx(i,j)=inival
686 tl_vfe(i,j)=inival
687 tl_vfx(i,j)=inival
688 tl_grad(i,j)=inival
689 tl_gzeta(i,j)=inival
690 tl_gzeta2(i,j)=inival
691# if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
692 tl_gzetasa(i,j)=inival
693# endif
694 tl_rhs_ubar(i,j)=inival
695 tl_rhs_vbar(i,j)=inival
696 tl_rhs_zeta(i,j)=inival
697 tl_zeta_new(i,j)=inival
698 tl_zwrk(i,j)=inival
699 END DO
700 END DO
701#endif
702
703
704
705
706
707#if defined DISTRIBUTE && !defined NESTING
708
709
710
711
712
713
714
715
716 DO j=jstrv-2,jendp2
717 DO i=istru-2,iendp2
718 dnew(i,j)=zeta(i,j,knew)+h(i,j)
719 drhs(i,j)=zeta(i,j,krhs)+h(i,j)
720 tl_drhs(i,j)=tl_zeta(i,j,krhs)+tl_h(i,j)
721 END DO
722 END DO
723 DO j=jstrv-2,jendp2
724 DO i=istru-1,iendp2
725 cff=0.5_r8*on_u(i,j)
726 cff1=cff*(drhs(i,j)+drhs(i-1,j))
727 tl_cff1=cff*(tl_drhs(i,j)+tl_drhs(i-1,j))
728 duon(i,j)=ubar(i,j,krhs)*cff1
729 tl_duon(i,j)=tl_ubar(i,j,krhs)*cff1+ &
730 & ubar(i,j,krhs)*tl_cff1
731# ifdef WEC_MELLOR
732 duson(i,j)=ubar_stokes(i,j)*cff1
733 tl_duson(i,j)=tl_ubar_stokes(i,j)*cff1+ &
734 & ubar_stokes(i,j)*tl_cff1
735 duon(i,j)=duon(i,j)+duson(i,j)
736 tl_duon(i,j)=tl_duon(i,j)+tl_duson(i,j)
737# endif
738 END DO
739 END DO
740 DO j=jstrv-1,jendp2
741 DO i=istru-2,iendp2
742 cff=0.5_r8*om_v(i,j)
743 cff1=cff*(drhs(i,j)+drhs(i,j-1))
744 tl_cff1=cff*(tl_drhs(i,j)+tl_drhs(i,j-1))
745 dvom(i,j)=vbar(i,j,krhs)*cff1
746 tl_dvom(i,j)=tl_vbar(i,j,krhs)*cff1+ &
747 & vbar(i,j,krhs)*tl_cff1
748# ifdef WEC_MELLOR
749 dvsom(i,j)=vbar_stokes(i,j)*cff1
750 tl_dvsom(i,j)=tl_vbar_stokes(i,j)*cff1+ &
751 & vbar_stokes(i,j)*tl_cff1
752 dvom(i,j)=dvom(i,j)+dvsom(i,j)
753 tl_dvom(i,j)=tl_dvom(i,j)+tl_dvsom(i,j)
754# endif
755 END DO
756 END DO
757
758#else
759
760 DO j=jstrvm2-1,jendp2
761 DO i=istrum2-1,iendp2
762 dnew(i,j)=zeta(i,j,knew)+h(i,j)
763 drhs(i,j)=zeta(i,j,krhs)+h(i,j)
764 tl_drhs(i,j)=tl_zeta(i,j,krhs)+tl_h(i,j)
765 END DO
766 END DO
767 DO j=jstrvm2-1,jendp2
768 DO i=istrum2,iendp2
769 cff=0.5_r8*on_u(i,j)
770 cff1=cff*(drhs(i,j)+drhs(i-1,j))
771 tl_cff1=cff*(tl_drhs(i,j)+tl_drhs(i-1,j))
772 duon(i,j)=ubar(i,j,krhs)*cff1
773 tl_duon(i,j)=tl_ubar(i,j,krhs)*cff1+ &
774 & ubar(i,j,krhs)*tl_cff1
775# ifdef WEC_MELLOR
776 duson(i,j)=ubar_stokes(i,j)*cff1
777 tl_duson(i,j)=tl_ubar_stokes(i,j)*cff1+ &
778 & ubar_stokes(i,j)*tl_cff1
779 duon(i,j)=duon(i,j)+duson(i,j)
780 tl_duon(i,j)=tl_duon(i,j)+tl_duson(i,j)
781# endif
782 END DO
783 END DO
784 DO j=jstrvm2,jendp2
785 DO i=istrum2-1,iendp2
786 cff=0.5_r8*om_v(i,j)
787 cff1=cff*(drhs(i,j)+drhs(i,j-1))
788 tl_cff1=cff*(tl_drhs(i,j)+tl_drhs(i,j-1))
789 dvom(i,j)=vbar(i,j,krhs)*cff1
790 tl_dvom(i,j)=tl_vbar(i,j,krhs)*cff1+ &
791 & vbar(i,j,krhs)*tl_cff1
792# ifdef WEC_MELLOR
793 dvsom(i,j)=vbar_stokes(i,j)*cff1
794 tl_dvsom(i,j)=tl_vbar_stokes(i,j)*cff1+ &
795 & vbar_stokes(i,j)*tl_cff1
796 dvom(i,j)=dvom(i,j)+dvsom(i,j)
797 tl_dvom(i,j)=tl_dvom(i,j)+tl_dvsom(i,j)
798# endif
799 END DO
800 END DO
801#endif
802#ifdef DISTRIBUTE
803
804 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
805 CALL exchange_u2d_tile (ng, tile, &
806 & imins, imaxs, jmins, jmaxs, &
807 & duon)
808 CALL exchange_u2d_tile (ng, tile, &
809 & imins, imaxs, jmins, jmaxs, &
810 & tl_duon)
811 CALL exchange_v2d_tile (ng, tile, &
812 & imins, imaxs, jmins, jmaxs, &
813 & dvom)
814 CALL exchange_v2d_tile (ng, tile, &
815 & imins, imaxs, jmins, jmaxs, &
816 & tl_dvom)
817 END IF
818
819 CALL mp_exchange2d (ng, tile, itlm, 2, &
820 & imins, imaxs, jmins, jmaxs, &
821 & nghostpoints, &
822 & ewperiodic(ng), nsperiodic(ng), &
823 & duon, dvom)
824 CALL mp_exchange2d (ng, tile, itlm, 2, &
825 & imins, imaxs, jmins, jmaxs, &
826 & nghostpoints, &
827 & ewperiodic(ng), nsperiodic(ng), &
828 & tl_duon, tl_dvom)
829#endif
830#if !defined FORWARD_RHS
831
832
833
834 DO j=jstr,jend
835 DO i=istr,iend
836 rhs_ubar(i,j)=0.0_r8
837 rhs_vbar(i,j)=0.0_r8
838 END DO
839 END DO
840#endif
841
842
843
844
845
846 IF (any(tl_volcons(:,ng))) THEN
847 CALL obc_flux_tile (ng, tile, &
848 & lbi, ubi, lbj, ubj, &
849 & imins, imaxs, jmins, jmaxs, &
850 & knew, &
851#ifdef MASKING
852 & umask, vmask, &
853#endif
854 & h, om_v, on_u, &
855 & ubar, vbar, zeta)
856
857
858
859
860 CALL set_duv_bc_tile (ng, tile, &
861 & lbi, ubi, lbj, ubj, &
862 & imins, imaxs, jmins, jmaxs, &
863 & krhs, &
864#ifdef MASKING
865 & umask, vmask, &
866#endif
867 & om_v, on_u, &
868 & ubar, vbar, &
869 & drhs, duon, dvom)
870 CALL tl_set_duv_bc_tile (ng, tile, &
871 & lbi, ubi, lbj, ubj, &
872 & imins, imaxs, jmins, jmaxs, &
873 & krhs, &
874#ifdef MASKING
875 & umask, vmask, &
876#endif
877 & om_v, on_u, ubar, vbar, &
878 & tl_ubar, tl_vbar, &
879 & drhs, duon, dvom, &
880 & tl_drhs, tl_duon, tl_dvom)
881 END IF
882#ifdef SOLVE3D
883
884
885
886
887
888 IF (predictor_2d_step(ng)) THEN
889 IF (first_2d_step) THEN
890
891
892
893 cff2=(-1.0_r8/12.0_r8)*weight(2,iif(ng)+1,ng)
894 DO j=jstrr,jendr
895 DO i=istrr,iendr
896
897
898 tl_zt_avg1(i,j)=0.0_r8
899 END DO
900 DO i=istr,iendr
901
902
903 tl_du_avg1(i,j)=0.0_r8
904
905
906 tl_du_avg2(i,j)=cff2*tl_duon(i,j)
907 END DO
908 END DO
909 DO j=jstr,jendr
910 DO i=istrr,iendr
911
912
913 tl_dv_avg1(i,j)=0.0_r8
914
915
916 tl_dv_avg2(i,j)=cff2*tl_dvom(i,j)
917 END DO
918 END DO
919 ELSE
920
921
922
923
924
925 cff1=weight(1,iif(ng)-1,ng)
926 cff2=(8.0_r8/12.0_r8)*weight(2,iif(ng) ,ng)- &
927 & (1.0_r8/12.0_r8)*weight(2,iif(ng)+1,ng)
928 DO j=jstrr,jendr
929 DO i=istrr,iendr
930
931
932 tl_zt_avg1(i,j)=tl_zt_avg1(i,j)+cff1*tl_zeta(i,j,krhs)
933 END DO
934 DO i=istr,iendr
935
936
937 tl_du_avg1(i,j)=tl_du_avg1(i,j)+cff1*tl_duon(i,j)
938# ifdef WEC_MELLOR
939
940
941 tl_du_avg1(i,j)=tl_du_avg1(i,j)-cff1*tl_duson(i,j)
942# endif
943
944
945 tl_du_avg2(i,j)=tl_du_avg2(i,j)+cff2*tl_duon(i,j)
946 END DO
947 END DO
948 DO j=jstr,jendr
949 DO i=istrr,iendr
950
951
952 tl_dv_avg1(i,j)=tl_dv_avg1(i,j)+cff1*tl_dvom(i,j)
953# ifdef WEC_MELLOR
954
955
956 tl_dv_avg1(i,j)=tl_dv_avg1(i,j)-cff1*tl_dvsom(i,j)
957# endif
958
959
960 tl_dv_avg2(i,j)=tl_dv_avg2(i,j)+cff2*tl_dvom(i,j)
961 END DO
962 END DO
963 END IF
964 ELSE
965 IF (first_2d_step) THEN
966 cff2=weight(2,iif(ng),ng)
967 ELSE
968 cff2=(5.0_r8/12.0_r8)*weight(2,iif(ng),ng)
969 END IF
970 DO j=jstrr,jendr
971 DO i=istr,iendr
972
973
974 tl_du_avg2(i,j)=tl_du_avg2(i,j)+cff2*tl_duon(i,j)
975 END DO
976 END DO
977 DO j=jstr,jendr
978 DO i=istrr,iendr
979
980
981 tl_dv_avg2(i,j)=tl_dv_avg2(i,j)+cff2*tl_dvom(i,j)
982 END DO
983 END DO
984 END IF
985
986
987
988# ifdef NESTING
989
990
991
992
993# endif
994
995 IF ((iif(ng).eq.(nfast(ng)+1)).and.predictor_2d_step(ng)) THEN
996 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
997
998
999
1000
1001 CALL exchange_r2d_tile (ng, tile, &
1002 & lbi, ubi, lbj, ubj, &
1003 & tl_zt_avg1)
1004
1005
1006
1007
1008 CALL exchange_u2d_tile (ng, tile, &
1009 & lbi, ubi, lbj, ubj, &
1010 & tl_du_avg1)
1011
1012
1013
1014
1015 CALL exchange_v2d_tile (ng, tile, &
1016 & lbi, ubi, lbj, ubj, &
1017 & tl_dv_avg1)
1018# ifdef NESTING
1019
1020
1021
1022
1023 CALL exchange_u2d_tile (ng, tile, &
1024 & lbi, ubi, lbj, ubj, &
1025 & tl_du_avg2)
1026
1027
1028
1029
1030 CALL exchange_v2d_tile (ng, tile, &
1031 & lbi, ubi, lbj, ubj, &
1032 & tl_dv_avg2)
1033# endif
1034 END IF
1035
1036# ifdef DISTRIBUTE
1037
1038
1039
1040
1041
1042
1043 CALL mp_exchange2d (ng, tile, itlm, 3, &
1044 & lbi, ubi, lbj, ubj, &
1045 & nghostpoints, &
1046 & ewperiodic(ng), nsperiodic(ng), &
1047 & tl_zt_avg1, tl_du_avg1, tl_dv_avg1)
1048# ifdef NESTING
1049
1050
1051
1052
1053
1054
1055 CALL mp_exchange2d (ng, tile, itlm, 2, &
1056 & lbi, ubi, lbj, ubj, &
1057 & nghostpoints, &
1058 & ewperiodic(ng), nsperiodic(ng), &
1059 & tl_du_avg2, tl_dv_avg2)
1060# endif
1061# endif
1062 END IF
1063#endif
1064#ifdef WET_DRY_NOT_YET
1065
1066
1067
1068
1069
1070
1071
1072
1073# ifdef MASKING
1074
1075# endif
1076
1077# ifdef SOLVE3D
1078
1079
1080# endif
1081
1082
1083
1084
1085
1086
1087
1088#endif
1089
1090
1091
1092
1093 IF (iif(ng).gt.nfast(ng)) RETURN
1094
1095
1096
1097
1098
1099
1100
1101
1102
1103#if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
1104 fac=1000.0_r8/rho0
1105#endif
1106#if defined STOCHASTIC_OPT && !defined STOCH_OPT_WHITE && \
1107
1108 IF (first_2d_step.and.soinitial(ng)) THEN
1109#else
1110 IF (first_2d_step) THEN
1111#endif
1112 cff1=dtfast(ng)
1113 DO j=jstrv-1,jend
1114 DO i=istru-1,iend
1115
1116
1117
1118 tl_rhs_zeta(i,j)=(tl_duon(i,j)-tl_duon(i+1,j))+ &
1119 & (tl_dvom(i,j)-tl_dvom(i,j+1))
1120
1121
1122
1123
1124 zeta_new(i,j)=zeta(i,j,knew)
1125 tl_zeta_new(i,j)=tl_zeta(i,j,kstp)+ &
1126 & pm(i,j)*pn(i,j)*cff1*tl_rhs_zeta(i,j)
1127#ifdef MASKING
1128 zeta_new(i,j)=zeta_new(i,j)*rmask(i,j)
1129 tl_zeta_new(i,j)=tl_zeta_new(i,j)*rmask(i,j)
1130#endif
1131
1132
1133 tl_dnew(i,j)=tl_zeta_new(i,j)+tl_h(i,j)
1134
1135 zwrk(i,j)=0.5_r8*(zeta(i,j,kstp)+zeta_new(i,j))
1136 tl_zwrk(i,j)=0.5_r8*(tl_zeta(i,j,kstp)+tl_zeta_new(i,j))
1137#if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
1138 gzeta(i,j)=(fac+rhos(i,j))*zwrk(i,j)
1139 tl_gzeta(i,j)=(fac+rhos(i,j))*tl_zwrk(i,j)+ &
1140 & tl_rhos(i,j)*zwrk(i,j)
1141 gzeta2(i,j)=gzeta(i,j)*zwrk(i,j)
1142 tl_gzeta2(i,j)=tl_gzeta(i,j)*zwrk(i,j)+ &
1143 & gzeta(i,j)*tl_zwrk(i,j)
1144 gzetasa(i,j)=zwrk(i,j)*(rhos(i,j)-rhoa(i,j))
1145 tl_gzetasa(i,j)=tl_zwrk(i,j)*(rhos(i,j)-rhoa(i,j))+ &
1146 & zwrk(i,j)*(tl_rhos(i,j)-tl_rhoa(i,j))
1147#else
1148 gzeta(i,j)=zwrk(i,j)
1149 tl_gzeta(i,j)=tl_zwrk(i,j)
1150 gzeta2(i,j)=zwrk(i,j)*zwrk(i,j)
1151 tl_gzeta2(i,j)=2.0_r8*tl_zwrk(i,j)*zwrk(i,j)
1152#endif
1153 END DO
1154 END DO
1155 ELSE IF (predictor_2d_step(ng)) THEN
1156 cff1=2.0_r8*dtfast(ng)
1157 cff4=4.0_r8/25.0_r8
1158 cff5=1.0_r8-2.0_r8*cff4
1159 DO j=jstrv-1,jend
1160 DO i=istru-1,iend
1161
1162
1163
1164 tl_rhs_zeta(i,j)=(tl_duon(i,j)-tl_duon(i+1,j))+ &
1165 & (tl_dvom(i,j)-tl_dvom(i,j+1))
1166
1167
1168
1169
1170 zeta_new(i,j)=zeta(i,j,knew)
1171 tl_zeta_new(i,j)=tl_zeta(i,j,kstp)+ &
1172 & pm(i,j)*pn(i,j)*cff1*tl_rhs_zeta(i,j)
1173#ifdef MASKING
1174 zeta_new(i,j)=zeta_new(i,j)*rmask(i,j)
1175 tl_zeta_new(i,j)=tl_zeta_new(i,j)*rmask(i,j)
1176#endif
1177
1178
1179 tl_dnew(i,j)=tl_zeta_new(i,j)+tl_h(i,j)
1180
1181 zwrk(i,j)=cff5*zeta(i,j,krhs)+ &
1182 & cff4*(zeta(i,j,kstp)+zeta_new(i,j))
1183 tl_zwrk(i,j)=cff5*tl_zeta(i,j,krhs)+ &
1184 & cff4*(tl_zeta(i,j,kstp)+tl_zeta_new(i,j))
1185#if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
1186 gzeta(i,j)=(fac+rhos(i,j))*zwrk(i,j)
1187 tl_gzeta(i,j)=(fac+rhos(i,j))*tl_zwrk(i,j)+ &
1188 & tl_rhos(i,j)*zwrk(i,j)
1189 gzeta2(i,j)=gzeta(i,j)*zwrk(i,j)
1190 tl_gzeta2(i,j)=tl_gzeta(i,j)*zwrk(i,j)+ &
1191 & gzeta(i,j)*tl_zwrk(i,j)
1192 gzetasa(i,j)=zwrk(i,j)*(rhos(i,j)-rhoa(i,j))
1193 tl_gzetasa(i,j)=tl_zwrk(i,j)*(rhos(i,j)-rhoa(i,j))+ &
1194 & zwrk(i,j)*(tl_rhos(i,j)-tl_rhoa(i,j))
1195#else
1196 gzeta(i,j)=zwrk(i,j)
1197 tl_gzeta(i,j)=tl_zwrk(i,j)
1198 gzeta2(i,j)=zwrk(i,j)*zwrk(i,j)
1199 tl_gzeta2(i,j)=2.0_r8*tl_zwrk(i,j)*zwrk(i,j)
1200#endif
1201 END DO
1202 END DO
1203 ELSE IF (corrector_2d_step) THEN
1204 cff1=dtfast(ng)*5.0_r8/12.0_r8
1205 cff2=dtfast(ng)*8.0_r8/12.0_r8
1206 cff3=dtfast(ng)*1.0_r8/12.0_r8
1207 cff4=2.0_r8/5.0_r8
1208 cff5=1.0_r8-cff4
1209 DO j=jstrv-1,jend
1210 DO i=istru-1,iend
1211
1212
1213
1214 tl_cff=cff1*((tl_duon(i,j)-tl_duon(i+1,j))+ &
1215 & (tl_dvom(i,j)-tl_dvom(i,j+1)))
1216
1217
1218
1219
1220
1221
1222 zeta_new(i,j)=zeta(i,j,knew)
1223 tl_zeta_new(i,j)=tl_zeta(i,j,kstp)+ &
1224 & pm(i,j)*pn(i,j)*(tl_cff+ &
1225 & cff2*tl_rzeta(i,j,kstp)- &
1226 & cff3*tl_rzeta(i,j,ptsk))
1227#ifdef MASKING
1228 zeta_new(i,j)=zeta_new(i,j)*rmask(i,j)
1229 tl_zeta_new(i,j)=tl_zeta_new(i,j)*rmask(i,j)
1230#endif
1231
1232
1233 tl_dnew(i,j)=tl_zeta_new(i,j)+tl_h(i,j)
1234
1235 zwrk(i,j)=cff5*zeta_new(i,j)+cff4*zeta(i,j,krhs)
1236 tl_zwrk(i,j)=cff5*tl_zeta_new(i,j)+cff4*tl_zeta(i,j,krhs)
1237#if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
1238 gzeta(i,j)=(fac+rhos(i,j))*zwrk(i,j)
1239 tl_gzeta(i,j)=(fac+rhos(i,j))*tl_zwrk(i,j)+ &
1240 & tl_rhos(i,j)*zwrk(i,j)
1241 gzeta2(i,j)=gzeta(i,j)*zwrk(i,j)
1242 tl_gzeta2(i,j)=tl_gzeta(i,j)*zwrk(i,j)+ &
1243 & gzeta(i,j)*tl_zwrk(i,j)
1244 gzetasa(i,j)=zwrk(i,j)*(rhos(i,j)-rhoa(i,j))
1245 tl_gzetasa(i,j)=tl_zwrk(i,j)*(rhos(i,j)-rhoa(i,j))+ &
1246 & zwrk(i,j)*(tl_rhos(i,j)-tl_rhoa(i,j))
1247#else
1248 gzeta(i,j)=zwrk(i,j)
1249 tl_gzeta(i,j)=tl_zwrk(i,j)
1250 gzeta2(i,j)=zwrk(i,j)*zwrk(i,j)
1251 tl_gzeta2(i,j)=2.0_r8*tl_zwrk(i,j)*zwrk(i,j)
1252#endif
1253 END DO
1254 END DO
1255 END IF
1256
1257
1258
1259#ifdef WET_DRY_NOT_YET
1260
1261
1262#endif
1263
1264 DO j=jstr,jend
1265 DO i=istr,iend
1266
1267
1268 tl_zeta(i,j,knew)=tl_zeta_new(i,j)
1269#if defined WET_DRY_NOT_YET && defined MASKING
1270
1271
1272
1273 tl_zeta(i,j,knew)=tl_zeta(i,j,knew)- &
1274 & tl_h(i,j)*(1.0_r8-rmask(i,j))
1275#endif
1276 END DO
1277 END DO
1278
1279
1280
1281 IF (predictor_2d_step(ng)) THEN
1282 DO j=jstr,jend
1283 DO i=istr,iend
1284
1285
1286 tl_rzeta(i,j,krhs)=tl_rhs_zeta(i,j)
1287 END DO
1288 END DO
1289
1290 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
1291
1292
1293
1294
1295 CALL exchange_r2d_tile (ng, tile, &
1296 & lbi, ubi, lbj, ubj, &
1297 & tl_rzeta(:,:,krhs))
1298 END IF
1299
1300#ifdef DISTRIBUTE
1301
1302
1303
1304
1305
1306
1307 CALL mp_exchange2d (ng, tile, itlm, 1, &
1308 & lbi, ubi, lbj, ubj, &
1309 & nghostpoints, &
1310 & ewperiodic(ng), nsperiodic(ng), &
1311 & tl_rzeta(:,:,krhs))
1312#endif
1313 END IF
1314
1315
1316
1317
1318
1319 IF (lwsrc(ng)) THEN
1320 DO is=1,nsrc(ng)
1321 IF (int(sources(ng)%Dsrc(is)).eq.2) THEN
1322 i=sources(ng)%Isrc(is)
1323 j=sources(ng)%Jsrc(is)
1324 IF (((istrr.le.i).and.(i.le.iendr)).and. &
1325 & ((jstrr.le.j).and.(j.le.jendr))) THEN
1326
1327
1328
1329
1330
1331 END IF
1332 END IF
1333 END DO
1334 END IF
1335
1336#if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
1337
1338
1339
1340
1341
1342 fac=0.5_r8*dtfast(ng)*ndtfast(ng)/(nfast(ng)*dt(ng))
1343 DO j=jstr,jend
1344 DO i=istr,iend
1345
1346
1347 tl_cff=fac*(tl_bed_thick(i,j,nstp)-tl_bed_thick(i,j,nnew))
1348
1349
1350 tl_h(i,j)=tl_h(i,j)-tl_cff
1351 END DO
1352 END DO
1353#endif
1354
1355
1356
1357
1358
1359
1360
1361
1362
1363 CALL tl_zetabc_tile (ng, tile, &
1364 & lbi, ubi, lbj, ubj, &
1365 & imins, imaxs, jmins, jmaxs, &
1366 & krhs, kstp, knew, &
1367 & zeta, tl_zeta)
1368
1369 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
1370
1371
1372
1373
1374 CALL exchange_r2d_tile (ng, tile, &
1375 & lbi, ubi, lbj, ubj, &
1376 & tl_zeta(:,:,knew))
1377 END IF
1378
1379#ifdef DISTRIBUTE
1380
1381
1382
1383
1384
1385
1386 CALL mp_exchange2d (ng, tile, itlm, 1, &
1387 & lbi, ubi, lbj, ubj, &
1388 & nghostpoints, &
1389 & ewperiodic(ng), nsperiodic(ng), &
1390 & tl_zeta(:,:,knew))
1391#endif
1392
1393
1394
1395
1396
1397
1398
1399
1400
1401 cff1=0.5_r8*g
1402 cff2=1.0_r8/3.0_r8
1403#if !defined SOLVE3D && defined ATM_PRESS
1404 fac3=0.5_r8*100.0_r8/rho0
1405#endif
1406 DO j=jstr,jend
1407 DO i=istru,iend
1408
1409
1410
1411
1412
1413#if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
1414
1415
1416
1417
1418
1419
1420
1421
1422#endif
1423
1424
1425
1426 tl_rhs_ubar(i,j)=cff1*on_u(i,j)* &
1427 & ((tl_h(i-1,j)+ &
1428 & tl_h(i ,j))* &
1429 & (gzeta(i-1,j)- &
1430 & gzeta(i ,j))+ &
1431 & (h(i-1,j)+ &
1432 & h(i ,j))* &
1433 & (tl_gzeta(i-1,j)- &
1434 & tl_gzeta(i ,j))+ &
1435#if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
1436 & (tl_h(i-1,j)- &
1437 & tl_h(i ,j))* &
1438 & (gzetasa(i-1,j)+ &
1439 & gzetasa(i ,j)+ &
1440 & cff2*(rhoa(i-1,j)- &
1441 & rhoa(i ,j))* &
1442 & (zwrk(i-1,j)- &
1443 & zwrk(i ,j)))+ &
1444 & (h(i-1,j)- &
1445 & h(i ,j))* &
1446 & (tl_gzetasa(i-1,j)+ &
1447 & tl_gzetasa(i ,j)+ &
1448 & cff2*((tl_rhoa(i-1,j)- &
1449 & tl_rhoa(i ,j))* &
1450 & (zwrk(i-1,j)- &
1451 & zwrk(i ,j))+ &
1452 & (rhoa(i-1,j)- &
1453 & rhoa(i ,j))* &
1454 & (tl_zwrk(i-1,j)- &
1455 & tl_zwrk(i ,j))))+ &
1456#endif
1457 & (tl_gzeta2(i-1,j)- &
1458 & tl_gzeta2(i ,j)))
1459#if defined ATM_PRESS && !defined SOLVE3D
1460
1461
1462
1463
1464
1465
1466 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)- &
1467 & fac3*on_u(i,j)* &
1468 & (tl_h(i-1,j)+tl_h(i,j)+ &
1469 & tl_gzeta(i-1,j)+tl_gzeta(i,j))* &
1470 & (pair(i,j)-pair(i-1,j))
1471#endif
1472#if defined TIDE_GENERATING_FORCES && !defined SOLVE3D
1473
1474
1475
1476
1477
1478
1479 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)- &
1480 & cff1*on_u(i,j)* &
1481 & ((tl_h(i-1,j)+tl_h(i,j)+ &
1482 & tl_gzeta(i-1,j)+tl_gzeta(i,j))* &
1483 & (eq_tide(i,j)-eq_tide(i-1,j))+ &
1484 & (h(i-1,j)+h(i,j)+ &
1485 & gzeta(i-1,j)+gzeta(i,j))* &
1486 & (tl_eq_tide(i,j)-tl_eq_tide(i-1,j)))
1487#endif
1488#ifdef DIAGNOSTICS_UV
1489
1490#endif
1491 END DO
1492 IF (j.ge.jstrv) THEN
1493 DO i=istr,iend
1494
1495
1496
1497
1498
1499#if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
1500
1501
1502
1503
1504
1505
1506
1507
1508#endif
1509
1510
1511
1512 tl_rhs_vbar(i,j)=cff1*om_v(i,j)* &
1513 & ((tl_h(i,j-1)+ &
1514 & tl_h(i,j ))* &
1515 & (gzeta(i,j-1)- &
1516 & gzeta(i,j ))+ &
1517 & (h(i,j-1)+ &
1518 & h(i,j ))* &
1519 & (tl_gzeta(i,j-1)- &
1520 & tl_gzeta(i,j ))+ &
1521#if defined VAR_RHO_2D_NOT_YET && defined SOLVE3D
1522 & (tl_h(i,j-1)- &
1523 & tl_h(i,j ))* &
1524 & (gzetasa(i,j-1)+ &
1525 & gzetasa(i,j )+ &
1526 & cff2*(rhoa(i,j-1)- &
1527 & rhoa(i,j ))* &
1528 & (zwrk(i,j-1)- &
1529 & zwrk(i,j )))+ &
1530 & (h(i,j-1)- &
1531 & h(i,j ))* &
1532 & (tl_gzetasa(i,j-1)+ &
1533 & tl_gzetasa(i,j )+ &
1534 & cff2*((tl_rhoa(i,j-1)- &
1535 & tl_rhoa(i,j ))* &
1536 & (zwrk(i,j-1)- &
1537 & zwrk(i,j ))+ &
1538 & (rhoa(i,j-1)- &
1539 & rhoa(i,j ))* &
1540 & (tl_zwrk(i,j-1)- &
1541 & tl_zwrk(i,j ))))+ &
1542#endif
1543 & (tl_gzeta2(i,j-1)- &
1544 & tl_gzeta2(i,j )))
1545#if defined ATM_PRESS && !defined SOLVE3D
1546
1547
1548
1549
1550
1551
1552 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)- &
1553 & fac3*om_v(i,j)* &
1554 & (tl_h(i,j-1)+tl_h(i,j)+ &
1555 & tl_gzeta(i,j-1)+tl_gzeta(i,j))* &
1556 & (pair(i,j)-pair(i,j-1))
1557#endif
1558#if defined TIDE_GENERATING_FORCES && !defined SOLVE3D
1559
1560
1561
1562
1563
1564
1565 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)- &
1566 & cff1*om_v(i,j)* &
1567 & ((tl_h(i,j-1)+tl_h(i,j)+ &
1568 & tl_gzeta(i,j-1)+tl_gzeta(i,j))* &
1569 & (eq_tide(i,j)-eq_tide(i,j-1))+ &
1570 & (h(i,j-1)+h(i,j)+ &
1571 & gzeta(i,j-1)+gzeta(i,j))* &
1572 & (tl_eq_tide(i,j)-tl_eq_tide(i,j-1)))
1573#endif
1574#ifdef DIAGNOSTICS_UV
1575
1576#endif
1577 END DO
1578 END IF
1579 END DO
1580#ifdef UV_ADV
1581
1582
1583
1584
1585
1586# ifdef UV_C2ADVECTION
1587
1588
1589
1590 DO j=jstr,jend
1591 DO i=istru-1,iend
1592
1593
1594# ifdef WEC_MELLOR
1595
1596
1597# endif
1598
1599
1600 tl_ufx(i,j)=0.25_r8* &
1601 & ((tl_duon(i,j)+tl_duon(i+1,j))* &
1602 & (ubar(i ,j,krhs)+ &
1603# ifdef WEC_MELLOR
1604 & ubar_stokes(i ,j)+ &
1605 & ubar_stokes(i+1,j)+ &
1606# endif
1607 & ubar(i+1,j,krhs))+ &
1608 & (duon(i,j)+duon(i+1,j))* &
1609 & (tl_ubar(i ,j,krhs)+ &
1610# ifdef WEC_MELLOR
1611 & tl_ubar_stokes(i ,j)+ &
1612 & tl_ubar_stokes(i+1,j)+ &
1613# endif
1614 & tl_ubar(i+1,j,krhs)))
1615 END DO
1616 END DO
1617
1618 DO j=jstr,jend+1
1619 DO i=istru,iend
1620
1621
1622# ifdef WEC_MELLOR
1623
1624
1625# endif
1626
1627
1628 tl_ufe(i,j)=0.25_r8* &
1629 & ((tl_dvom(i,j)+tl_dvom(i-1,j))* &
1630 & (ubar(i,j ,krhs)+ &
1631# ifdef WEC_MELLOR
1632 & ubar_stokes(i,j )+ &
1633 & ubar_stokes(i,j-1)+ &
1634# endif
1635 & ubar(i,j-1,krhs))+ &
1636 & (dvom(i,j)+dvom(i-1,j))* &
1637 & (tl_ubar(i,j ,krhs)+ &
1638# ifdef WEC_MELLOR
1639 & tl_ubar_stokes(i,j )+ &
1640 & tl_ubar_stokes(i,j-1)+ &
1641# endif
1642 & tl_ubar(i,j-1,krhs)))
1643 END DO
1644 END DO
1645
1646 DO j=jstrv,jend
1647 DO i=istr,iend+1
1648
1649
1650# ifdef WEC_MELLOR
1651
1652
1653# endif
1654
1655
1656 tl_vfx(i,j)=0.25_r8* &
1657 & ((tl_duon(i,j)+tl_duon(i,j-1))* &
1658 & (vbar(i ,j,krhs)+ &
1659# ifdef WEC_MELLOR
1660 & vbar_stokes(i ,j)+ &
1661 & vbar_stokes(i-1,j)+ &
1662# endif
1663 & vbar(i-1,j,krhs))+ &
1664 & (duon(i,j)+duon(i,j-1))* &
1665 & (tl_vbar(i ,j,krhs)+ &
1666# ifdef WEC_MELLOR
1667 & tl_vbar_stokes(i ,j)+ &
1668 & tl_vbar_stokes(i-1,j)+ &
1669# endif
1670 & tl_vbar(i-1,j,krhs)))
1671 END DO
1672 END DO
1673
1674 DO j=jstrv-1,jend
1675 DO i=istr,iend
1676
1677
1678# ifdef WEC_MELLOR
1679
1680
1681# endif
1682
1683
1684 tl_vfe(i,j)=0.25_r8* &
1685 & ((tl_dvom(i,j)+tl_dvom(i,j+1))* &
1686 & (vbar(i,j ,krhs)+ &
1687# ifdef WEC_MELLOR
1688 & vbar_stokes(i,j )+ &
1689 & vbar_stokes(i,j+1)+ &
1690# endif
1691 & vbar(i,j+1,krhs))+ &
1692 & (dvom(i,j)+dvom(i,j+1))* &
1693 & (tl_vbar(i,j ,krhs)+ &
1694# ifdef WEC_MELLOR
1695 & tl_vbar_stokes(i,j )+ &
1696 & tl_vbar_stokes(i,j+1)+ &
1697# endif
1698 & tl_vbar(i,j+1,krhs)))
1699 END DO
1700 END DO
1701# else
1702
1703
1704
1705 DO j=jstr,jend
1706 DO i=istrum1,iendp1
1707 grad(i,j)=ubar(i-1,j,krhs)-2.0_r8*ubar(i,j,krhs)+ &
1708# ifdef WEC_MELLOR
1709 & ubar_stokes(i-1,j)-2.0_r8*ubar_stokes(i,j)+ &
1710 & ubar_stokes(i+1,j)+ &
1711# endif
1712 & ubar(i+1,j,krhs)
1713 tl_grad(i,j)=tl_ubar(i-1,j,krhs)-2.0_r8*tl_ubar(i,j,krhs)+ &
1714# ifdef WEC_MELLOR
1715 & tl_ubar_stokes(i-1,j)-2.0_r8*tl_ubar_stokes(i,j)+&
1716 & tl_ubar_stokes(i+1,j)+ &
1717# endif
1718 & tl_ubar(i+1,j,krhs)
1719 dgrad(i,j)=duon(i-1,j)-2.0_r8*duon(i,j)+duon(i+1,j)
1720 tl_dgrad(i,j)=tl_duon(i-1,j)-2.0_r8*tl_duon(i,j)+ &
1721 & tl_duon(i+1,j)
1722 END DO
1723 END DO
1724 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
1725 IF (domain(ng)%Western_Edge(tile)) THEN
1726 DO j=jstr,jend
1727 grad(istr,j)=grad(istr+1,j)
1728 tl_grad(istr,j)=tl_grad(istr+1,j)
1729 dgrad(istr,j)=dgrad(istr+1,j)
1730 tl_dgrad(istr,j)=tl_dgrad(istr+1,j)
1731 END DO
1732 END IF
1733 END IF
1734 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
1735 IF (domain(ng)%Eastern_Edge(tile)) THEN
1736 DO j=jstr,jend
1737 grad(iend+1,j)=grad(iend,j)
1738 tl_grad(iend+1,j)=tl_grad(iend,j)
1739 dgrad(iend+1,j)=dgrad(iend,j)
1740 tl_dgrad(iend+1,j)=tl_dgrad(iend,j)
1741 END DO
1742 END IF
1743 END IF
1744
1745 cff=1.0_r8/6.0_r8
1746 DO j=jstr,jend
1747 DO i=istru-1,iend
1748
1749# ifdef WEC_MELLOR
1750
1751
1752# endif
1753
1754
1755
1756
1757
1758 tl_ufx(i,j)=0.25_r8* &
1759 & ((ubar(i ,j,krhs)+ &
1760# ifdef WEC_MELLOR
1761 & ubar_stokes(i ,j)+ &
1762 & ubar_stokes(i+1,j)+ &
1763# endif
1764 & ubar(i+1,j,krhs)- &
1765 & cff*(grad(i,j)+grad(i+1,j)))* &
1766 & (tl_duon(i,j)+tl_duon(i+1,j)- &
1767 & cff*(tl_dgrad(i,j)+tl_dgrad(i+1,j)))+ &
1768 & (tl_ubar(i ,j,krhs)+ &
1769# ifdef WEC_MELLOR
1770 & tl_ubar_stokes(i ,j)+ &
1771 & tl_ubar_stokes(i+1,j)+ &
1772# endif
1773 & tl_ubar(i+1,j,krhs)- &
1774 & cff*(tl_grad(i,j)+tl_grad(i+1,j)))* &
1775 & (duon(i,j)+duon(i+1,j)- &
1776 & cff*(dgrad(i,j)+dgrad(i+1,j))))
1777 END DO
1778 END DO
1779
1780 DO j=jstrm1,jendp1
1781 DO i=istru,iend
1782 grad(i,j)=ubar(i,j-1,krhs)-2.0_r8*ubar(i,j,krhs)+ &
1783# ifdef WEC_MELLOR
1784 & ubar_stokes(i,j-1)-2.0_r8*ubar_stokes(i,j)+ &
1785 & ubar_stokes(i,j+1)+ &
1786# endif
1787 & ubar(i,j+1,krhs)
1788 tl_grad(i,j)=tl_ubar(i,j-1,krhs)-2.0_r8*tl_ubar(i,j,krhs)+ &
1789# ifdef WEC_MELLOR
1790 & tl_ubar_stokes(i,j-1)-2.0_r8*tl_ubar_stokes(i,j)+&
1791 & tl_ubar_stokes(i,j+1)+ &
1792# endif
1793 & tl_ubar(i,j+1,krhs)
1794 END DO
1795 END DO
1796 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
1797 IF (domain(ng)%Southern_Edge(tile)) THEN
1798 DO i=istru,iend
1799 grad(i,jstr-1)=grad(i,jstr)
1800 tl_grad(i,jstr-1)=tl_grad(i,jstr)
1801 END DO
1802 END IF
1803 END IF
1804 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
1805 IF (domain(ng)%Northern_Edge(tile)) THEN
1806 DO i=istru,iend
1807 grad(i,jend+1)=grad(i,jend)
1808 tl_grad(i,jend+1)=tl_grad(i,jend)
1809 END DO
1810 END IF
1811 END IF
1812 DO j=jstr,jend+1
1813 DO i=istru-1,iend
1814 dgrad(i,j)=dvom(i-1,j)-2.0_r8*dvom(i,j)+dvom(i+1,j)
1815 tl_dgrad(i,j)=tl_dvom(i-1,j)-2.0_r8*tl_dvom(i,j)+ &
1816 & tl_dvom(i+1,j)
1817 END DO
1818 END DO
1819
1820 cff=1.0_r8/6.0_r8
1821 DO j=jstr,jend+1
1822 DO i=istru,iend
1823
1824# ifdef WEC_MELLOR
1825
1826
1827# endif
1828
1829
1830
1831
1832
1833 tl_ufe(i,j)=0.25_r8* &
1834 & ((tl_ubar(i,j ,krhs)+ &
1835# ifdef WEC_MELLOR
1836 & tl_ubar_stokes(i,j )+ &
1837 & tl_ubar_stokes(i,j-1)+ &
1838# endif
1839 & tl_ubar(i,j-1,krhs)- &
1840 & cff*(tl_grad(i,j)+tl_grad(i,j-1)))* &
1841 & (dvom(i,j)+dvom(i-1,j)- &
1842 & cff*(dgrad(i,j)+dgrad(i-1,j)))+ &
1843 & (ubar(i,j ,krhs)+ &
1844# ifdef WEC_MELLOR
1845 & ubar_stokes(i,j )+ &
1846 & ubar_stokes(i,j-1)+ &
1847# endif
1848 & ubar(i,j-1,krhs)- &
1849 & cff*(grad(i,j)+grad(i,j-1)))* &
1850 & (tl_dvom(i,j)+tl_dvom(i-1,j)- &
1851 & cff*(tl_dgrad(i,j)+tl_dgrad(i-1,j))))
1852 END DO
1853 END DO
1854
1855 DO j=jstrv,jend
1856 DO i=istrm1,iendp1
1857 grad(i,j)=vbar(i-1,j,krhs)-2.0_r8*vbar(i,j,krhs)+ &
1858# ifdef WEC_MELLOR
1859 & vbar_stokes(i-1,j)-2.0_r8*vbar_stokes(i,j)+ &
1860 & vbar_stokes(i+1,j)+ &
1861# endif
1862 & vbar(i+1,j,krhs)
1863 tl_grad(i,j)=tl_vbar(i-1,j,krhs)-2.0_r8*tl_vbar(i,j,krhs)+ &
1864# ifdef WEC_MELLOR
1865 & tl_vbar_stokes(i-1,j)-2.0_r8*tl_vbar_stokes(i,j)+&
1866 & tl_vbar_stokes(i+1,j)+ &
1867# endif
1868 & tl_vbar(i+1,j,krhs)
1869 END DO
1870 END DO
1871 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
1872 IF (domain(ng)%Western_Edge(tile)) THEN
1873 DO j=jstrv,jend
1874 grad(istr-1,j)=grad(istr,j)
1875 tl_grad(istr-1,j)=tl_grad(istr,j)
1876 END DO
1877 END IF
1878 END IF
1879 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
1880 IF (domain(ng)%Eastern_Edge(tile)) THEN
1881 DO j=jstrv,jend
1882 grad(iend+1,j)=grad(iend,j)
1883 tl_grad(iend+1,j)=tl_grad(iend,j)
1884 END DO
1885 END IF
1886 END IF
1887 DO j=jstrv-1,jend
1888 DO i=istr,iend+1
1889 dgrad(i,j)=duon(i,j-1)-2.0_r8*duon(i,j)+duon(i,j+1)
1890 tl_dgrad(i,j)=tl_duon(i,j-1)-2.0_r8*tl_duon(i,j)+ &
1891 & tl_duon(i,j+1)
1892 END DO
1893 END DO
1894
1895 cff=1.0_r8/6.0_r8
1896 DO j=jstrv,jend
1897 DO i=istr,iend+1
1898
1899# ifdef WEC_MELLOR
1900
1901
1902# endif
1903
1904
1905
1906
1907
1908 tl_vfx(i,j)=0.25_r8* &
1909 & ((tl_vbar(i ,j,krhs)+ &
1910# ifdef WEC_MELLOR
1911 & tl_vbar_stokes(i ,j)+ &
1912 & tl_vbar_stokes(i-1,j)+ &
1913# endif
1914 & tl_vbar(i-1,j,krhs)- &
1915 & cff*(tl_grad(i,j)+tl_grad(i-1,j)))* &
1916 & (duon(i,j)+duon(i,j-1)- &
1917 & cff*(dgrad(i,j)+dgrad(i,j-1)))+ &
1918 & (vbar(i ,j,krhs)+ &
1919# ifdef WEC_MELLOR
1920 & vbar_stokes(i ,j)+ &
1921 & vbar_stokes(i-1,j)+ &
1922# endif
1923 & vbar(i-1,j,krhs)- &
1924 & cff*(grad(i,j)+grad(i-1,j)))* &
1925 & (tl_duon(i,j)+tl_duon(i,j-1)- &
1926 & cff*(tl_dgrad(i,j)+tl_dgrad(i,j-1))))
1927 END DO
1928 END DO
1929
1930 DO j=jstrvm1,jendp1
1931 DO i=istr,iend
1932 grad(i,j)=vbar(i,j-1,krhs)-2.0_r8*vbar(i,j,krhs)+ &
1933# ifdef WEC_MELLOR
1934 & vbar_stokes(i,j-1)-2.0_r8*vbar_stokes(i,j)+ &
1935 & vbar_stokes(i,j+1)+ &
1936# endif
1937 & vbar(i,j+1,krhs)
1938 tl_grad(i,j)=tl_vbar(i,j-1,krhs)-2.0_r8*tl_vbar(i,j,krhs)+ &
1939# ifdef WEC_MELLOR
1940 & tl_vbar_stokes(i,j-1)-2.0_r8*tl_vbar_stokes(i,j)+&
1941 & tl_vbar_stokes(i,j+1)+ &
1942# endif
1943 & tl_vbar(i,j+1,krhs)
1944 dgrad(i,j)=dvom(i,j-1)-2.0_r8*dvom(i,j)+dvom(i,j+1)
1945 tl_dgrad(i,j)=tl_dvom(i,j-1)-2.0_r8*tl_dvom(i,j)+ &
1946 & tl_dvom(i,j+1)
1947 END DO
1948 END DO
1949 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
1950 IF (domain(ng)%Southern_Edge(tile)) THEN
1951 DO i=istr,iend
1952 grad(i,jstr)=grad(i,jstr+1)
1953 tl_grad(i,jstr)=tl_grad(i,jstr+1)
1954 dgrad(i,jstr)=dgrad(i,jstr+1)
1955 tl_dgrad(i,jstr)=tl_dgrad(i,jstr+1)
1956 END DO
1957 END IF
1958 END IF
1959 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
1960 IF (domain(ng)%Northern_Edge(tile)) THEN
1961 DO i=istr,iend
1962 grad(i,jend+1)=grad(i,jend)
1963 tl_grad(i,jend+1)=tl_grad(i,jend)
1964 dgrad(i,jend+1)=dgrad(i,jend)
1965 tl_dgrad(i,jend+1)=tl_dgrad(i,jend)
1966 END DO
1967 END IF
1968 END IF
1969
1970 cff=1.0_r8/6.0_r8
1971 DO j=jstrv-1,jend
1972 DO i=istr,iend
1973
1974# ifdef WEC_MELLOR
1975
1976
1977# endif
1978
1979
1980
1981
1982
1983 tl_vfe(i,j)=0.25_r8* &
1984 & ((tl_vbar(i,j ,krhs)+ &
1985# ifdef WEC_MELLOR
1986 & tl_vbar_stokes(i,j )+ &
1987 & tl_vbar_stokes(i,j+1)+ &
1988# endif
1989 & tl_vbar(i,j+1,krhs)- &
1990 & cff*(tl_grad(i,j)+tl_grad(i,j+1)))* &
1991 & (dvom(i,j)+dvom(i,j+1)- &
1992 & cff*(dgrad(i,j)+dgrad(i,j+1)))+ &
1993 & (vbar(i,j ,krhs)+ &
1994# ifdef WEC_MELLOR
1995 & vbar_stokes(i,j )+ &
1996 & vbar_stokes(i,j+1)+ &
1997# endif
1998 & vbar(i,j+1,krhs)- &
1999 & cff*(grad(i,j)+grad(i,j+1)))* &
2000 & (tl_dvom(i,j)+tl_dvom(i,j+1)- &
2001 & cff*(tl_dgrad(i,j)+tl_dgrad(i,j+1))))
2002 END DO
2003 END DO
2004# endif
2005
2006 DO j=jstr,jend
2007 DO i=istru,iend
2008
2009
2010 tl_cff1=tl_ufx(i,j)-tl_ufx(i-1,j)
2011
2012
2013 tl_cff2=tl_ufe(i,j+1)-tl_ufe(i,j)
2014
2015
2016 tl_fac=tl_cff1+tl_cff2
2017
2018
2019 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)-tl_fac
2020# if defined DIAGNOSTICS_UV
2021
2022
2023
2024# endif
2025 END DO
2026 END DO
2027 DO j=jstrv,jend
2028 DO i=istr,iend
2029
2030
2031 tl_cff1=tl_vfx(i+1,j)-tl_vfx(i,j)
2032
2033
2034 tl_cff2=tl_vfe(i,j)-tl_vfe(i,j-1)
2035
2036
2037 tl_fac=tl_cff1+tl_cff2
2038
2039
2040 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)-tl_fac
2041# if defined DIAGNOSTICS_UV
2042
2043
2044
2045# endif
2046 END DO
2047 END DO
2048#endif
2049#ifdef UV_COR
2050
2051
2052
2053
2054
2055 DO j=jstrv-1,jend
2056 DO i=istru-1,iend
2057 cff=0.5_r8*drhs(i,j)*fomn(i,j)
2058 tl_cff=0.5_r8*tl_drhs(i,j)*fomn(i,j)
2059
2060# ifdef WEC_MELLOR
2061
2062
2063# endif
2064
2065
2066 tl_ufx(i,j)=tl_cff*(vbar(i,j ,krhs)+ &
2067# ifdef WEC_MELLOR
2068 & vbar_stokes(i,j )+ &
2069 & vbar_stokes(i,j+1)+ &
2070# endif
2071 & vbar(i,j+1,krhs))+ &
2072 & cff*(tl_vbar(i,j ,krhs)+ &
2073# ifdef WEC_MELLOR
2074 & tl_vbar_stokes(i,j )+ &
2075 & tl_vbar_stokes(i,j+1)+ &
2076# endif
2077 & tl_vbar(i,j+1,krhs))
2078
2079# ifdef WEC_MELLOR
2080
2081
2082# endif
2083
2084
2085 tl_vfe(i,j)=tl_cff*(ubar(i ,j,krhs)+ &
2086# ifdef WEC_MELLOR
2087 & ubar_stokes(i ,j)+ &
2088 & ubar_stokes(i+1,j)+ &
2089# endif
2090 & ubar(i+1,j,krhs))+ &
2091 & cff*(tl_ubar(i ,j,krhs)+ &
2092# ifdef WEC_MELLOR
2093 & tl_ubar_stokes(i ,j)+ &
2094 & tl_ubar_stokes(i+1,j)+ &
2095# endif
2096 & tl_ubar(i+1,j,krhs))
2097 END DO
2098 END DO
2099 DO j=jstr,jend
2100 DO i=istru,iend
2101
2102
2103 tl_fac1=0.5_r8*(tl_ufx(i,j)+tl_ufx(i-1,j))
2104
2105
2106 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)+tl_fac1
2107# if defined DIAGNOSTICS_UV
2108
2109# endif
2110 END DO
2111 END DO
2112 DO j=jstrv,jend
2113 DO i=istr,iend
2114
2115
2116 tl_fac1=0.5_r8*(tl_vfe(i,j)+tl_vfe(i,j-1))
2117
2118
2119 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)-tl_fac1
2120# if defined DIAGNOSTICS_UV
2121
2122# endif
2123 END DO
2124 END DO
2125#endif
2126#if defined CURVGRID && defined UV_ADV
2127
2128
2129
2130
2131
2132 DO j=jstrv-1,jend
2133 DO i=istru-1,iend
2134 cff1=0.5_r8*(vbar(i,j ,krhs)+ &
2135# ifdef WEC_MELLOR
2136 & vbar_stokes(i,j )+ &
2137 & vbar_stokes(i,j+1)+ &
2138# endif
2139 & vbar(i,j+1,krhs))
2140 tl_cff1=0.5_r8*(tl_vbar(i,j ,krhs)+ &
2141# ifdef WEC_MELLOR
2142 & tl_vbar_stokes(i,j )+ &
2143 & tl_vbar_stokes(i,j+1)+ &
2144# endif
2145 & tl_vbar(i,j+1,krhs))
2146 cff2=0.5_r8*(ubar(i ,j,krhs)+ &
2147# ifdef WEC_MELLOR
2148 & ubar_stokes(i ,j)+ &
2149 & ubar_stokes(i+1,j)+ &
2150# endif
2151 & ubar(i+1,j,krhs))
2152 tl_cff2=0.5_r8*(tl_ubar(i ,j,krhs)+ &
2153# ifdef WEC_MELLOR
2154 & tl_ubar_stokes(i ,j)+ &
2155 & tl_ubar_stokes(i+1,j)+ &
2156# endif
2157 & tl_ubar(i+1,j,krhs))
2158 cff3=cff1*dndx(i,j)
2159 tl_cff3=tl_cff1*dndx(i,j)
2160 cff4=cff2*dmde(i,j)
2161 tl_cff4=tl_cff2*dmde(i,j)
2162 cff=drhs(i,j)*(cff3-cff4)
2163 tl_cff=tl_drhs(i,j)*(cff3-cff4)+ &
2164 & drhs(i,j)*(tl_cff3-tl_cff4)
2165
2166
2167 tl_ufx(i,j)=tl_cff*cff1+cff*tl_cff1
2168
2169
2170 tl_vfe(i,j)=tl_cff*cff2+cff*tl_cff2
2171# if defined DIAGNOSTICS_UV
2172
2173
2174
2175# endif
2176 END DO
2177 END DO
2178 DO j=jstr,jend
2179 DO i=istru,iend
2180
2181
2182 tl_fac1=0.5_r8*(tl_ufx(i,j)+tl_ufx(i-1,j))
2183
2184
2185 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)+tl_fac1
2186# if defined DIAGNOSTICS_UV
2187
2188
2189
2190
2191# endif
2192 END DO
2193 END DO
2194 DO j=jstrv,jend
2195 DO i=istr,iend
2196
2197
2198 tl_fac1=0.5_r8*(tl_vfe(i,j)+tl_vfe(i,j-1))
2199
2200
2201 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)-tl_fac1
2202# if defined DIAGNOSTICS_UV
2203
2204
2205
2206
2207# endif
2208 END DO
2209 END DO
2210#endif
2211#if defined UV_VIS2 || defined UV_VIS4
2212
2213
2214
2215
2216
2217# ifdef UV_VIS4
2218 DO j=jstrm1,jendp2
2219 DO i=istrm1,iendp2
2220# else
2221 DO j=jstr,jend+1
2222 DO i=istr,iend+1
2223# endif
2224 drhs_p(i,j)=0.25_r8*(drhs(i,j )+drhs(i-1,j )+ &
2225 & drhs(i,j-1)+drhs(i-1,j-1))
2226 tl_drhs_p(i,j)=0.25_r8*(tl_drhs(i,j )+tl_drhs(i-1,j )+ &
2227 & tl_drhs(i,j-1)+tl_drhs(i-1,j-1))
2228 END DO
2229 END DO
2230#endif
2231#ifdef UV_VIS2
2232
2233
2234
2235
2236
2237
2238
2239
2240 DO j=jstrv-1,jend
2241 DO i=istru-1,iend
2242
2243
2244
2245
2246
2247
2248
2249
2250 tl_cff=visc2_r(i,j)*0.5_r8* &
2251 & (tl_drhs(i,j)* &
2252 & (pmon_r(i,j)* &
2253 & ((pn(i ,j)+pn(i+1,j))*ubar(i+1,j,krhs)- &
2254 & (pn(i-1,j)+pn(i ,j))*ubar(i ,j,krhs))- &
2255 & pnom_r(i,j)* &
2256 & ((pm(i,j )+pm(i,j+1))*vbar(i,j+1,krhs)- &
2257 & (pm(i,j-1)+pm(i,j ))*vbar(i,j ,krhs)))+ &
2258 & drhs(i,j)* &
2259 & (pmon_r(i,j)* &
2260 & ((pn(i ,j)+pn(i+1,j))*tl_ubar(i+1,j,krhs)- &
2261 & (pn(i-1,j)+pn(i ,j))*tl_ubar(i ,j,krhs))- &
2262 & pnom_r(i,j)* &
2263 & ((pm(i,j )+pm(i,j+1))*tl_vbar(i,j+1,krhs)- &
2264 & (pm(i,j-1)+pm(i,j ))*tl_vbar(i,j ,krhs))))
2265
2266
2267 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*tl_cff
2268
2269
2270 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*tl_cff
2271 END DO
2272 END DO
2273 DO j=jstr,jend+1
2274 DO i=istr,iend+1
2275
2276
2277
2278
2279
2280
2281
2282
2283 tl_cff=visc2_p(i,j)*0.5_r8* &
2284 & (tl_drhs_p(i,j)* &
2285 & (pmon_p(i,j)* &
2286 & ((pn(i ,j-1)+pn(i ,j))*vbar(i ,j,krhs)- &
2287 & (pn(i-1,j-1)+pn(i-1,j))*vbar(i-1,j,krhs))+ &
2288 & pnom_p(i,j)* &
2289 & ((pm(i-1,j )+pm(i,j ))*ubar(i,j ,krhs)- &
2290 & (pm(i-1,j-1)+pm(i,j-1))*ubar(i,j-1,krhs)))+ &
2291 & drhs_p(i,j)* &
2292 & (pmon_p(i,j)* &
2293 & ((pn(i ,j-1)+pn(i ,j))*tl_vbar(i ,j,krhs)- &
2294 & (pn(i-1,j-1)+pn(i-1,j))*tl_vbar(i-1,j,krhs))+ &
2295 & pnom_p(i,j)* &
2296 & ((pm(i-1,j )+pm(i,j ))*tl_ubar(i,j ,krhs)- &
2297 & (pm(i-1,j-1)+pm(i,j-1))*tl_ubar(i,j-1,krhs))))
2298# ifdef MASKING
2299
2300
2301 tl_cff=tl_cff*pmask(i,j)
2302# endif
2303
2304
2305 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*tl_cff
2306
2307
2308 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*tl_cff
2309 END DO
2310 END DO
2311
2312
2313
2314 DO j=jstr,jend
2315 DO i=istru,iend
2316
2317
2318 tl_cff1=0.5_r8*(pn(i-1,j)+pn(i,j))* &
2319 & (tl_ufx(i,j )-tl_ufx(i-1,j))
2320
2321
2322 tl_cff2=0.5_r8*(pm(i-1,j)+pm(i,j))* &
2323 & (tl_ufe(i,j+1)-tl_ufe(i ,j))
2324
2325
2326 tl_fac=tl_cff1+tl_cff2
2327
2328
2329 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)+tl_fac
2330# if defined DIAGNOSTICS_UV
2331
2332
2333
2334# endif
2335 END DO
2336 END DO
2337 DO j=jstrv,jend
2338 DO i=istr,iend
2339
2340
2341 tl_cff1=0.5_r8*(pn(i,j-1)+pn(i,j))* &
2342 & (tl_vfx(i+1,j)-tl_vfx(i,j ))
2343
2344
2345 tl_cff2=0.5_r8*(pm(i,j-1)+pm(i,j))* &
2346 & (tl_vfe(i ,j)-tl_vfe(i,j-1))
2347
2348
2349 tl_fac=tl_cff1-tl_cff2
2350
2351
2352 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)+tl_fac
2353# if defined DIAGNOSTICS_UV
2354
2355
2356
2357# endif
2358 END DO
2359 END DO
2360#endif
2361#ifdef UV_VIS4
2362
2363
2364
2365
2366
2367
2368
2369
2370
2371
2372
2373
2374
2375 DO j=jstrvm2,jendp1
2376 DO i=istrum2,iendp1
2377 cff=visc4_r(i,j)*0.5_r8* &
2378 & (pmon_r(i,j)* &
2379 & ((pn(i ,j)+pn(i+1,j))*ubar(i+1,j,krhs)- &
2380 & (pn(i-1,j)+pn(i ,j))*ubar(i ,j,krhs))- &
2381 & pnom_r(i,j)* &
2382 & ((pm(i,j )+pm(i,j+1))*vbar(i,j+1,krhs)- &
2383 & (pm(i,j-1)+pm(i,j ))*vbar(i,j ,krhs)))
2384 tl_cff=visc4_r(i,j)*0.5_r8* &
2385 & (pmon_r(i,j)* &
2386 & ((pn(i ,j)+pn(i+1,j))*tl_ubar(i+1,j,krhs)- &
2387 & (pn(i-1,j)+pn(i ,j))*tl_ubar(i ,j,krhs))- &
2388 & pnom_r(i,j)* &
2389 & ((pm(i,j )+pm(i,j+1))*tl_vbar(i,j+1,krhs)- &
2390 & (pm(i,j-1)+pm(i,j ))*tl_vbar(i,j ,krhs)))
2391 ufx(i,j)=on_r(i,j)*on_r(i,j)*cff
2392 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*tl_cff
2393 vfe(i,j)=om_r(i,j)*om_r(i,j)*cff
2394 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*tl_cff
2395 END DO
2396 END DO
2397 DO j=jstrm1,jendp2
2398 DO i=istrm1,iendp2
2399 cff=visc4_p(i,j)*0.5_r8* &
2400 & (pmon_p(i,j)* &
2401 & ((pn(i ,j-1)+pn(i ,j))*vbar(i ,j,krhs)- &
2402 & (pn(i-1,j-1)+pn(i-1,j))*vbar(i-1,j,krhs))+ &
2403 & pnom_p(i,j)* &
2404 & ((pm(i-1,j )+pm(i,j ))*ubar(i,j ,krhs)- &
2405 & (pm(i-1,j-1)+pm(i,j-1))*ubar(i,j-1,krhs)))
2406 tl_cff=visc4_p(i,j)*0.5_r8* &
2407 & (pmon_p(i,j)* &
2408 & ((pn(i ,j-1)+pn(i ,j))*tl_vbar(i ,j,krhs)- &
2409 & (pn(i-1,j-1)+pn(i-1,j))*tl_vbar(i-1,j,krhs))+ &
2410 & pnom_p(i,j)* &
2411 & ((pm(i-1,j )+pm(i,j ))*tl_ubar(i,j ,krhs)- &
2412 & (pm(i-1,j-1)+pm(i,j-1))*tl_ubar(i,j-1,krhs)))
2413# ifdef MASKING
2414 cff=cff*pmask(i,j)
2415 tl_cff=tl_cff*pmask(i,j)
2416# endif
2417 ufe(i,j)=om_p(i,j)*om_p(i,j)*cff
2418 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*tl_cff
2419 vfx(i,j)=on_p(i,j)*on_p(i,j)*cff
2420 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*tl_cff
2421 END DO
2422 END DO
2423
2424
2425
2426 DO j=jstrm1,jendp1
2427 DO i=istrum1,iendp1
2428 lapu(i,j)=0.125_r8* &
2429 & (pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))* &
2430 & ((pn(i-1,j)+pn(i,j))* &
2431 & (ufx(i,j )-ufx(i-1,j))+ &
2432 & (pm(i-1,j)+pm(i,j))* &
2433 & (ufe(i,j+1)-ufe(i ,j)))
2434 tl_lapu(i,j)=0.125_r8* &
2435 & (pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))* &
2436 & ((pn(i-1,j)+pn(i,j))* &
2437 & (tl_ufx(i,j )-tl_ufx(i-1,j))+ &
2438 & (pm(i-1,j)+pm(i,j))* &
2439 & (tl_ufe(i,j+1)-tl_ufe(i ,j)))
2440 END DO
2441 END DO
2442 DO j=jstrvm1,jendp1
2443 DO i=istrm1,iendp1
2444 lapv(i,j)=0.125_r8* &
2445 & (pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))* &
2446 & ((pn(i,j-1)+pn(i,j))* &
2447 & (vfx(i+1,j)-vfx(i,j ))- &
2448 & (pm(i,j-1)+pm(i,j))* &
2449 & (vfe(i ,j)-vfe(i,j-1)))
2450 tl_lapv(i,j)=0.125_r8* &
2451 & (pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))* &
2452 & ((pn(i,j-1)+pn(i,j))* &
2453 & (tl_vfx(i+1,j)-tl_vfx(i,j ))- &
2454 & (pm(i,j-1)+pm(i,j))* &
2455 & (tl_vfe(i ,j)-tl_vfe(i,j-1)))
2456 END DO
2457 END DO
2458
2459
2460
2461
2462
2463 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
2464 IF (domain(ng)%Western_Edge(tile)) THEN
2465 IF (tl_lbc(iwest,isubar,ng)%closed) THEN
2466 DO j=jstrm1,jendp1
2467 lapu(istru-1,j)=0.0_r8
2468 tl_lapu(istru-1,j)=0.0_r8
2469 END DO
2470 ELSE
2471 DO j=jstrm1,jendp1
2472 lapu(istru-1,j)=lapu(istru,j)
2473 tl_lapu(istru-1,j)=tl_lapu(istru,j)
2474 END DO
2475 END IF
2476 IF (tl_lbc(iwest,isvbar,ng)%closed) THEN
2477 DO j=jstrvm1,jendp1
2478 lapv(istr-1,j)=gamma2(ng)*lapv(istr,j)
2479 tl_lapv(istr-1,j)=gamma2(ng)*tl_lapv(istr,j)
2480 END DO
2481 ELSE
2482 DO j=jstrvm1,jendp1
2483 lapv(istr-1,j)=0.0_r8
2484 tl_lapv(istr-1,j)=0.0_r8
2485 END DO
2486 END IF
2487 END IF
2488 END IF
2489
2490 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
2491 IF (domain(ng)%Eastern_Edge(tile)) THEN
2492 IF (tl_lbc(ieast,isubar,ng)%closed) THEN
2493 DO j=jstrm1,jendp1
2494 lapu(iend+1,j)=0.0_r8
2495 tl_lapu(iend+1,j)=0.0_r8
2496 END DO
2497 ELSE
2498 DO j=jstrm1,jendp1
2499 lapu(iend+1,j)=lapu(iend,j)
2500 tl_lapu(iend+1,j)=tl_lapu(iend,j)
2501 END DO
2502 END IF
2503 IF (tl_lbc(ieast,isvbar,ng)%closed) THEN
2504 DO j=jstrvm1,jendp1
2505 lapv(iend+1,j)=gamma2(ng)*lapv(iend,j)
2506 tl_lapv(iend+1,j)=gamma2(ng)*tl_lapv(iend,j)
2507 END DO
2508 ELSE
2509 DO j=jstrvm1,jendp1
2510 lapv(iend+1,j)=0.0_r8
2511 tl_lapv(iend+1,j)=0.0_r8
2512 END DO
2513 END IF
2514 END IF
2515 END IF
2516
2517 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
2518 IF (domain(ng)%Southern_Edge(tile)) THEN
2519 IF (tl_lbc(isouth,isubar,ng)%closed) THEN
2520 DO i=istrum1,iendp1
2521 lapu(i,jstr-1)=gamma2(ng)*lapu(i,jstr)
2522 tl_lapu(i,jstr-1)=gamma2(ng)*tl_lapu(i,jstr)
2523 END DO
2524 ELSE
2525 DO i=istrum1,iendp1
2526 lapu(i,jstr-1)=0.0_r8
2527 tl_lapu(i,jstr-1)=0.0_r8
2528 END DO
2529 END IF
2530 IF (tl_lbc(isouth,isvbar,ng)%closed) THEN
2531 DO i=istrm1,iendp1
2532 lapv(i,jstrv-1)=0.0_r8
2533 tl_lapv(i,jstrv-1)=0.0_r8
2534 END DO
2535 ELSE
2536 DO i=istrm1,iendp1
2537 lapv(i,jstrv-1)=lapv(i,jstrv)
2538 tl_lapv(i,jstrv-1)=tl_lapv(i,jstrv)
2539 END DO
2540 END IF
2541 END IF
2542 END IF
2543
2544 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
2545 IF (domain(ng)%Northern_Edge(tile)) THEN
2546 IF (tl_lbc(inorth,isubar,ng)%closed) THEN
2547 DO i=istrum1,iendp1
2548 lapu(i,jend+1)=gamma2(ng)*lapu(i,jend)
2549 tl_lapu(i,jend+1)=gamma2(ng)*tl_lapu(i,jend)
2550 END DO
2551 ELSE
2552 DO i=istrum1,iendp1
2553 lapu(i,jend+1)=0.0_r8
2554 tl_lapu(i,jend+1)=0.0_r8
2555 END DO
2556 END IF
2557 IF (tl_lbc(inorth,isvbar,ng)%closed) THEN
2558 DO i=istrm1,iendp1
2559 lapv(i,jend+1)=0.0_r8
2560 tl_lapv(i,jend+1)=0.0_r8
2561 END DO
2562 ELSE
2563 DO i=istrm1,iendp1
2564 lapv(i,jend+1)=lapv(i,jend)
2565 tl_lapv(i,jend+1)=tl_lapv(i,jend)
2566 END DO
2567 END IF
2568 END IF
2569 END IF
2570
2571 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
2572 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
2573 IF (domain(ng)%SouthWest_Corner(tile)) THEN
2574 lapu(istr ,jstr-1)=0.5_r8*(lapu(istr+1,jstr-1)+ &
2575 & lapu(istr ,jstr ))
2576 tl_lapu(istr ,jstr-1)=0.5_r8*(tl_lapu(istr+1,jstr-1)+ &
2577 & tl_lapu(istr ,jstr ))
2578 lapv(istr-1,jstr )=0.5_r8*(lapv(istr-1,jstr+1)+ &
2579 & lapv(istr ,jstr ))
2580 tl_lapv(istr-1,jstr )=0.5_r8*(tl_lapv(istr-1,jstr+1)+ &
2581 & tl_lapv(istr ,jstr ))
2582 END IF
2583 END IF
2584
2585 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
2586 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
2587 IF (domain(ng)%SouthEast_Corner(tile)) THEN
2588 lapu(iend+1,jstr-1)=0.5_r8*(lapu(iend ,jstr-1)+ &
2589 & lapu(iend+1,jstr ))
2590 tl_lapu(iend+1,jstr-1)=0.5_r8*(tl_lapu(iend ,jstr-1)+ &
2591 & tl_lapu(iend+1,jstr ))
2592 lapv(iend+1,jstr )=0.5_r8*(lapv(iend ,jstr )+ &
2593 & lapv(iend+1,jstr+1))
2594 tl_lapv(iend+1,jstr )=0.5_r8*(tl_lapv(iend ,jstr )+ &
2595 & tl_lapv(iend+1,jstr+1))
2596 END IF
2597 END IF
2598
2599 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
2600 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
2601 IF (domain(ng)%NorthWest_Corner(tile)) THEN
2602 lapu(istr ,jend+1)=0.5_r8*(lapu(istr+1,jend+1)+ &
2603 & lapu(istr ,jend ))
2604 tl_lapu(istr ,jend+1)=0.5_r8*(tl_lapu(istr+1,jend+1)+ &
2605 & tl_lapu(istr ,jend ))
2606 lapv(istr-1,jend+1)=0.5_r8*(lapv(istr ,jend+1)+ &
2607 & lapv(istr-1,jend ))
2608 tl_lapv(istr-1,jend+1)=0.5_r8*(tl_lapv(istr ,jend+1)+ &
2609 & tl_lapv(istr-1,jend ))
2610 END IF
2611 END IF
2612
2613 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
2614 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
2615 IF (domain(ng)%NorthEast_Corner(tile)) THEN
2616 lapu(iend+1,jend+1)=0.5_r8*(lapu(iend ,jend+1)+ &
2617 & lapu(iend+1,jend ))
2618 tl_lapu(iend+1,jend+1)=0.5_r8*(tl_lapu(iend ,jend+1)+ &
2619 & tl_lapu(iend+1,jend ))
2620 lapv(iend+1,jend+1)=0.5_r8*(lapv(iend ,jend+1)+ &
2621 & lapv(iend+1,jend ))
2622 tl_lapv(iend+1,jend+1)=0.5_r8*(tl_lapv(iend ,jend+1)+ &
2623 & tl_lapv(iend+1,jend ))
2624 END IF
2625 END IF
2626
2627
2628
2629
2630 DO j=jstrv-1,jend
2631 DO i=istru-1,iend
2632
2633
2634
2635
2636
2637
2638
2639
2640 tl_cff=visc4_r(i,j)*0.5_r8* &
2641 & (tl_drhs(i,j)* &
2642 & (pmon_r(i,j)* &
2643 & ((pn(i ,j)+pn(i+1,j))*lapu(i+1,j)- &
2644 & (pn(i-1,j)+pn(i ,j))*lapu(i ,j))- &
2645 & pnom_r(i,j)* &
2646 & ((pm(i,j )+pm(i,j+1))*lapv(i,j+1)- &
2647 & (pm(i,j-1)+pm(i,j ))*lapv(i,j )))+ &
2648 & drhs(i,j)* &
2649 & (pmon_r(i,j)* &
2650 & ((pn(i ,j)+pn(i+1,j))*tl_lapu(i+1,j)- &
2651 & (pn(i-1,j)+pn(i ,j))*tl_lapu(i ,j))- &
2652 & pnom_r(i,j)* &
2653 & ((pm(i,j )+pm(i,j+1))*tl_lapv(i,j+1)- &
2654 & (pm(i,j-1)+pm(i,j ))*tl_lapv(i,j ))))
2655
2656
2657 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*tl_cff
2658
2659
2660 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*tl_cff
2661 END DO
2662 END DO
2663 DO j=jstr,jend+1
2664 DO i=istr,iend+1
2665
2666
2667
2668
2669
2670
2671
2672
2673 tl_cff=visc4_p(i,j)*0.5_r8* &
2674 & (tl_drhs_p(i,j)* &
2675 & (pmon_p(i,j)* &
2676 & ((pn(i ,j-1)+pn(i ,j))*lapv(i ,j)- &
2677 & (pn(i-1,j-1)+pn(i-1,j))*lapv(i-1,j))+ &
2678 & pnom_p(i,j)* &
2679 & ((pm(i-1,j )+pm(i,j ))*lapu(i,j )- &
2680 & (pm(i-1,j-1)+pm(i,j-1))*lapu(i,j-1)))+ &
2681 & drhs_p(i,j)* &
2682 & (pmon_p(i,j)* &
2683 & ((pn(i ,j-1)+pn(i ,j))*tl_lapv(i ,j)- &
2684 & (pn(i-1,j-1)+pn(i-1,j))*tl_lapv(i-1,j))+ &
2685 & pnom_p(i,j)* &
2686 & ((pm(i-1,j )+pm(i,j ))*tl_lapu(i,j )- &
2687 & (pm(i-1,j-1)+pm(i,j-1))*tl_lapu(i,j-1))))
2688# ifdef MASKING
2689
2690
2691 tl_cff=tl_cff*pmask(i,j)
2692# endif
2693
2694
2695 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*tl_cff
2696
2697
2698 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*tl_cff
2699 END DO
2700 END DO
2701
2702
2703
2704 DO j=jstr,jend
2705 DO i=istru,iend
2706
2707
2708 tl_cff1=0.5_r8*(pn(i-1,j)+pn(i,j))* &
2709 & (tl_ufx(i,j )-tl_ufx(i-1,j))
2710
2711
2712 tl_cff2=0.5_r8*(pm(i-1,j)+pm(i,j))* &
2713 & (ufe(i,j+1)-ufe(i ,j))
2714
2715
2716 tl_fac=tl_cff1+tl_cff2
2717
2718
2719 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)+tl_fac
2720# if defined DIAGNOSTICS_UV
2721
2722
2723
2724# endif
2725 END DO
2726 END DO
2727 DO j=jstrv,jend
2728 DO i=istr,iend
2729
2730
2731 tl_cff1=0.5_r8*(pn(i,j-1)+pn(i,j))* &
2732 & (tl_vfx(i+1,j)-tl_vfx(i,j ))
2733
2734
2735 tl_cff2=0.5_r8*(pm(i,j-1)+pm(i,j))* &
2736 & (tl_vfe(i ,j)-tl_vfe(i,j-1))
2737
2738
2739 tl_fac=tl_cff1-tl_cff2
2740
2741
2742 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)+tl_fac
2743# if defined DIAGNOSTICS_UV
2744
2745
2746
2747# endif
2748 END DO
2749 END DO
2750#endif
2751#if defined WEC_MELLOR && \
2752 (
2753
2754
2755
2756
2757
2758 DO j=jstr,jend
2759 DO i=istru,iend
2760
2761
2762 tl_cff1=tl_rustr2d(i,j)*om_u(i,j)*on_u(i,j)
2763
2764
2765 tl_cff2=tl_rulag2d(i,j)
2766# ifndef SOLVE3D
2767
2768
2769 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)-tl_cff1-tl_cff2
2770# endif
2771# ifdef DIAGNOSTICS_UV
2772
2773# endif
2774 END DO
2775 END DO
2776 DO j=jstrv,jend
2777 DO i=istr,iend
2778
2779
2780 tl_cff1=tl_rvstr2d(i,j)*om_v(i,j)*on_v(i,j)
2781
2782
2783 tl_cff2=tl_rvlag2d(i,j)
2784# ifndef SOLVE3D
2785
2786
2787 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)-tl_cff1-tl_cff2
2788# endif
2789# ifdef DIAGNOSTICS_UV
2790
2791# endif
2792 END DO
2793 END DO
2794#endif
2795#ifndef SOLVE3D
2796
2797
2798
2799
2800
2801 DO j=jstr,jend
2802 DO i=istru,iend
2803
2804
2805 tl_fac=tl_bustr(i,j)*om_u(i,j)*on_u(i,j)
2806
2807
2808 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)-tl_fac
2809# ifdef DIAGNOSTICS_UV
2810
2811# endif
2812 END DO
2813 END DO
2814 DO j=jstrv,jend
2815 DO i=istr,iend
2816
2817
2818 tl_fac=tl_bvstr(i,j)*om_v(i,j)*on_v(i,j)
2819
2820
2821 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)-tl_fac
2822# ifdef DIAGNOSTICS_UV
2823
2824# endif
2825 END DO
2826 END DO
2827#else
2828# ifdef DIAGNOSTICS_UV
2829
2830
2831
2832
2833
2834
2835
2836
2837
2838
2839
2840
2841
2842# endif
2843#endif
2844
2845
2846
2847
2848
2849 IF (lnudgem2clm(ng)) THEN
2850 DO j=jstr,jend
2851 DO i=istru,iend
2852 cff=0.25_r8*(clima(ng)%M2nudgcof(i-1,j)+ &
2853 & clima(ng)%M2nudgcof(i ,j))* &
2854 & om_u(i,j)*on_u(i,j)
2855
2856
2857
2858
2859
2860 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)+ &
2861 & cff*((drhs(i-1,j)+drhs(i,j))* &
2862 & (-tl_ubar(i,j,krhs))+ &
2863 & (tl_drhs(i-1,j)+tl_drhs(i,j))* &
2864 & (clima(ng)%ubarclm(i,j)- &
2865 & ubar(i,j,krhs)))
2866 END DO
2867 END DO
2868 DO j=jstrv,jend
2869 DO i=istr,iend
2870 cff=0.25_r8*(clima(ng)%M2nudgcof(i,j-1)+ &
2871 & clima(ng)%M2nudgcof(i,j ))* &
2872 & om_v(i,j)*on_v(i,j)
2873
2874
2875
2876
2877
2878 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)+ &
2879 & cff*((drhs(i,j-1)+drhs(i,j))* &
2880 & (-tl_vbar(i,j,krhs))+ &
2881 & (tl_drhs(i,j-1)+tl_drhs(i,j))* &
2882 & (clima(ng)%vbarclm(i,j)- &
2883 & vbar(i,j,krhs)))
2884 END DO
2885 END DO
2886 END IF
2887
2888#ifdef SOLVE3D
2889
2890
2891
2892
2893
2894
2895
2896
2897
2898
2899
2900
2901
2902
2903
2904
2905
2906
2907 IF (first_2d_step.and.predictor_2d_step(ng)) THEN
2908#if defined STOCHASTIC_OPT && !defined STOCH_OPT_WHITE
2909 IF (iic(ng).eq.ntfirst(ng).and.soinitial(ng)) THEN
2910#else
2911 IF (iic(ng).eq.ntfirst(ng)) THEN
2912#endif
2913 DO j=jstr,jend
2914 DO i=istru,iend
2915
2916
2917 tl_rufrc(i,j)=tl_rufrc(i,j)-tl_rhs_ubar(i,j)
2918
2919
2920 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)+tl_rufrc(i,j)
2921
2922
2923 tl_ru(i,j,0,nstp)=tl_rufrc(i,j)
2924# ifdef DIAGNOSTICS_UV
2925
2926
2927
2928
2929
2930
2931
2932
2933
2934
2935
2936# endif
2937 END DO
2938 END DO
2939 DO j=jstrv,jend
2940 DO i=istr,iend
2941
2942
2943 tl_rvfrc(i,j)=tl_rvfrc(i,j)-tl_rhs_vbar(i,j)
2944
2945
2946 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)+tl_rvfrc(i,j)
2947
2948
2949 tl_rv(i,j,0,nstp)=tl_rvfrc(i,j)
2950# ifdef DIAGNOSTICS_UV
2951
2952
2953
2954
2955
2956
2957
2958
2959
2960
2961
2962# endif
2963 END DO
2964 END DO
2965#if defined STOCHASTIC_OPT && !defined STOCH_OPT_WHITE
2966 ELSE IF (iic(ng).eq.(ntfirst(ng)+1).and.soinitial(ng)) THEN
2967#else
2968 ELSE IF (iic(ng).eq.(ntfirst(ng)+1)) THEN
2969#endif
2970 DO j=jstr,jend
2971 DO i=istru,iend
2972
2973
2974 tl_rufrc(i,j)=tl_rufrc(i,j)-tl_rhs_ubar(i,j)
2975
2976
2977
2978 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)+ &
2979 & 1.5_r8*tl_rufrc(i,j)- &
2980 & 0.5_r8*tl_ru(i,j,0,nnew)
2981
2982
2983 tl_ru(i,j,0,nstp)=tl_rufrc(i,j)
2984# ifdef DIAGNOSTICS_UV
2985
2986
2987
2988
2989
2990
2991
2992
2993
2994
2995
2996
2997
2998
2999# endif
3000 END DO
3001 END DO
3002 DO j=jstrv,jend
3003 DO i=istr,iend
3004
3005
3006 tl_rvfrc(i,j)=tl_rvfrc(i,j)-tl_rhs_vbar(i,j)
3007
3008
3009
3010 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)+ &
3011 & 1.5_r8*tl_rvfrc(i,j)- &
3012 & 0.5_r8*tl_rv(i,j,0,nnew)
3013
3014
3015 tl_rv(i,j,0,nstp)=tl_rvfrc(i,j)
3016# ifdef DIAGNOSTICS_UV
3017
3018
3019
3020
3021
3022
3023
3024
3025
3026
3027
3028
3029
3030
3031# endif
3032 END DO
3033 END DO
3034 ELSE
3035 cff1=23.0_r8/12.0_r8
3036 cff2=16.0_r8/12.0_r8
3037 cff3= 5.0_r8/12.0_r8
3038 DO j=jstr,jend
3039 DO i=istru,iend
3040
3041
3042 tl_rufrc(i,j)=tl_rufrc(i,j)-tl_rhs_ubar(i,j)
3043
3044
3045
3046
3047
3048 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)+ &
3049 & cff1*tl_rufrc(i,j)- &
3050 & cff2*tl_ru(i,j,0,nnew)+ &
3051 & cff3*tl_ru(i,j,0,nstp)
3052
3053
3054 tl_ru(i,j,0,nstp)=tl_rufrc(i,j)
3055# ifdef DIAGNOSTICS_UV
3056
3057
3058
3059
3060
3061
3062
3063
3064
3065
3066
3067
3068
3069
3070
3071
3072
3073# endif
3074 END DO
3075 END DO
3076 DO j=jstrv,jend
3077 DO i=istr,iend
3078
3079
3080 tl_rvfrc(i,j)=tl_rvfrc(i,j)-tl_rhs_vbar(i,j)
3081
3082
3083
3084
3085
3086 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)+ &
3087 & cff1*tl_rvfrc(i,j)- &
3088 & cff2*tl_rv(i,j,0,nnew)+ &
3089 & cff3*tl_rv(i,j,0,nstp)
3090
3091
3092 tl_rv(i,j,0,nstp)=tl_rvfrc(i,j)
3093# ifdef DIAGNOSTICS_UV
3094
3095
3096
3097
3098
3099
3100
3101
3102
3103
3104
3105
3106
3107
3108
3109
3110
3111# endif
3112 END DO
3113 END DO
3114 END IF
3115 ELSE
3116 DO j=jstr,jend
3117 DO i=istru,iend
3118
3119
3120 tl_rhs_ubar(i,j)=tl_rhs_ubar(i,j)+tl_rufrc(i,j)
3121# ifdef DIAGNOSTICS_UV
3122
3123
3124
3125
3126
3127
3128# endif
3129 END DO
3130 END DO
3131 DO j=jstrv,jend
3132 DO i=istr,iend
3133
3134
3135 tl_rhs_vbar(i,j)=tl_rhs_vbar(i,j)+tl_rvfrc(i,j)
3136# ifdef DIAGNOSTICS_UV
3137
3138
3139
3140
3141
3142
3143# endif
3144 END DO
3145 END DO
3146 END IF
3147#else
3148
3149
3150
3151
3152
3153
3154
3155
3156
3157# ifdef DIAGNOSTICS_UV
3158
3159# endif
3160
3161
3162
3163
3164
3165
3166# ifdef DIAGNOSTICS_UV
3167
3168# endif
3169
3170
3171#endif
3172
3173
3174
3175
3176
3177
3178
3179 DO j=jstrv-1,jend
3180 DO i=istru-1,iend
3181 dstp(i,j)=zeta(i,j,kstp)+h(i,j)
3182 tl_dstp(i,j)=tl_zeta(i,j,kstp)+tl_h(i,j)
3183 END DO
3184 END DO
3185
3186
3187
3188
3189#ifdef WET_DRY_NOT_YET
3190
3191
3192#endif
3193
3194#if defined STOCHASTIC_OPT && !defined STOCH_OPT_WHITE && \
3195
3196 IF (first_2d_step.and.soinitial(ng)) THEN
3197#else
3198 IF (first_2d_step) THEN
3199#endif
3200 cff1=0.5_r8*dtfast(ng)
3201#ifdef WET_DRY_NOT_YET
3202 cff2=1.0_r8/cff1
3203#endif
3204 DO j=jstr,jend
3205 DO i=istru,iend
3206 cff=(pm(i,j)+pm(i-1,j))*(pn(i,j)+pn(i-1,j))
3207 fac=1.0_r8/(dnew(i,j)+dnew(i-1,j))
3208 tl_fac=-fac*fac*(tl_dnew(i,j)+tl_dnew(i-1,j))
3209
3210
3211
3212
3213 tl_ubar(i,j,knew)=(tl_ubar(i,j,kstp)* &
3214 & (dstp(i,j)+dstp(i-1,j))+ &
3215 & ubar(i,j,kstp)* &
3216 & (tl_dstp(i,j)+tl_dstp(i-1,j))+ &
3217 & cff*cff1*tl_rhs_ubar(i,j))*fac+ &
3218 & (ubar(i,j,kstp)* &
3219 & (dstp(i,j)+dstp(i-1,j))+ &
3220 & cff*cff1*rhs_ubar(i,j))*tl_fac
3221#ifdef MASKING
3222
3223
3224 tl_ubar(i,j,knew)=tl_ubar(i,j,knew)*umask(i,j)
3225#endif
3226#ifdef WET_DRY_NOT_YET
3227
3228
3229
3230
3231
3232
3233
3234 fac1=cff2/cff
3235
3236
3237
3238
3239 tl_rhs_ubar(i,j)=(tl_ubar(i,j,knew)* &
3240 & (dnew(i,j)+dnew(i-1,j))+ &
3241 & ubar(i,j,knew)* &
3242 & (tl_dnew(i,j)+tl_dnew(i-1,j))- &
3243 & tl_ubar(i,j,kstp)* &
3244 & (dstp(i,j)+dstp(i-1,j))- &
3245 & ubar(i,j,kstp)* &
3246 & (tl_dstp(i,j)+tl_dstp(i-1,j)))*fac1
3247#endif
3248 END DO
3249 END DO
3250 DO j=jstrv,jend
3251 DO i=istr,iend
3252 cff=(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
3253 fac=1.0_r8/(dnew(i,j)+dnew(i,j-1))
3254 tl_fac=-fac*fac*(tl_dnew(i,j)+tl_dnew(i,j-1))
3255
3256
3257
3258
3259 tl_vbar(i,j,knew)=(tl_vbar(i,j,kstp)* &
3260 & (dstp(i,j)+dstp(i,j-1))+ &
3261 & vbar(i,j,kstp)* &
3262 & (tl_dstp(i,j)+tl_dstp(i,j-1))+ &
3263 & cff*cff1*tl_rhs_vbar(i,j))*fac+ &
3264 & (vbar(i,j,kstp)* &
3265 & (dstp(i,j)+dstp(i,j-1))+ &
3266 & cff*cff1*rhs_vbar(i,j))*tl_fac
3267#ifdef MASKING
3268
3269
3270 tl_vbar(i,j,knew)=tl_vbar(i,j,knew)*vmask(i,j)
3271#endif
3272#ifdef WET_DRY_NOT_YET
3273
3274
3275
3276
3277
3278
3279
3280 fac1=cff2/cff
3281
3282
3283
3284
3285 tl_rhs_vbar(i,j)=(tl_vbar(i,j,knew)* &
3286 & (dnew(i,j)+dnew(i,j-1))+ &
3287 & vbar(i,j,knew)* &
3288 & (tl_dnew(i,j)+tl_dnew(i,j-1))- &
3289 & tl_vbar(i,j,kstp)* &
3290 & (dstp(i,j)+dstp(i,j-1))- &
3291 & vbar(i,j,kstp)* &
3292 & (tl_dstp(i,j)+tl_dstp(i,j-1)))*fac1
3293#endif
3294 END DO
3295 END DO
3296 ELSE IF (predictor_2d_step(ng)) THEN
3297 cff1=dtfast(ng)
3298#ifdef WET_DRY_NOT_YET
3299 cff2=1.0_r8/cff1
3300#endif
3301 DO j=jstr,jend
3302 DO i=istru,iend
3303 cff=(pm(i,j)+pm(i-1,j))*(pn(i,j)+pn(i-1,j))
3304 fac=1.0_r8/(dnew(i,j)+dnew(i-1,j))
3305 tl_fac=-fac*fac*(tl_dnew(i,j)+tl_dnew(i-1,j))
3306
3307
3308
3309
3310 tl_ubar(i,j,knew)=(tl_ubar(i,j,kstp)* &
3311 & (dstp(i,j)+dstp(i-1,j))+ &
3312 & ubar(i,j,kstp)* &
3313 & (tl_dstp(i,j)+tl_dstp(i-1,j))+ &
3314 & cff*cff1*tl_rhs_ubar(i,j))*fac+ &
3315 & (ubar(i,j,kstp)* &
3316 & (dstp(i,j)+dstp(i-1,j))+ &
3317 & cff*cff1*rhs_ubar(i,j))*tl_fac
3318#ifdef MASKING
3319
3320
3321 tl_ubar(i,j,knew)=tl_ubar(i,j,knew)*umask(i,j)
3322#endif
3323#ifdef WET_DRY_NOT_YET
3324
3325
3326
3327
3328
3329
3330
3331 fac1=cff2/cff
3332
3333
3334
3335
3336 tl_rhs_ubar(i,j)=(tl_ubar(i,j,knew)* &
3337 & (dnew(i,j)+dnew(i-1,j))+ &
3338 & ubar(i,j,knew)* &
3339 & (tl_dnew(i,j)+tl_dnew(i-1,j))- &
3340 & tl_ubar(i,j,kstp)* &
3341 & (dstp(i,j)+dstp(i-1,j))- &
3342 & ubar(i,j,kstp)* &
3343 & (tl_dstp(i,j)+tl_dstp(i-1,j)))*fac1
3344#endif
3345 END DO
3346 END DO
3347 DO j=jstrv,jend
3348 DO i=istr,iend
3349 cff=(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
3350 fac=1.0_r8/(dnew(i,j)+dnew(i,j-1))
3351 tl_fac=-fac*fac*(tl_dnew(i,j)+tl_dnew(i,j-1))
3352
3353
3354
3355
3356 tl_vbar(i,j,knew)=(tl_vbar(i,j,kstp)* &
3357 & (dstp(i,j)+dstp(i,j-1))+ &
3358 & vbar(i,j,kstp)* &
3359 & (tl_dstp(i,j)+tl_dstp(i,j-1))+ &
3360 & cff*cff1*tl_rhs_vbar(i,j))*fac+ &
3361 & (vbar(i,j,kstp)* &
3362 & (dstp(i,j)+dstp(i,j-1))+ &
3363 & cff*cff1*rhs_vbar(i,j))*tl_fac
3364#ifdef MASKING
3365
3366
3367 tl_vbar(i,j,knew)=tl_vbar(i,j,knew)*vmask(i,j)
3368#endif
3369#ifdef WET_DRY_NOT_YET
3370
3371
3372
3373
3374
3375
3376
3377 fac1=cff2/cff
3378
3379
3380
3381
3382 tl_rhs_vbar(i,j)=(tl_vbar(i,j,knew)* &
3383 & (dnew(i,j)+dnew(i,j-1))+ &
3384 & vbar(i,j,knew)* &
3385 & (tl_dnew(i,j)+tl_dnew(i,j-1))- &
3386 & tl_vbar(i,j,kstp)* &
3387 & (dstp(i,j)+dstp(i,j-1))- &
3388 & vbar(i,j,kstp)* &
3389 & (tl_dstp(i,j)+tl_dstp(i,j-1)))*fac1
3390#endif
3391 END DO
3392 END DO
3393 ELSE IF (corrector_2d_step) THEN
3394 cff1=0.5_r8*dtfast(ng)*5.0_r8/12.0_r8
3395 cff2=0.5_r8*dtfast(ng)*8.0_r8/12.0_r8
3396 cff3=0.5_r8*dtfast(ng)*1.0_r8/12.0_r8
3397#ifdef WET_DRY_NOT_YET
3398 cff4=1.0_r8/cff1
3399#endif
3400 DO j=jstr,jend
3401 DO i=istru,iend
3402 cff=(pm(i,j)+pm(i-1,j))*(pn(i,j)+pn(i-1,j))
3403 fac=1.0_r8/(dnew(i,j)+dnew(i-1,j))
3404 tl_fac=-fac*fac*(tl_dnew(i,j)+tl_dnew(i-1,j))
3405
3406
3407
3408
3409
3410
3411 tl_ubar(i,j,knew)=(tl_ubar(i,j,kstp)* &
3412 & (dstp(i,j)+dstp(i-1,j))+ &
3413 & ubar(i,j,kstp)* &
3414 & (tl_dstp(i,j)+tl_dstp(i-1,j))+ &
3415 & cff*(cff1*tl_rhs_ubar(i,j)+ &
3416 & cff2*tl_rubar(i,j,kstp)- &
3417 & cff3*tl_rubar(i,j,ptsk)))*fac+ &
3418 & (ubar(i,j,kstp)* &
3419 & (dstp(i,j)+dstp(i-1,j))+ &
3420 & cff*(cff1*rhs_ubar(i,j)+ &
3421 & cff2*rubar(i,j,kstp)- &
3422 & cff3*rubar(i,j,ptsk)))*tl_fac
3423#ifdef MASKING
3424
3425
3426 tl_ubar(i,j,knew)=tl_ubar(i,j,knew)*umask(i,j)
3427#endif
3428#ifdef WET_DRY_NOT_YET
3429
3430
3431
3432
3433
3434
3435
3436 fac1=1.0_r8/cff
3437
3438
3439
3440
3441
3442
3443 tl_rhs_ubar(i,j)=((tl_ubar(i,j,knew)* &
3444 & (dnew(i,j)+dnew(i-1,j))+ &
3445 & ubar(i,j,knew)* &
3446 & (tl_dnew(i,j)+tl_dnew(i-1,j))- &
3447 & tl_ubar(i,j,kstp)* &
3448 & (dstp(i,j)+dstp(i-1,j))- &
3449 & ubar(i,j,kstp)* &
3450 & (tl_dstp(i,j)+tl_dstp(i-1,j)))*fac1- &
3451 & cff2*tl_rubar(i,j,kstp)+ &
3452 & cff3*tl_rubar(i,j,ptsk))*cff4
3453#endif
3454 END DO
3455 END DO
3456 DO j=jstrv,jend
3457 DO i=istr,iend
3458 cff=(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
3459 fac=1.0_r8/(dnew(i,j)+dnew(i,j-1))
3460 tl_fac=-fac*fac*(tl_dnew(i,j)+tl_dnew(i,j-1))
3461
3462
3463
3464
3465
3466
3467 tl_vbar(i,j,knew)=(tl_vbar(i,j,kstp)* &
3468 & (dstp(i,j)+dstp(i,j-1))+ &
3469 & vbar(i,j,kstp)* &
3470 & (tl_dstp(i,j)+tl_dstp(i,j-1))+ &
3471 & cff*(cff1*tl_rhs_vbar(i,j)+ &
3472 & cff2*tl_rvbar(i,j,kstp)- &
3473 & cff3*tl_rvbar(i,j,ptsk)))*fac+ &
3474 & (vbar(i,j,kstp)* &
3475 & (dstp(i,j)+dstp(i,j-1))+ &
3476 & cff*(cff1*rhs_vbar(i,j)+ &
3477 & cff2*rvbar(i,j,kstp)- &
3478 & cff3*rvbar(i,j,ptsk)))*tl_fac
3479#ifdef MASKING
3480
3481
3482 tl_vbar(i,j,knew)=tl_vbar(i,j,knew)*vmask(i,j)
3483#endif
3484#ifdef WET_DRY_NOT_YET
3485
3486
3487
3488
3489
3490
3491
3492 fac1=1.0_r8/cff
3493
3494
3495
3496
3497
3498
3499 tl_rhs_vbar(i,j)=((tl_vbar(i,j,knew)* &
3500 & (dnew(i,j)+dnew(i,j-1))+ &
3501 & vbar(i,j,knew)* &
3502 & (tl_dnew(i,j)+tl_dnew(i,j-1))- &
3503 & tl_vbar(i,j,kstp)* &
3504 & (dstp(i,j)+dstp(i,j-1))- &
3505 & vbar(i,j,kstp)* &
3506 & (tl_dstp(i,j)+tl_dstp(i,j-1)))*fac1- &
3507 & cff2*tl_rvbar(i,j,kstp)+ &
3508 & cff3*tl_rvbar(i,j,ptsk))*cff4
3509#endif
3510 END DO
3511 END DO
3512 END IF
3513#ifdef DIAGNOSTICS_UV
3514
3515
3516
3517
3518
3519# ifdef MASKING
3520
3521
3522
3523
3524
3525
3526
3527
3528
3529
3530
3531
3532
3533
3534
3535# endif
3536# ifdef SOLVE3D
3537
3538
3539
3540
3541
3542
3543
3544
3545
3546
3547
3548
3549
3550
3551
3552
3553
3554
3555
3556
3557
3558
3559
3560
3561
3562
3563
3564
3565
3566
3567
3568
3569
3570
3571
3572
3573
3574
3575
3576
3577
3578
3579
3580
3581
3582
3583
3584
3585
3586
3587
3588
3589
3590
3591
3592# else
3593
3594
3595
3596
3597
3598
3599
3600
3601
3602
3603
3604
3605
3606
3607
3608
3609
3610
3611
3612
3613
3614
3615
3616
3617
3618
3619
3620
3621
3622
3623
3624
3625
3626
3627
3628
3629
3630
3631
3632
3633
3634
3635
3636
3637
3638
3639
3640
3641
3642
3643
3644
3645
3646
3647
3648
3649
3650
3651
3652
3653
3654
3655
3656
3657
3658
3659
3660
3661
3662
3663
3664
3665
3666
3667
3668
3669
3670
3671# endif
3672#endif
3673
3674
3675
3676
3677 IF (predictor_2d_step(ng)) THEN
3678 DO j=jstr,jend
3679 DO i=istru,iend
3680
3681
3682 tl_rubar(i,j,krhs)=tl_rhs_ubar(i,j)
3683 END DO
3684 END DO
3685 DO j=jstrv,jend
3686 DO i=istr,iend
3687
3688
3689 tl_rvbar(i,j,krhs)=tl_rhs_vbar(i,j)
3690 END DO
3691 END DO
3692#ifdef DIAGNOSTICS_UV
3693
3694
3695
3696
3697
3698
3699
3700
3701
3702
3703
3704
3705#endif
3706 END IF
3707
3708
3709
3710
3711
3712
3713
3714
3715
3716
3717
3718 CALL tl_u2dbc_tile (ng, tile, &
3719 & lbi, ubi, lbj, ubj, &
3720 & imins, imaxs, jmins, jmaxs, &
3721 & krhs, kstp, knew, &
3722 & ubar, vbar, zeta, &
3723 & tl_ubar, tl_vbar, tl_zeta)
3724
3725
3726
3727
3728
3729
3730 CALL tl_v2dbc_tile (ng, tile, &
3731 & lbi, ubi, lbj, ubj, &
3732 & imins, imaxs, jmins, jmaxs, &
3733 & krhs, kstp, knew, &
3734 & ubar, vbar, zeta, &
3735 & tl_ubar, tl_vbar, tl_zeta)
3736
3737
3738
3739
3740 IF (any(tl_volcons(:,ng))) THEN
3741 CALL tl_obc_flux_tile (ng, tile, &
3742 & lbi, ubi, lbj, ubj, &
3743 & imins, imaxs, jmins, jmaxs, &
3744 & knew, &
3745#ifdef MASKING
3746 & umask, vmask, &
3747#endif
3748 & h, tl_h, om_v, on_u, &
3749 & ubar, vbar, zeta, &
3750 & tl_ubar, tl_vbar, tl_zeta)
3751 END IF
3752
3753
3754
3755
3756
3757
3758
3759
3760 IF (luvsrc(ng)) THEN
3761 DO is=1,nsrc(ng)
3762 i=sources(ng)%Isrc(is)
3763 j=sources(ng)%Jsrc(is)
3764 IF (((istrr.le.i).and.(i.le.iendr)).and. &
3765 & ((jstrr.le.j).and.(j.le.jendr))) THEN
3766 IF (int(sources(ng)%Dsrc(is)).eq.0) THEN
3767 cff=1.0_r8/(on_u(i,j)* &
3768 & 0.5_r8*(zeta(i-1,j,knew)+h(i-1,j)+ &
3769 & zeta(i ,j,knew)+h(i ,j)))
3770 tl_cff=-cff*cff*on_u(i,j)* &
3771 & 0.5_r8*(tl_zeta(i-1,j,knew)+tl_h(i-1,j)+ &
3772 & tl_zeta(i ,j,knew)+tl_h(i ,j))
3773
3774
3775 tl_ubar(i,j,knew)=sources(ng)%tl_Qbar(is)*cff+ &
3776 & sources(ng)%Qbar(is)*tl_cff
3777 ELSE IF (int(sources(ng)%Dsrc(is)).eq.1) THEN
3778 cff=1.0_r8/(om_v(i,j)* &
3779 & 0.5_r8*(zeta(i,j-1,knew)+h(i,j-1)+ &
3780 & zeta(i,j ,knew)+h(i,j )))
3781 tl_cff=-cff*cff*om_v(i,j)* &
3782 & 0.5_r8*(tl_zeta(i,j-1,knew)+tl_h(i,j-1)+ &
3783 & tl_zeta(i,j ,knew)+tl_h(i,j ))
3784
3785
3786 tl_vbar(i,j,knew)=sources(ng)%tl_Qbar(is)*cff+ &
3787 & sources(ng)%Qbar(is)*tl_cff
3788 END IF
3789 END IF
3790 END DO
3791 END IF
3792
3793
3794
3795
3796
3797 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
3798
3799
3800
3801
3802 CALL exchange_u2d_tile (ng, tile, &
3803 & lbi, ubi, lbj, ubj, &
3804 & tl_ubar(:,:,knew))
3805
3806
3807
3808
3809 CALL exchange_v2d_tile (ng, tile, &
3810 & lbi, ubi, lbj, ubj, &
3811 & tl_vbar(:,:,knew))
3812 END IF
3813
3814#ifdef DISTRIBUTE
3815
3816
3817
3818
3819
3820
3821
3822 CALL mp_exchange2d (ng, tile, itlm, 2, &
3823 & lbi, ubi, lbj, ubj, &
3824 & nghostpoints, &
3825 & ewperiodic(ng), nsperiodic(ng), &
3826 & tl_ubar(:,:,knew), &
3827 & tl_vbar(:,:,knew))
3828#endif
3829
3830 RETURN