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.
ldfdyn_c3d.h90 in branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/LDF – NEMO

source: branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_c3d.h90 @ 4292

Last change on this file since 4292 was 4292, checked in by cetlod, 11 years ago

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