ROMS
Loading...
Searching...
No Matches
mod_grid Module Reference

Data Types

type  t_grid
 

Functions/Subroutines

subroutine, public allocate_grid (ng, extract_flag, lbi, ubi, lbj, ubj, lbij, ubij)
 
subroutine, public deallocate_grid (ng)
 
subroutine, public initialize_grid (ng, tile, model)
 

Variables

type(t_grid), dimension(:), allocatable grid
 

Function/Subroutine Documentation

◆ allocate_grid()

subroutine, public mod_grid::allocate_grid ( integer, intent(in) ng,
integer, intent(in) extract_flag,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) lbij,
integer, intent(in) ubij )

Definition at line 369 of file mod_grid.F.

371!
372!=======================================================================
373! !
374! This routine allocates all variables in the module for all nested !
375! grids. !
376! !
377!=======================================================================
378!
379 USE mod_param
380!
381! Local variable declarations.
382!
383 integer, intent(in) :: ng, Extract_Flag
384 integer, intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
385!
386! Local variable declarations.
387!
388 integer :: my_size
389!
390 real(r8) :: size2d
391!
392!-----------------------------------------------------------------------
393! Allocate and initialize module variables.
394!-----------------------------------------------------------------------
395!
396 IF (ng.eq.1) allocate ( grid(ngrids) )
397!
398! Set horizontal array size.
399!
400 size2d=real((ubi-lbi+1)*(ubj-lbj+1),r8)
401!
402! Nonlinear model state.
403!
404#if defined MASKING && defined PROPAGATOR
405 allocate ( grid(ng) % IJwaterR(lbi:ubi,lbj:ubj) )
406 dmem(ng)=dmem(ng)+size2d
407
408 allocate ( grid(ng) % IJwaterU(lbi:ubi,lbj:ubj) )
409 dmem(ng)=dmem(ng)+size2d
410
411 allocate ( grid(ng) % IJwaterV(lbi:ubi,lbj:ubj) )
412 dmem(ng)=dmem(ng)+size2d
413#endif
414
415 allocate ( grid(ng) % angler(lbi:ubi,lbj:ubj) )
416 dmem(ng)=dmem(ng)+size2d
417
418 allocate ( grid(ng) % CosAngler(lbi:ubi,lbj:ubj) )
419 dmem(ng)=dmem(ng)+size2d
420
421 allocate ( grid(ng) % SinAngler(lbi:ubi,lbj:ubj) )
422 dmem(ng)=dmem(ng)+size2d
423
424#if defined TIDE_GENERATING_FORCES
425 allocate ( grid(ng) % Cos2Lat(lbi:ubi,lbj:ubj) )
426 dmem(ng)=dmem(ng)+size2d
427
428 allocate ( grid(ng) % SinLat2(lbi:ubi,lbj:ubj) )
429 dmem(ng)=dmem(ng)+size2d
430#endif
431
432#if defined CURVGRID && defined UV_ADV
433 allocate ( grid(ng) % dmde(lbi:ubi,lbj:ubj) )
434 dmem(ng)=dmem(ng)+size2d
435
436 allocate ( grid(ng) % dndx(lbi:ubi,lbj:ubj) )
437 dmem(ng)=dmem(ng)+size2d
438#endif
439
440 allocate ( grid(ng) % f(lbi:ubi,lbj:ubj) )
441 dmem(ng)=dmem(ng)+size2d
442
443 allocate ( grid(ng) % fomn(lbi:ubi,lbj:ubj) )
444 dmem(ng)=dmem(ng)+size2d
445
446 allocate ( grid(ng) % grdscl(lbi:ubi,lbj:ubj) )
447 dmem(ng)=dmem(ng)+size2d
448
449 allocate ( grid(ng) % h(lbi:ubi,lbj:ubj) )
450 dmem(ng)=dmem(ng)+size2d
451
452 allocate ( grid(ng) % latp(lbi:ubi,lbj:ubj) )
453 dmem(ng)=dmem(ng)+size2d
454
455 allocate ( grid(ng) % latr(lbi:ubi,lbj:ubj) )
456 dmem(ng)=dmem(ng)+size2d
457
458 allocate ( grid(ng) % latu(lbi:ubi,lbj:ubj) )
459 dmem(ng)=dmem(ng)+size2d
460
461 allocate ( grid(ng) % latv(lbi:ubi,lbj:ubj) )
462 dmem(ng)=dmem(ng)+size2d
463
464 allocate ( grid(ng) % lonp(lbi:ubi,lbj:ubj))
465 dmem(ng)=dmem(ng)+size2d
466
467 allocate ( grid(ng) % lonr(lbi:ubi,lbj:ubj))
468 dmem(ng)=dmem(ng)+size2d
469
470 allocate ( grid(ng) % lonu(lbi:ubi,lbj:ubj))
471 dmem(ng)=dmem(ng)+size2d
472
473 allocate ( grid(ng) % lonv(lbi:ubi,lbj:ubj))
474 dmem(ng)=dmem(ng)+size2d
475
476 allocate ( grid(ng) % Mylon(lbi:ubi,lbj:ubj))
477 dmem(ng)=dmem(ng)+size2d
478
479 allocate ( grid(ng) % omn(lbi:ubi,lbj:ubj) )
480 dmem(ng)=dmem(ng)+size2d
481
482 allocate ( grid(ng) % om_p(lbi:ubi,lbj:ubj) )
483 dmem(ng)=dmem(ng)+size2d
484
485 allocate ( grid(ng) % om_r(lbi:ubi,lbj:ubj) )
486 dmem(ng)=dmem(ng)+size2d
487
488 allocate ( grid(ng) % om_u(lbi:ubi,lbj:ubj) )
489 dmem(ng)=dmem(ng)+size2d
490
491 allocate ( grid(ng) % om_v(lbi:ubi,lbj:ubj) )
492 dmem(ng)=dmem(ng)+size2d
493
494 allocate ( grid(ng) % on_p(lbi:ubi,lbj:ubj) )
495 dmem(ng)=dmem(ng)+size2d
496
497 allocate ( grid(ng) % on_r(lbi:ubi,lbj:ubj) )
498 dmem(ng)=dmem(ng)+size2d
499
500 allocate ( grid(ng) % on_u(lbi:ubi,lbj:ubj) )
501 dmem(ng)=dmem(ng)+size2d
502
503 allocate ( grid(ng) % on_v(lbi:ubi,lbj:ubj) )
504 dmem(ng)=dmem(ng)+size2d
505
506 allocate ( grid(ng) % pm(lbi:ubi,lbj:ubj) )
507 dmem(ng)=dmem(ng)+size2d
508
509 allocate ( grid(ng) % pn(lbi:ubi,lbj:ubj) )
510 dmem(ng)=dmem(ng)+size2d
511
512 allocate ( grid(ng) % pmon_p(lbi:ubi,lbj:ubj) )
513 dmem(ng)=dmem(ng)+size2d
514
515 allocate ( grid(ng) % pmon_r(lbi:ubi,lbj:ubj) )
516 dmem(ng)=dmem(ng)+size2d
517
518 allocate ( grid(ng) % pmon_u(lbi:ubi,lbj:ubj) )
519 dmem(ng)=dmem(ng)+size2d
520
521 allocate ( grid(ng) % pmon_v(lbi:ubi,lbj:ubj) )
522 dmem(ng)=dmem(ng)+size2d
523
524 allocate ( grid(ng) % pnom_p(lbi:ubi,lbj:ubj) )
525 dmem(ng)=dmem(ng)+size2d
526
527 allocate ( grid(ng) % pnom_r(lbi:ubi,lbj:ubj) )
528 dmem(ng)=dmem(ng)+size2d
529
530 allocate ( grid(ng) % pnom_u(lbi:ubi,lbj:ubj) )
531 dmem(ng)=dmem(ng)+size2d
532
533 allocate ( grid(ng) % pnom_v(lbi:ubi,lbj:ubj) )
534 dmem(ng)=dmem(ng)+size2d
535
536#if defined UV_LOGDRAG || defined GLS_MIXING || \
537 defined bbl_model || defined sediment
538 allocate ( grid(ng) % ZoBot(lbi:ubi,lbj:ubj) )
539 dmem(ng)=dmem(ng)+size2d
540#endif
541
542 allocate ( grid(ng) % rdrag(lbi:ubi,lbj:ubj) )
543 dmem(ng)=dmem(ng)+size2d
544
545#if defined UV_QDRAG
546 allocate ( grid(ng) % rdrag2(lbi:ubi,lbj:ubj) )
547 dmem(ng)=dmem(ng)+size2d
548#endif
549
550 allocate ( grid(ng) % xp(lbi:ubi,lbj:ubj) )
551 dmem(ng)=dmem(ng)+size2d
552
553 allocate ( grid(ng) % xr(lbi:ubi,lbj:ubj) )
554 dmem(ng)=dmem(ng)+size2d
555
556 allocate ( grid(ng) % xu(lbi:ubi,lbj:ubj) )
557 dmem(ng)=dmem(ng)+size2d
558
559 allocate ( grid(ng) % xv(lbi:ubi,lbj:ubj) )
560 dmem(ng)=dmem(ng)+size2d
561
562 allocate ( grid(ng) % yp(lbi:ubi,lbj:ubj) )
563 dmem(ng)=dmem(ng)+size2d
564
565 allocate ( grid(ng) % yr(lbi:ubi,lbj:ubj) )
566 dmem(ng)=dmem(ng)+size2d
567
568 allocate ( grid(ng) % yu(lbi:ubi,lbj:ubj) )
569 dmem(ng)=dmem(ng)+size2d
570
571 allocate ( grid(ng) % yv(lbi:ubi,lbj:ubj) )
572 dmem(ng)=dmem(ng)+size2d
573
574#ifdef SOLVE3D
575 allocate ( grid(ng) % Hz(lbi:ubi,lbj:ubj,n(ng)) )
576 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
577
578# ifdef ADJUST_BOUNDARY
579 allocate ( grid(ng) % Hz_bry(lbij:ubij,n(ng),4) )
580 dmem(ng)=dmem(ng)+4.0_r8*real((ubij-lbij)*n(ng),r8)
581# endif
582
583 allocate ( grid(ng) % Huon(lbi:ubi,lbj:ubj,n(ng)) )
584 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
585
586 allocate ( grid(ng) % Hvom(lbi:ubi,lbj:ubj,n(ng)) )
587 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
588
589 allocate ( grid(ng) % z0_r(lbi:ubi,lbj:ubj,n(ng)) )
590 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
591
592 allocate ( grid(ng) % z0_w(lbi:ubi,lbj:ubj,0:n(ng)) )
593 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
594
595 allocate ( grid(ng) % z_r(lbi:ubi,lbj:ubj,n(ng)) )
596 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
597
598 allocate ( grid(ng) % z_v(lbi:ubi,lbj:ubj,n(ng)) )
599 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
600
601 allocate ( grid(ng) % z_w(lbi:ubi,lbj:ubj,0:n(ng)) )
602 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
603
604# ifdef ICESHELF
605 allocate ( grid(ng) % IcePress(lbi:ubi,lbj:ubj) )
606 dmem(ng)=dmem(ng)+size2d
607
608 allocate ( grid(ng) % zice(lbi:ubi,lbj:ubj) )
609 dmem(ng)=dmem(ng)+size2d
610# endif
611
612#endif
613
614#ifdef MASKING
615 allocate ( grid(ng) % pmask(lbi:ubi,lbj:ubj) )
616 dmem(ng)=dmem(ng)+size2d
617
618 allocate ( grid(ng) % rmask(lbi:ubi,lbj:ubj) )
619 dmem(ng)=dmem(ng)+size2d
620
621 allocate ( grid(ng) % umask(lbi:ubi,lbj:ubj) )
622 dmem(ng)=dmem(ng)+size2d
623
624 allocate ( grid(ng) % vmask(lbi:ubi,lbj:ubj) )
625 dmem(ng)=dmem(ng)+size2d
626# if defined AVERAGES || \
627 (defined ad_averages && defined adjoint) || \
628 (defined rp_averages && defined tl_ioms) || \
629 (defined tl_averages && defined tangent)
630 allocate ( grid(ng) % pmask_avg(lbi:ubi,lbj:ubj) )
631 dmem(ng)=dmem(ng)+size2d
632
633 allocate ( grid(ng) % rmask_avg(lbi:ubi,lbj:ubj) )
634 dmem(ng)=dmem(ng)+size2d
635
636 allocate ( grid(ng) % umask_avg(lbi:ubi,lbj:ubj) )
637 dmem(ng)=dmem(ng)+size2d
638
639 allocate ( grid(ng) % vmask_avg(lbi:ubi,lbj:ubj) )
640 dmem(ng)=dmem(ng)+size2d
641# endif
642# ifdef DIAGNOSTICS
643 allocate ( grid(ng) % pmask_dia(lbi:ubi,lbj:ubj) )
644 dmem(ng)=dmem(ng)+size2d
645
646 allocate ( grid(ng) % rmask_dia(lbi:ubi,lbj:ubj) )
647 dmem(ng)=dmem(ng)+size2d
648
649 allocate ( grid(ng) % umask_dia(lbi:ubi,lbj:ubj) )
650 dmem(ng)=dmem(ng)+size2d
651
652 allocate ( grid(ng) % vmask_dia(lbi:ubi,lbj:ubj) )
653 dmem(ng)=dmem(ng)+size2d
654# endif
655 allocate ( grid(ng) % pmask_full(lbi:ubi,lbj:ubj) )
656 dmem(ng)=dmem(ng)+size2d
657
658 allocate ( grid(ng) % rmask_full(lbi:ubi,lbj:ubj) )
659 dmem(ng)=dmem(ng)+size2d
660
661 allocate ( grid(ng) % umask_full(lbi:ubi,lbj:ubj) )
662 dmem(ng)=dmem(ng)+size2d
663
664 allocate ( grid(ng) % vmask_full(lbi:ubi,lbj:ubj) )
665 dmem(ng)=dmem(ng)+size2d
666#endif
667
668#ifdef WET_DRY
669 allocate ( grid(ng) % pmask_wet(lbi:ubi,lbj:ubj) )
670 dmem(ng)=dmem(ng)+size2d
671
672 allocate ( grid(ng) % rmask_wet(lbi:ubi,lbj:ubj) )
673 dmem(ng)=dmem(ng)+size2d
674
675 allocate ( grid(ng) % umask_wet(lbi:ubi,lbj:ubj) )
676 dmem(ng)=dmem(ng)+size2d
677
678 allocate ( grid(ng) % vmask_wet(lbi:ubi,lbj:ubj) )
679 dmem(ng)=dmem(ng)+size2d
680# ifdef SOLVE3D
681 allocate ( grid(ng) % rmask_wet_avg(lbi:ubi,lbj:ubj) )
682 dmem(ng)=dmem(ng)+size2d
683# endif
684#endif
685
686#if defined AD_SENSITIVITY || defined I4DVAR_ANA_SENSITIVITY || \
687 defined opt_observations || defined sensitivity_4dvar || \
688 defined so_semi
689 allocate ( grid(ng) % Rscope(lbi:ubi,lbj:ubj) )
690 dmem(ng)=dmem(ng)+size2d
691
692 allocate ( grid(ng) % Uscope(lbi:ubi,lbj:ubj) )
693 dmem(ng)=dmem(ng)+size2d
694
695 allocate ( grid(ng) % Vscope(lbi:ubi,lbj:ubj) )
696 dmem(ng)=dmem(ng)+size2d
697#endif
698
699#if defined TANGENT || defined TL_IOMS
700!
701! Tangent linear model state.
702!
703 allocate ( grid(ng) % tl_h(lbi:ubi,lbj:ubj) )
704 dmem(ng)=dmem(ng)+size2d
705
706# ifdef SOLVE3D
707 allocate ( grid(ng) % tl_Hz(lbi:ubi,lbj:ubj,n(ng)) )
708 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
709
710# ifdef ADJUST_BOUNDARY
711 allocate ( grid(ng) % tl_Hz_bry(lbij:ubij,n(ng),4) )
712 dmem(ng)=dmem(ng)+4.0_r8*real((ubij-lbij)*n(ng),r8)*size2d
713# endif
714
715 allocate ( grid(ng) % tl_Huon(lbi:ubi,lbj:ubj,n(ng)) )
716 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
717
718 allocate ( grid(ng) % tl_Hvom(lbi:ubi,lbj:ubj,n(ng)) )
719 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
720
721 allocate ( grid(ng) % tl_z_r(lbi:ubi,lbj:ubj,n(ng)) )
722 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
723
724 allocate ( grid(ng) % tl_z_w(lbi:ubi,lbj:ubj,0:n(ng)) )
725 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
726# endif
727#endif
728
729#ifdef ADJOINT
730!
731! Adjoint model state.
732!
733 allocate ( grid(ng) % ad_h(lbi:ubi,lbj:ubj) )
734 dmem(ng)=dmem(ng)+size2d
735
736# ifdef SOLVE3D
737 allocate ( grid(ng) % ad_Hz(lbi:ubi,lbj:ubj,n(ng)) )
738 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
739
740# ifdef ADJUST_BOUNDARY
741 allocate ( grid(ng) % ad_Hz_bry(lbij:ubij,n(ng),4) )
742 dmem(ng)=dmem(ng)+4.0_r8*real((ubij-lbij)*n(ng),r8)*size2d
743# endif
744
745 allocate ( grid(ng) % ad_Huon(lbi:ubi,lbj:ubj,n(ng)) )
746 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
747
748 allocate ( grid(ng) % ad_Hvom(lbi:ubi,lbj:ubj,n(ng)) )
749 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
750
751 allocate ( grid(ng) % ad_z_r(lbi:ubi,lbj:ubj,n(ng)) )
752 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
753
754 allocate ( grid(ng) % ad_z_w(lbi:ubi,lbj:ubj,0:n(ng)) )
755 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
756# endif
757#endif
758
759#ifdef GRID_EXTRACT
760!
761! Allocate global arrays for extraction by interpolation.
762!
763 IF (extract_flag.ge.1) THEN
764 my_size=(iobounds(ng)%IUB_psi-iobounds(ng)%ILB_psi+1)* &
765 & (iobounds(ng)%JUB_psi-iobounds(ng)%JLB_psi+1)
766
767 allocate ( grid(ng) % Gangle_psi(my_size) )
768 dmem(ng)=dmem(ng)+real(my_size,r8)
769
770# ifdef MASKING
771 allocate ( grid(ng) % Gmask_psi(my_size) )
772 dmem(ng)=dmem(ng)+real(my_size,r8)
773# endif
774
775 allocate ( grid(ng) % Gx_psi(my_size) )
776 dmem(ng)=dmem(ng)+real(my_size,r8)
777
778 allocate ( grid(ng) % Gy_psi(my_size) )
779 dmem(ng)=dmem(ng)+real(my_size,r8)
780!
781 my_size=(iobounds(ng)%IUB_rho-iobounds(ng)%ILB_rho+1)* &
782 & (iobounds(ng)%JUB_rho-iobounds(ng)%JLB_rho+1)
783
784 allocate ( grid(ng) % Gangle_rho(my_size) )
785 dmem(ng)=dmem(ng)+real(my_size,r8)
786
787# ifdef MASKING
788 allocate ( grid(ng) % Gmask_rho(my_size) )
789 dmem(ng)=dmem(ng)+real(my_size,r8)
790# endif
791
792 allocate ( grid(ng) % Gx_rho(my_size) )
793 dmem(ng)=dmem(ng)+real(my_size,r8)
794
795 allocate ( grid(ng) % Gy_rho(my_size) )
796 dmem(ng)=dmem(ng)+real(my_size,r8)
797!
798 my_size=(iobounds(ng)%IUB_u-iobounds(ng)%ILB_u+1)* &
799 & (iobounds(ng)%JUB_u-iobounds(ng)%JLB_u+1)
800
801 allocate ( grid(ng) % Gangle_u(my_size) )
802 dmem(ng)=dmem(ng)+real(my_size,r8)
803
804# ifdef MASKING
805 allocate ( grid(ng) % Gmask_u(my_size) )
806 dmem(ng)=dmem(ng)+real(my_size,r8)
807# endif
808
809 allocate ( grid(ng) % Gx_u(my_size) )
810 dmem(ng)=dmem(ng)+real(my_size,r8)
811
812 allocate ( grid(ng) % Gy_u(my_size) )
813 dmem(ng)=dmem(ng)+real(my_size,r8)
814!
815 my_size=(iobounds(ng)%IUB_v-iobounds(ng)%ILB_v+1)* &
816 & (iobounds(ng)%JUB_v-iobounds(ng)%JLB_v+1)
817
818 allocate ( grid(ng) % Gangle_v(my_size) )
819 dmem(ng)=dmem(ng)+real(my_size,r8)
820
821# ifdef MASKING
822 allocate ( grid(ng) % Gmask_v(my_size) )
823 dmem(ng)=dmem(ng)+real(my_size,r8)
824# endif
825
826 allocate ( grid(ng) % Gx_v(my_size) )
827 dmem(ng)=dmem(ng)+real(my_size,r8)
828
829 allocate ( grid(ng) % Gy_v(my_size) )
830 dmem(ng)=dmem(ng)+real(my_size,r8)
831 END IF
832#endif
833!
834 RETURN
integer, dimension(:), allocatable n
Definition mod_param.F:479
real(r8), dimension(:), allocatable dmem
Definition mod_param.F:137
type(t_iobounds), dimension(:), allocatable iobounds
Definition mod_param.F:282
integer ngrids
Definition mod_param.F:113

