SUBROUTINE ice_bio_sms(nlay_i,kideb,kiut) ! !-----------------------------------------------------------------------------! ! *** ice_bio_sms *** ! Biological sources minus sinks ! (c) Martin Vancoppenolle, LOCEAN, October 2012 ! ! version : source_3.03 ! status : in permanent revision :-) ! !=============================================================================! ! INCLUDE 'type.com' INCLUDE 'para.com' INCLUDE 'const.com' INCLUDE 'ice.com' INCLUDE 'thermo.com' INCLUDE 'bio.com' INTEGER, DIMENSION(3) :: & zindex REAL, DIMENSION(nlay_bio) :: & zek, & zlim_nut LOGICAL :: & ln_write_bio ln_write_bio = .TRUE. zeps = 1.0e-14 ji = 1 ! !=============================================================================! ! !-----------------------------------------------------------------------------! ! 1) Control prints !-----------------------------------------------------------------------------! ! IF ( ln_write_bio ) THEN WRITE(numout,*) WRITE(numout,*) ' ** ice_bio_sms : ' WRITE(numout,*) ' ~~~~~~~~~~~~~~~~ ' WRITE(numout,*) CALL ice_brine WRITE(numout,*) ' Initial tracer values ' WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~ ' DO jn = 1, ntra_bio IF ( flag_active(jn) ) THEN WRITE(numout,*) ' Tracer : ', biotr_i_nam(jn) WRITE(numout,*) ' cbu_i_bio: ', ( cbu_i_bio(jn,layer), & layer = 1, nlay_bio ) ENDIF ! flag_active END DO ! jn WRITE(numout,*) WRITE(numout,*) ' Parameters ' WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~ ' WRITE(numout,*) ' si_c : ', si_c WRITE(numout,*) ' no3_c : ', no3_c WRITE(numout,*) ' po4_c : ', po4_c WRITE(numout,*) ' mumax_bio : ', mumax_bio WRITE(numout,*) ' klys_bio : ', klys_bio WRITE(numout,*) ' khs_si_bio: ', khs_si_bio WRITE(numout,*) ' ek_bio : ', ek_bio WRITE(numout,*) ' astar_alg : ', astar_alg WRITE(numout,*) ' rg_bio : ', rg_bio WRITE(numout,*) ' nn_lim_sal: ', nn_lim_sal WRITE(numout,*) ' lim_sal_wid :', lim_sal_wid WRITE(numout,*) ' lim_sal_smax:', lim_sal_smax WRITE(numout,*) ENDIF ! ln_write_bio ! !-----------------------------------------------------------------------------! ! 2) Source Minus Sink terms !-----------------------------------------------------------------------------! ! ! Reference of the upper layer IF ( ( c_grid .EQ. 'SL' ) .OR. ( c_grid .EQ. 'BA' ) ) THEN lim_dsi(:) = 0.; lim_no3(:) = 0.; lim_po4(:) = 0. lim_lig(:) = 0.; lim_tem(:) = 0.; lim_sal(:) = 0. lys_bio(:) = 0.; rem_bio(:) = 0.; syn_bio(:) = 0. rsp_bio(:) = 0.; exu_bio(:) = 0. lim_dsi(nlay_bio) = 1.; lim_no3(nlay_bio) = 1. lim_po4(nlay_bio) = 1. lim_lig(nlay_bio) = 1.; lim_tem(nlay_bio) = 1. lim_sal(nlay_bio) = 1. ELSE lim_dsi(:) = 1.; lim_no3(:) = 1.; lim_po4(:) = 1. lim_lig(:) = 1.; lim_tem(:) = 1.; lim_sal(:) = 1. lys_bio(:) = 0.; rem_bio(:) = 0.; syn_bio(:) = 0. rsp_bio(:) = 0.; exu_bio(:) = 0. ENDIF !------------------------------------------------------------------------ ! 2.1 Limitation factors !------------------------------------------------------------------------ ! 2.1.1 Nutrients & temperature !-------------------------------- DO layer = layer_00, nlay_bio IF ( ln_lim_dsi ) ! DSi limitation & lim_dsi(layer) = c_i_bio(1,layer) & / ( khs_si_bio + c_i_bio(1,layer) ) IF ( ln_lim_no3 ) ! NO3 limitation & lim_no3(layer) = c_i_bio(2,layer) & / ( khs_n_bio + c_i_bio(2,layer) ) IF ( ln_lim_po4 ) ! PO4 limitation & lim_po4(layer) = c_i_bio(3,layer) & / ( khs_p_bio + c_i_bio(3,layer) ) IF ( ln_lim_dsi .OR. ln_lim_no3 .OR. ln_lim_po4 ) ! Nut limitation ! & zlim_nut(layer) = lim_dsi(layer) * lim_no3(layer) ! & * lim_po4(layer) & zlim_nut(layer) = MIN ( lim_dsi(layer) , lim_no3(layer) & , lim_po4(layer) ) IF ( ln_lim_tem ) ! Temperature factor & lim_tem(layer) = EXP( rg_bio * tc_bio(layer) ) END DO ! layer ! 2.1.2 Brine salinity !---------------------- IF ( ln_lim_sal ) THEN SELECT CASE ( nn_lim_sal ) CASE (1) ! Arrigo and Sullivan (1992) DO layer = layer_00, nlay_bio zarg = 2.16 - 8.30e-5 * sbr_bio(layer)**2.11 & - 0.55 * LOG( sbr_bio(layer) ) lim_sal(layer) = EXP ( - zarg * zarg ) END DO CASE (2) ! Simple Gaussian ( Martouf and Grozny, 2010 ) DO layer = layer_00, nlay_bio zarg = ( sbr_bio(layer) - lim_sal_smax ) / lim_sal_wid lim_sal(layer) = EXP( - zarg * zarg ) END DO END SELECT ENDIF ! ln_lim_sal ! 2.1.3 Light !------------- IF ( ln_lim_lig ) THEN SELECT CASE ( nn_phs ) CASE (1) ! prescribed Ek zek(:) = ek_bio chlC_bio(:) = chla_c CASE (2) ! Ek depends on chla/C zek(:) = mumax_bio / ( alpha_bio * chla_c ) chlC_bio(:) = chla_c CASE (3) ! full account for limitation factors zek(:) = mumax_bio * lim_tem(:) * lim_sal(:) & * zlim_nut(:) / ( alpha_bio * chla_c ) chlC_bio(:) = chla_c CASE (4) ! varying chl-a carbon ratio zthmin = 0.01 ! min zthmax = 0.05 ! max DO layer = layer_00, nlay_bio chlC_bio (layer) = ( zthmax - ( zthmax - zthmin ) * & MIN( par_bio(layer) / Estar, 1.0 ) ) * zlim_nut(layer) END DO zek(:) = mumax_bio * lim_tem(:) * lim_sal(:) & * zlim_nut(:) / ( alpha_bio * chlC_bio(:) ) END SELECT DO layer = layer_00, nlay_bio ! Jassby & Platt 76 formulation lim_lig(layer) = TANH ( par_bio(layer) / zek(layer) ) END DO WRITE(numout,*) WRITE(numout,*) ' lim_lig : ', ( lim_lig(layer), & layer = layer_00, nlay_bio ) WRITE(numout,*) ' zlim_nut: ', ( zlim_nut(layer), & layer = layer_00, nlay_bio ) WRITE(numout,*) ' par_bio : ', ( par_bio(layer), & layer = layer_00, nlay_bio ) WRITE(numout,*) ' zek : ', ( zek(layer) , & layer = layer_00, nlay_bio ) WRITE(numout,*) ' chlC_bio: ', ( chlC_bio(layer), & layer = layer_00, nlay_bio ) WRITE(numout,*) ENDIF ! ln_lim_lig IF ( .NOT. ln_lim_lig ) chlC_bio(:) = chla_c !------------------------------------------------------------------------ ! 2.2 Loss & remineralization !------------------------------------------------------------------------ ! Exudation is zero exu_bio(:) = 0.0 DO layer = layer_00, nlay_bio ! Respiration (NP, NPDr) ! IF ( ln_rsp ) ... rsp_bio(layer) = krsp_bio * cbu_i_bio(jn_aoc,layer) * & lim_tem(layer) ! Lysis (NP, NPDr) IF ( ln_lys ) THEN IF ( nn_lys .GE. 1 ) ! no T-dependence & lys_bio(layer) = klys_bio * cbu_i_bio(jn_aoc,layer) IF ( nn_lys .EQ. 2 ) ! T-dependence & lys_bio(layer) = lys_bio(layer) * lim_tem(layer) ENDIF END DO ! Remineralization IF ( ln_rem ) THEN SELECT CASE ( nn_bio_opt ) CASE(0) ! NP rem_bio(layer_00:nlay_bio) = frem_bio * & lys_bio(layer_00:nlay_bio) CASE(1) ! NPDr IF ( nn_rem .GE. 1 ) ! no T-dependence & rem_bio(layer_00:nlay_bio) = krem_bio * & cbu_i_bio(jn_eoc,layer_00:nlay_bio) IF ( nn_rem .EQ. 2 ) ! T-dependence & rem_bio(layer_00:nlay_bio) = rem_bio(layer_00:nlay_bio) * & EXP( rg_bac * tc_bio(layer_00:nlay_bio) ) END SELECT ENDIF !------------------------------------------------------------------------ ! 2.3 Diatom Synthesis !------------------------------------------------------------------------ DO layer = layer_00, nlay_bio IF ( ln_syn ) THEN ! raw synthesis (in s-1) z_syn = mumax_bio & * lim_tem(layer) * lim_sal(layer) * lim_lig(layer) & * zlim_nut(layer) ! preclude excessive growth if insufficient nutrient stock zsyn1 = z_syn * cbu_i_bio(jn_aoc,layer) ! normal carbon production rate (mmolC/m3/s) zsyn2 = zsyn1; zsyn3 = zsyn1; zsyn4 = zsyn1 IF ( ln_lim_dsi ) ! prevent exhaustion of Si & zsyn2 = cbu_i_bio(jn_dsi,layer) / ddtb / si_c - zeps IF ( ln_lim_no3 ) ! prevent exhaustion of N & zsyn3 = cbu_i_bio(jn_din,layer) / ddtb / no3_c - zeps IF ( ln_lim_po4 ) ! prevent exhaustion of P & zsyn4 = cbu_i_bio(jn_dip,layer) / ddtb / po4_c - zeps ! nutrient stock-limited PP rate syn_bio(layer) = MIN( zsyn1, zsyn2, zsyn3, zsyn4 ) ELSE ! ln_syn is .FALSE. ! PP is prescribed to zero everywhere ! except in the bottom layer ! where it equals pp_presc ! if sun is shining IF ( ( layer .EQ. nlay_bio ) .AND. ( fsolgb(ji) .GT. 5.0 ) ) & syn_bio(layer) = pp_presc ENDIF ! ln_syn END DO ! layer !------------------------------------------------------------------------ ! 2.4 Update reservoirs !------------------------------------------------------------------------ DO layer = layer_00, nlay_bio ! DSi uptake zd_dsi1 = si_c * ( rem_bio(layer) + & rsp_bio(layer) - syn_bio(layer) ) * ddtb zd_dsi2 = - cbu_i_bio(jn_dsi,layer) ! Maximum DSi uptake zd_dsi = MAX ( zd_dsi1, zd_dsi2 ) ! NO3 uptake zd_no31 = no3_c * ( rem_bio(layer) + & rsp_bio(layer) - syn_bio(layer) ) * ddtb zd_no32 = - cbu_i_bio(jn_din,layer) ! Max NO3 uptake zd_no3 = MAX ( zd_no31, zd_no32 ) ! update PO4 zd_po41 = po4_c * ( rem_bio(layer) + & rsp_bio(layer) - syn_bio(layer) ) * ddtb zd_po42 = - cbu_i_bio(jn_dip,layer) ! Max PO4 uptake zd_po4 = MAX ( zd_po41, zd_po42 ) ! update algae zd_daf1 = ( syn_bio(layer) - lys_bio(layer) - & rsp_bio(layer) - exu_bio(layer) ) * ddtb zd_daf2 = - cbu_i_bio(jn_aoc,layer) ! Max algal carbon loss zd_daf = MAX ( zd_daf1, zd_daf2 ) ! update detritus IF ( nn_bio_opt .GE. 1 ) THEN zd_eoc1 = ( lys_bio(layer) + exu_bio(layer) - & rem_bio(layer) ) * ddtb zd_eoc2 = - cbu_i_bio(jn_eoc,layer) ! Max detrital carbon loss zd_eoc = MAX( zd_eoc1, zd_eoc2 ) ENDIF cbu_i_bio(jn_dsi,layer) = cbu_i_bio(jn_dsi,layer) + zd_dsi cbu_i_bio(jn_din,layer) = cbu_i_bio(jn_din,layer) + zd_no3 cbu_i_bio(jn_dip,layer) = cbu_i_bio(jn_dip,layer) + zd_po4 cbu_i_bio(jn_aoc,layer) = cbu_i_bio(jn_aoc,layer) + zd_daf IF ( nn_bio_opt .GE. 1 ) & cbu_i_bio(jn_eoc,layer) = cbu_i_bio(jn_eoc,layer) + zd_eoc IF ( ln_write_bio ) THEN WRITE(numout,*) ' Terms for layer: ', layer WRITE(numout,*) ' ~~~~~~~~~~~~~~~~ ' WRITE(numout,*) ' syn_bio : ', syn_bio(layer) WRITE(numout,*) ' lys_bio : ', lys_bio(layer) WRITE(numout,*) ' rem_bio : ', rem_bio(layer) WRITE(numout,*) ' exu_bio : ', exu_bio(layer) WRITE(numout,*) ' rsp_bio : ', rsp_bio(layer) WRITE(numout,*) ' lim_dsi : ', lim_dsi(layer) WRITE(numout,*) ' lim_no3 : ', lim_no3(layer) WRITE(numout,*) ' lim_po4 : ', lim_po4(layer) WRITE(numout,*) ' lim_lig : ', lim_lig(layer) WRITE(numout,*) ' lim_tem : ', lim_tem(layer) WRITE(numout,*) ' lim_sal : ', lim_sal(layer) WRITE(numout,*) WRITE(numout,*) ' tc_bio : ', tc_bio(layer) WRITE(numout,*) ' sbr_bio : ', sbr_bio(layer) WRITE(numout,*) ' par_bio : ', par_bio(layer) WRITE(numout,*) WRITE(numout,*) ' zd_dsi : ', zd_dsi WRITE(numout,*) ' zd_no3 : ', zd_no3 WRITE(numout,*) ' zd_po4 : ', zd_po4 WRITE(numout,*) ' zd_daf : ', zd_daf WRITE(numout,*) ' zd_eoc : ', zd_eoc WRITE(numout,*) WRITE(numout,*) ' dSi : ', cbu_i_bio(jn_dsi,layer) WRITE(numout,*) ' dIN : ', cbu_i_bio(jn_din,layer) WRITE(numout,*) ' dIP : ', cbu_i_bio(jn_dip,layer) WRITE(numout,*) ' AoC : ', cbu_i_bio(jn_aoc,layer) WRITE(numout,*) ' eoC : ', cbu_i_bio(jn_eoc,layer) WRITE(numout,*) ENDIF ! ln_write_bio END DO ! layer !------------------------------------------------------------------------ ! 2.5 Update DIC, Alk and Oxy !------------------------------------------------------------------------ DO layer = layer_00, nlay_bio zncp = - syn_bio(layer) + rsp_bio(layer) + rem_bio(layer) IF ( ln_carbon .AND. flag_active(jn_dic) & .AND. flag_active(jn_alk) ) THEN zd_dic1 = - zncp * ddtb zd_dic2 = - cbu_i_bio(jn_dic,layer) zd_dic = MAX( zd_dic1 , zd_dic2 ) zd_alk1 = no3_c * zncp * ddtb zd_alk2 = - cbu_i_bio(jn_alk,layer) zd_alk = MAX( zd_alk1 , zd_alk2 ) cbu_i_bio(jn_dic,layer) = cbu_i_bio(jn_dic,layer) + zd_dic cbu_i_bio(jn_alk,layer) = cbu_i_bio(jn_alk,layer) + zd_alk IF ( ln_write_bio ) THEN WRITE(numout,*) ' DIC : ', cbu_i_bio(jn_dic,layer) WRITE(numout,*) ' Alk : ', cbu_i_bio(jn_alk,layer) ENDIF ENDIF IF ( flag_active(jn_oxy) ) THEN zd_oxy1 = oxy_c * zncp * ddtb zd_oxy2 = - cbu_i_bio(jn_oxy,layer) zd_oxy = MAX( zd_oxy1, zd_oxy2 ) cbu_i_bio(jn_oxy,layer) = cbu_i_bio(jn_oxy,layer) + zd_oxy IF ( ln_write_bio ) THEN WRITE(numout,*) ' Oxy : ', cbu_i_bio(jn_oxy,layer) ENDIF ENDIF END DO ! !------------------------------------------------------------------------------! ! X) End of the routine !------------------------------------------------------------------------------! ! WRITE(numout,*) WRITE(numout,*) ' *** After biological sources and sinks *** ' WRITE(numout,*) ' model output ' DO jn = 1, ntra_bio IF ( flag_active(jn) ) THEN WRITE(numout,*) ' biotr_i_nam : ', biotr_i_nam(jn) WRITE(numout,*) ' cbu_i_bio : ', ( cbu_i_bio(jn, jk), jk = 1, & nlay_bio ) ENDIF END DO WRITE(numout,*) WRITE(numout,*) WRITE(numout,*) ' End of ice_bio_sms ' WRITE(numout,*) '~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~' !==============================================================================| ! Fin de la routine ice_bio_sms WRITE(numout,*) END