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.
p4zsms.F90 in branches/2013/dev_r3853_CNRS9_ConfSetting/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z – NEMO

source: branches/2013/dev_r3853_CNRS9_ConfSetting/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90 @ 3875

Last change on this file since 3875 was 3875, checked in by clevy, 11 years ago

Configuration Setting/Step? 1, see ticket:#1074

File size: 25.7 KB
Line 
1MODULE p4zsms
2   !!======================================================================
3   !!                         ***  MODULE p4zsms  ***
4   !! TOP :   PISCES Source Minus Sink manager
5   !!======================================================================
6   !! History :   1.0  !  2004-03 (O. Aumont) Original code
7   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
8   !!----------------------------------------------------------------------
9#if defined key_pisces
10   !!----------------------------------------------------------------------
11   !!   'key_pisces'                                       PISCES bio-model
12   !!----------------------------------------------------------------------
13   !!   p4zsms        :  Time loop of passive tracers sms
14   !!----------------------------------------------------------------------
15   USE oce_trc         !  shared variables between ocean and passive tracers
16   USE trc             !  passive tracers common variables
17   USE trcdta
18   USE sms_pisces      !  PISCES Source Minus Sink variables
19   USE p4zbio          !  Biological model
20   USE p4zche          !  Chemical model
21   USE p4zlys          !  Calcite saturation
22   USE p4zflx          !  Gas exchange
23   USE p4zsbc          !  External source of nutrients
24   USE p4zsed          !  Sedimentation
25   USE p4zint          !  time interpolation
26   USE iom             !  I/O manager
27   USE trdmod_oce      !  Ocean trends variables
28   USE trdmod_trc      !  TOP trends variables
29   USE sedmodel        !  Sediment model
30   USE prtctl_trc      !  print control for debugging
31
32   IMPLICIT NONE
33   PRIVATE
34
35   PUBLIC   p4z_sms_init    ! called in p4zsms.F90
36   PUBLIC   p4z_sms    ! called in p4zsms.F90
37
38   REAL(wp) :: alkbudget, no3budget, silbudget, ferbudget
39   INTEGER ::  numco2, numnut  !: logical unit for co2 budget
40
41   !!----------------------------------------------------------------------
42   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
43   !! $Id: p4zsms.F90 3320 2012-03-05 16:37:52Z cetlod $
44   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
45   !!----------------------------------------------------------------------
46
47CONTAINS
48
49   SUBROUTINE p4z_sms( kt )
50      !!---------------------------------------------------------------------
51      !!                     ***  ROUTINE p4z_sms  ***
52      !!
53      !! ** Purpose :   Managment of the call to Biological sources and sinks
54      !!              routines of PISCES bio-model
55      !!
56      !! ** Method  : - at each new day ...
57      !!              - several calls of bio and sed ???
58      !!              - ...
59      !!---------------------------------------------------------------------
60      !
61      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index     
62      !!
63      INTEGER ::   jnt, jn, jl
64      CHARACTER (len=25) :: charout
65      REAL(wp), POINTER, DIMENSION(:,:,:,:)  :: ztrdpis
66      !!---------------------------------------------------------------------
67      !
68      IF( nn_timing == 1 )  CALL timing_start('p4z_sms')
69      !
70      IF( l_trdtrc )  THEN
71         CALL wrk_alloc( jpi, jpj, jpk, jp_pisces, ztrdpis ) 
72         DO jn = 1, jp_pisces
73            jl = jn + jp_pcs0 - 1
74            ztrdpis(:,:,:,jn) = trn(:,:,:,jl)
75         ENDDO
76      ENDIF
77      !
78      IF( ln_rsttr .AND. kt == nittrc000 )                         CALL p4z_rst( nittrc000, 'READ' )  !* read or initialize all required fields
79      IF( ln_rsttr  .AND. ln_pisclo )                              CALL p4z_clo            ! damping on closed seas
80      IF( ln_pisdmp .AND. MOD( kt - nn_dttrc, nn_pisdmp ) == 0 )   CALL p4z_dmp( kt )      ! Relaxation of some tracers
81      !
82      IF( ndayflxtr /= nday_year ) THEN      ! New days
83         !
84         ndayflxtr = nday_year
85
86         IF(lwp) write(numout,*)
87         IF(lwp) write(numout,*) ' New chemical constants and various rates for biogeochemistry at new day : ', nday_year
88         IF(lwp) write(numout,*) '~~~~~~'
89
90         CALL p4z_che              ! computation of chemical constants
91         CALL p4z_int( kt )        ! computation of various rates for biogeochemistry
92         !
93      ENDIF
94
95      IF( ll_sbc ) CALL p4z_sbc( kt )   ! external sources of nutrients
96
97      DO jnt = 1, nrdttrc          ! Potential time splitting if requested
98         !
99         CALL p4z_bio (kt, jnt)    ! Biology
100         CALL p4z_sed (kt, jnt)    ! Sedimentation
101         !
102         DO jn = jp_pcs0, jp_pcs1
103            trb(:,:,:,jn) = trn(:,:,:,jn)
104         ENDDO
105         !
106      END DO
107
108      IF( l_trdtrc )  THEN
109         DO jn = 1, jp_pisces
110            jl = jn + jp_pcs0 - 1
111            ztrdpis(:,:,:,jn) = ( ztrdpis(:,:,:,jn) - trn(:,:,:,jl) ) * rfact2r
112         ENDDO
113      ENDIF
114
115      CALL p4z_lys( kt )             ! Compute CaCO3 saturation
116      CALL p4z_flx( kt )             ! Compute surface fluxes
117
118      DO jn = jp_pcs0, jp_pcs1
119        CALL lbc_lnk( trn(:,:,:,jn), 'T', 1. )
120        CALL lbc_lnk( trb(:,:,:,jn), 'T', 1. )
121        CALL lbc_lnk( tra(:,:,:,jn), 'T', 1. )
122      END DO
123      !
124      IF( lk_sed ) THEN 
125         !
126         CALL sed_model( kt )     !  Main program of Sediment model
127         !
128         DO jn = jp_pcs0, jp_pcs1
129           CALL lbc_lnk( trn(:,:,:,jn), 'T', 1. )
130         END DO
131         !
132      ENDIF
133      !
134      IF( lrst_trc )  CALL p4z_rst( kt, 'WRITE' )  !* Write PISCES informations in restart file
135      !
136      IF( l_trdtrc ) THEN
137         DO jn = 1, jp_pisces
138            jl = jn + jp_pcs0 - 1
139             ztrdpis(:,:,:,jn) = ztrdpis(:,:,:,jn) + tra(:,:,:,jl)
140             CALL trd_mod_trc( ztrdpis(:,:,:,jn), jn, jptra_trd_sms, kt )   ! save trends
141          END DO
142          CALL wrk_dealloc( jpi, jpj, jpk, jp_pisces, ztrdpis ) 
143      END IF
144      !
145      CALL p4z_chk_mass( kt ) ! Mass conservation checking
146
147      IF( nn_timing == 1 )  CALL timing_stop('p4z_sms')
148      !
149      !
150   END SUBROUTINE p4z_sms
151
152   SUBROUTINE p4z_sms_init
153      !!----------------------------------------------------------------------
154      !!                     ***  p4z_sms_init  *** 
155      !!
156      !! ** Purpose :   read PISCES namelist
157      !!
158      !! ** input   :   file 'namelist.trc.s' containing the following
159      !!             namelist: natext, natbio, natsms
160      !!                       natkriest ("key_kriest")
161      !!----------------------------------------------------------------------
162      NAMELIST/nampisbio/ nrdttrc, wsbio, xkmort, ferat3, wsbio2, niter1max, niter2max
163#if defined key_kriest
164      NAMELIST/nampiskrp/ xkr_eta, xkr_zeta, xkr_ncontent, xkr_mass_min, xkr_mass_max
165#endif
166      NAMELIST/nampisdmp/ ln_pisdmp, nn_pisdmp, ln_pisclo
167      NAMELIST/nampismass/ ln_check_mass
168      INTEGER :: ios                 ! Local integer output status for namelist read
169      !!----------------------------------------------------------------------
170
171      REWIND( numnatp_ref )              ! Namelist nampisbio in reference namelist : Pisces variables
172      READ  ( numnatp_ref, nampisbio, IOSTAT = ios, ERR = 901)
173901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisbio in reference namelist', lwp )
174
175      REWIND( numnatp_cfg )              ! Namelist nampisbio in configuration namelist : Pisces variables
176      READ  ( numnatp_cfg, nampisbio, IOSTAT = ios, ERR = 902 )
177902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisbio in configuration namelist', lwp )
178      WRITE ( numonp, nampisbio )
179
180      IF(lwp) THEN                         ! control print
181         WRITE(numout,*) ' Namelist : nampisbio'
182         WRITE(numout,*) '    frequence pour la biologie                nrdttrc   =', nrdttrc
183         WRITE(numout,*) '    POC sinking speed                         wsbio     =', wsbio
184         WRITE(numout,*) '    half saturation constant for mortality    xkmort    =', xkmort
185         WRITE(numout,*) '    Fe/C in zooplankton                       ferat3    =', ferat3
186         WRITE(numout,*) '    Big particles sinking speed               wsbio2    =', wsbio2
187         WRITE(numout,*) '    Maximum number of iterations for POC      niter1max =', niter1max
188         WRITE(numout,*) '    Maximum number of iterations for GOC      niter2max =', niter2max
189      ENDIF
190
191#if defined key_kriest
192
193      !                               ! nampiskrp : kriest parameters
194      !                               ! -----------------------------
195      xkr_eta      = 0.62
196      xkr_zeta     = 1.62
197      xkr_ncontent = 5.7E-6
198      xkr_mass_min = 0.0002
199      xkr_mass_max = 1.
200
201      REWIND( numnatp_ref )              ! Namelist nampiskrp in reference namelist : Pisces Kriest
202      READ  ( numnatp_ref, nampiskrp, IOSTAT = ios, ERR = 903)
203903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampiskrp in reference namelist', lwp )
204
205      REWIND( numnatp_cfg )              ! Namelist nampiskrp in configuration namelist : Pisces Kriest
206      READ  ( numnatp_cfg, nampiskrp, IOSTAT = ios, ERR = 904 )
207904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampiskrp in configuration namelist', lwp )
208      WRITE ( numonp, nampiskrp )
209
210      IF(lwp) THEN
211         WRITE(numout,*)
212         WRITE(numout,*) ' Namelist : nampiskrp'
213         WRITE(numout,*) '    Sinking  exponent                        xkr_eta      = ', xkr_eta
214         WRITE(numout,*) '    N content exponent                       xkr_zeta     = ', xkr_zeta
215         WRITE(numout,*) '    N content factor                         xkr_ncontent = ', xkr_ncontent
216         WRITE(numout,*) '    Minimum mass for Aggregates              xkr_mass_min = ', xkr_mass_min
217         WRITE(numout,*) '    Maximum mass for Aggregates              xkr_mass_max = ', xkr_mass_max
218         WRITE(numout,*)
219     ENDIF
220
221
222     ! Computation of some variables
223     xkr_massp = xkr_ncontent * 7.625 * xkr_mass_min**xkr_zeta
224
225#endif
226
227      ln_pisdmp = .true.
228      nn_pisdmp = 1
229      ln_pisclo = .false.
230
231      REWIND( numnatp_ref )              ! Namelist nampisdmp in reference namelist : Pisces damping
232      READ  ( numnatp_ref, nampisdmp, IOSTAT = ios, ERR = 905)
233905   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisdmp in reference namelist', lwp )
234
235      REWIND( numnatp_cfg )              ! Namelist nampisdmp in configuration namelist : Pisces damping
236      READ  ( numnatp_cfg, nampisdmp, IOSTAT = ios, ERR = 906 )
237906   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisdmp in configuration namelist', lwp )
238      WRITE ( numonp, nampisdmp )
239
240      IF(lwp) THEN                         ! control print
241         WRITE(numout,*)
242         WRITE(numout,*) ' Namelist : nampisdmp'
243         WRITE(numout,*) '    Relaxation of tracer to glodap mean value             ln_pisdmp      =', ln_pisdmp
244         WRITE(numout,*) '    Frequency of Relaxation                               nn_pisdmp      =', nn_pisdmp
245         WRITE(numout,*) '    Restoring of tracer to initial value  on closed seas  ln_pisclo      =', ln_pisclo
246         WRITE(numout,*) ' '
247      ENDIF
248
249      ln_check_mass = .false.
250      REWIND( numnatp_ref )              ! Namelist nampismass in reference namelist : Pisces mass conservation check
251      READ  ( numnatp_ref, nampismass, IOSTAT = ios, ERR = 907)
252907   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampismass in reference namelist', lwp )
253
254      REWIND( numnatp_cfg )              ! Namelist nampismass in configuration namelist : Pisces mass conservation check
255      READ  ( numnatp_cfg, nampismass, IOSTAT = ios, ERR = 908 )
256908   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampismass in configuration namelist', lwp )
257      WRITE ( numonp, nampismass )
258
259      IF(lwp) THEN                         ! control print
260         WRITE(numout,*) ' '
261         WRITE(numout,*) ' Namelist parameter for mass conservation checking'
262         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
263         WRITE(numout,*) '    Flag to check mass conservation of NO3/Si/TALK ln_check_mass = ', ln_check_mass
264      ENDIF
265
266   END SUBROUTINE p4z_sms_init
267
268   SUBROUTINE p4z_rst( kt, cdrw )
269      !!---------------------------------------------------------------------
270      !!                   ***  ROUTINE p4z_rst  ***
271      !!
272      !!  ** Purpose : Read or write variables in restart file:
273      !!
274      !!  WRITE(READ) mode:
275      !!       kt        : number of time step since the begining of the experiment at the
276      !!                   end of the current(previous) run
277      !!---------------------------------------------------------------------
278      INTEGER         , INTENT(in) ::   kt         ! ocean time-step
279      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag
280      !
281      INTEGER  ::  ji, jj, jk
282      REAL(wp) ::  zcaralk, zbicarb, zco3
283      REAL(wp) ::  ztmas, ztmas1
284      !!---------------------------------------------------------------------
285
286      IF( TRIM(cdrw) == 'READ' ) THEN
287         !
288         IF(lwp) WRITE(numout,*)
289         IF(lwp) WRITE(numout,*) ' p4z_rst : Read specific variables from pisces model '
290         IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~'
291         !
292         IF( iom_varid( numrtr, 'PH', ldstop = .FALSE. ) > 0 ) THEN
293            CALL iom_get( numrtr, jpdom_autoglo, 'PH' , hi(:,:,:)  )
294         ELSE
295!            hi(:,:,:) = 1.e-9
296            ! Set PH from  total alkalinity, borat (???), akb3 (???) and ak23 (???)
297            ! --------------------------------------------------------
298            DO jk = 1, jpk
299               DO jj = 1, jpj
300                  DO ji = 1, jpi
301                     ztmas   = tmask(ji,jj,jk)
302                     ztmas1  = 1. - tmask(ji,jj,jk)
303                     zcaralk = trn(ji,jj,jk,jptal) - borat(ji,jj,jk) / (  1. + 1.E-8 / ( rtrn + akb3(ji,jj,jk) )  )
304                     zco3    = ( zcaralk - trn(ji,jj,jk,jpdic) ) * ztmas + 0.5e-3 * ztmas1
305                     zbicarb = ( 2. * trn(ji,jj,jk,jpdic) - zcaralk )
306                     hi(ji,jj,jk) = ( ak23(ji,jj,jk) * zbicarb / zco3 ) * ztmas + 1.e-9 * ztmas1
307                  END DO
308               END DO
309            END DO
310         ENDIF
311         CALL iom_get( numrtr, jpdom_autoglo, 'Silicalim', xksi(:,:) )
312         IF( iom_varid( numrtr, 'Silicamax', ldstop = .FALSE. ) > 0 ) THEN
313            CALL iom_get( numrtr, jpdom_autoglo, 'Silicamax' , xksimax(:,:)  )
314         ELSE
315            xksimax(:,:) = xksi(:,:)
316         ENDIF
317         !
318      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN
319         IF( kt == nitrst ) THEN
320            IF(lwp) WRITE(numout,*)
321            IF(lwp) WRITE(numout,*) 'p4z_rst : write pisces restart file  kt =', kt
322            IF(lwp) WRITE(numout,*) '~~~~~~~'
323         ENDIF
324         CALL iom_rstput( kt, nitrst, numrtw, 'PH', hi(:,:,:) )
325         CALL iom_rstput( kt, nitrst, numrtw, 'Silicalim', xksi(:,:) )
326         CALL iom_rstput( kt, nitrst, numrtw, 'Silicamax', xksimax(:,:) )
327      ENDIF
328      !
329   END SUBROUTINE p4z_rst
330
331   SUBROUTINE p4z_dmp( kt )
332      !!----------------------------------------------------------------------
333      !!                    ***  p4z_dmp  ***
334      !!
335      !! ** purpose  : Relaxation of some tracers
336      !!----------------------------------------------------------------------
337      !
338      INTEGER, INTENT( in )  ::     kt ! time step
339      !
340      REAL(wp) ::  alkmean = 2426.     ! mean value of alkalinity ( Glodap ; for Goyet 2391. )
341      REAL(wp) ::  po4mean = 2.165     ! mean value of phosphates
342      REAL(wp) ::  no3mean = 30.90     ! mean value of nitrate
343      REAL(wp) ::  silmean = 91.51     ! mean value of silicate
344      !
345      REAL(wp) :: zarea, zalksum, zpo4sum, zno3sum, zsilsum
346      !!---------------------------------------------------------------------
347
348
349      IF(lwp)  WRITE(numout,*)
350      IF(lwp)  WRITE(numout,*) ' p4z_dmp : Restoring of nutrients at time-step kt = ', kt
351      IF(lwp)  WRITE(numout,*)
352
353      IF( cp_cfg == "orca" .AND. .NOT. lk_c1d ) THEN      ! ORCA configuration (not 1D) !
354         !                                                    ! --------------------------- !
355         ! set total alkalinity, phosphate, nitrate & silicate
356         zarea          = 1._wp / glob_sum( cvol(:,:,:) ) * 1e6             
357
358         zalksum = glob_sum( trn(:,:,:,jptal) * cvol(:,:,:)  ) * zarea
359         zpo4sum = glob_sum( trn(:,:,:,jppo4) * cvol(:,:,:)  ) * zarea * po4r
360         zno3sum = glob_sum( trn(:,:,:,jpno3) * cvol(:,:,:)  ) * zarea * rno3
361         zsilsum = glob_sum( trn(:,:,:,jpsil) * cvol(:,:,:)  ) * zarea
362 
363         IF(lwp) WRITE(numout,*) '       TALK mean : ', zalksum
364         trn(:,:,:,jptal) = trn(:,:,:,jptal) * alkmean / zalksum
365
366         IF(lwp) WRITE(numout,*) '       PO4  mean : ', zpo4sum
367         trn(:,:,:,jppo4) = trn(:,:,:,jppo4) * po4mean / zpo4sum
368
369         IF(lwp) WRITE(numout,*) '       NO3  mean : ', zno3sum
370         trn(:,:,:,jpno3) = trn(:,:,:,jpno3) * no3mean / zno3sum
371
372         IF(lwp) WRITE(numout,*) '       SiO3 mean : ', zsilsum
373         trn(:,:,:,jpsil) = MIN( 400.e-6,trn(:,:,:,jpsil) * silmean / zsilsum )
374         !
375      ENDIF
376
377   END SUBROUTINE p4z_dmp
378
379
380   SUBROUTINE p4z_chk_mass( kt )
381      !!----------------------------------------------------------------------
382      !!                  ***  ROUTINE p4z_chk_mass  ***
383      !!
384      !! ** Purpose :  Mass conservation check
385      !!
386      !!---------------------------------------------------------------------
387      !
388      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index     
389      !!
390      !!---------------------------------------------------------------------
391
392      IF( kt == nittrc000 ) THEN
393         IF( ln_check_mass .AND. lwp) THEN      !   Open budget file of NO3, ALK, Si, Fer
394            CALL ctl_opn( numco2, 'carbon.budget'  , 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
395            CALL ctl_opn( numnut, 'nutrient.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
396         ENDIF
397      ENDIF
398
399      IF( ln_check_mass .AND. kt == nitend ) THEN      !   Compute the budget of NO3, ALK, Si, Fer
400         no3budget = glob_sum( (   trn(:,:,:,jpno3) + trn(:,:,:,jpnh4)  &
401            &                    + trn(:,:,:,jpphy) + trn(:,:,:,jpdia)  &
402            &                    + trn(:,:,:,jpzoo) + trn(:,:,:,jpmes)  &
403            &                    + trn(:,:,:,jppoc)                     &
404#if ! defined key_kriest
405            &                    + trn(:,:,:,jpgoc)                     &
406#endif
407            &                    + trn(:,:,:,jpdoc)                     ) * cvol(:,:,:)  ) 
408         !
409         silbudget = glob_sum( (   trn(:,:,:,jpsil) + trn(:,:,:,jpgsi)  &
410            &                    + trn(:,:,:,jpdsi)                     ) * cvol(:,:,:)  )
411         !
412         alkbudget = glob_sum( (   trn(:,:,:,jpno3) * rno3              &
413            &                    + trn(:,:,:,jptal)                     &
414            &                    + trn(:,:,:,jpcal) * 2.                ) * cvol(:,:,:)  )
415         !
416         ferbudget = glob_sum( (   trn(:,:,:,jpfer) + trn(:,:,:,jpnfe)  &
417            &                    + trn(:,:,:,jpdfe)                     &
418#if ! defined key_kriest
419            &                    + trn(:,:,:,jpbfe)                     &
420#endif
421            &                    + trn(:,:,:,jpsfe)                     &
422            &                    + trn(:,:,:,jpzoo)                     &
423            &                    + trn(:,:,:,jpmes) * ferat3            ) * cvol(:,:,:)  )
424
425         !
426         t_atm_co2_flx  = t_atm_co2_flx / glob_sum( e1e2t(:,:) )
427         t_oce_co2_flx  = t_oce_co2_flx         * 12. / 1.e15 * (-1 )
428         tpp            = tpp           * 1000. * 12. / 1.E15
429         t_oce_co2_exp  = t_oce_co2_exp * 1000. * 12. / 1.E15
430         !
431         no3budget = no3budget / areatot
432         silbudget = silbudget / areatot
433         alkbudget = alkbudget / areatot
434         ferbudget = ferbudget / areatot
435         !
436         IF(lwp) THEN
437            WRITE(numco2,9000) ndastp, t_atm_co2_flx, t_oce_co2_flx, tpp, t_oce_co2_exp
438            WRITE(numnut,9500) ndastp, alkbudget, no3budget, silbudget, ferbudget
439         ENDIF
440         !
441      ENDIF
442       !
443 9000  FORMAT(i8,f10.5,e18.10,f10.5,f10.5)
444 9500  FORMAT(i8,4e18.10)     
445       !
446   END SUBROUTINE p4z_chk_mass
447
448   SUBROUTINE p4z_clo   
449      !!---------------------------------------------------------------------
450      !!                  ***  ROUTINE p4z_clo  ***
451      !!
452      !! ** Purpose :   Closed sea domain initialization
453      !!
454      !! ** Method  :   if a closed sea is located only in a model grid point
455      !!                we restore to initial data
456      !!
457      !! ** Action  :   ictsi1(), ictsj1() : south-west closed sea limits (i,j)
458      !!                ictsi2(), ictsj2() : north-east Closed sea limits (i,j)
459      !!----------------------------------------------------------------------
460      INTEGER, PARAMETER           ::   npicts   = 4        ! number of closed sea
461      INTEGER, DIMENSION(npicts)   ::   ictsi1, ictsj1      ! south-west closed sea limits (i,j)
462      INTEGER, DIMENSION(npicts)   ::   ictsi2, ictsj2      ! north-east closed sea limits (i,j)
463      INTEGER :: ji, jj, jk, jn, jl, jc                     ! dummy loop indices
464      INTEGER :: ierr                                       ! local integer
465      REAL(wp), POINTER, DIMENSION(:,:,:,:) ::  ztrcdta     ! 4D  workspace
466      !!----------------------------------------------------------------------
467
468      IF(lwp) WRITE(numout,*)
469      IF(lwp) WRITE(numout,*)' p4z_clo : closed seas '
470      IF(lwp) WRITE(numout,*)'~~~~~~~'
471
472      ! initial values
473      ictsi1(:) = 1  ;  ictsi2(:) = 1 
474      ictsj1(:) = 1  ;  ictsj2(:) = 1 
475
476      ! set the closed seas (in data domain indices)
477      ! -------------------
478
479      IF( cp_cfg == "orca" ) THEN
480         !
481         SELECT CASE ( jp_cfg )
482         !                                           ! =======================
483         CASE ( 2 )                                  !  ORCA_R2 configuration
484            !                                        ! =======================
485            !                                            ! Caspian Sea
486            ictsi1(1)   =  11  ;  ictsj1(1)   = 103
487            ictsi2(1)   =  17  ;  ictsj2(1)   = 112
488            !                                            ! Great North American Lakes
489            ictsi1(2)   =  97  ;  ictsj1(2)   = 107
490            ictsi2(2)   = 103  ;  ictsj2(2)   = 111
491            !                                            ! Black Sea 1 : west part of the Black Sea
492            ictsi1(3)   = 174  ; ictsj1(3)   = 107
493            ictsi2(3)   = 181  ; ictsj2(3)   = 112
494            !                                            ! Black Sea 2 : est part of the Black Sea
495            ictsi1(4)   =   2  ;  ictsj1(4)   = 107
496            ictsi2(4)   =   6  ;  ictsj2(4)   = 112
497            !                                        ! =======================
498         CASE ( 4 )                                  !  ORCA_R4 configuration
499            !                                        ! =======================
500            !                                            ! Caspian Sea
501            ictsi1(1)   =  4  ;  ictsj1(1)   = 53
502            ictsi2(1)   =  4  ;  ictsj2(1)   = 56
503            !                                            ! Great North American Lakes
504            ictsi1(2)   = 49  ;  ictsj1(2)   = 55
505            ictsi2(2)   = 51  ;  ictsj2(2)   = 56
506            !                                            ! Black Sea
507            ictsi1(3)   = 88  ;  ictsj1(3)   = 55
508            ictsi2(3)   = 91  ;  ictsj2(3)   = 56
509            !                                            ! Baltic Sea
510            ictsi1(4)   = 75  ;  ictsj1(4)   = 59
511            ictsi2(4)   = 76  ;  ictsj2(4)   = 61
512            !                                        ! =======================
513            !                                        ! =======================
514         CASE ( 025 )                                ! ORCA_R025 configuration
515            !                                        ! =======================
516                                                     ! Caspian + Aral sea
517            ictsi1(1)   = 1330 ; ictsj1(1)   = 645
518            ictsi2(1)   = 1400 ; ictsj2(1)   = 795
519            !                                        ! Azov Sea
520            ictsi1(2)   = 1284 ; ictsj1(2)   = 722
521            ictsi2(2)   = 1304 ; ictsj2(2)   = 747
522            !
523         END SELECT
524         !
525      ENDIF
526
527      ! convert the position in local domain indices
528      ! --------------------------------------------
529      DO jc = 1, npicts 
530         ictsi1(jc)   = mi0( ictsi1(jc) )
531         ictsj1(jc)   = mj0( ictsj1(jc) )
532
533         ictsi2(jc)   = mi1( ictsi2(jc) )
534         ictsj2(jc)   = mj1( ictsj2(jc) )
535      END DO
536
537      ! Restore close seas values to initial data
538      IF( ln_trcdta .AND. nb_trcdta > 0 )  THEN   ! Initialisation of tracer from a file that may also be used for damping
539         !
540         CALL wrk_alloc( jpi, jpj, jpk, nb_trcdta, ztrcdta )   ! Memory allocation
541         !
542         CALL trc_dta( nittrc000, ztrcdta )   ! read tracer data at nittrc000
543         !
544         DO jn = 1, jptra
545            IF( ln_trc_ini(jn) ) THEN      ! update passive tracers arrays with input data read from file
546                jl = n_trc_index(jn)
547                DO jc = 1, npicts
548                   DO jk = 1, jpkm1
549                      DO jj = ictsj1(jc), ictsj2(jc)
550                         DO ji = ictsi1(jc), ictsi2(jc)
551                            trn(ji,jj,jk,jn) = ztrcdta(ji,jj,jk,jl) * tmask(ji,jj,jk)
552                            trb(ji,jj,jk,jn) = trn(ji,jj,jk,jn)
553                         ENDDO
554                      ENDDO
555                   ENDDO
556                ENDDO
557             ENDIF
558          ENDDO
559          CALL wrk_dealloc( jpi, jpj, jpk, nb_trcdta, ztrcdta )
560      ENDIF
561      !
562   END SUBROUTINE p4z_clo
563#else
564   !!======================================================================
565   !!  Dummy module :                                   No PISCES bio-model
566   !!======================================================================
567CONTAINS
568   SUBROUTINE p4z_sms( kt )                   ! Empty routine
569      INTEGER, INTENT( in ) ::   kt
570      WRITE(*,*) 'p4z_sms: You should not have seen this print! error?', kt
571   END SUBROUTINE p4z_sms
572#endif 
573
574   !!======================================================================
575END MODULE p4zsms 
Note: See TracBrowser for help on using the repository browser.