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.
crsini.F90 in branches/UKMO/test_moci_test_suite_namelist_read/NEMOGCM/NEMO/OPA_SRC/CRS – NEMO

source: branches/UKMO/test_moci_test_suite_namelist_read/NEMOGCM/NEMO/OPA_SRC/CRS/crsini.F90 @ 9383

Last change on this file since 9383 was 9383, checked in by andmirek, 6 years ago

#2050 fixes and changes
File size: 12.1 KB
Line 
1MODULE crsini   
2   !!======================================================================
3   !!                         ***  MODULE crsini  ***
4   !!            Manage the grid coarsening module initialization
5   !!======================================================================
6   !!  History     2012-05   (J. Simeon, G. Madec, C. Ethe, C. Calone) Original code
7   !!----------------------------------------------------------------------
8
9   USE timing                   ! Timing
10   USE par_oce                  ! For parameter jpi,jpj,jphgr_msh
11   USE dom_oce                  ! For parameters in par_oce (jperio, lk_vvl)
12   USE crs                  ! Coarse grid domain
13   USE phycst, ONLY: omega, rad ! physical constants
14   USE wrk_nemo 
15   USE in_out_manager
16   USE par_kind, ONLY: wp
17   USE iom
18   USE crsdom
19   USE crsdomwri
20   USE crslbclnk
21   USE lib_mpp
22
23   IMPLICIT NONE
24   PRIVATE
25
26   PUBLIC  crs_init
27   PRIVATE crs_namelist
28
29   !! * Substitutions
30#  include "domzgr_substitute.h90"
31
32   !! $Id$
33CONTAINS
34   
35   SUBROUTINE crs_init 
36      !!-------------------------------------------------------------------
37      !!                     *** SUBROUTINE crs_init
38      !!  ** Purpose : Initialization of the grid coarsening module 
39      !!               1. Read namelist
40      !!               X2. MOVE TO crs_dom.F90 Set the domain definitions for coarse grid
41      !!                 a. Define the coarse grid starting/ending indices on parent grid
42      !!                    Here is where the T-pivot or F-pivot grids are discerned
43      !!                 b. TODO.  Include option for north-centric or equator-centric binning.
44      !!                 (centered only for odd reduction factors; even reduction bins bias equator to the south)
45      !!               3. Mask and mesh creation. => calls to crsfun
46      !!                  a. Use crsfun_mask to generate tmask,umask, vmask, fmask.
47      !!                  b. Use crsfun_coordinates to get coordinates
48      !!                  c. Use crsfun_UV to get horizontal scale factors
49      !!                  d. Use crsfun_TW to get initial vertical scale factors   
50      !!               4. Volumes and weights jes.... TODO. Updates for vvl? Where to do this? crsstp.F90?
51      !!                  a. Calculate initial coarse grid box volumes: pvol_T, pvol_W
52      !!                  b. Calculate initial coarse grid surface-averaging weights
53      !!                  c. Calculate initial coarse grid volume-averaging weights
54      !!                 
55      !!               X5. MOVE TO crs_dom_wri.F90 Using iom_rstput output the initial meshmask.
56      !!               ?. Another set of "masks" to generate
57      !!                  are the u- and v- surface areas for U- and V- area-weighted means.
58      !!                  Need to put this somewhere in section 3?
59      !! jes. What do to about the vvl?  GM.  could separate the weighting (denominator), so
60      !! output C*dA or C*dV as summation not mran, then do mean (division) at moment of output.
61      !! As is, crsfun takes into account vvl.   
62      !!      Talked about pre-setting the surface array to avoid IF/ENDIFS and division.
63      !!      But have then to make that preset array here and elsewhere.
64      !!      that is called every timestep...
65      !!
66      !!               - Read in pertinent data ?
67      !!-------------------------------------------------------------------
68      !! Local variables
69      INTEGER  :: ji,jj,jk      ! dummy indices
70      INTEGER  :: ierr                                ! allocation error status
71      INTEGER  ::   ios                 ! Local integer output status for namelist read
72      REAL(wp), DIMENSION(:,:,:), POINTER :: zfse3t, zfse3u, zfse3v, zfse3w
73
74      NAMELIST/namcrs/ nn_factx, nn_facty, nn_binref, nn_msh_crs, nn_crs_kz, ln_crs_wn
75      !!----------------------------------------------------------------------
76      !
77      IF( nn_timing == 1 )  CALL timing_start('crs_init')
78      !
79      IF(lwp) THEN
80         WRITE(numout,*)
81         WRITE(numout,*) 'crs_init : Initializing the grid coarsening module '
82      ENDIF
83
84     !---------------------------------------------------------
85     ! 1. Read Namelist file
86     !---------------------------------------------------------
87     !
88      IF(lwm) THEN
89         REWIND( numnam_ref )              ! Namelist namrun in reference namelist : Parameters of the run
90         READ  ( numnam_ref, namcrs, IOSTAT = ios, ERR = 901)
91901      CONTINUE
92      ENDIF
93      call mpp_bcast(ios)
94      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcrs in reference namelist', lwp )
95      IF(lwm) THEN
96         REWIND( numnam_cfg )              ! Namelist namrun in configuration namelist : Parameters of the run
97         READ  ( numnam_cfg, namcrs, IOSTAT = ios, ERR = 902 )
98902      CONTINUE
99      ENDIF
100      call mpp_bcast(ios)
101      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcrs in configuration namelist', lwp )
102
103      IF(lwm) WRITE ( numond, namcrs )
104
105     CALL crs_namelist()
106
107     IF(lwp) THEN
108        WRITE(numout,*)
109        WRITE(numout,*) 'crs_init: Namelist namcrs '
110        WRITE(numout,*) '   coarsening factor in i-direction      nn_factx   = ', nn_factx
111        WRITE(numout,*) '   coarsening factor in j-direction      nn_facty   = ', nn_facty
112        WRITE(numout,*) '   bin centering preference              nn_binref  = ', nn_binref
113        WRITE(numout,*) '   create (=1) a mesh file or not (=0)   nn_msh_crs = ', nn_msh_crs
114        WRITE(numout,*) '   type of Kz coarsening (0,1,2)         nn_crs_kz  = ', nn_crs_kz
115        WRITE(numout,*) '   wn coarsened or computed using hdivn  ln_crs_wn  = ', ln_crs_wn
116     ENDIF
117             
118     rfactx_r = 1. / nn_factx
119     rfacty_r = 1. / nn_facty
120
121     !---------------------------------------------------------
122     ! 2. Define Global Dimensions of the coarsened grid
123     !---------------------------------------------------------
124     CALL crs_dom_def     
125
126     !---------------------------------------------------------
127     ! 3. Mask and Mesh
128     !---------------------------------------------------------
129
130     !     Set up the masks and meshes     
131
132     !  3.a. Get the masks   
133 
134     CALL crs_dom_msk
135
136
137     !  3.b. Get the coordinates
138     !      Odd-numbered reduction factor, center coordinate on T-cell
139     !      Even-numbered reduction factor, center coordinate on U-,V- faces or f-corner.
140     !     
141     IF ( nresty /= 0 .AND. nrestx /= 0 ) THEN
142        CALL crs_dom_coordinates( gphit, glamt, 'T', gphit_crs, glamt_crs ) 
143        CALL crs_dom_coordinates( gphiu, glamu, 'U', gphiu_crs, glamu_crs )       
144        CALL crs_dom_coordinates( gphiv, glamv, 'V', gphiv_crs, glamv_crs ) 
145        CALL crs_dom_coordinates( gphif, glamf, 'F', gphif_crs, glamf_crs ) 
146     ELSEIF ( nresty /= 0 .AND. nrestx == 0 ) THEN
147        CALL crs_dom_coordinates( gphiu, glamu, 'T', gphit_crs, glamt_crs )
148        CALL crs_dom_coordinates( gphiu, glamu, 'U', gphiu_crs, glamu_crs )
149        CALL crs_dom_coordinates( gphif, glamf, 'V', gphiv_crs, glamv_crs )
150        CALL crs_dom_coordinates( gphif, glamf, 'F', gphif_crs, glamf_crs )
151     ELSEIF ( nresty == 0 .AND. nrestx /= 0 ) THEN
152        CALL crs_dom_coordinates( gphiv, glamv, 'T', gphit_crs, glamt_crs )
153        CALL crs_dom_coordinates( gphif, glamf, 'U', gphiu_crs, glamu_crs )
154        CALL crs_dom_coordinates( gphiv, glamv, 'V', gphiv_crs, glamv_crs )
155        CALL crs_dom_coordinates( gphif, glamf, 'F', gphif_crs, glamf_crs )
156     ELSE
157        CALL crs_dom_coordinates( gphif, glamf, 'T', gphit_crs, glamt_crs )
158        CALL crs_dom_coordinates( gphif, glamf, 'U', gphiu_crs, glamu_crs )
159        CALL crs_dom_coordinates( gphif, glamf, 'V', gphiv_crs, glamv_crs )
160        CALL crs_dom_coordinates( gphif, glamf, 'F', gphif_crs, glamf_crs )
161     ENDIF
162
163
164     !  3.c. Get the horizontal mesh
165
166     !      3.c.1 Horizontal scale factors
167
168     CALL crs_dom_hgr( e1t, e2t, 'T', e1t_crs, e2t_crs )
169     CALL crs_dom_hgr( e1u, e2u, 'U', e1u_crs, e2u_crs )
170     CALL crs_dom_hgr( e1v, e2v, 'V', e1v_crs, e2v_crs )
171     CALL crs_dom_hgr( e1f, e2f, 'F', e1f_crs, e2f_crs )
172
173     e1e2t_crs(:,:) = e1t_crs(:,:) * e2t_crs(:,:)
174     
175     
176     !      3.c.2 Coriolis factor 
177
178      SELECT CASE( jphgr_msh )   ! type of horizontal mesh
179
180      CASE ( 0, 1, 4 )           ! mesh on the sphere
181
182         ff_crs(:,:) = 2. * omega * SIN( rad * gphif_crs(:,:) )
183
184      CASE DEFAULT 
185
186       IF(lwp)    WRITE(numout,*) 'crsini.F90. crs_init. Only jphgr_msh = 0, 1 or 4 supported' 
187 
188      END SELECT 
189
190     !    3.d.1 mbathy ( vertical k-levels of bathymetry )     
191
192     CALL crs_dom_bat
193     
194     !
195     CALL wrk_alloc(jpi, jpj, jpk, zfse3t, zfse3u, zfse3v, zfse3w )
196     !
197     zfse3t(:,:,:) = fse3t(:,:,:)
198     zfse3u(:,:,:) = fse3u(:,:,:)
199     zfse3v(:,:,:) = fse3v(:,:,:)
200     zfse3w(:,:,:) = fse3w(:,:,:)
201
202     !    3.d.2   Surfaces
203     CALL crs_dom_sfc( tmask, 'W', e1e2w_crs, e1e2w_msk, p_e1=e1t, p_e2=e2t    )
204     CALL crs_dom_sfc( umask, 'U', e2e3u_crs, e2e3u_msk, p_e2=e2u, p_e3=zfse3u )
205     CALL crs_dom_sfc( vmask, 'V', e1e3v_crs, e1e3v_msk, p_e1=e1v, p_e3=zfse3v )
206   
207     facsurfu(:,:,:) = umask_crs(:,:,:) * e2e3u_msk(:,:,:) / e2e3u_crs(:,:,:)
208     facsurfv(:,:,:) = vmask_crs(:,:,:) * e1e3v_msk(:,:,:) / e1e3v_crs(:,:,:)
209
210     !    3.d.3   Vertical scale factors
211     !
212   
213 
214     CALL crs_dom_e3( e1t, e2t, zfse3t, e1e2w_crs, 'T', tmask, e3t_crs, e3t_max_crs)
215     CALL crs_dom_e3( e1u, e2u, zfse3u, e2e3u_crs, 'U', umask, e3u_crs, e3u_max_crs)
216     CALL crs_dom_e3( e1v, e2v, zfse3v, e1e3v_crs, 'V', vmask, e3v_crs, e3v_max_crs)
217     CALL crs_dom_e3( e1t, e2t, zfse3w, e1e2w_crs, 'W', tmask, e3w_crs, e3w_max_crs)
218
219     ! Reset 0 to e3t_0 or e3w_0
220     DO jk = 1, jpk
221        DO ji = 1, jpi_crs
222           DO jj = 1, jpj_crs
223              IF( e3t_crs(ji,jj,jk) == 0._wp ) e3t_crs(ji,jj,jk) = e3t_1d(jk)
224              IF( e3w_crs(ji,jj,jk) == 0._wp ) e3w_crs(ji,jj,jk) = e3w_1d(jk)
225              IF( e3u_crs(ji,jj,jk) == 0._wp ) e3u_crs(ji,jj,jk) = e3t_1d(jk)
226              IF( e3v_crs(ji,jj,jk) == 0._wp ) e3v_crs(ji,jj,jk) = e3t_1d(jk)
227           ENDDO
228        ENDDO
229     ENDDO
230
231     !    3.d.3   Vertical depth (meters)
232     CALL crs_dom_ope( gdept_0, 'MAX', 'T', tmask, gdept_crs, p_e3=zfse3t, psgn=1.0 ) 
233     CALL crs_dom_ope( gdepw_0, 'MAX', 'W', tmask, gdepw_crs, p_e3=zfse3w, psgn=1.0 )
234
235
236     !---------------------------------------------------------
237     ! 4.  Coarse grid ocean volume and averaging weights
238     !---------------------------------------------------------
239     ! 4.a. Ocean volume or area unmasked and masked
240     CALL crs_dom_facvol( tmask, 'T', e1t, e2t, zfse3t, ocean_volume_crs_t, facvol_t )
241     !
242     bt_crs(:,:,:) = ocean_volume_crs_t(:,:,:) * facvol_t(:,:,:)
243     !
244     r1_bt_crs(:,:,:) = 0._wp 
245     WHERE( bt_crs /= 0._wp ) r1_bt_crs(:,:,:) = 1._wp / bt_crs(:,:,:)
246
247     CALL crs_dom_facvol( tmask, 'W', e1t, e2t, zfse3w, ocean_volume_crs_w, facvol_w )
248     !
249     !---------------------------------------------------------
250     ! 5.  Write out coarse meshmask  (see OPA_SRC/DOM/domwri.F90 for ideas later)
251     !---------------------------------------------------------
252
253     IF( nn_msh_crs > 0 ) THEN
254        CALL dom_grid_crs   ! Save the parent grid information  & Switch to coarse grid domain
255        CALL crs_dom_wri     
256        CALL dom_grid_glo   ! Return to parent grid domain
257     ENDIF
258     
259     !---------------------------------------------------------
260     ! 7. Finish and clean-up
261     !---------------------------------------------------------
262     CALL wrk_dealloc(jpi, jpj, jpk, zfse3t, zfse3u, zfse3v, zfse3w )
263
264
265   END SUBROUTINE crs_init
266
267   SUBROUTINE crs_namelist()
268     !!---------------------------------------------------------------------
269     !!                   ***  ROUTINE crs_namelist  ***
270     !!                     
271     !! ** Purpose :   Broadcast namelist variables read by procesor lwm
272     !!
273     !! ** Method  :   use lib_mpp
274     !!----------------------------------------------------------------------
275#if defined key_mpp_mpi
276      CALL mpp_bcast(nn_factx)
277      CALL mpp_bcast(nn_facty)
278      CALL mpp_bcast(nn_binref)
279      CALL mpp_bcast(nn_msh_crs)
280      CALL mpp_bcast(nn_crs_kz)
281      CALL mpp_bcast(ln_crs_wn)
282#endif
283   END SUBROUTINE crs_namelist   
284   !!======================================================================
285
286END MODULE crsini
Note: See TracBrowser for help on using the repository browser.