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.
Changeset 7256 for branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC – NEMO

Ignore:
Timestamp:
2016-11-18T08:18:45+01:00 (7 years ago)
Author:
cbricaud
Message:

phaze NEMO routines in CRS branch with nemo_v3_6_STABLE branch at rev 7213 (09-09-2016) (merge -r 5519:7213 )

Location:
branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC
Files:
1 deleted
35 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zbio.F90

    r5602 r7256  
    599599 
    600600   !!====================================================================== 
    601 END MODULE  p2zbio 
     601END MODULE p2zbio 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zsms.F90

    r5602 r7256  
    8484 
    8585   !!====================================================================== 
    86 END MODULE  p2zsms 
     86END MODULE p2zsms 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zbio.F90

    r5602 r7256  
    109109 
    110110   !!====================================================================== 
    111 END MODULE  p4zbio 
    112  
     111END MODULE p4zbio 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zche.F90

    r5602 r7256  
    3232   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   fekeq    ! chemistry of Fe 
    3333   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   chemc    ! Solubilities of O2 and CO2 
    34    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   chemo2    ! Solubilities of O2 and CO2 
     34   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   chemo2   ! Solubilities of O2 and CO2 
     35   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   tempis   ! In situ temperature 
    3536 
    3637   REAL(wp), PUBLIC ::   atcox  = 0.20946         ! units atm 
     
    3940   REAL(wp) ::   o2atm  = 1. / ( 1000. * 0.20946 )   
    4041 
    41    REAL(wp) ::   akcc1  = -171.9065       ! coeff. for apparent solubility equilibrium 
    42    REAL(wp) ::   akcc2  =   -0.077993     ! Millero et al. 1995 from Mucci 1983 
    43    REAL(wp) ::   akcc3  = 2839.319         
    44    REAL(wp) ::   akcc4  =   71.595         
    45    REAL(wp) ::   akcc5  =   -0.77712       
    46    REAL(wp) ::   akcc6  =    0.00284263    
    47    REAL(wp) ::   akcc7  =  178.34         
    48    REAL(wp) ::   akcc8  =   -0.07711      
    49    REAL(wp) ::   akcc9  =    0.0041249    
    50  
    51    REAL(wp) ::   rgas   = 83.143         ! universal gas constants 
     42   REAL(wp) ::   rgas   = 83.14472       ! universal gas constants 
    5243   REAL(wp) ::   oxyco  = 1. / 22.4144   ! converts from liters of an ideal gas to moles 
    5344 
     
    5546   REAL(wp) ::   bor2   = 1. / 10.82 
    5647 
    57    REAL(wp) ::   ca0    = -162.8301      ! WEISS & PRICE 1980, units mol/(kg atm) 
    58    REAL(wp) ::   ca1    =  218.2968 
    59    REAL(wp) ::   ca2    =   90.9241 
    60    REAL(wp) ::   ca3    =   -1.47696 
    61    REAL(wp) ::   ca4    =    0.025695 
    62    REAL(wp) ::   ca5    =   -0.025225 
    63    REAL(wp) ::   ca6    =    0.0049867 
    64  
    65    REAL(wp) ::   c10    = -3670.7        ! Coeff. for 1. dissoc. of carbonic acid (Edmond and Gieskes, 1970)    
    66    REAL(wp) ::   c11    =    62.008      
    67    REAL(wp) ::   c12    =    -9.7944     
    68    REAL(wp) ::   c13    =     0.0118      
    69    REAL(wp) ::   c14    =    -0.000116 
    70  
    71    REAL(wp) ::   c20    = -1394.7       ! coeff. for 2. dissoc. of carbonic acid (Millero, 1995)    
    72    REAL(wp) ::   c21    =    -4.777    
    73    REAL(wp) ::   c22    =     0.0184    
    74    REAL(wp) ::   c23    =    -0.000118 
    75  
    7648   REAL(wp) ::   st1    =      0.14     ! constants for calculate concentrations for sulfate 
    7749   REAL(wp) ::   st2    =  1./96.062    !  (Morris & Riley 1966) 
    78    REAL(wp) ::   ks0    =    141.328  
    79    REAL(wp) ::   ks1    =  -4276.1   
    80    REAL(wp) ::   ks2    =    -23.093 
    81    REAL(wp) ::   ks3    = -13856.   
    82    REAL(wp) ::   ks4    =   324.57  
    83    REAL(wp) ::   ks5    =   -47.986 
    84    REAL(wp) ::   ks6    =  35474.  
    85    REAL(wp) ::   ks7    =   -771.54 
    86    REAL(wp) ::   ks8    =    114.723 
    87    REAL(wp) ::   ks9    =  -2698.   
    88    REAL(wp) ::   ks10   =   1776.  
    89    REAL(wp) ::   ks11   =      1. 
    90    REAL(wp) ::   ks12   =     -0.001005  
    9150 
    9251   REAL(wp) ::   ft1    =    0.000067   ! constants for calculate concentrations for fluorides 
    9352   REAL(wp) ::   ft2    = 1./18.9984    ! (Dickson & Riley 1979 ) 
    94    REAL(wp) ::   kf0    =  -12.641     
    95    REAL(wp) ::   kf1    = 1590.2     
    96    REAL(wp) ::   kf2    =    1.525     
    97    REAL(wp) ::   kf3    =    1.0      
    98    REAL(wp) ::   kf4    =   -0.001005 
    99  
    100    REAL(wp) ::   cb0    = -8966.90      ! Coeff. for 1. dissoc. of boric acid  
    101    REAL(wp) ::   cb1    = -2890.53      ! (Dickson and Goyet, 1994) 
    102    REAL(wp) ::   cb2    =   -77.942 
    103    REAL(wp) ::   cb3    =     1.728 
    104    REAL(wp) ::   cb4    =    -0.0996 
    105    REAL(wp) ::   cb5    =   148.0248 
    106    REAL(wp) ::   cb6    =   137.1942 
    107    REAL(wp) ::   cb7    =     1.62142 
    108    REAL(wp) ::   cb8    =   -24.4344 
    109    REAL(wp) ::   cb9    =   -25.085 
    110    REAL(wp) ::   cb10   =    -0.2474  
    111    REAL(wp) ::   cb11   =     0.053105 
    112  
    113    REAL(wp) ::   cw0    = -13847.26     ! Coeff. for dissoc. of water (Dickson and Riley, 1979 ) 
    114    REAL(wp) ::   cw1    =    148.9652   
    115    REAL(wp) ::   cw2    =    -23.6521 
    116    REAL(wp) ::   cw3    =    118.67  
    117    REAL(wp) ::   cw4    =     -5.977  
    118    REAL(wp) ::   cw5    =      1.0495   
    119    REAL(wp) ::   cw6    =     -0.01615 
    12053 
    12154   !                                    ! volumetric solubility constants for o2 in ml/L   
     
    185118      REAL(wp) ::   ztgg , ztgg2, ztgg3 , ztgg4 , ztgg5 
    186119      REAL(wp) ::   zpres, ztc  , zcl   , zcpexp, zoxy  , zcpexp2 
    187       REAL(wp) ::   zsqrt, ztr  , zlogt , zcek1 
    188       REAL(wp) ::   zis  , zis2 , zsal15, zisqrt 
     120      REAL(wp) ::   zsqrt, ztr  , zlogt , zcek1, zc1, zplat 
     121      REAL(wp) ::   zis  , zis2 , zsal15, zisqrt, za1  , za2 
    189122      REAL(wp) ::   zckb , zck1 , zck2  , zckw  , zak1 , zak2  , zakb , zaksp0, zakw 
    190123      REAL(wp) ::   zst  , zft  , zcks  , zckf  , zaksp1 
     
    193126      IF( nn_timing == 1 )  CALL timing_start('p4z_che') 
    194127      ! 
     128      ! Computations of chemical constants require in situ temperature 
     129      ! Here a quite simple formulation is used to convert  
     130      ! potential temperature to in situ temperature. The errors is less than  
     131      ! 0.04°C relative to an exact computation 
     132      ! --------------------------------------------------------------------- 
     133      DO jk = 1, jpk 
     134         DO jj = 1, jpj 
     135            DO ji = 1, jpi 
     136               zpres = fsdept(ji,jj,jk) / 1000. 
     137               za1 = 0.04 * ( 1.0 + 0.185 * tsn(ji,jj,jk,jp_tem) + 0.035 * (tsn(ji,jj,jk,jp_sal) - 35.0) ) 
     138               za2 = 0.0075 * ( 1.0 - tsn(ji,jj,jk,jp_tem) / 30.0 ) 
     139               tempis(ji,jj,jk) = tsn(ji,jj,jk,jp_tem) - za1 * zpres + za2 * zpres**2 
     140            END DO 
     141         END DO 
     142      END DO 
     143      ! 
    195144      ! CHEMICAL CONSTANTS - SURFACE LAYER 
    196145      ! ---------------------------------- 
     
    200149         DO ji = 1, jpi 
    201150            !                             ! SET ABSOLUTE TEMPERATURE 
    202             ztkel = tsn(ji,jj,1,jp_tem) + 273.16 
     151            ztkel = tempis(ji,jj,1) + 273.15 
    203152            zt    = ztkel * 0.01 
    204153            zt2   = zt * zt 
     
    208157            !                             ! LN(K0) OF SOLUBILITY OF CO2 (EQ. 12, WEISS, 1980) 
    209158            !                             !     AND FOR THE ATMOSPHERE FOR NON IDEAL GAS 
    210             zcek1 = ca0 + ca1 / zt + ca2 * zlogt + ca3 * zt2 + zsal * ( ca4 + ca5 * zt + ca6 * zt2 ) 
    211             !                             ! LN(K0) OF SOLUBILITY OF O2 and N2 in ml/L (EQ. 8, GARCIA AND GORDON, 1992) 
    212             ztgg  = LOG( ( 298.15 - tsn(ji,jj,1,jp_tem) ) / ztkel )  ! Set the GORDON & GARCIA scaled temperature 
    213             ztgg2 = ztgg  * ztgg 
    214             ztgg3 = ztgg2 * ztgg 
    215             ztgg4 = ztgg3 * ztgg 
    216             ztgg5 = ztgg4 * ztgg 
    217             zoxy  = ox0 + ox1 * ztgg + ox2 * ztgg2 + ox3 * ztgg3 + ox4 * ztgg4 + ox5 * ztgg5   & 
    218                    + zsal * ( ox6 + ox7 * ztgg + ox8 * ztgg2 + ox9 * ztgg3 ) +  ox10 * zsal2 
    219  
     159            zcek1 = 9345.17/ztkel - 60.2409 + 23.3585 * LOG(zt) + zsal*(0.023517 - 0.00023656*ztkel    & 
     160            &       + 0.0047036e-4*ztkel**2) 
    220161            !                             ! SET SOLUBILITIES OF O2 AND CO2  
    221             chemc(ji,jj,1) = EXP( zcek1 ) * 1.e-6 * rhop(ji,jj,1) / 1000.  ! mol/(L uatm) 
    222             chemc(ji,jj,2) = ( EXP( zoxy  ) * o2atm ) * oxyco              ! mol/(L atm) 
     162            chemc(ji,jj,1) = EXP( zcek1 ) * 1.e-6 * rhop(ji,jj,1) / 1000. ! mol/(kg uatm) 
     163            chemc(ji,jj,2) = -1636.75 + 12.0408*ztkel - 0.0327957*ztkel**2 + 0.0000316528*ztkel**3 
     164            chemc(ji,jj,3) = 57.7 - 0.118*ztkel 
    223165            ! 
    224166         END DO 
     
    233175!CDIR NOVERRCHK 
    234176            DO ji = 1, jpi 
    235               ztkel = tsn(ji,jj,jk,jp_tem) + 273.16 
     177              ztkel = tempis(ji,jj,jk) + 273.15 
    236178              zsal  = tsn(ji,jj,jk,jp_sal) + ( 1.- tmask(ji,jj,jk) ) * 35. 
    237179              zsal2 = zsal * zsal 
    238               ztgg  = LOG( ( 298.15 - tsn(ji,jj,jk,jp_tem) ) / ztkel )  ! Set the GORDON & GARCIA scaled temperature 
     180              ztgg  = LOG( ( 298.15 - tempis(ji,jj,jk) ) / ztkel )  ! Set the GORDON & GARCIA scaled temperature 
    239181              ztgg2 = ztgg  * ztgg 
    240182              ztgg3 = ztgg2 * ztgg 
     
    259201            DO ji = 1, jpi 
    260202 
    261                ! SET PRESSION 
    262                zpres   = 1.025e-1 * fsdept(ji,jj,jk) 
     203               ! SET PRESSION ACCORDING TO SAUNDER (1980) 
     204               zplat   = SIN ( ABS(gphit(ji,jj)*3.141592654/180.) ) 
     205               zc1 = 5.92E-3 + zplat**2 * 5.25E-3 
     206               zpres = ((1-zc1)-SQRT(((1-zc1)**2)-(8.84E-6*fsdept(ji,jj,jk)))) / 4.42E-6 
     207               zpres = zpres / 10.0 
    263208 
    264209               ! SET ABSOLUTE TEMPERATURE 
    265                ztkel   = tsn(ji,jj,jk,jp_tem) + 273.16 
     210               ztkel   = tempis(ji,jj,jk) + 273.15 
    266211               zsal    = tsn(ji,jj,jk,jp_sal) + ( 1.-tmask(ji,jj,jk) ) * 35. 
    267212               zsqrt  = SQRT( zsal ) 
     
    272217               zis2   = zis * zis 
    273218               zisqrt = SQRT( zis ) 
    274                ztc     = tsn(ji,jj,jk,jp_tem) + ( 1.- tmask(ji,jj,jk) ) * 20. 
     219               ztc     = tempis(ji,jj,jk) + ( 1.- tmask(ji,jj,jk) ) * 20. 
    275220 
    276221               ! CHLORINITY (WOOSTER ET AL., 1969) 
     
    284229 
    285230               ! DISSOCIATION CONSTANT FOR SULFATES on free H scale (Dickson 1990) 
    286                zcks    = EXP(  ks1 * ztr + ks0 + ks2 * zlogt                           & 
    287                   &                     + ( ks3 * ztr + ks4 + ks5 * zlogt ) * zisqrt   & 
    288                   &                     + ( ks6 * ztr + ks7 + ks8 * zlogt ) * zis      & 
    289                   &                     + ks9 * ztr * zis * zisqrt + ks10 * ztr *zis2 + LOG( ks11 + ks12 *zsal )  ) 
     231               zcks    = EXP(-4276.1 * ztr + 141.328 - 23.093 * zlogt         & 
     232               &         + (-13856. * ztr + 324.57 - 47.986 * zlogt) * zisqrt & 
     233               &         + (35474. * ztr - 771.54 + 114.723 * zlogt) * zis    & 
     234               &         - 2698. * ztr * zis**1.5 + 1776.* ztr * zis2         & 
     235               &         + LOG(1.0 - 0.001005 * zsal)) 
     236               ! 
     237               aphscale(ji,jj,jk) = ( 1. + zst / zcks ) 
    290238 
    291239               ! DISSOCIATION CONSTANT FOR FLUORIDES on free H scale (Dickson and Riley 79) 
    292                zckf    = EXP(  kf1 * ztr + kf0 + kf2 * zisqrt + LOG( kf3 + kf4 * zsal )  ) 
     240               zckf    = EXP( 1590.2*ztr - 12.641 + 1.525*zisqrt   & 
     241               &         + LOG(1.0d0 - 0.001005d0*zsal)            & 
     242               &         + LOG(1.0d0 + zst/zcks)) 
    293243 
    294244               ! DISSOCIATION CONSTANT FOR CARBONATE AND BORATE 
    295                zckb    = ( cb0 + cb1 * zsqrt + cb2  * zsal + cb3 * zsal15 + cb4 * zsal * zsal ) * ztr   & 
    296                   &    + ( cb5 + cb6 * zsqrt + cb7  * zsal )                                            & 
    297                   &    + ( cb8 + cb9 * zsqrt + cb10 * zsal ) * zlogt + cb11 * zsqrt * ztkel             & 
    298                   &    + LOG(  ( 1.+ zst / zcks + zft / zckf ) / ( 1.+ zst / zcks )  ) 
    299  
    300                zck1    = c10 * ztr + c11 + c12 * zlogt + c13 * zsal + c14 * zsal * zsal 
    301                zck2    = c20 * ztr + c21 + c22 * zsal   + c23 * zsal**2 
     245               zckb=  (-8966.90 - 2890.53*zsqrt - 77.942*zsal        & 
     246               &      + 1.728*zsal15 - 0.0996*zsal*zsal)*ztr         & 
     247               &      + (148.0248 + 137.1942*zsqrt + 1.62142*zsal)   & 
     248               &      + (-24.4344 - 25.085*zsqrt - 0.2474*zsal)      &  
     249               &      * zlogt + 0.053105*zsqrt*ztkel 
     250 
     251 
     252               ! DISSOCIATION COEFFICIENT FOR CARBONATE ACCORDING TO  
     253               ! MEHRBACH (1973) REFIT BY MILLERO (1995), seawater scale 
     254               zck1    = -1.0*(3633.86*ztr - 61.2172 + 9.6777*zlogt  & 
     255                  - 0.011555*zsal + 0.0001152*zsal*zsal) 
     256               zck2    = -1.0*(471.78*ztr + 25.9290 - 3.16967*zlogt      & 
     257                  - 0.01781*zsal + 0.0001122*zsal*zsal) 
    302258 
    303259               ! PKW (H2O) (DICKSON AND RILEY, 1979) 
    304                zckw    = cw0 * ztr + cw1 + cw2 * zlogt + ( cw3 * ztr + cw4 + cw5 * zlogt ) * zsqrt + cw6 * zsal 
    305  
     260               zckw = -13847.26*ztr + 148.9652 - 23.6521 * zlogt    &  
     261               &     + (118.67*ztr - 5.977 + 1.0495 * zlogt)        & 
     262               &     * zsqrt - 0.01615 * zsal 
    306263 
    307264               ! APPARENT SOLUBILITY PRODUCT K'SP OF CALCITE IN SEAWATER 
    308265               !       (S=27-43, T=2-25 DEG C) at pres =0 (atmos. pressure) (MUCCI 1983) 
    309                zaksp0  = akcc1 + akcc2 * ztkel + akcc3 * ztr + akcc4 * LOG10( ztkel )   & 
    310                   &   + ( akcc5 + akcc6 * ztkel + akcc7 * ztr ) * zsqrt + akcc8 * zsal + akcc9 * zsal15 
     266               zaksp0  = -171.9065 -0.077993*ztkel + 2839.319*ztr + 71.595*LOG10( ztkel )   & 
     267                  &      + (-0.77712 + 0.00284263*ztkel + 178.34*ztr) * zsqrt  & 
     268                  &      - 0.07711*zsal + 0.0041249*zsal15 
    311269 
    312270               ! K1, K2 OF CARBONIC ACID, KB OF BORIC ACID, KW (H2O) (LIT.?) 
     
    378336      !!                     ***  ROUTINE p4z_che_alloc  *** 
    379337      !!---------------------------------------------------------------------- 
    380       ALLOCATE( sio3eq(jpi,jpj,jpk), fekeq(jpi,jpj,jpk), chemc(jpi,jpj,2), chemo2(jpi,jpj,jpk), STAT=p4z_che_alloc ) 
     338      ALLOCATE( sio3eq(jpi,jpj,jpk), fekeq(jpi,jpj,jpk), chemc(jpi,jpj,3), chemo2(jpi,jpj,jpk),   & 
     339      &         tempis(jpi,jpj,jpk), STAT=p4z_che_alloc ) 
    381340      ! 
    382341      IF( p4z_che_alloc /= 0 )   CALL ctl_warn('p4z_che_alloc : failed to allocate arrays.') 
     
    396355 
    397356   !!====================================================================== 
    398 END MODULE  p4zche 
     357END MODULE p4zche 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zflx.F90

    r5602 r7256  
    8484      ! 
    8585      INTEGER  ::   ji, jj, jm, iind, iindm1 
    86       REAL(wp) ::   ztc, ztc2, ztc3, zws, zkgwan 
     86      REAL(wp) ::   ztc, ztc2, ztc3, ztc4, zws, zkgwan 
    8787      REAL(wp) ::   zfld, zflu, zfld16, zflu16, zfact 
     88      REAL(wp) ::   zvapsw, zsal, zfco2, zxc2, xCO2approx, ztkel, zfugcoeff 
    8889      REAL(wp) ::   zph, zah2, zbot, zdic, zalk, zsch_o2, zalka, zsch_co2 
    8990      REAL(wp) ::   zyr_dec, zdco2dt 
    9091      CHARACTER (len=25) :: charout 
    91       REAL(wp), POINTER, DIMENSION(:,:) :: zkgco2, zkgo2, zh2co3, zoflx, zw2d  
     92      REAL(wp), POINTER, DIMENSION(:,:) :: zkgco2, zkgo2, zh2co3, zoflx, zw2d, zpco2atm  
    9293      !!--------------------------------------------------------------------- 
    9394      ! 
    9495      IF( nn_timing == 1 )  CALL timing_start('p4z_flx') 
    9596      ! 
    96       CALL wrk_alloc( jpi, jpj, zkgco2, zkgo2, zh2co3, zoflx ) 
     97      CALL wrk_alloc( jpi, jpj, zkgco2, zkgo2, zh2co3, zoflx, zpco2atm ) 
    9798      ! 
    9899 
     
    135136 
    136137               ! CALCULATE [ALK]([CO3--], [HCO3-]) 
    137                zalk  = zalka - (  akw3(ji,jj,1) / zph - zph + zbot / ( 1.+ zph / akb3(ji,jj,1) )  ) 
     138               zalk  = zalka - (  akw3(ji,jj,1) / zph - zph / aphscale(ji,jj,1)    & 
     139               &       + zbot / ( 1.+ zph / akb3(ji,jj,1) )  ) 
    138140 
    139141               ! CALCULATE [H+] AND [H2CO3] 
     
    162164            ztc2 = ztc * ztc 
    163165            ztc3 = ztc * ztc2  
     166            ztc4 = ztc2 * ztc2  
    164167            ! Compute the schmidt Number both O2 and CO2 
    165             zsch_co2 = 2073.1 - 125.62 * ztc + 3.6276 * ztc2 - 0.043126 * ztc3 
    166             zsch_o2  = 1953.4 - 128.0  * ztc + 3.9918 * ztc2 - 0.050091 * ztc3 
     168            zsch_co2 = 2116.8 - 136.25 * ztc + 4.7353 * ztc2 - 0.092307 * ztc3 + 0.0007555 * ztc4 
     169            zsch_o2  = 1920.4 - 135.6  * ztc + 5.2122 * ztc2 - 0.109390 * ztc3 + 0.0009377 * ztc4 
    167170            !  wind speed  
    168171            zws  = wndm(ji,jj) * wndm(ji,jj) 
    169172            ! Compute the piston velocity for O2 and CO2 
    170             zkgwan = 0.3 * zws  + 2.5 * ( 0.5246 + 0.016256 * ztc + 0.00049946  * ztc2 ) 
     173            zkgwan = 0.251 * zws 
    171174            zkgwan = zkgwan * xconv * ( 1.- fr_i(ji,jj) ) * tmask(ji,jj,1) 
    172175# if defined key_degrad 
     
    181184      DO jj = 1, jpj 
    182185         DO ji = 1, jpi 
     186            ztkel  = tsn(ji,jj,1,jp_tem) + 273.15 
     187            zsal   = tsn(ji,jj,1,jp_sal) + ( 1.- tmask(ji,jj,1) ) * 35. 
     188            zvapsw = EXP(24.4543 - 67.4509*(100.0/ztkel) - 4.8489*LOG(ztkel/100) - 0.000544*zsal) 
     189            zpco2atm(ji,jj) = satmco2(ji,jj) * ( patm(ji,jj) - zvapsw ) 
     190            zxc2 = (1.0 - zpco2atm(ji,jj) * 1E-6 )**2 
     191            zfugcoeff = EXP(patm(ji,jj) * (chemc(ji,jj,2) + 2.0 * zxc2 * chemc(ji,jj,3) )   & 
     192            &           / (82.05736 * ztkel)) 
     193            zfco2 = zpco2atm(ji,jj) * zfugcoeff 
     194 
    183195            ! Compute CO2 flux for the sea and air 
    184             zfld = satmco2(ji,jj) * patm(ji,jj) * tmask(ji,jj,1) * chemc(ji,jj,1) * zkgco2(ji,jj)   ! (mol/L) * (m/s) 
    185             zflu = zh2co3(ji,jj) * tmask(ji,jj,1) * zkgco2(ji,jj)                                   ! (mol/L) (m/s) ? 
     196            zfld = zfco2 * chemc(ji,jj,1) * zkgco2(ji,jj)  ! (mol/L) * (m/s) 
     197            zflu = zh2co3(ji,jj) * zkgco2(ji,jj)                                   ! (mol/L) (m/s) ? 
    186198            oce_co2(ji,jj) = ( zfld - zflu ) * rfact2 * e1e2t(ji,jj) * tmask(ji,jj,1) * 1000. 
    187199            ! compute the trend 
    188             tra(ji,jj,1,jpdic) = tra(ji,jj,1,jpdic) + ( zfld - zflu ) * rfact2 / fse3t(ji,jj,1) 
     200            tra(ji,jj,1,jpdic) = tra(ji,jj,1,jpdic) + ( zfld - zflu ) * rfact2 / fse3t(ji,jj,1) * tmask(ji,jj,1) 
    189201 
    190202            ! Compute O2 flux  
    191             zfld16 = atcox * patm(ji,jj) * chemc(ji,jj,2) * tmask(ji,jj,1) * zkgo2(ji,jj)          ! (mol/L) * (m/s) 
    192             zflu16 = trb(ji,jj,1,jpoxy) * tmask(ji,jj,1) * zkgo2(ji,jj) 
    193             zoflx(ji,jj) = zfld16 - zflu16 
     203            zfld16 = patm(ji,jj) * chemo2(ji,jj,1) * zkgo2(ji,jj)          ! (mol/L) * (m/s) 
     204            zflu16 = trb(ji,jj,1,jpoxy) * zkgo2(ji,jj) 
     205            zoflx(ji,jj) = ( zfld16 - zflu16 ) * tmask(ji,jj,1) 
    194206            tra(ji,jj,1,jpoxy) = tra(ji,jj,1,jpoxy) + zoflx(ji,jj) * rfact2 / fse3t(ji,jj,1) 
    195207         END DO 
     
    222234         ENDIF 
    223235         IF( iom_use( "Dpco2" ) ) THEN 
    224            zw2d(:,:) = ( satmco2(:,:) * patm(:,:) - zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) 
     236           zw2d(:,:) = ( zpco2atm(:,:) - zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) 
    225237           CALL iom_put( "Dpco2" ,  zw2d ) 
    226238         ENDIF 
    227239         IF( iom_use( "Dpo2" ) )  THEN 
    228            zw2d(:,:) = ( atcox * patm(:,:) - trb(:,:,1,jpoxy) / ( chemc(:,:,2) + rtrn ) ) * tmask(:,:,1) 
     240           zw2d(:,:) = ( atcox * patm(:,:) - atcox * trn(:,:,1,jpoxy) / ( chemo2(:,:,1) + rtrn ) ) * tmask(:,:,1) 
    229241           CALL iom_put( "Dpo2"  , zw2d ) 
    230242         ENDIF 
     
    238250            trc2d(:,:,jp_pcs0_2d + 1) = zoflx(:,:) * 1000 * tmask(:,:,1)  
    239251            trc2d(:,:,jp_pcs0_2d + 2) = zkgco2(:,:) * tmask(:,:,1)  
    240             trc2d(:,:,jp_pcs0_2d + 3) = ( satmco2(:,:) * patm(:,:) - zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1)  
    241          ENDIF 
    242       ENDIF 
    243       ! 
    244       CALL wrk_dealloc( jpi, jpj, zkgco2, zkgo2, zh2co3, zoflx ) 
     252            trc2d(:,:,jp_pcs0_2d + 3) = ( zpco2atm(:,:) - zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) 
     253         ENDIF 
     254      ENDIF 
     255      ! 
     256      CALL wrk_dealloc( jpi, jpj, zkgco2, zkgo2, zh2co3, zoflx, zpco2atm ) 
    245257      ! 
    246258      IF( nn_timing == 1 )  CALL timing_stop('p4z_flx') 
     
    400412 
    401413   !!====================================================================== 
    402 END MODULE  p4zflx 
     414END MODULE p4zflx 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zint.F90

    r5602 r7256  
    8181 
    8282   !!====================================================================== 
    83 END MODULE  p4zint 
     83END MODULE p4zint 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zlim.F90

    r5602 r7256  
    4444   REAL(wp), PUBLIC ::  xkdoc       !:  2nd half-sat. of DOC remineralization   
    4545   REAL(wp), PUBLIC ::  concbfe     !:  Fe half saturation for bacteria  
     46   REAL(wp), PUBLIC ::  oxymin      !:  half saturation constant for anoxia 
    4647   REAL(wp), PUBLIC ::  qnfelim     !:  optimal Fe quota for nanophyto 
    4748   REAL(wp), PUBLIC ::  qdfelim     !:  optimal Fe quota for diatoms 
     
    121122               zlim1    = xnanono3(ji,jj,jk) + xnanonh4(ji,jj,jk) 
    122123               zlim2    = trb(ji,jj,jk,jppo4) / ( trb(ji,jj,jk,jppo4) + concbnh4 ) 
    123                zlim3    = trb(ji,jj,jk,jpfer) / ( concbfe + trb(ji,jj,jk,jpfer) ) 
     124               zlim3    = biron(ji,jj,jk)     / ( concbfe + biron(ji,jj,jk) ) 
    124125               zlim4    = trb(ji,jj,jk,jpdoc) / ( xkdoc   + trb(ji,jj,jk,jpdoc) ) 
    125126               xlimbacl(ji,jj,jk) = MIN( zlim1, zlim2, zlim3 ) 
     
    187188      END DO 
    188189      ! 
     190      DO jk = 1, jpkm1 
     191         DO jj = 1, jpj 
     192            DO ji = 1, jpi 
     193               ! denitrification factor computed from O2 levels 
     194               nitrfac(ji,jj,jk) = MAX(  0.e0, 0.4 * ( 6.e-6  - trb(ji,jj,jk,jpoxy) )    & 
     195                  &                                / ( oxymin + trb(ji,jj,jk,jpoxy) )  ) 
     196               nitrfac(ji,jj,jk) = MIN( 1., nitrfac(ji,jj,jk) ) 
     197            END DO 
     198         END DO 
     199      END DO 
    189200      ! 
    190201      IF( lk_iomput .AND. knt == nrdttrc ) THEN        ! save output diagnostics 
     
    216227      NAMELIST/nampislim/ concnno3, concdno3, concnnh4, concdnh4, concnfer, concdfer, concbfe,   & 
    217228         &                concbno3, concbnh4, xsizedia, xsizephy, xsizern, xsizerd,          &  
    218          &                xksi1, xksi2, xkdoc, qnfelim, qdfelim, caco3r 
     229         &                xksi1, xksi2, xkdoc, qnfelim, qdfelim, caco3r, oxymin 
    219230      INTEGER :: ios                 ! Local integer output status for namelist read 
    220231 
     
    249260         WRITE(numout,*) '    Minimum size criteria for nanophyto      xsizephy  = ', xsizephy 
    250261         WRITE(numout,*) '    Fe half saturation for bacteria          concbfe   = ', concbfe 
     262         WRITE(numout,*) '    halk saturation constant for anoxia       oxymin   =' , oxymin 
    251263         WRITE(numout,*) '    optimal Fe quota for nano.               qnfelim   = ', qnfelim 
    252264         WRITE(numout,*) '    Optimal Fe quota for diatoms             qdfelim   = ', qdfelim 
    253265      ENDIF 
    254  
     266      ! 
     267      nitrfac (:,:,:) = 0._wp 
     268      ! 
    255269   END SUBROUTINE p4z_lim_init 
    256270 
     
    265279 
    266280   !!====================================================================== 
    267 END MODULE  p4zlim 
     281END MODULE p4zlim 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zlys.F90

    r5602 r7256  
    6565      REAL(wp) ::   zomegaca, zexcess, zexcess0 
    6666      CHARACTER (len=25) :: charout 
    67       REAL(wp), POINTER, DIMENSION(:,:,:) :: zco3, zcaldiss    
     67      REAL(wp), POINTER, DIMENSION(:,:,:) :: zco3, zco3sat, zcaldiss    
    6868      !!--------------------------------------------------------------------- 
    6969      ! 
    7070      IF( nn_timing == 1 )  CALL timing_start('p4z_lys') 
    7171      ! 
    72       CALL wrk_alloc( jpi, jpj, jpk, zco3, zcaldiss ) 
     72      CALL wrk_alloc( jpi, jpj, jpk, zco3, zco3sat, zcaldiss ) 
    7373      ! 
    7474      zco3    (:,:,:) = 0. 
     
    9191                  zalka = trb(ji,jj,jk,jptal) / zfact 
    9292                  ! CALCULATE [ALK]([CO3--], [HCO3-]) 
    93                   zalk  = zalka - ( akw3(ji,jj,jk) / zph - zph + borat(ji,jj,jk) / ( 1. + zph / akb3(ji,jj,jk) ) ) 
     93                  zalk  = zalka - ( akw3(ji,jj,jk) / zph - zph / ( aphscale(ji,jj,jk) + rtrn )  & 
     94                  &       + borat(ji,jj,jk) / ( 1. + zph / akb3(ji,jj,jk) ) ) 
    9495                  ! CALCULATE [H+] and [CO3--] 
    9596                  zaldi = zdic - zalk 
     
    119120               zcalcon  = calcon * ( tsn(ji,jj,jk,jp_sal) / 35._wp ) 
    120121               zfact    = rhop(ji,jj,jk) / 1000._wp 
    121                zomegaca = ( zcalcon * zco3(ji,jj,jk) * zfact ) / aksp(ji,jj,jk)  
     122               zomegaca = ( zcalcon * zco3(ji,jj,jk) ) / ( aksp(ji,jj,jk) * zfact + rtrn ) 
     123               zco3sat(ji,jj,jk) = aksp(ji,jj,jk) * zfact / ( zcalcon + rtrn ) 
    122124 
    123125               ! SET DEGREE OF UNDER-/SUPERSATURATION 
     
    148150      IF( lk_iomput .AND. knt == nrdttrc ) THEN 
    149151         IF( iom_use( "PH"     ) ) CALL iom_put( "PH"    , -1. * LOG10( hi(:,:,:) )          * tmask(:,:,:) ) 
    150          IF( iom_use( "CO3"    ) ) CALL iom_put( "CO3"   , zco3(:,:,:) * 1.e+3               * tmask(:,:,:) ) 
    151          IF( iom_use( "CO3sat" ) ) CALL iom_put( "CO3sat", aksp(:,:,:) * 1.e+3 / calcon      * tmask(:,:,:) ) 
    152          IF( iom_use( "DCAL"   ) ) CALL iom_put( "DCAL"  , zcaldiss(:,:,:) * 1.e+3 * rfact2r   * tmask(:,:,:) ) 
     152         IF( iom_use( "CO3"    ) ) CALL iom_put( "CO3"   , zco3(:,:,:)    * 1.e+3            * tmask(:,:,:) ) 
     153         IF( iom_use( "CO3sat" ) ) CALL iom_put( "CO3sat", zco3sat(:,:,:) * 1.e+3            * tmask(:,:,:) ) 
     154         IF( iom_use( "DCAL"   ) ) CALL iom_put( "DCAL"  , zcaldiss(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:) ) 
    153155      ELSE 
    154          trc3d(:,:,:,jp_pcs0_3d    ) = -1. * LOG10( hi(:,:,:) ) * tmask(:,:,:) 
    155          trc3d(:,:,:,jp_pcs0_3d + 1) = zco3(:,:,:)              * tmask(:,:,:) 
    156          trc3d(:,:,:,jp_pcs0_3d + 2) = aksp(:,:,:) / calcon     * tmask(:,:,:) 
     156         IF( ln_diatrc ) THEN 
     157            trc3d(:,:,:,jp_pcs0_3d    ) = -1. * LOG10( hi(:,:,:) ) * tmask(:,:,:) 
     158            trc3d(:,:,:,jp_pcs0_3d + 1) = zco3(:,:,:)              * tmask(:,:,:) 
     159            trc3d(:,:,:,jp_pcs0_3d + 2) = zco3sat(:,:,:)           * tmask(:,:,:) 
     160         ENDIF 
    157161      ENDIF 
    158162      ! 
     
    163167      ENDIF 
    164168      ! 
    165       CALL wrk_dealloc( jpi, jpj, jpk, zco3, zcaldiss ) 
     169      CALL wrk_dealloc( jpi, jpj, jpk, zco3, zco3sat, zcaldiss ) 
    166170      ! 
    167171      IF( nn_timing == 1 )  CALL timing_stop('p4z_lys') 
     
    223227#endif  
    224228   !!====================================================================== 
    225 END MODULE  p4zlys 
     229END MODULE p4zlys 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zmeso.F90

    r5602 r7256  
    340340 
    341341   !!====================================================================== 
    342 END MODULE  p4zmeso 
     342END MODULE p4zmeso 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zmicro.F90

    r5602 r7256  
    273273 
    274274   !!====================================================================== 
    275 END MODULE  p4zmicro 
     275END MODULE p4zmicro 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zmort.F90

    r5602 r7256  
    277277 
    278278   !!====================================================================== 
    279 END MODULE  p4zmort 
     279END MODULE p4zmort 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zopt.F90

    r5602 r7256  
    7676      REAL(wp) ::   zchl 
    7777      REAL(wp) ::   zc0 , zc1 , zc2, zc3, z1_dep 
    78       REAL(wp), POINTER, DIMENSION(:,:  ) :: zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4, zqsr100 
     78      REAL(wp), POINTER, DIMENSION(:,:  ) :: zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 
     79      REAL(wp), POINTER, DIMENSION(:,:  ) :: zqsr100, zqsr_corr 
    7980      REAL(wp), POINTER, DIMENSION(:,:,:) :: zpar, ze0, ze1, ze2, ze3 
    8081      !!--------------------------------------------------------------------- 
     
    8384      ! 
    8485      ! Allocate temporary workspace 
    85       CALL wrk_alloc( jpi, jpj,      zqsr100, zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 ) 
     86      CALL wrk_alloc( jpi, jpj,      zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 ) 
     87      CALL wrk_alloc( jpi, jpj,      zqsr100, zqsr_corr ) 
    8688      CALL wrk_alloc( jpi, jpj, jpk, zpar, ze0, ze1, ze2, ze3 ) 
    8789 
     
    112114      !                                        !  -------------------------------------- 
    113115      IF( l_trcdm2dc ) THEN                     !  diurnal cycle 
    114          ! 1% of qsr to compute euphotic layer 
    115          zqsr100(:,:) = 0.01 * qsr_mean(:,:)     !  daily mean qsr 
    116          ! 
    117          CALL p4z_opt_par( kt, qsr_mean, ze1, ze2, ze3 )  
     116         ! 
     117         zqsr_corr(:,:) = qsr_mean(:,:) / ( 1. - fr_i(:,:) + rtrn ) 
     118         ! 
     119         CALL p4z_opt_par( kt, zqsr_corr, ze1, ze2, ze3, pqsr100 = zqsr100 )  
    118120         ! 
    119121         DO jk = 1, nksrp       
     
    123125         END DO 
    124126         ! 
    125          CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3 )  
     127         zqsr_corr(:,:) = qsr(:,:) / ( 1. - fr_i(:,:) + rtrn ) 
     128         ! 
     129         CALL p4z_opt_par( kt, zqsr_corr, ze1, ze2, ze3 )  
    126130         ! 
    127131         DO jk = 1, nksrp       
     
    130134         ! 
    131135      ELSE 
    132          ! 1% of qsr to compute euphotic layer 
    133          zqsr100(:,:) = 0.01 * qsr(:,:) 
    134          ! 
    135          CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3 )  
     136         ! 
     137         zqsr_corr(:,:) = qsr(:,:) / ( 1. - fr_i(:,:) + rtrn ) 
     138         ! 
     139         CALL p4z_opt_par( kt, zqsr_corr, ze1, ze2, ze3, pqsr100 = zqsr100 )  
    136140         ! 
    137141         DO jk = 1, nksrp       
     
    161165         DO jj = 1, jpj 
    162166           DO ji = 1, jpi 
    163               IF( etot_ndcy(ji,jj,jk) * tmask(ji,jj,jk) >= 0.43 * zqsr100(ji,jj) )  THEN 
     167              IF( etot_ndcy(ji,jj,jk) * tmask(ji,jj,jk) >= zqsr100(ji,jj) )  THEN 
    164168                 neln(ji,jj) = jk+1                    ! Euphotic level : 1rst T-level strictly below Euphotic layer 
    165169                 !                                     ! nb: ensure the compatibility with nmld_trc definition in trd_mld_trc_zint 
     
    226230      ENDIF 
    227231      ! 
    228       CALL wrk_dealloc( jpi, jpj,      zqsr100, zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 ) 
     232      CALL wrk_dealloc( jpi, jpj,      zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 ) 
     233      CALL wrk_dealloc( jpi, jpj,      zqsr100, zqsr_corr ) 
    229234      CALL wrk_dealloc( jpi, jpj, jpk, zpar,  ze0, ze1, ze2, ze3 ) 
    230235      ! 
     
    233238   END SUBROUTINE p4z_opt 
    234239 
    235    SUBROUTINE p4z_opt_par( kt, pqsr, pe1, pe2, pe3, pe0 )  
     240   SUBROUTINE p4z_opt_par( kt, pqsr, pe1, pe2, pe3, pe0, pqsr100 )  
    236241      !!---------------------------------------------------------------------- 
    237242      !!                  ***  routine p4z_opt_par  *** 
     
    242247      !!---------------------------------------------------------------------- 
    243248      !! * arguments 
    244       INTEGER, INTENT(in)                                       ::  kt            !   ocean time-step 
    245       REAL(wp), DIMENSION(jpi,jpj)    , INTENT(in)              ::  pqsr          !   shortwave 
    246       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout)           ::  pe1 , pe2 , pe3   !  PAR ( R-G-B) 
    247       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout), OPTIONAL ::  pe0   
     249      INTEGER, INTENT(in)                                        ::  kt            !   ocean time-step 
     250      REAL(wp), DIMENSION(jpi,jpj)    , INTENT(in)               ::  pqsr          !   shortwave 
     251      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout)            ::  pe1 , pe2 , pe3   !  PAR ( R-G-B) 
     252      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout), OPTIONAL  ::  pe0  
     253      REAL(wp), DIMENSION(jpi,jpj)    , INTENT(out)  , OPTIONAL  ::  pqsr100   
    248254      !! * local variables 
    249255      INTEGER    ::   ji, jj, jk     ! dummy loop indices 
     
    255261      ELSE                  ;  zqsr(:,:) = xparsw         * pqsr(:,:) 
    256262      ENDIF 
     263 
     264      !  Light at the euphotic depth  
     265      IF( PRESENT( pqsr100 ) )  pqsr100(:,:) = 0.01 * 3. * zqsr(:,:) 
    257266      ! 
    258267      IF( PRESENT( pe0 ) ) THEN     !  W-level 
     
    439448 
    440449   !!====================================================================== 
    441 END MODULE  p4zopt 
     450END MODULE p4zopt 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zprod.F90

    r5602 r7256  
    202202                      zdiattot    = ediat(ji,jj,jk) * zstrn(ji,jj) 
    203203                      ! 
    204                       zpislopead (ji,jj,jk) = pislope  * ( 1.+ zadap  * EXP( -znanotot ) ) 
    205                       zpislopead2(ji,jj,jk) = (pislope * zconctemp2 + pislope2 * zconctemp)  / ( trb(ji,jj,jk,jpdia) + rtrn ) 
    206  
    207                       zpislopen =  zpislopead(ji,jj,jk) * trb(ji,jj,jk,jpnch)                & 
    208                         &          / ( trb(ji,jj,jk,jpphy) * 12.                  + rtrn )   & 
    209                         &          / ( prmax(ji,jj,jk) * rday * xlimphy(ji,jj,jk) + rtrn ) 
    210  
    211                       zpislope2n = zpislopead2(ji,jj,jk) * trb(ji,jj,jk,jpdch)                & 
    212                         &          / ( trb(ji,jj,jk,jpdia) * 12.                  + rtrn )   & 
    213                         &          / ( prmax(ji,jj,jk) * rday * xlimdia(ji,jj,jk) + rtrn ) 
     204                      zpislopead (ji,jj,jk) = pislope  * ( 1.+ zadap  * EXP( -znanotot ) )           & 
     205                         &                   * trb(ji,jj,jk,jpnch) /( trb(ji,jj,jk,jpphy) * 12. + rtrn) 
     206                      zpislopead2(ji,jj,jk) = (pislope * zconctemp2 + pislope2 * zconctemp)  / ( trb(ji,jj,jk,jpdia) + rtrn )   & 
     207                         &                   * trb(ji,jj,jk,jpdch) /( trb(ji,jj,jk,jpdia) * 12. + rtrn) 
    214208 
    215209                      ! Computation of production function for Carbon 
    216210                      !  --------------------------------------------- 
     211                      zpislopen  =  zpislopead(ji,jj,jk)  / ( prmax(ji,jj,jk) * rday * xlimphy(ji,jj,jk) + rtrn ) 
     212                      zpislope2n =  zpislopead2(ji,jj,jk) / ( prmax(ji,jj,jk) * rday * xlimdia(ji,jj,jk) + rtrn ) 
    217213                      zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen  * znanotot ) ) 
    218214                      zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpislope2n * zdiattot ) ) 
     
    220216                      !  Computation of production function for Chlorophyll 
    221217                      !-------------------------------------------------- 
    222                       zprnch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislopen  * enano(ji,jj,jk) ) ) 
    223                       zprdch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislope2n * ediat(ji,jj,jk) ) ) 
     218                      zprnch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislopen  * znanotot ) ) 
     219                      zprdch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislope2n * zdiattot ) ) 
    224220                  ENDIF 
    225221               END DO 
     
    227223         END DO 
    228224      ENDIF 
    229  
    230  
     225       
    231226      !  Computation of a proxy of the N/C ratio 
    232227      !  --------------------------------------- 
     
    278273            zmxltst = MAX( 0.e0, hmld(ji,jj) - heup(ji,jj) ) 
    279274            zmxlday = zmxltst * zmxltst * r1_rday 
    280             zmixnano(ji,jj) = 1. - zmxlday / ( 2. + zmxlday ) 
    281             zmixdiat(ji,jj) = 1. - zmxlday / ( 4. + zmxlday ) 
     275            zmixnano(ji,jj) = 1. - zmxlday / ( 1. + zmxlday ) 
     276            zmixdiat(ji,jj) = 1. - zmxlday / ( 2. + zmxlday ) 
    282277         END DO 
    283278      END DO 
    284279  
    285       !  Mixed-layer effect on production                                                                                
     280      !  Mixed-layer effect on production  
     281      !  Sea-ice effect on production 
     282 
    286283      DO jk = 1, jpkm1 
    287284         DO jj = 1, jpj 
     
    291288                  zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * zmixdiat(ji,jj) 
    292289               ENDIF 
     290                  zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1. - fr_i(ji,jj) ) 
     291                  zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1. - fr_i(ji,jj) ) 
    293292            END DO 
    294293         END DO 
     
    330329      END DO 
    331330 
    332       IF( ln_newprod ) THEN 
    333 !CDIR NOVERRCHK 
    334          DO jk = 1, jpkm1 
    335 !CDIR NOVERRCHK 
    336             DO jj = 1, jpj 
    337 !CDIR NOVERRCHK 
    338                DO ji = 1, jpi 
    339                   IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN 
    340                      zprnch(ji,jj,jk) = zprnch(ji,jj,jk) * zmixnano(ji,jj) 
    341                      zprdch(ji,jj,jk) = zprdch(ji,jj,jk) * zmixdiat(ji,jj) 
    342                   ENDIF 
    343                   IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
    344                      !  production terms for nanophyto. ( chlorophyll ) 
    345                      znanotot = enano(ji,jj,jk) * zstrn(ji,jj) 
    346                      zprod    = rday * zprorca(ji,jj,jk) * zprnch(ji,jj,jk) * xlimphy(ji,jj,jk) 
    347                      zprochln(ji,jj,jk) = chlcmin * 12. * zprorca (ji,jj,jk) 
    348                      zprochln(ji,jj,jk) = zprochln(ji,jj,jk) + (chlcnm-chlcmin) * 12. * zprod / & 
    349                                         & (  zpislopead(ji,jj,jk) * znanotot +rtrn) 
    350                      !  production terms for diatomees ( chlorophyll ) 
    351                      zdiattot = ediat(ji,jj,jk) * zstrn(ji,jj) 
    352                      zprod = rday * zprorcad(ji,jj,jk) * zprdch(ji,jj,jk) * xlimdia(ji,jj,jk) 
    353                      zprochld(ji,jj,jk) = chlcmin * 12. * zprorcad(ji,jj,jk) 
    354                      zprochld(ji,jj,jk) = zprochld(ji,jj,jk) + (chlcdm-chlcmin) * 12. * zprod / & 
    355                                         & ( zpislopead2(ji,jj,jk) * zdiattot +rtrn ) 
    356                   ENDIF 
    357                END DO 
    358             END DO 
    359          END DO 
    360       ELSE 
    361 !CDIR NOVERRCHK 
    362          DO jk = 1, jpkm1 
    363 !CDIR NOVERRCHK 
    364             DO jj = 1, jpj 
    365 !CDIR NOVERRCHK 
    366                DO ji = 1, jpi 
    367                   IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
    368                      !  production terms for nanophyto. ( chlorophyll ) 
    369                      znanotot = enano(ji,jj,jk) 
    370                      zprod = rday * zprorca(ji,jj,jk) * zprnch(ji,jj,jk) * trb(ji,jj,jk,jpphy) * xlimphy(ji,jj,jk) 
    371                      zprochln(ji,jj,jk) = chlcmin * 12. * zprorca (ji,jj,jk) 
    372                      zprochln(ji,jj,jk) = zprochln(ji,jj,jk) + (chlcnm-chlcmin) * 144. * zprod            & 
    373                      &                    / ( zpislopead(ji,jj,jk) * trb(ji,jj,jk,jpnch) * znanotot +rtrn ) 
    374                      !  production terms for diatomees ( chlorophyll ) 
    375                      zdiattot = ediat(ji,jj,jk) 
    376                      zprod = rday * zprorcad(ji,jj,jk) * zprdch(ji,jj,jk) * trb(ji,jj,jk,jpdia) * xlimdia(ji,jj,jk) 
    377                      zprochld(ji,jj,jk) = chlcmin * 12. * zprorcad(ji,jj,jk) 
    378                      zprochld(ji,jj,jk) = zprochld(ji,jj,jk) + (chlcdm-chlcmin) * 144. * zprod             & 
    379                      &                    / ( zpislopead2(ji,jj,jk) * trb(ji,jj,jk,jpdch) * zdiattot +rtrn ) 
    380                   ENDIF 
    381                END DO 
    382             END DO 
    383          END DO 
    384       ENDIF 
     331!CDIR NOVERRCHK 
     332      DO jk = 1, jpkm1 
     333!CDIR NOVERRCHK 
     334         DO jj = 1, jpj 
     335!CDIR NOVERRCHK 
     336            DO ji = 1, jpi 
     337               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN 
     338                  zprnch(ji,jj,jk) = zprnch(ji,jj,jk) * zmixnano(ji,jj) 
     339                  zprdch(ji,jj,jk) = zprdch(ji,jj,jk) * zmixdiat(ji,jj) 
     340               ENDIF 
     341               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
     342                  !  production terms for nanophyto. ( chlorophyll ) 
     343                  znanotot = enano(ji,jj,jk) * zstrn(ji,jj) 
     344                  zprod    = rday * zprorca(ji,jj,jk) * zprnch(ji,jj,jk) * xlimphy(ji,jj,jk) 
     345                  zprochln(ji,jj,jk) = chlcmin * 12. * zprorca (ji,jj,jk) 
     346                  zprochln(ji,jj,jk) = zprochln(ji,jj,jk) + (chlcnm-chlcmin) * 12. * zprod / & 
     347                                     & (  zpislopead(ji,jj,jk) * znanotot +rtrn) 
     348                  !  production terms for diatomees ( chlorophyll ) 
     349                  zdiattot = ediat(ji,jj,jk) * zstrn(ji,jj) 
     350                  zprod = rday * zprorcad(ji,jj,jk) * zprdch(ji,jj,jk) * xlimdia(ji,jj,jk) 
     351                  zprochld(ji,jj,jk) = chlcmin * 12. * zprorcad(ji,jj,jk) 
     352                  zprochld(ji,jj,jk) = zprochld(ji,jj,jk) + (chlcdm-chlcmin) * 12. * zprod / & 
     353                                     & ( zpislopead2(ji,jj,jk) * zdiattot +rtrn ) 
     354               ENDIF 
     355            END DO 
     356         END DO 
     357      END DO 
    385358 
    386359      !   Update the arrays TRA which contain the biological sources and sinks 
     
    629602 
    630603   !!====================================================================== 
    631 END MODULE  p4zprod 
     604END MODULE p4zprod 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zrem.F90

    r5602 r7256  
    4444   REAL(wp), PUBLIC ::  xsiremlab  !: fast remineralisation rate of POC  
    4545   REAL(wp), PUBLIC ::  xsilab     !: fraction of labile biogenic silica  
    46    REAL(wp), PUBLIC ::  oxymin     !: halk saturation constant for anoxia  
    47  
    4846 
    4947   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   denitr     !: denitrification array 
     
    111109                  zdepprod(ji,jj,jk) = zdepmin**0.273 
    112110               ENDIF 
    113             END DO 
    114          END DO 
    115       END DO 
    116  
    117       DO jk = 1, jpkm1 
    118          DO jj = 1, jpj 
    119             DO ji = 1, jpi 
    120                ! denitrification factor computed from O2 levels 
    121                nitrfac(ji,jj,jk) = MAX(  0.e0, 0.4 * ( 6.e-6  - trb(ji,jj,jk,jpoxy) )    & 
    122                   &                                / ( oxymin + trb(ji,jj,jk,jpoxy) )  ) 
    123                nitrfac(ji,jj,jk) = MIN( 1., nitrfac(ji,jj,jk) ) 
    124111            END DO 
    125112         END DO 
     
    357344      !! 
    358345      !!---------------------------------------------------------------------- 
    359       NAMELIST/nampisrem/ xremik, xremip, nitrif, xsirem, xsiremlab, xsilab,   & 
    360       &                   oxymin 
     346      NAMELIST/nampisrem/ xremik, xremip, nitrif, xsirem, xsiremlab, xsilab 
    361347      INTEGER :: ios                 ! Local integer output status for namelist read 
    362348 
     
    380366         WRITE(numout,*) '    fraction of labile biogenic silica        xsilab    =', xsilab 
    381367         WRITE(numout,*) '    NH4 nitrification rate                    nitrif    =', nitrif 
    382          WRITE(numout,*) '    halk saturation constant for anoxia       oxymin    =', oxymin 
    383368      ENDIF 
    384369      ! 
    385       nitrfac (:,:,:) = 0._wp 
    386370      denitr  (:,:,:) = 0._wp 
    387371      denitnh4(:,:,:) = 0._wp 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsbc.F90

    r5602 r7256  
    159159      IF( ln_ndepo ) THEN 
    160160         IF( kt == nit000 .OR. ( kt /= nit000 .AND. ntimes_ndep > 1 ) ) THEN 
    161             CALL fld_read( kt, 1, sf_ndepo ) 
    162             DO jj = 1, jpj 
    163                DO ji = 1, jpi 
    164                   nitdep(ji,jj) = sf_ndepo(1)%fnow(ji,jj,1) / rno3 / ( 14E6 * ryyss * fse3t(ji,jj,1) + rtrn ) 
    165                END DO 
    166             END DO 
     161             zcoef = rno3 * 14E6 * ryyss 
     162             CALL fld_read( kt, 1, sf_ndepo ) 
     163             nitdep(:,:) = sf_ndepo(1)%fnow(:,:,1) / zcoef / fse3t(:,:,1)  
     164         ENDIF 
     165         IF( lk_vvl ) THEN 
     166           zcoef = rno3 * 14E6 * ryyss 
     167           nitdep(:,:) = sf_ndepo(1)%fnow(:,:,1) / zcoef / fse3t(:,:,1)  
    167168         ENDIF 
    168169      ENDIF 
     
    266267      IF( lk_offline ) THEN 
    267268        nk_rnf(:,:) = 1 
    268         h_rnf (:,:) = fsdept(:,:,1) 
     269        h_rnf (:,:) = e3t_0(:,:,1) 
    269270      ENDIF 
    270271 
     
    455456            DO jj = 1, jpj 
    456457               DO ji = 1, jpi 
    457                   zexpide   = MIN( 8.,( fsdept(ji,jj,jk) / 500. )**(-1.5) ) 
     458                  zexpide   = MIN( 8.,( gdept_0(ji,jj,jk) / 500. )**(-1.5) ) 
    458459                  zdenitide = -0.9543 + 0.7662 * LOG( zexpide ) - 0.235 * LOG( zexpide )**2 
    459460                  zcmask(ji,jj,jk) = zcmask(ji,jj,jk) * MIN( 1., EXP( zdenitide ) / 0.5 ) 
     
    465466         ironsed(:,:,jpk) = 0._wp 
    466467         DO jk = 1, jpkm1 
    467             ironsed(:,:,jk) = sedfeinput * zcmask(:,:,jk) / ( fse3t(:,:,jk) * rday ) 
     468            ironsed(:,:,jk) = sedfeinput * zcmask(:,:,jk) / ( e3t_0(:,:,jk) * rday ) 
    468469         END DO 
    469470         DEALLOCATE( zcmask) 
     
    483484         CALL iom_close( numhydro ) 
    484485         ! 
    485          hydrofe(:,:,:) = ( hydrofe(:,:,:) * hratio ) / ( cvol(:,:,:) * ryyss + rtrn ) / 1000._wp 
     486         DO jk = 1, jpk 
     487            hydrofe(:,:,jk) = ( hydrofe(:,:,jk) * hratio ) / ( e1e2t(:,:) * e3t_0(:,:,jk) * ryyss + rtrn ) / 1000._wp 
     488         ENDDO 
    486489         ! 
    487490      ENDIF 
     
    519522 
    520523   !!====================================================================== 
    521 END MODULE  p4zsbc 
     524END MODULE p4zsbc 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsed.F90

    r5602 r7256  
    7272      CHARACTER (len=25) :: charout 
    7373      REAL(wp), POINTER, DIMENSION(:,:  ) :: zpdep, zsidep, zwork1, zwork2, zwork3 
     74      REAL(wp), POINTER, DIMENSION(:,:)   :: zsedcal, zsedsi, zsedc 
    7475      REAL(wp), POINTER, DIMENSION(:,:  ) :: zdenit2d, zironice, zbureff 
    7576      REAL(wp), POINTER, DIMENSION(:,:  ) :: zwsbio3, zwsbio4, zwscal 
     
    8384      ! Allocate temporary workspace 
    8485      CALL wrk_alloc( jpi, jpj, zdenit2d, zwork1, zwork2, zwork3, zbureff ) 
     86      CALL wrk_alloc( jpi, jpj, zsedcal,  zsedsi, zsedc ) 
    8587      CALL wrk_alloc( jpi, jpj, zwsbio3, zwsbio4, zwscal ) 
    8688      CALL wrk_alloc( jpi, jpj, jpk, zsoufer ) 
     
    9193      zwork2  (:,:) = 0.e0 
    9294      zwork3  (:,:) = 0.e0 
     95      zsedsi   (:,:) = 0.e0 
     96      zsedcal  (:,:) = 0.e0 
     97      zsedc    (:,:) = 0.e0 
    9398 
    9499      ! Iron input/uptake due to sea ice : Crude parameterization based on Lancelot et al. 
     
    298303            tra(ji,jj,ikt,jptal) =  tra(ji,jj,ikt,jptal) + zcaloss * zrivalk * 2.0 
    299304            tra(ji,jj,ikt,jpdic) =  tra(ji,jj,ikt,jpdic) + zcaloss * zrivalk 
     305            zsedcal(ji,jj) = (1.0 - zrivalk) * zcaloss / zdep 
     306            zsedsi (ji,jj) = (1.0 - zrivsil) * zsiloss / zdep 
    300307#endif 
    301308         END DO 
     
    336343            tra(ji,jj,ikt,jptal) = tra(ji,jj,ikt,jptal) + rno3 * (zolimit + (1.+rdenit) * (zpdenit + zdenitt) ) 
    337344            tra(ji,jj,ikt,jpdic) = tra(ji,jj,ikt,jpdic) + zpdenit + zolimit + zdenitt 
    338             sdenit(ji,jj) = rdenit * zpdenit * fse3t(ji,jj,ikt) 
     345            sdenit(ji,jj) = rdenit * zpdenit / zdep 
     346            zsedc(ji,jj)   = (1. - zrivno3) * zwstpoc / zdep 
    339347#endif 
    340348         END DO 
     
    392400               CALL iom_put( "INTNFIX" , zwork1 )  
    393401            ENDIF 
     402            IF( iom_use("SedCal" ) ) CALL iom_put( "SedCal", zsedcal(:,:) * 1.e+3 ) 
     403            IF( iom_use("SedSi" ) )  CALL iom_put( "SedSi",  zsedsi (:,:) * 1.e+3 ) 
     404            IF( iom_use("SedC" ) )   CALL iom_put( "SedC",   zsedc  (:,:) * 1.e+3 ) 
     405            IF( iom_use("Sdenit" ) ) CALL iom_put( "Sdenit", sdenit (:,:) * 1.e+3 * rno3 ) 
    394406         ENDIF 
    395407      ELSE 
     
    405417      ! 
    406418      CALL wrk_dealloc( jpi, jpj, zdenit2d, zwork1, zwork2, zwork3, zbureff ) 
     419      CALL wrk_dealloc( jpi, jpj, zsedcal , zsedsi, zsedc ) 
    407420      CALL wrk_dealloc( jpi, jpj, zwsbio3, zwsbio4, zwscal ) 
    408421      CALL wrk_dealloc( jpi, jpj, jpk, zsoufer ) 
     
    436449 
    437450   !!====================================================================== 
    438 END MODULE  p4zsed 
     451END MODULE p4zsed 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsink.F90

    r5602 r7256  
    913913 
    914914   !!====================================================================== 
    915 END MODULE  p4zsink 
     915END MODULE p4zsink 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90

    r5602 r7256  
    3838 
    3939   REAL(wp) :: alkbudget, no3budget, silbudget, ferbudget, po4budget 
    40    REAL(wp) :: xfact1, xfact2 
     40   REAL(wp) :: xfact1, xfact2, xfact3 
    4141   INTEGER ::  numco2, numnut, numnit  !: logical unit for co2 budget 
    4242 
     
    133133         ! 
    134134         CALL p4z_bio( kt, jnt )   ! Biology 
    135          CALL p4z_sed( kt, jnt )   ! Sedimentation 
    136135         CALL p4z_lys( kt, jnt )   ! Compute CaCO3 saturation 
     136         CALL p4z_sed( kt, jnt )   ! Surface and Bottom boundary conditions 
    137137         CALL p4z_flx( kt, jnt )   ! Compute surface fluxes 
    138138         ! 
     
    474474      !!--------------------------------------------------------------------- 
    475475      ! 
    476       INTEGER , INTENT( in ) ::   kt      ! ocean time-step index       
    477       REAL(wp)               ::  zfact        
    478       REAL(wp) ::  zrdenittot, zsdenittot, znitrpottot 
     476      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index       
     477      REAL(wp)             ::  zrdenittot, zsdenittot, znitrpottot 
    479478      CHARACTER(LEN=100)   ::   cltxt 
    480479      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zvol 
     
    492491            xfact1 = rfact2r * 12. / 1.e15 * ryyss    ! conversion molC/kt --> PgC/yr 
    493492            xfact2 = 1.e+3 * rno3 * 14. / 1.e12 * ryyss   ! conversion molC/l/s ----> TgN/m3/yr 
     493            xfact3 = 1.e+3 * rfact2r * rno3   ! conversion molC/l/kt ----> molN/m3/s 
    494494            cltxt='time-step   Alkalinity        Nitrate        Phosphorus         Silicate           Iron' 
    495495            IF( lwp ) WRITE(numnut,*)  TRIM(cltxt) 
     
    574574      IF( iom_use( "Sdenit" ) .OR.  ( ln_check_mass .AND. kt == nitend )  ) THEN 
    575575         zsdenittot   = glob_sum ( sdenit(:,:) * e1e2t(:,:) ) 
    576          CALL iom_put( "Sdenit", sdenit(:,:) * zfact * tmask(:,:,1) )  ! Nitrate reduction in the sediments 
     576         CALL iom_put( "Sdenit", sdenit(:,:) * xfact3 * tmask(:,:,1) )  ! Nitrate reduction in the sediments 
    577577      ENDIF 
    578578 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/sms_pisces.F90

    r5602 r7256  
    101101   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hi         !: ??? 
    102102   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   excess     !: ??? 
     103   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   aphscale   !:  
     104 
    103105 
    104106   !!* Temperature dependancy of SMS terms 
     
    154156         &      ak23(jpi,jpj,jpk)    , aksp  (jpi,jpj,jpk) ,       & 
    155157         &      akw3(jpi,jpj,jpk)    , borat (jpi,jpj,jpk) ,       & 
    156          &      hi  (jpi,jpj,jpk)    , excess(jpi,jpj,jpk) ,     STAT=ierr(4) ) 
     158         &      hi  (jpi,jpj,jpk)    , excess(jpi,jpj,jpk) ,       & 
     159         &      aphscale(jpi,jpj,jpk),                           STAT=ierr(4) ) 
    157160         ! 
    158161      !* Temperature dependancy of SMS terms 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/trcice_pisces.F90

    r5602 r7256  
    2929CONTAINS 
    3030 
     31 
    3132   SUBROUTINE trc_ice_ini_pisces 
    3233      !!---------------------------------------------------------------------- 
    33       !!                   ***  ROUTINE trc_ice_ini_pisces *** 
     34      !!                   ***  ROUTINE trc_ini_pisces *** 
     35      !! 
     36      !! ** Purpose :   Initialisation of the PISCES biochemical model 
     37      !!---------------------------------------------------------------------- 
     38 
     39      IF( lk_p4z ) THEN  ;   CALL p4z_ice_ini   !  PISCES 
     40      ELSE               ;   CALL p2z_ice_ini   !  LOBSTER 
     41      ENDIF 
     42 
     43   END SUBROUTINE trc_ice_ini_pisces 
     44 
     45 
     46   SUBROUTINE p4z_ice_ini 
     47 
     48#if defined key_pisces  
     49      !!---------------------------------------------------------------------- 
     50      !!                   ***  ROUTINE p4z_ice_ini *** 
    3451      !! 
    3552      !! ** Purpose :   PISCES fake sea ice model setting 
     
    5875 
    5976                                        !--- Dummy variables 
    60       REAL(wp), DIMENSION(jptra,2) & 
    61                ::  zratio            ! effective ice-ocean tracer cc ratio 
     77      REAL(wp), DIMENSION(jp_pisces,2)  :: zratio  ! effective ice-ocean tracer cc ratio 
     78      REAL(wp), DIMENSION(jp_pisces,4)  :: zpisc   ! prescribes concentration  
     79      !                                            !  1:global, 2:Arctic, 3:Antarctic, 4:Baltic 
     80 
    6281      REAL(wp), DIMENSION(2) :: zrs  ! ice-ocean salinity ratio, 1 - global, 2- Baltic 
    6382      REAL(wp) :: zsice_bal          ! prescribed ice salinity in the Baltic 
     
    8099      ! fluxes 
    81100 
    82       !--- Global case  
    83       IF ( cn_trc_o(jpdic) == 'GL ' ) trc_o(:,:,jpdic) =  1.99e-3_wp  
    84       IF ( cn_trc_o(jpdoc) == 'GL ' ) trc_o(:,:,jpdoc) =  2.04e-5_wp  
    85       IF ( cn_trc_o(jptal) == 'GL ' ) trc_o(:,:,jptal) =  2.31e-3_wp  
    86       IF ( cn_trc_o(jpoxy) == 'GL ' ) trc_o(:,:,jpoxy) =  2.47e-4_wp 
    87       IF ( cn_trc_o(jpcal) == 'GL ' ) trc_o(:,:,jpcal) =  1.04e-8_wp 
    88       IF ( cn_trc_o(jppo4) == 'GL ' ) trc_o(:,:,jppo4) =  5.77e-7_wp / po4r  
    89       IF ( cn_trc_o(jppoc) == 'GL ' ) trc_o(:,:,jppoc) =  1.27e-6_wp   
     101      !--- Global values 
     102      zpisc(jpdic,1) =  1.99e-3_wp  
     103      zpisc(jpdoc,1) =  2.04e-5_wp  
     104      zpisc(jptal,1) =  2.31e-3_wp  
     105      zpisc(jpoxy,1) =  2.47e-4_wp 
     106      zpisc(jpcal,1) =  1.04e-8_wp 
     107      zpisc(jppo4,1) =  5.77e-7_wp / po4r  
     108      zpisc(jppoc,1) =  1.27e-6_wp   
    90109#  if ! defined key_kriest 
    91       IF ( cn_trc_o(jpgoc) == 'GL ' ) trc_o(:,:,jpgoc) =  5.23e-8_wp   
    92       IF ( cn_trc_o(jpbfe) == 'GL ' ) trc_o(:,:,jpbfe) =  9.84e-13_wp  
     110      zpisc(jpgoc,1) =  5.23e-8_wp   
     111      zpisc(jpbfe,1) =  9.84e-13_wp  
    93112#  else 
    94       IF ( cn_trc_o(jpnum) == 'GL ' ) trc_o(:,:,jpnum) = 0. ! could not get this value since did not use it 
     113      zpisc(jpnum,1) = 0. ! could not get this value since did not use it 
    95114#  endif 
    96       IF ( cn_trc_o(jpsil) == 'GL ' ) trc_o(:,:,jpsil) =  7.36e-6_wp   
    97       IF ( cn_trc_o(jpdsi) == 'GL ' ) trc_o(:,:,jpdsi) =  1.07e-7_wp  
    98       IF ( cn_trc_o(jpgsi) == 'GL ' ) trc_o(:,:,jpgsi) =  1.53e-8_wp 
    99       IF ( cn_trc_o(jpphy) == 'GL ' ) trc_o(:,:,jpphy) =  9.57e-8_wp 
    100       IF ( cn_trc_o(jpdia) == 'GL ' ) trc_o(:,:,jpdia) =  4.24e-7_wp 
    101       IF ( cn_trc_o(jpzoo) == 'GL ' ) trc_o(:,:,jpzoo) =  6.07e-7_wp 
    102       IF ( cn_trc_o(jpmes) == 'GL ' ) trc_o(:,:,jpmes) =  3.44e-7_wp 
    103       IF ( cn_trc_o(jpfer) == 'GL ' ) trc_o(:,:,jpfer) =  4.06e-10_wp 
    104       IF ( cn_trc_o(jpsfe) == 'GL ' ) trc_o(:,:,jpsfe) =  2.51e-11_wp 
    105       IF ( cn_trc_o(jpdfe) == 'GL ' ) trc_o(:,:,jpdfe) =  6.57e-12_wp 
    106       IF ( cn_trc_o(jpnfe) == 'GL ' ) trc_o(:,:,jpnfe) =  1.76e-11_wp 
    107       IF ( cn_trc_o(jpnch) == 'GL ' ) trc_o(:,:,jpnch) =  1.67e-7_wp 
    108       IF ( cn_trc_o(jpdch) == 'GL ' ) trc_o(:,:,jpdch) =  1.02e-7_wp 
    109       IF ( cn_trc_o(jpno3) == 'GL ' ) trc_o(:,:,jpno3) =  5.79e-6_wp / rno3  
    110       IF ( cn_trc_o(jpnh4) == 'GL ' ) trc_o(:,:,jpnh4) =  3.22e-7_wp / rno3 
     115      zpisc(jpsil,1) =  7.36e-6_wp   
     116      zpisc(jpdsi,1) =  1.07e-7_wp  
     117      zpisc(jpgsi,1) =  1.53e-8_wp 
     118      zpisc(jpphy,1) =  9.57e-8_wp 
     119      zpisc(jpdia,1) =  4.24e-7_wp 
     120      zpisc(jpzoo,1) =  6.07e-7_wp 
     121      zpisc(jpmes,1) =  3.44e-7_wp 
     122      zpisc(jpfer,1) =  4.06e-10_wp 
     123      zpisc(jpsfe,1) =  2.51e-11_wp 
     124      zpisc(jpdfe,1) =  6.57e-12_wp 
     125      zpisc(jpnfe,1) =  1.76e-11_wp 
     126      zpisc(jpnch,1) =  1.67e-7_wp 
     127      zpisc(jpdch,1) =  1.02e-7_wp 
     128      zpisc(jpno3,1) =  5.79e-6_wp / rno3  
     129      zpisc(jpnh4,1) =  3.22e-7_wp / rno3 
    111130 
    112131      !--- Arctic specificities (dissolved inorganic & DOM) 
    113       IF ( cn_trc_o(jpdic) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpdic) =  1.98e-3_wp  ; END WHERE ; ENDIF 
    114       IF ( cn_trc_o(jpdoc) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpdoc) =  6.00e-6_wp  ; END WHERE ; ENDIF 
    115       IF ( cn_trc_o(jptal) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jptal) =  2.13e-3_wp  ; END WHERE ; ENDIF 
    116       IF ( cn_trc_o(jpoxy) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpoxy) =  3.65e-4_wp  ; END WHERE ; ENDIF 
    117       IF ( cn_trc_o(jpcal) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpcal) =  1.50e-9_wp  ; END WHERE ; ENDIF 
    118       IF ( cn_trc_o(jppo4) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jppo4) =  4.09e-7_wp / po4r ; END WHERE ; ENDIF 
    119       IF ( cn_trc_o(jppoc) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jppoc) =  4.05e-7_wp  ; END WHERE ; ENDIF 
     132      zpisc(jpdic,2) =  1.98e-3_wp  
     133      zpisc(jpdoc,2) =  6.00e-6_wp  
     134      zpisc(jptal,2) =  2.13e-3_wp  
     135      zpisc(jpoxy,2) =  3.65e-4_wp   
     136      zpisc(jpcal,2) =  1.50e-9_wp   
     137      zpisc(jppo4,2) =  4.09e-7_wp / po4r  
     138      zpisc(jppoc,2) =  4.05e-7_wp   
    120139#  if ! defined key_kriest 
    121       IF ( cn_trc_o(jpgoc) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpgoc) =  2.84e-8_wp  ; END WHERE ; ENDIF 
    122       IF ( cn_trc_o(jpbfe) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpbfe) =  7.03e-13_wp ; END WHERE ; ENDIF 
     140      zpisc(jpgoc,2) =  2.84e-8_wp   
     141      zpisc(jpbfe,2) =  7.03e-13_wp  
    123142#  else 
    124       IF ( cn_trc_o(jpnum) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpnum) =  0.00e-00_wp ; END WHERE ; ENDIF 
     143      zpisc(jpnum,2) =  0.00e-00_wp  
    125144#  endif 
    126       IF ( cn_trc_o(jpsil) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpsil) =  6.87e-6_wp  ; END WHERE ; ENDIF 
    127       IF ( cn_trc_o(jpdsi) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpdsi) =  1.73e-7_wp  ; END WHERE ; ENDIF 
    128       IF ( cn_trc_o(jpgsi) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpgsi) =  7.93e-9_wp  ; END WHERE ; ENDIF 
    129       IF ( cn_trc_o(jpphy) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpphy) =  5.25e-7_wp  ; END WHERE ; ENDIF 
    130       IF ( cn_trc_o(jpdia) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpdia) =  7.75e-7_wp  ; END WHERE ; ENDIF 
    131       IF ( cn_trc_o(jpzoo) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpzoo) =  3.34e-7_wp  ; END WHERE ; ENDIF 
    132       IF ( cn_trc_o(jpmes) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpmes) =  2.49e-7_wp  ; END WHERE ; ENDIF 
    133       IF ( cn_trc_o(jpfer) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpfer) =  1.43e-9_wp  ; END WHERE ; ENDIF 
    134       IF ( cn_trc_o(jpsfe) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpsfe) =  2.21e-11_wp ; END WHERE ; ENDIF 
    135       IF ( cn_trc_o(jpdfe) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpdfe) =  2.04e-11_wp ; END WHERE ; ENDIF 
    136       IF ( cn_trc_o(jpnfe) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpnfe) =  1.75e-11_wp ; END WHERE ; ENDIF 
    137       IF ( cn_trc_o(jpnch) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpnch) =  1.46e-07_wp ; END WHERE ; ENDIF 
    138       IF ( cn_trc_o(jpdch) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpdch) =  2.36e-07_wp ; END WHERE ; ENDIF 
    139       IF ( cn_trc_o(jpno3) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpno3) =  3.51e-06_wp / rno3 ; END WHERE ; ENDIF 
    140       IF ( cn_trc_o(jpnh4) == 'AA ' ) THEN ; WHERE( gphit(:,:) >= 00._wp ) ; trc_o(:,:,jpnh4) =  6.15e-08_wp / rno3 ; END WHERE ; ENDIF 
     145      zpisc(jpsil,2) =  6.87e-6_wp   
     146      zpisc(jpdsi,2) =  1.73e-7_wp  
     147      zpisc(jpgsi,2) =  7.93e-9_wp 
     148      zpisc(jpphy,2) =  5.25e-7_wp   
     149      zpisc(jpdia,2) =  7.75e-7_wp  
     150      zpisc(jpzoo,2) =  3.34e-7_wp 
     151      zpisc(jpmes,2) =  2.49e-7_wp   
     152      zpisc(jpfer,2) =  1.43e-9_wp  
     153      zpisc(jpsfe,2) =  2.21e-11_wp  
     154      zpisc(jpdfe,2) =  2.04e-11_wp  
     155      zpisc(jpnfe,2) =  1.75e-11_wp  
     156      zpisc(jpnch,2) =  1.46e-07_wp  
     157      zpisc(jpdch,2) =  2.36e-07_wp  
     158      zpisc(jpno3,2) =  3.51e-06_wp / rno3  
     159      zpisc(jpnh4,2) =  6.15e-08_wp / rno3  
    141160 
    142161      !--- Antarctic specificities (dissolved inorganic & DOM) 
    143       IF ( cn_trc_o(jpdic) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpdic) =  2.20e-3_wp  ; END WHERE ; ENDIF 
    144       IF ( cn_trc_o(jpdoc) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpdoc) =  7.02e-6_wp  ; END WHERE ; ENDIF 
    145       IF ( cn_trc_o(jptal) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jptal) =  2.37e-3_wp  ; END WHERE ; ENDIF 
    146       IF ( cn_trc_o(jpoxy) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpoxy) =  3.42e-4_wp  ; END WHERE ; ENDIF 
    147       IF ( cn_trc_o(jpcal) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpcal) =  3.17e-9_wp  ; END WHERE ; ENDIF 
    148       IF ( cn_trc_o(jppo4) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jppo4) =  1.88e-6_wp / po4r  ; END WHERE ; ENDIF 
    149       IF ( cn_trc_o(jppoc) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jppoc) =  1.13e-6_wp  ; END WHERE ; ENDIF 
     162      zpisc(jpdic,3) =  2.20e-3_wp   
     163      zpisc(jpdoc,3) =  7.02e-6_wp   
     164      zpisc(jptal,3) =  2.37e-3_wp   
     165      zpisc(jpoxy,3) =  3.42e-4_wp   
     166      zpisc(jpcal,3) =  3.17e-9_wp   
     167      zpisc(jppo4,3) =  1.88e-6_wp / po4r   
     168      zpisc(jppoc,3) =  1.13e-6_wp   
    150169#  if ! defined key_kriest 
    151       IF ( cn_trc_o(jpgoc) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpgoc) =  2.89e-8_wp  ; END WHERE ; ENDIF 
    152       IF ( cn_trc_o(jpbfe) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpbfe) =  5.63e-13_wp ; END WHERE ; ENDIF 
     170      zpisc(jpgoc,3) =  2.89e-8_wp   
     171      zpisc(jpbfe,3) =  5.63e-13_wp  
    153172#  else 
    154       IF ( cn_trc_o(jpnum) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpnum) =  0.00e-00_wp ; END WHERE ; ENDIF 
     173      zpisc(jpnum,3) =  0.00e-00_wp  
    155174#  endif 
    156       IF ( cn_trc_o(jpsil) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpsil) =  4.96e-5_wp  ; END WHERE ; ENDIF 
    157       IF ( cn_trc_o(jpdsi) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpdsi) =  5.63e-7_wp  ; END WHERE ; ENDIF 
    158       IF ( cn_trc_o(jpgsi) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpgsi) =  5.35e-8_wp  ; END WHERE ; ENDIF 
    159       IF ( cn_trc_o(jpphy) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpphy) =  8.10e-7_wp  ; END WHERE ; ENDIF 
    160       IF ( cn_trc_o(jpdia) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpdia) =  5.77e-7_wp  ; END WHERE ; ENDIF 
    161       IF ( cn_trc_o(jpzoo) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpzoo) =  6.68e-7_wp  ; END WHERE ; ENDIF 
    162       IF ( cn_trc_o(jpmes) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpmes) =  3.55e-7_wp  ; END WHERE ; ENDIF 
    163       IF ( cn_trc_o(jpfer) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpfer) =  1.62e-10_wp ; END WHERE ; ENDIF 
    164       IF ( cn_trc_o(jpsfe) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpsfe) =  2.29e-11_wp ; END WHERE ; ENDIF 
    165       IF ( cn_trc_o(jpdfe) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpdfe) =  8.75e-12_wp ; END WHERE ; ENDIF 
    166       IF ( cn_trc_o(jpnfe) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpnfe) =  1.48e-11_wp ; END WHERE ; ENDIF 
    167       IF ( cn_trc_o(jpnch) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpnch) =  2.02e-7_wp  ; END WHERE ; ENDIF 
    168       IF ( cn_trc_o(jpdch) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpdch) =  1.60e-7_wp  ; END WHERE ; ENDIF 
    169       IF ( cn_trc_o(jpno3) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpno3) =  2.64e-5_wp / rno3  ; END WHERE ; ENDIF 
    170       IF ( cn_trc_o(jpnh4) == 'AA ' ) THEN ; WHERE( gphit(:,:) <  00._wp ) ; trc_o(:,:,jpnh4) =  3.39e-7_wp / rno3  ; END WHERE ; ENDIF 
     175      zpisc(jpsil,3) =  4.96e-5_wp   
     176      zpisc(jpdsi,3) =  5.63e-7_wp  
     177      zpisc(jpgsi,3) =  5.35e-8_wp 
     178      zpisc(jpphy,3) =  8.10e-7_wp   
     179      zpisc(jpdia,3) =  5.77e-7_wp  
     180      zpisc(jpzoo,3) =  6.68e-7_wp 
     181      zpisc(jpmes,3) =  3.55e-7_wp   
     182      zpisc(jpfer,3) =  1.62e-10_wp 
     183      zpisc(jpsfe,3) =  2.29e-11_wp  
     184      zpisc(jpdfe,3) =  8.75e-12_wp 
     185      zpisc(jpnfe,3) =  1.48e-11_wp  
     186      zpisc(jpnch,3) =  2.02e-7_wp   
     187      zpisc(jpdch,3) =  1.60e-7_wp   
     188      zpisc(jpno3,3) =  2.64e-5_wp / rno3   
     189      zpisc(jpnh4,3) =  3.39e-7_wp / rno3   
    171190 
    172191      !--- Baltic Sea particular case for ORCA configurations 
    173       IF( cp_cfg == "orca" ) THEN            ! Baltic mask 
    174          WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND.    & 
    175                 54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp ) 
    176          trc_o(:,:,jpdic) = 1.14e-3_wp 
    177          trc_o(:,:,jpdoc) = 1.06e-5_wp 
    178          trc_o(:,:,jptal) = 1.16e-3_wp 
    179          trc_o(:,:,jpoxy) = 3.71e-4_wp 
    180          trc_o(:,:,jpcal) = 1.51e-9_wp 
    181          trc_o(:,:,jppo4) = 2.85e-9_wp / po4r 
    182          trc_o(:,:,jppoc) = 4.84e-7_wp 
     192      zpisc(jpdic,4) = 1.14e-3_wp 
     193      zpisc(jpdoc,4) = 1.06e-5_wp 
     194      zpisc(jptal,4) = 1.16e-3_wp 
     195      zpisc(jpoxy,4) = 3.71e-4_wp 
     196      zpisc(jpcal,4) = 1.51e-9_wp 
     197      zpisc(jppo4,4) = 2.85e-9_wp / po4r 
     198      zpisc(jppoc,4) = 4.84e-7_wp 
    183199#  if ! defined key_kriest 
    184          trc_o(:,:,jpgoc) = 1.05e-8_wp 
    185          trc_o(:,:,jpbfe) = 4.97e-13_wp 
     200      zpisc(jpgoc,4) = 1.05e-8_wp 
     201      zpisc(jpbfe,4) = 4.97e-13_wp 
    186202#  else 
    187          trc_o(:,:,jpnum) = 0. ! could not get this value 
     203      zpisc(jpnum,4) = 0. ! could not get this value 
    188204#  endif 
    189          trc_o(:,:,jpsil) = 4.91e-5_wp 
    190          trc_o(:,:,jpdsi) = 3.25e-7_wp 
    191          trc_o(:,:,jpgsi) = 1.93e-8_wp 
    192          trc_o(:,:,jpphy) = 6.64e-7_wp 
    193          trc_o(:,:,jpdia) = 3.41e-7_wp 
    194          trc_o(:,:,jpzoo) = 3.83e-7_wp 
    195          trc_o(:,:,jpmes) = 0.225e-6_wp 
    196          trc_o(:,:,jpfer) = 2.45e-9_wp 
    197          trc_o(:,:,jpsfe) = 3.89e-11_wp 
    198          trc_o(:,:,jpdfe) = 1.33e-11_wp 
    199          trc_o(:,:,jpnfe) = 2.62e-11_wp 
    200          trc_o(:,:,jpnch) = 1.17e-7_wp 
    201          trc_o(:,:,jpdch) = 9.69e-8_wp 
    202          trc_o(:,:,jpno3) = 5.36e-5_wp / rno3 
    203          trc_o(:,:,jpnh4) = 7.18e-7_wp / rno3 
    204          END WHERE 
    205       ENDIF ! cfg 
     205      zpisc(jpsil,4) = 4.91e-5_wp 
     206      zpisc(jpdsi,4) = 3.25e-7_wp 
     207      zpisc(jpgsi,4) = 1.93e-8_wp 
     208      zpisc(jpphy,4) = 6.64e-7_wp 
     209      zpisc(jpdia,4) = 3.41e-7_wp 
     210      zpisc(jpzoo,4) = 3.83e-7_wp 
     211      zpisc(jpmes,4) = 0.225e-6_wp 
     212      zpisc(jpfer,4) = 2.45e-9_wp 
     213      zpisc(jpsfe,4) = 3.89e-11_wp 
     214      zpisc(jpdfe,4) = 1.33e-11_wp 
     215      zpisc(jpnfe,4) = 2.62e-11_wp 
     216      zpisc(jpnch,4) = 1.17e-7_wp 
     217      zpisc(jpdch,4) = 9.69e-8_wp 
     218      zpisc(jpno3,4) = 5.36e-5_wp / rno3 
     219      zpisc(jpnh4,4) = 7.18e-7_wp / rno3 
     220  
     221      DO jn = jp_pcs0, jp_pcs1 
     222         IF( cn_trc_o(jn) == 'GL ' ) trc_o(:,:,jn) = zpisc(jn,1)  ! Global case 
     223         IF( cn_trc_o(jn) == 'AA ' ) THEN  
     224            WHERE( gphit(:,:) >= 0._wp ) ; trc_o(:,:,jn) = zpisc(jn,2) ; END WHERE ! Arctic  
     225            WHERE( gphit(:,:) <  0._wp ) ; trc_o(:,:,jn) = zpisc(jn,3) ; END WHERE ! Antarctic  
     226         ENDIF 
     227         IF( cp_cfg == "orca" ) THEN     !  Baltic Sea particular case for ORCA configurations 
     228             WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND.    & 
     229                    54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp ) 
     230                    trc_o(:,:,jn) = zpisc(jn,4) 
     231            END WHERE 
     232         ENDIF  
     233      ENDDO 
     234 
     235 
    206236 
    207237      !----------------------------- 
     
    217247 
    218248      DO jn = jp_pcs0, jp_pcs1 
    219          IF ( trc_ice_ratio(jn) >= 0._wp )  zratio(jn,:) = trc_ice_ratio(jn) 
    220          IF ( trc_ice_ratio(jn) == -1._wp ) zratio(jn,:) = zrs(:) 
    221          IF ( trc_ice_ratio(jn) == -2._wp ) zratio(jn,:) = -9999.99_wp 
     249         IF( trc_ice_ratio(jn) >= 0._wp )  zratio(jn,:) = trc_ice_ratio(jn) 
     250         IF( trc_ice_ratio(jn) == -1._wp ) zratio(jn,:) = zrs(:) 
     251         IF( trc_ice_ratio(jn) == -2._wp ) zratio(jn,:) = -9999.99_wp 
    222252      END DO 
    223253 
     
    227257      DO jn = jp_pcs0, jp_pcs1 
    228258         !-- Everywhere but in the Baltic 
    229          IF ( trc_ice_ratio(jn) >= -1._wp ) THEN !! no prescribed concentration 
    230                                               !! (typically everything but iron)  
     259         IF ( trc_ice_ratio(jn) >= -1._wp ) THEN ! no prescribed conc. ; typically everything but iron)  
    231260            trc_i(:,:,jn) = zratio(jn,1) * trc_o(:,:,jn)  
    232          ELSE                                 !! prescribed concentration 
     261         ELSE                                    ! prescribed concentration 
    233262            trc_i(:,:,jn) = trc_ice_prescr(jn) 
    234263         ENDIF 
    235264        
    236265         !-- Baltic 
    237          IF( cp_cfg == "orca" ) THEN !! Baltic treated seperately for ORCA configs 
    238             IF ( trc_ice_ratio(jn) >= - 1._wp ) THEN !! no prescribed concentration 
    239                                                  !! (typically everything but iron)  
     266         IF( cp_cfg == "orca" ) THEN  ! Baltic treated seperately for ORCA configs 
     267            IF ( trc_ice_ratio(jn) >= - 1._wp ) THEN ! no prescribed conc. ; typically everything but iron)  
    240268               WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND.    & 
    241269                      54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp ) 
    242270                     trc_i(:,:,jn) = zratio(jn,2) * trc_o(:,:,jn)  
    243271               END WHERE 
    244             ELSE                                 !! prescribed tracer concentration in ice 
     272            ELSE                                 ! prescribed tracer concentration in ice 
    245273               WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND.    & 
    246274                   54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp ) 
     
    251279      ! 
    252280      END DO ! jn 
    253  
    254    END SUBROUTINE trc_ice_ini_pisces 
     281#endif 
     282 
     283   END SUBROUTINE p4z_ice_ini 
     284 
     285   SUBROUTINE p2z_ice_ini 
     286#if defined key_pisces_reduced  
     287      !!---------------------------------------------------------------------- 
     288      !!                   ***  ROUTINE p2z_ice_ini *** 
     289      !! 
     290      !! ** Purpose :   Initialisation of the LOBSTER biochemical model 
     291      !!---------------------------------------------------------------------- 
     292#endif 
     293   END SUBROUTINE p2z_ice_ini 
     294 
    255295 
    256296#else 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/PISCES/trcini_pisces.F90

    r5602 r7256  
    115115      po4r    =   1._wp / 122._wp 
    116116      o2nit   =  32._wp / 122._wp 
    117       rdenit  = 105._wp /  16._wp 
     117      o2ut    = 133._wp / 122._wp 
     118      rdenit  =  ( ( o2ut + o2nit ) * 0.80 - rno3 - rno3 * 0.60 ) / rno3 
    118119      rdenita =   3._wp /  5._wp 
    119       o2ut    = 133._wp / 122._wp 
     120 
    120121 
    121122      ! Initialization of tracer concentration in case of  no restart  
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/TRP/trcdmp.F90

    r6101 r7256  
    3535   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   restotr   ! restoring coeff. on tracers (s-1) 
    3636 
    37    INTEGER, PARAMETER           ::   npncts   = 5        ! number of closed sea 
     37   INTEGER, PARAMETER           ::   npncts   = 8        ! number of closed sea 
    3838   INTEGER, DIMENSION(npncts)   ::   nctsi1, nctsj1      ! south-west closed sea limits (i,j) 
    3939   INTEGER, DIMENSION(npncts)   ::   nctsi2, nctsj2      ! north-east closed sea limits (i,j) 
     
    107107                
    108108               jl = n_trc_index(jn)  
    109                CALL trc_dta( kt, sf_trcdta(jl),rf_trfac(jl) )   ! read tracer data at nit000 
    110                ztrcdta(:,:,:) = sf_trcdta(jl)%fnow(:,:,:) 
     109               CALL trc_dta( kt, sf_trcdta(jl), rf_trfac(jl), ztrcdta )   ! read tracer data at nit000 
    111110 
    112111               SELECT CASE ( nn_zdmp_tr ) 
     
    208207            ! 
    209208                                                        ! Caspian Sea 
    210             nctsi1(1)   = 332  ; nctsj1(1)   = 243 - isrow 
    211             nctsi2(1)   = 344  ; nctsj2(1)   = 275 - isrow 
     209            nctsi1(1)   = 333  ; nctsj1(1)   = 243 - isrow 
     210            nctsi2(1)   = 342  ; nctsj2(1)   = 274 - isrow 
     211            !                                           ! Lake Superior 
     212            nctsi1(2)   = 198  ; nctsj1(2)   = 258 - isrow 
     213            nctsi2(2)   = 204  ; nctsj2(2)   = 262 - isrow 
     214            !                                           ! Lake Michigan 
     215            nctsi1(3)   = 201  ; nctsj1(3)   = 250 - isrow 
     216            nctsi2(3)   = 203  ; nctsj2(3)   = 256 - isrow 
     217            !                                           ! Lake Huron 
     218            nctsi1(4)   = 204  ; nctsj1(4)   = 252 - isrow 
     219            nctsi2(4)   = 209  ; nctsj2(4)   = 256 - isrow 
     220            !                                           ! Lake Erie 
     221            nctsi1(5)   = 206  ; nctsj1(5)   = 249 - isrow 
     222            nctsi2(5)   = 209  ; nctsj2(5)   = 251 - isrow 
     223            !                                           ! Lake Ontario 
     224            nctsi1(6)   = 210  ; nctsj1(6)   = 252 - isrow 
     225            nctsi2(6)   = 212  ; nctsj2(6)   = 252 - isrow 
     226            !                                           ! Victoria Lake 
     227            nctsi1(7)   = 321  ; nctsj1(7)   = 180 - isrow 
     228            nctsi2(7)   = 322  ; nctsj2(7)   = 189 - isrow 
     229            !                                           ! Baltic Sea 
     230            nctsi1(8)   = 297  ; nctsj1(8)   = 270 - isrow 
     231            nctsi2(8)   = 308  ; nctsj2(8)   = 293 - isrow 
    212232            !                                         
    213233            !                                           ! ======================= 
     
    283303            IF( ln_trc_ini(jn) ) THEN      ! update passive tracers arrays with input data read from file 
    284304                jl = n_trc_index(jn) 
    285                 CALL trc_dta( kt, sf_trcdta(jl),rf_trfac(jl) )   ! read tracer data at nit000 
    286                 ztrcdta(:,:,:) = sf_trcdta(jl)%fnow(:,:,:) 
     305                CALL trc_dta( kt, sf_trcdta(jl), rf_trfac(jl), ztrcdta )   ! read tracer data at nit000 
    287306                DO jc = 1, npncts 
    288307                   DO jk = 1, jpkm1 
    289308                      DO jj = nctsj1(jc), nctsj2(jc) 
    290309                         DO ji = nctsi1(jc), nctsi2(jc) 
    291                             trn(ji,jj,jk,jn) = ztrcdta(ji,jj,jk) * tmask(ji,jj,jk) 
     310                            trn(ji,jj,jk,jn) = ztrcdta(ji,jj,jk) 
    292311                            trb(ji,jj,jk,jn) = trn(ji,jj,jk,jn) 
    293312                         ENDDO 
     
    317336      IF( nn_timing == 1 )  CALL timing_start('trc_dmp_init') 
    318337      ! 
     338      !Allocate arrays 
     339      IF( trc_dmp_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'trc_dmp_init: unable to allocate arrays' ) 
    319340 
    320341      IF( lzoom )   nn_zdmp_tr = 0           ! restoring to climatology at closed north or south boundaries 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90

    r5602 r7256  
    1818   USE trc             ! ocean passive tracers variables 
    1919   USE trcnam_trp      ! passive tracers transport namelist variables 
    20    USE ldftra_oce      ! lateral diffusion coefficient on tracers 
     20   USE ldftra_oce,ONLY: ln_traldf_grif,rn_aht_0,rn_ahtb_0,lk_traldf_eiv     ! lateral diffusion coefficient on tracers 
    2121   USE ldfslp          ! ??? 
    2222   USE traldf_bilapg   ! lateral mixing            (tra_ldf_bilapg routine) 
     
    5656      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index 
    5757      !! 
    58       INTEGER            :: jn 
     58      INTEGER            :: ji, jj, jk, jn 
     59      REAL(wp)           :: zdep 
    5960      CHARACTER (len=22) :: charout 
    6061      REAL(wp), POINTER, DIMENSION(:,:,:,:) ::   ztrtrd 
     
    6667 
    6768      rldf = rldf_rat 
    68  
     69      ! 
     70      r_fact_lap(:,:,:) = 1. 
     71      DO jk= 1, jpk 
     72         DO jj = 1, jpj 
     73            DO ji = 1, jpi 
     74               IF( fsdept(ji,jj,jk) > 200. .AND. gphit(ji,jj) < 5. .AND. gphit(ji,jj) > -5. ) THEN 
     75                  zdep = MAX( fsdept(ji,jj,jk) - 1000., 0. ) / 1000. 
     76                  r_fact_lap(ji,jj,jk) = MAX( 1., rn_fact_lap * EXP( -zdep ) ) 
     77               ENDIF 
     78            END DO 
     79         END DO 
     80      END DO 
     81      ! 
    6982      IF( l_trdtrc )  THEN 
    7083         CALL wrk_alloc( jpi, jpj, jpk, jptra, ztrtrd ) 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/TRP/trcnam_trp.F90

    r5602 r7256  
    4040   REAL(wp), PUBLIC ::   rn_ahtrc_0          !: diffusivity coefficient for passive tracer (m2/s) 
    4141   REAL(wp), PUBLIC ::   rn_ahtrb_0          !: background diffusivity coefficient for passive tracer (m2/s) 
     42   REAL(wp), PUBLIC ::   rn_fact_lap         !: Enhanced zonal diffusivity coefficent in the equatorial domain 
    4243 
    4344   !                                        !!: ** Treatment of Negative concentrations ( nam_trcrad ) 
     
    7475      NAMELIST/namtrc_ldf/ ln_trcldf_lap  ,     & 
    7576         &                 ln_trcldf_bilap, ln_trcldf_level,     & 
    76          &                 ln_trcldf_hor  , ln_trcldf_iso  , rn_ahtrc_0, rn_ahtrb_0 
     77         &                 ln_trcldf_hor  , ln_trcldf_iso  , rn_ahtrc_0, rn_ahtrb_0,   & 
     78         &                 rn_fact_lap 
     79 
    7780      NAMELIST/namtrc_zdf/ ln_trczdf_exp  , nn_trczdf_exp 
    7881      NAMELIST/namtrc_rad/ ln_trcrad 
     
    127130         WRITE(numout,*) '      diffusivity coefficient                                 rn_ahtrc_0 = ', rn_ahtrc_0 
    128131         WRITE(numout,*) '      background hor. diffusivity                             rn_ahtrb_0 = ', rn_ahtrb_0 
     132         WRITE(numout,*) '      enhanced zonal diffusivity                             rn_fact_lap = ', rn_fact_lap 
    129133      ENDIF 
    130134 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/TRP/trcnxt.F90

    r7210 r7256  
    104104      ENDIF 
    105105 
     106#if defined key_agrif 
     107      CALL Agrif_trc                   ! AGRIF zoom boundaries 
     108#endif 
    106109      ! Update after tracer on domain lateral boundaries 
    107110      DO jn = 1, jptra 
     
    112115#if defined key_bdy 
    113116!!      CALL bdy_trc( kt )               ! BDY open boundaries 
    114 #endif 
    115 #if defined key_agrif 
    116       CALL Agrif_trc                   ! AGRIF zoom boundaries 
    117117#endif 
    118118 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/TRP/trcsbc.F90

    r7210 r7256  
    103103         IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
    104104 
    105          IF( ln_rsttr .AND.    &                     ! Restart: read in restart  file 
     105         IF( ln_rsttr .AND. .NOT.ln_top_euler .AND.   &                     ! Restart: read in restart  file 
    106106            iom_varid( numrtr, 'sbc_'//TRIM(ctrcnm(1))//'_b', ldstop = .FALSE. ) > 0 ) THEN 
    107107            IF(lwp) WRITE(numout,*) '          nittrc000-nn_dttrc surface tracer content forcing fields red in the restart file' 
     
    172172            END DO 
    173173         ENDIF 
     174         ! 
     175         CALL lbc_lnk( sbc_trc(:,:,jn), 'T', 1. ) 
    174176         !                                       Concentration dilution effect on tracers due to evaporation & precipitation  
    175177         DO jj = 2, jpj 
     
    190192      !                                           Write in the tracer restar  file 
    191193      !                                          ******************************* 
    192       IF( lrst_trc ) THEN 
     194      IF( lrst_trc .AND. .NOT.ln_top_euler ) THEN 
    193195         IF(lwp) WRITE(numout,*) 
    194196         IF(lwp) WRITE(numout,*) 'sbc : ocean surface tracer content forcing fields written in tracer restart file ',   & 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/TRP/trctrp.F90

    r7217 r7256  
    1818   USE trcnam_trp      ! passive tracers transport namelist variables 
    1919   USE trabbl          ! bottom boundary layer               (trc_bbl routine) 
    20    USE trabbl_crs      ! bottom boundary layer               (trc_bbl routine) 
    2120   USE trcbbl          ! bottom boundary layer               (trc_bbl routine) 
    22    USE trcbbl_crs      ! bottom boundary layer               (trc_bbl routine) 
    2321   USE zdfkpp          ! KPP non-local tracer fluxes         (trc_kpp routine) 
    2422   USE trcdmp          ! internal damping                    (trc_dmp routine) 
     
    7674      IF( .NOT. lk_c1d ) THEN 
    7775         ! 
    78                                CALL trc_sbc( kstp ) 
    79          IF( ln_crs_top ) THEN ;    CALL trc_bbl_crs( kstp ) 
    80          ELSE              ;    CALL trc_bbl( kstp ) 
    81          ENDIF 
    82          IF( ln_trcdmp )        CALL trc_dmp( kstp )            ! internal damping trends 
    83  
    84          IF( ln_crs_top ) THEN ;    CALL trc_adv_crs( kstp ) 
    85          ELSE              ;    CALL trc_adv( kstp ) 
     76                                       CALL trc_sbc( kstp ) 
     77         IF( lk_trabbl )               CALL trc_bbl( kstp ) 
     78         IF( ln_trcdmp )               CALL trc_dmp( kstp )            ! internal damping trends 
     79         IF( ln_crs_top ) THEN    ;    CALL trc_adv_crs( kstp ) 
     80         ELSE                     ;    CALL trc_adv( kstp ) 
    8681         ENDIF 
    8782 
    88          IF( ln_trcdmp_clo )    CALL trc_dmp_clo( kstp )        ! internal damping trends on closed seas only 
    89          IF( ln_crs_top ) THEN ;    CALL trc_ldf_crs( kstp ) 
    90          ELSE              ;    CALL trc_ldf( kstp ) 
     83         IF( ln_zps ) THEN 
     84           IF( ln_crs_top ) THEN 
     85              CALL zps_hde_crs( kstp, jptra, trn, gtru, gtrv ) 
     86           ELSE 
     87              IF( ln_isfcav ) THEN ; CALL zps_hde_isf( kstp, jptra, trb, pgtu=gtru, pgtv=gtrv, pgtui=gtrui, pgtvi=gtrvi )  ! both top & bottom 
     88              ELSE                 ; CALL zps_hde    ( kstp, jptra, trb, gtru, gtrv )                                      !  only bottom 
     89              ENDIF 
     90           ENDIF 
     91         ENDIF 
     92 
     93         IF( ln_crs_top ) THEN    ;    CALL trc_ldf_crs( kstp ) 
     94         ELSE                     ;    CALL trc_ldf( kstp ) 
    9195         ENDIF 
    9296         IF( .NOT. lk_offline .AND. lk_zdfkpp )    & 
    93             &                   CALL trc_kpp( kstp )            ! KPP non-local tracer fluxes 
     97            &                          CALL trc_kpp( kstp )            ! KPP non-local tracer fluxes 
    9498#if defined key_agrif 
    9599         IF(.NOT. Agrif_Root()) CALL Agrif_Sponge_trc           ! tracers sponge 
    96100#endif 
    97          IF( ln_crs_top ) THEN ;    CALL trc_zdf_crs( kstp ) 
    98          ELSE              ;    CALL trc_zdf( kstp ) 
     101         IF( ln_crs_top ) THEN    ;    CALL trc_zdf_crs( kstp ) 
     102         ELSE                     ;    CALL trc_zdf( kstp ) 
    99103         ENDIF 
    100  
    101                                 CALL trc_nxt( kstp )            ! tracer fields at next time step      
    102          IF( ln_trcrad )        CALL trc_rad( kstp )            ! Correct artificial negative concentrations 
     104         ! 
     105                                       CALL trc_nxt( kstp )            ! tracer fields at next time step      
     106         IF( ln_trcrad )               CALL trc_rad( kstp )            ! Correct artificial negative concentrations 
     107         IF( ln_trcdmp_clo )           CALL trc_dmp_clo( kstp )        ! internal damping trends on closed seas only 
    103108 
    104109#if defined key_agrif 
    105110      IF( .NOT. Agrif_Root())   CALL Agrif_Update_Trc( kstp )   ! Update tracer at AGRIF zoom boundaries : children only 
    106111#endif 
    107           ! Partial steps: now horizontal gradient of passive 
    108          IF( ln_zps    )THEN 
    109             IF( ln_crs_top ) THEN   
    110                CALL zps_hde_crs( kstp, jptra, trn, gtru, gtrv ) 
    111             ELSE 
    112                IF( ln_isfcav )THEN 
    113                   CALL zps_hde_isf( kstp, jptra, trn, pgtu=gtru, pgtv=gtrv, pgtui=gtrui, pgtvi=gtrvi )  ! Partial steps: now horizontal gradient of passive 
    114                ELSE 
    115                   CALL zps_hde    ( kstp, jptra, trn, gtru, gtrv )   ! Partial steps: now horizontal gradient of passive 
    116                ENDIF 
    117             ENDIF 
    118          ENDIF 
    119                                                                 ! tracers at the bottom ocean level 
    120          ! 
     112 
    121113      ELSE                                               ! 1D vertical configuration 
    122114                                CALL trc_sbc( kstp )            ! surface boundary condition 
     
    130122      ! 
    131123      IF( nn_timing == 1 )   CALL timing_stop('trc_trp') 
     124      ! 
     1259400  FORMAT(a25,i4,D23.16) 
    132126      ! 
    133127   END SUBROUTINE trc_trp 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/oce_trc.F90

    r7210 r7256  
    219219   USE crs , ONLY :  ahtw     =>   ahtw_crs        !: lateral diffusivity coef. at w-points  
    220220   USE crs , ONLY :  ahtt     =>   ahtt_crs        !: lateral diffusivity coef. at t-points 
     221   USE crs , ONLY :  r_fact_lap     =>  r_fact_lap_crs        !: enhanced zonal diffusivity coefficient 
    221222   USE ldftra_oce , ONLY :  rldf     =>   rldf 
    222223   USE crs , ONLY :  trc_i => trc_i_crs 
     
    459460   USE ldftra_oce , ONLY :  aeiw     =>   aeiw        !: eddy induced velocity coef. at w-points (m2/s)  
    460461   USE ldftra_oce , ONLY :  lk_traldf_eiv  =>  lk_traldf_eiv     !: eddy induced velocity flag 
     462   USE ldftra_oce , ONLY :  r_fact_lap     =>  r_fact_lap        !: enhanced zonal diffusivity coefficient 
    461463 
    462464   !* vertical diffusion * 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/trcdta.F90

    r6772 r7256  
    7777      ALLOCATE( n_trc_index(ntrc), slf_i(ntrc), STAT=ierr0 ) 
    7878      IF( ierr0 > 0 ) THEN 
    79          CALL ctl_stop( 'trc_nam: unable to allocate n_trc_index' )   ;   RETURN 
     79         CALL ctl_stop( 'trc_dta_init: unable to allocate n_trc_index' )   ;   RETURN 
    8080      ENDIF 
    8181      nb_trcdta      = 0 
     
    9191      IF(lwp) THEN 
    9292         WRITE(numout,*) ' ' 
     93         WRITE(numout,*) 'trc_dta_init : Passive tracers Initial Conditions ' 
     94         WRITE(numout,*) '~~~~~~~~~~~~~~ ' 
    9395         WRITE(numout,*) ' number of passive tracers to be initialize by data :', ntra 
    9496         WRITE(numout,*) ' ' 
     
    107109         DO jn = 1, ntrc 
    108110            IF( ln_trc_ini(jn) )  THEN    ! open input file only if ln_trc_ini(jn) is true 
    109                clndta = TRIM( sn_trcdta(jn)%clvar )  
    110                clntrc = TRIM( ctrcnm   (jn)       )  
     111               clndta = TRIM( sn_trcdta(jn)%clvar ) 
     112               if (jn > jptra) then 
     113                  clntrc='Dummy' ! By pass weird formats in ocean.output if ntrc > jptra 
     114               else 
     115                  clntrc = TRIM( ctrcnm   (jn)       ) 
     116               endif 
    111117               zfact  = rn_trfac(jn) 
    112                IF( clndta /=  clntrc ) THEN  
    113                   CALL ctl_warn( 'trc_dta_init: passive tracer data initialisation :  ',   & 
    114                   &              'the variable name in the data file : '//clndta//   &  
    115                   &              '  must be the same than the name of the passive tracer : '//clntrc//' ') 
     118               IF( clndta /=  clntrc ) THEN 
     119                  CALL ctl_warn( 'trc_dta_init: passive tracer data initialisation    ',   & 
     120                  &              'Input name of data file : '//TRIM(clndta)//   & 
     121                  &              ' differs from that of tracer : '//TRIM(clntrc)//' ') 
    116122               ENDIF 
    117                WRITE(numout,*) ' read an initial file for passive tracer number :', jn, ' name : ', clndta, &  
    118                &               ' multiplicative factor : ', zfact 
     123               WRITE(numout,'(a, i4,3a,e11.3)') ' Read IC file for tracer number :', & 
     124               &            jn, ', name : ', TRIM(clndta), ', Multiplicative Scaling factor : ', zfact 
    119125            ENDIF 
    120126         END DO 
     
    124130         ALLOCATE( sf_trcdta(nb_trcdta), rf_trfac(nb_trcdta), STAT=ierr1 ) 
    125131         IF( ierr1 > 0 ) THEN 
    126             CALL ctl_stop( 'trc_dta_ini: unable to allocate  sf_trcdta structure' )   ;   RETURN 
     132            CALL ctl_stop( 'trc_dta_init: unable to allocate  sf_trcdta structure' )   ;   RETURN 
    127133         ENDIF 
    128134         ! 
     
    135141               IF( sn_trcdta(jn)%ln_tint )  ALLOCATE( sf_trcdta(jl)%fdta(jpi,jpj,jpk,2) , STAT=ierr3 ) 
    136142               IF( ierr2 + ierr3 > 0 ) THEN 
    137                  CALL ctl_stop( 'trc_dta : unable to allocate passive tracer data arrays' )   ;   RETURN 
     143                 CALL ctl_stop( 'trc_dta_init : unable to allocate passive tracer data arrays' )   ;   RETURN 
    138144               ENDIF 
    139145            ENDIF 
     
    141147         ENDDO 
    142148         !                         ! fill sf_trcdta with slf_i and control print 
    143          CALL fld_fill( sf_trcdta, slf_i, cn_dir, 'trc_dta', 'Passive tracer data', 'namtrc' ) 
     149         CALL fld_fill( sf_trcdta, slf_i, cn_dir, 'trc_dta_init', 'Passive tracer data', 'namtrc' ) 
    144150         ! 
    145151      ENDIF 
     
    151157 
    152158 
    153    SUBROUTINE trc_dta( kt, sf_dta, zrf_trfac ) 
     159   SUBROUTINE trc_dta( kt, sf_dta, ptrfac, ptrc) 
    154160      !!---------------------------------------------------------------------- 
    155161      !!                   ***  ROUTINE trc_dta  *** 
     
    164170      !!---------------------------------------------------------------------- 
    165171      INTEGER                     , INTENT(in   ) ::   kt     ! ocean time-step 
    166       TYPE(FLD), DIMENSION(1)   , INTENT(inout) ::   sf_dta     ! array of information on the field to read 
    167       REAL(wp)                  , INTENT(in   ) ::   zrf_trfac  ! multiplication factor 
     172      TYPE(FLD), DIMENSION(1)     , INTENT(inout) ::   sf_dta     ! array of information on the field to read 
     173      REAL(wp)                    , INTENT(in   ) ::   ptrfac  ! multiplication factor 
     174      REAL(wp), DIMENSION(jpi,jpj,jpk), OPTIONAL  , INTENT(out  ) ::   ptrc 
    168175      ! 
    169176      INTEGER ::   ji, jj, jk, jl, jkk, ik    ! dummy loop indices 
    170177      REAL(wp)::   zl, zi 
    171178      REAL(wp), DIMENSION(jpk) ::  ztp                ! 1D workspace 
     179      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrcdta   ! 3D  workspace 
    172180      CHARACTER(len=100) :: clndta 
    173181      !!---------------------------------------------------------------------- 
     
    177185      IF( nb_trcdta > 0 ) THEN 
    178186         ! 
     187         CALL wrk_alloc( jpi, jpj, jpk, ztrcdta )    ! Memory allocation 
     188         ! 
    179189         CALL fld_read( kt, 1, sf_dta )      !==   read data at kt time step   ==! 
     190         ztrcdta(:,:,:) = sf_dta(1)%fnow(:,:,:) * tmask(:,:,:)    ! Mask 
    180191         ! 
    181192         IF( ln_sco ) THEN                   !==   s- or mixed s-zps-coordinate   ==! 
     
    186197            ENDIF 
    187198            ! 
    188                DO jj = 1, jpj                         ! vertical interpolation of T & S 
     199            DO jj = 1, jpj                         ! vertical interpolation of T & S 
     200               DO ji = 1, jpi 
     201                  DO jk = 1, jpk                        ! determines the intepolated T-S profiles at each (i,j) points 
     202                     zl = fsdept_n(ji,jj,jk) 
     203                     IF(     zl < gdept_1d(1  ) ) THEN         ! above the first level of data 
     204                        ztp(jk) = ztrcdta(ji,jj,1) 
     205                     ELSEIF( zl > gdept_1d(jpk) ) THEN         ! below the last level of data 
     206                        ztp(jk) =  ztrcdta(ji,jj,jpkm1) 
     207                     ELSE                                      ! inbetween : vertical interpolation between jkk & jkk+1 
     208                        DO jkk = 1, jpkm1                                  ! when  gdept(jkk) < zl < gdept(jkk+1) 
     209                           IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN 
     210                              zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk)) 
     211                              ztp(jk) = ztrcdta(ji,jj,jkk) + ( ztrcdta(ji,jj,jkk+1) - & 
     212                                        ztrcdta(ji,jj,jkk) ) * zi  
     213                           ENDIF 
     214                        END DO 
     215                     ENDIF 
     216                  END DO 
     217                  DO jk = 1, jpkm1 
     218                    ztrcdta(ji,jj,jk) = ztp(jk) * tmask(ji,jj,jk)     ! mask required for mixed zps-s-coord 
     219                  END DO 
     220                  ztrcdta(ji,jj,jpk) = 0._wp 
     221                END DO 
     222            END DO 
     223            !  
     224         ELSE                                !==   z- or zps- coordinate   ==! 
     225            ! 
     226            IF( ln_zps ) THEN                      ! zps-coordinate (partial steps) interpolation at the last ocean level 
     227               DO jj = 1, jpj 
    189228                  DO ji = 1, jpi 
    190                      DO jk = 1, jpk                        ! determines the intepolated T-S profiles at each (i,j) points 
    191                         zl = fsdept_n(ji,jj,jk) 
    192                         IF(     zl < gdept_1d(1  ) ) THEN         ! above the first level of data 
    193                            ztp(jk) =  sf_dta(1)%fnow(ji,jj,1) 
    194                         ELSEIF( zl > gdept_1d(jpk) ) THEN         ! below the last level of data 
    195                            ztp(jk) =  sf_dta(1)%fnow(ji,jj,jpkm1) 
    196                         ELSE                                      ! inbetween : vertical interpolation between jkk & jkk+1 
    197                            DO jkk = 1, jpkm1                                  ! when  gdept(jkk) < zl < gdept(jkk+1) 
    198                               IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN 
    199                                  zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk)) 
    200                                  ztp(jk) = sf_dta(1)%fnow(ji,jj,jkk) + ( sf_dta(1)%fnow(ji,jj,jkk+1) - & 
    201                                            sf_dta(1)%fnow(ji,jj,jkk) ) * zi  
    202                               ENDIF 
    203                            END DO 
    204                         ENDIF 
    205                      END DO 
    206                      DO jk = 1, jpkm1 
    207                         sf_dta(1)%fnow(ji,jj,jk) = ztp(jk) * tmask(ji,jj,jk)     ! mask required for mixed zps-s-coord 
    208                      END DO 
    209                      sf_dta(1)%fnow(ji,jj,jpk) = 0._wp 
     229                     ik = mbkt(ji,jj)  
     230                     IF( ik > 1 ) THEN 
     231                        zl = ( gdept_1d(ik) - fsdept_n(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) ) 
     232                        ztrcdta(ji,jj,ik) = (1.-zl) * ztrcdta(ji,jj,ik) + zl * ztrcdta(ji,jj,ik-1) 
     233                     ENDIF 
     234                     ik = mikt(ji,jj) 
     235                     IF( ik > 1 ) THEN 
     236                        zl = ( fsdept_n(ji,jj,ik) - gdept_1d(ik) ) / ( gdept_1d(ik+1) - gdept_1d(ik) ) 
     237                        ztrcdta(ji,jj,ik) = (1.-zl) * ztrcdta(ji,jj,ik) + zl * ztrcdta(ji,jj,ik+1) 
     238                     ENDIF 
    210239                  END DO 
    211240               END DO 
    212             !  
    213          ELSE                                !==   z- or zps- coordinate   ==! 
    214             !                              
    215                sf_dta(1)%fnow(:,:,:) = sf_dta(1)%fnow(:,:,:) * tmask(:,:,:)    ! Mask 
    216                ! 
    217                IF( ln_zps ) THEN                      ! zps-coordinate (partial steps) interpolation at the last ocean level 
    218                   DO jj = 1, jpj 
    219                      DO ji = 1, jpi 
    220                         ik = mbkt(ji,jj)  
    221                         IF( ik > 1 ) THEN 
    222                            zl = ( gdept_1d(ik) - fsdept_n(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) ) 
    223                            sf_dta(1)%fnow(ji,jj,ik) = (1.-zl) * sf_dta(1)%fnow(ji,jj,ik) + zl * sf_dta(1)%fnow(ji,jj,ik-1) 
    224                         ENDIF 
    225                         ik = mikt(ji,jj) 
    226                         IF( ik > 1 ) THEN 
    227                            zl = ( gdept_0(ji,jj,ik) - gdept_1d(ik) ) / ( gdept_1d(ik+1) - gdept_1d(ik) ) 
    228                            sf_dta(1)%fnow(ji,jj,ik) = (1.-zl) * sf_dta(1)%fnow(ji,jj,ik) + zl * sf_dta(1)%fnow(ji,jj,ik+1) 
    229                         ENDIF 
    230                      END DO 
    231                   END DO 
    232                ENDIF 
    233             ! 
    234          ENDIF 
    235          ! 
    236          sf_dta(1)%fnow(:,:,:) = sf_dta(1)%fnow(:,:,:) * zrf_trfac   !  multiplicative factor 
     241            ENDIF 
     242            ! 
     243         ENDIF 
     244         ! 
     245         ! Add multiplicative factor 
     246         ztrcdta(:,:,:) = ztrcdta(:,:,:) * ptrfac 
     247         ! 
     248         ! Data structure for trc_ini (and BFMv5.1 coupling) 
     249         IF( .NOT. PRESENT(ptrc) ) sf_dta(1)%fnow(:,:,:) = ztrcdta(:,:,:) 
     250         ! 
     251         ! Data structure for trc_dmp 
     252         IF( PRESENT(ptrc) )  ptrc(:,:,:) = ztrcdta(:,:,:) 
    237253         ! 
    238254         IF( lwp .AND. kt == nit000 ) THEN 
     
    241257               WRITE(numout,*) 
    242258               WRITE(numout,*)'  level = 1' 
    243                CALL prihre( sf_dta(1)%fnow(:,:,1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout ) 
     259               CALL prihre( ztrcdta(:,:,1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout ) 
    244260               WRITE(numout,*)'  level = ', jpk/2 
    245                CALL prihre( sf_dta(1)%fnow(:,:,jpk/2), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout ) 
     261               CALL prihre( ztrcdta(:,:,jpk/2), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout ) 
    246262               WRITE(numout,*)'  level = ', jpkm1 
    247                CALL prihre( sf_dta(1)%fnow(:,:,jpkm1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout ) 
     263               CALL prihre( ztrcdta(:,:,jpkm1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout ) 
    248264               WRITE(numout,*) 
    249265         ENDIF 
     266         ! 
     267         CALL wrk_dealloc( jpi, jpj, jpk, ztrcdta ) 
     268         ! 
    250269      ENDIF 
    251270      ! 
     
    258277   !!---------------------------------------------------------------------- 
    259278CONTAINS 
    260    SUBROUTINE trc_dta( kt, sf_dta, zrf_trfac )        ! Empty routine 
     279   SUBROUTINE trc_dta( kt, sf_dta, ptrfac, ptrc)        ! Empty routine 
    261280      WRITE(*,*) 'trc_dta: You should not have seen this print! error?', kt 
    262281   END SUBROUTINE trc_dta 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/trcini.F90

    r7217 r7256  
    2525   USE trcini_my_trc   ! MY_TRC   initialisation 
    2626   USE trcdta          ! initialisation from files 
    27    USE zpshde,ONLY: zps_hde, zps_hde_isf    ! partial step: hor. derivative   (zps_hde routine) 
    28    USE zpshde_crs      ! partial step: hor. derivative   (zps_hde routine) 
     27   USE daymod          ! calendar manager 
    2928   USE prtctl_trc      ! Print control passive tracers (prt_ctl_trc_init routine) 
    3029   USE trcsub          ! variables to substep passive tracers 
     
    6362      INTEGER ::   jk, jn, jl    ! dummy loop indices 
    6463      CHARACTER (len=25) :: charout 
    65       REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrcdta   ! 4D  workspace 
    6664      !!--------------------------------------------------------------------- 
    6765      ! 
     
    123121        IF( ln_trcdta .AND. nb_trcdta > 0 ) THEN  ! Initialisation of tracer from a file that may also be used for damping 
    124122            ! 
    125             CALL wrk_alloc( jpi, jpj, jpk, ztrcdta )    ! Memory allocation 
    126             ! 
    127123            DO jn = 1, jptra 
    128124               IF( ln_trc_ini(jn) ) THEN      ! update passive tracers arrays with input data read from file 
    129125                  jl = n_trc_index(jn)  
    130                   CALL trc_dta( nit000, sf_trcdta(jl),rf_trfac(jl) )   ! read tracer data at nit000 
    131                   ztrcdta(:,:,:) = sf_trcdta(jl)%fnow(:,:,:) 
    132                   trn(:,:,:,jn) = ztrcdta(:,:,:) * tmask(:,:,:)   
     126                  CALL trc_dta( nit000, sf_trcdta(jl), rf_trfac(jl) )   ! read tracer data at nit000 
     127                  trn(:,:,:,jn) = sf_trcdta(jl)%fnow(:,:,:)  
    133128                  IF( .NOT.ln_trcdmp .AND. .NOT.ln_trcdmp_clo ) THEN      !== deallocate data structure   ==! 
    134129                     !                                                    (data used only for initialisation) 
     
    140135               ENDIF 
    141136            ENDDO 
    142             CALL wrk_dealloc( jpi, jpj, jpk, ztrcdta ) 
     137            ! 
    143138        ENDIF 
    144139        ! 
     
    148143  
    149144      tra(:,:,:,:) = 0._wp 
    150       IF( ln_crs_top)  THEN 
    151          CALL zps_hde_crs( nit000, jptra, trn, gtru, gtrv ) 
    152       ELSE 
    153          IF( ln_zps .AND. .NOT. lk_c1d .AND. .NOT. ln_isfcav )   &              ! Partial steps: before horizontal gradient of passive 
    154          &    CALL zps_hde    ( nit000, jptra, trn, gtru, gtrv  )  ! Partial steps: before horizontal gradient 
    155          IF( ln_zps .AND. .NOT. lk_c1d .AND.       ln_isfcav )   & 
    156          &    CALL zps_hde_isf( nit000, jptra, trn, pgtu=gtru, pgtv=gtrv, pgtui=gtrui, pgtvi=gtrvi )       ! tracers at the bottom ocean level 
    157       ENDIF 
    158  
    159145      ! 
    160146      IF( nn_dttrc /= 1 )        CALL trc_sub_ini      ! Initialize variables for substepping passive tracers 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/trcnam.F90

    r5602 r7256  
    397397   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    398398   !!====================================================================== 
    399 END MODULE  trcnam 
     399END MODULE trcnam 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/trcrst.F90

    r5602 r7256  
    307307         IF(lwp) WRITE(numout,9000) jn, TRIM( ctrcnm(jn) ), zmean, zmin, zmax, zdrift 
    308308      END DO 
    309       WRITE(numout,*)  
     309      IF(lwp) WRITE(numout,*)  
    3103109000  FORMAT(' tracer nb :',i2,'    name :',a10,'    mean :',e18.10,'    min :',e18.10, & 
    311311      &      '    max :',e18.10,'    drift :',e18.10, ' %') 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/trcsms.F90

    r3680 r7256  
    7575 
    7676   !!====================================================================== 
    77 END MODULE  trcsms 
     77END MODULE trcsms 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/trcstp.F90

    r7222 r7256  
    2020   USE trdtrc_oce 
    2121   USE trdmxl_trc 
    22    USE iom, ONLY : lk_iomput , iom_close 
     22   USE iom, ONLY : lk_iomput , iom_close, iom_varid, jpdom_autoglo, iom_get, iom_rstput 
    2323   USE in_out_manager 
    2424   USE trcsub 
     
    3333   REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE ::   qsr_arr ! save qsr during TOP time-step 
    3434   REAL(wp) :: rdt_sampl 
    35    INTEGER  :: nb_rec_per_days 
    36    INTEGER  :: isecfst, iseclast 
     35   INTEGER  :: nb_rec_per_day 
     36   REAL(wp) :: rsecfst, rseclast 
    3737   LOGICAL  :: llnew 
    3838 
     
    6060      REAL(wp)              ::  ztrai 
    6161      CHARACTER (len=25)    ::  charout  
    62  
    6362      !!------------------------------------------------------------------- 
    6463      ! 
     
    9594                                   CALL trc_sms      ( kt )       ! tracers: sinks and sources 
    9695                                   CALL trc_trp      ( kt )       ! transport of passive tracers 
     96 
    9797         IF( kt == nittrc000 ) THEN 
    9898            CALL iom_close( numrtr )       ! close input tracer restart file 
     
    106106      ENDIF 
    107107      ! 
     108 
    108109      ztrai = 0._wp                                                   !  content of all tracers 
    109110      DO jn = 1, jptra 
     
    111112      END DO 
    112113      IF( lwp ) WRITE(numstr,9300) kt,  ztrai / areatot 
    113 9300  FORMAT(i10,e18.10) 
     1149300  FORMAT(i10,D23.16) 
    114115      ! 
    115116      IF( nn_timing == 1 )   CALL timing_stop('trc_stp') 
     
    124125      !!               of diurnal cycle 
    125126      !! 
    126       !! ** Method  : store in TOP the qsr every hour ( or every time-step the latter  
     127      !! ** Method  : store in TOP the qsr every hour ( or every time-step if the latter  
    127128      !!              is greater than 1 hour ) and then, compute the  mean with  
    128129      !!              a moving average over 24 hours.  
     
    131132      INTEGER, INTENT(in) ::   kt 
    132133      INTEGER  :: jn 
     134      REAL(wp) :: zkt 
     135      CHARACTER(len=1)               ::   cl1                      ! 1 character 
     136      CHARACTER(len=2)               ::   cl2                      ! 2 characters 
    133137 
    134138      IF( kt == nit000 ) THEN 
    135139         IF( ln_cpl )  THEN   
    136             rdt_sampl = 86400. / ncpl_qsr_freq 
    137             nb_rec_per_days = ncpl_qsr_freq 
     140            rdt_sampl = rday / ncpl_qsr_freq 
     141            nb_rec_per_day = ncpl_qsr_freq 
    138142         ELSE   
    139             rdt_sampl = MAX( 3600., rdt * nn_dttrc ) 
    140             nb_rec_per_days = INT( 86400 / rdt_sampl ) 
     143            rdt_sampl = MAX( 3600., rdttrc(1) ) 
     144            nb_rec_per_day = INT( rday / rdt_sampl ) 
    141145         ENDIF 
    142146         ! 
    143147         IF( lwp ) THEN 
    144148            WRITE(numout,*)  
    145             WRITE(numout,*) ' Sampling frequency dt = ', rdt_sampl, 's','   Number of sampling per day  nrec = ', nb_rec_per_days 
     149            WRITE(numout,*) ' Sampling frequency dt = ', rdt_sampl, 's','   Number of sampling per day  nrec = ', nb_rec_per_day 
    146150            WRITE(numout,*)  
    147151         ENDIF 
    148152         ! 
    149          ALLOCATE( qsr_arr(jpi,jpj,nb_rec_per_days ) ) 
    150          DO jn = 1, nb_rec_per_days 
    151             qsr_arr(:,:,jn) = qsr(:,:) 
    152          ENDDO 
    153          qsr_mean(:,:) = qsr(:,:) 
    154          ! 
    155          isecfst  = nsec_year + nsec1jan000   !   number of seconds between Jan. 1st 00h of nit000 year and the middle of time step 
    156          iseclast = isecfst 
    157          ! 
    158       ENDIF 
    159       ! 
    160       iseclast = nsec_year + nsec1jan000 
    161       llnew   = ( iseclast - isecfst )  > INT( rdt_sampl )   !   new shortwave to store 
    162       IF( kt /= nittrc000 .AND. llnew ) THEN 
     153         ALLOCATE( qsr_arr(jpi,jpj,nb_rec_per_day ) ) 
     154         ! 
     155         !                                            !* Restart: read in restart file 
     156         IF( ln_rsttr .AND. iom_varid( numrtr, 'qsr_mean' , ldstop = .FALSE. ) > 0 .AND. & 
     157                            iom_varid( numrtr, 'qsr_arr_1', ldstop = .FALSE. ) > 0 .AND. & 
     158                            iom_varid( numrtr, 'ktdcy'    , ldstop = .FALSE. ) > 0 ) THEN  
     159            CALL iom_get( numrtr, 'ktdcy', zkt )   !  A mean of qsr 
     160            rsecfst = INT( zkt ) * rdttrc(1) 
     161            IF(lwp) WRITE(numout,*) 'trc_qsr_mean:   qsr_mean read in the restart file at time-step rsecfst =', rsecfst, ' s ' 
     162            CALL iom_get( numrtr, jpdom_autoglo, 'qsr_mean', qsr_mean )   !  A mean of qsr 
     163            DO jn = 1, nb_rec_per_day  
     164             IF( jn <= 9 )  THEN 
     165               WRITE(cl1,'(i1)') jn 
     166               CALL iom_get( numrtr, jpdom_autoglo, 'qsr_arr_'//cl1, qsr_arr(:,:,jn) )   !  A mean of qsr 
     167             ELSE 
     168               WRITE(cl2,'(i2.2)') jn 
     169               CALL iom_get( numrtr, jpdom_autoglo, 'qsr_arr_'//cl2, qsr_arr(:,:,jn) )   !  A mean of qsr 
     170             ENDIF 
     171           ENDDO 
     172         ELSE                                         !* no restart: set from nit000 values 
     173            IF(lwp) WRITE(numout,*) 'trc_qsr_mean:   qsr_mean set to nit000 values' 
     174            rsecfst  = kt * rdttrc(1) 
     175            ! 
     176            qsr_mean(:,:) = qsr(:,:) 
     177            DO jn = 1, nb_rec_per_day 
     178               qsr_arr(:,:,jn) = qsr_mean(:,:) 
     179            ENDDO 
     180         ENDIF 
     181         ! 
     182      ENDIF 
     183      ! 
     184      rseclast = kt * rdttrc(1) 
     185      ! 
     186      llnew   = ( rseclast - rsecfst ) .ge.  rdt_sampl    !   new shortwave to store 
     187      IF( llnew ) THEN 
    163188          IF( lwp ) WRITE(numout,*) ' New shortwave to sample for TOP at time kt = ', kt, & 
    164              &                      ' time = ', (iseclast+rdt*nn_dttrc/2.)/3600.,'hours ' 
    165           isecfst = iseclast 
    166           DO jn = 1, nb_rec_per_days - 1 
     189             &                      ' time = ', rseclast/3600.,'hours ' 
     190          rsecfst = rseclast 
     191          DO jn = 1, nb_rec_per_day - 1 
    167192             qsr_arr(:,:,jn) = qsr_arr(:,:,jn+1) 
    168193          ENDDO 
    169           qsr_arr (:,:,nb_rec_per_days) = qsr(:,:) 
    170           qsr_mean(:,:                ) = SUM( qsr_arr(:,:,:), 3 ) / nb_rec_per_days 
     194          qsr_arr (:,:,nb_rec_per_day) = qsr(:,:) 
     195          qsr_mean(:,:                ) = SUM( qsr_arr(:,:,:), 3 ) / nb_rec_per_day 
     196      ENDIF 
     197      ! 
     198      IF( lrst_trc ) THEN    !* Write the mean of qsr in restart file  
     199         IF(lwp) WRITE(numout,*) 
     200         IF(lwp) WRITE(numout,*) 'trc_mean_qsr : write qsr_mean in restart file  kt =', kt 
     201         IF(lwp) WRITE(numout,*) '~~~~~~~' 
     202         zkt = REAL( kt, wp ) 
     203         CALL iom_rstput( kt, nitrst, numrtw, 'ktdcy', zkt ) 
     204          DO jn = 1, nb_rec_per_day  
     205             IF( jn <= 9 )  THEN 
     206               WRITE(cl1,'(i1)') jn 
     207               CALL iom_rstput( kt, nitrst, numrtw, 'qsr_arr_'//cl1, qsr_arr(:,:,jn) ) 
     208             ELSE 
     209               WRITE(cl2,'(i2.2)') jn 
     210               CALL iom_rstput( kt, nitrst, numrtw, 'qsr_arr_'//cl2, qsr_arr(:,:,jn) ) 
     211             ENDIF 
     212         ENDDO 
     213         CALL iom_rstput( kt, nitrst, numrtw, 'qsr_mean', qsr_mean(:,:) ) 
    171214      ENDIF 
    172215      ! 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/trcsub.F90

    r6772 r7256  
    1616   USE iom_def, ONLY : jprstlib 
    1717   USE lbclnk 
    18 !#if defined key_zdftke 
    19 !   USE zdftke          ! twice TKE (en) 
    20 !#endif 
    21 #if defined key_zdfgls 
    22    USE zdfgls, ONLY: en 
    23 #endif 
    24 !   USE trabbl 
    25 !   USE zdf_oce 
    26 !   USE domvvl 
    27    USE divcur, ONLY : div_cur           ! hor. divergence and curl      (div & cur routines) 
     18   !cbr USE trabbl 
     19   !cbr USE zdf_oce 
     20   !cbr USE domvvl 
     21   USE divcur, ONLY : div_cur        ! hor. divergence and curl      (div & cur routines) 
    2822   USE sbcrnf, ONLY: h_rnf, nk_rnf   ! River runoff  
    2923   USE bdy_oce 
Note: See TracChangeset for help on using the changeset viewer.