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/2013/dev_LOCEAN_2013/NEMOGCM/NEMO/OPA_SRC/CRS – NEMO

source: branches/2013/dev_LOCEAN_2013/NEMOGCM/NEMO/OPA_SRC/CRS/crsini.F90 @ 4154

Last change on this file since 4154 was 4154, checked in by cetlod, 11 years ago

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