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 7403 for branches/2016/dev_merge_2016/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90 – NEMO

Ignore:
Timestamp:
2016-11-30T17:56:53+01:00 (8 years ago)
Author:
timgraham
Message:

Merge dev_INGV_METO_merge_2016 into branch

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90

    r6421 r7403  
    66   !! History :   1.0  !  2004-03 (O. Aumont) Original code 
    77   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90 
    8    !!---------------------------------------------------------------------- 
    9 #if defined key_pisces 
    10    !!---------------------------------------------------------------------- 
    11    !!   'key_pisces'                                       PISCES bio-model 
    128   !!---------------------------------------------------------------------- 
    139   !!   p4zsms         :  Time loop of passive tracers sms 
     
    6965      INTEGER ::   ji, jj, jk, jnt, jn, jl 
    7066      REAL(wp) ::  ztra 
    71 #if defined key_kriest 
    72       REAL(wp) ::  zcoef1, zcoef2 
    73 #endif 
    7467      CHARACTER (len=25) :: charout 
    7568      !!--------------------------------------------------------------------- 
     
    8376        CALL p4z_che                              ! initialize the chemical constants 
    8477        ! 
    85         IF( .NOT. ln_rsttr ) THEN  ;   CALL p4z_ph_ini   !  set PH at kt=nit000  
     78        IF( .NOT. ln_rsttr ) THEN  ;   CALL ahini_for_at(hi)   !  set PH at kt=nit000  
    8679        ELSE                       ;   CALL p4z_rst( nittrc000, 'READ' )  !* read or initialize all required fields  
    8780        ENDIF 
     
    9184      IF( ln_pisdmp .AND. MOD( kt - nn_dttrc, nn_pisdmp ) == 0 )   CALL p4z_dmp( kt )      ! Relaxation of some tracers 
    9285      ! 
    93       !                                                                    !   set time step size (Euler/Leapfrog) 
    94       IF( ( neuler == 0 .AND. kt == nittrc000 ) .OR. ln_top_euler ) THEN   ;    rfact = rdttrc     !  at nittrc000 
    95       ELSEIF( kt <= nittrc000 + nn_dttrc )                          THEN   ;    rfact = 2. * rdttrc   ! at nittrc000 or nittrc000+nn_dttrc (Leapfrog) 
    96       ENDIF 
     86      rfact = r2dttrc 
    9787      ! 
    9888      IF( ( ln_top_euler .AND. kt == nittrc000 )  .OR. ( .NOT.ln_top_euler .AND. kt <= nittrc000 + nn_dttrc ) ) THEN 
    9989         rfactr  = 1. / rfact 
    100          rfact2  = rfact / FLOAT( nrdttrc ) 
     90         rfact2  = rfact / REAL( nrdttrc, wp ) 
    10191         rfact2r = 1. / rfact2 
    10292         xstep = rfact2 / rday         ! Time step duration for biology 
     
    165155      END DO 
    166156 
    167 #if defined key_kriest 
    168       !  
    169       zcoef1 = 1.e0 / xkr_massp  
    170       zcoef2 = 1.e0 / xkr_massp / 1.1 
    171       DO jk = 1,jpkm1 
    172          trb(:,:,jk,jpnum) = MAX(  trb(:,:,jk,jpnum), trb(:,:,jk,jppoc) * zcoef1 / xnumm(jk)  ) 
    173          trb(:,:,jk,jpnum) = MIN(  trb(:,:,jk,jpnum), trb(:,:,jk,jppoc) * zcoef2              ) 
    174       END DO 
    175       ! 
    176 #endif 
    177       ! 
    178157      ! 
    179158      IF( l_trdtrc ) THEN 
     
    212191      !! ** input   :   file 'namelist.trc.s' containing the following 
    213192      !!             namelist: natext, natbio, natsms 
    214       !!                       natkriest ("key_kriest") 
    215       !!---------------------------------------------------------------------- 
    216       NAMELIST/nampisbio/ nrdttrc, wsbio, xkmort, ferat3, wsbio2, niter1max, niter2max 
    217 #if defined key_kriest 
    218       NAMELIST/nampiskrp/ xkr_eta, xkr_zeta, xkr_ncontent, xkr_mass_min, xkr_mass_max 
    219 #endif 
     193      !!---------------------------------------------------------------------- 
     194      NAMELIST/nampisbio/ nrdttrc, wsbio, xkmort, ferat3, wsbio2, wsbio2max, wsbio2scale,    & 
     195         &                   niter1max, niter2max, wfep, ldocp, ldocz, lthet,  & 
     196         &                   no3rat3, po4rat3 
     197 
    220198      NAMELIST/nampisdmp/ ln_pisdmp, nn_pisdmp 
    221199      NAMELIST/nampismass/ ln_check_mass 
     
    234212      IF(lwp) THEN                         ! control print 
    235213         WRITE(numout,*) ' Namelist : nampisbio' 
    236          WRITE(numout,*) '    frequence pour la biologie                nrdttrc   =', nrdttrc 
    237          WRITE(numout,*) '    POC sinking speed                         wsbio     =', wsbio 
    238          WRITE(numout,*) '    half saturation constant for mortality    xkmort    =', xkmort 
    239          WRITE(numout,*) '    Fe/C in zooplankton                       ferat3    =', ferat3 
    240          WRITE(numout,*) '    Big particles sinking speed               wsbio2    =', wsbio2 
     214         WRITE(numout,*) '    frequence pour la biologie                nrdttrc    =', nrdttrc 
     215         WRITE(numout,*) '    POC sinking speed                         wsbio      =', wsbio 
     216         WRITE(numout,*) '    half saturation constant for mortality    xkmort     =', xkmort  
     217         IF( ln_p5z ) THEN 
     218            WRITE(numout,*) '    N/C in zooplankton                        no3rat3    =', no3rat3 
     219            WRITE(numout,*) '    P/C in zooplankton                        po4rat3    =', po4rat3 
     220         ENDIF 
     221         WRITE(numout,*) '    Fe/C in zooplankton                       ferat3     =', ferat3 
     222         WRITE(numout,*) '    Big particles sinking speed               wsbio2     =', wsbio2 
     223         WRITE(numout,*) '    Big particles maximum sinking speed       wsbio2max  =', wsbio2max 
     224         WRITE(numout,*) '    Big particles sinking speed length scale  wsbio2scale =', wsbio2scale 
    241225         WRITE(numout,*) '    Maximum number of iterations for POC      niter1max =', niter1max 
    242226         WRITE(numout,*) '    Maximum number of iterations for GOC      niter2max =', niter2max 
    243       ENDIF 
    244  
    245 #if defined key_kriest 
    246  
    247       !                               ! nampiskrp : kriest parameters 
    248       !                               ! ----------------------------- 
    249       REWIND( numnatp_ref )              ! Namelist nampiskrp in reference namelist : Pisces Kriest 
    250       READ  ( numnatp_ref, nampiskrp, IOSTAT = ios, ERR = 903) 
    251 903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampiskrp in reference namelist', lwp ) 
    252  
    253       REWIND( numnatp_cfg )              ! Namelist nampiskrp in configuration namelist : Pisces Kriest 
    254       READ  ( numnatp_cfg, nampiskrp, IOSTAT = ios, ERR = 904 ) 
    255 904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampiskrp in configuration namelist', lwp ) 
    256       IF(lwm) WRITE ( numonp, nampiskrp ) 
    257  
    258       IF(lwp) THEN 
    259          WRITE(numout,*) 
    260          WRITE(numout,*) ' Namelist : nampiskrp' 
    261          WRITE(numout,*) '    Sinking  exponent                        xkr_eta      = ', xkr_eta 
    262          WRITE(numout,*) '    N content exponent                       xkr_zeta     = ', xkr_zeta 
    263          WRITE(numout,*) '    N content factor                         xkr_ncontent = ', xkr_ncontent 
    264          WRITE(numout,*) '    Minimum mass for Aggregates              xkr_mass_min = ', xkr_mass_min 
    265          WRITE(numout,*) '    Maximum mass for Aggregates              xkr_mass_max = ', xkr_mass_max 
    266          WRITE(numout,*) 
    267      ENDIF 
    268  
    269  
    270      ! Computation of some variables 
    271      xkr_massp = xkr_ncontent * 7.625 * xkr_mass_min**xkr_zeta 
    272  
    273 #endif 
     227         IF( ln_ligand ) THEN 
     228            WRITE(numout,*) '    FeP sinking speed                             wfep   =', wfep 
     229            IF( ln_p4z ) THEN 
     230              WRITE(numout,*) '    Phyto ligand production per unit doc          ldocp  =', ldocp 
     231              WRITE(numout,*) '    Zoo ligand production per unit doc            ldocz  =', ldocz 
     232              WRITE(numout,*) '    Proportional loss of ligands due to Fe uptake lthet  =', lthet 
     233            ENDIF 
     234         ENDIF 
     235      ENDIF 
     236 
    274237 
    275238      REWIND( numnatp_ref )              ! Namelist nampisdmp in reference namelist : Pisces damping 
     
    308271   END SUBROUTINE p4z_sms_init 
    309272 
    310    SUBROUTINE p4z_ph_ini 
    311       !!--------------------------------------------------------------------- 
    312       !!                   ***  ROUTINE p4z_ini_ph  *** 
    313       !! 
    314       !!  ** Purpose : Initialization of chemical variables of the carbon cycle 
    315       !!--------------------------------------------------------------------- 
    316       INTEGER  ::  ji, jj, jk 
    317       REAL(wp) ::  zcaralk, zbicarb, zco3 
    318       REAL(wp) ::  ztmas, ztmas1 
    319       !!--------------------------------------------------------------------- 
    320  
    321       ! Set PH from  total alkalinity, borat (???), akb3 (???) and ak23 (???) 
    322       ! -------------------------------------------------------- 
    323       DO jk = 1, jpk 
    324          DO jj = 1, jpj 
    325             DO ji = 1, jpi 
    326                ztmas   = tmask(ji,jj,jk) 
    327                ztmas1  = 1. - tmask(ji,jj,jk) 
    328                zcaralk = trb(ji,jj,jk,jptal) - borat(ji,jj,jk) / (  1. + 1.E-8 / ( rtrn + akb3(ji,jj,jk) )  ) 
    329                zco3    = ( zcaralk - trb(ji,jj,jk,jpdic) ) * ztmas + 0.5e-3 * ztmas1 
    330                zbicarb = ( 2. * trb(ji,jj,jk,jpdic) - zcaralk ) 
    331                hi(ji,jj,jk) = ( ak23(ji,jj,jk) * zbicarb / zco3 ) * ztmas + 1.e-9 * ztmas1 
    332             END DO 
    333          END DO 
    334      END DO 
    335      ! 
    336    END SUBROUTINE p4z_ph_ini 
    337  
    338273   SUBROUTINE p4z_rst( kt, cdrw ) 
    339274      !!--------------------------------------------------------------------- 
     
    349284      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag 
    350285      ! 
    351       INTEGER  ::  ji, jj, jk 
    352       REAL(wp) ::  zcaralk, zbicarb, zco3 
    353       REAL(wp) ::  ztmas, ztmas1 
    354286      !!--------------------------------------------------------------------- 
    355287 
     
    363295            CALL iom_get( numrtr, jpdom_autoglo, 'PH' , hi(:,:,:)  ) 
    364296         ELSE 
    365 !            hi(:,:,:) = 1.e-9  
    366             CALL p4z_ph_ini 
     297            CALL ahini_for_at(hi) 
    367298         ENDIF 
    368299         CALL iom_get( numrtr, jpdom_autoglo, 'Silicalim', xksi(:,:) ) 
     
    379310         ENDIF 
    380311         ! 
     312         IF( ln_p5z ) THEN 
     313            IF( iom_varid( numrtr, 'sized', ldstop = .FALSE. ) > 0 ) THEN 
     314               CALL iom_get( numrtr, jpdom_autoglo, 'sizep' , sized(:,:,:)  ) 
     315               CALL iom_get( numrtr, jpdom_autoglo, 'sizen' , sized(:,:,:)  ) 
     316               CALL iom_get( numrtr, jpdom_autoglo, 'sized' , sized(:,:,:)  ) 
     317            ELSE 
     318               sizep(:,:,:) = 1. 
     319               sizen(:,:,:) = 1. 
     320               sized(:,:,:) = 1. 
     321            ENDIF 
     322        ENDIF 
     323        ! 
    381324      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN 
    382325         IF( kt == nitrst ) THEN 
     
    389332         CALL iom_rstput( kt, nitrst, numrtw, 'Silicamax', xksimax(:,:) ) 
    390333         CALL iom_rstput( kt, nitrst, numrtw, 'tcflxcum', t_oce_co2_flx_cum ) 
     334         IF( ln_p5z ) THEN 
     335            CALL iom_rstput( kt, nitrst, numrtw, 'sizep', sized(:,:,:) ) 
     336            CALL iom_rstput( kt, nitrst, numrtw, 'sizen', sized(:,:,:) ) 
     337            CALL iom_rstput( kt, nitrst, numrtw, 'sized', sized(:,:,:) ) 
     338         ENDIF 
    391339      ENDIF 
    392340      ! 
     
    475423      REAL(wp)             ::  zrdenittot, zsdenittot, znitrpottot 
    476424      CHARACTER(LEN=100)   ::   cltxt 
    477       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zvol 
    478425      INTEGER :: jk 
     426      REAL(wp), POINTER, DIMENSION(:,:,:) :: zwork 
    479427      !!---------------------------------------------------------------------- 
    480428 
     
    496444      ENDIF 
    497445 
     446      CALL wrk_alloc( jpi, jpj, jpk, zwork ) 
    498447      ! 
    499448      IF( iom_use( "pno3tot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN 
    500449         !   Compute the budget of NO3, ALK, Si, Fer 
    501          no3budget = glob_sum( (   trn(:,:,:,jpno3) + trn(:,:,:,jpnh4)  & 
    502             &                    + trn(:,:,:,jpphy) + trn(:,:,:,jpdia)  & 
    503             &                    + trn(:,:,:,jpzoo) + trn(:,:,:,jpmes)  & 
    504             &                    + trn(:,:,:,jppoc)                     & 
    505 #if ! defined key_kriest 
    506             &                    + trn(:,:,:,jpgoc)                     & 
    507 #endif 
    508             &                    + trn(:,:,:,jpdoc)                     ) * cvol(:,:,:)  ) 
    509          ! 
    510          no3budget = no3budget / areatot 
    511          CALL iom_put( "pno3tot", no3budget ) 
     450         IF( ln_p4z ) THEN 
     451            zwork(:,:,:) =    trn(:,:,:,jpno3) + trn(:,:,:,jpnh4)                      & 
     452               &          +   trn(:,:,:,jpphy) + trn(:,:,:,jpdia)                      & 
     453               &          +   trn(:,:,:,jppoc) + trn(:,:,:,jpgoc)  + trn(:,:,:,jpdoc)  &         
     454               &          +   trn(:,:,:,jpzoo) + trn(:,:,:,jpmes)  
     455        ELSE 
     456            zwork(:,:,:) =    trn(:,:,:,jpno3) + trn(:,:,:,jpnh4) + trn(:,:,:,jpnph)   & 
     457               &          +   trn(:,:,:,jpndi) + trn(:,:,:,jpnpi)                      &  
     458               &          +   trn(:,:,:,jppon) + trn(:,:,:,jpgon) + trn(:,:,:,jpdon)   & 
     459               &          + ( trn(:,:,:,jpzoo) + trn(:,:,:,jpmes) ) * no3rat3  
     460        ENDIF 
     461        ! 
     462        no3budget = glob_sum( zwork(:,:,:) * cvol(:,:,:)  )   
     463        no3budget = no3budget / areatot 
     464        CALL iom_put( "pno3tot", no3budget ) 
    512465      ENDIF 
    513466      ! 
    514467      IF( iom_use( "ppo4tot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN 
    515          po4budget = glob_sum( (   trn(:,:,:,jppo4)                     & 
    516             &                    + trn(:,:,:,jpphy) + trn(:,:,:,jpdia)  & 
    517             &                    + trn(:,:,:,jpzoo) + trn(:,:,:,jpmes)  & 
    518             &                    + trn(:,:,:,jppoc)                     & 
    519 #if ! defined key_kriest 
    520             &                    + trn(:,:,:,jpgoc)                     & 
    521 #endif 
    522             &                    + trn(:,:,:,jpdoc)                     ) * cvol(:,:,:)  ) 
    523          po4budget = po4budget / areatot 
    524          CALL iom_put( "ppo4tot", po4budget ) 
     468         IF( ln_p4z ) THEN 
     469            zwork(:,:,:) =    trn(:,:,:,jppo4)                                         & 
     470               &          +   trn(:,:,:,jpphy) + trn(:,:,:,jpdia)                      & 
     471               &          +   trn(:,:,:,jppoc) + trn(:,:,:,jpgoc)  + trn(:,:,:,jpdoc)  &         
     472               &          +   trn(:,:,:,jpzoo) + trn(:,:,:,jpmes)  
     473        ELSE 
     474            zwork(:,:,:) =    trn(:,:,:,jppo4) + trn(:,:,:,jppph)                      & 
     475               &          +   trn(:,:,:,jppdi) + trn(:,:,:,jpppi)                      &  
     476               &          +   trn(:,:,:,jppop) + trn(:,:,:,jpgop) + trn(:,:,:,jpdop)   & 
     477               &          + ( trn(:,:,:,jpzoo) + trn(:,:,:,jpmes) ) * po4rat3  
     478        ENDIF 
     479        ! 
     480        po4budget = glob_sum( zwork(:,:,:) * cvol(:,:,:)  )   
     481        po4budget = po4budget / areatot 
     482        CALL iom_put( "ppo4tot", po4budget ) 
    525483      ENDIF 
    526484      ! 
    527485      IF( iom_use( "psiltot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN 
    528          silbudget = glob_sum( (   trn(:,:,:,jpsil) + trn(:,:,:,jpgsi)  & 
    529             &                    + trn(:,:,:,jpdsi)                     ) * cvol(:,:,:)  ) 
    530          ! 
     486         zwork(:,:,:) =  trn(:,:,:,jpsil) + trn(:,:,:,jpgsi) + trn(:,:,:,jpdsi)  
     487         ! 
     488         silbudget = glob_sum( zwork(:,:,:) * cvol(:,:,:)  )   
    531489         silbudget = silbudget / areatot 
    532490         CALL iom_put( "psiltot", silbudget ) 
     
    534492      ! 
    535493      IF( iom_use( "palktot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN 
    536          alkbudget = glob_sum( (   trn(:,:,:,jpno3) * rno3              & 
    537             &                    + trn(:,:,:,jptal)                     & 
    538             &                    + trn(:,:,:,jpcal) * 2.                ) * cvol(:,:,:)  ) 
    539          ! 
     494         zwork(:,:,:) =  trn(:,:,:,jpno3) * rno3 + trn(:,:,:,jptal) + trn(:,:,:,jpcal) * 2.               
     495         ! 
     496         alkbudget = glob_sum( zwork(:,:,:) * cvol(:,:,:)  )         ! 
    540497         alkbudget = alkbudget / areatot 
    541498         CALL iom_put( "palktot", alkbudget ) 
     
    543500      ! 
    544501      IF( iom_use( "pfertot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN 
    545          ferbudget = glob_sum( (   trn(:,:,:,jpfer) + trn(:,:,:,jpnfe)  & 
    546             &                    + trn(:,:,:,jpdfe)                     & 
    547 #if ! defined key_kriest 
    548             &                    + trn(:,:,:,jpbfe)                     & 
    549 #endif 
    550             &                    + trn(:,:,:,jpsfe)                     & 
    551             &                    + trn(:,:,:,jpzoo) * ferat3            & 
    552             &                    + trn(:,:,:,jpmes) * ferat3            ) * cvol(:,:,:)  ) 
    553          ! 
     502         zwork(:,:,:) =   trn(:,:,:,jpfer) + trn(:,:,:,jpnfe) + trn(:,:,:,jpdfe)   & 
     503            &         +   trn(:,:,:,jpbfe) + trn(:,:,:,jpsfe)                      & 
     504            &         + ( trn(:,:,:,jpzoo) + trn(:,:,:,jpmes) )  * ferat3     
     505         IF( ln_ligand)  zwork(:,:,:) = zwork(:,:,:) + trn(:,:,:,jpfep)                 
     506         ! 
     507         ferbudget = glob_sum( zwork(:,:,:) * cvol(:,:,:)  )   
    554508         ferbudget = ferbudget / areatot 
    555509         CALL iom_put( "pfertot", ferbudget ) 
    556510      ENDIF 
    557511      ! 
    558  
     512      CALL wrk_dealloc( jpi, jpj, jpk, zwork ) 
     513      ! 
    559514      ! Global budget of N SMS : denitrification in the water column and in the sediment 
    560515      !                          nitrogen fixation by the diazotrophs 
     
    600555   END SUBROUTINE p4z_chk_mass 
    601556 
    602 #else 
    603    !!====================================================================== 
    604    !!  Dummy module :                                   No PISCES bio-model 
    605    !!====================================================================== 
    606 CONTAINS 
    607    SUBROUTINE p4z_sms( kt )                   ! Empty routine 
    608       INTEGER, INTENT( in ) ::   kt 
    609       WRITE(*,*) 'p4z_sms: You should not have seen this print! error?', kt 
    610    END SUBROUTINE p4z_sms 
    611 #endif  
    612  
    613557   !!====================================================================== 
    614558END MODULE p4zsms  
Note: See TracChangeset for help on using the changeset viewer.