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 3294 for trunk/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zlys.F90 – NEMO

Ignore:
Timestamp:
2012-01-28T17:44:18+01:00 (12 years ago)
Author:
rblod
Message:

Merge of 3.4beta into the trunk

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zlys.F90

    r2715 r3294  
    99   !!             1.0  !  2004     (O. Aumont) modifications 
    1010   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90 
     11   !!                  !  2011-02  (J. Simeon, J. Orr)  Calcon salinity dependence 
     12   !!             3.4  !  2011-06  (O. Aumont, C. Ethe) Improvment of calcite dissolution 
    1113   !!---------------------------------------------------------------------- 
    1214#if defined key_pisces 
     
    1719   !!   p4z_lys_init   :   Read the namelist parameters 
    1820   !!---------------------------------------------------------------------- 
    19    USE trc 
    20    USE oce_trc         ! 
    21    USE trc 
    22    USE sms_pisces 
    23    USE prtctl_trc 
    24    USE iom 
     21   USE oce_trc         !  shared variables between ocean and passive tracers 
     22   USE trc             !  passive tracers common variables  
     23   USE sms_pisces      !  PISCES Source Minus Sink variables 
     24   USE prtctl_trc      !  print control for debugging 
     25   USE iom             !  I/O manager 
    2526 
    2627   IMPLICIT NONE 
     
    5758      !! ** Method  : - ??? 
    5859      !!--------------------------------------------------------------------- 
    59       USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 
    60       USE wrk_nemo, ONLY: zco3 => wrk_3d_2, zcaldiss => wrk_3d_3  
    6160      ! 
    6261      INTEGER, INTENT(in) ::   kt ! ocean time step 
    6362      INTEGER  ::   ji, jj, jk, jn 
    64       REAL(wp) ::   zbot, zalk, zdic, zph, zremco3, zah2 
    65       REAL(wp) ::   zdispot, zfact, zalka 
     63      REAL(wp) ::   zalk, zdic, zph, zah2 
     64      REAL(wp) ::   zdispot, zfact, zcalcon, zalka, zaldi 
    6665      REAL(wp) ::   zomegaca, zexcess, zexcess0 
    67 #if defined key_diatrc && defined key_iomput 
    6866      REAL(wp) ::   zrfact2 
    69 #endif 
    7067      CHARACTER (len=25) :: charout 
     68      REAL(wp), POINTER, DIMENSION(:,:,:) :: zco3, zcaldiss    
    7169      !!--------------------------------------------------------------------- 
    72  
    73       IF(  wrk_in_use(3, 2,3) ) THEN 
    74          CALL ctl_stop('p4z_lys: requested workspace arrays unavailable')  ;  RETURN 
    75       END IF 
    76  
    77       zco3(:,:,:) = 0. 
    78 # if defined key_diatrc && defined key_iomput 
     70      ! 
     71      IF( nn_timing == 1 )  CALL timing_start('p4z_lys') 
     72      ! 
     73      CALL wrk_alloc( jpi, jpj, jpk, zco3, zcaldiss ) 
     74      ! 
     75      zco3    (:,:,:) = 0. 
    7976      zcaldiss(:,:,:) = 0. 
    80 # endif 
    8177      !     ------------------------------------------- 
    8278      !     COMPUTE [CO3--] and [H+] CONCENTRATIONS 
     
    9187!CDIR NOVERRCHK 
    9288               DO ji = 1, jpi 
    93  
    94                   ! SET DUMMY VARIABLE FOR TOTAL BORATE 
    95                   zbot  = borat(ji,jj,jk) 
    96  
    97                   ! SET DUMMY VARIABLE FOR TOTAL BORATE 
    98                   zbot  = borat(ji,jj,jk) 
    99                   zfact = rhop (ji,jj,jk) / 1000. + rtrn 
    100  
    101                   ! SET DUMMY VARIABLE FOR [H+] 
    102                   zph   = hi(ji,jj,jk) * tmask(ji,jj,jk) / zfact + ( 1.-tmask(ji,jj,jk) ) * 1.e-9 
    103  
    104                   ! SET DUMMY VARIABLE FOR [SUM(CO2)]GIVEN  
     89                  zfact = rhop(ji,jj,jk) / 1000. + rtrn 
     90                  zph  = hi(ji,jj,jk) * tmask(ji,jj,jk) / zfact + ( 1.-tmask(ji,jj,jk) ) * 1.e-9 ! [H+] 
    10591                  zdic  = trn(ji,jj,jk,jpdic) / zfact 
    10692                  zalka = trn(ji,jj,jk,jptal) / zfact 
    107  
    10893                  ! CALCULATE [ALK]([CO3--], [HCO3-]) 
    109                   zalk  = zalka - (  akw3(ji,jj,jk) / zph - zph   & 
    110                      &             + zbot / (1.+ zph / akb3(ji,jj,jk) )  ) 
    111  
     94                  zalk  = zalka - ( akw3(ji,jj,jk) / zph - zph + borat(ji,jj,jk) / ( 1. + zph / akb3(ji,jj,jk) ) ) 
    11295                  ! CALCULATE [H+] and [CO3--] 
    113                   zah2 = SQRT( (zdic-zalk)*(zdic-zalk)+   & 
    114                      &     4.*(zalk*ak23(ji,jj,jk)/ak13(ji,jj,jk))   & 
    115                      &     *(2*zdic-zalk)) 
    116  
    117                   zah2=0.5*ak13(ji,jj,jk)/zalk*((zdic-zalk)+zah2) 
    118                   zco3(ji,jj,jk) = zalk/(2.+zah2/ak23(ji,jj,jk))*zfact 
    119  
    120                   hi(ji,jj,jk)  = zah2*zfact 
    121  
     96                  zaldi = zdic - zalk 
     97                  zah2  = SQRT( zaldi * zaldi + 4.* ( zalk * ak23(ji,jj,jk) / ak13(ji,jj,jk) ) * ( zdic + zaldi ) ) 
     98                  zah2  = 0.5 * ak13(ji,jj,jk) / zalk * ( zaldi + zah2 ) 
     99                  ! 
     100                  zco3(ji,jj,jk) = zalk / ( 2. + zah2 / ak23(ji,jj,jk) ) * zfact 
     101                  hi(ji,jj,jk)   = zah2 * zfact 
    122102               END DO 
    123103            END DO 
     
    137117 
    138118               ! DEVIATION OF [CO3--] FROM SATURATION VALUE 
    139                zomegaca = ( calcon * zco3(ji,jj,jk) ) / aksp(ji,jj,jk) 
     119               ! Salinity dependance in zomegaca and divide by rhop/1000 to have good units 
     120               zcalcon  = calcon * ( tsn(ji,jj,jk,jp_sal) / 35._wp ) 
     121               zfact    = rhop(ji,jj,jk) / 1000._wp 
     122               zomegaca = ( zcalcon * zco3(ji,jj,jk) * zfact ) / aksp(ji,jj,jk)  
    140123 
    141124               ! SET DEGREE OF UNDER-/SUPERSATURATION 
    142                zexcess0 = MAX( 0., ( 1.- zomegaca ) ) 
     125               excess(ji,jj,jk) = 1._wp - zomegaca 
     126               zexcess0 = MAX( 0., excess(ji,jj,jk) ) 
    143127               zexcess  = zexcess0**nca 
    144128 
     
    146130               !       (ACCORDING TO THIS FORMULATION ALSO SOME PARTICULATE 
    147131               !       CACO3 GETS DISSOLVED EVEN IN THE CASE OF OVERSATURATION) 
     132               zdispot = kdca * zexcess * trn(ji,jj,jk,jpcal) 
    148133# if defined key_degrad 
    149               zdispot = kdca * zexcess * trn(ji,jj,jk,jpcal) * facvol(ji,jj,jk) 
    150 # else 
    151               zdispot = kdca * zexcess * trn(ji,jj,jk,jpcal) 
     134               zdispot = zdispot * facvol(ji,jj,jk) 
    152135# endif 
    153  
    154136              !  CHANGE OF [CO3--] , [ALK], PARTICULATE [CACO3], 
    155137              !       AND [SUM(CO2)] DUE TO CACO3 DISSOLUTION/PRECIPITATION 
    156               zremco3 = zdispot / rmtss 
    157               zco3(ji,jj,jk) = zco3(ji,jj,jk) + zremco3 * rfact 
    158               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + 2. * zremco3 
    159               tra(ji,jj,jk,jpcal) = tra(ji,jj,jk,jpcal) -      zremco3 
    160               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) +      zremco3 
    161  
    162 # if defined key_diatrc && defined key_iomput 
    163               zcaldiss(ji,jj,jk) = zremco3  ! calcite dissolution 
    164 # endif 
     138              zcaldiss(ji,jj,jk)  = zdispot / rmtss  ! calcite dissolution 
     139              zco3(ji,jj,jk)      = zco3(ji,jj,jk) + zcaldiss(ji,jj,jk) * rfact 
     140              ! 
     141              tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + 2. * zcaldiss(ji,jj,jk) 
     142              tra(ji,jj,jk,jpcal) = tra(ji,jj,jk,jpcal) -      zcaldiss(ji,jj,jk) 
     143              tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) +      zcaldiss(ji,jj,jk) 
    165144            END DO 
    166145         END DO 
    167146      END DO 
    168  
    169 # if defined key_diatrc 
    170 #  if ! defined key_iomput 
    171       trc3d(:,:,:,jp_pcs0_3d    ) = hi  (:,:,:)          * tmask(:,:,:) 
    172       trc3d(:,:,:,jp_pcs0_3d + 1) = zco3(:,:,:)          * tmask(:,:,:) 
    173       trc3d(:,:,:,jp_pcs0_3d + 2) = aksp(:,:,:) / calcon * tmask(:,:,:) 
    174 #  else 
    175       zrfact2 = 1.e3 * rfact2r 
    176       CALL iom_put( "PH"    , hi      (:,:,:)           * tmask(:,:,:) ) 
    177       CALL iom_put( "CO3"   , zco3    (:,:,:)           * tmask(:,:,:) ) 
    178       CALL iom_put( "CO3sat", aksp    (:,:,:) / calcon  * tmask(:,:,:) ) 
    179       CALL iom_put( "DCAL"  , zcaldiss(:,:,:) * zrfact2 * tmask(:,:,:) ) 
    180 #  endif 
    181 # endif 
    182       ! 
    183        IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
    184          WRITE(charout, FMT="('lys ')") 
    185          CALL prt_ctl_trc_info(charout) 
    186          CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm) 
    187        ENDIF 
    188  
    189       IF( wrk_not_released(3, 2,3) ) CALL ctl_stop('p4z_lys: failed to release workspace arrays') 
     147      ! 
     148      IF( ln_diatrc )  THEN 
     149         ! 
     150         IF( lk_iomput ) THEN 
     151            zrfact2 = 1.e3 * rfact2r 
     152            CALL iom_put( "PH"    , hi      (:,:,:)           * tmask(:,:,:) ) 
     153            CALL iom_put( "CO3"   , zco3    (:,:,:)           * tmask(:,:,:) ) 
     154            CALL iom_put( "CO3sat", aksp    (:,:,:) / calcon  * tmask(:,:,:) ) 
     155            CALL iom_put( "DCAL"  , zcaldiss(:,:,:) * zrfact2 * tmask(:,:,:) ) 
     156         ELSE 
     157            trc3d(:,:,:,jp_pcs0_3d    ) = hi  (:,:,:)          * tmask(:,:,:) 
     158            trc3d(:,:,:,jp_pcs0_3d + 1) = zco3(:,:,:)          * tmask(:,:,:) 
     159            trc3d(:,:,:,jp_pcs0_3d + 2) = aksp(:,:,:) / calcon * tmask(:,:,:) 
     160         ENDIF 
     161         ! 
     162      ENDIF 
     163      ! 
     164      IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
     165        WRITE(charout, FMT="('lys ')") 
     166        CALL prt_ctl_trc_info(charout) 
     167        CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm) 
     168      ENDIF 
     169      ! 
     170      CALL wrk_dealloc( jpi, jpj, jpk, zco3, zcaldiss ) 
     171      ! 
     172      IF( nn_timing == 1 )  CALL timing_stop('p4z_lys') 
    190173      ! 
    191174   END SUBROUTINE p4z_lys 
     
    199182      !! 
    200183      !! ** Method  :   Read the nampiscal namelist and check the parameters 
    201       !!      called at the first timestep (nit000) 
     184      !!      called at the first timestep (nittrc000) 
    202185      !! 
    203186      !! ** input   :   Namelist nampiscal 
     
    207190      NAMELIST/nampiscal/ kdca, nca 
    208191 
    209       REWIND( numnat )                     ! read numnat 
    210       READ  ( numnat, nampiscal ) 
     192      REWIND( numnatp )                     ! read numnatp 
     193      READ  ( numnatp, nampiscal ) 
    211194 
    212195      IF(lwp) THEN                         ! control print 
Note: See TracChangeset for help on using the changeset viewer.