Philippe 07/03/2019: IO bugfix: io_set_mnhversion must be called by all the processes
[MNH-git_open_source-lfs.git] / src / SURFEX / read_oceann.F90
1 !SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
2 !SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
3 !SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
4 !SFX_LIC for details. version 1.
5 !     #########
6       SUBROUTINE READ_OCEAN_n (DTCO, O, OR, U, HPROGRAM)
7 !     #########################################
8 !
9 !!****  *READ_OCEAN_n* - read oceanic variables
10 !!
11 !!
12 !!    PURPOSE
13 !!    -------
14 !!
15 !!**  METHOD
16 !!    ------
17 !!
18 !!    EXTERNAL
19 !!    --------
20 !!
21 !!
22 !!    IMPLICIT ARGUMENTS
23 !!    ------------------
24 !!
25 !!    REFERENCE
26 !!    ---------
27 !!
28 !!
29 !!    AUTHOR
30 !!    ------
31 !!      C. Lebeaupin Brossier   *Meteo France*
32 !!
33 !!    MODIFICATIONS
34 !!    -------------
35 !!      Original    04/2007 
36 !!      Modofied    07/2012, P. Le Moigne : CMO1D phasing
37 !-------------------------------------------------------------------------------
38 !
39 !*       0.    DECLARATIONS
40 !              ------------
41 !
42 !
43 !
44 !
45 USE MODD_DATA_COVER_n, ONLY : DATA_COVER_t
46 USE MODD_OCEAN_n, ONLY : OCEAN_t
47 USE MODD_OCEAN_REL_n, ONLY : OCEAN_REL_t
48 USE MODD_SURF_ATM_n, ONLY : SURF_ATM_t
49 !
50 USE MODD_SURF_PAR, ONLY : XUNDEF
51 USE MODD_OCEAN_GRID, ONLY : NOCKMIN,NOCKMAX,XZHOC
52 !
53 !
54 USE MODI_READ_SURF
55 USE MODI_OCEAN_MERCATORVERGRID
56 USE MODI_PREP_OCEAN_MERCATORVERGRID
57 !
58 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
59 USE PARKIND1  ,ONLY : JPRB
60 !
61 USE MODI_GET_TYPE_DIM_n
62 !
63 IMPLICIT NONE
64 !
65 !*       0.1   Declarations of arguments
66 !              -------------------------
67 !
68 !
69 TYPE(DATA_COVER_t), INTENT(INOUT) :: DTCO
70 TYPE(OCEAN_t), INTENT(INOUT) :: O
71 TYPE(OCEAN_REL_t), INTENT(INOUT) :: OR
72 TYPE(SURF_ATM_t), INTENT(INOUT) :: U
73 !
74  CHARACTER(LEN=6),  INTENT(IN)  :: HPROGRAM ! calling program
75 !
76 !*       0.2   Declarations of local variables
77 !              -------------------------------
78 !
79 INTEGER           :: ILU          ! 1D physical dimension
80 !
81 INTEGER           :: IRESP          ! Error code after redding
82 !
83  CHARACTER(LEN=4)  :: YLVL
84 !
85  CHARACTER(LEN=12) :: YRECFM         ! Name of the article to be read
86  CHARACTER(LEN=14) :: YFORM          ! Writing format
87 REAL, DIMENSION(:),ALLOCATABLE  :: ZWORK      ! 1D array to write data in file
88 !
89 INTEGER :: JLEVEL ! loop counter on oceanic levels
90 INTEGER :: J      ! loop counter on sea grid points
91 !
92 INTEGER           :: IVERSION       ! surface version
93 REAL(KIND=JPRB) :: ZHOOK_HANDLE
94 !
95 !-------------------------------------------------------------------------------
96 IF (LHOOK) CALL DR_HOOK('READ_OCEAN_N',0,ZHOOK_HANDLE)
97 O%NOCTCOUNT=0
98 !
99 YRECFM='VERSION'
100  CALL READ_SURF(HPROGRAM,YRECFM,IVERSION,IRESP)
101 !
102 !* flag to use or not Ocean model
103 !
104 IF (IVERSION<=3) THEN
105    O%LMERCATOR=.FALSE.
106 ELSE
107    YRECFM='SEA_OCEAN'
108    CALL READ_SURF(HPROGRAM,YRECFM,O%LMERCATOR,IRESP)
109 ENDIF
110 !
111 IF (.NOT. O%LMERCATOR) THEN
112   ALLOCATE(O%XSEAT(0,0))
113   ALLOCATE(O%XSEAS(0,0))  
114   ALLOCATE(O%XSEAU(0,0))
115   ALLOCATE(O%XSEAV(0,0))  
116   ALLOCATE(OR%XSEAT_REL(0,0))
117   ALLOCATE(OR%XSEAS_REL(0,0))
118   ALLOCATE(OR%XSEAU_REL(0,0))
119   ALLOCATE(OR%XSEAV_REL(0,0))
120   !
121   ALLOCATE(O%XSEAE(0,0))
122   ALLOCATE(O%XSEABATH(0,0))
123   ALLOCATE(O%XSEAHMO(0))
124   ALLOCATE(O%XLE        (0,0))
125   ALLOCATE(O%XLK        (0,0))
126   ALLOCATE(O%XKMEL      (0,0))
127   ALLOCATE(O%XKMELM     (0,0))
128   ALLOCATE(O%XSEATEND   (0))
129   ALLOCATE(O%XDTFSOL(0,0))
130   ALLOCATE(O%XDTFNSOL(0))
131   IF (LHOOK) CALL DR_HOOK('READ_OCEAN_N',1,ZHOOK_HANDLE)
132   RETURN
133 ENDIF
134 !
135 !-------------------------------------------------------------------------------
136 !
137 IF (IVERSION>=8) THEN
138
139   NOCKMIN = 0
140   YRECFM='SEA_NBLEVEL'
141   CALL READ_SURF(HPROGRAM,YRECFM,NOCKMAX,IRESP)
142   !
143   ALLOCATE(XZHOC(NOCKMIN:NOCKMAX))
144   XZHOC(NOCKMIN) = 0.
145   ! Read vertical grid
146   DO JLEVEL = NOCKMIN+1,NOCKMAX
147     WRITE(YLVL,'(I4)') JLEVEL
148     YRECFM='LEVL_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
149     CALL READ_SURF(HPROGRAM,YRECFM,XZHOC(JLEVEL),IRESP)
150   END DO
151   !
152 ELSE
153   !
154   NOCKMIN = 0
155   CALL PREP_OCEAN_MERCATORVERGRID(HPROGRAM,.TRUE.)
156   !
157 ENDIF
158 !
159  CALL OCEAN_MERCATORVERGRID
160 !
161 ! Relaxation time and logical
162 YRECFM='TAU_REL_OC'
163  CALL READ_SURF(HPROGRAM,YRECFM,OR%XTAU_REL,IRESP)
164 !
165 YRECFM='LREL_CUR_OC'
166  CALL READ_SURF(HPROGRAM,YRECFM,OR%LREL_CUR,IRESP)
167
168 YRECFM='LREL_TS_OC'
169  CALL READ_SURF(HPROGRAM,YRECFM,OR%LREL_TS,IRESP)
170 YRECFM='LFLX_NULL_OC'
171  CALL READ_SURF(HPROGRAM,YRECFM,OR%LFLUX_NULL,IRESP)
172 YRECFM='LFLX_CORR_OC'
173  CALL READ_SURF(HPROGRAM,YRECFM,OR%LFLX_CORR,IRESP)
174 YRECFM='CORR_FLX_OC'
175  CALL READ_SURF(HPROGRAM,YRECFM,OR%XQCORR,IRESP)
176 YRECFM='LDIAPYC_OC'
177  CALL READ_SURF(HPROGRAM,YRECFM,OR%LDIAPYCNAL,IRESP)
178 !
179 !* 1D physical dimension
180 !
181 YRECFM='SIZE_SEA'
182  CALL GET_TYPE_DIM_n(DTCO, U, 'SEA   ',ILU)
183 !
184 !*       2.     Prognostic fields:
185 !               -----------------
186 !
187 ALLOCATE(ZWORK(ILU))
188 !* oceanic temperature
189 !
190 ALLOCATE(O%XSEAT(ILU,NOCKMIN:NOCKMAX))
191 !
192 DO JLEVEL=NOCKMIN+1,NOCKMAX
193   WRITE(YLVL,'(I4)') JLEVEL
194   YRECFM='TEMP_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
195   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
196   O%XSEAT(:,JLEVEL)=ZWORK(:)
197 END DO
198 O%XSEAT(:,NOCKMIN)=O%XSEAT(:,NOCKMIN+1)
199 !
200 !* relaxation profile for oceanic temperature
201 !
202 ALLOCATE(OR%XSEAT_REL(ILU,NOCKMIN:NOCKMAX))
203 !
204 DO JLEVEL=NOCKMIN+1,NOCKMAX
205   WRITE(YLVL,'(I4)') JLEVEL
206   YRECFM='T_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
207   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
208   OR%XSEAT_REL(:,JLEVEL)=ZWORK(:)
209 END DO
210 OR%XSEAT_REL(:,NOCKMIN)=OR%XSEAT_REL(:,NOCKMIN+1)
211 !
212 !* oceanic salinity
213 !
214 ALLOCATE(O%XSEAS(ILU,NOCKMIN:NOCKMAX))
215 !
216 DO JLEVEL=NOCKMIN+1,NOCKMAX
217   WRITE(YLVL,'(I4)') JLEVEL
218   YRECFM='SALT_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
219   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
220   O%XSEAS(:,JLEVEL)=ZWORK(:)
221 END DO
222 O%XSEAS(:,NOCKMIN)=O%XSEAS(:,NOCKMIN+1)
223 !
224 !* oceanic salinity profile of relaxation
225 !
226 ALLOCATE(OR%XSEAS_REL(ILU,NOCKMIN:NOCKMAX))
227 !
228 DO JLEVEL=NOCKMIN+1,NOCKMAX
229   WRITE(YLVL,'(I4)') JLEVEL
230   YRECFM='S_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
231   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
232   OR%XSEAS_REL(:,JLEVEL)=ZWORK(:)
233 END DO
234 OR%XSEAS_REL(:,NOCKMIN)=OR%XSEAS_REL(:,NOCKMIN+1)
235 !
236 !* oceanic current
237 !
238 ALLOCATE(OR%XSEAU_REL(ILU,NOCKMIN:NOCKMAX))
239 ALLOCATE(OR%XSEAV_REL(ILU,NOCKMIN:NOCKMAX))
240 !
241 DO JLEVEL=NOCKMIN+1,NOCKMAX
242   WRITE(YLVL,'(I4)') JLEVEL
243   YRECFM='U_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
244   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
245   OR%XSEAU_REL(:,JLEVEL)=ZWORK(:)
246 END DO
247 OR%XSEAU_REL(:,NOCKMIN)=OR%XSEAU_REL(:,NOCKMIN+1)
248 !
249 DO JLEVEL=NOCKMIN+1,NOCKMAX
250   WRITE(YLVL,'(I4)') JLEVEL
251   YRECFM='V_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
252   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
253   OR%XSEAV_REL(:,JLEVEL)=ZWORK(:)
254 END DO
255 OR%XSEAV_REL(:,NOCKMIN)=OR%XSEAV_REL(:,NOCKMIN+1)
256 !
257 ALLOCATE(O%XSEAU(ILU,NOCKMIN:NOCKMAX))
258 ALLOCATE(O%XSEAV(ILU,NOCKMIN:NOCKMAX))
259 !
260 DO JLEVEL=NOCKMIN+1,NOCKMAX
261   WRITE(YLVL,'(I4)') JLEVEL
262   YRECFM='UCUR_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
263   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
264   O%XSEAU(:,JLEVEL)=ZWORK(:)
265 END DO
266 DO JLEVEL=NOCKMIN+1,NOCKMAX
267   WRITE(YLVL,'(I4)') JLEVEL
268   YRECFM='VCUR_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
269   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
270   O%XSEAV(:,JLEVEL)=ZWORK(:)
271 END DO
272 O%XSEAU(:,NOCKMIN)=O%XSEAU(:,NOCKMIN+1)
273 O%XSEAV(:,NOCKMIN)=O%XSEAV(:,NOCKMIN+1)
274 !
275 !* oceanic turbulent kinetic energy
276 !
277 ALLOCATE(O%XSEAE(ILU,NOCKMIN:NOCKMAX))
278 !
279 DO JLEVEL=NOCKMIN+1,NOCKMAX
280   WRITE(YLVL,'(I4)') JLEVEL
281   YRECFM='TKE_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
282   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
283   O%XSEAE(:,JLEVEL)=ZWORK(:)
284 END DO
285 O%XSEAE(:,NOCKMIN)=O%XSEAE(:,NOCKMIN+1)
286 !
287 !
288 !-------------------------------------------------------------------------------
289 !
290 !*       4.     Semi-prognostic fields:
291 !               ----------------------
292 !
293 !* bathymetry indice
294 !
295 ALLOCATE(O%XSEABATH(ILU,NOCKMIN:NOCKMAX))
296 !
297 DO JLEVEL=NOCKMIN+1,NOCKMAX
298   WRITE(YLVL,'(I4)') JLEVEL
299   YRECFM='SEAINDBATH'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
300   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
301   O%XSEABATH(:,JLEVEL)=ZWORK(:)
302 END DO
303 O%XSEABATH(:,NOCKMIN)=1.
304 !
305 !-------------------------------------------------------------------------------
306 !Complete undefined values for the oceanic 1D model convergence
307 DO J=1,ILU
308   DO JLEVEL=NOCKMIN+2,NOCKMAX
309     IF (O%XSEABATH(J,JLEVEL)==0.) THEN
310       O%XSEAT(J,JLEVEL)=O%XSEAT(J,JLEVEL-1)
311       O%XSEAS(J,JLEVEL)=O%XSEAS(J,JLEVEL-1)
312       O%XSEAU(J,JLEVEL)=O%XSEAU(J,JLEVEL-1)
313       O%XSEAV(J,JLEVEL)=O%XSEAV(J,JLEVEL-1)
314       O%XSEAE(J,JLEVEL)=O%XSEAE(J,JLEVEL-1)
315       !
316       OR%XSEAT_REL(J,JLEVEL)=OR%XSEAT_REL(J,JLEVEL-1)
317       OR%XSEAS_REL(J,JLEVEL)=OR%XSEAS_REL(J,JLEVEL-1)
318       OR%XSEAU_REL(J,JLEVEL)=OR%XSEAU_REL(J,JLEVEL-1)
319       OR%XSEAV_REL(J,JLEVEL)=OR%XSEAV_REL(J,JLEVEL-1)      
320     ENDIF
321   ENDDO
322 ENDDO
323 !
324 DEALLOCATE(ZWORK)
325 !-------------------------------------------------------------------------------
326 ALLOCATE(O%XSEAHMO(ILU))
327 YRECFM='SEA_HMO'
328  CALL READ_SURF(HPROGRAM,YRECFM,O%XSEAHMO(:),IRESP)
329 !
330 !-------------------------------------------------------------------------------
331 ALLOCATE(O%XLE        (SIZE(O%XSEAT,1),NOCKMIN:NOCKMAX))
332 ALLOCATE(O%XLK        (SIZE(O%XSEAT,1),NOCKMIN:NOCKMAX))
333 ALLOCATE(O%XKMEL      (SIZE(O%XSEAT,1),NOCKMIN:NOCKMAX))
334 ALLOCATE(O%XKMELM     (SIZE(O%XSEAT,1),NOCKMIN:NOCKMAX))
335 O%XLE(:,:)    =XUNDEF
336 O%XLK(:,:)    =XUNDEF
337 O%XKMEL(:,:)  =XUNDEF
338 O%XKMELM(:,:) =XUNDEF
339 !
340 ALLOCATE(O%XSEATEND   (SIZE(O%XSEAT,1)))
341 O%XSEATEND(:) =XUNDEF
342 !
343 ALLOCATE(O%XDTFSOL(ILU,NOCKMIN:NOCKMAX))
344 ALLOCATE(O%XDTFNSOL(ILU))
345 !
346 O%XDTFSOL(:,:) = XUNDEF 
347 O%XDTFNSOL(:) = XUNDEF 
348 !
349 IF (LHOOK) CALL DR_HOOK('READ_OCEAN_N',1,ZHOOK_HANDLE)
350 !
351 !------------------------------------------------------------------------------
352 END SUBROUTINE READ_OCEAN_n