References mod_param::dmem, grid, mod_param::iobounds, mod_param::n, mod_param::ngrids, and mod_kinds::r8.

Referenced by mod_arrays::roms_allocate_arrays().

Here is the caller graph for this function:

◆ deallocate_grid()

subroutine, public mod_grid::deallocate_grid ( integer, intent(in) ng)

Definition at line 837 of file mod_grid.F.

838!
839!=======================================================================
840! !
841! This routine deallocates all variables in the module for all nested !
842! grids. !
843! !
844!=======================================================================
845!
846 USE mod_param, ONLY : ngrids
847#ifdef SUBOBJECT_DEALLOCATION
848 USE destroy_mod, ONLY : destroy
849#endif
850!
851! Imported variable declarations.
852!
853 integer, intent(in) :: ng
854!
855! Local variable declarations.
856!
857 character (len=*), parameter :: MyFile = &
858 & __FILE__//", deallocate_grid"
859
860#ifdef SUBOBJECT_DEALLOCATION
861!
862!-----------------------------------------------------------------------
863! Deallocate each variable in the derived-type T_GRID structure
864! separately.
865!-----------------------------------------------------------------------
866!
867! Nonlinear model state.
868!
869# if defined MASKING && defined PROPAGATOR
870 IF (.not.destroy(ng, grid(ng)%IJwaterR, myfile, &
871 & __line__, 'GRID(ng)%IJwaterR')) RETURN
872
873 IF (.not.destroy(ng, grid(ng)%IJwaterU, myfile, &
874 & __line__, 'GRID(ng)%IJwaterU')) RETURN
875
876 IF (.not.destroy(ng, grid(ng)%IJwaterV, myfile, &
877 & __line__, 'GRID(ng)%IJwaterV')) RETURN
878# endif
879
880 IF (.not.destroy(ng, grid(ng)%angler, myfile, &
881 & __line__, 'GRID(ng)%angler')) RETURN
882
883 IF (.not.destroy(ng, grid(ng)%CosAngler, myfile, &
884 & __line__, 'GRID(ng)%CosAngler')) RETURN
885
886 IF (.not.destroy(ng, grid(ng)%SinAngler, myfile, &
887 & __line__, 'GRID(ng)%SinAngler')) RETURN
888
889# if defined TIDE_GENERATING_FORCES
890 IF (.not.destroy(ng, grid(ng)%Cos2Lat, myfile, &
891 & __line__, 'GRID(ng)%Cos2Lat')) RETURN
892
893 IF (.not.destroy(ng, grid(ng)%SinLat2, myfile, &
894 & __line__, 'GRID(ng)%SinLat2')) RETURN
895# endif
896
897# if defined CURVGRID && defined UV_ADV
898 IF (.not.destroy(ng, grid(ng)%dmde, myfile, &
899 & __line__, 'GRID(ng)%dmde')) RETURN
900
901 IF (.not.destroy(ng, grid(ng)%dndx, myfile, &
902 & __line__, 'GRID(ng)%dndx')) RETURN
903# endif
904
905 IF (.not.destroy(ng, grid(ng)%f, myfile, &
906 & __line__, 'GRID(ng)%f')) RETURN
907
908 IF (.not.destroy(ng, grid(ng)%fomn, myfile, &
909 & __line__, 'GRID(ng)%fomn')) RETURN
910
911 IF (.not.destroy(ng, grid(ng)%grdscl, myfile, &
912 & __line__, 'GRID(ng)%grdscl')) RETURN
913
914 IF (.not.destroy(ng, grid(ng)%h, myfile, &
915 & __line__, 'GRID(ng)%h')) RETURN
916
917 IF (.not.destroy(ng, grid(ng)%latp, myfile, &
918 & __line__, 'GRID(ng)%latp')) RETURN
919
920 IF (.not.destroy(ng, grid(ng)%latr, myfile, &
921 & __line__, 'GRID(ng)%latr')) RETURN
922
923 IF (.not.destroy(ng, grid(ng)%latu, myfile, &
924 & __line__, 'GRID(ng)%latu')) RETURN
925
926 IF (.not.destroy(ng, grid(ng)%latv, myfile, &
927 & __line__, 'GRID(ng)%latv')) RETURN
928
929 IF (.not.destroy(ng, grid(ng)%lonp, myfile, &
930 & __line__, 'GRID(ng)%lonp')) RETURN
931
932 IF (.not.destroy(ng, grid(ng)%lonr, myfile, &
933 & __line__, 'GRID(ng)%lonr')) RETURN
934
935 IF (.not.destroy(ng, grid(ng)%lonu, myfile, &
936 & __line__, 'GRID(ng)%lonu')) RETURN
937
938 IF (.not.destroy(ng, grid(ng)%lonv, myfile, &
939 & __line__, 'GRID(ng)%lonv')) RETURN
940
941 IF (.not.destroy(ng, grid(ng)%Mylon, myfile, &
942 & __line__, 'GRID(ng)%Mylon')) RETURN
943
944 IF (.not.destroy(ng, grid(ng)%omn, myfile, &
945 & __line__, 'GRID(ng)%omn')) RETURN
946
947 IF (.not.destroy(ng, grid(ng)%om_p, myfile, &
948 & __line__, 'GRID(ng)%om_p')) RETURN
949
950 IF (.not.destroy(ng, grid(ng)%om_r, myfile, &
951 & __line__, 'GRID(ng)%om_r')) RETURN
952
953 IF (.not.destroy(ng, grid(ng)%om_u, myfile, &
954 & __line__, 'GRID(ng)%om_u')) RETURN
955
956 IF (.not.destroy(ng, grid(ng)%om_v, myfile, &
957 & __line__, 'GRID(ng)%om_v')) RETURN
958
959 IF (.not.destroy(ng, grid(ng)%on_p, myfile, &
960 & __line__, 'GRID(ng)%on_p')) RETURN
961
962 IF (.not.destroy(ng, grid(ng)%on_r, myfile, &
963 & __line__, 'GRID(ng)%on_r')) RETURN
964
965 IF (.not.destroy(ng, grid(ng)%on_u, myfile, &
966 & __line__, 'GRID(ng)%on_u')) RETURN
967
968 IF (.not.destroy(ng, grid(ng)%on_v, myfile, &
969 & __line__, 'GRID(ng)%on_v')) RETURN
970
971 IF (.not.destroy(ng, grid(ng)%pm, myfile, &
972 & __line__, 'GRID(ng)%pm')) RETURN
973
974 IF (.not.destroy(ng, grid(ng)%pn, myfile, &
975 & __line__, 'GRID(ng)%pn')) RETURN
976
977 IF (.not.destroy(ng, grid(ng)%pmon_p, myfile, &
978 & __line__, 'GRID(ng)%pmon_p')) RETURN
979
980 IF (.not.destroy(ng, grid(ng)%pmon_r, myfile, &
981 & __line__, 'GRID(ng)%pmon_r')) RETURN
982
983 IF (.not.destroy(ng, grid(ng)%pmon_u, myfile, &
984 & __line__, 'GRID(ng)%pmon_u')) RETURN
985
986 IF (.not.destroy(ng, grid(ng)%pmon_v, myfile, &
987 & __line__, 'GRID(ng)%pmon_v')) RETURN
988
989 IF (.not.destroy(ng, grid(ng)%pnom_p, myfile, &
990 & __line__, 'GRID(ng)%pnom_p')) RETURN
991
992 IF (.not.destroy(ng, grid(ng)%pnom_r, myfile, &
993 & __line__, 'GRID(ng)%pnom_r')) RETURN
994
995 IF (.not.destroy(ng, grid(ng)%pnom_u, myfile, &
996 & __line__, 'GRID(ng)%pnom_u')) RETURN
997
998 IF (.not.destroy(ng, grid(ng)%pnom_v, myfile, &
999 & __line__, 'GRID(ng)%pnom_v')) RETURN
1000
1001# if defined UV_LOGDRAG || defined GLS_MIXING || \
1002 defined bbl_model || defined sediment
1003 IF (.not.destroy(ng, grid(ng)%ZoBot, myfile, &
1004 & __line__, 'GRID(ng)%ZoBot')) RETURN
1005# endif
1006
1007 IF (.not.destroy(ng, grid(ng)%rdrag, myfile, &
1008 & __line__, 'GRID(ng)%rdrag')) RETURN
1009
1010# if defined UV_QDRAG
1011 IF (.not.destroy(ng, grid(ng)%rdrag2, myfile, &
1012 & __line__, 'GRID(ng)%rdrag2')) RETURN
1013# endif
1014
1015 IF (.not.destroy(ng, grid(ng)%xp, myfile, &
1016 & __line__, 'GRID(ng)%xp')) RETURN
1017
1018 IF (.not.destroy(ng, grid(ng)%xr, myfile, &
1019 & __line__, 'GRID(ng)%xr')) RETURN
1020
1021 IF (.not.destroy(ng, grid(ng)%xu, myfile, &
1022 & __line__, 'GRID(ng)%xu')) RETURN
1023
1024 IF (.not.destroy(ng, grid(ng)%xv, myfile, &
1025 & __line__, 'GRID(ng)%xv')) RETURN
1026
1027 IF (.not.destroy(ng, grid(ng)%yp, myfile, &
1028 & __line__, 'GRID(ng)%yp')) RETURN
1029
1030 IF (.not.destroy(ng, grid(ng)%yr, myfile, &
1031 & __line__, 'GRID(ng)%yr')) RETURN
1032
1033 IF (.not.destroy(ng, grid(ng)%yu, myfile, &
1034 & __line__, 'GRID(ng)%yu')) RETURN
1035
1036 IF (.not.destroy(ng, grid(ng)%yv, myfile, &
1037 & __line__, 'GRID(ng)%yv')) RETURN
1038
1039# ifdef SOLVE3D
1040 IF (.not.destroy(ng, grid(ng)%Hz, myfile, &
1041 & __line__, 'GRID(ng)%Hz')) RETURN
1042
1043# ifdef ADJUST_BOUNDARY
1044 IF (.not.destroy(ng, grid(ng)%Hz_bry, myfile, &
1045 & __line__, 'GRID(ng)%Hz_bry')) RETURN
1046# endif
1047
1048 IF (.not.destroy(ng, grid(ng)%Huon, myfile, &
1049 & __line__, 'GRID(ng)%Huon')) RETURN
1050
1051 IF (.not.destroy(ng, grid(ng)%Hvom, myfile, &
1052 & __line__, 'GRID(ng)%Hvom')) RETURN
1053
1054 IF (.not.destroy(ng, grid(ng)%z0_r, myfile, &
1055 & __line__, 'GRID(ng)%z0_r')) RETURN
1056
1057 IF (.not.destroy(ng, grid(ng)%z0_w, myfile, &
1058 & __line__, 'GRID(ng)%z0_w')) RETURN
1059
1060 IF (.not.destroy(ng, grid(ng)%z_r, myfile, &
1061 & __line__, 'GRID(ng)%z_r')) RETURN
1062
1063 IF (.not.destroy(ng, grid(ng)%z_v, myfile, &
1064 & __line__, 'GRID(ng)%z_v')) RETURN
1065
1066 IF (.not.destroy(ng, grid(ng)%z_w, myfile, &
1067 & __line__, 'GRID(ng)%z_w')) RETURN
1068
1069# ifdef ICESHELF
1070 IF (.not.destroy(ng, grid(ng)%IcePress, myfile, &
1071 & __line__, 'GRID(ng)%IcePress')) RETURN
1072
1073 IF (.not.destroy(ng, grid(ng)%zice, myfile, &
1074 & __line__, 'GRID(ng)%zice')) RETURN
1075# endif
1076
1077# endif
1078
1079# ifdef MASKING
1080 IF (.not.destroy(ng, grid(ng)%pmask, myfile, &
1081 & __line__, 'GRID(ng)%pmask')) RETURN
1082
1083 IF (.not.destroy(ng, grid(ng)%rmask, myfile, &
1084 & __line__, 'GRID(ng)%rmask')) RETURN
1085
1086 IF (.not.destroy(ng, grid(ng)%umask, myfile, &
1087 & __line__, 'GRID(ng)%umask')) RETURN
1088
1089 IF (.not.destroy(ng, grid(ng)%vmask, myfile, &
1090 & __line__, 'GRID(ng)%vmask')) RETURN
1091
1092# if defined AVERAGES || \
1093 (defined ad_averages && defined adjoint) || \
1094 (defined rp_averages && defined tl_ioms) || \
1095 (defined tl_averages && defined tangent)
1096 IF (.not.destroy(ng, grid(ng)%pmask_avg, myfile, &
1097 & __line__, 'GRID(ng)%pmask_avg')) RETURN
1098
1099 IF (.not.destroy(ng, grid(ng)%rmask_avg, myfile, &
1100 & __line__, 'GRID(ng)%rmask_avg')) RETURN
1101
1102 IF (.not.destroy(ng, grid(ng)%umask_avg, myfile, &
1103 & __line__, 'GRID(ng)%umask_avg')) RETURN
1104
1105 IF (.not.destroy(ng, grid(ng)%vmask_avg, myfile, &
1106 & __line__, 'GRID(ng)%vmask_avg')) RETURN
1107# endif
1108
1109# ifdef DIAGNOSTICS
1110 IF (.not.destroy(ng, grid(ng)%pmask_dia, myfile, &
1111 & __line__, 'GRID(ng)%pmask_dia')) RETURN
1112
1113 IF (.not.destroy(ng, grid(ng)%rmask_dia, myfile, &
1114 & __line__, 'GRID(ng)%rmask_dia')) RETURN
1115
1116 IF (.not.destroy(ng, grid(ng)%umask_dia, myfile, &
1117 & __line__, 'GRID(ng)%umask_dia')) RETURN
1118
1119 IF (.not.destroy(ng, grid(ng)%vmask_dia, myfile, &
1120 & __line__, 'GRID(ng)%vmask_dia')) RETURN
1121# endif
1122
1123 IF (.not.destroy(ng, grid(ng)%pmask_full, myfile, &
1124 & __line__, 'GRID(ng)%pmask_full')) RETURN
1125
1126 IF (.not.destroy(ng, grid(ng)%rmask_full, myfile, &
1127 & __line__, 'GRID(ng)%rmask_full')) RETURN
1128
1129 IF (.not.destroy(ng, grid(ng)%umask_full, myfile, &
1130 & __line__, 'GRID(ng)%umask_full')) RETURN
1131
1132 IF (.not.destroy(ng, grid(ng)%vmask_full, myfile, &
1133 & __line__, 'GRID(ng)%vmask_full')) RETURN
1134# endif
1135
1136# ifdef WET_DRY
1137 IF (.not.destroy(ng, grid(ng)%pmask_wet, myfile, &
1138 & __line__, 'GRID(ng)%pmask_wet')) RETURN
1139
1140 IF (.not.destroy(ng, grid(ng)%rmask_wet, myfile, &
1141 & __line__, 'GRID(ng)%rmask_wet')) RETURN
1142
1143 IF (.not.destroy(ng, grid(ng)%umask_wet, myfile, &
1144 & __line__, 'GRID(ng)%umask_wet')) RETURN
1145
1146 IF (.not.destroy(ng, grid(ng)%vmask_wet, myfile, &
1147 & __line__, 'GRID(ng)%vmask_wet')) RETURN
1148
1149# ifdef SOLVE3D
1150 IF (.not.destroy(ng, grid(ng)%rmask_wet_avg, myfile, &
1151 & __line__, 'GRID(ng)%rmask_wet_avg')) RETURN
1152# endif
1153# endif
1154
1155# if defined AD_SENSITIVITY || defined I4DVAR_ANA_SENSITIVITY || \
1156 defined opt_observations || defined sensitivity_4dvar || \
1157 defined so_semi
1158 IF (.not.destroy(ng, grid(ng)%Rscope, myfile, &
1159 & __line__, 'GRID(ng)%Rscope')) RETURN
1160
1161 IF (.not.destroy(ng, grid(ng)%Uscope, myfile, &
1162 & __line__, 'GRID(ng)%Uscope')) RETURN
1163
1164 IF (.not.destroy(ng, grid(ng)%Vscope, myfile, &
1165 & __line__, 'GRID(ng)%Vscope')) RETURN
1166# endif
1167
1168# if defined TANGENT || defined TL_IOMS
1169!
1170! Tangent linear model state.
1171!
1172 IF (.not.destroy(ng, grid(ng)%tl_h, myfile, &
1173 & __line__, 'GRID(ng)%tl_h')) RETURN
1174
1175# ifdef SOLVE3D
1176 IF (.not.destroy(ng, grid(ng)%tl_Hz, myfile, &
1177 & __line__, 'GRID(ng)%tl_Hz')) RETURN
1178
1179# ifdef ADJUST_BOUNDARY
1180 IF (.not.destroy(ng, grid(ng)%tl_Hz_bry, myfile, &
1181 & __line__, 'GRID(ng)%tl_Hz_bry')) RETURN
1182# endif
1183
1184 IF (.not.destroy(ng, grid(ng)%tl_Huon, myfile, &
1185 & __line__, 'GRID(ng)%tl_Huon')) RETURN
1186
1187 IF (.not.destroy(ng, grid(ng)%tl_Hvom, myfile, &
1188 & __line__, 'GRID(ng)%tl_Hvom')) RETURN
1189
1190 IF (.not.destroy(ng, grid(ng)%tl_z_r, myfile, &
1191 & __line__, 'GRID(ng)%tl_z_r')) RETURN
1192
1193 IF (.not.destroy(ng, grid(ng)%tl_z_w, myfile, &
1194 & __line__, 'GRID(ng)%tl_z_w')) RETURN
1195# endif
1196# endif
1197
1198# ifdef ADJOINT
1199!
1200! Adjoint model state.
1201!
1202 IF (.not.destroy(ng, grid(ng)%ad_h, myfile, &
1203 & __line__, 'GRID(ng)%ad_h')) RETURN
1204
1205# ifdef SOLVE3D
1206 IF (.not.destroy(ng, grid(ng)%ad_Hz, myfile, &
1207 & __line__, 'GRID(ng)%ad_Hz')) RETURN
1208
1209# ifdef ADJUST_BOUNDARY
1210 IF (.not.destroy(ng, grid(ng)%ad_Hz_bry, myfile, &
1211 & __line__, 'GRID(ng)%ad_Hz_bry')) RETURN
1212# endif
1213
1214 IF (.not.destroy(ng, grid(ng)%ad_Huon, myfile, &
1215 & __line__, 'GRID(ng)%ad_Huon')) RETURN
1216
1217 IF (.not.destroy(ng, grid(ng)%ad_Hvom, myfile, &
1218 & __line__, 'GRID(ng)%ad_Hvom')) RETURN
1219
1220 IF (.not.destroy(ng, grid(ng)%ad_z_r, myfile, &
1221 & __line__, 'GRID(ng)%ad_z_r')) RETURN
1222
1223 IF (.not.destroy(ng, grid(ng)%ad_z_w, myfile, &
1224 & __line__, 'GRID(ng)%ad_z_w')) RETURN
1225# endif
1226# endif
1227#endif
1228!
1229!-----------------------------------------------------------------------
1230! Deallocate derived-type GRID structure.
1231!-----------------------------------------------------------------------
1232!
1233 IF (ng.eq.ngrids) THEN
1234 IF (allocated(grid)) deallocate ( grid )
1235 END IF
1236!
1237 RETURN

