FMS  2025.04
Flexible Modeling System
gradient.F90
1 !***********************************************************************
2 !* Apache License 2.0
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* Licensed under the Apache License, Version 2.0 (the "License");
7 !* you may not use this file except in compliance with the License.
8 !* You may obtain a copy of the License at
9 !*
10 !* http://www.apache.org/licenses/LICENSE-2.0
11 !*
12 !* FMS is distributed in the hope that it will be useful, but WITHOUT
13 !* WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied;
14 !* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
15 !* PARTICULAR PURPOSE. See the License for the specific language
16 !* governing permissions and limitations under the License.
17 !***********************************************************************
18 !> @defgroup gradient_mod gradient_mod
19 !> @ingroup mosaic
20 !> @brief Implements some utility routines to calculate gradient.
21 !> @author Zhi Liang
22 !!
23 !! Currently only gradient on cubic grid is implemented. Also a public interface
24 !! is provided to calculate grid information needed to calculate gradient.
25 
26 !> @addtogroup gradient_mod
27 !> @{
28 module gradient_mod
29 
30 use mpp_mod, only : mpp_error, fatal
31 use constants_mod, only : radius
32 use platform_mod
33 
34 implicit none
35 private
36 
37 
38 public :: gradient_cubic
39 public :: calc_cubic_grid_info
40 
41 ! Include variable "version" to be written to log file.
42 #include<file_version.h>
43 
44 contains
45 
46 
47 !#####################################################################
48 !> Pin has halo size = 1.
49 !! @param pin the size of pin will be (nx+2,ny+2), T-cell center, with halo = 1
50 !! @param the size of dx will be (nx, ny+1), N-cell center
51 !! @param the size of dy will be (nx+1, ny), E-cell center
52 !! @param the size of area will be (nx, ny), T-cell center.
53 !! @param The size of edge_w will be (ny+1), C-cell center
54 !! @param The size of edge_e will be (ny+1), C-cell center
55 !! @param The size of edge_s will be (nx+1), C-cell center
56 !! @param The size of edge_n will be (nx+1), C-cell center
57 !! @param The size of en_n will be (3,nx,ny+1), N-cell center
58 !! @param The size of en_e will be (3,nx+1,ny), E-cell center
59 !! @param The size of vlon will be (3,nx, ny) T-cell center
60 !! @param The size of vlat will be (3,nx, ny), T-cell center
61 subroutine gradient_cubic(pin, dx, dy, area, edge_w, edge_e, edge_s, edge_n, &
62  en_n, en_e, vlon, vlat, grad_x, grad_y, on_west_edge, &
63  on_east_edge, on_south_edge, on_north_edge)
64 
65  real(r8_kind), dimension(:,: ), intent(in ) :: pin, dx, dy, area
66  real(r8_kind), dimension(: ), intent(in ) :: edge_w, edge_e, edge_s, edge_n
67  real(r8_kind), dimension(:,:,:), intent(in ) :: en_n, en_e
68  real(r8_kind), dimension(:,:,:), intent(in ) :: vlon, vlat
69  real(r8_kind), dimension(:,: ), intent(out) :: grad_x, grad_y
70  logical, intent(in ) :: on_west_edge, on_east_edge, on_south_edge, on_north_edge
71  integer :: nx, ny
72 
73 
74  nx = size(grad_x,1)
75  ny = size(grad_x,2)
76 
77  if(size(pin,1) .NE. nx+2 .OR. size(pin,2) .NE. ny+2)call mpp_error(fatal, &
78  & "gradient_mod:size of pin should be (nx+2, ny+2)")
79  if(size(dx,1) .NE. nx .OR. size(dx,2) .NE. ny+1 ) call mpp_error(fatal, &
80  & "gradient_mod: size of dx should be (nx,ny+1)")
81  if(size(dy,1) .NE. nx+1 .OR. size(dy,2) .NE. ny ) call mpp_error(fatal, &
82  & "gradient_mod: size of dy should be (nx+1,ny)")
83  if(size(area,1) .NE. nx .OR. size(area,2) .NE. ny ) call mpp_error(fatal, &
84  & "gradient_mod: size of area should be (nx,ny)")
85  if(size(vlon,1) .NE. 3 .OR. size(vlon,2) .NE. nx .OR. size(vlon,3) .NE. ny) &
86  call mpp_error(fatal, "gradient_mod: size of vlon should be (3,nx,ny)")
87  if(size(vlat,1) .NE. 3 .OR. size(vlat,2) .NE. nx .OR. size(vlat,3) .NE. ny) &
88  call mpp_error(fatal, "gradient_mod: size of vlat should be (3,nx,ny)")
89  if(size(edge_w) .NE. ny+1) call mpp_error(fatal, "gradient_mod: size of edge_w should be (ny+1)")
90  if(size(edge_e) .NE. ny+1) call mpp_error(fatal, "gradient_mod: size of edge_e should be (ny+1)")
91  if(size(edge_s) .NE. nx+1) call mpp_error(fatal, "gradient_mod: size of edge_s should be (nx+1)")
92  if(size(edge_n) .NE. nx+1) call mpp_error(fatal, "gradient_mod: size of edge_n should be (nx+1)")
93  if(size(en_n,1) .NE. 3 .OR. size(en_n,2) .NE. nx .OR. size(en_n,3) .NE. ny+1 ) &
94  call mpp_error(fatal, "gradient_mod:size of en_n should be (3, nx, ny+1)")
95  if(size(en_e,1) .NE. 3 .OR. size(en_e,2) .NE. nx+1 .OR. size(en_e,3) .NE. ny ) &
96  call mpp_error(fatal, "gradient_mod:size of en_e should be (3, nx+1, ny)")
97 
98  call grad_c2l(nx, ny, pin, dx, dy, area, edge_w, edge_e, edge_s, edge_n, en_n, en_e, vlon, vlat, &
99  grad_x, grad_y, on_west_edge, on_east_edge, on_south_edge, on_north_edge)
100 
101  return
102 
103 end subroutine gradient_cubic
104 
105 
106 subroutine calc_cubic_grid_info(xt, yt, xc, yc, dx, dy, area, edge_w, edge_e, edge_s, edge_n, &
107  en_n, en_e, vlon, vlat, on_west_edge, on_east_edge, on_south_edge, on_north_edge )
108  real(r8_kind), dimension(:,: ), intent(in ) :: xt, yt, xc, yc
109  real(r8_kind), dimension(:,: ), intent(out) :: dx, dy, area
110  real(r8_kind), dimension(: ), intent(out) :: edge_w, edge_e, edge_s, edge_n
111  real(r8_kind), dimension(:,:,:), intent(out) :: en_n, en_e
112  real(r8_kind), dimension(:,:,:), intent(out) :: vlon, vlat
113  logical, intent(in ) :: on_west_edge, on_east_edge, on_south_edge, on_north_edge
114  integer :: nx, ny, nxp, nyp
115 
116 
117  nx = size(area,1)
118  ny = size(area,2)
119  nxp = nx+1
120  nyp = ny+1
121 
122  if(size(xt,1) .NE. nx+2 .OR. size(xt,2) .NE. ny+2 ) call mpp_error(fatal, &
123  & "gradient_mod: size of xt should be (nx+2,ny+2)")
124  if(size(yt,1) .NE. nx+2 .OR. size(yt,2) .NE. ny+2 ) call mpp_error(fatal, &
125  & "gradient_mod: size of yt should be (nx+2,ny+2)")
126  if(size(xc,1) .NE. nxp .OR. size(xc,2) .NE. nyp ) call mpp_error(fatal, &
127  & "gradient_mod: size of xc should be (nx+1,ny+1)")
128  if(size(yc,1) .NE. nxp .OR. size(yc,2) .NE. nyp ) call mpp_error(fatal, &
129  & "gradient_mod: size of yc should be (nx+1,ny+1)")
130  if(size(dx,1) .NE. nx .OR. size(dx,2) .NE. nyp ) call mpp_error(fatal, &
131  & "gradient_mod: size of dx should be (nx,ny+1)")
132  if(size(dy,1) .NE. nxp .OR. size(dy,2) .NE. ny ) call mpp_error(fatal, &
133  & "gradient_mod: size of dy should be (nx+1,ny)")
134  if(size(area,1) .NE. nx .OR. size(area,2) .NE. ny ) call mpp_error(fatal, &
135  & "gradient_mod: size of area should be (nx,ny)")
136  if(size(vlon,1) .NE. 3 .OR. size(vlon,2) .NE. nx .OR. size(vlon,3) .NE. ny) &
137  call mpp_error(fatal, "gradient_mod: size of vlon should be (3,nx,ny)")
138  if(size(vlat,1) .NE. 3 .OR. size(vlat,2) .NE. nx .OR. size(vlat,3) .NE. ny) &
139  call mpp_error(fatal, "gradient_mod: size of vlat should be (3,nx,ny)")
140  if(size(edge_w) .NE. ny+1) call mpp_error(fatal, "gradient_mod: size of edge_w should be (ny-1)")
141  if(size(edge_e) .NE. ny+1) call mpp_error(fatal, "gradient_mod: size of edge_e should be (ny-1)")
142  if(size(edge_s) .NE. nx+1) call mpp_error(fatal, "gradient_mod: size of edge_s should be (nx-1)")
143  if(size(edge_n) .NE. nx+1) call mpp_error(fatal, "gradient_mod: size of edge_n should be (nx-1)")
144  if(size(en_n,1) .NE. 3 .OR. size(en_n,2) .NE. nx .OR. size(en_n,3) .NE. nyp ) &
145  call mpp_error(fatal, "gradient_mod:size of en_n should be (3, nx, ny+1)")
146  if(size(en_e,1) .NE. 3 .OR. size(en_e,2) .NE. nxp .OR. size(en_e,3) .NE. ny ) &
147  call mpp_error(fatal, "gradient_mod:size of en_e should be (3, nx+1, ny)")
148 
149 
150  call calc_c2l_grid_info(nx, ny, xt, yt, xc, yc, dx, dy, area, edge_w, edge_e, edge_s, edge_n, &
151  en_n, en_e, vlon, vlat, on_west_edge, on_east_edge, on_south_edge, on_north_edge )
152 
153 
154  return
155 
156 end subroutine calc_cubic_grid_info
157 
158 end module gradient_mod
159 !> @}
160 ! close documentation grouping
mpp_mod::mpp_error
Error handler.
Definition: mpp.F90:381