1 !***********************************************************************
2 !* GNU Lesser General Public License
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
6 !* FMS is free software: you can redistribute it and/or modify it under
7 !* the terms of the GNU Lesser General Public License as published by
8 !* the Free Software Foundation, either version 3 of the License, or (at
9 !* your option) any later version.
11 !* FMS is distributed in the hope that it will be useful, but WITHOUT
12 !* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 !* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
16 !* You should have received a copy of the GNU Lesser General Public
17 !* License along with FMS. If not, see <http://www.gnu.org/licenses/>.
18 !***********************************************************************
19 !> @defgroup gradient_mod gradient_mod
21 !> @brief Implements some utility routines to calculate gradient.
24 !! Currently only gradient on cubic grid is implemented. Also a public interface
25 !! is provided to calculate grid information needed to calculate gradient.
27 !> @addtogroup gradient_mod
31 use mpp_mod, only : mpp_error, FATAL
32 use constants_mod, only : RADIUS
39 public :: gradient_cubic
40 public :: calc_cubic_grid_info
42 ! Include variable "version" to be written to log file.
43 #include<file_version.h>
48 !#####################################################################
49 !> Pin has halo size = 1.
50 !! @param pin the size of pin will be (nx+2,ny+2), T-cell center, with halo = 1
51 !! @param the size of dx will be (nx, ny+1), N-cell center
52 !! @param the size of dy will be (nx+1, ny), E-cell center
53 !! @param the size of area will be (nx, ny), T-cell center.
54 !! @param The size of edge_w will be (ny+1), C-cell center
55 !! @param The size of edge_e will be (ny+1), C-cell center
56 !! @param The size of edge_s will be (nx+1), C-cell center
57 !! @param The size of edge_n will be (nx+1), C-cell center
58 !! @param The size of en_n will be (3,nx,ny+1), N-cell center
59 !! @param The size of en_e will be (3,nx+1,ny), E-cell center
60 !! @param The size of vlon will be (3,nx, ny) T-cell center
61 !! @param The size of vlat will be (3,nx, ny), T-cell center
62 subroutine gradient_cubic(pin, dx, dy, area, edge_w, edge_e, edge_s, edge_n, &
63 en_n, en_e, vlon, vlat, grad_x, grad_y, on_west_edge, &
64 on_east_edge, on_south_edge, on_north_edge)
66 real(r8_kind), dimension(:,: ), intent(in ) :: pin, dx, dy, area
67 real(r8_kind), dimension(: ), intent(in ) :: edge_w, edge_e, edge_s, edge_n
68 real(r8_kind), dimension(:,:,:), intent(in ) :: en_n, en_e
69 real(r8_kind), dimension(:,:,:), intent(in ) :: vlon, vlat
70 real(r8_kind), dimension(:,: ), intent(out) :: grad_x, grad_y
71 logical, intent(in ) :: on_west_edge, on_east_edge, on_south_edge, on_north_edge
78 if(size(pin,1) .NE. nx+2 .OR. size(pin,2) .NE. ny+2)call mpp_error(FATAL, &
79 & "gradient_mod:size of pin should be (nx+2, ny+2)")
80 if(size(dx,1) .NE. nx .OR. size(dx,2) .NE. ny+1 ) call mpp_error(FATAL, &
81 & "gradient_mod: size of dx should be (nx,ny+1)")
82 if(size(dy,1) .NE. nx+1 .OR. size(dy,2) .NE. ny ) call mpp_error(FATAL, &
83 & "gradient_mod: size of dy should be (nx+1,ny)")
84 if(size(area,1) .NE. nx .OR. size(area,2) .NE. ny ) call mpp_error(FATAL, &
85 & "gradient_mod: size of area should be (nx,ny)")
86 if(size(vlon,1) .NE. 3 .OR. size(vlon,2) .NE. nx .OR. size(vlon,3) .NE. ny) &
87 call mpp_error(FATAL, "gradient_mod: size of vlon should be (3,nx,ny)")
88 if(size(vlat,1) .NE. 3 .OR. size(vlat,2) .NE. nx .OR. size(vlat,3) .NE. ny) &
89 call mpp_error(FATAL, "gradient_mod: size of vlat should be (3,nx,ny)")
90 if(size(edge_w) .NE. ny+1) call mpp_error(FATAL, "gradient_mod: size of edge_w should be (ny+1)")
91 if(size(edge_e) .NE. ny+1) call mpp_error(FATAL, "gradient_mod: size of edge_e should be (ny+1)")
92 if(size(edge_s) .NE. nx+1) call mpp_error(FATAL, "gradient_mod: size of edge_s should be (nx+1)")
93 if(size(edge_n) .NE. nx+1) call mpp_error(FATAL, "gradient_mod: size of edge_n should be (nx+1)")
94 if(size(en_n,1) .NE. 3 .OR. size(en_n,2) .NE. nx .OR. size(en_n,3) .NE. ny+1 ) &
95 call mpp_error(FATAL, "gradient_mod:size of en_n should be (3, nx, ny+1)")
96 if(size(en_e,1) .NE. 3 .OR. size(en_e,2) .NE. nx+1 .OR. size(en_e,3) .NE. ny ) &
97 call mpp_error(FATAL, "gradient_mod:size of en_e should be (3, nx+1, ny)")
99 call grad_c2l(nx, ny, pin, dx, dy, area, edge_w, edge_e, edge_s, edge_n, en_n, en_e, vlon, vlat, &
100 grad_x, grad_y, on_west_edge, on_east_edge, on_south_edge, on_north_edge)
104 end subroutine gradient_cubic
107 subroutine calc_cubic_grid_info(xt, yt, xc, yc, dx, dy, area, edge_w, edge_e, edge_s, edge_n, &
108 en_n, en_e, vlon, vlat, on_west_edge, on_east_edge, on_south_edge, on_north_edge )
109 real(r8_kind), dimension(:,: ), intent(in ) :: xt, yt, xc, yc
110 real(r8_kind), dimension(:,: ), intent(out) :: dx, dy, area
111 real(r8_kind), dimension(: ), intent(out) :: edge_w, edge_e, edge_s, edge_n
112 real(r8_kind), dimension(:,:,:), intent(out) :: en_n, en_e
113 real(r8_kind), dimension(:,:,:), intent(out) :: vlon, vlat
114 logical, intent(in ) :: on_west_edge, on_east_edge, on_south_edge, on_north_edge
115 integer :: nx, ny, nxp, nyp
123 if(size(xt,1) .NE. nx+2 .OR. size(xt,2) .NE. ny+2 ) call mpp_error(FATAL, &
124 & "gradient_mod: size of xt should be (nx+2,ny+2)")
125 if(size(yt,1) .NE. nx+2 .OR. size(yt,2) .NE. ny+2 ) call mpp_error(FATAL, &
126 & "gradient_mod: size of yt should be (nx+2,ny+2)")
127 if(size(xc,1) .NE. nxp .OR. size(xc,2) .NE. nyp ) call mpp_error(FATAL, &
128 & "gradient_mod: size of xc should be (nx+1,ny+1)")
129 if(size(yc,1) .NE. nxp .OR. size(yc,2) .NE. nyp ) call mpp_error(FATAL, &
130 & "gradient_mod: size of yc should be (nx+1,ny+1)")
131 if(size(dx,1) .NE. nx .OR. size(dx,2) .NE. nyp ) call mpp_error(FATAL, &
132 & "gradient_mod: size of dx should be (nx,ny+1)")
133 if(size(dy,1) .NE. nxp .OR. size(dy,2) .NE. ny ) call mpp_error(FATAL, &
134 & "gradient_mod: size of dy should be (nx+1,ny)")
135 if(size(area,1) .NE. nx .OR. size(area,2) .NE. ny ) call mpp_error(FATAL, &
136 & "gradient_mod: size of area should be (nx,ny)")
137 if(size(vlon,1) .NE. 3 .OR. size(vlon,2) .NE. nx .OR. size(vlon,3) .NE. ny) &
138 call mpp_error(FATAL, "gradient_mod: size of vlon should be (3,nx,ny)")
139 if(size(vlat,1) .NE. 3 .OR. size(vlat,2) .NE. nx .OR. size(vlat,3) .NE. ny) &
140 call mpp_error(FATAL, "gradient_mod: size of vlat should be (3,nx,ny)")
141 if(size(edge_w) .NE. ny+1) call mpp_error(FATAL, "gradient_mod: size of edge_w should be (ny-1)")
142 if(size(edge_e) .NE. ny+1) call mpp_error(FATAL, "gradient_mod: size of edge_e should be (ny-1)")
143 if(size(edge_s) .NE. nx+1) call mpp_error(FATAL, "gradient_mod: size of edge_s should be (nx-1)")
144 if(size(edge_n) .NE. nx+1) call mpp_error(FATAL, "gradient_mod: size of edge_n should be (nx-1)")
145 if(size(en_n,1) .NE. 3 .OR. size(en_n,2) .NE. nx .OR. size(en_n,3) .NE. nyp ) &
146 call mpp_error(FATAL, "gradient_mod:size of en_n should be (3, nx, ny+1)")
147 if(size(en_e,1) .NE. 3 .OR. size(en_e,2) .NE. nxp .OR. size(en_e,3) .NE. ny ) &
148 call mpp_error(FATAL, "gradient_mod:size of en_e should be (3, nx+1, ny)")
151 call calc_c2l_grid_info(nx, ny, xt, yt, xc, yc, dx, dy, area, edge_w, edge_e, edge_s, edge_n, &
152 en_n, en_e, vlon, vlat, on_west_edge, on_east_edge, on_south_edge, on_north_edge )
157 end subroutine calc_cubic_grid_info
159 end module gradient_mod
161 ! close documentation grouping