Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

zgrloc.F90 @ 15733

Last change on this file since 15733 was 15278, checked in by dbruciaferri, 3 years ago
reorganising code in more general structure and adding partial-steps for MEs
File size: 11.9 KB

Line
1	MODULE zgrloc
2	!!==============================================================================
3	!! * MODULE zgrmes *
4	!! Ocean initialization : Multiple Enveloped s coordinate (MES)
5	!!==============================================================================
6	!! NEMO 4.0 ! 2021-07 (D. Bruciaferri)
7	!!----------------------------------------------------------------------
8	!
9	USE oce ! ocean variables
10	USE dom_oce ! ocean domain
11	USE depth_e3 ! depth <=> e3
12	USE closea ! closed seas
13	USE mes ! MEs-coordinates
14	!
15	USE in_out_manager ! I/O manager
16	USE iom ! I/O library
17	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
18	USE lib_mpp ! distributed memory computing library
19	USE wrk_nemo ! Memory allocation
20	USE timing ! Timing
21
22	IMPLICIT NONE
23	PRIVATE
24
25	PUBLIC zgr_loc ! called by zgrmes.F90
26	!
27	! INPUT GLOBAL VERTICAL COORDINATE SYSTEM
28	!
29	REAL(wp), POINTER, DIMENSION(:) :: e3t_1d_in , e3w_1d_in !: ref. scale factors (m)
30	REAL(wp), POINTER, DIMENSION(:) :: gdept_1d_in, gdepw_1d_in !: ref. depth (m)
31	!
32	REAL(wp), POINTER, DIMENSION(:,:) :: mbathy_in
33	!
34	REAL(wp), POINTER, DIMENSION(:,:,:) :: e3t_in, e3u_in , e3v_in , e3f_in !: scale factors [m]
35	REAL(wp), POINTER, DIMENSION(:,:,:) :: e3w_in, e3uw_in, e3vw_in !: - -
36	!
37	REAL(wp), POINTER, DIMENSION(:,:,:) :: gdept_in, gdepw_in !: depths [m]
38	!
39	!! * Substitutions
40	# include "vectopt_loop_substitute.h90"
41
42	CONTAINS
43
44	! =====================================================================================================
45
46	SUBROUTINE zgr_loc
47	!!---------------------------------------------------------------------
48	!! * ROUTINE zgr_loc *
49	!!
50	!! ** Purpose : Wrap the steps to generate a localised vertical
51	!! coordinate systems
52	!!----------------------------------------------------------------------
53	INTEGER :: ji, jj, jk ! dummy loop index
54	!-----------------------------------------------------------------------
55	!
56	CALL wrk_alloc( jpk, gdept_1d_in, gdepw_1d_in, e3t_1d_in , e3w_1d_in)
57	CALL wrk_alloc( jpi, jpj, mbathy_in)
58	CALL wrk_alloc( jpi, jpj, jpk, gdept_in, gdepw_in, e3t_in, e3w_in)
59	CALL wrk_alloc( jpi, jpj, jpk, e3u_in , e3v_in , e3f_in, e3uw_in, e3vw_in)
60	!
61	! Reading the input vertical grid that will be used globally
62	CALL zgr_dom_read
63	! Creating the local cooridnate system within the global input vertical grid
64	IF ( ln_mes ) CALL zgr_mes_local
65	!
66	CALL wrk_dealloc( jpk, gdept_1d_in, gdepw_1d_in, e3t_1d_in , e3w_1d_in)
67	CALL wrk_dealloc( jpi, jpj, mbathy_in)
68	CALL wrk_dealloc( jpi, jpj, jpk, gdept_in, gdepw_in, e3t_in, e3w_in)
69	CALL wrk_dealloc( jpi, jpj, jpk, e3u_in , e3v_in , e3f_in, e3uw_in, e3vw_in)
70
71	END SUBROUTINE zgr_loc
72
73	! =====================================================================================================
74
75	SUBROUTINE zgr_dom_read
76	!!---------------------------------------------------------------------
77	!! * ROUTINE zgr_dom_read *
78	!!
79	!! ** Purpose : Read the vertical information in the domain configuration file
80	!!
81	!!----------------------------------------------------------------------
82	INTEGER :: jk ! dummy loop index
83	INTEGER :: inum ! local logical unit
84	REAL(WP) :: z_cav
85	REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace
86	!!----------------------------------------------------------------------
87	!
88	IF(lwp) THEN
89	WRITE(numout,*)
90	WRITE(numout,*) ' zgr_dom_read : read the vertical coordinates in domain_cfg_global.nc file'
91	WRITE(numout,*) ' ~~~~~~~~'
92	ENDIF
93	!
94	CALL iom_open( 'domain_cfg_global.nc', inum )
95	!
96	! !* ocean cavities under iceshelves
97	CALL iom_get( inum, 'ln_isfcav', z_cav )
98	IF( z_cav == 0._wp ) THEN ; ln_isfcav = .false. ; ELSE ; ln_isfcav = .true. ; ENDIF
99	!
100	! !* vertical scale factors
101	CALL iom_get( inum, jpdom_unknown, 'e3t_1d' , e3t_1d_in ) ! 1D reference coordinate
102	CALL iom_get( inum, jpdom_unknown, 'e3w_1d' , e3w_1d_in )
103	!
104	CALL iom_get( inum, jpdom_global, 'bottom_level' , mbathy_in ) ! 2D mbathy
105	!
106	CALL iom_get( inum, jpdom_global, 'e3t_0' , e3t_in ) ! 3D coordinate
107	CALL iom_get( inum, jpdom_global, 'e3u_0' , e3u_in )
108	CALL iom_get( inum, jpdom_global, 'e3v_0' , e3v_in )
109	CALL iom_get( inum, jpdom_global, 'e3f_0' , e3f_in )
110	CALL iom_get( inum, jpdom_global, 'e3w_0' , e3w_in )
111	CALL iom_get( inum, jpdom_global, 'e3uw_0' , e3uw_in )
112	CALL iom_get( inum, jpdom_global, 'e3vw_0' , e3vw_in )
113	!
114	! !* depths
115	! !- old depth definition (obsolescent feature)
116	IF( iom_varid( inum, 'gdept_1d', ldstop = .FALSE. ) > 0 .AND. &
117	& iom_varid( inum, 'gdepw_1d', ldstop = .FALSE. ) > 0 .AND. &
118	& iom_varid( inum, 'gdept_0' , ldstop = .FALSE. ) > 0 .AND. &
119	& iom_varid( inum, 'gdepw_0' , ldstop = .FALSE. ) > 0 ) THEN
120	CALL ctl_warn( 'zgr_dom_read : old definition of depths and scale factors used ', &
121	& ' depths at t- and w-points read in the domain configuration file')
122	CALL iom_get( inum, jpdom_unknown, 'gdept_1d', gdept_1d_in )
123	CALL iom_get( inum, jpdom_unknown, 'gdepw_1d', gdepw_1d_in )
124	CALL iom_get( inum, jpdom_global , 'gdept_0' , gdept_in )
125	CALL iom_get( inum, jpdom_global , 'gdepw_0' , gdepw_in )
126	!
127	ELSE !- depths computed from e3. scale factors
128	CALL e3_to_depth( e3t_1d_in, e3w_1d_in, gdept_1d_in, gdepw_1d_in ) ! 1D reference depth
129	CALL e3_to_depth( e3t_in , e3w_in , gdept_in , gdepw_in ) ! 3D depths
130	IF(lwp) THEN
131	WRITE(numout,*)
132	WRITE(numout,*) ' GLOBAL reference 1D z-coordinate depth and scale factors:'
133	WRITE(numout, "(9x,' level gdept_1d gdepw_1d e3t_1d e3w_1d ')" )
134	WRITE(numout, "(10x, i4, 4f9.2)" ) ( jk, gdept_1d_in(jk), gdepw_1d_in(jk), e3t_1d_in(jk), e3w_1d_in(jk), jk = 1, jpk )
135	ENDIF
136	ENDIF
137	!
138	CALL iom_close( inum )
139	!
140	END SUBROUTINE zgr_dom_read
141
142	! =====================================================================================================
143
144	SUBROUTINE zgr_mes_local
145	!!---------------------------------------------------------------------
146	!! * ROUTINE zgr_dom_merge *
147	!!
148	!! ** Purpose : Create a local MEs grid within a gloabal grid
149	!! using different vertical coordinates.
150	!!
151	!!----------------------------------------------------------------------
152	INTEGER :: ji, jj, jk ! dummy loop index
153	INTEGER :: inum ! local logical unit
154	REAL(wp), DIMENSION(jpi,jpj) :: l2g_wgt ! weigths for computing model
155	! levels in transition area
156	!!----------------------------------------------------------------------
157
158	IF(lwp) THEN
159	WRITE(numout,*)
160	WRITE(numout,*) ' zgr_mes_merge : read the vertical coordinates in domain_cfg_global.nc file'
161	WRITE(numout,*) ' ~~~~~~~~'
162	ENDIF
163	!
164	CALL iom_open( 'bathy_meter.nc', inum )
165	CALL iom_get( inum, jpdom_data, 's2z_msk', l2g_msk)
166	CALL iom_get( inum, jpdom_data, 's2z_wgt', l2g_wgt)
167
168	DO jj = 1,jpj
169	DO ji = 1,jpi
170	SELECT CASE (INT(l2g_msk(ji,jj)))
171	CASE (0) ! global zps area
172	mbathy (ji,jj ) = mbathy_in(ji,jj)
173	gdept_0(ji,jj,:) = gdept_in(ji,jj,:)
174	gdepw_0(ji,jj,:) = gdepw_in(ji,jj,:)
175	e3t_0 (ji,jj,:) = e3t_in (ji,jj,:)
176	e3w_0 (ji,jj,:) = e3w_in (ji,jj,:)
177	e3u_0 (ji,jj,:) = e3u_in (ji,jj,:)
178	e3v_0 (ji,jj,:) = e3v_in (ji,jj,:)
179	e3f_0 (ji,jj,:) = e3f_in (ji,jj,:)
180	e3uw_0 (ji,jj,:) = e3uw_in (ji,jj,:)
181	e3vw_0 (ji,jj,:) = e3vw_in (ji,jj,:)
182	CASE (1) ! MEs to zps transition area
183	gdept_0(ji,jj,:) = l2g_wgt(ji,jj) * gdept_0(ji,jj,:) + &
184	& ( 1._wp - l2g_wgt(ji,jj) ) * gdept_in(ji,jj,:)
185	gdepw_0(ji,jj,:) = l2g_wgt(ji,jj) * gdepw_0(ji,jj,:) + &
186	& ( 1._wp - l2g_wgt(ji,jj) ) * gdepw_in(ji,jj,:)
187	CASE (2) ! MEs area
188	CYCLE
189	END SELECT
190	END DO
191	END DO
192	!
193	! e3t, e3w for transition zone
194	! as finite differences
195	DO jj = 1, jpj
196	DO ji = 1, jpi
197	IF ( l2g_msk(ji,jj) == 1._wp ) THEN
198	DO jk = 1,jpkm1
199	e3t_0(ji,jj,jk) = gdepw_0(ji,jj,jk+1) - gdepw_0(ji,jj,jk)
200	e3w_0(ji,jj,jk+1) = gdept_0(ji,jj,jk+1) - gdept_0(ji,jj,jk)
201	ENDDO
202	! Surface
203	jk = 1
204	e3w_0(ji,jj,jk) = 2.0_wp * (gdept_0(ji,jj,1) - gdepw_0(ji,jj,1))
205	!
206	! Bottom
207	jk = jpk
208	e3t_0(ji,jj,jk) = 2.0_wp * (gdept_0(ji,jj,jk) - gdepw_0(ji,jj,jk))
209	END IF
210	END DO
211	END DO
212	!
213	! MBATHY transition zone
214	DO jj = 1, jpj
215	DO ji = 1, jpi
216	IF ( l2g_msk(ji,jj) == 1._wp ) THEN
217	DO jk = 1, jpkm1
218	IF( scobot(ji,jj) >= gdept_0(ji,jj,jk) ) mbathy(ji,jj) = MAX( 2, jk )
219	IF( scobot(ji,jj) == 0.e0 ) mbathy(ji,jj) = 0
220	IF( scobot(ji,jj) < 0.e0 ) mbathy(ji,jj) = INT( scobot(ji,jj)) ! do we need it?
221	END DO
222	END IF
223	END DO
224	END DO
225	!
226	! Computing e3u_0, e3v_0, e3f_0, e3uw_0, e3vw_0
227	! for transition zone
228	!
229	DO jj = 1, jpjm1
230	DO ji = 1, jpim1
231	IF ( l2g_msk(ji,jj) == 1._wp ) THEN
232	DO jk = 1, jpk
233	e3u_0(ji,jj,jk)=(MIN(1,mbathy(ji,jj))* e3t_0(ji,jj,jk) + &
234	MIN(1,mbathy(ji+1,jj))*e3t_0(ji+1,jj,jk) ) / &
235	MAX( 1, MIN(1,mbathy(ji,jj))+MIN(1,mbathy(ji+1,jj)) )
236
237	e3v_0(ji,jj,jk)=(MIN(1,mbathy(ji,jj))* e3t_0(ji,jj,jk) + &
238	MIN(1,mbathy(ji,jj+1))*e3t_0(ji,jj+1,jk) ) / &
239	MAX( 1, MIN(1,mbathy(ji,jj))+MIN(1,mbathy(ji,jj+1)) )
240
241	e3uw_0(ji,jj,jk)=(MIN(1,mbathy(ji,jj))* e3w_0(ji,jj,jk) + &
242	MIN(1,mbathy(ji+1,jj))*e3w_0(ji+1,jj,jk) ) / &
243	MAX( 1, MIN(1,mbathy(ji,jj))+MIN(1,mbathy(ji+1,jj)) )
244
245	e3vw_0(ji,jj,jk)=(MIN(1,mbathy(ji,jj))* e3w_0(ji,jj,jk) + &
246	MIN(1,mbathy(ji,jj+1))*e3w_0(ji,jj+1,jk) ) / &
247	MAX( 1, MIN(1,mbathy(ji,jj))+MIN(1,mbathy(ji,jj+1)) )
248
249	e3f_0(ji,jj,jk)=(MIN(1,mbathy(ji,jj))* e3t_0(ji,jj,jk) + &
250	MIN(1,mbathy(ji+1,jj))*e3t_0(ji+1,jj,jk) + &
251	MIN(1,mbathy(ji+1,jj+1))* e3t_0(ji+1,jj+1,jk) + &
252	MIN(1,mbathy(ji,jj+1))*e3t_0(ji,jj+1,jk) ) / &
253	MAX( 1, MIN(1,mbathy(ji,jj))+MIN(1,mbathy(ji,jj+1)) &
254	+ MIN(1,mbathy(ji+1,jj))+MIN(1,mbathy(ji+1,jj+1)) )
255	END DO
256	END IF
257	END DO
258	END DO
259
260	END SUBROUTINE zgr_mes_local
261
262	END MODULE zgrloc

Note: See TracBrowser for help on using the repository browser.

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

Context Navigation

source: NEMO/branches/UKMO/tools_r4.0-HEAD_dev_MEs/DOMAINcfg/src/zgrloc.F90 @ 15733

Download in other formats: