Juan 24/05/2017: add MPPDB_CHECK on SURFEX routine for reprod check with MNH_PARALLEL key
[MNH-git_open_source-lfs.git] / src / SURFEX / pgd_grid.F90
1 !SURFEX_LIC Copyright 1994-2014 Meteo-France 
2 !SURFEX_LIC This is part of the SURFEX software governed by the CeCILL-C  licence
3 !SURFEX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
4 !SURFEX_LIC for details. version 1.
5 !     #########
6       SUBROUTINE PGD_GRID(HPROGRAM,HFILE,HFILETYPE,OGRID,HGRID,KGRID_PAR,PGRID_PAR)
7 !     ##########################################################
8 !!
9 !!    PURPOSE
10 !!    -------
11 !!   Reads in namelist the grid type and parameters.
12 !!
13 !!    METHOD
14 !!    ------
15 !!   
16 !!    EXTERNAL
17 !!    --------
18 !!
19 !!
20 !!    IMPLICIT ARGUMENTS
21 !!    ------------------
22 !!
23 !!
24 !!    REFERENCE
25 !!    ---------
26 !!
27 !!    AUTHOR
28 !!    ------
29 !!
30 !!    V. Masson                   Meteo-France
31 !!
32 !!    MODIFICATION
33 !!    ------------
34 !!
35 !!    Original     01/2004
36 !!    E. Martin    10/2007 IGN grid
37 !!    M. Moge      05/02/2015 parallelization (using local sizes, GET_MEAN_OF_COORD_SQRT_ll, SET_NAM_GRID_CONF_PROJ_LOCAL) + MPPDB_CHECK
38 !!    M. Moge      01/03/2015 call SPLIT_GRID if CPROGRAM == 'PGD   ' + remove SET_NAM_GRID_CONF_PROJ_LOCAL
39 !!    M. Moge      01/03/2015 change in the input arguments of PGD_GRID_IO_INIT : passing IDXRATIO, IDYRATIO
40 !----------------------------------------------------------------------------
41 !
42 !*    0.     DECLARATION
43 !            -----------
44 !
45 USE MODD_SURFEX_MPI,     ONLY : NSIZE, NINDEX, NPIO, NRANK
46 USE MODD_SURFEX_OMP,     ONLY : NINDX2, NWORK, XWORK, XWORK2, XWORK3, &
47                                 NWORK_FULL, XWORK_FULL, XWORK2_FULL
48 !
49 USE MODD_PGD_GRID,       ONLY : NL, XGRID_PAR, NGRID_PAR, XMESHLENGTH
50 USE MODN_PGD_GRID
51 USE MODD_SURF_ATM_GRID_n, ONLY : XLAT, XLON, XMESH_SIZE, XJPDIR
52 USE MODD_SURF_ATM_n,     ONLY : NDIM_FULL, NSIZE_FULL
53 USE MODD_CSTS,           ONLY : XPI, XRADIUS
54 !
55 USE MODI_DEFAULT_GRID
56 USE MODI_GRID_FROM_FILE
57 USE MODI_OPEN_NAMELIST
58 USE MODI_TEST_NAM_VAR_SURF
59 USE MODI_CLOSE_NAMELIST
60 USE MODI_GET_LUOUT
61 USE MODI_READ_NAM_GRIDTYPE
62 USE MODI_LATLON_GRID
63 !
64 USE MODE_POS_SURF
65 !
66 !
67 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
68 USE PARKIND1  ,ONLY : JPRB
69 !
70 USE MODI_ABOR1_SFX
71 !
72 USE MODI_PGD_GRID_IO_INIT
73 #ifdef MNH_PARALLEL
74 USE MODE_TOOLS_ll, ONLY : GET_MEAN_OF_COORD_SQRT_ll
75 !
76 USE MODI_GET_SIZE_FULL_n
77 #ifdef MNH_PARALLEL
78 USE MODE_MPPDB
79 #endif
80 USE MODI_SPLIT_GRID
81 USE MODD_CONF, ONLY : CPROGRAM
82 #endif
83 !
84 IMPLICIT NONE
85 !
86 !*    0.1    Declaration of dummy arguments
87 !            ------------------------------
88 !
89  CHARACTER(LEN=6),  INTENT(IN)   :: HPROGRAM   ! program calling the surface
90  CHARACTER(LEN=28), INTENT(IN)   :: HFILE      ! atmospheric file name
91  CHARACTER(LEN=6),  INTENT(IN)   :: HFILETYPE  ! atmospheric file type
92 LOGICAL,           INTENT(IN)   :: OGRID      ! .true. if grid is imposed by atm. model
93  CHARACTER(LEN=10), INTENT(OUT)  :: HGRID      ! grid type
94 INTEGER,           INTENT(OUT)  :: KGRID_PAR  ! size of PGRID_PAR
95 REAL, DIMENSION(:), POINTER     :: PGRID_PAR  ! parameters defining this grid
96 !
97 !
98 !*    0.2    Declaration of local variables
99 !            ------------------------------
100 !
101 INTEGER           :: ILUOUT ! output listing logical unit
102 INTEGER           :: ILUNAM ! namelist file  logical unit
103 LOGICAL           :: GFOUND ! Flag true if namelist is present
104 REAL(KIND=JPRB) :: ZHOOK_HANDLE
105 INTEGER :: IIMAX_ll, IJMAX_ll ! global size of son model
106 LOGICAL :: GRECT
107 !
108  INTEGER :: IXOR = 1            ! position of modified bottom left point
109  INTEGER :: IYOR = 1            ! according to initial grid
110  INTEGER :: IXSIZE = -999       ! number of grid meshes in initial grid to be
111  INTEGER :: IYSIZE = -999       ! covered by the modified grid
112  INTEGER :: IDXRATIO = 1        ! resolution ratio between modified grid
113  INTEGER :: IDYRATIO = 1        ! and initial grid
114 NAMELIST/NAM_INIFILE_CONF_PROJ/IXOR,IYOR,IXSIZE,IYSIZE,IDXRATIO,IDYRATIO
115 !
116 !*    0.3    Declaration of namelists
117 !            ------------------------
118 !
119 !------------------------------------------------------------------------------
120 !
121 !*    1.      Defaults
122 !             --------
123 !
124 IF (LHOOK) CALL DR_HOOK('PGD_GRID',0,ZHOOK_HANDLE)
125  CALL DEFAULT_GRID(HPROGRAM,CGRID)
126 !
127 YINIFILE  = '                         '
128 YFILETYPE = '      '
129 !
130 IF (OGRID) THEN
131   YINIFILE  = HFILE
132   YFILETYPE = HFILETYPE
133 END IF
134 !
135  CALL GET_LUOUT(HPROGRAM,ILUOUT)
136 !------------------------------------------------------------------------------
137 !
138 !*    2.      Open namelist
139 !             -------------
140 !
141 IF (.NOT. OGRID) THEN
142   CALL OPEN_NAMELIST(HPROGRAM,ILUNAM)
143 !
144 !------------------------------------------------------------------------------
145 !
146 !*    3.      Read grid type
147 !             --------------
148 !
149   CALL POSNAM(ILUNAM,'NAM_PGD_GRID',GFOUND,ILUOUT)
150   IF (GFOUND) READ(UNIT=ILUNAM,NML=NAM_PGD_GRID)
151 !
152 !------------------------------------------------------------------------------
153 !
154 !*    5.      Close namelist
155 !             --------------
156 !
157   CALL CLOSE_NAMELIST(HPROGRAM,ILUNAM)
158 !
159 END IF
160 !-------------------------------------------------------------------------------
161 !
162 !*       4.    check of grid and input file types
163 !              ----------------------------------
164
165  CALL TEST_NAM_VAR_SURF(ILUOUT,'CGRID',CGRID,'CONF PROJ ','NONE      ','LONLAT REG','CARTESIAN ','GAUSS     ',&
166           'IGN       ','LONLATVAL ')  
167  CALL TEST_NAM_VAR_SURF(ILUOUT,'YFILETYPE',YFILETYPE,'      ','MESONH','LFI   ','ASCII ')
168 !
169 !
170 !------------------------------------------------------------------------------
171 !
172 !*    5.      Initializes grid characteristics
173 !             --------------------------------
174 !
175 !*    5.1     From another file
176 !             -----------------
177 !
178 IF (LEN_TRIM(YFILETYPE)>0 .AND. LEN_TRIM(YINIFILE)>0 ) THEN
179   IF (YFILETYPE=='MESONH' .OR. YFILETYPE=='LFI   ' .OR. YFILETYPE=='ASCII ') THEN
180     CALL GRID_FROM_FILE(HPROGRAM,YINIFILE,YFILETYPE,OGRID,CGRID,NGRID_PAR,XGRID_PAR,NL)
181     HGRID     = CGRID
182     IF ( HGRID == "IGN       " .OR. HGRID == "GAUSS     " .OR. HGRID == "NONE      " ) THEN
183       GRECT = .FALSE.
184     ELSE
185       GRECT = .TRUE.
186     ENDIF
187     ! on lit la taille globale du modele fils dans la namelist
188     CALL OPEN_NAMELIST(HPROGRAM,ILUNAM)
189     CALL POSNAM(ILUNAM,'NAM_INIFILE_CONF_PROJ',GFOUND,ILUOUT)
190     IF (GFOUND) THEN 
191       READ(UNIT=ILUNAM,NML=NAM_INIFILE_CONF_PROJ)
192       IIMAX_ll = IXSIZE*IDXRATIO
193       IJMAX_ll = IYSIZE*IDYRATIO
194     ENDIF
195     CALL CLOSE_NAMELIST(HPROGRAM,ILUNAM)
196
197     !*    3.      Additional actions for I/O
198     !
199     IF (GFOUND) THEN 
200 #ifdef MNH_PARALLEL
201       CALL PGD_GRID_IO_INIT(HPROGRAM,NGRID_PAR,XGRID_PAR, HGRID, GRECT, IIMAX_ll, IJMAX_ll, IDXRATIO, IDYRATIO)
202 #else
203       CALL PGD_GRID_IO_INIT(HPROGRAM)
204 #endif
205       NDIM_FULL = NL
206     ELSE
207 #ifdef MNH_PARALLEL
208       CALL PGD_GRID_IO_INIT(HPROGRAM,NGRID_PAR,XGRID_PAR, HGRID, GRECT)
209 #else
210       CALL PGD_GRID_IO_INIT(HPROGRAM)
211 #endif
212     ENDIF
213     NSIZE     = NDIM_FULL
214 #ifdef MNH_PARALLEL
215     CALL GET_SIZE_FULL_n(HPROGRAM,NDIM_FULL,NSIZE_FULL)
216     NL = NSIZE_FULL
217 #else
218     NSIZE_FULL = NL
219 #endif
220   ELSE
221     CALL ABOR1_SFX('PGD_GRID: FILE TYPE NOT SUPPORTED '//HFILETYPE//' FOR FILE '//HFILE)
222   END IF
223   !we don't need to call SPLIT_GRID, the grid has been splitted in GRID_FROM_FILE
224 !
225 ELSE
226 !
227 !*    5.2     Grid not initialized
228 !             --------------------
229 !
230   IF (CGRID=='NONE      ' .OR. CGRID=='          ') THEN
231     CALL ABOR1_SFX('PGD_GRID: GRID TYPE NOT INITIALIZED, CGRID='//CGRID)
232
233 !
234 !*    5.3     Grid initialized
235 !             ----------------
236 !
237   ELSE
238 !
239     CALL READ_NAM_GRIDTYPE(HPROGRAM,CGRID,NGRID_PAR,XGRID_PAR,NL)
240     HGRID     = CGRID
241     !*    3.      Additional actions for I/O
242     !
243 #ifdef MNH_PARALLEL
244     CALL PGD_GRID_IO_INIT(HPROGRAM,NGRID_PAR,XGRID_PAR)
245 #else
246     CALL PGD_GRID_IO_INIT(HPROGRAM)
247 #endif
248     NDIM_FULL = NL
249     NSIZE     = NDIM_FULL
250 #ifdef MNH_PARALLEL
251     CALL GET_SIZE_FULL_n(HPROGRAM,NDIM_FULL,NSIZE_FULL)
252     NL = NSIZE_FULL
253 #else
254     NSIZE_FULL = NL
255 #endif
256   END IF
257 #ifdef MNH_PARALLEL
258   ! IF we are in PREP_PGD, we need to split the grid. Otherwise, the grid was read in parallel and is already splitted
259   IF ( CPROGRAM == 'PGD   ') THEN
260     CALL SPLIT_GRID('MESONH',NGRID_PAR,XGRID_PAR)
261   ENDIF
262 #endif
263
264 END IF
265 !
266 !
267 IF (.NOT.ALLOCATED(NINDEX)) THEN
268   ALLOCATE(NINDEX(NDIM_FULL))
269   NINDEX(:) = 0
270 ENDIF
271 NINDX2    = NDIM_FULL
272 ALLOCATE(NWORK(NDIM_FULL))
273 ALLOCATE(XWORK(NDIM_FULL))
274 ALLOCATE(XWORK2(NDIM_FULL,10))
275 ALLOCATE(XWORK3(NDIM_FULL,10,10))
276 IF (NRANK==NPIO) THEN
277   ALLOCATE(NWORK_FULL(NDIM_FULL))
278   ALLOCATE(XWORK_FULL(NDIM_FULL))
279   ALLOCATE(XWORK2_FULL(NDIM_FULL,10))
280 ELSE
281   ALLOCATE(NWORK_FULL(0))
282   ALLOCATE(XWORK_FULL(0))
283   ALLOCATE(XWORK2_FULL(0,0))
284 ENDIF
285 !
286 KGRID_PAR = NGRID_PAR
287 ALLOCATE(PGRID_PAR(KGRID_PAR))
288 PGRID_PAR = XGRID_PAR
289 !
290 !------------------------------------------------------------------------------
291 !
292 !*    6.      Latitude and longitude
293 !             ----------------------
294 !
295 ALLOCATE(XLAT       (NSIZE_FULL))
296 ALLOCATE(XLON       (NSIZE_FULL))
297 ALLOCATE(XMESH_SIZE (NSIZE_FULL))
298 ALLOCATE(XJPDIR     (NSIZE_FULL))
299  CALL LATLON_GRID(CGRID,NGRID_PAR,NSIZE_FULL,ILUOUT,XGRID_PAR,XLAT,XLON,XMESH_SIZE,XJPDIR)
300 #ifdef MNH_PARALLEL
301  CALL MPPDB_CHECK_SURFEX2D(XLAT,"PGD_GRID after LATLON_GRID:XLAT",PRECISION,ILUOUT)
302  CALL MPPDB_CHECK_SURFEX2D(XLON,"PGD_GRID after LATLON_GRID:XLON",PRECISION,ILUOUT)
303  CALL MPPDB_CHECK_SURFEX2D(XMESH_SIZE,"PGD_GRID after LATLON_GRID:XMESH_SIZE",PRECISION,ILUOUT)
304 #endif
305 !
306 !------------------------------------------------------------------------------
307 !
308 !*    7.      Average grid length (in degrees)
309 !             --------------------------------
310 !
311 !* in meters
312 #ifdef MNH_PARALLEL
313 CALL GET_MEAN_OF_COORD_SQRT_ll(XMESH_SIZE,NSIZE_FULL,NDIM_FULL,XMESHLENGTH)
314 #else
315 XMESHLENGTH = SUM ( SQRT(XMESH_SIZE) ) / NL
316 #endif
317 !
318 !* in degrees (of latitude)
319 XMESHLENGTH = XMESHLENGTH *180. / XPI / XRADIUS
320 IF (LHOOK) CALL DR_HOOK('PGD_GRID',1,ZHOOK_HANDLE)
321 !
322 !-------------------------------------------------------------------------------
323 !
324 END SUBROUTINE PGD_GRID