References grid, and mod_param::ngrids.

Referenced by mod_arrays::roms_deallocate_arrays().

Here is the caller graph for this function:

◆ initialize_grid()

subroutine, public mod_grid::initialize_grid ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model )

Definition at line 1240 of file mod_grid.F.

1241!
1242!=======================================================================
1243! !
1244! This routine initialize all variables in the module using first !
1245! touch distribution policy. In shared-memory configuration, this !
1246! operation actually performs propagation of the "shared arrays" !
1247! across the cluster, unless another policy is specified to !
1248! override the default. !
1249! !
1250!=======================================================================
1251!
1252 USE mod_param
1253 USE mod_scalars
1254!
1255! Imported variable declarations.
1256!
1257 integer, intent(in) :: ng, tile, model
1258!
1259! Local variable declarations.
1260!
1261 integer :: Imin, Imax, Jmin, Jmax
1262 integer :: i, j
1263#ifdef SOLVE3D
1264 integer :: k
1265#endif
1266
1267 real(r8), parameter :: IniVal = 0.0_r8
1268 real(r8) :: IniMetricVal
1269
1270#include "set_bounds.h"
1271!
1272! Set array initialization range.
1273!
1274#ifdef DISTRIBUTE
1275 imin=bounds(ng)%LBi(tile)
1276 imax=bounds(ng)%UBi(tile)
1277 jmin=bounds(ng)%LBj(tile)
1278 jmax=bounds(ng)%UBj(tile)
1279#else
1280 IF (domain(ng)%Western_Edge(tile)) THEN
1281 imin=bounds(ng)%LBi(tile)
1282 ELSE
1283 imin=istr
1284 END IF
1285 IF (domain(ng)%Eastern_Edge(tile)) THEN
1286 imax=bounds(ng)%UBi(tile)
1287 ELSE
1288 imax=iend
1289 END IF
1290 IF (domain(ng)%Southern_Edge(tile)) THEN
1291 jmin=bounds(ng)%LBj(tile)
1292 ELSE
1293 jmin=jstr
1294 END IF
1295 IF (domain(ng)%Northern_Edge(tile)) THEN
1296 jmax=bounds(ng)%UBj(tile)
1297 ELSE
1298 jmax=jend
1299 END IF
1300#endif
1301!
1302! Set initialization value that it is special in nexting to just
1303! load contact points that have not been initialized from the
1304! regular physical grid. This is done to make sure that all these
1305! important metric values have been set-up correctly.
1306!
1307#ifdef NESTING
1308 inimetricval=spval ! very large value
1309#else
1310 inimetricval=inival
1311#endif
1312!
1313!-----------------------------------------------------------------------
1314! Initialize module variables.
1315!-----------------------------------------------------------------------
1316!
1317! Nonlinear model state.
1318!
1319 IF ((model.eq.0).or.(model.eq.inlm)) THEN
1320 DO j=jmin,jmax
1321 DO i=imin,imax
1322#if defined MASKING && defined PROPAGATOR
1323 grid(ng) % IJwaterR(i,j) = 0
1324 grid(ng) % IJwaterU(i,j) = 0
1325 grid(ng) % IJwaterV(i,j) = 0
1326#endif
1327 grid(ng) % angler(i,j) = inimetricval
1328 grid(ng) % CosAngler(i,j) = inival
1329 grid(ng) % SinAngler(i,j) = inival
1330
1331#if defined TIDE_GENERATING_FORCES
1332 grid(ng) % Cos2Lat(i,j) = inival
1333 grid(ng) % SinLat2(i,j) = inival
1334#endif
1335
1336#if defined CURVGRID && defined UV_ADV
1337 grid(ng) % dmde(i,j) = inimetricval
1338 grid(ng) % dndx(i,j) = inimetricval
1339#endif
1340 grid(ng) % f(i,j) = inimetricval
1341 grid(ng) % fomn(i,j) = inival
1342 grid(ng) % grdscl(i,j) = inival
1343
1344 grid(ng) % h(i,j) = inimetricval
1345
1346 grid(ng) % latp(i,j) = inival
1347 grid(ng) % latr(i,j) = inimetricval
1348 grid(ng) % latu(i,j) = inimetricval
1349 grid(ng) % latv(i,j) = inimetricval
1350 grid(ng) % lonp(i,j) = inival
1351 grid(ng) % lonr(i,j) = inimetricval
1352 grid(ng) % lonu(i,j) = inimetricval
1353 grid(ng) % lonv(i,j) = inimetricval
1354 grid(ng) % MyLon(i,j) = inimetricval
1355
1356 grid(ng) % omn(i,j) = inival
1357 grid(ng) % om_p(i,j) = inival
1358 grid(ng) % om_r(i,j) = inival
1359 grid(ng) % om_u(i,j) = inival
1360 grid(ng) % om_v(i,j) = inival
1361 grid(ng) % on_p(i,j) = inival
1362 grid(ng) % on_r(i,j) = inival
1363 grid(ng) % on_u(i,j) = inival
1364 grid(ng) % on_v(i,j) = inival
1365
1366 grid(ng) % pm(i,j) = inimetricval
1367 grid(ng) % pn(i,j) = inimetricval
1368
1369 grid(ng) % pmon_p(i,j) = inival
1370 grid(ng) % pmon_r(i,j) = inival
1371 grid(ng) % pmon_u(i,j) = inival
1372 grid(ng) % pmon_v(i,j) = inival
1373 grid(ng) % pnom_p(i,j) = inival
1374 grid(ng) % pnom_r(i,j) = inival
1375 grid(ng) % pnom_u(i,j) = inival
1376 grid(ng) % pnom_v(i,j) = inival
1377
1378#if defined UV_LOGDRAG || defined GLS_MIXING || \
1379 defined bbl_model || defined sediment
1380 grid(ng) % ZoBot(i,j) = zob(ng)
1381#endif
1382 grid(ng) % rdrag(i,j) = rdrg(ng)
1383#if defined UV_QDRAG
1384 grid(ng) % rdrag2(i,j) = rdrg2(ng)
1385#endif
1386
1387 grid(ng) % xp(i,j) = inival
1388 grid(ng) % xr(i,j) = inimetricval
1389 grid(ng) % xu(i,j) = inimetricval
1390 grid(ng) % xv(i,j) = inimetricval
1391 grid(ng) % yp(i,j) = inival
1392 grid(ng) % yr(i,j) = inimetricval
1393 grid(ng) % yu(i,j) = inimetricval
1394 grid(ng) % yv(i,j) = inimetricval
1395
1396#if defined ICESHELF && defined SOLVE3D
1397 grid(ng) % IcePress(i,j) = inival
1398 grid(ng) % zice(i,j) = inival
1399#endif
1400
1401#ifdef MASKING
1402 grid(ng) % pmask(i,j) = inival
1403 grid(ng) % rmask(i,j) = inimetricval
1404 grid(ng) % umask(i,j) = inimetricval
1405 grid(ng) % vmask(i,j) = inimetricval
1406# if defined AVERAGES || \
1407 (defined ad_averages && defined adjoint) || \
1408 (defined rp_averages && defined tl_ioms) || \
1409 (defined tl_averages && defined tangent)
1410 grid(ng) % pmask_avg(i,j) = inival
1411 grid(ng) % rmask_avg(i,j) = inival
1412 grid(ng) % umask_avg(i,j) = inival
1413 grid(ng) % vmask_avg(i,j) = inival
1414# endif
1415# ifdef DIAGNOSTICS
1416 grid(ng) % pmask_dia(i,j) = inival
1417 grid(ng) % rmask_dia(i,j) = inival
1418 grid(ng) % umask_dia(i,j) = inival
1419 grid(ng) % vmask_dia(i,j) = inival
1420# endif
1421 grid(ng) % pmask_full(i,j) = inival
1422 grid(ng) % rmask_full(i,j) = inival
1423 grid(ng) % umask_full(i,j) = inival
1424 grid(ng) % vmask_full(i,j) = inival
1425#endif
1426
1427#ifdef WET_DRY
1428 grid(ng) % pmask_wet(i,j) = inival
1429 grid(ng) % rmask_wet(i,j) = inival
1430 grid(ng) % umask_wet(i,j) = inival
1431 grid(ng) % vmask_wet(i,j) = inival
1432# ifdef SOLVE3D
1433 grid(ng) % rmask_wet_avg(i,j) = inival
1434# endif
1435#endif
1436
1437#if defined AD_SENSITIVITY || defined I4DVAR_ANA_SENSITIVITY || \
1438 defined opt_observations || defined sensitivity_4dvar || \
1439 defined so_semi
1440 grid(ng) % Rscope(i,j) = inival
1441 grid(ng) % Uscope(i,j) = inival
1442 grid(ng) % Vscope(i,j) = inival
1443#endif
1444 END DO
1445
1446#ifdef SOLVE3D
1447 DO k=1,n(ng)
1448 DO i=imin,imax
1449 grid(ng) % Hz(i,j,k) = inival
1450 grid(ng) % Huon(i,j,k) = inival
1451 grid(ng) % Hvom(i,j,k) = inival
1452 grid(ng) % z0_r(i,j,k) = inival
1453 grid(ng) % z_r(i,j,k) = inival
1454 grid(ng) % z_v(i,j,k) = inival
1455 END DO
1456 END DO
1457 DO k=0,n(ng)
1458 DO i=imin,imax
1459 grid(ng) % z0_w(i,j,k) = inival
1460 grid(ng) % z_w(i,j,k) = inival
1461 END DO
1462 END DO
1463#endif
1464 END DO
1465#if defined ADJUST_BOUNDARY && defined SOLVE3D
1466 grid(ng) % Hz_bry = inival
1467#endif
1468 END IF
1469
1470#if defined TANGENT || defined TL_IOMS
1471!
1472! Tangent linear model state.
1473!
1474 IF ((model.eq.0).or.(model.eq.itlm).or.(model.eq.irpm)) THEN
1475 DO j=jmin,jmax
1476 DO i=imin,imax
1477 grid(ng) % tl_h(i,j) = inival
1478 END DO
1479# ifdef SOLVE3D
1480 DO k=1,n(ng)
1481 DO i=imin,imax
1482 grid(ng) % tl_Hz(i,j,k) = inival
1483 grid(ng) % tl_Huon(i,j,k) = inival
1484 grid(ng) % tl_Hvom(i,j,k) = inival
1485 grid(ng) % tl_z_r(i,j,k) = inival
1486 END DO
1487 END DO
1488 DO k=0,n(ng)
1489 DO i=imin,imax
1490 grid(ng) % tl_z_w(i,j,k) = inival
1491 END DO
1492 END DO
1493# endif
1494 END DO
1495# if defined ADJUST_BOUNDARY && defined SOLVE3D
1496 grid(ng) % tl_Hz_bry = inival
1497# endif
1498 END IF
1499#endif
1500
1501#ifdef ADJOINT
1502!
1503! Adjoint model state.
1504!
1505 IF ((model.eq.0).or.(model.eq.iadm)) THEN
1506 DO j=jmin,jmax
1507 DO i=imin,imax
1508 grid(ng) % ad_h(i,j) = inival
1509 END DO
1510# ifdef SOLVE3D
1511 DO k=1,n(ng)
1512 DO i=imin,imax
1513 grid(ng) % ad_Hz(i,j,k) = inival
1514 grid(ng) % ad_Huon(i,j,k) = inival
1515 grid(ng) % ad_Hvom(i,j,k) = inival
1516 grid(ng) % ad_z_r(i,j,k) = inival
1517 END DO
1518 END DO
1519 DO k=0,n(ng)
1520 DO i=imin,imax
1521 grid(ng) % ad_z_w(i,j,k) = inival
1522 END DO
1523 END DO
1524# endif
1525 END DO
1526# if defined ADJUST_BOUNDARY && defined SOLVE3D
1527 grid(ng) % ad_Hz_bry = inival
1528# endif
1529 END IF
1530#endif
1531!
1532 RETURN
integer, parameter inlm
Definition mod_param.F:662
integer, parameter irpm
Definition mod_param.F:664
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer, parameter iadm
Definition mod_param.F:665
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, parameter itlm
Definition mod_param.F:663
real(dp), parameter spval
real(r8), dimension(:), allocatable rdrg
real(r8), dimension(:), allocatable zob
real(r8), dimension(:), allocatable rdrg2

