New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
modvariables.F90 in branches/UKMO/r8727_WAVE-2_Clementi_add_coupling/NEMOGCM/EXTERNAL/AGRIF/AGRIF_FILES – NEMO

source: branches/UKMO/r8727_WAVE-2_Clementi_add_coupling/NEMOGCM/EXTERNAL/AGRIF/AGRIF_FILES/modvariables.F90 @ 8749

Last change on this file since 8749 was 8749, checked in by jcastill, 6 years ago

Remove svn keywords

File size: 5.8 KB
Line 
1module Agrif_Variables
2!
3    use Agrif_CurgridFunctions
4!
5    implicit none
6!
7contains
8!
9!===================================================================================================
10!  subroutine Agrif_Declare_Variable
11!
12!> Declare a new variable profile
13!---------------------------------------------------------------------------------------------------
14subroutine Agrif_Declare_Variable ( posvar, firstpoint, raf, lb, ub, varid, torestore )
15!---------------------------------------------------------------------------------------------------
16    integer,      dimension(:), intent(in)  :: posvar     !< position of the variable on the cell
17                                                          !! (1 for the border of the edge, 2 for the center)
18    integer,      dimension(:), intent(in)  :: firstpoint !< index of the first point in the real domain
19    character(1), dimension(:), intent(in)  :: raf        !< Array indicating the type of dimension (space or not)
20                                                          !!   for each of them
21    integer,      dimension(:), intent(in)  :: lb         !< Lower bounds of the array
22    integer,      dimension(:), intent(in)  :: ub         !< Upper bounds of the array
23    integer,                    intent(out) :: varid      !< Id number of the newly created variable
24    logical,      optional,     intent(in)  :: torestore  !< Indicates if the array restore is used
25!---------------------------------------------------------------------------------------------------
26    type(Agrif_Variables_List), pointer :: new_varlist
27    type(Agrif_Variable),       pointer :: var
28    integer                             :: nbdim, i
29    logical                             :: restore
30
31    restore = .FALSE.
32    if ( Agrif_Mygrid % ngridstep /= 0 ) then
33        if (present(torestore)) restore = torestore
34    endif
35!
36    nbdim = SIZE(posvar)
37!
38    allocate(new_varlist)
39    allocate(new_varlist % var)
40
41    var => new_varlist % var
42
43    allocate(var % posvar(nbdim))
44    allocate(var % interptab(nbdim))
45    allocate(var % coords(nbdim))
46!
47    var % nbdim          = nbdim
48    var % interptab      = raf(1:nbdim)
49    var % posvar         = posvar(1:nbdim)
50    var % point(1:nbdim) = firstpoint(1:nbdim)
51    var % restore = restore
52!
53    do i = 1,nbdim
54        select case( raf(i) )
55            case('x')    ; var % coords(i) = 1
56            case('y')    ; var % coords(i) = 2
57            case('z')    ; var % coords(i) = 3
58            case('N')    ; var % coords(i) = 0
59            case default ; var % coords(i) = 0
60        end select
61    enddo
62!
63    var % lb(1:nbdim) = lb(1:nbdim)
64    var % ub(1:nbdim) = ub(1:nbdim)
65
66    if ( restore ) then
67        select case(nbdim)
68        case(1)
69            allocate(var % Restore1D(lb(1):ub(1)))
70            var % Restore1D = 0
71        case(2)
72            allocate(var % Restore2D(lb(1):ub(1),   &
73                                     lb(2):ub(2)))
74            var % Restore2D = 0
75        case(3)
76            allocate(var % Restore3D(lb(1):ub(1),   &
77                                     lb(2):ub(2),   &
78                                     lb(3):ub(3)))
79            var % Restore3D = 0
80        case(4)
81            allocate(var % Restore4D(lb(1):ub(1),   &
82                                     lb(2):ub(2),   &
83                                     lb(3):ub(3),   &
84                                     lb(4):ub(4)))
85            var % Restore4D = 0
86        case(5)
87            allocate(var % Restore5D(lb(1):ub(1),   &
88                                     lb(2):ub(2),   &
89                                     lb(3):ub(3),   &
90                                     lb(4):ub(4),   &
91                                     lb(5):ub(5)))
92            var % Restore5D = 0
93        end select
94    endif
95
96    new_varlist % next => Agrif_Curgrid % variables
97
98    Agrif_Curgrid % variables => new_varlist
99    Agrif_Curgrid % Nbvariables = Agrif_Curgrid % Nbvariables + 1
100
101    varid = -Agrif_Curgrid % Nbvariables
102
103    var % parent_var => Agrif_Search_Variable(Agrif_Curgrid % parent, Agrif_Curgrid % nbvariables)
104    var % root_var   => Agrif_Search_Variable(Agrif_Mygrid, Agrif_Curgrid % nbvariables)
105!---------------------------------------------------------------------------------------------------
106end subroutine Agrif_Declare_Variable
107!===================================================================================================
108!
109!===================================================================================================
110!  function Agrif_Search_Variable
111!
112!> Returns a pointer to the variable varid for the grid grid.
113!---------------------------------------------------------------------------------------------------
114function Agrif_Search_Variable ( grid, varid ) result(outvar)
115!---------------------------------------------------------------------------------------------------
116    type(Agrif_Grid), pointer       :: grid     !< Pointer on the current grid.
117    integer,          intent(in)    :: varid    !< ID number of the variable we are looking for.
118!
119    type(Agrif_Variable),       pointer :: outvar
120    type(Agrif_Variables_List), pointer :: parcours
121    integer :: nb, varidinv
122!
123    if ( .not.associated(grid) ) then
124        outvar => NULL()
125        return
126    endif
127!
128    parcours => grid % variables
129
130    if (.not. associated(parcours)) then   ! can occur on the grand mother grid
131           outvar => NULL()                ! during the first call by agrif_mygrid
132           return
133    endif
134
135    varidinv = 1 + grid % nbvariables - varid
136
137    do nb = 1,varidinv-1
138        parcours => parcours % next
139    enddo
140
141    outvar => parcours % var
142!---------------------------------------------------------------------------------------------------
143end function Agrif_Search_variable
144!===================================================================================================
145!
146end module Agrif_Variables
Note: See TracBrowser for help on using the repository browser.