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/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC/CRS – NEMO

source: branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC/CRS/crsini.F90 @ 6772

Last change on this file since 6772 was 6772, checked in by cbricaud, 8 years ago

clean in coarsening branch
  • Property svn:keywords set to Id
File size: 21.7 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   USE ldftra_crs
23   USE ieee_arithmetic
24
25   IMPLICIT NONE
26   PRIVATE
27
28   PUBLIC  crs_init
29
30   !! * Substitutions
31#  include "domzgr_substitute.h90"
32
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) :: zmin,zmax
73      REAL(wp), DIMENSION(:,:,:), POINTER :: zfse3t, zfse3u, zfse3v, zfse3w
74
75      NAMELIST/namcrs/ nn_factx, nn_facty, nn_binref, nn_msh_crs, nn_crs_kz, ln_crs_wn, ln_crs_top
76      !!----------------------------------------------------------------------
77      !
78      IF( nn_timing == 1 )  CALL timing_start('crs_init')
79      !
80      IF(lwp) THEN
81         WRITE(numout,*)
82         WRITE(numout,*) 'crs_init : Initializing the grid coarsening module '
83      ENDIF
84
85     !---------------------------------------------------------
86     ! 1. Read Namelist file
87     !---------------------------------------------------------
88     !
89
90      REWIND( numnam_ref )              ! Namelist namrun in reference namelist : Parameters of the run
91      READ  ( numnam_ref, namcrs, IOSTAT = ios, ERR = 901)
92901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcrs in reference namelist', lwp )
93
94      REWIND( numnam_cfg )              ! Namelist namrun in configuration namelist : Parameters of the run
95      READ  ( numnam_cfg, namcrs, IOSTAT = ios, ERR = 902 )
96902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcrs in configuration namelist', lwp )
97      IF(lwm) WRITE ( numond, namcrs )
98
99     IF(lwp) THEN
100        WRITE(numout,*)
101        WRITE(numout,*) 'crs_init: Namelist namcrs '
102        WRITE(numout,*) '   coarsening factor in i-direction      nn_factx   = ', nn_factx
103        WRITE(numout,*) '   coarsening factor in j-direction      nn_facty   = ', nn_facty
104        WRITE(numout,*) '   bin centering preference              nn_binref  = ', nn_binref
105        WRITE(numout,*) '   create (=1) a mesh file or not (=0)   nn_msh_crs = ', nn_msh_crs
106        WRITE(numout,*) '   type of Kz coarsening (0,1,2)         nn_crs_kz  = ', nn_crs_kz
107        WRITE(numout,*) '   wn coarsened or computed using hdivn  ln_crs_wn  = ', ln_crs_wn
108     ENDIF
109             
110     rfactx_r = 1. / nn_factx
111     rfacty_r = 1. / nn_facty
112
113write(narea+200,*)"crsini0",nstop; call flush(narea+200)
114
115     !---------------------------------------------------------
116     ! 2. Define Global Dimensions of the coarsened grid
117     !---------------------------------------------------------
118     CALL crs_dom_def     
119write(narea+200,*)"crsini1",nstop; call flush(narea+200)
120
121     !---------------------------------------------------------
122     ! 3. Mask and Mesh
123     !---------------------------------------------------------
124
125     !     Set up the masks and meshes     
126
127     !  3.a. Get the masks   
128 
129     CALL crs_dom_msk
130write(narea+200,*)"crsini2",nstop; call flush(narea+200)
131      CALL mppsync()
132
133     !IF( narea==279 )THEN
134     !WRITE(narea+200,*)"tutu1 ",jpi,jpj,nldi,nlei,nldj,nlej
135     !DO jj=1,jpj
136     !   WRITE(narea+200,*)"tutu2 ",jj,MINVAL(tmask(:,jj,:)),MAXVAL(tmask(:,jj,:))
137     !ENDDO
138     !ENDIF
139
140     !  3.b. Get the coordinates
141     !      Odd-numbered reduction factor, center coordinate on T-cell
142     !      Even-numbered reduction factor, center coordinate on U-,V- faces or f-corner.
143     !     
144     !IF ( nresty /= 0 .AND. nrestx /= 0 ) THEN
145        CALL crs_dom_coordinates( gphit, glamt, 'T', gphit_crs, glamt_crs ) 
146        CALL crs_dom_coordinates( gphiu, glamu, '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     !ELSEIF ( nresty /= 0 .AND. nrestx == 0 ) THEN
150     !   CALL crs_dom_coordinates( gphiu, glamu, 'T', gphit_crs, glamt_crs )
151     !   CALL crs_dom_coordinates( gphiu, glamu, '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     !ELSEIF ( nresty == 0 .AND. nrestx /= 0 ) THEN
155     !   CALL crs_dom_coordinates( gphiv, glamv, 'T', gphit_crs, glamt_crs )
156     !   CALL crs_dom_coordinates( gphif, glamf, 'U', gphiu_crs, glamu_crs )
157     !   CALL crs_dom_coordinates( gphiv, glamv, 'V', gphiv_crs, glamv_crs )
158     !   CALL crs_dom_coordinates( gphif, glamf, 'F', gphif_crs, glamf_crs )
159     !ELSE
160     !   CALL crs_dom_coordinates( gphif, glamf, 'T', gphit_crs, glamt_crs )
161     !   CALL crs_dom_coordinates( gphif, glamf, 'U', gphiu_crs, glamu_crs )
162     !   CALL crs_dom_coordinates( gphif, glamf, 'V', gphiv_crs, glamv_crs )
163     !   CALL crs_dom_coordinates( gphif, glamf, 'F', gphif_crs, glamf_crs )
164     !ENDIF
165      CALL mppsync()
166
167write(narea+200,*)"crsini3",nstop; call flush(narea+200)
168
169     !  3.c. Get the horizontal mesh
170
171     !      3.c.1 Horizontal scale factors
172
173     CALL crs_dom_hgr( e1t, e2t, 'T', e1t_crs, e2t_crs )
174     CALL crs_dom_hgr( e1u, e2u, 'U', e1u_crs, e2u_crs )
175     CALL crs_dom_hgr( e1v, e2v, 'V', e1v_crs, e2v_crs )
176     CALL crs_dom_hgr( e1f, e2f, 'F', e1f_crs, e2f_crs )
177
178     DO ji=nldi_crs,nlei_crs
179     DO jj=nldj_crs,nlej_crs
180        IF( e1t_crs(ji,jj)==0._wp .AND. tmask_crs(ji,jj,1) .NE. 0._wp )WRITE(narea+8000-1,*)"e1t_crs=0",ji,jj;CALL FLUSH(narea+8000-1) 
181        IF( e1u_crs(ji,jj)==0._wp .AND. umask_crs(ji,jj,1) .NE. 0._wp )WRITE(narea+8000-1,*)"e1u_crs=0",ji,jj;CALL FLUSH(narea+8000-1) 
182        IF( e1v_crs(ji,jj)==0._wp .AND. vmask_crs(ji,jj,1) .NE. 0._wp )WRITE(narea+8000-1,*)"e1v_crs=0",ji,jj;CALL FLUSH(narea+8000-1) 
183        IF( e1f_crs(ji,jj)==0._wp .AND. fmask_crs(ji,jj,1) .NE. 0._wp )WRITE(narea+8000-1,*)"e1f_crs=0",ji,jj;CALL FLUSH(narea+8000-1) 
184        IF( e2t_crs(ji,jj)==0._wp .AND. tmask_crs(ji,jj,1) .NE. 0._wp )WRITE(narea+8000-1,*)"e2t_crs=0",ji,jj;CALL FLUSH(narea+8000-1) 
185        IF( e2u_crs(ji,jj)==0._wp .AND. umask_crs(ji,jj,1) .NE. 0._wp )WRITE(narea+8000-1,*)"e2u_crs=0",ji,jj;CALL FLUSH(narea+8000-1) 
186        IF( e2v_crs(ji,jj)==0._wp .AND. vmask_crs(ji,jj,1) .NE. 0._wp )WRITE(narea+8000-1,*)"e2v_crs=0",ji,jj;CALL FLUSH(narea+8000-1) 
187        IF( e2f_crs(ji,jj)==0._wp .AND. fmask_crs(ji,jj,1) .NE. 0._wp )WRITE(narea+8000-1,*)"e2f_crs=0",ji,jj;CALL FLUSH(narea+8000-1) 
188     ENDDO
189     ENDDO
190
191
192     WHERE(e1t_crs == 0._wp) e1t_crs=r_inf
193     WHERE(e1u_crs == 0._wp) e1u_crs=r_inf
194     WHERE(e1v_crs == 0._wp) e1v_crs=r_inf
195     WHERE(e1f_crs == 0._wp) e1f_crs=r_inf
196     WHERE(e2t_crs == 0._wp) e2t_crs=r_inf
197     WHERE(e2u_crs == 0._wp) e2u_crs=r_inf
198     WHERE(e2v_crs == 0._wp) e2v_crs=r_inf
199     WHERE(e2f_crs == 0._wp) e2f_crs=r_inf
200
201     zmin=MINVAL(e1t_crs,mask=(tmask_crs(:,:,1)==1));CALL mpp_min(zmin);zmax=MAXVAL(e1t_crs,mask=(tmask_crs(:,:,1)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e1t_crs ",zmin,zmax
202     zmin=MINVAL(e1u_crs,mask=(umask_crs(:,:,1)==1));CALL mpp_min(zmin);zmax=MAXVAL(e1u_crs,mask=(umask_crs(:,:,1)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e1u_crs ",zmin,zmax
203     zmin=MINVAL(e1v_crs,mask=(vmask_crs(:,:,1)==1));CALL mpp_min(zmin);zmax=MAXVAL(e1v_crs,mask=(vmask_crs(:,:,1)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e1v_crs ",zmin,zmax
204     zmin=MINVAL(e1f_crs,mask=(fmask_crs(:,:,1)==1));CALL mpp_min(zmin);zmax=MAXVAL(e1f_crs,mask=(fmask_crs(:,:,1)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e1f_crs ",zmin,zmax
205     zmin=MINVAL(e2t_crs,mask=(tmask_crs(:,:,1)==1));CALL mpp_min(zmin);zmax=MAXVAL(e2t_crs,mask=(tmask_crs(:,:,1)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e2t_crs ",zmin,zmax
206     zmin=MINVAL(e2u_crs,mask=(umask_crs(:,:,1)==1));CALL mpp_min(zmin);zmax=MAXVAL(e2u_crs,mask=(umask_crs(:,:,1)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e2u_crs ",zmin,zmax
207     zmin=MINVAL(e2v_crs,mask=(vmask_crs(:,:,1)==1));CALL mpp_min(zmin);zmax=MAXVAL(e2v_crs,mask=(vmask_crs(:,:,1)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e2v_crs ",zmin,zmax
208     zmin=MINVAL(e2f_crs,mask=(fmask_crs(:,:,1)==1));CALL mpp_min(zmin);zmax=MAXVAL(e2f_crs,mask=(fmask_crs(:,:,1)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e2f_crs ",zmin,zmax
209
210
211     e1e2t_crs(:,:) = e1t_crs(:,:) * e2t_crs(:,:)
212     
213write(narea+200,*)"crsini4",nstop; call flush(narea+200)
214     
215     !      3.c.2 Coriolis factor 
216
217      SELECT CASE( jphgr_msh )   ! type of horizontal mesh
218
219      CASE ( 0, 1, 4 )           ! mesh on the sphere
220
221         ff_crs(:,:) = 2. * omega * SIN( rad * gphif_crs(:,:) )
222
223      CASE DEFAULT 
224
225       IF(lwp)    WRITE(numout,*) 'crsini.F90. crs_init. Only jphgr_msh = 0, 1 or 4 supported' 
226 
227      END SELECT 
228
229     !    3.d.1 mbathy ( vertical k-levels of bathymetry )     
230
231     CALL crs_dom_bat
232     
233     !
234     CALL wrk_alloc(jpi, jpj, jpk, zfse3t, zfse3u, zfse3v, zfse3w )
235     !
236     zfse3t(:,:,:) = e3t_0(:,:,:) !fse3t(:,:,:)
237     zfse3u(:,:,:) = e3u_0(:,:,:) !fse3u(:,:,:)
238     zfse3v(:,:,:) = e3v_0(:,:,:) !fse3v(:,:,:)
239     zfse3w(:,:,:) = e3w_0(:,:,:) !fse3w(:,:,:)
240
241     !    3.d.2   Surfaces
242     e2e3u_crs(:,:,:)=0._wp
243     e2e3u_msk(:,:,:)=0._wp
244     e1e3v_crs(:,:,:)=0._wp
245     e1e3v_msk(:,:,:)=0._wp
246     CALL crs_dom_sfc( tmask, 'W', e1e2w_crs, e1e2w_msk, p_e1=e1t, p_e2=e2t    )
247     CALL crs_dom_sfc( umask, 'U', e2e3u_crs, e2e3u_msk, p_e2=e2u, p_e3=zfse3u )
248     CALL crs_dom_sfc( vmask, 'V', e1e3v_crs, e1e3v_msk, p_e1=e1v, p_e3=zfse3v )
249
250     DO ji=nldi_crs,nlei_crs
251     DO jj=nldj_crs,nlej_crs
252     DO jk=1,jpk   
253        IF( e1e2w_crs(ji,jj,jk)==0._wp .AND. tmask_crs(ji,jj,jk)==1._wp )WRITE(narea+8000-1,*)"e1e2w_crs=0",ji,jj,jk;CALL FLUSH(narea+8000-1)
254        IF( e1e2w_msk(ji,jj,jk)==0._wp .AND. tmask_crs(ji,jj,jk)==1._wp )WRITE(narea+8000-1,*)"e1e2w_msk=0",ji,jj,jk;CALL FLUSH(narea+8000-1)
255        IF( e2e3u_crs(ji,jj,jk)==0._wp .AND. umask_crs(ji,jj,jk)==1._wp )WRITE(narea+8000-1,*)"e2e3u_crs=0",ji,jj,jk;CALL FLUSH(narea+8000-1)
256        IF( e2e3u_msk(ji,jj,jk)==0._wp .AND. umask_crs(ji,jj,jk)==1._wp )WRITE(narea+8000-1,*)"e2e3u_msk=0",ji,jj,jk;CALL FLUSH(narea+8000-1)
257        IF( e1e3v_crs(ji,jj,jk)==0._wp .AND. vmask_crs(ji,jj,jk)==1._wp )WRITE(narea+8000-1,*)"e1e3v_crs=0",ji,jj,jk;CALL FLUSH(narea+8000-1)
258        IF( e1e3v_msk(ji,jj,jk)==0._wp .AND. vmask_crs(ji,jj,jk)==1._wp )WRITE(narea+8000-1,*)"e1e3v_msk=0",ji,jj,jk;CALL FLUSH(narea+8000-1)
259     ENDDO
260     ENDDO
261     ENDDO
262write(narea+200,*)"crsini5",nstop; call flush(narea+200)
263
264!     WHERE(e1e2w_crs == 0._wp) e1e2w_crs=r_inf
265!     WHERE(e2e3u_crs == 0._wp) e2e3u_crs=r_inf
266!     WHERE(e1e3v_crs == 0._wp) e1e3v_crs=r_inf
267!     WHERE(e1e2w_msk == 0._wp) e1e2w_msk=r_inf
268!     WHERE(e2e3u_msk == 0._wp) e2e3u_msk=r_inf
269!     WHERE(e1e3v_msk == 0._wp) e1e3v_msk=r_inf
270     zmin=MINVAL(e1e2w_crs,mask=(tmask_crs(:,:,:)==1));CALL mpp_min(zmin);zmax=MAXVAL(e1e2w_crs,mask=(tmask_crs(:,:,:)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e1e2w_crs ",zmin,zmax
271     zmin=MINVAL(e2e3u_crs,mask=(umask_crs(:,:,:)==1));CALL mpp_min(zmin);zmax=MAXVAL(e2e3u_crs,mask=(umask_crs(:,:,:)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e2e3u_crs ",zmin,zmax
272     zmin=MINVAL(e1e3v_crs,mask=(vmask_crs(:,:,:)==1));CALL mpp_min(zmin);zmax=MAXVAL(e1e3v_crs,mask=(vmask_crs(:,:,:)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e1e3v_crs ",zmin,zmax
273     zmin=MINVAL(e1e2w_msk,mask=(tmask_crs(:,:,:)==1));CALL mpp_min(zmin);zmax=MAXVAL(e1e2w_msk,mask=(tmask_crs(:,:,:)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e1e2w_msk ",zmin,zmax
274     zmin=MINVAL(e2e3u_msk,mask=(umask_crs(:,:,:)==1));CALL mpp_min(zmin);zmax=MAXVAL(e2e3u_msk,mask=(umask_crs(:,:,:)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e2e3u_msk ",zmin,zmax
275     zmin=MINVAL(e1e3v_msk,mask=(vmask_crs(:,:,:)==1));CALL mpp_min(zmin);zmax=MAXVAL(e1e3v_msk,mask=(vmask_crs(:,:,:)==1));CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"crs e1e3v_msk ",zmin,zmax
276   
277     !cbr facsurfu(:,:,:) = umask_crs(:,:,:) * e2e3u_msk(:,:,:) / e2e3u_crs(:,:,:)
278     !cbr facsurfv(:,:,:) = vmask_crs(:,:,:) * e1e3v_msk(:,:,:) / e1e3v_crs(:,:,:)
279     DO jk=1,jpk
280        DO ji=1,jpi_crs
281           DO jj=1,jpj_crs
282
283              facsurfu(ji,jj,jk) = umask_crs(ji,jj,jk) * e2e3u_msk(ji,jj,jk) 
284              IF( e2e3u_crs(ji,jj,jk) .NE. 0._wp ) facsurfu(ji,jj,jk) = facsurfu(ji,jj,jk) / e2e3u_crs(ji,jj,jk)
285
286              facsurfv(ji,jj,jk) = vmask_crs(ji,jj,jk) * e1e3v_msk(ji,jj,jk) 
287              IF( e1e3v_crs(ji,jj,jk) .NE. 0._wp ) facsurfv(ji,jj,jk) = facsurfv(ji,jj,jk) / e1e3v_crs(ji,jj,jk)
288
289           ENDDO
290        ENDDO
291     ENDDO
292
293        DO ji=nldi_crs,nlei_crs
294           DO jj=nldj_crs,nlej_crs
295        IF( ABS(e2u_crs(ji,jj)) .LE. 1.e-5 )WRITE(narea+8000-1,*)"UNDERFLOW e2u_crs",ji,jj,e2u_crs(ji,jj),umask_crs(ji,jj,1) ; CALL FLUSH(narea+8000-1)
296        IF( ABS(e1v_crs(ji,jj)) .LE. 1.e-5 )WRITE(narea+8000-1,*)"UNDERFLOW e1v_crs",ji,jj,e1v_crs(ji,jj),vmask_crs(ji,jj,1) ; CALL FLUSH(narea+8000-1)
297           ENDDO
298        ENDDO
299
300
301     !    3.d.3   Vertical scale factors
302     !
303     zmin=MINVAL(e2u_crs(nldi_crs:nlei_crs,nldj_crs:nlej_crs));zmax=MAXVAL(e2u_crs(nldi_crs:nlei_crs,nldj_crs:nlej_crs));CALL mpp_min(zmin);CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"e2u_crs",zmin,zmax;CALL FLUSH(numout) 
304     zmin=MINVAL(e1v_crs(nldi_crs:nlei_crs,nldj_crs:nlej_crs));zmax=MAXVAL(e1v_crs(nldi_crs:nlei_crs,nldj_crs:nlej_crs));CALL mpp_min(zmin);CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"e1v_crs",zmin,zmax;CALL FLUSH(numout) 
305     zmin=MINVAL(e1e2w_crs(nldi_crs:nlei_crs,nldj_crs:nlej_crs,:));zmax=MAXVAL(e1e2w_crs(nldi_crs:nlei_crs,nldj_crs:nlej_crs,:));CALL mpp_min(zmin);CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"e1e2w_crs",zmin,zmax;CALL FLUSH(numout) 
306     zmin=MINVAL(zfse3u);zmax=MAXVAL(zfse3u);CALL mpp_min(zmin);CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"zfse3u",zmin,zmax;CALL FLUSH(numout) 
307     zmin=MINVAL(zfse3v);zmax=MAXVAL(zfse3v);CALL mpp_min(zmin);CALL mpp_max(zmax);IF(lwp)WRITE(numout,*)"zfse3v",zmin,zmax;CALL FLUSH(numout) 
308
309     CALL crs_dom_e3( e1t, e2t, zfse3t, p_sfc_3d_crs=e1e2w_crs, cd_type='T', p_mask=tmask, p_e3_crs=e3t_0_crs, p_e3_max_crs=e3t_max_0_crs)
310     CALL crs_dom_e3( e1t, e2t, zfse3w, p_sfc_3d_crs=e1e2w_crs, cd_type='W', p_mask=tmask, p_e3_crs=e3w_0_crs, p_e3_max_crs=e3w_max_0_crs)
311     CALL crs_dom_e3( e1u, e2u, zfse3u, p_sfc_2d_crs=e2u_crs  , cd_type='U', p_mask=umask, p_e3_crs=e3u_0_crs, p_e3_max_crs=e3u_max_0_crs)
312     CALL crs_dom_e3( e1v, e2v, zfse3v, p_sfc_2d_crs=e1v_crs  , cd_type='V', p_mask=vmask, p_e3_crs=e3v_0_crs, p_e3_max_crs=e3v_max_0_crs)
313     WHERE(e3t_max_0_crs == 0._wp) e3t_max_0_crs=r_inf
314     WHERE(e3u_max_0_crs == 0._wp) e3u_max_0_crs=r_inf
315     WHERE(e3v_max_0_crs == 0._wp) e3v_max_0_crs=r_inf
316     WHERE(e3w_max_0_crs == 0._wp) e3w_max_0_crs=r_inf
317
318write(narea+200,*)"crsini6",nstop; call flush(narea+200)
319#if defined key_vvl
320     e3t_max_n_crs=e3t_max_0_crs
321     e3u_max_n_crs=e3u_max_0_crs
322     e3v_max_n_crs=e3v_max_0_crs
323     e3w_max_n_crs=e3w_max_0_crs
324#endif
325
326     ht_0_crs(:,:)=0._wp
327     DO jk = 1, jpk
328        ht_0_crs(:,:)=ht_0_crs(:,:)+e3t_0_crs(:,:,jk)*tmask_crs(:,:,jk)
329     ENDDO
330
331#if defined key_vvl
332     e3t_0_crs(:,:,:) = e3t_0_crs(:,:,:) * tmask_crs(:,:,:)
333     e3u_0_crs(:,:,:) = e3u_0_crs(:,:,:) * umask_crs(:,:,:)
334     e3v_0_crs(:,:,:) = e3v_0_crs(:,:,:) * vmask_crs(:,:,:)
335     e3w_0_crs(:,:,:) = e3w_0_crs(:,:,:) * tmask_crs(:,:,:)
336#endif
337
338write(narea+200,*)"crsini7",nstop; call flush(narea+200)
339     ! Reset 0 to e3t_0 or e3w_0
340     DO jk = 1, jpk
341        DO ji = 1, jpi_crs
342           DO jj = 1, jpj_crs
343              IF( e3t_0_crs(ji,jj,jk) == 0._wp ) e3t_0_crs(ji,jj,jk) = e3t_1d(jk)
344              IF( e3w_0_crs(ji,jj,jk) == 0._wp ) e3w_0_crs(ji,jj,jk) = e3w_1d(jk)
345              IF( e3u_0_crs(ji,jj,jk) == 0._wp ) e3u_0_crs(ji,jj,jk) = e3t_1d(jk)
346              IF( e3v_0_crs(ji,jj,jk) == 0._wp ) e3v_0_crs(ji,jj,jk) = e3t_1d(jk)
347           ENDDO
348        ENDDO
349     ENDDO
350
351#if defined key_vvl
352     e3t_b_crs(:,:,:) = e3t_0_crs(:,:,:)
353     e3u_b_crs(:,:,:) = e3u_0_crs(:,:,:)
354     e3v_b_crs(:,:,:) = e3v_0_crs(:,:,:)
355     e3w_b_crs(:,:,:) = e3w_0_crs(:,:,:)
356
357     e3t_n_crs(:,:,:) = e3t_0_crs(:,:,:)
358     e3u_n_crs(:,:,:) = e3u_0_crs(:,:,:)
359     e3v_n_crs(:,:,:) = e3v_0_crs(:,:,:)
360     e3w_n_crs(:,:,:) = e3w_0_crs(:,:,:)
361
362     e3t_a_crs(:,:,:) = e3t_0_crs(:,:,:)
363     e3u_a_crs(:,:,:) = e3u_0_crs(:,:,:)
364     e3v_a_crs(:,:,:) = e3v_0_crs(:,:,:)
365     e3w_a_crs(:,:,:) = e3w_0_crs(:,:,:)
366#endif
367
368     !    3.d.3   Vertical depth (meters)
369     !cbr: il semblerait que p_e3=... ne soit pas utile ici !!!!!!!!!
370     CALL crs_dom_ope( gdept_0, 'MAX', 'T', tmask, gdept_0_crs, p_e3=zfse3t, psgn=1.0 ) 
371     CALL crs_dom_ope( gdepw_0, 'MAX', 'W', tmask, gdepw_0_crs, p_e3=zfse3w, psgn=1.0 )
372#if defined key_vvl
373     gdept_n_crs(:,:,:) = gdept_0_crs(:,:,:)
374     gdepw_n_crs(:,:,:) = gdepw_0_crs(:,:,:)
375#endif
376
377write(narea+200,*)"crsini8",nstop; call flush(narea+200)
378
379     !---------------------------------------------------------
380     ! 4.  Coarse grid ocean volume and averaging weights
381     !---------------------------------------------------------
382     CALL crs_dom_facvol( tmask, 'T', e1t, e2t, zfse3t, ocean_volume_crs_t, facvol_t )
383     !
384     bt_crs(:,:,:) = ocean_volume_crs_t(:,:,:) * facvol_t(:,:,:)
385     !
386     r1_bt_crs(:,:,:) = 0._wp 
387     WHERE( bt_crs /= 0._wp ) r1_bt_crs(:,:,:) = 1._wp / bt_crs(:,:,:)
388
389     CALL crs_dom_facvol( tmask, 'W', e1t, e2t, zfse3w, ocean_volume_crs_w, facvol_w )
390
391
392     !---------------------------------------------------------
393     ! 5.  Coarse grid ocean volume and averaging weights
394     !---------------------------------------------------------
395     !CALL dom_grid_crs   ! Save the parent grid information  & Switch to coarse grid domain
396     !CALL ldf_tra_crs_init
397     !CALL dom_grid_glo   ! Return to parent grid domain
398
399write(narea+200,*)"crsini9",nstop; call flush(narea+200)
400
401     !
402     !---------------------------------------------------------
403     ! 6.  Write out coarse meshmask  (see OPA_SRC/DOM/domwri.F90 for ideas later)
404     !---------------------------------------------------------
405
406     IF( nn_msh_crs > 0 ) THEN
407        CALL dom_grid_crs   ! Save the parent grid information  & Switch to coarse grid domain
408        CALL crs_dom_wri     
409        CALL dom_grid_glo   ! Return to parent grid domain
410     ENDIF
411   
412
413      rhop_crs(:,:,:)=0._wp ; rhd_crs(:,:,:)=0._wp ; rb2_crs(:,:,:)=0._wp
414
415write(narea+200,*)"crsini10",nstop; call flush(narea+200)
416 
417     !---------------------------------------------------------
418     ! 7. Finish and clean-up
419     !---------------------------------------------------------
420     CALL wrk_dealloc(jpi, jpj, jpk, zfse3t, zfse3u, zfse3v, zfse3w )
421
422      IF( nn_timing == 1 )  CALL timing_stop('crs_init')
423
424   END SUBROUTINE crs_init
425   
426   !!======================================================================
427
428END MODULE crsini
Note: See TracBrowser for help on using the repository browser.