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 13463 for NEMO/branches/2019/dev_r11351_fldread_with_XIOS/src/TOP/PISCES/P4Z/p4zflx.F90 – NEMO

Ignore:
Timestamp:
2020-09-14T17:40:34+02:00 (4 years ago)
Author:
andmirek
Message:

Ticket #2195:update to trunk 13461

Location:
NEMO/branches/2019/dev_r11351_fldread_with_XIOS
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r11351_fldread_with_XIOS

    • Property svn:externals
      •  

        old new  
        33^/utils/build/mk@HEAD         mk 
        44^/utils/tools@HEAD            tools 
        5 ^/vendors/AGRIF/dev@HEAD      ext/AGRIF 
         5^/vendors/AGRIF/dev_r12970_AGRIF_CMEMS      ext/AGRIF 
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 
         8 
         9# SETTE 
         10^/utils/CI/sette@13382        sette 
  • NEMO/branches/2019/dev_r11351_fldread_with_XIOS/src/TOP/PISCES/P4Z/p4zflx.F90

    r10425 r13463  
    1919   USE sms_pisces     !  PISCES Source Minus Sink variables 
    2020   USE p4zche         !  Chemical model 
    21    USE prtctl_trc     !  print control for debugging 
     21   USE prtctl         !  print control for debugging 
    2222   USE iom            !  I/O manager 
    2323   USE fldread        !  read input fields 
     
    5252   REAL(wp) ::   xconv  = 0.01_wp / 3600._wp   !: coefficients for conversion  
    5353 
     54   !! * Substitutions 
     55#  include "do_loop_substitute.h90" 
     56#  include "domzgr_substitute.h90" 
    5457   !!---------------------------------------------------------------------- 
    5558   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
     
    5962CONTAINS 
    6063 
    61    SUBROUTINE p4z_flx ( kt, knt ) 
     64   SUBROUTINE p4z_flx ( kt, knt, Kbb, Kmm, Krhs ) 
    6265      !!--------------------------------------------------------------------- 
    6366      !!                     ***  ROUTINE p4z_flx  *** 
     
    7174      !!--------------------------------------------------------------------- 
    7275      INTEGER, INTENT(in) ::   kt, knt   ! 
     76      INTEGER, INTENT(in) ::   Kbb, Kmm, Krhs      ! time level indices 
    7377      ! 
    7478      INTEGER  ::   ji, jj, jm, iind, iindm1 
     
    8084      CHARACTER (len=25) ::   charout 
    8185      REAL(wp), DIMENSION(jpi,jpj) ::   zkgco2, zkgo2, zh2co3, zoflx,  zpco2atm   
    82       REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zw2d 
    8386      !!--------------------------------------------------------------------- 
    8487      ! 
     
    107110      IF( l_co2cpl )   satmco2(:,:) = atm_co2(:,:) 
    108111 
    109       DO jj = 1, jpj 
    110          DO ji = 1, jpi 
    111             ! DUMMY VARIABLES FOR DIC, H+, AND BORATE 
    112             zfact = rhop(ji,jj,1) / 1000. + rtrn 
    113             zdic  = trb(ji,jj,1,jpdic) 
    114             zph   = MAX( hi(ji,jj,1), 1.e-10 ) / zfact 
    115             ! CALCULATE [H2CO3] 
    116             zh2co3(ji,jj) = zdic/(1. + ak13(ji,jj,1)/zph + ak13(ji,jj,1)*ak23(ji,jj,1)/zph**2) 
    117          END DO 
    118       END DO 
     112      DO_2D( 1, 1, 1, 1 ) 
     113         ! DUMMY VARIABLES FOR DIC, H+, AND BORATE 
     114         zfact = rhop(ji,jj,1) / 1000. + rtrn 
     115         zdic  = tr(ji,jj,1,jpdic,Kbb) 
     116         zph   = MAX( hi(ji,jj,1), 1.e-10 ) / zfact 
     117         ! CALCULATE [H2CO3] 
     118         zh2co3(ji,jj) = zdic/(1. + ak13(ji,jj,1)/zph + ak13(ji,jj,1)*ak23(ji,jj,1)/zph**2) 
     119      END_2D 
    119120 
    120121      ! -------------- 
     
    125126      ! ------------------------------------------- 
    126127 
    127       DO jj = 1, jpj 
    128          DO ji = 1, jpi 
    129             ztc  = MIN( 35., tsn(ji,jj,1,jp_tem) ) 
    130             ztc2 = ztc * ztc 
    131             ztc3 = ztc * ztc2  
    132             ztc4 = ztc2 * ztc2  
    133             ! Compute the schmidt Number both O2 and CO2 
    134             zsch_co2 = 2116.8 - 136.25 * ztc + 4.7353 * ztc2 - 0.092307 * ztc3 + 0.0007555 * ztc4 
    135             zsch_o2  = 1920.4 - 135.6  * ztc + 5.2122 * ztc2 - 0.109390 * ztc3 + 0.0009377 * ztc4 
    136             !  wind speed  
    137             zws  = wndm(ji,jj) * wndm(ji,jj) 
    138             ! Compute the piston velocity for O2 and CO2 
    139             zkgwan = 0.251 * zws 
    140             zkgwan = zkgwan * xconv * ( 1.- fr_i(ji,jj) ) * tmask(ji,jj,1) 
    141             ! compute gas exchange for CO2 and O2 
    142             zkgco2(ji,jj) = zkgwan * SQRT( 660./ zsch_co2 ) 
    143             zkgo2 (ji,jj) = zkgwan * SQRT( 660./ zsch_o2 ) 
    144          END DO 
    145       END DO 
    146  
    147  
    148       DO jj = 1, jpj 
    149          DO ji = 1, jpi 
    150             ztkel = tempis(ji,jj,1) + 273.15 
    151             zsal  = salinprac(ji,jj,1) + ( 1.- tmask(ji,jj,1) ) * 35. 
    152             zvapsw    = EXP(24.4543 - 67.4509*(100.0/ztkel) - 4.8489*LOG(ztkel/100) - 0.000544*zsal) 
    153             zpco2atm(ji,jj) = satmco2(ji,jj) * ( patm(ji,jj) - zvapsw ) 
    154             zxc2      = ( 1.0 - zpco2atm(ji,jj) * 1E-6 )**2 
    155             zfugcoeff = EXP( patm(ji,jj) * (chemc(ji,jj,2) + 2.0 * zxc2 * chemc(ji,jj,3) )   & 
    156             &           / ( 82.05736 * ztkel )) 
    157             zfco2 = zpco2atm(ji,jj) * zfugcoeff 
    158  
    159             ! Compute CO2 flux for the sea and air 
    160             zfld = zfco2 * chemc(ji,jj,1) * zkgco2(ji,jj)  ! (mol/L) * (m/s) 
    161             zflu = zh2co3(ji,jj) * zkgco2(ji,jj)                                   ! (mol/L) (m/s) ? 
    162             oce_co2(ji,jj) = ( zfld - zflu ) * rfact2 * e1e2t(ji,jj) * tmask(ji,jj,1) * 1000. 
    163             ! compute the trend 
    164             tra(ji,jj,1,jpdic) = tra(ji,jj,1,jpdic) + ( zfld - zflu ) * rfact2 / e3t_n(ji,jj,1) * tmask(ji,jj,1) 
    165  
    166             ! Compute O2 flux  
    167             zfld16 = patm(ji,jj) * chemo2(ji,jj,1) * zkgo2(ji,jj)          ! (mol/L) * (m/s) 
    168             zflu16 = trb(ji,jj,1,jpoxy) * zkgo2(ji,jj) 
    169             zoflx(ji,jj) = ( zfld16 - zflu16 ) * tmask(ji,jj,1) 
    170             tra(ji,jj,1,jpoxy) = tra(ji,jj,1,jpoxy) + zoflx(ji,jj) * rfact2 / e3t_n(ji,jj,1) 
    171          END DO 
    172       END DO 
     128      DO_2D( 1, 1, 1, 1 ) 
     129         ztc  = MIN( 35., ts(ji,jj,1,jp_tem,Kmm) ) 
     130         ztc2 = ztc * ztc 
     131         ztc3 = ztc * ztc2  
     132         ztc4 = ztc2 * ztc2  
     133         ! Compute the schmidt Number both O2 and CO2 
     134         zsch_co2 = 2116.8 - 136.25 * ztc + 4.7353 * ztc2 - 0.092307 * ztc3 + 0.0007555 * ztc4 
     135         zsch_o2  = 1920.4 - 135.6  * ztc + 5.2122 * ztc2 - 0.109390 * ztc3 + 0.0009377 * ztc4 
     136         !  wind speed  
     137         zws  = wndm(ji,jj) * wndm(ji,jj) 
     138         ! Compute the piston velocity for O2 and CO2 
     139         zkgwan = 0.251 * zws 
     140         zkgwan = zkgwan * xconv * ( 1.- fr_i(ji,jj) ) * tmask(ji,jj,1) 
     141         ! compute gas exchange for CO2 and O2 
     142         zkgco2(ji,jj) = zkgwan * SQRT( 660./ zsch_co2 ) 
     143         zkgo2 (ji,jj) = zkgwan * SQRT( 660./ zsch_o2 ) 
     144      END_2D 
     145 
     146 
     147      DO_2D( 1, 1, 1, 1 ) 
     148         ztkel = tempis(ji,jj,1) + 273.15 
     149         zsal  = salinprac(ji,jj,1) + ( 1.- tmask(ji,jj,1) ) * 35. 
     150         zvapsw    = EXP(24.4543 - 67.4509*(100.0/ztkel) - 4.8489*LOG(ztkel/100) - 0.000544*zsal) 
     151         zpco2atm(ji,jj) = satmco2(ji,jj) * ( patm(ji,jj) - zvapsw ) 
     152         zxc2      = ( 1.0 - zpco2atm(ji,jj) * 1E-6 )**2 
     153         zfugcoeff = EXP( patm(ji,jj) * (chemc(ji,jj,2) + 2.0 * zxc2 * chemc(ji,jj,3) )   & 
     154         &           / ( 82.05736 * ztkel )) 
     155         zfco2 = zpco2atm(ji,jj) * zfugcoeff 
     156 
     157         ! Compute CO2 flux for the sea and air 
     158         zfld = zfco2 * chemc(ji,jj,1) * zkgco2(ji,jj)  ! (mol/L) * (m/s) 
     159         zflu = zh2co3(ji,jj) * zkgco2(ji,jj)                                   ! (mol/L) (m/s) ? 
     160         oce_co2(ji,jj) = ( zfld - zflu ) * tmask(ji,jj,1)  
     161         ! compute the trend 
     162         tr(ji,jj,1,jpdic,Krhs) = tr(ji,jj,1,jpdic,Krhs) + oce_co2(ji,jj) * rfact2 / e3t(ji,jj,1,Kmm) 
     163 
     164         ! Compute O2 flux  
     165         zfld16 = patm(ji,jj) * chemo2(ji,jj,1) * zkgo2(ji,jj)          ! (mol/L) * (m/s) 
     166         zflu16 = tr(ji,jj,1,jpoxy,Kbb) * zkgo2(ji,jj) 
     167         zoflx(ji,jj) = ( zfld16 - zflu16 ) * tmask(ji,jj,1) 
     168         tr(ji,jj,1,jpoxy,Krhs) = tr(ji,jj,1,jpoxy,Krhs) + zoflx(ji,jj) * rfact2 / e3t(ji,jj,1,Kmm) 
     169      END_2D 
    173170 
    174171      IF( iom_use("tcflx") .OR. iom_use("tcflxcum") .OR. kt == nitrst   & 
    175172         &                 .OR. (ln_check_mass .AND. kt == nitend) )    & 
    176          t_oce_co2_flx  = glob_sum( 'p4zflx', oce_co2(:,:) )                    !  Total Flux of Carbon 
     173         t_oce_co2_flx  = glob_sum( 'p4zflx', oce_co2(:,:) * e1e2t(:,:) * 1000. )                    !  Total Flux of Carbon 
    177174      t_oce_co2_flx_cum = t_oce_co2_flx_cum + t_oce_co2_flx       !  Cumulative Total Flux of Carbon 
    178175!      t_atm_co2_flx     = glob_sum( 'p4zflx', satmco2(:,:) * e1e2t(:,:) )       ! Total atmospheric pCO2 
    179176      t_atm_co2_flx     =  atcco2      ! Total atmospheric pCO2 
    180177  
    181       IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
     178      IF(sn_cfctl%l_prttrc)   THEN  ! print mean trends (used for debugging) 
    182179         WRITE(charout, FMT="('flx ')") 
    183          CALL prt_ctl_trc_info(charout) 
    184          CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm) 
     180         CALL prt_ctl_info( charout, cdcomp = 'top' ) 
     181         CALL prt_ctl(tab4d_1=tr(:,:,:,:,Krhs), mask1=tmask, clinfo=ctrcnm) 
    185182      ENDIF 
    186183 
    187184      IF( lk_iomput .AND. knt == nrdttrc ) THEN 
    188          ALLOCATE( zw2d(jpi,jpj) )   
    189          IF( iom_use( "Cflx"  ) )  THEN 
    190             zw2d(:,:) = oce_co2(:,:) / e1e2t(:,:) * rfact2r 
    191             CALL iom_put( "Cflx"     , zw2d )  
    192          ENDIF 
    193          IF( iom_use( "Oflx"  ) )  THEN 
    194             zw2d(:,:) =  zoflx(:,:) * 1000 * tmask(:,:,1) 
    195             CALL iom_put( "Oflx" , zw2d ) 
    196          ENDIF 
    197          IF( iom_use( "Kg"    ) )  THEN 
    198             zw2d(:,:) =  zkgco2(:,:) * tmask(:,:,1) 
    199             CALL iom_put( "Kg"   , zw2d ) 
    200          ENDIF 
    201          IF( iom_use( "Dpco2" ) ) THEN 
    202            zw2d(:,:) = ( zpco2atm(:,:) - zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) 
    203            CALL iom_put( "Dpco2" ,  zw2d ) 
    204          ENDIF 
    205          IF( iom_use( "Dpo2" ) )  THEN 
    206            zw2d(:,:) = ( atcox * patm(:,:) - atcox * trb(:,:,1,jpoxy) / ( chemo2(:,:,1) + rtrn ) ) * tmask(:,:,1) 
    207            CALL iom_put( "Dpo2"  , zw2d ) 
    208          ENDIF 
    209          CALL iom_put( "tcflx"    , t_oce_co2_flx * rfact2r )   ! molC/s 
    210          CALL iom_put( "tcflxcum" , t_oce_co2_flx_cum       )   ! molC 
    211          ! 
    212          DEALLOCATE( zw2d ) 
     185         CALL iom_put( "AtmCo2"  , satmco2(:,:) * tmask(:,:,1) )   ! Atmospheric CO2 concentration 
     186         CALL iom_put( "Cflx"    , oce_co2(:,:) * 1000. )  
     187         CALL iom_put( "Oflx"    , zoflx(:,:) * 1000.  ) 
     188         CALL iom_put( "Kg"      , zkgco2(:,:) * tmask(:,:,1)  ) 
     189         CALL iom_put( "Dpco2"   , ( zpco2atm(:,:) - zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) ) 
     190         CALL iom_put( "pCO2sea" , ( zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) ) 
     191         CALL iom_put( "Dpo2"    , ( atcox * patm(:,:) - atcox * tr(:,:,1,jpoxy,Kbb) / ( chemo2(:,:,1) + rtrn ) ) * tmask(:,:,1) ) 
     192         CALL iom_put( "tcflx"   , t_oce_co2_flx     )   ! molC/s 
     193         CALL iom_put( "tcflxcum", t_oce_co2_flx_cum )   ! molC 
    213194      ENDIF 
    214195      ! 
     
    239220      ENDIF 
    240221      ! 
    241       REWIND( numnatp_ref )              ! Namelist nampisext in reference namelist : Pisces atm. conditions 
    242222      READ  ( numnatp_ref, nampisext, IOSTAT = ios, ERR = 901) 
    243 901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nampisext in reference namelist', lwp ) 
    244       REWIND( numnatp_cfg )              ! Namelist nampisext in configuration namelist : Pisces atm. conditions 
     223901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nampisext in reference namelist' ) 
    245224      READ  ( numnatp_cfg, nampisext, IOSTAT = ios, ERR = 902 ) 
    246 902   IF( ios >  0 )   CALL ctl_nam ( ios , 'nampisext in configuration namelist', lwp ) 
     225902   IF( ios >  0 )   CALL ctl_nam ( ios , 'nampisext in configuration namelist' ) 
    247226      IF(lwm) WRITE ( numonp, nampisext ) 
    248227      ! 
     
    320299         ENDIF 
    321300         ! 
    322          REWIND( numnatp_ref )              ! Namelist nampisatm in reference namelist : Pisces atm. sea level pressure file 
    323301         READ  ( numnatp_ref, nampisatm, IOSTAT = ios, ERR = 901) 
    324 901      IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisatm in reference namelist', lwp ) 
    325          REWIND( numnatp_cfg )              ! Namelist nampisatm in configuration namelist : Pisces atm. sea level pressure file  
     302901      IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisatm in reference namelist' ) 
    326303         READ  ( numnatp_cfg, nampisatm, IOSTAT = ios, ERR = 902 ) 
    327 902      IF( ios >  0 )   CALL ctl_nam ( ios , 'nampisatm in configuration namelist', lwp ) 
     304902      IF( ios >  0 )   CALL ctl_nam ( ios , 'nampisatm in configuration namelist' ) 
    328305         IF(lwm) WRITE ( numonp, nampisatm ) 
    329306         ! 
Note: See TracChangeset for help on using the changeset viewer.