References mod_param::bounds, mod_param::domain, grid, mod_param::iadm, mod_param::inlm, mod_param::irpm, mod_param::itlm, mod_param::n, mod_scalars::rdrg, mod_scalars::rdrg2, mod_scalars::spval, and mod_scalars::zob.

Referenced by ad_initial(), and mod_arrays::roms_initialize_arrays().

Here is the caller graph for this function:

Variable Documentation

◆ grid

type (t_grid), dimension(:), allocatable mod_grid::grid

Definition at line 365 of file mod_grid.F.

365 TYPE (T_GRID), allocatable :: GRID(:)

Referenced by ad_balance_mod::ad_balance(), ad_bc_2d_mod::ad_bc_u2d_tile(), ad_bc_3d_mod::ad_bc_u3d_tile(), ad_bc_2d_mod::ad_bc_v2d_tile(), ad_bc_3d_mod::ad_bc_v3d_tile(), ad_biology_mod::ad_biology(), ad_bulk_flux_mod::ad_bulk_flux(), ad_convolution_mod::ad_convolution(), ad_nesting_mod::ad_correct_tracer_tile(), ad_diag_mod::ad_diag(), dotproduct_mod::ad_dotproduct(), ad_nesting_mod::ad_fine2coarse(), ad_nesting_mod::ad_get_composite(), ad_get_data(), ad_get_idata(), ad_htobs_mod::ad_htobs(), ad_ini_fields_mod::ad_ini_fields(), ini_adjust_mod::ad_ini_perturb(), ad_ini_fields_mod::ad_ini_zeta(), ad_misfit_mod::ad_misfit(), ad_obc_adjust_mod::ad_obc2d_adjust(), ad_obc_volcons_mod::ad_obc_flux(), ad_omega_mod::ad_omega(), ad_ini_fields_mod::ad_out_fields(), ad_ini_fields_mod::ad_out_zeta(), ad_pack(), ad_pack_tile(), ad_pre_step3d_mod::ad_pre_step3d(), ad_prsgrd_mod::ad_prsgrd(), ad_nesting_mod::ad_put_composite(), ad_nesting_mod::ad_put_refine2d(), ad_nesting_mod::ad_put_refine3d(), ad_rho_eos_mod::ad_rho_eos(), ad_rhs3d_mod::ad_rhs3d(), ad_set_avg_mod::ad_set_avg_tile(), ad_set_data_tile(), ad_set_depth_mod::ad_set_depth(), ad_set_depth_mod::ad_set_depth_bry(), ad_set_massflux_mod::ad_set_massflux(), ad_set_vbc_mod::ad_set_vbc(), ad_set_vbc_mod::ad_set_vbc_tile(), dotproduct_mod::ad_statenorm(), ad_step2d_mod::ad_step2d(), ad_step3d_t_mod::ad_step3d_t(), ad_step3d_uv_mod::ad_step3d_uv(), ad_t3dbc_mod::ad_t3dbc_tile(), ad_t3dmix2_mod::ad_t3dmix2(), ad_t3dmix4_mod::ad_t3dmix4(), ad_t3drelax_mod::ad_t3drelax(), ad_u2dbc_mod::ad_u2dbc_tile(), ad_u3dbc_mod::ad_u3dbc_tile(), ad_unpack(), ad_unpack_tile(), ad_uv3dmix2_mod::ad_uv3dmix2(), ad_uv3dmix4_mod::ad_uv3dmix4(), ad_uv3drelax_mod::ad_uv3drelax(), uv_var_change_mod::ad_uv_a2c_grid(), uv_var_change_mod::ad_uv_c2a_grid(), ad_v2dbc_mod::ad_v2dbc_tile(), ad_v3dbc_mod::ad_v3dbc_tile(), ad_wrt_his_mod::ad_wrt_his_nf90(), ad_wrt_his_mod::ad_wrt_his_pio(), ad_wvelocity_mod::ad_wvelocity(), ad_nesting_mod::ad_z_weights(), ad_zetabc_mod::ad_zetabc_tile(), adfromtl_mod::adfromtl(), adsen_force_mod::adsen_force(), adsen_initial_mod::adsen_initial(), allocate_grid(), analytical_mod::ana_dqdsst(), analytical_mod::ana_drag(), analytical_mod::ana_drag_tile(), analytical_mod::ana_fsobc_tile(), ana_grid(), analytical_mod::ana_initial(), analytical_mod::ana_m2obc(), analytical_mod::ana_m2obc_tile(), analytical_mod::ana_m3obc_tile(), analytical_mod::ana_mask(), analytical_mod::ana_nlminitial_tile(), analytical_mod::ana_psource(), analytical_mod::ana_scope(), analytical_mod::ana_sediment(), analytical_mod::ana_smflux(), analytical_mod::ana_specir(), analytical_mod::ana_spinning(), analytical_mod::ana_srflux(), analytical_mod::ana_tclima_tile(), analytical_mod::ana_tobc(), analytical_mod::ana_vmix(), analytical_mod::ana_winds(), analytical_mod::ana_wtype(), analytical_mod::ana_wwave(), back_cost_mod::back_cost(), zeta_balance_mod::balance_ref(), bbl_output_mod::bbl_wrt_nf90(), bbl_output_mod::bbl_wrt_pio(), bbl_mod::bblm(), bc_2d_mod::bc_u2d_tile(), bc_3d_mod::bc_u3d_tile(), bc_4d_mod::bc_u4d_tile(), bc_2d_mod::bc_v2d_tile(), bc_3d_mod::bc_v3d_tile(), bc_4d_mod::bc_v4d_tile(), zeta_balance_mod::biconj(), biology_mod::biology(), bulk_flux_mod::bulk_flux(), cgradient_mod::cgradient(), comp_jb0_mod::comp_jb0(), nesting_mod::correct_tracer_tile(), cmeps_roms_mod::createscalarfield(), deallocate_grid(), diag_mod::diag(), extract_sta_mod::extract_sta2d(), extract_sta_mod::extract_sta3d(), nesting_mod::fine2coarse(), nesting_mod::get_composite(), get_data(), get_grid_mod::get_grid_nf90(), get_grid_mod::get_grid_pio(), get_idata(), nesting_mod::get_metrics(), get_nudgcoef_mod::get_nudgcoef_nf90(), get_nudgcoef_mod::get_nudgcoef_pio(), get_state_mod::get_state_nf90(), get_state_mod::get_state_pio(), get_wetdry_mod::get_wetdry_nf90(), get_wetdry_mod::get_wetdry_pio(), gls_corstep_mod::gls_corstep(), gls_prestep_mod::gls_prestep(), grid_coords(), ice_advect_mod::ice_advect_tile(), ice_thermo_mod::ice_thermo(), ini_adjust_mod::ini_adjust(), ini_adjust_mod::ini_adjust_tile(), ini_fields_mod::ini_fields(), ini_hmixcoef_mod::ini_hmixcoef(), ini_lanczos_mod::ini_lanczos(), ini_adjust_mod::ini_perturb(), ini_fields_mod::ini_zeta(), initialize_grid(), lmd_bkpp(), lmd_vmix_mod::lmd_finish(), lmd_skpp(), lmd_vmix_mod::lmd_vmix(), ini_adjust_mod::load_adtotl(), ini_adjust_mod::load_tltoad(), nesting_mod::mask_hweights(), metrics_mod::metrics(), my25_corstep_mod::my25_corstep(), my25_prestep_mod::my25_prestep(), nf_fread2d_mod::nf_fread2d::nf90_fread2d(), dotproduct_mod::nl_dotproduct(), normalization_mod::normalization(), obc_volcons_mod::obc_flux(), obs_write_mod::obs_operator(), obs_write_mod::obs_write_nf90(), obs_write_mod::obs_write_pio(), omega_mod::omega(), biology_floats_mod::oyster_floats_tile(), nf_fread2d_mod::nf_fread2d::pio_fread2d(), posterior_mod::posterior(), posterior_var_mod::posterior_var(), pre_step3d_mod::pre_step3d(), prsgrd_mod::prsgrd(), nesting_mod::put_composite(), nesting_mod::put_refine2d(), nesting_mod::put_refine3d(), random_ic_mod::random_ic(), set_massflux_mod::reset_massflux(), rho_eos_mod::rho_eos(), rhs3d_mod::rhs3d(), cmeps_roms_mod::roms_export(), esmf_roms_mod::roms_export(), roms_interp_create(), cmeps_roms_mod::roms_rotate(), esmf_roms_mod::roms_rotate(), cmeps_roms_mod::roms_setgridarrays(), esmf_roms_mod::roms_setgridarrays(), rp_set_depth_mod::rp_bath(), rp_biology_mod::rp_biology(), rp_bulk_flux_mod::rp_bulk_flux(), rp_diag_mod::rp_diag(), rp_get_data(), rp_get_idata(), ini_adjust_mod::rp_ini_adjust(), rp_ini_fields_mod::rp_ini_fields(), rp_ini_fields_mod::rp_ini_zeta(), rp_initial(), rp_obc_adjust_mod::rp_obc2d_adjust(), rp_obc_volcons_mod::rp_obc_flux(), rp_omega_mod::rp_omega(), rp_pre_step3d_mod::rp_pre_step3d(), rp_prsgrd_mod::rp_prsgrd(), rp_rho_eos_mod::rp_rho_eos(), rp_rhs3d_mod::rp_rhs3d(), rp_set_data_tile(), rp_set_depth_mod::rp_set_depth(), rp_set_depth_mod::rp_set_depth_bry(), rp_set_massflux_mod::rp_set_massflux(), rp_set_vbc_mod::rp_set_vbc(), rp_step2d_mod::rp_step2d(), rp_step3d_t_mod::rp_step3d_t(), rp_step3d_uv_mod::rp_step3d_uv(), rp_t3dbc_mod::rp_t3dbc_tile(), rp_t3dmix2_mod::rp_t3dmix2(), rp_t3dmix4_mod::rp_t3dmix4(), rp_t3drelax_mod::rp_t3drelax(), rp_u2dbc_mod::rp_u2dbc_tile(), rp_u3dbc_mod::rp_u3dbc_tile(), rp_uv3dmix2_mod::rp_uv3dmix2(), rp_uv3dmix4_mod::rp_uv3dmix4(), rp_uv3drelax_mod::rp_uv3drelax(), rp_v2dbc_mod::rp_v2dbc_tile(), rp_v3dbc_mod::rp_v3dbc_tile(), rp_wrt_ini_mod::rp_wrt_ini_nf90(), rp_wrt_ini_mod::rp_wrt_ini_pio(), rp_zetabc_mod::rp_zetabc_tile(), sed_bed_mod::sed_bed(), sed_bedload(), sed_fluxes_mod::sed_fluxes(), sed_settling_mod::sed_settling(), sediment_output_mod::sediment_wrt_nf90(), sediment_output_mod::sediment_wrt_pio(), sediment_output_mod::sediment_wrt_station_nf90(), sediment_output_mod::sediment_wrt_station_pio(), set_avg_mod::set_avg(), set_avg_mod::set_avg_tile(), set_data_tile(), set_depth_mod::set_depth(), set_depth_mod::set_depth0(), set_depth_mod::set_depth_bry(), set_avg_mod::set_detide_tile(), set_diags_tile(), set_masks_mod::set_masks(), set_massflux_mod::set_massflux(), set_tides_mod::set_tides(), set_vbc_mod::set_vbc(), set_vbc_mod::set_vbc_tile(), state_join_mod::state_join_nf90(), state_join_mod::state_join_pio(), stats_modobs_mod::stats_modobs_nf90(), stats_modobs_mod::stats_modobs_pio(), step2d_mod::step2d(), step3d_t_mod::step3d_t(), step3d_uv_mod::step3d_uv(), step_floats_mod::step_floats_tile(), stiffness_mod::stiffness(), t3dbc_mod::t3dbc_tile(), t3dmix2_mod::t3dmix2(), t3dmix4_mod::t3dmix4(), time_corr_mod::time_corr_nf90(), time_corr_mod::time_corr_pio(), tkebc_mod::tkebc_tile(), tl_balance_mod::tl_balance(), tl_set_depth_mod::tl_bath(), tl_biology_mod::tl_biology(), tl_bulk_flux_mod::tl_bulk_flux(), tl_convolution_mod::tl_convolution(), tl_nesting_mod::tl_correct_tracer_tile(), tl_diag_mod::tl_diag(), dotproduct_mod::tl_dotproduct(), tl_nesting_mod::tl_fine2coarse(), tl_nesting_mod::tl_get_composite(), tl_get_data(), tl_get_idata(), tl_ini_fields_mod::tl_ini_fields(), ini_adjust_mod::tl_ini_perturb(), tl_ini_fields_mod::tl_ini_zeta(), inner2state_mod::tl_inner2state(), tl_obc_adjust_mod::tl_obc2d_adjust(), tl_obc_volcons_mod::tl_obc_flux(), tl_omega_mod::tl_omega(), tl_pack(), tl_pre_step3d_mod::tl_pre_step3d(), tl_prsgrd_mod::tl_prsgrd(), tl_nesting_mod::tl_put_composite(), tl_nesting_mod::tl_put_refine2d(), tl_nesting_mod::tl_put_refine3d(), tl_rho_eos_mod::tl_rho_eos(), tl_rhs3d_mod::tl_rhs3d(), tl_set_avg_mod::tl_set_avg_tile(), tl_set_data_tile(), tl_set_depth_mod::tl_set_depth(), tl_set_depth_mod::tl_set_depth_bry(), tl_set_massflux_mod::tl_set_massflux(), tl_set_vbc_mod::tl_set_vbc(), tl_set_vbc_mod::tl_set_vbc_tile(), dotproduct_mod::tl_statenorm(), tl_step2d_mod::tl_step2d(), tl_step3d_t_mod::tl_step3d_t(), tl_step3d_uv_mod::tl_step3d_uv(), tl_t3dbc_mod::tl_t3dbc_tile(), tl_t3dmix2_mod::tl_t3dmix2(), tl_t3dmix4_mod::tl_t3dmix4(), tl_t3drelax_mod::tl_t3drelax(), tl_u2dbc_mod::tl_u2dbc_tile(), tl_u3dbc_mod::tl_u3dbc_tile(), tl_unpack(), tl_unpack_tile(), tl_uv3dmix2_mod::tl_uv3dmix2(), tl_uv3dmix4_mod::tl_uv3dmix4(), tl_uv3drelax_mod::tl_uv3drelax(), uv_var_change_mod::tl_uv_a2c_grid(), uv_var_change_mod::tl_uv_c2a_grid(), tl_v2dbc_mod::tl_v2dbc_tile(), tl_v3dbc_mod::tl_v3dbc_tile(), tl_wrt_his_mod::tl_wrt_his_nf90(), tl_wrt_his_mod::tl_wrt_his_pio(), tl_wrt_ini_mod::tl_wrt_ini_nf90(), tl_wrt_ini_mod::tl_wrt_ini_pio(), tl_nesting_mod::tl_z_weights(), tl_zetabc_mod::tl_zetabc_tile(), u2dbc_mod::u2dbc_tile(), u3dbc_mod::u3dbc_tile(), uv3dmix2_mod::uv3dmix2(), uv3dmix4_mod::uv3dmix4(), uv_var_change_mod::uv_a2c_grid(), uv_var_change_mod::uv_c2a_grid(), v2dbc_mod::v2dbc_tile(), v3dbc_mod::v3dbc_tile(), vorticity_mod::vorticity(), vwalk_floats_mod::vwalk_floats_tile(), wetdry_mod::wetdry(), wpoints(), wrt_aug_imp_mod::wrt_aug_imp_nf90(), wrt_aug_imp_mod::wrt_aug_imp_pio(), wrt_avg_mod::wrt_avg_nf90(), wrt_avg_mod::wrt_avg_pio(), wrt_dai_mod::wrt_dai_nf90(), wrt_dai_mod::wrt_dai_pio(), wrt_diags_mod::wrt_diags_nf90(), wrt_diags_mod::wrt_diags_pio(), wrt_error_mod::wrt_error_nf90(), wrt_error_mod::wrt_error_pio(), wrt_evolved_mod::wrt_evolved_nf90(), wrt_evolved_mod::wrt_evolved_pio(), wrt_ini_mod::wrt_frc_ad_nf90(), wrt_ini_mod::wrt_frc_ad_pio(), wrt_ini_mod::wrt_frc_nf90(), wrt_ini_mod::wrt_frc_pio(), wrt_hessian_mod::wrt_hessian_nf90(), wrt_hessian_mod::wrt_hessian_pio(), wrt_his_mod::wrt_his_nf90(), wrt_his_mod::wrt_his_pio(), wrt_impulse_mod::wrt_impulse_nf90(), wrt_impulse_mod::wrt_impulse_pio(), wrt_info_mod::wrt_info::wrt_info_nf90(), wrt_info_mod::wrt_info::wrt_info_pio(), wrt_ini_mod::wrt_ini_nf90(), wrt_ini_mod::wrt_ini_pio(), wrt_quick_mod::wrt_quick_nf90(), wrt_quick_mod::wrt_quick_pio(), wrt_rst_mod::wrt_rst_nf90(), wrt_rst_mod::wrt_rst_pio(), wrt_state_mod::wrt_state_nf90(), wrt_state_mod::wrt_state_pio(), wrt_station_mod::wrt_station_nf90(), wrt_station_mod::wrt_station_pio(), wrt_tides_mod::wrt_tides_nf90(), wrt_tides_mod::wrt_tides_pio(), wvelocity_mod::wvelocity(), nesting_mod::z_weights(), and zetabc_mod::zetabc_tile().