source: branches/UKMO/dev_isf_eorca1_dynldfc3d_UKESM_GO6package_r9321/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_c3d.h90 @ 11257

Last change on this file since 11257 was 11257, checked in by mathiot, 15 months ago

add scaling depending on grid size for eorca1 south of -75N

File size: 18.9 KB
Line 
1   !!----------------------------------------------------------------------
2   !!                        ***  ldfdyn_c3d.h90  ***
3   !!----------------------------------------------------------------------
4
5   !!----------------------------------------------------------------------
6   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
7   !! $Id$
8   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
9   !!----------------------------------------------------------------------
10
11   !!----------------------------------------------------------------------
12   !!   'key_dynldf_c3d'             3D lateral eddy viscosity coefficients
13   !!----------------------------------------------------------------------
14
15   SUBROUTINE ldf_dyn_c3d( ld_print )
16      !!----------------------------------------------------------------------
17      !!                  ***  ROUTINE ldf_dyn_c3d  ***
18      !!                   
19      !! ** Purpose :   initializations of the horizontal ocean physics
20      !!
21      !! ** Method  :   3D eddy viscosity coef. ( longitude, latitude, depth )
22      !!       laplacian operator   : ahm1, ahm2 defined at T- and F-points
23      !!                              ahm2, ahm4 never used
24      !!       bilaplacian operator : ahm1, ahm2 never used
25      !!                           :  ahm3, ahm4 defined at U- and V-points
26      !!       ??? explanation of the default is missing
27      !!----------------------------------------------------------------------
28      USE ldftra_oce, ONLY :   aht0
29      USE iom
30      !!
31      LOGICAL, INTENT (in) ::   ld_print   ! If true, output arrays on numout
32      !!
33      INTEGER  ::   ji, jj, jk   ! dummy loop indices
34      REAL(wp) ::   zr = 0.2     ! maximum of the reduction factor at the bottom ocean ( 0 < zr < 1 )
35      REAL(wp) ::   zh = 500.    ! depth of at which start the reduction ( > dept(1) )
36      REAL(wp) ::   zd_max       ! maximum grid spacing over the global domain
37      REAL(wp) ::   za00, zc, zd, zetmax, zefmax, zeumax, zevmax   ! local scalars
38      REAL(wp), POINTER, DIMENSION(:) :: zcoef   
39      !!----------------------------------------------------------------------
40      !
41      CALL wrk_alloc( jpk, zcoef )
42      !
43      IF(lwp) WRITE(numout,*)
44      IF(lwp) WRITE(numout,*) 'ldf_dyn_c3d : 3D lateral eddy viscosity coefficient'
45      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
46
47     
48      ! Set ahm1 and ahm2  ( T- and F- points) (used for laplacian operators
49      ! =================                       whatever its orientation is)
50      IF( ln_dynldf_lap ) THEN
51         ! define ahm1 and ahm2 at the right grid point position
52         ! (USER: modify ahm1 and ahm2 following your desiderata)
53
54         zd_max = MAX( MAXVAL( e1t(:,:) ), MAXVAL( e2t(:,:) ) )
55         IF( lk_mpp )   CALL mpp_max( zd_max )   ! max over the global domain
56
57         IF(lwp) WRITE(numout,*) '              laplacian operator: ahm proportional to e1'
58         IF(lwp) WRITE(numout,*) '              maximum grid-spacing = ', zd_max, ' maximum value for ahm = ', ahm0
59
60         za00 = ahm0 / zd_max
61
62         IF( ln_dynldf_iso ) THEN
63            IF(lwp) WRITE(numout,*) '              Caution, as implemented now, the isopycnal part of momentum'
64            IF(lwp) WRITE(numout,*) '                 mixing use aht0 as eddy viscosity coefficient. Thus, it is'
65            IF(lwp) WRITE(numout,*) '                 uniform and you must be sure that your ahm is greater than'
66            IF(lwp) WRITE(numout,*) '                 aht0 everywhere in the model domain.'
67         ENDIF
68
69         CALL ldf_zpf( .TRUE. , 1000., 500., 0.25, fsdept(:,:,:), ahm1 )   ! vertical profile
70         CALL ldf_zpf( .TRUE. , 1000., 500., 0.25, fsdept(:,:,:), ahm2 )   ! vertical profile
71         DO jk = 1,jpk
72            DO jj = 1, jpj
73               DO ji = 1, jpi
74                  zetmax = MAX( e1t(ji,jj), e2t(ji,jj) )
75                  zefmax = MAX( e1f(ji,jj), e2f(ji,jj) )
76                  ahm1(ji,jj,jk) = za00 * zetmax * ahm1(ji,jj,jk)
77                  ahm2(ji,jj,jk) = za00 * zefmax * ahm2(ji,jj,jk)
78               END DO
79            END DO
80         END DO
81
82
83         ! Special case for ORCA R1, R2 and R4 configurations (overwrite the value of ahm1 ahm2)
84         ! ==============================================
85         IF( cp_cfg == "orca" .AND. ( jp_cfg == 1 .OR. jp_cfg == 2 .OR. jp_cfg == 4 ) ) THEN
86            IF(lwp) WRITE(numout,*)
87            IF(lwp) WRITE(numout,*) '              ORCA R1, R2 or R4: overwrite the previous definition of ahm'
88            IF(lwp) WRITE(numout,*) '              ================='
89            CALL ldf_dyn_c3d_orca( ld_print )
90         ENDIF
91
92      ENDIF
93     
94      ! Control print
95      IF(lwp .AND. ld_print ) THEN
96         WRITE(numout,*)
97         WRITE(numout,*) '         3D ahm1 array (k=1)'
98         CALL prihre( ahm1(:,:,1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1.e-3, numout )
99         WRITE(numout,*)
100         WRITE(numout,*) '         3D ahm2 array (k=1)'
101         CALL prihre( ahm2(:,:,1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1.e-3, numout )
102      ENDIF
103
104
105      ! ahm3 and ahm4 at U- and V-points (used for bilaplacian operator
106      ! ================================  whatever its orientation is)
107      ! (USER: modify ahm3 and ahm4 following your desiderata)
108      ! Here: ahm is proportional to the cube of the maximum of the gridspacing
109      !       in the to horizontal direction
110
111      IF( ln_dynldf_bilap ) THEN
112
113         zd_max = MAX( MAXVAL( e1u(:,:) ), MAXVAL( e2u(:,:) ) )
114         IF( lk_mpp )   CALL mpp_max( zd_max )   ! max over the global domain
115
116         IF(lwp) WRITE(numout,*) '              bi-laplacian operator: ahm proportional to e1**3 '
117         IF(lwp) WRITE(numout,*) '              maximum grid-spacing = ', zd_max, ' maximum value for ahm = ', ahm0
118
119         za00 = ahm0_blp / ( zd_max * zd_max * zd_max )
120         DO jj = 1, jpj
121            DO ji = 1, jpi
122               zeumax = MAX( e1u(ji,jj), e2u(ji,jj) )
123               zevmax = MAX( e1v(ji,jj), e2v(ji,jj) )
124               ahm3(ji,jj,1) = za00 * zeumax * zeumax * zeumax
125               ahm4(ji,jj,1) = za00 * zevmax * zevmax * zevmax
126            END DO
127         END DO
128
129         zh = MAX( zh, fsdept(1,1,1) )   ! at least the first reach ahm0
130         IF( ln_zco ) THEN               ! z-coordinate, same profile everywhere
131            IF(lwp) WRITE(numout,'(36x," ahm ", 7x)')
132            DO jk = 1, jpk
133               IF( fsdept(1,1,jk) <= zh ) THEN
134                  zcoef(jk) = 1.e0
135               ELSE
136                  zcoef(jk) = 1.e0 + ( zr - 1.e0 )   &
137                     &               * (  1. - EXP( ( fsdept(1,1,jk   ) - zh ) / zh )  )   &
138                     &               / (  1. - EXP( ( fsdept(1,1,jpkm1) - zh ) / zh )  )
139               ENDIF
140               ahm3(:,:,jk) = ahm3(:,:,1) * zcoef(jk)
141               ahm4(:,:,jk) = ahm4(:,:,1) * zcoef(jk)
142               IF(lwp) WRITE(numout,'(34x,E7.2,8x,i3)') zcoef(jk) * ahm0, jk
143            END DO
144         ELSE                            ! partial steps or s-ccordinate
145            zc = MAXVAL( fsdept(:,:,jpkm1) )
146            IF( lk_mpp )   CALL mpp_max( zc )   ! max over the global domain
147
148            zc = 1. / (  1. - EXP( ( zc - zh ) / zh )  )
149            DO jk = 2, jpkm1
150               DO jj = 1, jpj
151                  DO ji = 1, jpi
152                     IF( fsdept(ji,jj,jk) <= zh ) THEN
153                        ahm3(ji,jj,jk) = ahm3(ji,jj,1)
154                        ahm4(ji,jj,jk) = ahm4(ji,jj,1)
155                     ELSE
156                        zd = 1.e0 + ( zr - 1.e0 ) * (  1. - EXP( ( fsdept(ji,jj,jk) - zh ) / zh )  ) * zc
157                        ahm3(ji,jj,jk) = ahm3(ji,jj,1) * zd
158                        ahm4(ji,jj,jk) = ahm4(ji,jj,1) * zd
159                     ENDIF
160                  END DO
161               END DO
162            END DO
163            ahm3(:,:,jpk) = ahm3(:,:,jpkm1)
164            ahm4(:,:,jpk) = ahm4(:,:,jpkm1)
165            IF(lwp) WRITE(numout,'(36x," ahm ", 7x)')
166            DO jk = 1, jpk
167               IF(lwp) WRITE(numout,'(30x,E10.2,8x,i3)') ahm3(1,1,jk), jk
168            END DO
169         ENDIF
170
171         ! Control print
172         IF( lwp .AND. ld_print ) THEN
173            WRITE(numout,*)
174            WRITE(numout,*) 'inildf: ahm3 array at level 1'
175            CALL prihre(ahm3(:,:,1  ),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)
176            WRITE(numout,*)
177            WRITE(numout,*) 'inildf: ahm4 array at level 1'
178            CALL prihre(ahm4(:,:,jpk),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)
179         ENDIF
180      ENDIF
181      !
182      CALL wrk_dealloc( jpk, zcoef )
183      !
184   END SUBROUTINE ldf_dyn_c3d
185
186
187   SUBROUTINE ldf_dyn_c3d_orca( ld_print )
188      !!----------------------------------------------------------------------
189      !!                  ***  ROUTINE ldf_dyn_c3d  ***
190      !!                   
191      !! ** Purpose :   ORCA R1, R2 and R4 only
192      !!
193      !! ** Method  :   blah blah blah ....
194      !!----------------------------------------------------------------------
195      USE ldftra_oce, ONLY:   aht0
196      USE iom
197      !!
198      LOGICAL, INTENT(in) ::   ld_print   ! If true, output arrays on numout
199      !!
200      INTEGER ::   ji, jj, jk, jn      ! dummy loop indices
201      INTEGER ::   ii0, ii1, ij0, ij1  ! local integers
202      INTEGER ::   inum, iim, ijm      !
203      INTEGER ::   ifreq, il1, il2, ij, ii
204      REAL(wp) ::   zahmeq, zcoff, zcoft, zmsk   ! local scalars
205      REAL(wp) ::   zemax , zemin, zetmax, zefmax, zeref, zahmm, zemax75
206      CHARACTER (len=15) ::   clexp
207      INTEGER , POINTER, DIMENSION(:,:)  :: icof
208      INTEGER , POINTER, DIMENSION(:,:)  :: imsk
209      REAL(wp), POINTER, DIMENSION(:  )  :: zcoef   
210      REAL(wp), POINTER, DIMENSION(:,:)  :: zahm0
211      !
212      REAL(wp), DIMENSION(:,:), ALLOCATABLE ::   ztemp2d  ! temporary array to read ahmcoef file
213      !!----------------------------------------------------------------------
214      !
215      CALL wrk_alloc( jpi   , jpj   , icof  )
216      CALL wrk_alloc( jpi   , jpj   , imsk  )
217      CALL wrk_alloc( jpk   ,         zcoef )
218      CALL wrk_alloc( jpi   , jpj   , zahm0 )
219      !
220      IF(lwp) WRITE(numout,*)
221      IF(lwp) WRITE(numout,*) 'ldfdyn_c3d_orca : 3D eddy viscosity coefficient'
222      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~'
223      IF(lwp) WRITE(numout,*) '        orca R1, R2 or R4 configuration: reduced in the surface Eq. strip '
224
225      ! Read 2d integer array to specify western boundary increase in the
226      ! ===================== equatorial strip (20N-20S) defined at t-points
227      ALLOCATE( ztemp2d(jpi,jpj) )
228      ztemp2d(:,:) = 0.
229      CALL iom_open ( 'ahmcoef.nc', inum )
230      CALL iom_get  ( inum, jpdom_data, 'icof', ztemp2d)
231      icof(:,:)  = NINT(ztemp2d(:,:))
232      CALL iom_close( inum )
233      DEALLOCATE(ztemp2d)
234
235      ! Set ahm1 and ahm2
236      ! =================
237     
238      ! define ahm1 and ahm2 at the right grid point position
239      ! (USER: modify ahm1 and ahm2 following your desiderata)
240      ! biharmonic : ahm1 (ahm2) defined at u- (v-) point
241      ! harmonic   : ahm1 (ahm2) defined at t- (f-) point
242     
243      ! first level : as for 2D coefficients
244     
245      ! Decrease ahm to zahmeq m2/s in the tropics
246      ! (from 90 to 20 degre: ahm = constant
247      ! from 20 to  2.5 degre: ahm = decrease in (1-cos)/2
248      ! from  2.5 to  0 degre: ahm = constant
249      ! symmetric in the south hemisphere)
250     
251      IF( jp_cfg == 4 )   THEN
252         zahmeq = 5.0 * aht0
253         zahmm  = min( 160000.0, ahm0)
254         zemax = MAXVAL ( e1t(:,:) * e2t(:,:), ssmask(:,:) .GE. 0.5 )
255         zemin = MINVAL ( e1t(:,:) * e2t(:,:), ssmask(:,:) .GE. 0.5 )
256         zeref = MAXVAL ( e1t(:,:) * e2t(:,:),   &
257             &   ssmask(:,:) .GE. 0.5 .AND. ABS(gphit(:,:)) .GT. 50. )
258 
259         DO jj = 1, jpj
260           DO ji = 1, jpi
261              zmsk = e1t(ji,jj) * e2t(ji,jj)
262              IF( abs(gphit(ji,jj)) .LE. 15 ) THEN
263                 zahm0(ji,jj) = ahm0
264              ELSE
265                 IF( zmsk .GE. zeref ) THEN
266                    zahm0(ji,jj) = ahm0
267                 ELSE
268                    zahm0(ji,jj) = zahmm + (ahm0-zahmm)*(1.0 -   &
269                        &          cos((rpi*0.5*(zmsk-zemin)/(zeref-zemin))))
270                 ENDIF
271              ENDIF
272           END DO
273         END DO
274      ENDIF
275
276      IF( jp_cfg == 2 )   THEN
277         zahmeq     = aht0
278         zahmm      = ahm0
279         zahm0(:,:) = ahm0
280      ENDIF
281
282      IF( jp_cfg == 1 )   THEN
283         zahmeq     = aht0  ! reduced to aht0 on equator; set to ahm0 if no tropical reduction is required
284         zahmm      = ahm0
285         zahm0(:,:) = ahm0
286      ENDIF
287
288      DO jj = 1, jpj
289         DO ji = 1, jpi
290            IF( ABS(gphif(ji,jj)) >= 20.) THEN
291               ahm2(ji,jj,1) =  zahm0(ji,jj)
292            ELSEIF( ABS(gphif(ji,jj)) <= 2.5) THEN
293               ahm2(ji,jj,1) =  zahmeq
294            ELSE
295               ahm2(ji,jj,1) = zahmeq + (zahm0(ji,jj)-zahmeq)/2.   &
296                  &            *(1.-COS( rad*(ABS(gphif(ji,jj))-2.5)*180./17.5 ) )
297            ENDIF
298            IF( ABS(gphit(ji,jj)) >= 20.) THEN
299               ahm1(ji,jj,1) =  zahm0(ji,jj)
300            ELSEIF( ABS(gphit(ji,jj)) <= 2.5) THEN
301               ahm1(ji,jj,1) =  zahmeq
302            ELSE
303               ahm1(ji,jj,1) = zahmeq + (zahm0(ji,jj)-zahmeq)/2.   &
304                  &            *(1.-COS( rad*(ABS(gphit(ji,jj))-2.5)*180./17.5 ) )
305            ENDIF
306         END DO
307      END DO
308     
309      ! increase along western boundaries of equatorial strip
310      ! t-point
311      DO jj = 1, jpjm1
312         DO ji = 1, jpim1
313            zcoft = float( icof(ji,jj) ) / 100.
314            ahm1(ji,jj,1) = zcoft * zahm0(ji,jj) + (1.-zcoft) * ahm1(ji,jj,1)
315         END DO
316      END DO
317      ! f-point
318      icof(:,:) = icof(:,:) * ssmask(:,:)
319      DO jj = 1, jpjm1
320         DO ji = 1, jpim1   ! NO vector opt.
321            zmsk = ssmask(ji,jj+1) + ssmask(ji+1,jj+1) + ssmask(ji,jj) + ssmask(ji,jj+1)
322            IF( zmsk == 0. ) THEN
323               zcoff = 1.
324            ELSE
325               zcoff = FLOAT( icof(ji,jj+1) + icof(ji+1,jj+1) + icof(ji,jj) + icof(ji,jj+1) )   &
326                     / (zmsk * 100.)
327            ENDIF
328            ahm2(ji,jj,1) = zcoff * zahm0(ji,jj) + (1.-zcoff) * ahm2(ji,jj,1)
329         END DO
330      END DO
331
332      ! other level: re-increase the coef in the deep ocean
333      !==================================================================
334      ! Prior to v3.3, zcoeff was hardwired according to k-index jk.
335      !
336      ! From v3.3 onwards this has been generalised to a function of
337      ! depth so that it can be used with any number of levels.
338      !
339      ! The function has been chosen to match the original values (shown
340      ! in the following comments) when using the standard 31 ORCA levels. 
341      !     DO jk = 1, 21
342      !        zcoef(jk) = 1._wp
343      !     END DO
344      !     zcoef(22) = 2._wp
345      !     zcoef(23) = 3._wp
346      !     zcoef(24) = 5._wp
347      !     zcoef(25) = 7._wp
348      !     zcoef(26) = 9._wp
349      !     DO jk = 27, jpk
350      !        zcoef(jk) = 10._wp
351      !     END DO
352      !==================================================================
353
354       IF(lwp) THEN
355          WRITE(numout,*)
356          WRITE(numout,*) '         1D zcoef array '
357          WRITE(numout,*) '         ~~~~~~~~~~~~~~ '
358          WRITE(numout,*)
359          WRITE(numout,*) '    jk        zcoef '
360       ENDIF
361
362      DO jk=1, jpk
363         zcoef(jk) = 1.0_wp + NINT(9.0_wp*(gdept_1d(jk)-800.0_wp)/(3000.0_wp-800.0_wp))
364         zcoef(jk) = MIN(10.0_wp, MAX(1.0_wp, zcoef(jk)))
365         IF(lwp) WRITE(numout,'(4x,i3,6x,f7.3)') jk,zcoef(jk)
366      END DO
367
368      DO jk = 2, jpk
369         ahm1(:,:,jk) = MIN( zahm0(:,:), zcoef(jk) * ahm1(:,:,1) )
370         ahm2(:,:,jk) = MIN( zahm0(:,:), zcoef(jk) * ahm2(:,:,1) )
371      END DO
372
373      IF( jp_cfg == 4 )   THEN         ! Limit AHM in Gibraltar strait
374         ij0 = 50   ;   ij1 = 53
375         ii0 = 69   ;   ii1 = 71
376         DO jk = 1, jpk
377            ahm1(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),jk) = MIN( zahmm, ahm1(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),jk) )
378            ahm2(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),jk) = MIN( zahmm, ahm2(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),jk) )
379         END DO
380      ENDIF
381      CALL lbc_lnk( ahm1, 'T', 1. )   ! Lateral boundary conditions (unchanged sign)
382      CALL lbc_lnk( ahm2, 'F', 1. )
383
384      IF( jp_cfg == 1 )   THEN         ! Limit AHM south of -75 (critical for Giant ice shelves with small e1/e2)
385         ! special orca1 treatment remove the grid size scaling. Need to restore it south of -75N
386         ! define max grid size south of -75
387         imsk(:,:) = 0
388         WHERE( gphit(:,:) < -75.0 ) imsk(:,:) = 1 
389         zemax75 = MAX( MAXVAL( e1t(:,:)*imsk(:,:) ), MAXVAL( e2t(:,:)*imsk(:,:) ) )
390         IF( lk_mpp )   CALL mpp_max(zemax75)
391         !
392         ! apply grid size scaling of ahm south of -75 (no change north of it)
393         DO jj = 1, jpj
394            DO ji = 1, jpi
395               IF ( gphit(ji,jj) < -75.0 ) THEN
396                  zetmax = MAX( e1t(ji,jj), e2t(ji,jj) )
397                  zefmax = MAX( e1f(ji,jj), e2f(ji,jj) )
398                  ahm1(ji,jj,:) = zetmax / zemax75 * ahm1(ji,jj,:)
399                  ahm2(ji,jj,:) = zefmax / zemax75 * ahm2(ji,jj,:)
400               END IF
401            END DO
402         END DO
403      END IF
404
405
406      IF(lwp) THEN                    ! Control print
407         WRITE(numout,*)
408         WRITE(numout,*) '         3D ahm1 array (k=1)'
409         CALL prihre( ahm1(:,:,1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1.e-3, numout )
410         WRITE(numout,*)
411         WRITE(numout,*) '         3D ahm2 array (k=1)'
412         CALL prihre( ahm2(:,:,1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1.e-3, numout )
413         WRITE(numout,*)
414         WRITE(numout,*) '         3D ahm2 array (k=jpk)'
415         CALL prihre( ahm2(:,:,jpk), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1.e-3, numout )
416      ENDIF
417
418
419      ! Set ahm3 and ahm4
420      ! =================
421
422      ! define ahm3 and ahm4 at the right grid point position
423      ! initialization to a constant value
424      !     (USER: modify ahm3 and ahm4 following your desiderata)
425      !     harmonic isopycnal or geopotential:
426      !                          ahm3 (ahm4) defined at u- (v-) point
427      DO jk = 1, jpk
428         DO jj = 2, jpj
429            DO ji = 2, jpi
430               ahm3(ji,jj,jk) = 0.5 * ( ahm2(ji,jj,jk) + ahm2(ji  ,jj-1,jk) )
431               ahm4(ji,jj,jk) = 0.5 * ( ahm2(ji,jj,jk) + ahm2(ji-1,jj  ,jk) )
432            END DO
433         END DO
434      END DO
435      ahm3 ( :, 1, :) = ahm3 ( :, 2, :)
436      ahm4 ( :, 1, :) = ahm4 ( :, 2, :)
437     
438      CALL lbc_lnk( ahm3, 'U', 1. )    ! Lateral boundary conditions (unchanged sign)
439      CALL lbc_lnk( ahm4, 'V', 1. )
440
441      ! Control print
442
443      IF( lwp .AND. ld_print ) THEN
444         WRITE(numout,*)
445         WRITE(numout,*) '         ahm3 array level 1'
446         CALL prihre(ahm3(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)
447         WRITE(numout,*)
448         WRITE(numout,*) '         ahm4 array level 1'
449         CALL prihre(ahm4(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)
450      ENDIF
451      !
452      CALL wrk_dealloc( jpi   , jpj   , icof  )
453      CALL wrk_dealloc( jpi   , jpj   , imsk  )
454      CALL wrk_dealloc( jpk   ,         zcoef )
455      CALL wrk_dealloc( jpi   , jpj   , zahm0 )
456      !
457   END SUBROUTINE ldf_dyn_c3d_orca
Note: See TracBrowser for help on using the repository browser.