MOM6
MOM_transcribe_grid.F90
Go to the documentation of this file.
1 !> Module with routines for copying information from a shared dynamic horizontal
2 !! grid to an ocean-specific horizontal grid and the reverse.
4 
5 ! This file is part of MOM6. See LICENSE.md for the license.
6 
8 use mom_domains, only : to_all, scalar_pair, cgrid_ne, agrid, bgrid_ne, corner
10 use mom_error_handler, only : mom_error, mom_mesg, fatal, warning
13 
14 implicit none ; private
15 
17 
18 contains
19 
20 !> Copies information from a dynamic (shared) horizontal grid type into an
21 !! ocean_grid_type.
22 subroutine copy_dyngrid_to_mom_grid(dG, oG, US)
23  type(dyn_horgrid_type), intent(in) :: dg !< Common horizontal grid type
24  type(ocean_grid_type), intent(inout) :: og !< Ocean grid type
25  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
26 
27  integer :: isd, ied, jsd, jed ! Common data domains.
28  integer :: isdb, iedb, jsdb, jedb ! Common data domains.
29  integer :: ido, jdo, ido2, jdo2 ! Indexing offsets between the grids.
30  integer :: igst, jgst ! Global starting indices.
31  integer :: i, j
32 
33  ! MOM_grid_init and create_dyn_horgrid are called outside of this routine.
34  ! This routine copies over the fields that were set by MOM_initialized_fixed.
35 
36  ! Determine the indexing offsets between the grids.
37  ido = dg%idg_offset - og%idg_offset
38  jdo = dg%jdg_offset - og%jdg_offset
39 
40  isd = max(og%isd, dg%isd+ido) ; jsd = max(og%jsd, dg%jsd+jdo)
41  ied = min(og%ied, dg%ied+ido) ; jed = min(og%jed, dg%jed+jdo)
42  isdb = max(og%IsdB, dg%IsdB+ido) ; jsdb = max(og%JsdB, dg%JsdB+jdo)
43  iedb = min(og%IedB, dg%IedB+ido) ; jedb = min(og%JedB, dg%JedB+jdo)
44 
45  ! Check that the grids conform.
46  if ((isd > og%isc) .or. (ied < og%ied) .or. (jsd > og%jsc) .or. (jed > og%jed)) &
47  call mom_error(fatal, "copy_dyngrid_to_MOM_grid called with incompatible grids.")
48 
49  do i=isd,ied ; do j=jsd,jed
50  og%geoLonT(i,j) = dg%geoLonT(i+ido,j+jdo)
51  og%geoLatT(i,j) = dg%geoLatT(i+ido,j+jdo)
52  og%dxT(i,j) = dg%dxT(i+ido,j+jdo)
53  og%dyT(i,j) = dg%dyT(i+ido,j+jdo)
54  og%areaT(i,j) = dg%areaT(i+ido,j+jdo)
55  og%bathyT(i,j) = dg%bathyT(i+ido,j+jdo)
56 
57  og%dF_dx(i,j) = dg%dF_dx(i+ido,j+jdo)
58  og%dF_dy(i,j) = dg%dF_dy(i+ido,j+jdo)
59  og%sin_rot(i,j) = dg%sin_rot(i+ido,j+jdo)
60  og%cos_rot(i,j) = dg%cos_rot(i+ido,j+jdo)
61  og%mask2dT(i,j) = dg%mask2dT(i+ido,j+jdo)
62  enddo ; enddo
63 
64  do i=isdb,iedb ; do j=jsd,jed
65  og%geoLonCu(i,j) = dg%geoLonCu(i+ido,j+jdo)
66  og%geoLatCu(i,j) = dg%geoLatCu(i+ido,j+jdo)
67  og%dxCu(i,j) = dg%dxCu(i+ido,j+jdo)
68  og%dyCu(i,j) = dg%dyCu(i+ido,j+jdo)
69  og%dy_Cu(i,j) = dg%dy_Cu(i+ido,j+jdo)
70 
71  og%mask2dCu(i,j) = dg%mask2dCu(i+ido,j+jdo)
72  og%areaCu(i,j) = dg%areaCu(i+ido,j+jdo)
73  og%IareaCu(i,j) = dg%IareaCu(i+ido,j+jdo)
74  enddo ; enddo
75 
76  do i=isd,ied ; do j=jsdb,jedb
77  og%geoLonCv(i,j) = dg%geoLonCv(i+ido,j+jdo)
78  og%geoLatCv(i,j) = dg%geoLatCv(i+ido,j+jdo)
79  og%dxCv(i,j) = dg%dxCv(i+ido,j+jdo)
80  og%dyCv(i,j) = dg%dyCv(i+ido,j+jdo)
81  og%dx_Cv(i,j) = dg%dx_Cv(i+ido,j+jdo)
82 
83  og%mask2dCv(i,j) = dg%mask2dCv(i+ido,j+jdo)
84  og%areaCv(i,j) = dg%areaCv(i+ido,j+jdo)
85  og%IareaCv(i,j) = dg%IareaCv(i+ido,j+jdo)
86  enddo ; enddo
87 
88  do i=isdb,iedb ; do j=jsdb,jedb
89  og%geoLonBu(i,j) = dg%geoLonBu(i+ido,j+jdo)
90  og%geoLatBu(i,j) = dg%geoLatBu(i+ido,j+jdo)
91  og%dxBu(i,j) = dg%dxBu(i+ido,j+jdo)
92  og%dyBu(i,j) = dg%dyBu(i+ido,j+jdo)
93  og%areaBu(i,j) = dg%areaBu(i+ido,j+jdo)
94  og%CoriolisBu(i,j) = dg%CoriolisBu(i+ido,j+jdo)
95  og%mask2dBu(i,j) = dg%mask2dBu(i+ido,j+jdo)
96  enddo ; enddo
97 
98  og%bathymetry_at_vel = dg%bathymetry_at_vel
99  if (og%bathymetry_at_vel) then
100  do i=isdb,iedb ; do j=jsd,jed
101  og%Dblock_u(i,j) = dg%Dblock_u(i+ido,j+jdo)
102  og%Dopen_u(i,j) = dg%Dopen_u(i+ido,j+jdo)
103  enddo ; enddo
104  do i=isd,ied ; do j=jsdb,jedb
105  og%Dblock_v(i,j) = dg%Dblock_v(i+ido,j+jdo)
106  og%Dopen_v(i,j) = dg%Dopen_v(i+ido,j+jdo)
107  enddo ; enddo
108  endif
109 
110  og%gridLonT(og%isg:og%ieg) = dg%gridLonT(dg%isg:dg%ieg)
111  og%gridLatT(og%jsg:og%jeg) = dg%gridLatT(dg%jsg:dg%jeg)
112  ! The more complicated logic here avoids segmentation faults if one grid uses
113  ! global symmetric memory while the other does not. Because a northeast grid
114  ! convention is being used, the upper bounds for each array correspond.
115  ! Note that the dynamic grid always uses symmetric memory.
116  ido2 = dg%IegB-og%IegB ; igst = max(og%IsgB, (dg%isg-1)-ido2)
117  jdo2 = dg%JegB-og%JegB ; jgst = max(og%JsgB, (dg%jsg-1)-jdo2)
118  do i=igst,og%IegB ; og%gridLonB(i) = dg%gridLonB(i+ido2) ; enddo
119  do j=jgst,og%JegB ; og%gridLatB(j) = dg%gridLatB(j+jdo2) ; enddo
120 
121  ! Copy various scalar variables and strings.
122  og%x_axis_units = dg%x_axis_units ; og%y_axis_units = dg%y_axis_units
123  og%areaT_global = dg%areaT_global ; og%IareaT_global = dg%IareaT_global
124  og%south_lat = dg%south_lat ; og%west_lon = dg%west_lon
125  og%len_lat = dg%len_lat ; og%len_lon = dg%len_lon
126  og%Rad_Earth = dg%Rad_Earth ; og%max_depth = dg%max_depth
127 
128 ! Update the halos in case the dynamic grid has smaller halos than the ocean grid.
129  call pass_var(og%areaT, og%Domain)
130  call pass_var(og%bathyT, og%Domain)
131  call pass_var(og%geoLonT, og%Domain)
132  call pass_var(og%geoLatT, og%Domain)
133  call pass_vector(og%dxT, og%dyT, og%Domain, to_all+scalar_pair, agrid)
134  call pass_vector(og%dF_dx, og%dF_dy, og%Domain, to_all, agrid)
135  call pass_vector(og%cos_rot, og%sin_rot, og%Domain, to_all, agrid)
136  call pass_var(og%mask2dT, og%Domain)
137 
138  call pass_vector(og%areaCu, og%areaCv, og%Domain, to_all+scalar_pair, cgrid_ne)
139  call pass_vector(og%dyCu, og%dxCv, og%Domain, to_all+scalar_pair, cgrid_ne)
140  call pass_vector(og%dxCu, og%dyCv, og%Domain, to_all+scalar_pair, cgrid_ne)
141  call pass_vector(og%dy_Cu, og%dx_Cv, og%Domain, to_all+scalar_pair, cgrid_ne)
142  call pass_vector(og%mask2dCu, og%mask2dCv, og%Domain, to_all+scalar_pair, cgrid_ne)
143  call pass_vector(og%IareaCu, og%IareaCv, og%Domain, to_all+scalar_pair, cgrid_ne)
144  call pass_vector(og%IareaCu, og%IareaCv, og%Domain, to_all+scalar_pair, cgrid_ne)
145  call pass_vector(og%geoLatCu, og%geoLatCv, og%Domain, to_all+scalar_pair, cgrid_ne)
146 
147  call pass_var(og%areaBu, og%Domain, position=corner)
148  call pass_var(og%geoLonBu, og%Domain, position=corner, inner_halo=og%isc-isd)
149  call pass_var(og%geoLatBu, og%Domain, position=corner)
150  call pass_vector(og%dxBu, og%dyBu, og%Domain, to_all+scalar_pair, bgrid_ne)
151  call pass_var(og%CoriolisBu, og%Domain, position=corner)
152  call pass_var(og%mask2dBu, og%Domain, position=corner)
153 
154  if (og%bathymetry_at_vel) then
155  call pass_vector(og%Dblock_u, og%Dblock_v, og%Domain, to_all+scalar_pair, cgrid_ne)
156  call pass_vector(og%Dopen_u, og%Dopen_v, og%Domain, to_all+scalar_pair, cgrid_ne)
157  endif
158 
159  call set_derived_metrics(og, us)
160 
161 end subroutine copy_dyngrid_to_mom_grid
162 
163 
164 !> Copies information from an ocean_grid_type into a dynamic (shared)
165 !! horizontal grid type.
166 subroutine copy_mom_grid_to_dyngrid(oG, dG, US)
167  type(ocean_grid_type), intent(in) :: og !< Ocean grid type
168  type(dyn_horgrid_type), intent(inout) :: dg !< Common horizontal grid type
169  type(unit_scale_type), optional, intent(in) :: us !< A dimensional unit scaling type
170 
171  integer :: isd, ied, jsd, jed ! Common data domains.
172  integer :: isdb, iedb, jsdb, jedb ! Common data domains.
173  integer :: ido, jdo, ido2, jdo2 ! Indexing offsets between the grids.
174  integer :: igst, jgst ! Global starting indices.
175  integer :: i, j
176 
177  ! MOM_grid_init and create_dyn_horgrid are called outside of this routine.
178  ! This routine copies over the fields that were set by MOM_initialized_fixed.
179 
180  ! Determine the indexing offsets between the grids.
181  ido = og%idG_offset - dg%idG_offset
182  jdo = og%jdG_offset - dg%jdG_offset
183 
184  isd = max(dg%isd, og%isd+ido) ; jsd = max(dg%jsd, og%jsd+jdo)
185  ied = min(dg%ied, og%ied+ido) ; jed = min(dg%jed, og%jed+jdo)
186  isdb = max(dg%IsdB, og%IsdB+ido) ; jsdb = max(dg%JsdB, og%JsdB+jdo)
187  iedb = min(dg%IedB, og%IedB+ido) ; jedb = min(dg%JedB, og%JedB+jdo)
188 
189  ! Check that the grids conform.
190  if ((isd > dg%isc) .or. (ied < dg%ied) .or. (jsd > dg%jsc) .or. (jed > dg%jed)) &
191  call mom_error(fatal, "copy_dyngrid_to_MOM_grid called with incompatible grids.")
192 
193  do i=isd,ied ; do j=jsd,jed
194  dg%geoLonT(i,j) = og%geoLonT(i+ido,j+jdo)
195  dg%geoLatT(i,j) = og%geoLatT(i+ido,j+jdo)
196  dg%dxT(i,j) = og%dxT(i+ido,j+jdo)
197  dg%dyT(i,j) = og%dyT(i+ido,j+jdo)
198  dg%areaT(i,j) = og%areaT(i+ido,j+jdo)
199  dg%bathyT(i,j) = og%bathyT(i+ido,j+jdo)
200 
201  dg%dF_dx(i,j) = og%dF_dx(i+ido,j+jdo)
202  dg%dF_dy(i,j) = og%dF_dy(i+ido,j+jdo)
203  dg%sin_rot(i,j) = og%sin_rot(i+ido,j+jdo)
204  dg%cos_rot(i,j) = og%cos_rot(i+ido,j+jdo)
205  dg%mask2dT(i,j) = og%mask2dT(i+ido,j+jdo)
206  enddo ; enddo
207 
208  do i=isdb,iedb ; do j=jsd,jed
209  dg%geoLonCu(i,j) = og%geoLonCu(i+ido,j+jdo)
210  dg%geoLatCu(i,j) = og%geoLatCu(i+ido,j+jdo)
211  dg%dxCu(i,j) = og%dxCu(i+ido,j+jdo)
212  dg%dyCu(i,j) = og%dyCu(i+ido,j+jdo)
213  dg%dy_Cu(i,j) = og%dy_Cu(i+ido,j+jdo)
214 
215  dg%mask2dCu(i,j) = og%mask2dCu(i+ido,j+jdo)
216  dg%areaCu(i,j) = og%areaCu(i+ido,j+jdo)
217  dg%IareaCu(i,j) = og%IareaCu(i+ido,j+jdo)
218  enddo ; enddo
219 
220  do i=isd,ied ; do j=jsdb,jedb
221  dg%geoLonCv(i,j) = og%geoLonCv(i+ido,j+jdo)
222  dg%geoLatCv(i,j) = og%geoLatCv(i+ido,j+jdo)
223  dg%dxCv(i,j) = og%dxCv(i+ido,j+jdo)
224  dg%dyCv(i,j) = og%dyCv(i+ido,j+jdo)
225  dg%dx_Cv(i,j) = og%dx_Cv(i+ido,j+jdo)
226 
227  dg%mask2dCv(i,j) = og%mask2dCv(i+ido,j+jdo)
228  dg%areaCv(i,j) = og%areaCv(i+ido,j+jdo)
229  dg%IareaCv(i,j) = og%IareaCv(i+ido,j+jdo)
230  enddo ; enddo
231 
232  do i=isdb,iedb ; do j=jsdb,jedb
233  dg%geoLonBu(i,j) = og%geoLonBu(i+ido,j+jdo)
234  dg%geoLatBu(i,j) = og%geoLatBu(i+ido,j+jdo)
235  dg%dxBu(i,j) = og%dxBu(i+ido,j+jdo)
236  dg%dyBu(i,j) = og%dyBu(i+ido,j+jdo)
237  dg%areaBu(i,j) = og%areaBu(i+ido,j+jdo)
238  dg%CoriolisBu(i,j) = og%CoriolisBu(i+ido,j+jdo)
239  dg%mask2dBu(i,j) = og%mask2dBu(i+ido,j+jdo)
240  enddo ; enddo
241 
242  dg%bathymetry_at_vel = og%bathymetry_at_vel
243  if (dg%bathymetry_at_vel) then
244  do i=isdb,iedb ; do j=jsd,jed
245  dg%Dblock_u(i,j) = og%Dblock_u(i+ido,j+jdo)
246  dg%Dopen_u(i,j) = og%Dopen_u(i+ido,j+jdo)
247  enddo ; enddo
248  do i=isd,ied ; do j=jsdb,jedb
249  dg%Dblock_v(i,j) = og%Dblock_v(i+ido,j+jdo)
250  dg%Dopen_v(i,j) = og%Dopen_v(i+ido,j+jdo)
251  enddo ; enddo
252  endif
253 
254  dg%gridLonT(dg%isg:dg%ieg) = og%gridLonT(og%isg:og%ieg)
255  dg%gridLatT(dg%jsg:dg%jeg) = og%gridLatT(og%jsg:og%jeg)
256 
257  ! The more complicated logic here avoids segmentation faults if one grid uses
258  ! global symmetric memory while the other does not. Because a northeast grid
259  ! convention is being used, the upper bounds for each array correspond.
260  ! Note that the dynamic grid always uses symmetric memory.
261  ido2 = og%IegB-dg%IegB ; igst = max(dg%isg-1, og%IsgB-ido2)
262  jdo2 = og%JegB-dg%JegB ; jgst = max(dg%jsg-1, og%JsgB-jdo2)
263  do i=igst,dg%IegB ; dg%gridLonB(i) = og%gridLonB(i+ido2) ; enddo
264  do j=jgst,dg%JegB ; dg%gridLatB(j) = og%gridLatB(j+jdo2) ; enddo
265 
266  ! Copy various scalar variables and strings.
267  dg%x_axis_units = og%x_axis_units ; dg%y_axis_units = og%y_axis_units
268  dg%areaT_global = og%areaT_global ; dg%IareaT_global = og%IareaT_global
269  dg%south_lat = og%south_lat ; dg%west_lon = og%west_lon
270  dg%len_lat = og%len_lat ; dg%len_lon = og%len_lon
271  dg%Rad_Earth = og%Rad_Earth ; dg%max_depth = og%max_depth
272 
273 ! Update the halos in case the dynamic grid has smaller halos than the ocean grid.
274  call pass_var(dg%areaT, dg%Domain)
275  call pass_var(dg%bathyT, dg%Domain)
276  call pass_var(dg%geoLonT, dg%Domain)
277  call pass_var(dg%geoLatT, dg%Domain)
278  call pass_vector(dg%dxT, dg%dyT, dg%Domain, to_all+scalar_pair, agrid)
279  call pass_vector(dg%dF_dx, dg%dF_dy, dg%Domain, to_all, agrid)
280  call pass_vector(dg%cos_rot, dg%sin_rot, dg%Domain, to_all, agrid)
281  call pass_var(dg%mask2dT, dg%Domain)
282 
283  call pass_vector(dg%areaCu, dg%areaCv, dg%Domain, to_all+scalar_pair, cgrid_ne)
284  call pass_vector(dg%dyCu, dg%dxCv, dg%Domain, to_all+scalar_pair, cgrid_ne)
285  call pass_vector(dg%dxCu, dg%dyCv, dg%Domain, to_all+scalar_pair, cgrid_ne)
286  call pass_vector(dg%dy_Cu, dg%dx_Cv, dg%Domain, to_all+scalar_pair, cgrid_ne)
287  call pass_vector(dg%mask2dCu, dg%mask2dCv, dg%Domain, to_all+scalar_pair, cgrid_ne)
288  call pass_vector(dg%IareaCu, dg%IareaCv, dg%Domain, to_all+scalar_pair, cgrid_ne)
289  call pass_vector(dg%IareaCu, dg%IareaCv, dg%Domain, to_all+scalar_pair, cgrid_ne)
290  call pass_vector(dg%geoLatCu, dg%geoLatCv, dg%Domain, to_all+scalar_pair, cgrid_ne)
291 
292  call pass_var(dg%areaBu, dg%Domain, position=corner)
293  call pass_var(dg%geoLonBu, dg%Domain, position=corner, inner_halo=dg%isc-isd)
294  call pass_var(dg%geoLatBu, dg%Domain, position=corner)
295  call pass_vector(dg%dxBu, dg%dyBu, dg%Domain, to_all+scalar_pair, bgrid_ne)
296  call pass_var(dg%CoriolisBu, dg%Domain, position=corner)
297  call pass_var(dg%mask2dBu, dg%Domain, position=corner)
298 
299  if (dg%bathymetry_at_vel) then
300  call pass_vector(dg%Dblock_u, dg%Dblock_v, dg%Domain, to_all+scalar_pair, cgrid_ne)
301  call pass_vector(dg%Dopen_u, dg%Dopen_v, dg%Domain, to_all+scalar_pair, cgrid_ne)
302  endif
303 
304  call set_derived_dyn_horgrid(dg, us)
305 
306 end subroutine copy_mom_grid_to_dyngrid
307 
308 end module mom_transcribe_grid
mom_dyn_horgrid
Contains a shareable dynamic type for describing horizontal grids and metric data and utilty routines...
Definition: MOM_dyn_horgrid.F90:3
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
mom_transcribe_grid::copy_mom_grid_to_dyngrid
subroutine, public copy_mom_grid_to_dyngrid(oG, dG, US)
Copies information from an ocean_grid_type into a dynamic (shared) horizontal grid type.
Definition: MOM_transcribe_grid.F90:167
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_dyn_horgrid::set_derived_dyn_horgrid
subroutine, public set_derived_dyn_horgrid(G, US)
set_derived_dyn_horgrid calculates metric terms that are derived from other metrics.
Definition: MOM_dyn_horgrid.F90:317
mom_domains::pass_vector
Do a halo update on a pair of arrays representing the two components of a vector.
Definition: MOM_domains.F90:54
mom_transcribe_grid
Module with routines for copying information from a shared dynamic horizontal grid to an ocean-specif...
Definition: MOM_transcribe_grid.F90:3
mom_domains::to_all
integer, parameter, public to_all
A flag for passing in all directions.
Definition: MOM_domains.F90:132
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_transcribe_grid::copy_dyngrid_to_mom_grid
subroutine, public copy_dyngrid_to_mom_grid(dG, oG, US)
Copies information from a dynamic (shared) horizontal grid type into an ocean_grid_type.
Definition: MOM_transcribe_grid.F90:23
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_grid::set_derived_metrics
subroutine, public set_derived_metrics(G, US)
set_derived_metrics calculates metric terms that are derived from other metrics.
Definition: MOM_grid.F90:413
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_dyn_horgrid::dyn_horgrid_type
Describes the horizontal ocean grid with only dynamic memory arrays.
Definition: MOM_dyn_horgrid.F90:23
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26