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

Ignore:
Timestamp:
2011-03-30T17:58:35+02:00 (13 years ago)
Author:
rblod
Message:

First attempt to put dynamic allocation on the trunk

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zsink.F90

    r2528 r2715  
    2121   PUBLIC   p4z_sink         ! called in p4zbio.F90 
    2222   PUBLIC   p4z_sink_init    ! called in trcsms_pisces.F90 
    23  
    24    !! * Shared module variables 
    25    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   &   !: 
    26      wsbio3, wsbio4,      &    !: POC and GOC sinking speeds 
    27      wscal                     !: Calcite and BSi sinking speeds 
    28  
    29    !! * Module variables 
    30    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   &   !: 
    31      sinking, sinking2,   &    !: POC sinking fluxes (different meanings depending on the parameterization 
    32      sinkcal, sinksil,    &    !: CaCO3 and BSi sinking fluxes 
    33      sinkfer                   !: Small BFe sinking flux 
    34  
    35    INTEGER  :: & 
    36       iksed  = 10              ! 
     23   PUBLIC   p4z_sink_alloc 
     24 
     25   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wsbio3   !: POC sinking speed  
     26   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wsbio4   !: GOC sinking speed 
     27   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wscal    !: Calcite and BSi sinking speeds 
     28 
     29   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sinking, sinking2  !: POC sinking fluxes  
     30   !                                                          !  (different meanings depending on the parameterization) 
     31   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sinkcal, sinksil   !: CaCO3 and BSi sinking fluxes 
     32   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sinkfer            !: Small BFe sinking fluxes 
     33#if ! defined key_kriest 
     34   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sinkfer2           !: Big iron sinking fluxes 
     35#endif 
     36 
     37   INTEGER  :: iksed  = 10 
    3738 
    3839#if  defined key_kriest 
    39    REAL(wp)          ::       &    
    40       xkr_sfact    = 250.  ,  &   !: Sinking factor 
    41       xkr_stick    = 0.2   ,  &   !: Stickiness 
    42       xkr_nnano    = 2.337 ,  &   !: Nbr of cell in nano size class 
    43       xkr_ndiat    = 3.718 ,  &   !: Nbr of cell in diatoms size class 
    44       xkr_nmeso    = 7.147 ,  &   !: Nbr of cell in mesozoo  size class 
    45       xkr_naggr    = 9.877        !: Nbr of cell in aggregates  size class 
    46  
    47    REAL(wp)          ::       &    
    48       xkr_frac 
    49  
    50    REAL(wp), PUBLIC ::        & 
    51       xkr_dnano            ,  &   !: Size of particles in nano pool 
    52       xkr_ddiat            ,  &   !: Size of particles in diatoms pool 
    53       xkr_dmeso            ,  &   !: Size of particles in mesozoo pool 
    54       xkr_daggr            ,  &   !: Size of particles in aggregates pool 
    55       xkr_wsbio_min        ,  &   !: min vertical particle speed 
    56       xkr_wsbio_max               !: max vertical particle speed 
    57  
    58    REAL(wp), PUBLIC, DIMENSION(jpk) ::   &   !: 
    59       xnumm                       !:     maximum number of particles in aggregates 
    60  
    61 #endif 
    62  
    63 #if ! defined key_kriest 
    64    REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &   !: 
    65      sinkfer2                  !: Big Fe sinking flux 
    66 #endif  
     40   REAL(wp) ::  xkr_sfact    = 250.     !: Sinking factor 
     41   REAL(wp) ::  xkr_stick    = 0.2      !: Stickiness 
     42   REAL(wp) ::  xkr_nnano    = 2.337    !: Nbr of cell in nano size class 
     43   REAL(wp) ::  xkr_ndiat    = 3.718    !: Nbr of cell in diatoms size class 
     44   REAL(wp) ::  xkr_nmeso    = 7.147    !: Nbr of cell in mesozoo  size class 
     45   REAL(wp) ::  xkr_naggr    = 9.877    !: Nbr of cell in aggregates  size class 
     46 
     47   REAL(wp) ::  xkr_frac  
     48 
     49   REAL(wp), PUBLIC ::  xkr_dnano       !: Size of particles in nano pool 
     50   REAL(wp), PUBLIC ::  xkr_ddiat       !: Size of particles in diatoms pool 
     51   REAL(wp), PUBLIC ::  xkr_dmeso       !: Size of particles in mesozoo pool 
     52   REAL(wp), PUBLIC ::  xkr_daggr       !: Size of particles in aggregates pool 
     53   REAL(wp), PUBLIC ::  xkr_wsbio_min   !: min vertical particle speed 
     54   REAL(wp), PUBLIC ::  xkr_wsbio_max   !: max vertical particle speed 
     55 
     56   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   xnumm   !:  maximum number of particles in aggregates 
     57#endif 
    6758 
    6859   !!* Substitution 
    … …  
    7162   !! NEMO/TOP 3.3 , NEMO Consortium (2010) 
    7263   !! $Id$  
    73    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
     64   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    7465   !!---------------------------------------------------------------------- 
    75  
    7666CONTAINS 
    7767 
    7868#if defined key_kriest 
     69   !!---------------------------------------------------------------------- 
     70   !!   'key_kriest'                                                    ??? 
     71   !!---------------------------------------------------------------------- 
    7972 
    8073   SUBROUTINE p4z_sink ( kt, jnt ) 
    … …  
    8780      !! ** Method  : - ??? 
    8881      !!--------------------------------------------------------------------- 
    89  
     82      USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
     83      USE wrk_nemo, ONLY:   znum3d => wrk_3d_2 
     84      ! 
    9085      INTEGER, INTENT(in) :: kt, jnt 
     86      ! 
    9187      INTEGER  :: ji, jj, jk 
    9288      REAL(wp) :: zagg1, zagg2, zagg3, zagg4, zagg5, zaggsi, zaggsh 
    … …  
    9995      INTEGER  :: ik1 
    10096#endif 
    101       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   znum3d 
    10297      CHARACTER (len=25) :: charout 
    103  
    104       !!--------------------------------------------------------------------- 
    105  
     98      !!--------------------------------------------------------------------- 
     99      ! 
     100      IF( wrk_in_use(3, 2 ) ) THEN 
     101         CALL ctl_stop('p4z_sink: requested workspace arrays unavailable')   ;   RETURN 
     102      ENDIF 
     103       
    106104      !     Initialisation of variables used to compute Sinking Speed 
    107105      !     --------------------------------------------------------- 
    108106 
    109        znum3d(:,:,:) = 0.e0 
    110        zval1 = 1. + xkr_zeta 
    111        zval2 = 1. + xkr_zeta + xkr_eta 
    112        zval3 = 1. + xkr_eta 
    113  
    114      !     Computation of the vertical sinking speed : Kriest et Evans, 2000 
    115      !     ----------------------------------------------------------------- 
     107      znum3d(:,:,:) = 0.e0 
     108      zval1 = 1. + xkr_zeta 
     109      zval2 = 1. + xkr_zeta + xkr_eta 
     110      zval3 = 1. + xkr_eta 
     111 
     112      !     Computation of the vertical sinking speed : Kriest et Evans, 2000 
     113      !     ----------------------------------------------------------------- 
    116114 
    117115      DO jk = 1, jpkm1 
    … …  
    131129                  zdiv1 = zeps - zval3 
    132130                  wsbio3(ji,jj,jk) = xkr_wsbio_min * ( zeps - zval1 ) / zdiv    & 
    133      &                             - xkr_wsbio_max *   zgm * xkr_eta  / zdiv 
     131                     &             - xkr_wsbio_max *   zgm * xkr_eta  / zdiv 
    134132                  wsbio4(ji,jj,jk) = xkr_wsbio_min *   ( zeps-1. )    / zdiv1   & 
    135      &                             - xkr_wsbio_max *   zfm * xkr_eta  / zdiv1 
     133                     &             - xkr_wsbio_max *   zfm * xkr_eta  / zdiv1 
    136134                  IF( znum == 1.1)   wsbio3(ji,jj,jk) = wsbio4(ji,jj,jk) 
    137135               ENDIF 
    … …  
    140138      END DO 
    141139 
    142       wscal(:,:,:) = MAX( wsbio3(:,:,:), 50. ) 
     140      wscal(:,:,:) = MAX( wsbio3(:,:,:), 50._wp ) 
    143141 
    144142      !   INITIALIZE TO ZERO ALL THE SINKING ARRAYS 
    … …  
    305303#endif 
    306304      ! 
    307        IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
     305      IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
    308306         WRITE(charout, FMT="('sink')") 
    309307         CALL prt_ctl_trc_info(charout) 
    310308         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm) 
    311        ENDIF 
    312  
     309      ENDIF 
     310      ! 
     311      IF( wrk_not_released(3, 2 ) )   CALL ctl_stop('p4z_sink: failed to release workspace arrays') 
     312      ! 
    313313   END SUBROUTINE p4z_sink 
     314 
    314315 
    315316   SUBROUTINE p4z_sink_init 
    … …  
    324325      !! 
    325326      !! ** input   :   Namelist nampiskrs 
    326       !! 
    327327      !!---------------------------------------------------------------------- 
    328328      INTEGER  ::   jk, jn, kiter 
    … …  
    330330      REAL(wp) ::   zws, zwr, zwl,wmax, znummax 
    331331      REAL(wp) ::   zmin, zmax, zl, zr, xacc 
    332  
     332      ! 
    333333      NAMELIST/nampiskrs/ xkr_sfact, xkr_stick ,  & 
    334334         &                xkr_nnano, xkr_ndiat, xkr_nmeso, xkr_naggr 
    335  
    336335      !!---------------------------------------------------------------------- 
     336      ! 
    337337      REWIND( numnat )                     ! read nampiskrs 
    338338      READ  ( numnat, nampiskrs ) 
    … …  
    347347         WRITE(numout,*) '    Nbr of cell in mesozoo size class        xkr_nmeso    = ', xkr_nmeso 
    348348         WRITE(numout,*) '    Nbr of cell in aggregates size class     xkr_naggr    = ', xkr_naggr 
    349      ENDIF 
    350  
    351  
    352      ! max and min vertical particle speed 
    353      xkr_wsbio_min = xkr_sfact * xkr_mass_min**xkr_eta 
    354      xkr_wsbio_max = xkr_sfact * xkr_mass_max**xkr_eta 
    355      WRITE(numout,*) ' max and min vertical particle speed ', xkr_wsbio_min, xkr_wsbio_max 
    356  
    357      ! 
    358      !    effect of the sizes of the different living pools on particle numbers 
    359      !    nano = 2um-20um -> mean size=6.32 um -> ws=2.596 -> xnum=xnnano=2.337 
    360      !    diat and microzoo = 10um-200um -> 44.7 -> 8.732 -> xnum=xndiat=3.718 
    361      !    mesozoo = 200um-2mm -> 632.45 -> 45.14 -> xnum=xnmeso=7.147 
    362      !    aggregates = 200um-10mm -> 1414 -> 74.34 -> xnum=xnaggr=9.877 
    363      !    doc aggregates = 1um 
    364      ! ---------------------------------------------------------- 
    365  
    366      xkr_dnano = 1. / ( xkr_massp * xkr_nnano ) 
    367      xkr_ddiat = 1. / ( xkr_massp * xkr_ndiat ) 
    368      xkr_dmeso = 1. / ( xkr_massp * xkr_nmeso ) 
    369      xkr_daggr = 1. / ( xkr_massp * xkr_naggr ) 
     349      ENDIF 
     350 
     351 
     352      ! max and min vertical particle speed 
     353      xkr_wsbio_min = xkr_sfact * xkr_mass_min**xkr_eta 
     354      xkr_wsbio_max = xkr_sfact * xkr_mass_max**xkr_eta 
     355      WRITE(numout,*) ' max and min vertical particle speed ', xkr_wsbio_min, xkr_wsbio_max 
     356 
     357      ! 
     358      !    effect of the sizes of the different living pools on particle numbers 
     359      !    nano = 2um-20um -> mean size=6.32 um -> ws=2.596 -> xnum=xnnano=2.337 
     360      !    diat and microzoo = 10um-200um -> 44.7 -> 8.732 -> xnum=xndiat=3.718 
     361      !    mesozoo = 200um-2mm -> 632.45 -> 45.14 -> xnum=xnmeso=7.147 
     362      !    aggregates = 200um-10mm -> 1414 -> 74.34 -> xnum=xnaggr=9.877 
     363      !    doc aggregates = 1um 
     364      ! ---------------------------------------------------------- 
     365 
     366      xkr_dnano = 1. / ( xkr_massp * xkr_nnano ) 
     367      xkr_ddiat = 1. / ( xkr_massp * xkr_ndiat ) 
     368      xkr_dmeso = 1. / ( xkr_massp * xkr_nmeso ) 
     369      xkr_daggr = 1. / ( xkr_massp * xkr_naggr ) 
    370370 
    371371      !!--------------------------------------------------------------------- 
    … …  
    379379      WRITE(numout,*)'    kriest : Compute maximum number of particles in aggregates' 
    380380 
    381       xacc     =  0.001 
     381      xacc     =  0.001_wp 
    382382      kiter    = 50 
    383       zmin     =  1.10 
     383      zmin     =  1.10_wp 
    384384      zmax     = xkr_mass_max / xkr_mass_min 
    385385      xkr_frac = zmax 
    … …  
    402402            &     xkr_frac**( -xkr_zeta / znum ) / zdiv ) & 
    403403            & - wmax 
    404 iflag:  DO jn = 1, kiter 
    405            IF( zwl == 0.e0 ) THEN 
    406               znummax = zl 
    407            ELSE IF ( zwr == 0.e0 ) THEN 
    408               znummax = zr 
    409            ELSE 
    410               znummax = ( zr + zl ) / 2. 
    411               zdiv = xkr_zeta + xkr_eta - xkr_eta * znummax 
    412               znum = znummax - 1. 
    413               zws =  xkr_wsbio_min * xkr_zeta / zdiv & 
    414                  & - ( xkr_wsbio_max * xkr_eta * znum * & 
    415                  &     xkr_frac**( -xkr_zeta / znum ) / zdiv ) & 
    416                  & - wmax 
    417               IF( zws * zwl < 0. ) THEN 
    418                  zr = znummax 
    419               ELSE 
    420                  zl = znummax 
    421               ENDIF 
    422               zdiv = xkr_zeta + xkr_eta - xkr_eta * zl 
    423               znum = zl - 1. 
    424               zwl =  xkr_wsbio_min * xkr_zeta / zdiv & 
    425                  & - ( xkr_wsbio_max * xkr_eta * znum * & 
    426                  &     xkr_frac**( -xkr_zeta / znum ) / zdiv ) & 
    427                  & - wmax 
    428  
    429               zdiv = xkr_zeta + xkr_eta - xkr_eta * zr 
    430               znum = zr - 1. 
    431               zwr =  xkr_wsbio_min * xkr_zeta / zdiv & 
    432                  & - ( xkr_wsbio_max * xkr_eta * znum * & 
    433                  &     xkr_frac**( -xkr_zeta / znum ) / zdiv ) & 
    434                  & - wmax 
    435  
    436               IF ( ABS ( zws )  <= xacc ) EXIT iflag 
    437  
    438            ENDIF 
    439  
    440         END DO iflag 
    441  
    442         xnumm(jk) = znummax 
    443         WRITE(numout,*) '       jk = ', jk, ' wmax = ', wmax,' xnum max = ', xnumm(jk) 
    444  
    445      END DO 
    446  
     404iflag:   DO jn = 1, kiter 
     405            IF    ( zwl == 0._wp ) THEN   ;   znummax = zl 
     406            ELSEIF( zwr == 0._wp ) THEN   ;   znummax = zr 
     407            ELSE 
     408               znummax = ( zr + zl ) / 2. 
     409               zdiv = xkr_zeta + xkr_eta - xkr_eta * znummax 
     410               znum = znummax - 1. 
     411               zws =  xkr_wsbio_min * xkr_zeta / zdiv & 
     412                  & - ( xkr_wsbio_max * xkr_eta * znum * & 
     413                  &     xkr_frac**( -xkr_zeta / znum ) / zdiv ) & 
     414                  & - wmax 
     415               IF( zws * zwl < 0. ) THEN   ;   zr = znummax 
     416               ELSE                        ;   zl = znummax 
     417               ENDIF 
     418               zdiv = xkr_zeta + xkr_eta - xkr_eta * zl 
     419               znum = zl - 1. 
     420               zwl =  xkr_wsbio_min * xkr_zeta / zdiv & 
     421                  & - ( xkr_wsbio_max * xkr_eta * znum * & 
     422                  &     xkr_frac**( -xkr_zeta / znum ) / zdiv ) & 
     423                  & - wmax 
     424 
     425               zdiv = xkr_zeta + xkr_eta - xkr_eta * zr 
     426               znum = zr - 1. 
     427               zwr =  xkr_wsbio_min * xkr_zeta / zdiv & 
     428                  & - ( xkr_wsbio_max * xkr_eta * znum * & 
     429                  &     xkr_frac**( -xkr_zeta / znum ) / zdiv ) & 
     430                  & - wmax 
     431               ! 
     432               IF ( ABS ( zws )  <= xacc ) EXIT iflag 
     433               ! 
     434            ENDIF 
     435            ! 
     436         END DO iflag 
     437 
     438         xnumm(jk) = znummax 
     439         WRITE(numout,*) '       jk = ', jk, ' wmax = ', wmax,' xnum max = ', xnumm(jk) 
     440         ! 
     441      END DO 
     442      ! 
    447443  END SUBROUTINE p4z_sink_init 
    448444 
    … …  
    476472         DO jj = 1, jpj 
    477473            DO ji=1,jpi 
    478                zfact = MAX( 0., fsdepw(ji,jj,jk+1) - hmld(ji,jj) ) / 4000. 
     474               zfact = MAX( 0., fsdepw(ji,jj,jk+1) - hmld(ji,jj) ) / 4000._wp 
    479475               wsbio4(ji,jj,jk) = wsbio2 + ( 200.- wsbio2 ) * zfact 
    480476            END DO 
    … …  
    584580#endif 
    585581      ! 
    586        IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
     582      IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
    587583         WRITE(charout, FMT="('sink')") 
    588584         CALL prt_ctl_trc_info(charout) 
    589585         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm) 
    590        ENDIF 
    591  
     586      ENDIF 
     587      ! 
    592588   END SUBROUTINE p4z_sink 
     589 
    593590 
    594591   SUBROUTINE p4z_sink_init 
    … …  
    611608      !!      transport term, i.e.  div(u*tra). 
    612609      !!--------------------------------------------------------------------- 
     610      USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 
     611      USE wrk_nemo, ONLY: ztraz => wrk_3d_2, zakz => wrk_3d_3, zwsink2 => wrk_3d_4 
     612      ! 
    613613      INTEGER , INTENT(in   )                         ::   jp_tra    ! tracer index index       
    614614      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pwsink    ! sinking speed 
    … …  
    617617      INTEGER  ::   ji, jj, jk, jn 
    618618      REAL(wp) ::   zigma,zew,zign, zflx, zstep 
    619       REAL(wp), DIMENSION(jpi,jpj,jpk) ::  ztraz, zakz 
    620       REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zwsink2 
    621       !!--------------------------------------------------------------------- 
    622  
     619      !!--------------------------------------------------------------------- 
     620 
     621      IF(  wrk_in_use(3, 2,3,4 ) ) THEN 
     622         CALL ctl_stop('p4z_sink2: requested workspace arrays unavailable') 
     623         RETURN 
     624      END IF 
    623625 
    624626      zstep = rfact2 / 2. 
    … …  
    701703      END DO 
    702704 
    703       trn(:,:,:,jp_tra) = trb(:,:,:,jp_tra) 
    704       psinkflx(:,:,:)   = 2. * psinkflx(:,:,:) 
    705  
     705      trn     (:,:,:,jp_tra) = trb(:,:,:,jp_tra) 
     706      psinkflx(:,:,:)        = 2. * psinkflx(:,:,:) 
     707      ! 
     708      IF( wrk_not_released(3, 2,3,4) )   CALL ctl_stop('p4z_sink2: failed to release workspace arrays') 
    706709      ! 
    707710   END SUBROUTINE p4z_sink2 
    708711 
     712 
     713   INTEGER FUNCTION p4z_sink_alloc() 
     714      !!---------------------------------------------------------------------- 
     715      !!                     ***  ROUTINE p4z_sink_alloc  *** 
     716      !!---------------------------------------------------------------------- 
     717      ALLOCATE( wsbio3 (jpi,jpj,jpk) , wsbio4  (jpi,jpj,jpk) , wscal(jpi,jpj,jpk) ,     & 
     718         &      sinking(jpi,jpj,jpk) , sinking2(jpi,jpj,jpk)                      ,     &                 
     719         &      sinkcal(jpi,jpj,jpk) , sinksil (jpi,jpj,jpk)                      ,     &                 
     720#if defined key_kriest 
     721         &      xnumm(jpk)                                                        ,     &                 
     722#else 
     723         &      sinkfer2(jpi,jpj,jpk)                                             ,     &                 
     724#endif 
     725         &      sinkfer(jpi,jpj,jpk)                                              , STAT=p4z_sink_alloc )                 
     726         ! 
     727      IF( p4z_sink_alloc /= 0 ) CALL ctl_warn('p4z_sink_alloc : failed to allocate arrays.') 
     728      ! 
     729   END FUNCTION p4z_sink_alloc 
     730    
    709731#else 
    710732   !!====================================================================== 
Note: See TracChangeset for help on using the changeset viewer.