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

source: branches/UKMO/dev_r10171_test_crs_AMM7/NEMOGCM/NEMO/OPA_SRC/CRS/crsdomwri.F90 @ 10207

Last change on this file since 10207 was 10207, checked in by cmao, 6 years ago

remove svn keyword
File size: 16.2 KB
Line 
1MODULE crsdomwri
2   !!======================================================================
3   !! Coarse Ocean initialization : write the coarse ocean domain mesh and mask files
4   !!======================================================================
5   !! History :  OPA  ! 1997-02  (G. Madec)  Original code
6   !!            8.1  ! 1999-11  (M. Imbard)  NetCDF FORMAT with IOIPSL
7   !!   NEMO     1.0  ! 2002-08  (G. Madec)  F90 and several file
8   !!            3.0  ! 2008-01  (S. Masson) add dom_uniq_crs
9   !!            4.0  ! 2011-01  (A. R. Porter, STFC Daresbury) dynamical allocation
10   !!                 ! 2012-06  (J. Simeon, C. Calone, C Ethe )  Reduced and modified for coarse grid
11   !!----------------------------------------------------------------------
12
13   !!----------------------------------------------------------------------
14   !!   crs_dom_wri    : create and write mesh and mask file(s)
15   !!----------------------------------------------------------------------
16   USE timing          ! Timing
17   USE dom_oce         ! ocean space and time domain
18   USE in_out_manager  ! I/O manager
19   USE par_kind, ONLY: wp
20   USE lib_mpp         ! MPP library
21   USE iom_def
22   USE iom
23   USE crs         ! coarse grid domain
24   USE crsdom         ! coarse grid domain
25   USE crslbclnk       ! crs mediator to lbclnk
26   USE wrk_nemo        ! Working array
27
28
29
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC crs_dom_wri        ! routine called by crsini.F90
34
35CONTAINS
36
37   SUBROUTINE crs_dom_wri
38      !!----------------------------------------------------------------------
39      !!                  ***  ROUTINE crs_dom_wri  ***
40      !!
41      !! ** Purpose :   Create the NetCDF file(s) which contain(s) all the
42      !!      ocean domain informations (mesh and mask arrays). This (these)
43      !!      file(s) is (are) used for visualisation (SAXO software) and
44      !!      diagnostic computation.
45      !!
46      !! ** Method  :   Write in a file all the arrays generated in routines
47      !!      crsini for meshes and mask. In three separate files:
48      !!      domain size, horizontal grid-point position,
49      !!      masks, depth and vertical scale factors
50      !!     
51      !! ** Output files :   mesh_hgr_crs.nc, mesh_zgr_crs.nc, mesh_mask.nc
52      !!----------------------------------------------------------------------
53      !!
54      INTEGER           ::   inum0    ! temprary units for 'mesh_mask.nc' file
55      INTEGER           ::   inum1    ! temprary units for 'mesh.nc'      file
56      INTEGER           ::   inum2    ! temprary units for 'mask.nc'      file
57      INTEGER           ::   inum3    ! temprary units for 'mesh_hgr.nc'  file
58      INTEGER           ::   inum4    ! temprary units for 'mesh_zgr.nc'  file
59      CHARACTER(len=21) ::   clnam0   ! filename (mesh and mask informations)
60      CHARACTER(len=21) ::   clnam1   ! filename (mesh informations)
61      CHARACTER(len=21) ::   clnam2   ! filename (mask informations)
62      CHARACTER(len=21) ::   clnam3   ! filename (horizontal mesh informations)
63      CHARACTER(len=21) ::   clnam4   ! filename (vertical   mesh informations)
64      INTEGER           ::   ji, jj, jk   ! dummy loop indices
65      INTEGER           ::   iji,ijj
66      !                                   !  workspaces
67      REAL(wp), POINTER, DIMENSION(:,:  ) :: zprt, zprw 
68      REAL(wp), POINTER, DIMENSION(:,:,:) :: zdepu, zdepv
69      REAL(wp), POINTER, DIMENSION(:,:  ) :: ze3tp, ze3wp
70      !!----------------------------------------------------------------------
71      !
72      IF( nn_timing == 1 )  CALL timing_start('crs_dom_wri')
73      !
74      CALL wrk_alloc( jpi_crs, jpj_crs,      zprt , zprw  )
75      CALL wrk_alloc( jpi_crs, jpj_crs,      ze3tp, ze3wp )
76      CALL wrk_alloc( jpi_crs, jpj_crs, jpk, zdepu, zdepv )
77
78      ze3tp(:,:) = 0.0
79      ze3wp(:,:) = 0.0
80
81      !
82      IF(lwp) WRITE(numout,*)
83      IF(lwp) WRITE(numout,*) 'crs_dom_wri : create NetCDF mesh and mask information file(s)'
84      IF(lwp) WRITE(numout,*) '~~~~~~~'
85     
86      clnam0 = 'mesh_mask_crs'  ! filename (mesh and mask informations)
87      clnam1 = 'mesh_crs'       ! filename (mesh informations)
88      clnam2 = 'mask_crs'       ! filename (mask informations)
89      clnam3 = 'mesh_hgr_crs'   ! filename (horizontal mesh informations)
90      clnam4 = 'mesh_zgr_crs'   ! filename (vertical   mesh informations)
91     
92
93      SELECT CASE ( MOD(nn_msh_crs, 3) )
94         !                                  ! ============================
95      CASE ( 1 )                            !  create 'mesh_mask.nc' file
96         !                                  ! ============================
97         CALL iom_open( TRIM(clnam0), inum0, ldwrt = .TRUE., kiolib = jprstlib )
98         inum2 = inum0                                            ! put all the informations
99         inum3 = inum0                                            ! in unit inum0
100         inum4 = inum0
101         
102         !                                  ! ============================
103      CASE ( 2 )                            !  create 'mesh.nc' and
104         !                                  !         'mask.nc' files
105         !                                  ! ============================
106         CALL iom_open( TRIM(clnam1), inum1, ldwrt = .TRUE., kiolib = jprstlib )
107         CALL iom_open( TRIM(clnam2), inum2, ldwrt = .TRUE., kiolib = jprstlib )
108         inum3 = inum1                                            ! put mesh informations
109         inum4 = inum1                                            ! in unit inum1
110         !                                  ! ============================
111      CASE ( 0 )                            !  create 'mesh_hgr.nc'
112         !                                  !         'mesh_zgr.nc' and
113         !                                  !         'mask.nc'     files
114         !                                  ! ============================
115         CALL iom_open( TRIM(clnam2), inum2, ldwrt = .TRUE., kiolib = jprstlib )
116         CALL iom_open( TRIM(clnam3), inum3, ldwrt = .TRUE., kiolib = jprstlib )
117         CALL iom_open( TRIM(clnam4), inum4, ldwrt = .TRUE., kiolib = jprstlib )
118         !
119      END SELECT
120 
121      !========================================================
122      !                                                         ! masks (inum2)
123
124      CALL iom_rstput( 0, 0, inum2, 'tmask', tmask_crs, ktype = jp_i1 )     !    ! land-sea mask
125      CALL iom_rstput( 0, 0, inum2, 'umask', umask_crs, ktype = jp_i1 )
126      CALL iom_rstput( 0, 0, inum2, 'vmask', vmask_crs, ktype = jp_i1 )
127      CALL iom_rstput( 0, 0, inum2, 'fmask', fmask_crs, ktype = jp_i1 )
128      CALL iom_rstput( 0, 0, inum2, 'tmaskutil', tmask_i_crs, ktype = jp_i1 )
129                                   !    ! unique point mask
130      CALL dom_uniq_crs( zprw, 'U' )
131      zprt = umask_crs(:,:,1) * zprw
132      CALL iom_rstput( 0, 0, inum2, 'umaskutil', zprt, ktype = jp_i1 ) 
133      CALL dom_uniq_crs( zprw, 'V' )
134      zprt = vmask_crs(:,:,1) * zprw
135      CALL iom_rstput( 0, 0, inum2, 'vmaskutil', zprt, ktype = jp_i1 ) 
136      CALL dom_uniq_crs( zprw, 'F' )
137      zprt = fmask_crs(:,:,1) * zprw
138      CALL iom_rstput( 0, 0, inum2, 'fmaskutil', zprt, ktype = jp_i1 ) 
139      !========================================================
140      !                                                         ! horizontal mesh (inum3)
141      CALL iom_rstput( 0, 0, inum3, 'glamt', glamt_crs, ktype = jp_r4 )     !    ! latitude
142      CALL iom_rstput( 0, 0, inum3, 'glamu', glamu_crs, ktype = jp_r4 )
143      CALL iom_rstput( 0, 0, inum3, 'glamv', glamv_crs, ktype = jp_r4 )
144      CALL iom_rstput( 0, 0, inum3, 'glamf', glamf_crs, ktype = jp_r4 )
145     
146      CALL iom_rstput( 0, 0, inum3, 'gphit', gphit_crs, ktype = jp_r4 )     !    ! longitude
147      CALL iom_rstput( 0, 0, inum3, 'gphiu', gphiu_crs, ktype = jp_r4 )
148      CALL iom_rstput( 0, 0, inum3, 'gphiv', gphiv_crs, ktype = jp_r4 )
149      CALL iom_rstput( 0, 0, inum3, 'gphif', gphif_crs, ktype = jp_r4 )
150     
151      CALL iom_rstput( 0, 0, inum3, 'e1t', e1t_crs, ktype = jp_r8 )         !    ! e1 scale factors
152      CALL iom_rstput( 0, 0, inum3, 'e1u', e1u_crs, ktype = jp_r8 )
153      CALL iom_rstput( 0, 0, inum3, 'e1v', e1v_crs, ktype = jp_r8 )
154      CALL iom_rstput( 0, 0, inum3, 'e1f', e1f_crs, ktype = jp_r8 )
155     
156      CALL iom_rstput( 0, 0, inum3, 'e2t', e2t_crs, ktype = jp_r8 )         !    ! e2 scale factors
157      CALL iom_rstput( 0, 0, inum3, 'e2u', e2u_crs, ktype = jp_r8 )
158      CALL iom_rstput( 0, 0, inum3, 'e2v', e2v_crs, ktype = jp_r8 )
159      CALL iom_rstput( 0, 0, inum3, 'e2f', e2f_crs, ktype = jp_r8 )
160     
161      CALL iom_rstput( 0, 0, inum3, 'ff', ff_crs, ktype = jp_r8 )           !    ! coriolis factor
162
163      !========================================================
164      !                                                         ! vertical mesh (inum4)
165!     ! note that mbkt is set to 1 over land ==> use surface tmask_crs
166      zprt(:,:) = tmask_crs(:,:,1) * REAL( mbkt_crs(:,:) , wp )
167      CALL iom_rstput( 0, 0, inum4, 'mbathy', zprt, ktype = jp_i2 )     !    ! nb of ocean T-points
168
169      IF( ln_zps ) THEN                       ! z-coordinate - partial steps
170
171           
172         IF ( nn_msh_crs <= 6 ) THEN
173            CALL iom_rstput( 0, 0, inum4, 'e3t', e3t_0_crs )     
174            CALL iom_rstput( 0, 0, inum4, 'e3w', e3w_0_crs )     
175            CALL iom_rstput( 0, 0, inum4, 'e3u', e3u_0_crs )     
176            CALL iom_rstput( 0, 0, inum4, 'e3v', e3v_0_crs )     
177            CALL iom_rstput( 0, 0, inum4, 'e3t_max_crs', e3t_max_0_crs )     
178            CALL iom_rstput( 0, 0, inum4, 'e3w_max_crs', e3w_max_0_crs )     
179         ELSE
180            DO jj = 1,jpj_crs   
181               DO ji = 1,jpi_crs
182                  ze3tp(ji,jj) = e3t_0_crs(ji,jj,mbkt_crs(ji,jj)) * tmask_crs(ji,jj,1)
183                  ze3wp(ji,jj) = e3w_0_crs(ji,jj,mbkt_crs(ji,jj)) * tmask_crs(ji,jj,1)
184               END DO
185            END DO
186
187            CALL crs_lbc_lnk( ze3tp,'T', 1.0 )
188            CALL crs_lbc_lnk( ze3wp,'W', 1.0 )
189 
190            CALL iom_rstput( 0, 0, inum4, 'e3t_ps', ze3tp )     
191            CALL iom_rstput( 0, 0, inum4, 'e3w_ps', ze3wp )
192         ENDIF
193
194         IF ( nn_msh_crs <= 3 ) THEN
195            CALL iom_rstput( 0, 0, inum4, 'gdept', gdept_0_crs, ktype = jp_r4 ) 
196            DO jk = 1,jpk   
197               DO jj = 1, jpj_crsm1   
198                  DO ji = 1, jpi_crsm1  ! jes what to do for fs_jpim1??vector opt.
199                     zdepu(ji,jj,jk) = MIN( gdept_0_crs(ji,jj,jk) , gdept_0_crs(ji+1,jj  ,jk) ) * umask_crs(ji,jj,jk)
200                     zdepv(ji,jj,jk) = MIN( gdept_0_crs(ji,jj,jk) , gdept_0_crs(ji  ,jj+1,jk) ) * vmask_crs(ji,jj,jk)
201                  END DO   
202               END DO   
203            END DO
204
205            CALL crs_lbc_lnk( zdepu,'U', 1. )   ;   CALL crs_lbc_lnk( zdepv,'V', 1. ) 
206            CALL iom_rstput( 0, 0, inum4, 'gdepu', zdepu, ktype = jp_r4 )
207            CALL iom_rstput( 0, 0, inum4, 'gdepv', zdepv, ktype = jp_r4 )
208            CALL iom_rstput( 0, 0, inum4, 'gdepw', gdepw_0_crs, ktype = jp_r4 )
209         ELSE
210            DO jj = 1,jpj_crs   
211               DO ji = 1,jpi_crs
212                  zprt(ji,jj) = gdept_0(ji,jj,mbkt(ji,jj)  ) * tmask(ji,jj,1)
213                  zprw(ji,jj) = gdepw_0(ji,jj,mbkt(ji,jj)+1) * tmask(ji,jj,1)
214               END DO
215            END DO
216            CALL iom_rstput( 0, 0, inum4, 'hdept', zprt, ktype = jp_r4 )     
217            CALL iom_rstput( 0, 0, inum4, 'hdepw', zprw, ktype = jp_r4 ) 
218         ENDIF
219
220         CALL iom_rstput( 0, 0, inum4, 'gdept_1d', gdept_1d )     !    ! reference z-coord.
221         CALL iom_rstput( 0, 0, inum4, 'gdepw_1d', gdepw_1d )
222         CALL iom_rstput( 0, 0, inum4, 'e3t_1d'  , e3t_1d   )
223         CALL iom_rstput( 0, 0, inum4, 'e3w_1d'  , e3w_1d   )
224
225         CALL iom_rstput(  0, 0, inum4, 'ocean_volume_t', ocean_volume_crs_t ) 
226         CALL iom_rstput(  0, 0, inum4, 'facvol_t' , facvol_t  ) 
227         CALL iom_rstput(  0, 0, inum4, 'facvol_w' , facvol_w  ) 
228         CALL iom_rstput(  0, 0, inum4, 'facsurfu' , facsurfu  ) 
229         CALL iom_rstput(  0, 0, inum4, 'facsurfv' , facsurfv  ) 
230         CALL iom_rstput(  0, 0, inum4, 'e1e2w_msk', e1e2w_msk ) 
231         CALL iom_rstput(  0, 0, inum4, 'e2e3u_msk', e2e3u_msk ) 
232         CALL iom_rstput(  0, 0, inum4, 'e1e3v_msk', e1e3v_msk )
233         CALL iom_rstput(  0, 0, inum4, 'e1e2w'    , e1e2w_crs ) 
234         CALL iom_rstput(  0, 0, inum4, 'e2e3u'    , e2e3u_crs ) 
235         CALL iom_rstput(  0, 0, inum4, 'e1e3v'    , e1e3v_crs )
236         CALL iom_rstput(  0, 0, inum4, 'bt'       , bt_crs    )
237         CALL iom_rstput(  0, 0, inum4, 'r1_bt'    , r1_bt_crs )
238
239         CALL iom_rstput(  0, 0, inum4, 'crs_surfu_wgt', crs_surfu_wgt ) 
240         CALL iom_rstput(  0, 0, inum4, 'crs_surfv_wgt', crs_surfv_wgt ) 
241         CALL iom_rstput(  0, 0, inum4, 'crs_volt_wgt' , crs_volt_wgt  ) 
242
243      ENDIF
244     
245     IF( ln_zco ) THEN
246         !                                                      ! z-coordinate - full steps
247        CALL iom_rstput( 0, 0, inum4, 'gdept_1d', gdept_1d )     !    ! depth
248        CALL iom_rstput( 0, 0, inum4, 'gdepw_1d', gdepw_1d )
249        CALL iom_rstput( 0, 0, inum4, 'e3t_1d'  , e3t_1d   )     !    ! scale factors
250        CALL iom_rstput( 0, 0, inum4, 'e3w_1d'  , e3w_1d   )
251     ENDIF
252      !                                     ! ============================
253      !                                     !        close the files
254      !                                     ! ============================
255      SELECT CASE ( MOD(nn_msh_crs, 3) )
256      CASE ( 1 )               
257         CALL iom_close( inum0 )
258      CASE ( 2 )
259         CALL iom_close( inum1 )
260         CALL iom_close( inum2 )
261      CASE ( 0 )
262         CALL iom_close( inum2 )
263         CALL iom_close( inum3 )
264         CALL iom_close( inum4 )
265      END SELECT
266      !
267      CALL wrk_dealloc( jpi_crs, jpj_crs,      zprt , zprw  )
268      CALL wrk_dealloc( jpi_crs, jpj_crs,      ze3tp, ze3wp )
269      CALL wrk_dealloc( jpi_crs, jpj_crs, jpk, zdepu, zdepv )
270      !
271      IF( nn_timing == 1 )  CALL timing_stop('crs_dom_wri')
272      !
273     
274   END SUBROUTINE crs_dom_wri
275
276
277   SUBROUTINE dom_uniq_crs( puniq, cdgrd )
278      !!----------------------------------------------------------------------
279      !!                  ***  ROUTINE crs_dom_uniq_crs  ***
280      !!                   
281      !! ** Purpose :   identify unique point of a grid (TUVF)
282      !!
283      !! ** Method  :   1) apply crs_lbc_lnk on an array with different values for each element
284      !!                2) check which elements have been changed
285      !!----------------------------------------------------------------------
286      !
287      CHARACTER(len=1)        , INTENT(in   ) ::   cdgrd   !
288      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   puniq   !
289      !
290      REAL(wp) ::  zshift   ! shift value link to the process number
291      INTEGER  ::  ji       ! dummy loop indices
292      LOGICAL, DIMENSION(SIZE(puniq,1),SIZE(puniq,2),1) ::  lldbl  ! store whether each point is unique or not
293      REAL(wp), POINTER, DIMENSION(:,:) :: ztstref
294      !!----------------------------------------------------------------------
295      !
296      IF( nn_timing == 1 )  CALL timing_start('crs_dom_uniq_crs')
297      !
298      CALL wrk_alloc( jpi_crs, jpj_crs, ztstref )
299      !
300      ! build an array with different values for each element
301      ! in mpp: make sure that these values are different even between process
302      ! -> apply a shift value according to the process number
303      zshift = jpi_crs * jpj_crs * ( narea - 1 )
304      ztstref(:,:) = RESHAPE( (/ (zshift + REAL(ji,wp), ji = 1, jpi_crs*jpj_crs) /), (/ jpi_crs, jpj_crs /) )
305      !
306      puniq(:,:) = ztstref(:,:)                   ! default definition
307      CALL crs_lbc_lnk( puniq,cdgrd, 1. )            ! apply boundary conditions
308      lldbl(:,:,1) = puniq(:,:) == ztstref(:,:)   ! check which values have been changed
309      !
310      puniq(:,:) = 1.                             ! default definition
311      ! fill only the inner part of the cpu with llbl converted into real
312      puniq(nldi_crs:nlei_crs,nldj_crs:nlej_crs) = REAL( COUNT( lldbl(nldi_crs:nlei_crs,nldj_crs:nlej_crs,:), dim = 3 ) , wp )
313      !
314      CALL wrk_dealloc( jpi_crs, jpj_crs, ztstref )
315      !
316      IF( nn_timing == 1 )  CALL timing_stop('crs_dom_uniq_crs')
317      !
318     
319   END SUBROUTINE dom_uniq_crs
320
321   !!======================================================================
322
323END MODULE crsdomwri
324
325
Note: See TracBrowser for help on using the repository browser.