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 @ 9366

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

#2050 first version. Compiled OK in moci test suite

File size: 12.0 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      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcrs in reference namelist', lwm )
92         REWIND( numnam_cfg )              ! Namelist namrun in configuration namelist : Parameters of the run
93         READ  ( numnam_cfg, namcrs, IOSTAT = ios, ERR = 902 )
94902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcrs in configuration namelist', lwm )
95      ENDIF
96      IF(lwm) WRITE ( numond, namcrs )
97
98     CALL crs_namelist()
99
100     IF(lwp) THEN
101        WRITE(numout,*)
102        WRITE(numout,*) 'crs_init: Namelist namcrs '
103        WRITE(numout,*) '   coarsening factor in i-direction      nn_factx   = ', nn_factx
104        WRITE(numout,*) '   coarsening factor in j-direction      nn_facty   = ', nn_facty
105        WRITE(numout,*) '   bin centering preference              nn_binref  = ', nn_binref
106        WRITE(numout,*) '   create (=1) a mesh file or not (=0)   nn_msh_crs = ', nn_msh_crs
107        WRITE(numout,*) '   type of Kz coarsening (0,1,2)         nn_crs_kz  = ', nn_crs_kz
108        WRITE(numout,*) '   wn coarsened or computed using hdivn  ln_crs_wn  = ', ln_crs_wn
109     ENDIF
110             
111     rfactx_r = 1. / nn_factx
112     rfacty_r = 1. / nn_facty
113
114     !---------------------------------------------------------
115     ! 2. Define Global Dimensions of the coarsened grid
116     !---------------------------------------------------------
117     CALL crs_dom_def     
118
119     !---------------------------------------------------------
120     ! 3. Mask and Mesh
121     !---------------------------------------------------------
122
123     !     Set up the masks and meshes     
124
125     !  3.a. Get the masks   
126 
127     CALL crs_dom_msk
128
129
130     !  3.b. Get the coordinates
131     !      Odd-numbered reduction factor, center coordinate on T-cell
132     !      Even-numbered reduction factor, center coordinate on U-,V- faces or f-corner.
133     !     
134     IF ( nresty /= 0 .AND. nrestx /= 0 ) THEN
135        CALL crs_dom_coordinates( gphit, glamt, 'T', gphit_crs, glamt_crs ) 
136        CALL crs_dom_coordinates( gphiu, glamu, 'U', gphiu_crs, glamu_crs )       
137        CALL crs_dom_coordinates( gphiv, glamv, 'V', gphiv_crs, glamv_crs ) 
138        CALL crs_dom_coordinates( gphif, glamf, 'F', gphif_crs, glamf_crs ) 
139     ELSEIF ( nresty /= 0 .AND. nrestx == 0 ) THEN
140        CALL crs_dom_coordinates( gphiu, glamu, 'T', gphit_crs, glamt_crs )
141        CALL crs_dom_coordinates( gphiu, glamu, 'U', gphiu_crs, glamu_crs )
142        CALL crs_dom_coordinates( gphif, glamf, 'V', gphiv_crs, glamv_crs )
143        CALL crs_dom_coordinates( gphif, glamf, 'F', gphif_crs, glamf_crs )
144     ELSEIF ( nresty == 0 .AND. nrestx /= 0 ) THEN
145        CALL crs_dom_coordinates( gphiv, glamv, 'T', gphit_crs, glamt_crs )
146        CALL crs_dom_coordinates( gphif, glamf, 'U', gphiu_crs, glamu_crs )
147        CALL crs_dom_coordinates( gphiv, glamv, 'V', gphiv_crs, glamv_crs )
148        CALL crs_dom_coordinates( gphif, glamf, 'F', gphif_crs, glamf_crs )
149     ELSE
150        CALL crs_dom_coordinates( gphif, glamf, 'T', gphit_crs, glamt_crs )
151        CALL crs_dom_coordinates( gphif, glamf, 'U', gphiu_crs, glamu_crs )
152        CALL crs_dom_coordinates( gphif, glamf, 'V', gphiv_crs, glamv_crs )
153        CALL crs_dom_coordinates( gphif, glamf, 'F', gphif_crs, glamf_crs )
154     ENDIF
155
156
157     !  3.c. Get the horizontal mesh
158
159     !      3.c.1 Horizontal scale factors
160
161     CALL crs_dom_hgr( e1t, e2t, 'T', e1t_crs, e2t_crs )
162     CALL crs_dom_hgr( e1u, e2u, 'U', e1u_crs, e2u_crs )
163     CALL crs_dom_hgr( e1v, e2v, 'V', e1v_crs, e2v_crs )
164     CALL crs_dom_hgr( e1f, e2f, 'F', e1f_crs, e2f_crs )
165
166     e1e2t_crs(:,:) = e1t_crs(:,:) * e2t_crs(:,:)
167     
168     
169     !      3.c.2 Coriolis factor 
170
171      SELECT CASE( jphgr_msh )   ! type of horizontal mesh
172
173      CASE ( 0, 1, 4 )           ! mesh on the sphere
174
175         ff_crs(:,:) = 2. * omega * SIN( rad * gphif_crs(:,:) )
176
177      CASE DEFAULT 
178
179       IF(lwp)    WRITE(numout,*) 'crsini.F90. crs_init. Only jphgr_msh = 0, 1 or 4 supported' 
180 
181      END SELECT 
182
183     !    3.d.1 mbathy ( vertical k-levels of bathymetry )     
184
185     CALL crs_dom_bat
186     
187     !
188     CALL wrk_alloc(jpi, jpj, jpk, zfse3t, zfse3u, zfse3v, zfse3w )
189     !
190     zfse3t(:,:,:) = fse3t(:,:,:)
191     zfse3u(:,:,:) = fse3u(:,:,:)
192     zfse3v(:,:,:) = fse3v(:,:,:)
193     zfse3w(:,:,:) = fse3w(:,:,:)
194
195     !    3.d.2   Surfaces
196     CALL crs_dom_sfc( tmask, 'W', e1e2w_crs, e1e2w_msk, p_e1=e1t, p_e2=e2t    )
197     CALL crs_dom_sfc( umask, 'U', e2e3u_crs, e2e3u_msk, p_e2=e2u, p_e3=zfse3u )
198     CALL crs_dom_sfc( vmask, 'V', e1e3v_crs, e1e3v_msk, p_e1=e1v, p_e3=zfse3v )
199   
200     facsurfu(:,:,:) = umask_crs(:,:,:) * e2e3u_msk(:,:,:) / e2e3u_crs(:,:,:)
201     facsurfv(:,:,:) = vmask_crs(:,:,:) * e1e3v_msk(:,:,:) / e1e3v_crs(:,:,:)
202
203     !    3.d.3   Vertical scale factors
204     !
205   
206 
207     CALL crs_dom_e3( e1t, e2t, zfse3t, e1e2w_crs, 'T', tmask, e3t_crs, e3t_max_crs)
208     CALL crs_dom_e3( e1u, e2u, zfse3u, e2e3u_crs, 'U', umask, e3u_crs, e3u_max_crs)
209     CALL crs_dom_e3( e1v, e2v, zfse3v, e1e3v_crs, 'V', vmask, e3v_crs, e3v_max_crs)
210     CALL crs_dom_e3( e1t, e2t, zfse3w, e1e2w_crs, 'W', tmask, e3w_crs, e3w_max_crs)
211
212     ! Reset 0 to e3t_0 or e3w_0
213     DO jk = 1, jpk
214        DO ji = 1, jpi_crs
215           DO jj = 1, jpj_crs
216              IF( e3t_crs(ji,jj,jk) == 0._wp ) e3t_crs(ji,jj,jk) = e3t_1d(jk)
217              IF( e3w_crs(ji,jj,jk) == 0._wp ) e3w_crs(ji,jj,jk) = e3w_1d(jk)
218              IF( e3u_crs(ji,jj,jk) == 0._wp ) e3u_crs(ji,jj,jk) = e3t_1d(jk)
219              IF( e3v_crs(ji,jj,jk) == 0._wp ) e3v_crs(ji,jj,jk) = e3t_1d(jk)
220           ENDDO
221        ENDDO
222     ENDDO
223
224     !    3.d.3   Vertical depth (meters)
225     CALL crs_dom_ope( gdept_0, 'MAX', 'T', tmask, gdept_crs, p_e3=zfse3t, psgn=1.0 ) 
226     CALL crs_dom_ope( gdepw_0, 'MAX', 'W', tmask, gdepw_crs, p_e3=zfse3w, psgn=1.0 )
227
228
229     !---------------------------------------------------------
230     ! 4.  Coarse grid ocean volume and averaging weights
231     !---------------------------------------------------------
232     ! 4.a. Ocean volume or area unmasked and masked
233     CALL crs_dom_facvol( tmask, 'T', e1t, e2t, zfse3t, ocean_volume_crs_t, facvol_t )
234     !
235     bt_crs(:,:,:) = ocean_volume_crs_t(:,:,:) * facvol_t(:,:,:)
236     !
237     r1_bt_crs(:,:,:) = 0._wp 
238     WHERE( bt_crs /= 0._wp ) r1_bt_crs(:,:,:) = 1._wp / bt_crs(:,:,:)
239
240     CALL crs_dom_facvol( tmask, 'W', e1t, e2t, zfse3w, ocean_volume_crs_w, facvol_w )
241     !
242     !---------------------------------------------------------
243     ! 5.  Write out coarse meshmask  (see OPA_SRC/DOM/domwri.F90 for ideas later)
244     !---------------------------------------------------------
245
246     IF( nn_msh_crs > 0 ) THEN
247        CALL dom_grid_crs   ! Save the parent grid information  & Switch to coarse grid domain
248        CALL crs_dom_wri     
249        CALL dom_grid_glo   ! Return to parent grid domain
250     ENDIF
251     
252     !---------------------------------------------------------
253     ! 7. Finish and clean-up
254     !---------------------------------------------------------
255     CALL wrk_dealloc(jpi, jpj, jpk, zfse3t, zfse3u, zfse3v, zfse3w )
256
257
258   END SUBROUTINE crs_init
259
260   SUBROUTINE crs_namelist()
261     !!---------------------------------------------------------------------
262     !!                   ***  ROUTINE crs_namelist  ***
263     !!                     
264     !! ** Purpose :   Broadcast namelist variables read by procesor lwm
265     !!
266     !! ** Method  :   use lib_mpp
267     !!----------------------------------------------------------------------
268#if defined key_mpp_mpi
269      CALL mpp_bcast(nn_factx)
270      CALL mpp_bcast(nn_facty)
271      CALL mpp_bcast(nn_binref)
272      CALL mpp_bcast(nn_msh_crs)
273      CALL mpp_bcast(nn_crs_kz)
274      CALL mpp_bcast(ln_crs_wn)
275#endif
276   END SUBROUTINE crs_namelist   
277   !!======================================================================
278
279END MODULE crsini
Note: See TracBrowser for help on using the repository browser.