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.
trcini.pisces.h90 in trunk/NEMO/TOP_SRC/SMS – NEMO

source: trunk/NEMO/TOP_SRC/SMS/trcini.pisces.h90 @ 336

Last change on this file since 336 was 336, checked in by opalod, 19 years ago

nemo_v1_update_024 : CE + RB + CT : new evolution of modules
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 14.9 KB
Line 
1   !!----------------------------------------------------------------------
2   !!                    ***  trcini.pisces.h90 ***
3   !!----------------------------------------------------------------------
4#  include "domzgr_substitute.h90"
5#  include "passivetrc_substitute.h90"
6CONTAINS
7
8   SUBROUTINE trc_ini
9      !!-----------------------------------------------------------------
10      !!
11      !!                   ***  ROUTINE trc_ini ***
12      !!                     
13      !!
14      !!  Purpose :
15      !!  ---------
16      !!     Initialisation of PISCES biological and chemical variables
17      !!
18      !!   INPUT :
19      !!   -----
20      !!      common
21      !!              all the common defined in opa
22      !!
23      !!
24      !!   OUTPUT :                   : no
25      !!   ------
26      !!
27      !!   EXTERNAL :
28      !!   ----------
29      !!         p4zche
30      !!
31      !!   MODIFICATIONS:
32      !!   --------------
33      !!      original  : 1988-07  E. MAIER-REIMER      MPI HAMBURG
34      !!      additions : 1999-10  O. Aumont and C. Le Quere
35      !!      additions : 2002     O. Aumont (PISCES)
36      !!     03-2005 O. Aumont and A. El Moussaoui F90
37      !!---------------------------------------------------------------------
38      !! local declarations
39      !! ==================
40      INTEGER :: ji,jj,jk
41      INTEGER :: ichl,iband,mo
42      INTEGER , PARAMETER :: jpmois = 12,      &
43         jpan   = 1
44
45      REAL(wp) :: xtoto,expide,denitide,ztra,zmaskt
46      REAL(wp) , DIMENSION (jpi,jpj) :: riverdoc,river,ndepo
47      CHARACTER (len=34) :: clname
48
49      INTEGER :: ipi,ipj,ipk,itime
50      INTEGER , DIMENSION (jpmois) :: istep
51      INTEGER , DIMENSION (jpan) :: istep0
52      REAL(wp) :: zsecond, zdate0
53      REAL(wp) , DIMENSION (jpi,jpj) :: zlon,zlat
54      REAL(wp), DIMENSION (jpk) :: zlev
55      INTEGER :: numriv,numdust,numbath,numdep
56
57      !! 1. initialization
58      !! -----------------
59
60      !! computation of the record length for direct access FILE
61      !! this length depend of 512 for the t3d machine
62      !!
63      rfact = rdttra(1) * float(ndttrc)
64      rfactr = 1./rfact
65      IF(lwp) WRITE(numout,*) ' Tracer time step=',rfact,' rdt=',rdt
66      rfact2= rfact / float(nrdttrc)
67      rfact2r = 1./rfact2
68      IF(lwp) write(numout,*) ' Biology time step=',rfact2
69
70      !!    INITIALISE DUST INPUT FROM ATMOSPHERE
71      !!    -------------------------------------
72
73      IF (bdustfer) THEN
74         clname='dust.orca.nc'
75         CALL flinopen(clname,mig(1),nlci,mjg(1),nlcj,.false.,ipi,ipj,0        &
76            &      ,zlon,zlat,zlev,itime,istep,zdate0,zsecond,numdust)
77         CALL flinget(numdust,'dust',jpidta,jpjdta,0,jpmois,1,                 &
78            &        12,mig(1),nlci,mjg(1),nlcj,dustmo(1:nlci,1:nlcj,:) )
79         CALL flinclo(numdust)
80
81         ! Extra-halo initialization in MPP
82         IF( lk_mpp ) THEN
83            DO ji = nlci+1, jpi
84               dustmo(ji,:,:) = dustmo(1,:,:)
85            ENDDO
86            DO jj = nlcj+1, jpj
87               dustmo(:,jj,:)=dustmo(:,1,:)
88            ENDDO
89         ENDIF
90      ELSE
91         dustmo(:,:,:)=0.
92      ENDIF
93
94      !!    INITIALISE THE NUTRIENT INPUT BY RIVERS
95      !!    ---------------------------------------
96
97      IF (briver) THEN
98         clname='river.orca.nc'
99         CALL flinopen(clname,mig(1),nlci,mjg(1),nlcj,.false.,ipi,ipj,0        &
100            &      ,zlon,zlat,zlev,itime,istep0,zdate0,zsecond,numriv)
101         CALL flinget(numriv,'riverdic',jpidta,jpjdta,0,jpan,1,                &
102            &        1,mig(1),nlci,mjg(1),nlcj,river(1:nlci,1:nlcj) )
103         CALL flinget(numriv,'riverdoc',jpidta,jpjdta,0,jpan,1,                &
104            &        1,mig(1),nlci,mjg(1),nlcj,riverdoc(1:nlci,1:nlcj) )
105         CALL flinclo(numriv)
106
107         ! Extra-halo initialization in MPP
108         IF( lk_mpp ) THEN
109            DO ji = nlci+1, jpi
110               river(ji,:) = river(1,:)
111               riverdoc(ji,:) = riverdoc(1,:)
112            ENDDO
113            DO jj = nlcj+1, jpj
114               river(:,jj)=river(:,1)
115               riverdoc(:,jj) = riverdoc(:,1)
116            ENDDO
117         ENDIF
118
119      ELSE
120         river(:,:)=0.
121         riverdoc(:,:)=0.
122      endif
123
124      !!    INITIALISE THE N INPUT BY DUST
125      !!  ---------------------------------------
126
127      IF (bndepo) THEN
128         clname='ndeposition.orca.nc'
129         CALL flinopen(clname,mig(1),nlci,mjg(1),nlcj,.false.,ipi,ipj,0        &
130            &      ,zlon,zlat,zlev,itime,istep0,zdate0,zsecond,numdep)
131         CALL flinget(numdep,'ndep',jpidta,jpjdta,0,jpan,1,                   &
132            &        1,mig(1),nlci,mjg(1),nlcj,ndepo(1:nlci,1:nlcj) )
133         CALL flinclo(numdep)
134
135         ! Extra-halo initialization in MPP
136         IF( lk_mpp ) THEN
137            DO ji = nlci+1, jpi
138               ndepo(ji,:) = ndepo(1,:)
139            ENDDO
140            DO jj = nlcj+1, jpj
141               ndepo(:,jj)=ndepo(:,1)
142            ENDDO
143         ENDIF
144
145      ELSE
146         ndepo(:,:)=0.
147      ENDIF
148
149      !!    Computation of the coastal mask.
150      !!    Computation of an island mask to enhance coastal supply
151      !!    of iron
152      !!    -------------------------------------------------------
153
154      IF (bsedinput) THEN
155         clname='bathy.orca.nc'
156         CALL flinopen(clname,mig(1),nlci,mjg(1),nlcj,.false.,ipi,ipj,ipk      &
157            &      ,zlon,zlat,zlev,itime,istep0,zdate0,zsecond,numbath)
158         CALL flinget(numbath,'bathy',jpidta,jpjdta,jpk,jpan,1,                &
159            &        1,mig(1),nlci,mjg(1),nlcj,cmask(1:nlci,1:nlcj,1:jpk) )
160         CALL flinclo(numbath)
161
162         do jk=1,5
163            do jj=2,jpj-1
164               do ji=2,jpi-1
165                  if (tmask(ji,jj,jk).ne.0) then
166                     zmaskt=tmask(ji+1,jj,jk)*tmask(ji-1,jj,jk)*tmask(ji,jj+1,jk)    &
167                        &          *tmask(ji,jj-1,jk)*tmask(ji,jj,jk+1)
168                     if (zmaskt.eq.0) then
169                        cmask(ji,jj,jk)=0.1
170                     endif
171                  endif
172               end do
173            end do
174         end do
175
176
177         ! Extra-halo initialization in MPP
178         IF( lk_mpp ) THEN
179            DO ji = nlci+1, jpi
180               cmask(ji,:,:) = cmask(1,:,:)
181            ENDDO
182            DO jj = nlcj+1, jpj
183               cmask(:,jj,:)=cmask(:,1,:)
184            ENDDO
185         ENDIF
186
187         DO jk = 1, jpk
188            DO jj = 1, jpj
189               DO ji = 1, jpi
190                  expide=min(8.,(fsdept(ji,jj,jk)/500.)**(-1.5))
191                  denitide=-0.9543+0.7662*log(expide)-0.235*log(expide)**2
192                  cmask(ji,jj,jk)=cmask(ji,jj,jk)*min(1.,exp(denitide)/0.5)
193               END DO
194            END DO
195         END DO
196
197      ELSE
198         cmask(:,:,:)=0.
199      ENDIF
200
201      ! Lateral boundary conditions on ( avt, en )   (sign unchanged)
202      CALL lbc_lnk( cmask , 'T', 1. )
203
204      !!     Computation of the total atmospheric supply of Si
205      !!     -------------------------------------------------
206
207      sumdepsi=0.
208      DO mo=1,12
209         DO jj=2,jpjm1
210            DO ji=2,jpim1
211               sumdepsi=sumdepsi+dustmo(ji,jj,mo)/(12.*rmoss)*8.8        &
212                  *0.075/28.1*e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,1)
213            END DO
214         END DO
215      END DO
216
217      IF( lk_mpp )   CALL mpp_sum( sumdepsi )  ! sum over the global domain
218
219      !!    COMPUTATION OF THE N/P RELEASE DUE TO COASTAL RIVERS
220      !!    COMPUTATION OF THE Si RELEASE DUE TO COASTAL RIVERS
221      !!    ---------------------------------------------------
222
223      DO jj=1,jpj
224         DO ji=1,jpi
225            cotdep(ji,jj)=river(ji,jj)*1E9/(12.*raass                          &
226               *e1t(ji,jj)*e2t(ji,jj)*fse3t(ji,jj,1)+rtrn)*tmask(ji,jj,1)
227            rivinp(ji,jj)=(river(ji,jj)+riverdoc(ji,jj))*1E9                   &
228               /(31.6*raass*e1t(ji,jj)*e2t(ji,jj)*fse3t(ji,jj,1)+rtrn)   &
229               *tmask(ji,jj,1)
230            nitdep(ji,jj)=7.6*ndepo(ji,jj)*tmask(ji,jj,1)/(14E6*raass          &
231               *fse3t(ji,jj,1)+rtrn)
232         END DO
233      END DO
234      ! Lateral boundary conditions on ( cotdep, rivinp, nitdep )   (sign unchanged)
235      CALL lbc_lnk( cotdep , 'T', 1. )  ;  CALL lbc_lnk( rivinp , 'T', 1. )  ;  CALL lbc_lnk( nitdep , 'T', 1. )
236
237      rivpo4input=0.
238      rivalkinput=0.
239      nitdepinput=0.
240      DO jj=2,jpjm1
241         DO ji=2,jpim1
242            rivpo4input=rivpo4input+rivinp(ji,jj)*(e1t(ji,jj)*e2t(ji,jj)    &
243               *fse3t(ji,jj,1))*tmask(ji,jj,1)*raass
244            rivalkinput=rivalkinput+cotdep(ji,jj)*(e1t(ji,jj)*e2t(ji,jj)    &
245               *fse3t(ji,jj,1))*tmask(ji,jj,1)*raass
246            nitdepinput=nitdepinput+nitdep(ji,jj)*(e1t(ji,jj)*e2t(ji,jj)    &
247               *fse3t(ji,jj,1))*tmask(ji,jj,1)*raass
248         END DO
249      END DO
250
251      IF( lk_mpp ) THEN
252         CALL mpp_sum( rivpo4input )  ! sum over the global domain
253         CALL mpp_sum( rivalkinput )  ! sum over the global domain
254         CALL mpp_sum( nitdepinput )  ! sum over the global domain
255      ENDIF
256
257
258      !!    Coastal supply of iron
259      !!    ----------------------
260
261      DO jk=1,jpkm1
262         ironsed(:,:,jk)=sedfeinput*cmask(:,:,jk)         &
263            /(fse3t(:,:,jk)*rjjss)
264      END DO
265
266      ! Lateral boundary conditions on ( ironsed )   (sign unchanged)
267      CALL lbc_lnk( ironsed , 'T', 1. )
268      !!----------------------------------------------------------------------
269      !!
270      !! Initialize biological variables
271      !!
272      !!----------------------------------------------------------------------
273      !! Set biological ratios
274      !! ---------------------
275
276      rno3   = (16.+2.)/122.
277      po4r   = 1./122.
278      o2nit  = 32./122.
279      rdenit = 97.6/16.
280      o2ut   = 140./122.
281
282      !!----------------------------------------------------------------------
283      !!
284      !! Initialize chemical variables
285      !!
286      !!----------------------------------------------------------------------
287
288      !! set pre-industrial atmospheric [co2] (ppm) and o2/n2 ratio
289      !! ----------------------------------------------------------
290
291      atcox = 0.20946
292
293      !! Set lower/upper limits for temperature and salinity
294      !! ---------------------------------------------------
295
296      salchl = 1./1.80655
297      calcon = 1.03E-2
298
299      !! Set coefficients for apparent solubility equilibrium
300      !!   of calcite (Ingle, 1800, eq. 6)
301      !! ----------------------------------------------------
302
303      akcc1 = -34.452
304      akcc2 = -39.866
305      akcc3 = 110.21
306      akcc4 = -7.5752E-6
307
308
309      !! Set coefficients for seawater pressure correction
310      !! -------------------------------------------------
311
312      devk1(1) = -25.5
313      devk2(1) = 0.1271
314      devk3(1) = 0.
315      devk4(1) = -3.08E-3
316      devk5(1) = 0.0877E-3
317
318      devk1(2) = -15.82
319      devk2(2) = -0.0219
320      devk3(2) = 0.
321      devk4(2) = 1.13E-3
322      devk5(2) = -0.1475E-3
323
324      devk1(3) = -29.48
325      devk2(3) = 0.1622
326      devk3(3) = 2.608E-3
327      devk4(3) = -2.84E-3
328      devk5(3) = 0.
329
330      devk1(4) = -25.60
331      devk2(4) = 0.2324
332      devk3(4) = -3.6246E-3
333      devk4(4) = -5.13E-3
334      devk5(4) = 0.0794E-3
335
336      devkst = 0.23
337      devks  = 35.4
338
339      !! Set universal gas constants
340      !! ---------------------------
341
342      rgas = 83.143
343      oxyco = 1./22.4144
344
345      !! Set boron constants
346      !! -------------------
347
348      bor1 = 0.00023
349      bor2 = 1./10.82
350
351      !! Set volumetric solubility constants for co2 in ml/l (Weiss, 1974)
352      !! -----------------------------------------------------------------
353
354      c00 = -60.2409
355      c01 = 93.4517
356      c02 = 23.3585
357      c03 = 0.023517
358      c04 = -0.023656
359      c05 = 0.0047036
360
361      ca0 = -162.8301
362      ca1 = 218.2968
363      ca2 = 90.9241
364      ca3 = -1.47696
365      ca4 = 0.025695
366      ca5 = -0.025225
367      ca6 = 0.0049867
368
369      !! Set coeff. for 1. dissoc. of carbonic acid (Edmond and Gieskes, 1970)
370      !! ---------------------------------------------------------------------
371
372      c10 = -3670.7
373      c11 =  62.008
374      c12 = -9.7944
375      c13 = 0.0118
376      c14 = -0.000116
377
378      !! Set coeff. for 2. dissoc. of carbonic acid (Edmond and Gieskes, 1970)
379      !! ---------------------------------------------------------------------
380
381      c20 = -1394.7
382      c21 = -4.777
383      c22 = 0.0184
384      c23 = -0.000118
385
386      !! Set coeff. for 1. dissoc. of boric acid (Edmond and Gieskes, 1970)
387      !! ------------------------------------------------------------------
388
389      cb0  = -8966.90
390      cb1  = -2890.53
391      cb2  = -77.942
392      cb3  = 1.728
393      cb4  = -0.0996
394      cb5  = 148.0248
395      cb6  = 137.1942
396      cb7  = 1.62142
397      cb8  = -24.4344
398      cb9  = -25.085
399      cb10 = -0.2474
400      cb11 = 0.053105
401
402      !! Set coeff. for dissoc. of water (Dickson and Riley, 1979,
403      !!   eq. 7, coefficient cw2 corrected from 0.9415 to 0.09415
404      !!   after pers. commun. to B. Bacastow, 1988)
405      !! ---------------------------------------------------------
406
407      cw0 = -13847.26
408      cw1 = 148.9652
409      cw2 = -23.6521
410      cw3 = 118.67
411      cw4 = -5.977
412      cw5 = 1.0495
413      cw6 = -0.01615
414
415      !
416      ! Set coeff. for dissoc. of phosphate (Millero (1974)
417      ! ---------------------------------------------------
418      !
419      cp10 = 115.525
420      cp11 = -4576.752
421      cp12 = -18.453
422      cp13 = -106.736
423      cp14 = 0.69171
424      cp15 = -0.65643
425      cp16 = -0.01844
426
427      cp20 = 172.0883
428      cp21 = -8814.715
429      cp22 = -27.927
430      cp23 = -160.340
431      cp24 = 1.3566
432      cp25 = 0.37335
433      cp26 = -0.05778
434
435
436      cp30 = -18.141
437      cp31 = -3070.75
438      cp32 = 17.27039
439      cp33 = 2.81197
440      cp34 = -44.99486
441      cp35 = -0.09984
442      !
443      ! Set coeff. for dissoc. of phosphate (Millero (1974)
444      ! ---------------------------------------------------
445      !
446      cs10 = 117.385
447      cs11 = -8904.2
448      cs12 = -19.334
449      cs13 = -458.79
450      cs14 =  3.5913
451      cs15 = 188.74
452      cs16 = -1.5998
453      cs17 = -12.1652
454      cs18 = 0.07871
455      cs19 = -0.001005
456
457      !! Set volumetric solubility constants for o2 in ml/l (Weiss, 1970)
458      !! ----------------------------------------------------------------
459
460      ox0 = -58.3877
461      ox1 = 85.8079
462      ox2 = 23.8439
463      ox3 = -0.034892
464      ox4 = 0.015568
465      ox5 = -0.0019387
466
467      !!  FROM THE NEW BIOOPTIC MODEL PROPOSED JM ANDRE, WE READ HERE
468      !!  A PRECOMPUTED ARRAY CORRESPONDING TO THE ATTENUATION COEFFICIENT
469
470      open(49,file='kRGB61.txt',form='formatted')
471      do ichl=1,61
472         READ(49,*) xtoto,(xkrgb(iband,ichl),iband = 1,3)
473      end do
474      close(49)
475
476#if defined key_off_degrad
477
478      !! Read volume for degraded regions (DEGINIT)
479      !! ------------------------------------------
480
481#    if defined key_vpp
482      CALL READ3S(902,facvol,jpi,jpj,jpk)
483#    else
484      READ (902) facvol
485#    endif
486#endif
487
488
489      !!  Call p4zche to initialize the chemical constants
490      !!  ------------------------------------------------
491
492      CALL p4zche
493      !!
494      !!  Initialize a counter for the computation of chemistry
495      !!
496      ndayflxtr=0
497
498      IF(lwp) WRITE(numout,*) ' Initialisation of PISCES done'
499
500   END SUBROUTINE trc_ini
Note: See TracBrowser for help on using the repository browser.