- Timestamp:
- 2015-07-10T13:28:53+02:00 (9 years ago)
- Location:
- branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z
- Files:
-
- 17 edited
-
p4zbio.F90 (modified) (3 diffs)
-
p4zche.F90 (modified) (1 diff, 1 prop)
-
p4zfechem.F90 (modified) (13 diffs)
-
p4zflx.F90 (modified) (8 diffs, 1 prop)
-
p4zint.F90 (modified) (2 diffs, 1 prop)
-
p4zlim.F90 (modified) (7 diffs)
-
p4zlys.F90 (modified) (6 diffs)
-
p4zmeso.F90 (modified) (16 diffs)
-
p4zmicro.F90 (modified) (13 diffs)
-
p4zmort.F90 (modified) (8 diffs, 1 prop)
-
p4zopt.F90 (modified) (17 diffs)
-
p4zprod.F90 (modified) (14 diffs)
-
p4zrem.F90 (modified) (12 diffs)
-
p4zsbc.F90 (modified) (11 diffs, 1 prop)
-
p4zsed.F90 (modified) (23 diffs, 1 prop)
-
p4zsink.F90 (modified) (27 diffs)
-
p4zsms.F90 (modified) (17 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zbio.F90
r4793 r5581 44 44 CONTAINS 45 45 46 SUBROUTINE p4z_bio ( kt, jnt )46 SUBROUTINE p4z_bio ( kt, knt ) 47 47 !!--------------------------------------------------------------------- 48 48 !! *** ROUTINE p4z_bio *** … … 54 54 !! ** Method : - ??? 55 55 !!--------------------------------------------------------------------- 56 INTEGER, INTENT(in) :: kt, jnt 57 INTEGER :: ji, jj, jk, jn 58 REAL(wp) :: ztra 59 #if defined key_kriest 60 REAL(wp) :: zcoef1, zcoef2 61 #endif 56 INTEGER, INTENT(in) :: kt, knt 57 INTEGER :: ji, jj, jk, jn 62 58 CHARACTER (len=25) :: charout 63 59 … … 80 76 81 77 82 CALL p4z_opt ( kt, jnt ) ! Optic: PAR in the water column83 CALL p4z_sink ( kt, jnt ) ! vertical flux of particulate organic matter84 CALL p4z_fechem(kt, jnt ) ! Iron chemistry/scavenging85 CALL p4z_lim ( kt, jnt ) ! co-limitations by the various nutrients86 CALL p4z_prod ( kt, jnt ) ! phytoplankton growth rate over the global ocean.78 CALL p4z_opt ( kt, knt ) ! Optic: PAR in the water column 79 CALL p4z_sink ( kt, knt ) ! vertical flux of particulate organic matter 80 CALL p4z_fechem(kt, knt ) ! Iron chemistry/scavenging 81 CALL p4z_lim ( kt, knt ) ! co-limitations by the various nutrients 82 CALL p4z_prod ( kt, knt ) ! phytoplankton growth rate over the global ocean. 87 83 ! ! (for each element : C, Si, Fe, Chl ) 88 84 CALL p4z_mort ( kt ) ! phytoplankton mortality 89 90 CALL p4z_micro( kt, jnt ) ! microzooplankton91 CALL p4z_meso ( kt, jnt ) ! mesozooplankton92 CALL p4z_rem ( kt, jnt ) ! remineralization terms of organic matter+scavenging of Fe85 ! ! zooplankton sources/sinks routines 86 CALL p4z_micro( kt, knt ) ! microzooplankton 87 CALL p4z_meso ( kt, knt ) ! mesozooplankton 88 CALL p4z_rem ( kt, knt ) ! remineralization terms of organic matter+scavenging of Fe 93 89 ! ! test if tracers concentrations fall below 0. 94 xnegtr(:,:,:) = 1.e0 95 DO jn = jp_pcs0, jp_pcs1 96 DO jk = 1, jpk 97 DO jj = 1, jpj 98 DO ji = 1, jpi 99 IF( ( trn(ji,jj,jk,jn) + tra(ji,jj,jk,jn) ) < 0.e0 ) THEN 100 ztra = ABS( trn(ji,jj,jk,jn) ) / ( ABS( tra(ji,jj,jk,jn) ) + rtrn ) 101 102 xnegtr(ji,jj,jk) = MIN( xnegtr(ji,jj,jk), ztra ) 103 ENDIF 104 END DO 105 END DO 106 END DO 107 END DO 108 ! ! where at least 1 tracer concentration becomes negative 109 ! ! 110 DO jn = jp_pcs0, jp_pcs1 111 trn(:,:,:,jn) = trn(:,:,:,jn) + xnegtr(:,:,:) * tra(:,:,:,jn) 112 END DO 113 114 115 tra(:,:,:,:) = 0.e0 116 117 #if defined key_kriest 118 ! 119 zcoef1 = 1.e0 / xkr_massp 120 zcoef2 = 1.e0 / xkr_massp / 1.1 121 DO jk = 1,jpkm1 122 trn(:,:,jk,jpnum) = MAX( trn(:,:,jk,jpnum), trn(:,:,jk,jppoc) * zcoef1 / xnumm(jk) ) 123 trn(:,:,jk,jpnum) = MIN( trn(:,:,jk,jpnum), trn(:,:,jk,jppoc) * zcoef2 ) 124 END DO 125 #endif 126 127 ! 90 ! ! 128 91 IF(ln_ctl) THEN ! print mean trends (used for debugging) 129 92 WRITE(charout, FMT="('bio ')") 130 93 CALL prt_ctl_trc_info(charout) 131 CALL prt_ctl_trc(tab4d=tr n, mask=tmask, clinfo=ctrcnm)94 CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm) 132 95 ENDIF 133 96 ! -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zche.F90
- Property svn:keywords set to Id
r4793 r5581 168 168 !!---------------------------------------------------------------------- 169 169 !! NEMO/TOP 3.3 , NEMO Consortium (2010) 170 !! $Id : p4zche.F90 3294 2012-01-28 16:44:18Z rblod$170 !! $Id$ 171 171 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 172 172 !!---------------------------------------------------------------------- -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zfechem.F90
r4624 r5581 48 48 CONTAINS 49 49 50 SUBROUTINE p4z_fechem( kt, jnt )50 SUBROUTINE p4z_fechem( kt, knt ) 51 51 !!--------------------------------------------------------------------- 52 52 !! *** ROUTINE p4z_fechem *** … … 62 62 !!--------------------------------------------------------------------- 63 63 ! 64 INTEGER, INTENT(in) :: kt, jnt ! ocean time step64 INTEGER, INTENT(in) :: kt, knt ! ocean time step 65 65 ! 66 66 INTEGER :: ji, jj, jk, jic … … 101 101 ! ------------------------------------------------- 102 102 IF( ln_ligvar ) THEN 103 ztotlig(:,:,:) = 0.09 * tr n(:,:,:,jpdoc) * 1E6 + ligand * 1E9103 ztotlig(:,:,:) = 0.09 * trb(:,:,:,jpdoc) * 1E6 + ligand * 1E9 104 104 ztotlig(:,:,:) = MIN( ztotlig(:,:,:), 10. ) 105 105 ELSE … … 127 127 zionic = 19.9201 * tsn(ji,jj,jk,jp_sal) / ( 1000. - 1.00488 * tsn(ji,jj,jk,jp_sal) + rtrn ) 128 128 zph = -LOG10( MAX( hi(ji,jj,jk), rtrn) ) 129 zoxy = tr n(ji,jj,jk,jpoxy) * ( rhop(ji,jj,jk) / 1.e3 )129 zoxy = trb(ji,jj,jk,jpoxy) * ( rhop(ji,jj,jk) / 1.e3 ) 130 130 ! Fe2+ oxydation rate from Santana-Casiano et al. (2005) 131 131 zkox = 35.407 - 6.7109 * zph + 0.5342 * zph * zph - 5362.6 / ( tsn(ji,jj,jk,jp_tem) + 273.15 ) & … … 137 137 zkph1 = zkph2 / 5. 138 138 ! pass the dfe concentration from PISCES 139 ztfe = tr n(ji,jj,jk,jpfer) * 1e9139 ztfe = trb(ji,jj,jk,jpfer) * 1e9 140 140 ! ---------------------------------------------------------- 141 141 ! ANALYTICAL SOLUTION OF ROOTS OF THE FE3+ EQUATION … … 204 204 zkeq = fekeq(ji,jj,jk) 205 205 zfesatur = zTL1(ji,jj,jk) * 1E-9 206 ztfe = tr n(ji,jj,jk,jpfer)206 ztfe = trb(ji,jj,jk,jpfer) 207 207 ! Fe' is the root of a 2nd order polynom 208 208 zFe3 (ji,jj,jk) = ( -( 1. + zfesatur * zkeq - zkeq * ztfe ) & … … 210 210 & + 4. * ztfe * zkeq) ) / ( 2. * zkeq ) 211 211 zFe3 (ji,jj,jk) = zFe3(ji,jj,jk) * 1E9 212 zFeL1(ji,jj,jk) = MAX( 0., tr n(ji,jj,jk,jpfer) * 1E9 - zFe3(ji,jj,jk) )212 zFeL1(ji,jj,jk) = MAX( 0., trb(ji,jj,jk,jpfer) * 1E9 - zFe3(ji,jj,jk) ) 213 213 END DO 214 214 END DO … … 240 240 ENDIF 241 241 #if defined key_kriest 242 ztrc = ( tr n(ji,jj,jk,jppoc) + trn(ji,jj,jk,jpcal) + trn(ji,jj,jk,jpgsi) ) * 1.e6242 ztrc = ( trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpcal) + trb(ji,jj,jk,jpgsi) ) * 1.e6 243 243 #else 244 ztrc = ( tr n(ji,jj,jk,jppoc) + trn(ji,jj,jk,jpgoc) + trn(ji,jj,jk,jpcal) + trn(ji,jj,jk,jpgsi) ) * 1.e6244 ztrc = ( trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + trb(ji,jj,jk,jpcal) + trb(ji,jj,jk,jpgsi) ) * 1.e6 245 245 #endif 246 IF( ln_dust ) zdust = dust(ji,jj) / ( wdust *rday ) * tmask(ji,jj,jk) ! dust in kg/m2/s246 IF( ln_dust ) zdust = dust(ji,jj) / ( wdust / rday ) * tmask(ji,jj,jk) ! dust in kg/m2/s 247 247 zlam1b = 3.e-5 + xlamdust * zdust + xlam1 * ztrc 248 248 zscave = zfeequi * zlam1b * zstep … … 251 251 ! to later allocate scavenged iron to the different organic pools 252 252 ! --------------------------------------------------------- 253 zdenom1 = xlam1 * tr n(ji,jj,jk,jppoc) / zlam1b253 zdenom1 = xlam1 * trb(ji,jj,jk,jppoc) / zlam1b 254 254 #if ! defined key_kriest 255 zdenom2 = xlam1 * tr n(ji,jj,jk,jpgoc) / zlam1b255 zdenom2 = xlam1 * trb(ji,jj,jk,jpgoc) / zlam1b 256 256 #endif 257 257 … … 262 262 zlamfac = MIN( 1. , zlamfac ) 263 263 zdep = MIN( 1., 1000. / fsdept(ji,jj,jk) ) 264 zlam1b = xlam1 * MAX( 0.e0, ( tr n(ji,jj,jk,jpfer) * 1.e9 - ztotlig(ji,jj,jk) ) )265 zcoag = zfeequi * zlam1b * zstep + 1E-4 * ( 1. - zlamfac ) * zdep * zstep * tr n(ji,jj,jk,jpfer)264 zlam1b = xlam1 * MAX( 0.e0, ( trb(ji,jj,jk,jpfer) * 1.e9 - ztotlig(ji,jj,jk) ) ) 265 zcoag = zfeequi * zlam1b * zstep + 1E-4 * ( 1. - zlamfac ) * zdep * zstep * trb(ji,jj,jk,jpfer) 266 266 267 267 ! Compute the coagulation of colloidal iron. This parameterization … … 269 269 ! It requires certainly some more work as it is very poorly constrained. 270 270 ! ---------------------------------------------------------------- 271 zlam1a = ( 0.369 * 0.3 * tr n(ji,jj,jk,jpdoc) + 102.4 * trn(ji,jj,jk,jppoc) ) * xdiss(ji,jj,jk) &272 & + ( 114. * 0.3 * tr n(ji,jj,jk,jpdoc) + 5.09E3 * trn(ji,jj,jk,jppoc) )271 zlam1a = ( 0.369 * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4 * trb(ji,jj,jk,jppoc) ) * xdiss(ji,jj,jk) & 272 & + ( 114. * 0.3 * trb(ji,jj,jk,jpdoc) + 5.09E3 * trb(ji,jj,jk,jppoc) ) 273 273 zaggdfea = zlam1a * zstep * zfecoll 274 274 #if defined key_kriest … … 278 278 tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zscave * zdenom1 + zaggdfea + zaggdfeb 279 279 #else 280 zlam1b = 3.53E3 * tr n(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk)280 zlam1b = 3.53E3 * trb(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk) 281 281 zaggdfeb = zlam1b * zstep * zfecoll 282 282 ! … … 292 292 ! ---------------------------------------- 293 293 IF( ln_fechem ) THEN 294 biron(:,:,:) = MAX( 0., tr n(:,:,:,jpfer) - zFeP(:,:,:) * 1E-9 )294 biron(:,:,:) = MAX( 0., trb(:,:,:,jpfer) - zFeP(:,:,:) * 1E-9 ) 295 295 ELSE 296 biron(:,:,:) = tr n(:,:,:,jpfer)296 biron(:,:,:) = trb(:,:,:,jpfer) 297 297 ENDIF 298 298 299 299 ! Output of some diagnostics variables 300 300 ! --------------------------------- 301 IF( ln_diatrc .AND. lk_iomput ) THEN 302 IF( jnt == nrdttrc ) THEN 303 CALL iom_put("Fe3" , zFe3 (:,:,:) * tmask(:,:,:) ) ! Fe3+ 304 CALL iom_put("FeL1" , zFeL1 (:,:,:) * tmask(:,:,:) ) ! FeL1 305 CALL iom_put("TL1" , zTL1 (:,:,:) * tmask(:,:,:) ) ! TL1 306 CALL iom_put("Totlig" , ztotlig(:,:,:) * tmask(:,:,:) ) ! TL 307 CALL iom_put("Biron" , biron (:,:,:) * 1e9 * tmask(:,:,:) ) ! biron 308 IF( ln_fechem ) THEN 309 CALL iom_put("Fe2" , zFe2 (:,:,:) * tmask(:,:,:) ) ! Fe2+ 310 CALL iom_put("FeL2", zFeL2 (:,:,:) * tmask(:,:,:) ) ! FeL2 311 CALL iom_put("FeP" , zFeP (:,:,:) * tmask(:,:,:) ) ! FeP 312 CALL iom_put("TL2" , zTL2 (:,:,:) * tmask(:,:,:) ) ! TL2 313 ENDIF 301 IF( lk_iomput .AND. knt == nrdttrc ) THEN 302 IF( iom_use("Fe3") ) CALL iom_put("Fe3" , zFe3 (:,:,:) * tmask(:,:,:) ) ! Fe3+ 303 IF( iom_use("FeL1") ) CALL iom_put("FeL1" , zFeL1 (:,:,:) * tmask(:,:,:) ) ! FeL1 304 IF( iom_use("TL1") ) CALL iom_put("TL1" , zTL1 (:,:,:) * tmask(:,:,:) ) ! TL1 305 IF( iom_use("Totlig") ) CALL iom_put("Totlig" , ztotlig(:,:,:) * tmask(:,:,:) ) ! TL 306 IF( iom_use("Biron") ) CALL iom_put("Biron" , biron (:,:,:) * 1e9 * tmask(:,:,:) ) ! biron 307 IF( ln_fechem ) THEN 308 IF( iom_use("Fe2") ) CALL iom_put("Fe2" , zFe2 (:,:,:) * tmask(:,:,:) ) ! Fe2+ 309 IF( iom_use("FeL2") ) CALL iom_put("FeL2" , zFeL2 (:,:,:) * tmask(:,:,:) ) ! FeL2 310 IF( iom_use("FeP") ) CALL iom_put("FeP" , zFeP (:,:,:) * tmask(:,:,:) ) ! FeP 311 IF( iom_use("TL2") ) CALL iom_put("TL2" , zTL2 (:,:,:) * tmask(:,:,:) ) ! TL2 314 312 ENDIF 315 313 ENDIF -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zflx.F90
- Property svn:keywords set to Id
r4793 r5581 63 63 !!---------------------------------------------------------------------- 64 64 !! NEMO/TOP 3.3 , NEMO Consortium (2010) 65 !! $Id : p4zflx.F90 3294 2012-01-28 16:44:18Z rblod$65 !! $Id$ 66 66 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 67 67 !!---------------------------------------------------------------------- 68 68 CONTAINS 69 69 70 SUBROUTINE p4z_flx ( kt )70 SUBROUTINE p4z_flx ( kt, knt ) 71 71 !!--------------------------------------------------------------------- 72 72 !! *** ROUTINE p4z_flx *** … … 81 81 !!--------------------------------------------------------------------- 82 82 ! 83 INTEGER, INTENT(in) :: kt !83 INTEGER, INTENT(in) :: kt, knt ! 84 84 ! 85 85 INTEGER :: ji, jj, jm, iind, iindm1 … … 89 89 REAL(wp) :: zyr_dec, zdco2dt 90 90 CHARACTER (len=25) :: charout 91 REAL(wp), POINTER, DIMENSION(:,:) :: zkgco2, zkgo2, zh2co3, zoflx 91 REAL(wp), POINTER, DIMENSION(:,:) :: zkgco2, zkgo2, zh2co3, zoflx, zw2d 92 92 !!--------------------------------------------------------------------- 93 93 ! … … 101 101 ! IS USED TO COMPUTE AIR-SEA FLUX OF CO2 102 102 103 IF( kt /= nit000 ) CALL p4z_patm( kt ) ! Get sea-level pressure (E&K [1981] climatology) for use in flux calcs103 IF( kt /= nit000 .AND. knt == 1 ) CALL p4z_patm( kt ) ! Get sea-level pressure (E&K [1981] climatology) for use in flux calcs 104 104 105 105 IF( ln_co2int ) THEN … … 130 130 zbot = borat(ji,jj,1) 131 131 zfact = rhop(ji,jj,1) / 1000. + rtrn 132 zdic = tr n(ji,jj,1,jpdic) / zfact132 zdic = trb(ji,jj,1,jpdic) / zfact 133 133 zph = MAX( hi(ji,jj,1), 1.e-10 ) / zfact 134 zalka = tr n(ji,jj,1,jptal) / zfact134 zalka = trb(ji,jj,1,jptal) / zfact 135 135 136 136 ! CALCULATE [ALK]([CO3--], [HCO3-]) … … 184 184 zfld = satmco2(ji,jj) * patm(ji,jj) * tmask(ji,jj,1) * chemc(ji,jj,1) * zkgco2(ji,jj) ! (mol/L) * (m/s) 185 185 zflu = zh2co3(ji,jj) * tmask(ji,jj,1) * zkgco2(ji,jj) ! (mol/L) (m/s) ? 186 oce_co2(ji,jj) = ( zfld - zflu ) * rfact * e1e2t(ji,jj) * tmask(ji,jj,1) * 1000.186 oce_co2(ji,jj) = ( zfld - zflu ) * rfact2 * e1e2t(ji,jj) * tmask(ji,jj,1) * 1000. 187 187 ! compute the trend 188 tra(ji,jj,1,jpdic) = tra(ji,jj,1,jpdic) + ( zfld - zflu ) / fse3t(ji,jj,1)188 tra(ji,jj,1,jpdic) = tra(ji,jj,1,jpdic) + ( zfld - zflu ) * rfact2 / fse3t(ji,jj,1) 189 189 190 190 ! Compute O2 flux 191 191 zfld16 = atcox * patm(ji,jj) * chemc(ji,jj,2) * tmask(ji,jj,1) * zkgo2(ji,jj) ! (mol/L) * (m/s) 192 zflu16 = tr n(ji,jj,1,jpoxy) * tmask(ji,jj,1) * zkgo2(ji,jj)192 zflu16 = trb(ji,jj,1,jpoxy) * tmask(ji,jj,1) * zkgo2(ji,jj) 193 193 zoflx(ji,jj) = zfld16 - zflu16 194 tra(ji,jj,1,jpoxy) = tra(ji,jj,1,jpoxy) + zoflx(ji,jj) / fse3t(ji,jj,1)194 tra(ji,jj,1,jpoxy) = tra(ji,jj,1,jpoxy) + zoflx(ji,jj) * rfact2 / fse3t(ji,jj,1) 195 195 END DO 196 196 END DO 197 197 198 t_oce_co2_flx = t_oce_co2_flx + glob_sum( oce_co2(:,:) ) ! Cumulative Total Flux of Carbon 199 t_atm_co2_flx = glob_sum( satmco2(:,:) * e1e2t(:,:) ) ! Total atmospheric pCO2 200 198 t_oce_co2_flx = glob_sum( oce_co2(:,:) ) ! Total Flux of Carbon 199 t_oce_co2_flx_cum = t_oce_co2_flx_cum + t_oce_co2_flx ! Cumulative Total Flux of Carbon 200 ! t_atm_co2_flx = glob_sum( satmco2(:,:) * e1e2t(:,:) ) ! Total atmospheric pCO2 201 t_atm_co2_flx = atcco2 ! Total atmospheric pCO2 202 201 203 IF(ln_ctl) THEN ! print mean trends (used for debugging) 202 204 WRITE(charout, FMT="('flx ')") … … 205 207 ENDIF 206 208 207 IF( ln_diatrc ) THEN 208 IF( lk_iomput ) THEN 209 CALL iom_put( "Cflx" , oce_co2(:,:) / e1e2t(:,:) / rfact ) 210 CALL iom_put( "Oflx" , zoflx(:,:) * 1000 * tmask(:,:,1) ) 211 CALL iom_put( "Kg" , zkgco2(:,:) * tmask(:,:,1) ) 212 CALL iom_put( "Dpco2", ( satmco2(:,:) * patm(:,:) - zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) ) 213 CALL iom_put( "Dpo2" , ( atcox * patm(:,:) - trn(:,:,1,jpoxy) / ( chemc(:,:,2) + rtrn ) ) * tmask(:,:,1) ) 214 ELSE 215 trc2d(:,:,jp_pcs0_2d ) = oce_co2(:,:) / e1e2t(:,:) / rfact 209 IF( lk_iomput .AND. knt == nrdttrc ) THEN 210 CALL wrk_alloc( jpi, jpj, zw2d ) 211 IF( iom_use( "Cflx" ) ) THEN 212 zw2d(:,:) = oce_co2(:,:) / e1e2t(:,:) * rfact2r 213 CALL iom_put( "Cflx" , zw2d ) 214 ENDIF 215 IF( iom_use( "Oflx" ) ) THEN 216 zw2d(:,:) = zoflx(:,:) * 1000 * tmask(:,:,1) 217 CALL iom_put( "Oflx" , zw2d ) 218 ENDIF 219 IF( iom_use( "Kg" ) ) THEN 220 zw2d(:,:) = zkgco2(:,:) * tmask(:,:,1) 221 CALL iom_put( "Kg" , zw2d ) 222 ENDIF 223 IF( iom_use( "Dpco2" ) ) THEN 224 zw2d(:,:) = ( satmco2(:,:) * patm(:,:) - zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) 225 CALL iom_put( "Dpco2" , zw2d ) 226 ENDIF 227 IF( iom_use( "Dpo2" ) ) THEN 228 zw2d(:,:) = ( atcox * patm(:,:) - trb(:,:,1,jpoxy) / ( chemc(:,:,2) + rtrn ) ) * tmask(:,:,1) 229 CALL iom_put( "Dpo2" , zw2d ) 230 ENDIF 231 IF( iom_use( "tcflx" ) ) CALL iom_put( "tcflx" , t_oce_co2_flx * rfact2r ) ! molC/s 232 CALL iom_put( "tcflxcum" , t_oce_co2_flx_cum ) ! molC 233 ! 234 CALL wrk_dealloc( jpi, jpj, zw2d ) 235 ELSE 236 IF( ln_diatrc ) THEN 237 trc2d(:,:,jp_pcs0_2d ) = oce_co2(:,:) / e1e2t(:,:) * rfact2r 216 238 trc2d(:,:,jp_pcs0_2d + 1) = zoflx(:,:) * 1000 * tmask(:,:,1) 217 239 trc2d(:,:,jp_pcs0_2d + 2) = zkgco2(:,:) * tmask(:,:,1) … … 290 312 ! 291 313 oce_co2(:,:) = 0._wp ! Initialization of Flux of Carbon 314 t_oce_co2_flx = 0._wp 292 315 t_atm_co2_flx = 0._wp 293 t_oce_co2_flx = 0._wp294 316 ! 295 317 CALL p4z_patm( nit000 ) -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zint.F90
- Property svn:keywords set to Id
r4793 r5581 26 26 !!---------------------------------------------------------------------- 27 27 !! NEMO/TOP 3.3 , NEMO Consortium (2010) 28 !! $Id : p4zint.F90 3294 2012-01-28 16:44:18Z rblod$28 !! $Id$ 29 29 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 30 30 !!---------------------------------------------------------------------- … … 56 56 DO ji = 1, jpi 57 57 DO jj = 1, jpj 58 zvar = tr n(ji,jj,1,jpsil) * trn(ji,jj,1,jpsil)58 zvar = trb(ji,jj,1,jpsil) * trb(ji,jj,1,jpsil) 59 59 xksimax(ji,jj) = MAX( xksimax(ji,jj), ( 1.+ 7.* zvar / ( xksilim * xksilim + zvar ) ) * 1e-6 ) 60 60 END DO -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zlim.F90
r4793 r5581 62 62 CONTAINS 63 63 64 SUBROUTINE p4z_lim( kt, jnt )64 SUBROUTINE p4z_lim( kt, knt ) 65 65 !!--------------------------------------------------------------------- 66 66 !! *** ROUTINE p4z_lim *** … … 72 72 !!--------------------------------------------------------------------- 73 73 ! 74 INTEGER, INTENT(in) :: kt, jnt74 INTEGER, INTENT(in) :: kt, knt 75 75 ! 76 76 INTEGER :: ji, jj, jk 77 77 REAL(wp) :: zlim1, zlim2, zlim3, zlim4, zno3, zferlim 78 78 REAL(wp) :: zconcd, zconcd2, zconcn, zconcn2 79 REAL(wp) :: z1_tr ndia, z1_trnphy, ztem1, ztem2, zetot1, zetot279 REAL(wp) :: z1_trbdia, z1_trbphy, ztem1, ztem2, zetot1, zetot2 80 80 REAL(wp) :: zdenom, zratio, zironmin 81 81 REAL(wp) :: zconc1d, zconc1dnh4, zconc0n, zconc0nnh4 … … 90 90 ! Tuning of the iron concentration to a minimum level that is set to the detection limit 91 91 !------------------------------------- 92 zno3 = tr n(ji,jj,jk,jpno3) / 40.e-692 zno3 = trb(ji,jj,jk,jpno3) / 40.e-6 93 93 zferlim = MAX( 3e-11 * zno3 * zno3, 5e-12 ) 94 94 zferlim = MIN( zferlim, 7e-11 ) 95 tr n(ji,jj,jk,jpfer) = MAX( trn(ji,jj,jk,jpfer), zferlim )95 trb(ji,jj,jk,jpfer) = MAX( trb(ji,jj,jk,jpfer), zferlim ) 96 96 97 97 ! Computation of a variable Ks for iron on diatoms taking into account 98 98 ! that increasing biomass is made of generally bigger cells 99 99 !------------------------------------------------ 100 zconcd = MAX( 0.e0 , tr n(ji,jj,jk,jpdia) - xsizedia )101 zconcd2 = tr n(ji,jj,jk,jpdia) - zconcd102 zconcn = MAX( 0.e0 , tr n(ji,jj,jk,jpphy) - xsizephy )103 zconcn2 = tr n(ji,jj,jk,jpphy) - zconcn104 z1_tr nphy = 1. / ( trn(ji,jj,jk,jpphy) + rtrn )105 z1_tr ndia = 1. / ( trn(ji,jj,jk,jpdia) + rtrn )106 107 concdfe(ji,jj,jk) = MAX( concdfer, ( zconcd2 * concdfer + concdfer * xsizerd * zconcd ) * z1_tr ndia )108 zconc1d = MAX( concdno3, ( zconcd2 * concdno3 + concdno3 * xsizerd * zconcd ) * z1_tr ndia )109 zconc1dnh4 = MAX( concdnh4, ( zconcd2 * concdnh4 + concdnh4 * xsizerd * zconcd ) * z1_tr ndia )110 111 concnfe(ji,jj,jk) = MAX( concnfer, ( zconcn2 * concnfer + concnfer * xsizern * zconcn ) * z1_tr nphy )112 zconc0n = MAX( concnno3, ( zconcn2 * concnno3 + concnno3 * xsizern * zconcn ) * z1_tr nphy )113 zconc0nnh4 = MAX( concnnh4, ( zconcn2 * concnnh4 + concnnh4 * xsizern * zconcn ) * z1_tr nphy )100 zconcd = MAX( 0.e0 , trb(ji,jj,jk,jpdia) - xsizedia ) 101 zconcd2 = trb(ji,jj,jk,jpdia) - zconcd 102 zconcn = MAX( 0.e0 , trb(ji,jj,jk,jpphy) - xsizephy ) 103 zconcn2 = trb(ji,jj,jk,jpphy) - zconcn 104 z1_trbphy = 1. / ( trb(ji,jj,jk,jpphy) + rtrn ) 105 z1_trbdia = 1. / ( trb(ji,jj,jk,jpdia) + rtrn ) 106 107 concdfe(ji,jj,jk) = MAX( concdfer, ( zconcd2 * concdfer + concdfer * xsizerd * zconcd ) * z1_trbdia ) 108 zconc1d = MAX( concdno3, ( zconcd2 * concdno3 + concdno3 * xsizerd * zconcd ) * z1_trbdia ) 109 zconc1dnh4 = MAX( concdnh4, ( zconcd2 * concdnh4 + concdnh4 * xsizerd * zconcd ) * z1_trbdia ) 110 111 concnfe(ji,jj,jk) = MAX( concnfer, ( zconcn2 * concnfer + concnfer * xsizern * zconcn ) * z1_trbphy ) 112 zconc0n = MAX( concnno3, ( zconcn2 * concnno3 + concnno3 * xsizern * zconcn ) * z1_trbphy ) 113 zconc0nnh4 = MAX( concnnh4, ( zconcn2 * concnnh4 + concnnh4 * xsizern * zconcn ) * z1_trbphy ) 114 114 115 115 ! Michaelis-Menten Limitation term for nutrients Small bacteria 116 116 ! ------------------------------------------------------------- 117 zdenom = 1. / ( concbno3 * concbnh4 + concbnh4 * tr n(ji,jj,jk,jpno3) + concbno3 * trn(ji,jj,jk,jpnh4) )118 xnanono3(ji,jj,jk) = tr n(ji,jj,jk,jpno3) * concbnh4 * zdenom119 xnanonh4(ji,jj,jk) = tr n(ji,jj,jk,jpnh4) * concbno3 * zdenom117 zdenom = 1. / ( concbno3 * concbnh4 + concbnh4 * trb(ji,jj,jk,jpno3) + concbno3 * trb(ji,jj,jk,jpnh4) ) 118 xnanono3(ji,jj,jk) = trb(ji,jj,jk,jpno3) * concbnh4 * zdenom 119 xnanonh4(ji,jj,jk) = trb(ji,jj,jk,jpnh4) * concbno3 * zdenom 120 120 ! 121 121 zlim1 = xnanono3(ji,jj,jk) + xnanonh4(ji,jj,jk) 122 zlim2 = tr n(ji,jj,jk,jppo4) / ( trn(ji,jj,jk,jppo4) + concbnh4 )123 zlim3 = tr n(ji,jj,jk,jpfer) / ( concbfe + trn(ji,jj,jk,jpfer) )124 zlim4 = tr n(ji,jj,jk,jpdoc) / ( xkdoc + trn(ji,jj,jk,jpdoc) )122 zlim2 = trb(ji,jj,jk,jppo4) / ( trb(ji,jj,jk,jppo4) + concbnh4 ) 123 zlim3 = trb(ji,jj,jk,jpfer) / ( concbfe + trb(ji,jj,jk,jpfer) ) 124 zlim4 = trb(ji,jj,jk,jpdoc) / ( xkdoc + trb(ji,jj,jk,jpdoc) ) 125 125 xlimbacl(ji,jj,jk) = MIN( zlim1, zlim2, zlim3 ) 126 126 xlimbac (ji,jj,jk) = MIN( zlim1, zlim2, zlim3 ) * zlim4 … … 128 128 ! Michaelis-Menten Limitation term for nutrients Small flagellates 129 129 ! ----------------------------------------------- 130 zdenom = 1. / ( zconc0n * zconc0nnh4 + zconc0nnh4 * tr n(ji,jj,jk,jpno3) + zconc0n * trn(ji,jj,jk,jpnh4) )131 xnanono3(ji,jj,jk) = tr n(ji,jj,jk,jpno3) * zconc0nnh4 * zdenom132 xnanonh4(ji,jj,jk) = tr n(ji,jj,jk,jpnh4) * zconc0n * zdenom130 zdenom = 1. / ( zconc0n * zconc0nnh4 + zconc0nnh4 * trb(ji,jj,jk,jpno3) + zconc0n * trb(ji,jj,jk,jpnh4) ) 131 xnanono3(ji,jj,jk) = trb(ji,jj,jk,jpno3) * zconc0nnh4 * zdenom 132 xnanonh4(ji,jj,jk) = trb(ji,jj,jk,jpnh4) * zconc0n * zdenom 133 133 ! 134 134 zlim1 = xnanono3(ji,jj,jk) + xnanonh4(ji,jj,jk) 135 zlim2 = tr n(ji,jj,jk,jppo4) / ( trn(ji,jj,jk,jppo4) + zconc0nnh4 )136 zratio = tr n(ji,jj,jk,jpnfe) * z1_trnphy137 zironmin = xcoef1 * tr n(ji,jj,jk,jpnch) * z1_trnphy + xcoef2 * zlim1 + xcoef3 * xnanono3(ji,jj,jk)135 zlim2 = trb(ji,jj,jk,jppo4) / ( trb(ji,jj,jk,jppo4) + zconc0nnh4 ) 136 zratio = trb(ji,jj,jk,jpnfe) * z1_trbphy 137 zironmin = xcoef1 * trb(ji,jj,jk,jpnch) * z1_trbphy + xcoef2 * zlim1 + xcoef3 * xnanono3(ji,jj,jk) 138 138 zlim3 = MAX( 0.,( zratio - zironmin ) / qnfelim ) 139 139 xnanopo4(ji,jj,jk) = zlim2 … … 143 143 ! Michaelis-Menten Limitation term for nutrients Diatoms 144 144 ! ---------------------------------------------- 145 zdenom = 1. / ( zconc1d * zconc1dnh4 + zconc1dnh4 * tr n(ji,jj,jk,jpno3) + zconc1d * trn(ji,jj,jk,jpnh4) )146 xdiatno3(ji,jj,jk) = tr n(ji,jj,jk,jpno3) * zconc1dnh4 * zdenom147 xdiatnh4(ji,jj,jk) = tr n(ji,jj,jk,jpnh4) * zconc1d * zdenom145 zdenom = 1. / ( zconc1d * zconc1dnh4 + zconc1dnh4 * trb(ji,jj,jk,jpno3) + zconc1d * trb(ji,jj,jk,jpnh4) ) 146 xdiatno3(ji,jj,jk) = trb(ji,jj,jk,jpno3) * zconc1dnh4 * zdenom 147 xdiatnh4(ji,jj,jk) = trb(ji,jj,jk,jpnh4) * zconc1d * zdenom 148 148 ! 149 149 zlim1 = xdiatno3(ji,jj,jk) + xdiatnh4(ji,jj,jk) 150 zlim2 = tr n(ji,jj,jk,jppo4) / ( trn(ji,jj,jk,jppo4) + zconc1dnh4 )151 zlim3 = tr n(ji,jj,jk,jpsil) / ( trn(ji,jj,jk,jpsil) + xksi(ji,jj) )152 zratio = tr n(ji,jj,jk,jpdfe) * z1_trndia153 zironmin = xcoef1 * tr n(ji,jj,jk,jpdch) * z1_trndia + xcoef2 * zlim1 + xcoef3 * xdiatno3(ji,jj,jk)150 zlim2 = trb(ji,jj,jk,jppo4) / ( trb(ji,jj,jk,jppo4) + zconc1dnh4 ) 151 zlim3 = trb(ji,jj,jk,jpsil) / ( trb(ji,jj,jk,jpsil) + xksi(ji,jj) ) 152 zratio = trb(ji,jj,jk,jpdfe) * z1_trbdia 153 zironmin = xcoef1 * trb(ji,jj,jk,jpdch) * z1_trbdia + xcoef2 * zlim1 + xcoef3 * xdiatno3(ji,jj,jk) 154 154 zlim4 = MAX( 0., ( zratio - zironmin ) / qdfelim ) 155 155 xdiatpo4(ji,jj,jk) = zlim2 … … 166 166 DO jj = 1, jpj 167 167 DO ji = 1, jpi 168 zlim1 = ( tr n(ji,jj,jk,jpno3) * concnnh4 + trn(ji,jj,jk,jpnh4) * concnno3 ) &169 & / ( concnno3 * concnnh4 + concnnh4 * tr n(ji,jj,jk,jpno3) + concnno3 * trn(ji,jj,jk,jpnh4) )170 zlim2 = tr n(ji,jj,jk,jppo4) / ( trn(ji,jj,jk,jppo4) + concnnh4 )171 zlim3 = tr n(ji,jj,jk,jpfer) / ( trn(ji,jj,jk,jpfer) + 5.E-11 )168 zlim1 = ( trb(ji,jj,jk,jpno3) * concnnh4 + trb(ji,jj,jk,jpnh4) * concnno3 ) & 169 & / ( concnno3 * concnnh4 + concnnh4 * trb(ji,jj,jk,jpno3) + concnno3 * trb(ji,jj,jk,jpnh4) ) 170 zlim2 = trb(ji,jj,jk,jppo4) / ( trb(ji,jj,jk,jppo4) + concnnh4 ) 171 zlim3 = trb(ji,jj,jk,jpfer) / ( trb(ji,jj,jk,jpfer) + 5.E-11 ) 172 172 ztem1 = MAX( 0., tsn(ji,jj,jk,jp_tem) ) 173 173 ztem2 = tsn(ji,jj,jk,jp_tem) - 10. 174 zetot1 = MAX( 0., etot (ji,jj,jk) - 1.) / ( 4. + etot(ji,jj,jk) )175 zetot2 = 30. / ( 30. + etot (ji,jj,jk) )174 zetot1 = MAX( 0., etot_ndcy(ji,jj,jk) - 1.) / ( 4. + etot_ndcy(ji,jj,jk) ) 175 zetot2 = 30. / ( 30. + etot_ndcy(ji,jj,jk) ) 176 176 177 177 xfracal(ji,jj,jk) = caco3r * MIN( zlim1, zlim2, zlim3 ) & 178 178 & * ztem1 / ( 0.1 + ztem1 ) & 179 & * MAX( 1., tr n(ji,jj,jk,jpphy) * 1.e6 / 2. ) &179 & * MAX( 1., trb(ji,jj,jk,jpphy) * 1.e6 / 2. ) & 180 180 & * zetot1 * zetot2 & 181 181 & * ( 1. + EXP(-ztem2 * ztem2 / 25. ) ) & … … 187 187 END DO 188 188 ! 189 IF( ln_diatrc .AND. lk_iomput .AND. jnt == nrdttrc ) THEN ! save output diagnostics 190 ! 191 CALL iom_put( "xfracal", xfracal(:,:,:) * tmask(:,:,:) ) ! euphotic layer deptht 192 CALL iom_put( "LNnut" , xlimphy(:,:,:) * tmask(:,:,:) ) ! Nutrient limitation term 193 CALL iom_put( "LDnut" , xlimdia(:,:,:) * tmask(:,:,:) ) ! Nutrient limitation term 194 CALL iom_put( "LNFe" , xlimnfe(:,:,:) * tmask(:,:,:) ) ! Iron limitation term 195 CALL iom_put( "LDFe" , xlimdfe(:,:,:) * tmask(:,:,:) ) ! Iron limitation term 196 ! 189 ! 190 IF( lk_iomput .AND. knt == nrdttrc ) THEN ! save output diagnostics 191 IF( iom_use( "xfracal" ) ) CALL iom_put( "xfracal", xfracal(:,:,:) * tmask(:,:,:) ) ! euphotic layer deptht 192 IF( iom_use( "LNnut" ) ) CALL iom_put( "LNnut" , xlimphy(:,:,:) * tmask(:,:,:) ) ! Nutrient limitation term 193 IF( iom_use( "LDnut" ) ) CALL iom_put( "LDnut" , xlimdia(:,:,:) * tmask(:,:,:) ) ! Nutrient limitation term 194 IF( iom_use( "LNFe" ) ) CALL iom_put( "LNFe" , xlimnfe(:,:,:) * tmask(:,:,:) ) ! Iron limitation term 195 IF( iom_use( "LDFe" ) ) CALL iom_put( "LDFe" , xlimdfe(:,:,:) * tmask(:,:,:) ) ! Iron limitation term 197 196 ENDIF 198 199 197 ! 200 198 IF( nn_timing == 1 ) CALL timing_stop('p4z_lim') -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zlys.F90
r4793 r5581 48 48 CONTAINS 49 49 50 SUBROUTINE p4z_lys( kt )50 SUBROUTINE p4z_lys( kt, knt ) 51 51 !!--------------------------------------------------------------------- 52 52 !! *** ROUTINE p4z_lys *** … … 59 59 !!--------------------------------------------------------------------- 60 60 ! 61 INTEGER, INTENT(in) :: kt ! ocean time step61 INTEGER, INTENT(in) :: kt, knt ! ocean time step 62 62 INTEGER :: ji, jj, jk, jn 63 63 REAL(wp) :: zalk, zdic, zph, zah2 64 64 REAL(wp) :: zdispot, zfact, zcalcon, zalka, zaldi 65 65 REAL(wp) :: zomegaca, zexcess, zexcess0 66 REAL(wp) :: zrfact267 66 CHARACTER (len=25) :: charout 68 67 REAL(wp), POINTER, DIMENSION(:,:,:) :: zco3, zcaldiss … … 89 88 zfact = rhop(ji,jj,jk) / 1000. + rtrn 90 89 zph = hi(ji,jj,jk) * tmask(ji,jj,jk) / zfact + ( 1.-tmask(ji,jj,jk) ) * 1.e-9 ! [H+] 91 zdic = tr n(ji,jj,jk,jpdic) / zfact92 zalka = tr n(ji,jj,jk,jptal) / zfact90 zdic = trb(ji,jj,jk,jpdic) / zfact 91 zalka = trb(ji,jj,jk,jptal) / zfact 93 92 ! CALCULATE [ALK]([CO3--], [HCO3-]) 94 93 zalk = zalka - ( akw3(ji,jj,jk) / zph - zph + borat(ji,jj,jk) / ( 1. + zph / akb3(ji,jj,jk) ) ) … … 130 129 ! (ACCORDING TO THIS FORMULATION ALSO SOME PARTICULATE 131 130 ! CACO3 GETS DISSOLVED EVEN IN THE CASE OF OVERSATURATION) 132 zdispot = kdca * zexcess * tr n(ji,jj,jk,jpcal)131 zdispot = kdca * zexcess * trb(ji,jj,jk,jpcal) 133 132 # if defined key_degrad 134 133 zdispot = zdispot * facvol(ji,jj,jk) … … 136 135 ! CHANGE OF [CO3--] , [ALK], PARTICULATE [CACO3], 137 136 ! AND [SUM(CO2)] DUE TO CACO3 DISSOLUTION/PRECIPITATION 138 zcaldiss(ji,jj,jk) = zdispot / rmtss! calcite dissolution139 zco3(ji,jj,jk) = zco3(ji,jj,jk) + zcaldiss(ji,jj,jk) * rfact137 zcaldiss(ji,jj,jk) = zdispot * rfact2 / rmtss ! calcite dissolution 138 zco3(ji,jj,jk) = zco3(ji,jj,jk) + zcaldiss(ji,jj,jk) 140 139 ! 141 140 tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + 2. * zcaldiss(ji,jj,jk) … … 146 145 END DO 147 146 ! 148 IF( ln_diatrc ) THEN 149 ! 150 IF( lk_iomput ) THEN 151 zrfact2 = 1.e3 * rfact2r 152 CALL iom_put( "PH" , -1. * LOG10( hi(:,:,:) ) * tmask(:,:,:) ) 153 CALL iom_put( "CO3" , zco3 (:,:,:) * 1e+3 * tmask(:,:,:) ) 154 CALL iom_put( "CO3sat", aksp (:,:,:) * 1e+3 / calcon * tmask(:,:,:) ) 155 CALL iom_put( "DCAL" , zcaldiss(:,:,:) * zrfact2 * tmask(:,:,:) ) 156 ELSE 157 trc3d(:,:,:,jp_pcs0_3d ) = -1. * LOG10( hi(:,:,:) ) * tmask(:,:,:) 158 trc3d(:,:,:,jp_pcs0_3d + 1) = zco3(:,:,:) * tmask(:,:,:) 159 trc3d(:,:,:,jp_pcs0_3d + 2) = aksp(:,:,:) / calcon * tmask(:,:,:) 160 ENDIF 161 ! 147 148 IF( lk_iomput .AND. knt == nrdttrc ) THEN 149 IF( iom_use( "PH" ) ) CALL iom_put( "PH" , -1. * LOG10( hi(:,:,:) ) * tmask(:,:,:) ) 150 IF( iom_use( "CO3" ) ) CALL iom_put( "CO3" , zco3(:,:,:) * 1.e+3 * tmask(:,:,:) ) 151 IF( iom_use( "CO3sat" ) ) CALL iom_put( "CO3sat", aksp(:,:,:) * 1.e+3 / calcon * tmask(:,:,:) ) 152 IF( iom_use( "DCAL" ) ) CALL iom_put( "DCAL" , zcaldiss(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:) ) 153 ELSE 154 trc3d(:,:,:,jp_pcs0_3d ) = -1. * LOG10( hi(:,:,:) ) * tmask(:,:,:) 155 trc3d(:,:,:,jp_pcs0_3d + 1) = zco3(:,:,:) * tmask(:,:,:) 156 trc3d(:,:,:,jp_pcs0_3d + 2) = aksp(:,:,:) / calcon * tmask(:,:,:) 162 157 ENDIF 163 158 ! -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zmeso.F90
r4793 r5581 60 60 CONTAINS 61 61 62 SUBROUTINE p4z_meso( kt, jnt )62 SUBROUTINE p4z_meso( kt, knt ) 63 63 !!--------------------------------------------------------------------- 64 64 !! *** ROUTINE p4z_meso *** … … 68 68 !! ** Method : - ??? 69 69 !!--------------------------------------------------------------------- 70 INTEGER, INTENT(in) :: kt, jnt ! ocean time step70 INTEGER, INTENT(in) :: kt, knt ! ocean time step 71 71 INTEGER :: ji, jj, jk 72 72 REAL(wp) :: zcompadi, zcompaph, zcompapoc, zcompaz, zcompam … … 83 83 REAL(wp) :: zgrazfffp, zgrazfffg, zgrazffep, zgrazffeg 84 84 CHARACTER (len=25) :: charout 85 REAL(wp) :: zrfact2 86 REAL(wp), POINTER, DIMENSION(:,:,:) :: zgrazing 85 REAL(wp), POINTER, DIMENSION(:,:,:) :: zgrazing, zw3d 87 86 88 87 !!--------------------------------------------------------------------- … … 90 89 IF( nn_timing == 1 ) CALL timing_start('p4z_meso') 91 90 ! 92 IF( l n_diatrc .AND. lk_iomput ) THEN91 IF( lk_iomput ) THEN 93 92 CALL wrk_alloc( jpi, jpj, jpk, zgrazing ) 94 93 zgrazing(:,:,:) = 0._wp … … 98 97 DO jj = 1, jpj 99 98 DO ji = 1, jpi 100 zcompam = MAX( ( tr n(ji,jj,jk,jpmes) - 1.e-9 ), 0.e0 )99 zcompam = MAX( ( trb(ji,jj,jk,jpmes) - 1.e-9 ), 0.e0 ) 101 100 # if defined key_degrad 102 101 zstep = xstep * facvol(ji,jj,jk) … … 108 107 ! Respiration rates of both zooplankton 109 108 ! ------------------------------------- 110 zrespz2 = resrat2 * zfact * tr n(ji,jj,jk,jpmes) / ( xkmort + trn(ji,jj,jk,jpmes) ) &109 zrespz2 = resrat2 * zfact * trb(ji,jj,jk,jpmes) / ( xkmort + trb(ji,jj,jk,jpmes) ) & 111 110 & + resrat2 * zfact * 3. * nitrfac(ji,jj,jk) 112 111 … … 114 113 ! no real reason except that it seems to be more stable and may mimic predation 115 114 ! --------------------------------------------------------------- 116 ztortz2 = mzrat2 * 1.e6 * zfact * tr n(ji,jj,jk,jpmes)115 ztortz2 = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes) 117 116 ! 118 zcompadi = MAX( ( tr n(ji,jj,jk,jpdia) - xthresh2dia ), 0.e0 )119 zcompaz = MAX( ( tr n(ji,jj,jk,jpzoo) - xthresh2zoo ), 0.e0 )117 zcompadi = MAX( ( trb(ji,jj,jk,jpdia) - xthresh2dia ), 0.e0 ) 118 zcompaz = MAX( ( trb(ji,jj,jk,jpzoo) - xthresh2zoo ), 0.e0 ) 120 119 ! Size effect of nanophytoplankton on grazing : the smaller it is, the less prone 121 120 ! it is to predation by mesozooplankton 122 121 ! ------------------------------------------------------------------------------- 123 zcompaph = MAX( ( tr n(ji,jj,jk,jpphy) - xthresh2phy ), 0.e0 ) &122 zcompaph = MAX( ( trb(ji,jj,jk,jpphy) - xthresh2phy ), 0.e0 ) & 124 123 & * MIN(1., MAX( 0., ( quotan(ji,jj,jk) - 0.2) / 0.3 ) ) 125 zcompapoc = MAX( ( tr n(ji,jj,jk,jppoc) - xthresh2poc ), 0.e0 )124 zcompapoc = MAX( ( trb(ji,jj,jk,jppoc) - xthresh2poc ), 0.e0 ) 126 125 127 126 zfood = xprefc * zcompadi + xprefz * zcompaz + xprefp * zcompaph + xprefpoc * zcompapoc … … 129 128 zdenom = zfoodlim / ( xkgraz2 + zfoodlim ) 130 129 zdenom2 = zdenom / ( zfood + rtrn ) 131 zgraze2 = grazrat2 * zstep * tgfunc2(ji,jj,jk) * tr n(ji,jj,jk,jpmes)130 zgraze2 = grazrat2 * zstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpmes) 132 131 133 132 zgrazd = zgraze2 * xprefc * zcompadi * zdenom2 … … 136 135 zgrazpoc = zgraze2 * xprefpoc * zcompapoc * zdenom2 137 136 138 zgraznf = zgrazn * tr n(ji,jj,jk,jpnfe) / ( trn(ji,jj,jk,jpphy) + rtrn)139 zgrazf = zgrazd * tr n(ji,jj,jk,jpdfe) / ( trn(ji,jj,jk,jpdia) + rtrn)140 zgrazpof = zgrazpoc * tr n(ji,jj,jk,jpsfe) / ( trn(ji,jj,jk,jppoc) + rtrn)137 zgraznf = zgrazn * trb(ji,jj,jk,jpnfe) / ( trb(ji,jj,jk,jpphy) + rtrn) 138 zgrazf = zgrazd * trb(ji,jj,jk,jpdfe) / ( trb(ji,jj,jk,jpdia) + rtrn) 139 zgrazpof = zgrazpoc * trb(ji,jj,jk,jpsfe) / ( trb(ji,jj,jk,jppoc) + rtrn) 141 140 142 141 ! Mesozooplankton flux feeding on GOC … … 145 144 # if ! defined key_kriest 146 145 zgrazffeg = grazflux * zstep * wsbio4(ji,jj,jk) & 147 & * tgfunc2(ji,jj,jk) * tr n(ji,jj,jk,jpgoc) * trn(ji,jj,jk,jpmes)148 zgrazfffg = zgrazffeg * tr n(ji,jj,jk,jpbfe) / (trn(ji,jj,jk,jpgoc) + rtrn)146 & * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpgoc) * trb(ji,jj,jk,jpmes) 147 zgrazfffg = zgrazffeg * trb(ji,jj,jk,jpbfe) / (trb(ji,jj,jk,jpgoc) + rtrn) 149 148 # endif 150 149 zgrazffep = grazflux * zstep * wsbio3(ji,jj,jk) & 151 & * tgfunc2(ji,jj,jk) * tr n(ji,jj,jk,jppoc) * trn(ji,jj,jk,jpmes)152 zgrazfffp = zgrazffep * tr n(ji,jj,jk,jpsfe) / (trn(ji,jj,jk,jppoc) + rtrn)150 & * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpmes) 151 zgrazfffp = zgrazffep * trb(ji,jj,jk,jpsfe) / (trb(ji,jj,jk,jppoc) + rtrn) 153 152 ! 154 153 # if ! defined key_kriest … … 159 158 ! diatoms based aggregates are more prone to fractionation 160 159 ! since they are more porous (marine snow instead of fecal pellets) 161 zratio = tr n(ji,jj,jk,jpgsi) / ( trn(ji,jj,jk,jpgoc) + rtrn )160 zratio = trb(ji,jj,jk,jpgsi) / ( trb(ji,jj,jk,jpgoc) + rtrn ) 162 161 zratio2 = zratio * zratio 163 162 zfrac = zproport * grazflux * zstep * wsbio4(ji,jj,jk) & 164 & * tr n(ji,jj,jk,jpgoc) * trn(ji,jj,jk,jpmes) &165 & * ( 0. 1 + 3.9* zratio2 / ( 1.**2 + zratio2 ) )166 zfracfe = zfrac * tr n(ji,jj,jk,jpbfe) / (trn(ji,jj,jk,jpgoc) + rtrn)163 & * trb(ji,jj,jk,jpgoc) * trb(ji,jj,jk,jpmes) & 164 & * ( 0.2 + 3.8 * zratio2 / ( 1.**2 + zratio2 ) ) 165 zfracfe = zfrac * trb(ji,jj,jk,jpbfe) / (trb(ji,jj,jk,jpgoc) + rtrn) 167 166 168 167 zgrazffep = zproport * zgrazffep … … 186 185 187 186 ! Total grazing ( grazing by microzoo is already computed in p4zmicro ) 188 IF( l n_diatrc .AND. lk_iomput ) zgrazing(ji,jj,jk) = zgraztot187 IF( lk_iomput ) zgrazing(ji,jj,jk) = zgraztot 189 188 190 189 ! Mesozooplankton efficiency … … 216 215 tra(ji,jj,jk,jpzoo) = tra(ji,jj,jk,jpzoo) - zgrazz 217 216 tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) - zgrazn 218 tra(ji,jj,jk,jpnch) = tra(ji,jj,jk,jpnch) - zgrazn * tr n(ji,jj,jk,jpnch) / ( trn(ji,jj,jk,jpphy) + rtrn )219 tra(ji,jj,jk,jpdch) = tra(ji,jj,jk,jpdch) - zgrazd * tr n(ji,jj,jk,jpdch) / ( trn(ji,jj,jk,jpdia) + rtrn )220 tra(ji,jj,jk,jpdsi) = tra(ji,jj,jk,jpdsi) - zgrazd * tr n(ji,jj,jk,jpdsi) / ( trn(ji,jj,jk,jpdia) + rtrn )221 tra(ji,jj,jk,jpgsi) = tra(ji,jj,jk,jpgsi) + zgrazd * tr n(ji,jj,jk,jpdsi) / ( trn(ji,jj,jk,jpdia) + rtrn )217 tra(ji,jj,jk,jpnch) = tra(ji,jj,jk,jpnch) - zgrazn * trb(ji,jj,jk,jpnch) / ( trb(ji,jj,jk,jpphy) + rtrn ) 218 tra(ji,jj,jk,jpdch) = tra(ji,jj,jk,jpdch) - zgrazd * trb(ji,jj,jk,jpdch) / ( trb(ji,jj,jk,jpdia) + rtrn ) 219 tra(ji,jj,jk,jpdsi) = tra(ji,jj,jk,jpdsi) - zgrazd * trb(ji,jj,jk,jpdsi) / ( trb(ji,jj,jk,jpdia) + rtrn ) 220 tra(ji,jj,jk,jpgsi) = tra(ji,jj,jk,jpgsi) + zgrazd * trb(ji,jj,jk,jpdsi) / ( trb(ji,jj,jk,jpdia) + rtrn ) 222 221 tra(ji,jj,jk,jpnfe) = tra(ji,jj,jk,jpnfe) - zgraznf 223 222 tra(ji,jj,jk,jpdfe) = tra(ji,jj,jk,jpdfe) - zgrazf … … 232 231 tra(ji,jj,jk,jpcal) = tra(ji,jj,jk,jpcal) + zprcaca 233 232 #if defined key_kriest 234 znumpoc = tr n(ji,jj,jk,jpnum) / ( trn(ji,jj,jk,jppoc) + rtrn )233 znumpoc = trb(ji,jj,jk,jpnum) / ( trb(ji,jj,jk,jppoc) + rtrn ) 235 234 tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zmortzgoc - zgrazpoc - zgrazffep + zgrapoc2 236 235 tra(ji,jj,jk,jpnum) = tra(ji,jj,jk,jpnum) - zgrazpoc * znumpoc + zgrapoc2 * xkr_dmeso & … … 249 248 END DO 250 249 ! 251 IF( ln_diatrc .AND. lk_iomput .AND. jnt == nrdttrc ) THEN 252 zrfact2 = 1.e3 * rfact2r 253 CALL iom_put( "GRAZ2", zgrazing(:,:,:) * zrfact2 * tmask(:,:,:) ) ! Total grazing of phyto by zooplankton 254 CALL iom_put( "PCAL" , prodcal(:,:,:) * zrfact2 * tmask(:,:,:) ) ! Calcite production 250 IF( lk_iomput .AND. knt == nrdttrc ) THEN 251 CALL wrk_alloc( jpi, jpj, jpk, zw3d ) 252 IF( iom_use( "GRAZ2" ) ) THEN 253 zw3d(:,:,:) = zgrazing(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:) ! Total grazing of phyto by zooplankton 254 CALL iom_put( "GRAZ2", zw3d ) 255 ENDIF 256 IF( iom_use( "PCAL" ) ) THEN 257 zw3d(:,:,:) = prodcal(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:) ! Calcite production 258 CALL iom_put( "PCAL", zw3d ) 259 ENDIF 260 CALL wrk_dealloc( jpi, jpj, jpk, zw3d ) 255 261 ENDIF 256 262 ! … … 261 267 ENDIF 262 268 ! 263 IF( l n_diatrc .AND. lk_iomput ) CALL wrk_dealloc( jpi, jpj, jpk, zgrazing )269 IF( lk_iomput ) CALL wrk_dealloc( jpi, jpj, jpk, zgrazing ) 264 270 ! 265 271 IF( nn_timing == 1 ) CALL timing_stop('p4z_meso') -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zmicro.F90
r4793 r5581 59 59 CONTAINS 60 60 61 SUBROUTINE p4z_micro( kt, jnt )61 SUBROUTINE p4z_micro( kt, knt ) 62 62 !!--------------------------------------------------------------------- 63 63 !! *** ROUTINE p4z_micro *** … … 68 68 !!--------------------------------------------------------------------- 69 69 INTEGER, INTENT(in) :: kt ! ocean time step 70 INTEGER, INTENT(in) :: jnt70 INTEGER, INTENT(in) :: knt 71 71 ! 72 72 INTEGER :: ji, jj, jk … … 79 79 REAL(wp) :: zgrazp, zgrazm, zgrazsd 80 80 REAL(wp) :: zgrazmf, zgrazsf, zgrazpf 81 REAL(wp) :: zrfact2 82 REAL(wp), POINTER, DIMENSION(:,:,:) :: zgrazing 81 REAL(wp), POINTER, DIMENSION(:,:,:) :: zgrazing, zw3d 83 82 CHARACTER (len=25) :: charout 84 83 !!--------------------------------------------------------------------- … … 86 85 IF( nn_timing == 1 ) CALL timing_start('p4z_micro') 87 86 ! 88 IF( l n_diatrc .AND. lk_iomput ) CALL wrk_alloc( jpi, jpj, jpk, zgrazing )87 IF( lk_iomput ) CALL wrk_alloc( jpi, jpj, jpk, zgrazing ) 89 88 ! 90 89 DO jk = 1, jpkm1 91 90 DO jj = 1, jpj 92 91 DO ji = 1, jpi 93 zcompaz = MAX( ( tr n(ji,jj,jk,jpzoo) - 1.e-9 ), 0.e0 )92 zcompaz = MAX( ( trb(ji,jj,jk,jpzoo) - 1.e-9 ), 0.e0 ) 94 93 zstep = xstep 95 94 # if defined key_degrad … … 100 99 ! Respiration rates of both zooplankton 101 100 ! ------------------------------------- 102 zrespz = resrat * zfact * tr n(ji,jj,jk,jpzoo) / ( xkmort + trn(ji,jj,jk,jpzoo) ) &101 zrespz = resrat * zfact * trb(ji,jj,jk,jpzoo) / ( xkmort + trb(ji,jj,jk,jpzoo) ) & 103 102 & + resrat * zfact * 3. * nitrfac(ji,jj,jk) 104 103 … … 106 105 ! no real reason except that it seems to be more stable and may mimic predation. 107 106 ! --------------------------------------------------------------- 108 ztortz = mzrat * 1.e6 * zfact * tr n(ji,jj,jk,jpzoo)109 110 zcompadi = MIN( MAX( ( tr n(ji,jj,jk,jpdia) - xthreshdia ), 0.e0 ), xsizedia )111 zcompaph = MAX( ( tr n(ji,jj,jk,jpphy) - xthreshphy ), 0.e0 )112 zcompapoc = MAX( ( tr n(ji,jj,jk,jppoc) - xthreshpoc ), 0.e0 )107 ztortz = mzrat * 1.e6 * zfact * trb(ji,jj,jk,jpzoo) 108 109 zcompadi = MIN( MAX( ( trb(ji,jj,jk,jpdia) - xthreshdia ), 0.e0 ), xsizedia ) 110 zcompaph = MAX( ( trb(ji,jj,jk,jpphy) - xthreshphy ), 0.e0 ) 111 zcompapoc = MAX( ( trb(ji,jj,jk,jppoc) - xthreshpoc ), 0.e0 ) 113 112 114 113 ! Microzooplankton grazing … … 118 117 zdenom = zfoodlim / ( xkgraz + zfoodlim ) 119 118 zdenom2 = zdenom / ( zfood + rtrn ) 120 zgraze = grazrat * zstep * tgfunc2(ji,jj,jk) * tr n(ji,jj,jk,jpzoo)119 zgraze = grazrat * zstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpzoo) 121 120 122 121 zgrazp = zgraze * xpref2p * zcompaph * zdenom2 … … 124 123 zgrazsd = zgraze * xpref2d * zcompadi * zdenom2 125 124 126 zgrazpf = zgrazp * tr n(ji,jj,jk,jpnfe) / (trn(ji,jj,jk,jpphy) + rtrn)127 zgrazmf = zgrazm * tr n(ji,jj,jk,jpsfe) / (trn(ji,jj,jk,jppoc) + rtrn)128 zgrazsf = zgrazsd * tr n(ji,jj,jk,jpdfe) / (trn(ji,jj,jk,jpdia) + rtrn)125 zgrazpf = zgrazp * trb(ji,jj,jk,jpnfe) / (trb(ji,jj,jk,jpphy) + rtrn) 126 zgrazmf = zgrazm * trb(ji,jj,jk,jpsfe) / (trb(ji,jj,jk,jppoc) + rtrn) 127 zgrazsf = zgrazsd * trb(ji,jj,jk,jpdfe) / (trb(ji,jj,jk,jpdia) + rtrn) 129 128 ! 130 129 zgraztot = zgrazp + zgrazm + zgrazsd … … 137 136 ! Various remineralization and excretion terms 138 137 ! -------------------------------------------- 139 zgrasrat = zgraztotf/ ( zgraztot + rtrn )140 zgrasratn = zgraztotn/ ( zgraztot + rtrn )138 zgrasrat = ( zgraztotf + rtrn ) / ( zgraztot + rtrn ) 139 zgrasratn = ( zgraztotn + rtrn ) / ( zgraztot + rtrn ) 141 140 zepshert = MIN( 1., zgrasratn, zgrasrat / ferat3) 142 141 zepsherv = zepshert * MIN( epsher, (1. - unass) * zgrasrat / ferat3, (1. - unass) * zgrasratn ) … … 166 165 tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) - zgrazp 167 166 tra(ji,jj,jk,jpdia) = tra(ji,jj,jk,jpdia) - zgrazsd 168 tra(ji,jj,jk,jpnch) = tra(ji,jj,jk,jpnch) - zgrazp * tr n(ji,jj,jk,jpnch)/(trn(ji,jj,jk,jpphy)+rtrn)169 tra(ji,jj,jk,jpdch) = tra(ji,jj,jk,jpdch) - zgrazsd * tr n(ji,jj,jk,jpdch)/(trn(ji,jj,jk,jpdia)+rtrn)170 tra(ji,jj,jk,jpdsi) = tra(ji,jj,jk,jpdsi) - zgrazsd * tr n(ji,jj,jk,jpdsi)/(trn(ji,jj,jk,jpdia)+rtrn)171 tra(ji,jj,jk,jpgsi) = tra(ji,jj,jk,jpgsi) + zgrazsd * tr n(ji,jj,jk,jpdsi)/(trn(ji,jj,jk,jpdia)+rtrn)167 tra(ji,jj,jk,jpnch) = tra(ji,jj,jk,jpnch) - zgrazp * trb(ji,jj,jk,jpnch)/(trb(ji,jj,jk,jpphy)+rtrn) 168 tra(ji,jj,jk,jpdch) = tra(ji,jj,jk,jpdch) - zgrazsd * trb(ji,jj,jk,jpdch)/(trb(ji,jj,jk,jpdia)+rtrn) 169 tra(ji,jj,jk,jpdsi) = tra(ji,jj,jk,jpdsi) - zgrazsd * trb(ji,jj,jk,jpdsi)/(trb(ji,jj,jk,jpdia)+rtrn) 170 tra(ji,jj,jk,jpgsi) = tra(ji,jj,jk,jpgsi) + zgrazsd * trb(ji,jj,jk,jpdsi)/(trb(ji,jj,jk,jpdia)+rtrn) 172 171 tra(ji,jj,jk,jpnfe) = tra(ji,jj,jk,jpnfe) - zgrazpf 173 172 tra(ji,jj,jk,jpdfe) = tra(ji,jj,jk,jpdfe) - zgrazsf … … 185 184 #if defined key_kriest 186 185 tra(ji,jj,jk,jpnum) = tra(ji,jj,jk,jpnum) + zmortz * xkr_dmicro & 187 - zgrazm * tr n(ji,jj,jk,jpnum) / ( trn(ji,jj,jk,jppoc) + rtrn )186 - zgrazm * trb(ji,jj,jk,jpnum) / ( trb(ji,jj,jk,jppoc) + rtrn ) 188 187 #endif 189 188 END DO … … 191 190 END DO 192 191 ! 193 IF( ln_diatrc .AND. lk_iomput .AND. jnt == nrdttrc ) THEN 194 zrfact2 = 1.e3 * rfact2r 195 CALL iom_put( "GRAZ1" , zgrazing(:,:,:) * zrfact2 * tmask(:,:,:) ) ! Total grazing of phyto by zooplankton 192 IF( lk_iomput .AND. knt == nrdttrc ) THEN 193 CALL wrk_alloc( jpi, jpj, jpk, zw3d ) 194 IF( iom_use( "GRAZ1" ) ) THEN 195 zw3d(:,:,:) = zgrazing(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:) ! Total grazing of phyto by zooplankton 196 CALL iom_put( "GRAZ1", zw3d ) 197 ENDIF 198 CALL wrk_dealloc( jpi, jpj, jpk, zw3d ) 196 199 ENDIF 197 200 ! … … 202 205 ENDIF 203 206 ! 204 IF( l n_diatrc .AND. lk_iomput ) CALL wrk_dealloc( jpi, jpj, jpk, zgrazing )207 IF( lk_iomput ) CALL wrk_dealloc( jpi, jpj, jpk, zgrazing ) 205 208 ! 206 209 IF( nn_timing == 1 ) CALL timing_stop('p4z_micro') -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zmort.F90
- Property svn:keywords set to Id
r4793 r5581 39 39 !!---------------------------------------------------------------------- 40 40 !! NEMO/TOP 3.3 , NEMO Consortium (2010) 41 !! $Id : p4zmort.F90 3160 2011-11-20 14:27:18Z cetlod$41 !! $Id$ 42 42 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 43 43 !!---------------------------------------------------------------------- … … 85 85 DO jj = 1, jpj 86 86 DO ji = 1, jpi 87 zcompaph = MAX( ( tr n(ji,jj,jk,jpphy) - 1e-8 ), 0.e0 )87 zcompaph = MAX( ( trb(ji,jj,jk,jpphy) - 1e-8 ), 0.e0 ) 88 88 zstep = xstep 89 89 # if defined key_degrad … … 94 94 ! due to turbulence is negligible. Mortality is also set 95 95 ! to 0 96 zsizerat = MIN(1., MAX( 0., (quotan(ji,jj,jk) - 0.2) / 0.3) ) * tr n(ji,jj,jk,jpphy)96 zsizerat = MIN(1., MAX( 0., (quotan(ji,jj,jk) - 0.2) / 0.3) ) * trb(ji,jj,jk,jpphy) 97 97 ! Squared mortality of Phyto similar to a sedimentation term during 98 98 ! blooms (Doney et al. 1996) … … 102 102 ! increased when nutrients are limiting phytoplankton growth 103 103 ! as observed for instance in case of iron limitation. 104 ztortp = mprat * xstep * zcompaph / ( xkmort + tr n(ji,jj,jk,jpphy) ) * zsizerat104 ztortp = mprat * xstep * zcompaph / ( xkmort + trb(ji,jj,jk,jpphy) ) * zsizerat 105 105 106 106 zmortp = zrespp + ztortp … … 108 108 ! Update the arrays TRA which contains the biological sources and sinks 109 109 110 zfactfe = tr n(ji,jj,jk,jpnfe)/(trn(ji,jj,jk,jpphy)+rtrn)111 zfactch = tr n(ji,jj,jk,jpnch)/(trn(ji,jj,jk,jpphy)+rtrn)110 zfactfe = trb(ji,jj,jk,jpnfe)/(trb(ji,jj,jk,jpphy)+rtrn) 111 zfactch = trb(ji,jj,jk,jpnch)/(trb(ji,jj,jk,jpphy)+rtrn) 112 112 tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) - zmortp 113 113 tra(ji,jj,jk,jpnch) = tra(ji,jj,jk,jpnch) - zmortp * zfactch … … 172 172 DO ji = 1, jpi 173 173 174 zcompadi = MAX( ( tr n(ji,jj,jk,jpdia) - 1e-9), 0. )174 zcompadi = MAX( ( trb(ji,jj,jk,jpdia) - 1e-9), 0. ) 175 175 176 176 ! Aggregation term for diatoms is increased in case of nutrient … … 186 186 zlim2 = xlimdia(ji,jj,jk) * xlimdia(ji,jj,jk) 187 187 zlim1 = 0.25 * ( 1. - zlim2 ) / ( 0.25 + zlim2 ) 188 zrespp2 = 1.e6 * zstep * ( wchld + wchldm * zlim1 ) * xdiss(ji,jj,jk) * zcompadi * tr n(ji,jj,jk,jpdia)188 zrespp2 = 1.e6 * zstep * ( wchld + wchldm * zlim1 ) * xdiss(ji,jj,jk) * zcompadi * trb(ji,jj,jk,jpdia) 189 189 190 190 ! Phytoplankton mortality. 191 191 ! ------------------------ 192 ztortp2 = mprat2 * zstep * tr n(ji,jj,jk,jpdia) / ( xkmort + trn(ji,jj,jk,jpdia) ) * zcompadi192 ztortp2 = mprat2 * zstep * trb(ji,jj,jk,jpdia) / ( xkmort + trb(ji,jj,jk,jpdia) ) * zcompadi 193 193 194 194 zmortp2 = zrespp2 + ztortp2 … … 196 196 ! Update the arrays tra which contains the biological sources and sinks 197 197 ! --------------------------------------------------------------------- 198 zfactch = tr n(ji,jj,jk,jpdch) / ( trn(ji,jj,jk,jpdia) + rtrn )199 zfactfe = tr n(ji,jj,jk,jpdfe) / ( trn(ji,jj,jk,jpdia) + rtrn )200 zfactsi = tr n(ji,jj,jk,jpdsi) / ( trn(ji,jj,jk,jpdia) + rtrn )198 zfactch = trb(ji,jj,jk,jpdch) / ( trb(ji,jj,jk,jpdia) + rtrn ) 199 zfactfe = trb(ji,jj,jk,jpdfe) / ( trb(ji,jj,jk,jpdia) + rtrn ) 200 zfactsi = trb(ji,jj,jk,jpdsi) / ( trb(ji,jj,jk,jpdia) + rtrn ) 201 201 tra(ji,jj,jk,jpdia) = tra(ji,jj,jk,jpdia) - zmortp2 202 202 tra(ji,jj,jk,jpdch) = tra(ji,jj,jk,jpdch) - zmortp2 * zfactch -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zopt.F90
r4793 r5581 35 35 REAL(wp) :: parlux !: Fraction of shortwave as PAR 36 36 REAL(wp) :: xparsw !: parlux/3 37 REAL(wp) :: xsi0r !: 1. /rn_si0 37 38 38 39 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_par ! structure of input par … … 42 43 43 44 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: enano, ediat !: PAR for phyto, nano and diat 45 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: etot_ndcy !: PAR over 24h in case of diurnal cycle 44 46 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: emoy !: averaged PAR in the mixed layer 47 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ekb, ekg, ekr !: wavelength (Red-Green-Blue) 45 48 46 49 INTEGER :: nksrp ! levels below which the light cannot penetrate ( depth larger than 391 m) … … 57 60 CONTAINS 58 61 59 SUBROUTINE p4z_opt( kt, jnt )62 SUBROUTINE p4z_opt( kt, knt ) 60 63 !!--------------------------------------------------------------------- 61 64 !! *** ROUTINE p4z_opt *** … … 67 70 !!--------------------------------------------------------------------- 68 71 ! 69 INTEGER, INTENT(in) :: kt, jnt ! ocean time step72 INTEGER, INTENT(in) :: kt, knt ! ocean time step 70 73 ! 71 74 INTEGER :: ji, jj, jk 72 75 INTEGER :: irgb 73 REAL(wp) :: zchl , zxsi0r76 REAL(wp) :: zchl 74 77 REAL(wp) :: zc0 , zc1 , zc2, zc3, z1_dep 75 REAL(wp), POINTER, DIMENSION(:,: ) :: zdepmoy, zetmp , zetmp1, zetmp276 REAL(wp), POINTER, DIMENSION(:,:,:) :: z ekg, zekr, zekb, ze0, ze1, ze2, ze378 REAL(wp), POINTER, DIMENSION(:,: ) :: zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4, zqsr100 79 REAL(wp), POINTER, DIMENSION(:,:,:) :: zpar, ze0, ze1, ze2, ze3 77 80 !!--------------------------------------------------------------------- 78 81 ! … … 80 83 ! 81 84 ! Allocate temporary workspace 82 CALL wrk_alloc( jpi, jpj, z depmoy, zetmp, zetmp1, zetmp2 )83 CALL wrk_alloc( jpi, jpj, jpk, z ekg, zekr, zekb, ze0, ze1, ze2, ze3 )84 85 IF( jnt == 1 .AND. ln_varpar ) CALL p4z_optsbc( kt )85 CALL wrk_alloc( jpi, jpj, zqsr100, zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 ) 86 CALL wrk_alloc( jpi, jpj, jpk, zpar, ze0, ze1, ze2, ze3 ) 87 88 IF( knt == 1 .AND. ln_varpar ) CALL p4z_opt_sbc( kt ) 86 89 87 90 ! Initialisation of variables used to compute PAR 88 91 ! ----------------------------------------------- 89 ze1(:,:,jpk) = 0._wp 90 ze2(:,:,jpk) = 0._wp 91 ze3(:,:,jpk) = 0._wp 92 92 ze1(:,:,:) = 0._wp 93 ze2(:,:,:) = 0._wp 94 ze3(:,:,:) = 0._wp 93 95 ! !* attenuation coef. function of Chlorophyll and wavelength (Red-Green-Blue) 94 96 DO jk = 1, jpkm1 ! -------------------------------------------------------- … … 97 99 !CDIR NOVERRCHK 98 100 DO ji = 1, jpi 99 zchl = ( tr n(ji,jj,jk,jpnch) + trn(ji,jj,jk,jpdch) + rtrn ) * 1.e6101 zchl = ( trb(ji,jj,jk,jpnch) + trb(ji,jj,jk,jpdch) + rtrn ) * 1.e6 100 102 zchl = MIN( 10. , MAX( 0.05, zchl ) ) 101 103 irgb = NINT( 41 + 20.* LOG10( zchl ) + rtrn ) 102 104 ! 103 zekb(ji,jj,jk) = xkrgb(1,irgb) * fse3t(ji,jj,jk)104 zekg(ji,jj,jk) = xkrgb(2,irgb) * fse3t(ji,jj,jk)105 zekr(ji,jj,jk) = xkrgb(3,irgb) * fse3t(ji,jj,jk)105 ekb(ji,jj,jk) = xkrgb(1,irgb) * fse3t(ji,jj,jk) 106 ekg(ji,jj,jk) = xkrgb(2,irgb) * fse3t(ji,jj,jk) 107 ekr(ji,jj,jk) = xkrgb(3,irgb) * fse3t(ji,jj,jk) 106 108 END DO 107 109 END DO 108 110 END DO 109 110 111 111 ! !* Photosynthetically Available Radiation (PAR) 112 112 ! ! -------------------------------------- 113 114 IF( ln_varpar ) THEN 115 ze1(:,:,1) = par_varsw(:,:) * qsr(:,:) * EXP( -0.5 * zekb(:,:,1) ) 116 ze2(:,:,1) = par_varsw(:,:) * qsr(:,:) * EXP( -0.5 * zekg(:,:,1) ) 117 ze3(:,:,1) = par_varsw(:,:) * qsr(:,:) * EXP( -0.5 * zekr(:,:,1) ) 113 IF( l_trcdm2dc ) THEN ! diurnal cycle 114 ! 1% of qsr to compute euphotic layer 115 zqsr100(:,:) = 0.01 * qsr_mean(:,:) ! daily mean qsr 116 ! 117 CALL p4z_opt_par( kt, qsr_mean, ze1, ze2, ze3 ) 118 ! 119 DO jk = 1, nksrp 120 etot_ndcy(:,:,jk) = ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk) 121 enano (:,:,jk) = 2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 122 ediat (:,:,jk) = 1.6 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.7 * ze3(:,:,jk) 123 END DO 124 ! 125 CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3 ) 126 ! 127 DO jk = 1, nksrp 128 etot(:,:,jk) = ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk) 129 END DO 130 ! 118 131 ELSE 119 ze1(:,:,1) = xparsw * qsr(:,:) * EXP( -0.5 * zekb(:,:,1) ) 120 ze2(:,:,1) = xparsw * qsr(:,:) * EXP( -0.5 * zekg(:,:,1) ) 121 ze3(:,:,1) = xparsw * qsr(:,:) * EXP( -0.5 * zekr(:,:,1) ) 122 ENDIF 123 124 !CDIR NOVERRCHK 125 DO jj = 1, jpj 126 !CDIR NOVERRCHK 127 DO ji = 1, jpi 128 zc1 = ze1(ji,jj,1) 129 zc2 = ze2(ji,jj,1) 130 zc3 = ze3(ji,jj,1) 131 etot (ji,jj,1) = ( zc1 + zc2 + zc3 ) 132 enano(ji,jj,1) = ( 2.1 * zc1 + 0.42 * zc2 + 0.4 * zc3 ) 133 ediat(ji,jj,1) = ( 1.6 * zc1 + 0.69 * zc2 + 0.7 * zc3 ) 134 END DO 135 END DO 136 137 138 DO jk = 2, nksrp 139 !CDIR NOVERRCHK 140 DO jj = 1, jpj 141 !CDIR NOVERRCHK 142 DO ji = 1, jpi 143 zc1 = ze1(ji,jj,jk-1) * EXP( -0.5 * ( zekb(ji,jj,jk-1) + zekb(ji,jj,jk) ) ) 144 zc2 = ze2(ji,jj,jk-1) * EXP( -0.5 * ( zekg(ji,jj,jk-1) + zekg(ji,jj,jk) ) ) 145 zc3 = ze3(ji,jj,jk-1) * EXP( -0.5 * ( zekr(ji,jj,jk-1) + zekr(ji,jj,jk) ) ) 146 ze1 (ji,jj,jk) = zc1 147 ze2 (ji,jj,jk) = zc2 148 ze3 (ji,jj,jk) = zc3 149 etot (ji,jj,jk) = ( zc1 + zc2 + zc3 ) 150 enano(ji,jj,jk) = ( 2.1 * zc1 + 0.42 * zc2 + 0.4 * zc3 ) 151 ediat(ji,jj,jk) = ( 1.6 * zc1 + 0.69 * zc2 + 0.7 * zc3 ) 152 END DO 153 END DO 154 END DO 132 ! 1% of qsr to compute euphotic layer 133 zqsr100(:,:) = 0.01 * qsr(:,:) 134 ! 135 CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3 ) 136 ! 137 DO jk = 1, nksrp 138 etot (:,:,jk) = ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk) 139 enano(:,:,jk) = 2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 140 ediat(:,:,jk) = 1.6 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.7 * ze3(:,:,jk) 141 END DO 142 etot_ndcy(:,:,:) = etot(:,:,:) 143 ENDIF 144 155 145 156 146 IF( ln_qsr_bio ) THEN !* heat flux accros w-level (used in the dynamics) 157 147 ! ! ------------------------ 158 zxsi0r = 1.e0 / rn_si0 159 ! 160 ze0(:,:,1) = rn_abs * qsr(:,:) 161 ! ! surface value : separation in R-G-B + near surface 162 IF( ln_varpar ) THEN 163 ze0(:,:,1) = ( 1. - 3. * par_varsw(:,:) ) * qsr(:,:) 164 ze1(:,:,1) = par_varsw(:,:) * qsr(:,:) 165 ze2(:,:,1) = par_varsw(:,:) * qsr(:,:) 166 ze3(:,:,1) = par_varsw(:,:) * qsr(:,:) 167 ELSE 168 ze0(:,:,1) = ( 1. - 3. * xparsw ) * qsr(:,:) 169 ze1(:,:,1) = xparsw * qsr(:,:) 170 ze2(:,:,1) = xparsw * qsr(:,:) 171 ze3(:,:,1) = xparsw * qsr(:,:) 172 ENDIF 148 CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3, pe0=ze0 ) 149 ! 173 150 etot3(:,:,1) = qsr(:,:) * tmask(:,:,1) 174 !175 !176 151 DO jk = 2, nksrp + 1 177 !CDIR NOVERRCHK 178 DO jj = 1, jpj 179 !CDIR NOVERRCHK 180 DO ji = 1, jpi 181 zc0 = ze0(ji,jj,jk-1) * EXP( -fse3t(ji,jj,jk-1) * zxsi0r ) 182 zc1 = ze1(ji,jj,jk-1) * EXP( -zekb(ji,jj,jk-1 ) ) 183 zc2 = ze2(ji,jj,jk-1) * EXP( -zekg(ji,jj,jk-1 ) ) 184 zc3 = ze3(ji,jj,jk-1) * EXP( -zekr(ji,jj,jk-1 ) ) 185 ze0(ji,jj,jk) = zc0 186 ze1(ji,jj,jk) = zc1 187 ze2(ji,jj,jk) = zc2 188 ze3(ji,jj,jk) = zc3 189 etot3(ji,jj,jk) = ( zc0 + zc1 + zc2 + zc3 ) * tmask(ji,jj,jk) 190 END DO 191 ! 192 END DO 193 ! 194 END DO 195 ! 196 ENDIF 197 152 etot3(:,:,jk) = ( ze0(:,:,jk) + ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk) ) * tmask(:,:,jk) 153 END DO 154 ! ! ------------------------ 155 ENDIF 198 156 ! !* Euphotic depth and level 199 157 neln(:,:) = 1 ! ------------------------ … … 203 161 DO jj = 1, jpj 204 162 DO ji = 1, jpi 205 IF( etot (ji,jj,jk) * tmask(ji,jj,jk) >= 0.0043 * qsr(ji,jj) ) THEN163 IF( etot_ndcy(ji,jj,jk) * tmask(ji,jj,jk) >= 0.43 * zqsr100(ji,jj) ) THEN 206 164 neln(ji,jj) = jk+1 ! Euphotic level : 1rst T-level strictly below Euphotic layer 207 165 ! ! nb: ensure the compatibility with nmld_trc definition in trd_mld_trc_zint … … 211 169 END DO 212 170 END DO 213 171 ! 214 172 heup(:,:) = MIN( 300., heup(:,:) ) 215 216 173 ! !* mean light over the mixed layer 217 174 zdepmoy(:,:) = 0.e0 ! ------------------------------- 218 zetmp (:,:) = 0.e0219 175 zetmp1 (:,:) = 0.e0 220 176 zetmp2 (:,:) = 0.e0 177 zetmp3 (:,:) = 0.e0 178 zetmp4 (:,:) = 0.e0 221 179 222 180 DO jk = 1, nksrp … … 226 184 DO ji = 1, jpi 227 185 IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN 228 zetmp (ji,jj) = zetmp (ji,jj) + etot (ji,jj,jk) * fse3t(ji,jj,jk) 229 zetmp1 (ji,jj) = zetmp1 (ji,jj) + enano(ji,jj,jk) * fse3t(ji,jj,jk) 230 zetmp2 (ji,jj) = zetmp2 (ji,jj) + ediat(ji,jj,jk) * fse3t(ji,jj,jk) 186 zetmp1 (ji,jj) = zetmp1 (ji,jj) + etot (ji,jj,jk) * fse3t(ji,jj,jk) ! remineralisation 187 zetmp2 (ji,jj) = zetmp2 (ji,jj) + etot_ndcy(ji,jj,jk) * fse3t(ji,jj,jk) ! production 188 zetmp3 (ji,jj) = zetmp3 (ji,jj) + enano (ji,jj,jk) * fse3t(ji,jj,jk) ! production 189 zetmp4 (ji,jj) = zetmp4 (ji,jj) + ediat (ji,jj,jk) * fse3t(ji,jj,jk) ! production 231 190 zdepmoy(ji,jj) = zdepmoy(ji,jj) + fse3t(ji,jj,jk) 232 191 ENDIF … … 235 194 END DO 236 195 ! 237 emoy(:,:,:) = etot(:,:,:) 196 emoy(:,:,:) = etot(:,:,:) ! remineralisation 197 zpar(:,:,:) = etot_ndcy(:,:,:) ! diagnostic : PAR with no diurnal cycle 238 198 ! 239 199 DO jk = 1, nksrp … … 244 204 IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN 245 205 z1_dep = 1. / ( zdepmoy(ji,jj) + rtrn ) 246 emoy (ji,jj,jk) = zetmp (ji,jj) * z1_dep 247 enano(ji,jj,jk) = zetmp1(ji,jj) * z1_dep 248 ediat(ji,jj,jk) = zetmp2(ji,jj) * z1_dep 206 emoy (ji,jj,jk) = zetmp1(ji,jj) * z1_dep 207 zpar (ji,jj,jk) = zetmp2(ji,jj) * z1_dep 208 enano(ji,jj,jk) = zetmp3(ji,jj) * z1_dep 209 ediat(ji,jj,jk) = zetmp4(ji,jj) * z1_dep 249 210 ENDIF 250 211 END DO 251 212 END DO 252 213 END DO 253 254 IF( ln_diatrc ) THEN ! save output diagnostics 214 ! 215 IF( lk_iomput ) THEN 216 IF( knt == nrdttrc ) THEN 217 IF( iom_use( "Heup" ) ) CALL iom_put( "Heup" , heup(:,: ) * tmask(:,:,1) ) ! euphotic layer deptht 218 IF( iom_use( "PARDM" ) ) CALL iom_put( "PARDM", zpar(:,:,:) * tmask(:,:,:) ) ! Photosynthetically Available Radiation 219 IF( iom_use( "PAR" ) ) CALL iom_put( "PAR" , emoy(:,:,:) * tmask(:,:,:) ) ! Photosynthetically Available Radiation 220 ENDIF 221 ELSE 222 IF( ln_diatrc ) THEN ! save output diagnostics 223 trc2d(:,:, jp_pcs0_2d + 10) = heup(:,: ) * tmask(:,:,1) 224 trc3d(:,:,:,jp_pcs0_3d + 3) = etot(:,:,:) * tmask(:,:,:) 225 ENDIF 226 ENDIF 227 ! 228 CALL wrk_dealloc( jpi, jpj, zqsr100, zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 ) 229 CALL wrk_dealloc( jpi, jpj, jpk, zpar, ze0, ze1, ze2, ze3 ) 230 ! 231 IF( nn_timing == 1 ) CALL timing_stop('p4z_opt') 232 ! 233 END SUBROUTINE p4z_opt 234 235 SUBROUTINE p4z_opt_par( kt, pqsr, pe1, pe2, pe3, pe0 ) 236 !!---------------------------------------------------------------------- 237 !! *** routine p4z_opt_par *** 238 !! 239 !! ** purpose : compute PAR of each wavelength (Red-Green-Blue) 240 !! for a given shortwave radiation 241 !! 242 !!---------------------------------------------------------------------- 243 !! * arguments 244 INTEGER, INTENT(in) :: kt ! ocean time-step 245 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in) :: pqsr ! shortwave 246 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pe1 , pe2 , pe3 ! PAR ( R-G-B) 247 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout), OPTIONAL :: pe0 248 !! * local variables 249 INTEGER :: ji, jj, jk ! dummy loop indices 250 REAL(wp), DIMENSION(jpi,jpj) :: zqsr ! shortwave 251 !!---------------------------------------------------------------------- 252 253 ! Real shortwave 254 IF( ln_varpar ) THEN ; zqsr(:,:) = par_varsw(:,:) * pqsr(:,:) 255 ELSE ; zqsr(:,:) = xparsw * pqsr(:,:) 256 ENDIF 257 ! 258 IF( PRESENT( pe0 ) ) THEN ! W-level 259 ! 260 pe0(:,:,1) = pqsr(:,:) - 3. * zqsr(:,:) ! ( 1 - 3 * alpha ) * q 261 pe1(:,:,1) = zqsr(:,:) 262 pe2(:,:,1) = zqsr(:,:) 263 pe3(:,:,1) = zqsr(:,:) 264 ! 265 DO jk = 2, nksrp + 1 266 !CDIR NOVERRCHK 267 DO jj = 1, jpj 268 !CDIR NOVERRCHK 269 DO ji = 1, jpi 270 pe0(ji,jj,jk) = pe0(ji,jj,jk-1) * EXP( -fse3t(ji,jj,jk-1) * xsi0r ) 271 pe1(ji,jj,jk) = pe1(ji,jj,jk-1) * EXP( -ekb(ji,jj,jk-1 ) ) 272 pe2(ji,jj,jk) = pe2(ji,jj,jk-1) * EXP( -ekg(ji,jj,jk-1 ) ) 273 pe3(ji,jj,jk) = pe3(ji,jj,jk-1) * EXP( -ekr(ji,jj,jk-1 ) ) 274 END DO 275 ! 276 END DO 277 ! 278 END DO 255 279 ! 256 IF( lk_iomput ) THEN 257 IF( jnt == nrdttrc ) THEN 258 CALL iom_put( "Heup", heup(:,: ) * tmask(:,:,1) ) ! euphotic layer deptht 259 CALL iom_put( "PAR" , emoy(:,:,:) * tmask(:,:,:) ) ! Photosynthetically Available Radiation 260 ENDIF 261 ELSE 262 trc2d(:,:, jp_pcs0_2d + 10) = heup(:,: ) * tmask(:,:,1) 263 trc3d(:,:,:,jp_pcs0_3d + 3) = etot(:,:,:) * tmask(:,:,:) 264 ENDIF 280 ELSE ! T- level 265 281 ! 266 ENDIF 267 ! 268 CALL wrk_dealloc( jpi, jpj, zdepmoy, zetmp, zetmp1, zetmp2 ) 269 CALL wrk_dealloc( jpi, jpj, jpk, zekg, zekr, zekb, ze0, ze1, ze2, ze3 ) 270 ! 271 IF( nn_timing == 1 ) CALL timing_stop('p4z_opt') 272 ! 273 END SUBROUTINE p4z_opt 274 275 SUBROUTINE p4z_optsbc( kt ) 276 !!---------------------------------------------------------------------- 277 !! *** routine p4z_optsbc *** 282 pe1(:,:,1) = zqsr(:,:) * EXP( -0.5 * ekb(:,:,1) ) 283 pe2(:,:,1) = zqsr(:,:) * EXP( -0.5 * ekg(:,:,1) ) 284 pe3(:,:,1) = zqsr(:,:) * EXP( -0.5 * ekr(:,:,1) ) 285 ! 286 DO jk = 2, nksrp 287 !CDIR NOVERRCHK 288 DO jj = 1, jpj 289 !CDIR NOVERRCHK 290 DO ji = 1, jpi 291 pe1(ji,jj,jk) = pe1(ji,jj,jk-1) * EXP( -0.5 * ( ekb(ji,jj,jk-1) + ekb(ji,jj,jk) ) ) 292 pe2(ji,jj,jk) = pe2(ji,jj,jk-1) * EXP( -0.5 * ( ekg(ji,jj,jk-1) + ekg(ji,jj,jk) ) ) 293 pe3(ji,jj,jk) = pe3(ji,jj,jk-1) * EXP( -0.5 * ( ekr(ji,jj,jk-1) + ekr(ji,jj,jk) ) ) 294 END DO 295 END DO 296 END DO 297 ! 298 ENDIF 299 ! 300 END SUBROUTINE p4z_opt_par 301 302 303 SUBROUTINE p4z_opt_sbc( kt ) 304 !!---------------------------------------------------------------------- 305 !! *** routine p4z_opt_sbc *** 278 306 !! 279 307 !! ** purpose : read and interpolate the variable PAR fraction … … 286 314 !!---------------------------------------------------------------------- 287 315 !! * arguments 288 INTEGER , INTENT( in ) :: kt! ocean time step316 INTEGER , INTENT(in) :: kt ! ocean time step 289 317 290 318 !! * local declarations … … 299 327 IF( kt == nit000 .OR. ( kt /= nit000 .AND. ntimes_par > 1 ) ) THEN 300 328 CALL fld_read( kt, 1, sf_par ) 301 par_varsw(:,:) = ( sf_par(1)%fnow(:,:,1) ) /3.0329 par_varsw(:,:) = ( sf_par(1)%fnow(:,:,1) ) / 3.0 302 330 ENDIF 303 331 ENDIF … … 305 333 IF( nn_timing == 1 ) CALL timing_stop('p4z_optsbc') 306 334 ! 307 END SUBROUTINE p4z_opt sbc335 END SUBROUTINE p4z_opt_sbc 308 336 309 337 SUBROUTINE p4z_opt_init … … 349 377 ! 350 378 xparsw = parlux / 3.0 379 xsi0r = 1.e0 / rn_si0 351 380 ! 352 381 ! Variable PAR at the surface of the ocean … … 374 403 IF(lwp) WRITE(numout,*) ' level of light extinction = ', nksrp, ' ref depth = ', gdepw_1d(nksrp+1), ' m' 375 404 ! 376 etot (:,:,:) = 0._wp 377 enano(:,:,:) = 0._wp 378 ediat(:,:,:) = 0._wp 379 IF( ln_qsr_bio ) etot3(:,:,:) = 0._wp 405 ekr (:,:,:) = 0._wp 406 ekb (:,:,:) = 0._wp 407 ekg (:,:,:) = 0._wp 408 etot (:,:,:) = 0._wp 409 etot_ndcy(:,:,:) = 0._wp 410 enano (:,:,:) = 0._wp 411 ediat (:,:,:) = 0._wp 412 IF( ln_qsr_bio ) etot3 (:,:,:) = 0._wp 380 413 ! 381 414 IF( nn_timing == 1 ) CALL timing_stop('p4z_opt_init') … … 388 421 !! *** ROUTINE p4z_opt_alloc *** 389 422 !!---------------------------------------------------------------------- 390 ALLOCATE( enano(jpi,jpj,jpk), ediat(jpi,jpj,jpk), emoy (jpi,jpj,jpk), STAT=p4z_opt_alloc ) 423 ALLOCATE( ekb(jpi,jpj,jpk) , ekr(jpi,jpj,jpk), ekg(jpi,jpj,jpk), & 424 & enano(jpi,jpj,jpk) , ediat(jpi,jpj,jpk), & 425 & etot_ndcy(jpi,jpj,jpk), emoy (jpi,jpj,jpk), STAT=p4z_opt_alloc ) 391 426 ! 392 427 IF( p4z_opt_alloc /= 0 ) CALL ctl_warn('p4z_opt_alloc : failed to allocate arrays.') -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zprod.F90
r4793 r5581 64 64 CONTAINS 65 65 66 SUBROUTINE p4z_prod( kt , jnt )66 SUBROUTINE p4z_prod( kt , knt ) 67 67 !!--------------------------------------------------------------------- 68 68 !! *** ROUTINE p4z_prod *** … … 74 74 !!--------------------------------------------------------------------- 75 75 ! 76 INTEGER, INTENT(in) :: kt, jnt76 INTEGER, INTENT(in) :: kt, knt 77 77 ! 78 78 INTEGER :: ji, jj, jk … … 83 83 REAL(wp) :: zpislopen , zpislope2n 84 84 REAL(wp) :: zrum, zcodel, zargu, zval 85 REAL(wp) :: z rfact285 REAL(wp) :: zfact 86 86 CHARACTER (len=25) :: charout 87 REAL(wp), POINTER, DIMENSION(:,: ) :: zmixnano, zmixdiat, zstrn 88 REAL(wp), POINTER, DIMENSION(:,:,:) :: zpislopead, zpislopead2, zprdia, zprbio, zprdch, zprnch, zysopt 87 REAL(wp), POINTER, DIMENSION(:,: ) :: zmixnano, zmixdiat, zstrn, zw2d 88 REAL(wp), POINTER, DIMENSION(:,:,:) :: zpislopead, zpislopead2, zprdia, zprbio, zprdch, zprnch, zysopt, zw3d 89 89 REAL(wp), POINTER, DIMENSION(:,:,:) :: zprorca, zprorcad, zprofed, zprofen, zprochln, zprochld, zpronew, zpronewd 90 90 !!--------------------------------------------------------------------- … … 129 129 END DO 130 130 131 IF( ln_newprod ) THEN 132 ! Impact of the day duration on phytoplankton growth 133 DO jk = 1, jpkm1 134 DO jj = 1 ,jpj 135 DO ji = 1, jpi 136 IF( etot(ji,jj,jk) > 1.E-3 ) THEN 137 zval = MAX( 1., zstrn(ji,jj) ) 138 zval = 1.5 * zval / ( 12. + zval ) 139 zprbio(ji,jj,jk) = prmax(ji,jj,jk) * zval 140 zprdia(ji,jj,jk) = zprbio(ji,jj,jk) 141 ENDIF 142 END DO 143 END DO 144 END DO 145 ENDIF 131 ! Impact of the day duration on phytoplankton growth 132 DO jk = 1, jpkm1 133 DO jj = 1 ,jpj 134 DO ji = 1, jpi 135 IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 136 zval = MAX( 1., zstrn(ji,jj) ) 137 zval = 1.5 * zval / ( 12. + zval ) 138 zprbio(ji,jj,jk) = prmax(ji,jj,jk) * zval 139 zprdia(ji,jj,jk) = zprbio(ji,jj,jk) 140 ENDIF 141 END DO 142 END DO 143 END DO 146 144 147 145 ! Maximum light intensity … … 157 155 DO ji = 1, jpi 158 156 ! Computation of the P-I slope for nanos and diatoms 159 IF( etot (ji,jj,jk) > 1.E-3 ) THEN157 IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 160 158 ztn = MAX( 0., tsn(ji,jj,jk,jp_tem) - 15. ) 161 159 zadap = xadap * ztn / ( 2.+ ztn ) 162 zconctemp = MAX( 0.e0 , tr n(ji,jj,jk,jpdia) - xsizedia )163 zconctemp2 = tr n(ji,jj,jk,jpdia) - zconctemp160 zconctemp = MAX( 0.e0 , trb(ji,jj,jk,jpdia) - xsizedia ) 161 zconctemp2 = trb(ji,jj,jk,jpdia) - zconctemp 164 162 znanotot = enano(ji,jj,jk) * zstrn(ji,jj) 165 163 zdiattot = ediat(ji,jj,jk) * zstrn(ji,jj) 166 164 ! 167 165 zpislopead (ji,jj,jk) = pislope * ( 1.+ zadap * EXP( -znanotot ) ) & 168 & * tr n(ji,jj,jk,jpnch) /( trn(ji,jj,jk,jpphy) * 12. + rtrn)166 & * trb(ji,jj,jk,jpnch) /( trb(ji,jj,jk,jpphy) * 12. + rtrn) 169 167 ! 170 zpislopead2(ji,jj,jk) = (pislope * zconctemp2 + pislope2 * zconctemp) / ( tr n(ji,jj,jk,jpdia) + rtrn ) &171 & * tr n(ji,jj,jk,jpdch) /( trn(ji,jj,jk,jpdia) * 12. + rtrn)168 zpislopead2(ji,jj,jk) = (pislope * zconctemp2 + pislope2 * zconctemp) / ( trb(ji,jj,jk,jpdia) + rtrn ) & 169 & * trb(ji,jj,jk,jpdch) /( trb(ji,jj,jk,jpdia) * 12. + rtrn) 172 170 173 171 ! Computation of production function for Carbon … … 196 194 197 195 ! Computation of the P-I slope for nanos and diatoms 198 IF( etot (ji,jj,jk) > 1.E-3 ) THEN196 IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 199 197 ztn = MAX( 0., tsn(ji,jj,jk,jp_tem) - 15. ) 200 198 zadap = ztn / ( 2.+ ztn ) 201 zconctemp = MAX( 0.e0 , trn(ji,jj,jk,jpdia) - xsizedia ) 202 zconctemp2 = trn(ji,jj,jk,jpdia) - zconctemp 199 zconctemp = MAX( 0.e0 , trb(ji,jj,jk,jpdia) - xsizedia ) 200 zconctemp2 = trb(ji,jj,jk,jpdia) - zconctemp 201 znanotot = enano(ji,jj,jk) * zstrn(ji,jj) 202 zdiattot = ediat(ji,jj,jk) * zstrn(ji,jj) 203 203 ! 204 zpislopead (ji,jj,jk) = pislope * ( 1.+ zadap * EXP( - 0.21 * enano(ji,jj,jk)) )205 zpislopead2(ji,jj,jk) = (pislope * zconctemp2 + pislope2 * zconctemp) / ( tr n(ji,jj,jk,jpdia) + rtrn )206 207 zpislopen = zpislopead(ji,jj,jk) * tr n(ji,jj,jk,jpnch) &208 & / ( tr n(ji,jj,jk,jpphy) * 12. + rtrn ) &204 zpislopead (ji,jj,jk) = pislope * ( 1.+ zadap * EXP( -znanotot ) ) 205 zpislopead2(ji,jj,jk) = (pislope * zconctemp2 + pislope2 * zconctemp) / ( trb(ji,jj,jk,jpdia) + rtrn ) 206 207 zpislopen = zpislopead(ji,jj,jk) * trb(ji,jj,jk,jpnch) & 208 & / ( trb(ji,jj,jk,jpphy) * 12. + rtrn ) & 209 209 & / ( prmax(ji,jj,jk) * rday * xlimphy(ji,jj,jk) + rtrn ) 210 210 211 zpislope2n = zpislopead2(ji,jj,jk) * tr n(ji,jj,jk,jpdch) &212 & / ( tr n(ji,jj,jk,jpdia) * 12. + rtrn ) &211 zpislope2n = zpislopead2(ji,jj,jk) * trb(ji,jj,jk,jpdch) & 212 & / ( trb(ji,jj,jk,jpdia) * 12. + rtrn ) & 213 213 & / ( prmax(ji,jj,jk) * rday * xlimdia(ji,jj,jk) + rtrn ) 214 214 215 215 ! Computation of production function for Carbon 216 216 ! --------------------------------------------- 217 zprbio(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk)) )218 zprdia(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislope2n * ediat(ji,jj,jk)) )217 zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen * znanotot ) ) 218 zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpislope2n * zdiattot ) ) 219 219 220 220 ! Computation of production function for Chlorophyll 221 221 !-------------------------------------------------- 222 zprnch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) * zstrn(ji,jj)) )223 zprdch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislope2n * ediat(ji,jj,jk) * zstrn(ji,jj)) )222 zprnch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) ) ) 223 zprdch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislope2n * ediat(ji,jj,jk) ) ) 224 224 ENDIF 225 225 END DO … … 252 252 DO ji = 1, jpi 253 253 254 IF( etot (ji,jj,jk) > 1.E-3 ) THEN254 IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 255 255 ! Si/C of diatoms 256 256 ! ------------------------ … … 258 258 ! Si/C is arbitrariliy increased for very high Si concentrations 259 259 ! to mimic the very high ratios observed in the Southern Ocean (silpot2) 260 zlim = tr n(ji,jj,jk,jpsil) / ( trn(ji,jj,jk,jpsil) + xksi1 )260 zlim = trb(ji,jj,jk,jpsil) / ( trb(ji,jj,jk,jpsil) + xksi1 ) 261 261 zsilim = MIN( zprdia(ji,jj,jk) / ( prmax(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) ) 262 262 zsilfac = 4.4 * EXP( -4.23 * zsilim ) * MAX( 0.e0, MIN( 1., 2.2 * ( zlim - 0.5 ) ) ) + 1.e0 263 zsiborn = tr n(ji,jj,jk,jpsil) * trn(ji,jj,jk,jpsil) * trn(ji,jj,jk,jpsil)263 zsiborn = trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) 264 264 IF (gphit(ji,jj) < -30 ) THEN 265 265 zsilfac2 = 1. + 2. * zsiborn / ( zsiborn + xksi2**3 ) … … 302 302 !CDIR NOVERRCHK 303 303 DO ji = 1, jpi 304 IF( etot (ji,jj,jk) > 1.E-3 ) THEN304 IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 305 305 ! production terms for nanophyto. 306 zprorca(ji,jj,jk) = zprbio(ji,jj,jk) * xlimphy(ji,jj,jk) * tr n(ji,jj,jk,jpphy) * rfact2306 zprorca(ji,jj,jk) = zprbio(ji,jj,jk) * xlimphy(ji,jj,jk) * trb(ji,jj,jk,jpphy) * rfact2 307 307 zpronew(ji,jj,jk) = zprorca(ji,jj,jk) * xnanono3(ji,jj,jk) / ( xnanono3(ji,jj,jk) + xnanonh4(ji,jj,jk) + rtrn ) 308 308 ! 309 zratio = tr n(ji,jj,jk,jpnfe) / ( trn(ji,jj,jk,jpphy) + rtrn )309 zratio = trb(ji,jj,jk,jpnfe) / ( trb(ji,jj,jk,jpphy) + rtrn ) 310 310 zratio = zratio / fecnm 311 311 zmax = MAX( 0., ( 1. - zratio ) / ABS( 1.05 - zratio ) ) … … 313 313 & * ( 4. - 4.5 * xlimnfe(ji,jj,jk) / ( xlimnfe(ji,jj,jk) + 0.5 ) ) & 314 314 & * biron(ji,jj,jk) / ( biron(ji,jj,jk) + concnfe(ji,jj,jk) ) & 315 & * zmax * tr n(ji,jj,jk,jpphy) * rfact2315 & * zmax * trb(ji,jj,jk,jpphy) * rfact2 316 316 ! production terms for diatomees 317 zprorcad(ji,jj,jk) = zprdia(ji,jj,jk) * xlimdia(ji,jj,jk) * tr n(ji,jj,jk,jpdia) * rfact2317 zprorcad(ji,jj,jk) = zprdia(ji,jj,jk) * xlimdia(ji,jj,jk) * trb(ji,jj,jk,jpdia) * rfact2 318 318 zpronewd(ji,jj,jk) = zprorcad(ji,jj,jk) * xdiatno3(ji,jj,jk) / ( xdiatno3(ji,jj,jk) + xdiatnh4(ji,jj,jk) + rtrn ) 319 319 ! 320 zratio = tr n(ji,jj,jk,jpdfe) / ( trn(ji,jj,jk,jpdia) + rtrn )320 zratio = trb(ji,jj,jk,jpdfe) / ( trb(ji,jj,jk,jpdia) + rtrn ) 321 321 zratio = zratio / fecdm 322 322 zmax = MAX( 0., ( 1. - zratio ) / ABS( 1.05 - zratio ) ) … … 324 324 & * ( 4. - 4.5 * xlimdfe(ji,jj,jk) / ( xlimdfe(ji,jj,jk) + 0.5 ) ) & 325 325 & * biron(ji,jj,jk) / ( biron(ji,jj,jk) + concdfe(ji,jj,jk) ) & 326 & * zmax * tr n(ji,jj,jk,jpdia) * rfact2326 & * zmax * trb(ji,jj,jk,jpdia) * rfact2 327 327 ENDIF 328 328 END DO … … 341 341 zprdch(ji,jj,jk) = zprdch(ji,jj,jk) * zmixdiat(ji,jj) 342 342 ENDIF 343 IF( etot (ji,jj,jk) > 1.E-3 ) THEN343 IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 344 344 ! production terms for nanophyto. ( chlorophyll ) 345 345 znanotot = enano(ji,jj,jk) * zstrn(ji,jj) … … 365 365 !CDIR NOVERRCHK 366 366 DO ji = 1, jpi 367 IF( etot (ji,jj,jk) > 1.E-3 ) THEN367 IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 368 368 ! production terms for nanophyto. ( chlorophyll ) 369 znanotot = enano(ji,jj,jk) * zstrn(ji,jj)370 zprod = rday * zprorca(ji,jj,jk) * zprnch(ji,jj,jk) * tr n(ji,jj,jk,jpphy) * xlimphy(ji,jj,jk)369 znanotot = enano(ji,jj,jk) 370 zprod = rday * zprorca(ji,jj,jk) * zprnch(ji,jj,jk) * trb(ji,jj,jk,jpphy) * xlimphy(ji,jj,jk) 371 371 zprochln(ji,jj,jk) = chlcmin * 12. * zprorca (ji,jj,jk) 372 372 zprochln(ji,jj,jk) = zprochln(ji,jj,jk) + (chlcnm-chlcmin) * 144. * zprod & 373 & / ( zpislopead(ji,jj,jk) * tr n(ji,jj,jk,jpnch) * znanotot +rtrn )373 & / ( zpislopead(ji,jj,jk) * trb(ji,jj,jk,jpnch) * znanotot +rtrn ) 374 374 ! production terms for diatomees ( chlorophyll ) 375 zdiattot = ediat(ji,jj,jk) * zstrn(ji,jj)376 zprod = rday * zprorcad(ji,jj,jk) * zprdch(ji,jj,jk) * tr n(ji,jj,jk,jpdia) * xlimdia(ji,jj,jk)375 zdiattot = ediat(ji,jj,jk) 376 zprod = rday * zprorcad(ji,jj,jk) * zprdch(ji,jj,jk) * trb(ji,jj,jk,jpdia) * xlimdia(ji,jj,jk) 377 377 zprochld(ji,jj,jk) = chlcmin * 12. * zprorcad(ji,jj,jk) 378 378 zprochld(ji,jj,jk) = zprochld(ji,jj,jk) + (chlcdm-chlcmin) * 144. * zprod & 379 & / ( zpislopead2(ji,jj,jk) * tr n(ji,jj,jk,jpdch) * zdiattot +rtrn )379 & / ( zpislopead2(ji,jj,jk) * trb(ji,jj,jk,jpdch) * zdiattot +rtrn ) 380 380 ENDIF 381 381 END DO … … 412 412 END DO 413 413 414 ! Total primary production per year 415 tpp = tpp + glob_sum( ( zprorca(:,:,:) + zprorcad(:,:,:) ) * cvol(:,:,:) ) 416 417 IF( ln_diatrc ) THEN 418 ! 419 zrfact2 = 1.e3 * rfact2r ! conversion from mol/L/timestep into mol/m3/s 420 IF( lk_iomput ) THEN 421 IF( jnt == nrdttrc ) THEN 422 CALL iom_put( "PPPHY" , zprorca (:,:,:) * zrfact2 * tmask(:,:,:) ) ! primary production by nanophyto 423 CALL iom_put( "PPPHY2" , zprorcad(:,:,:) * zrfact2 * tmask(:,:,:) ) ! primary production by diatom 424 CALL iom_put( "PPNEWN" , zpronew (:,:,:) * zrfact2 * tmask(:,:,:) ) ! new primary production by nanophyto 425 CALL iom_put( "PPNEWD" , zpronewd(:,:,:) * zrfact2 * tmask(:,:,:) ) ! new primary production by diatom 426 CALL iom_put( "PBSi" , zprorcad(:,:,:) * zrfact2 * tmask(:,:,:) * zysopt(:,:,:) ) ! biogenic silica production 427 CALL iom_put( "PFeD" , zprofed (:,:,:) * zrfact2 * tmask(:,:,:) ) ! biogenic iron production by diatom 428 CALL iom_put( "PFeN" , zprofen (:,:,:) * zrfact2 * tmask(:,:,:) ) ! biogenic iron production by nanophyto 429 CALL iom_put( "Mumax" , prmax(:,:,:) * tmask(:,:,:) ) ! Maximum growth rate 430 CALL iom_put( "MuN" , zprbio(:,:,:) * xlimphy(:,:,:) * tmask(:,:,:) ) ! Realized growth rate for nanophyto 431 CALL iom_put( "MuD" , zprdia(:,:,:) * xlimdia(:,:,:) * tmask(:,:,:) ) ! Realized growth rate for diatoms 432 CALL iom_put( "LNlight", zprbio (:,:,:) / (prmax(:,:,:) + rtrn) * tmask(:,:,:) ) ! light limitation term 433 CALL iom_put( "LDlight", zprdia (:,:,:) / (prmax(:,:,:) + rtrn) * tmask(:,:,:) ) ! light limitation term 434 ENDIF 435 ELSE 436 trc3d(:,:,:,jp_pcs0_3d + 4) = zprorca (:,:,:) * zrfact2 * tmask(:,:,:) 437 trc3d(:,:,:,jp_pcs0_3d + 5) = zprorcad(:,:,:) * zrfact2 * tmask(:,:,:) 438 trc3d(:,:,:,jp_pcs0_3d + 6) = zpronew (:,:,:) * zrfact2 * tmask(:,:,:) 439 trc3d(:,:,:,jp_pcs0_3d + 7) = zpronewd(:,:,:) * zrfact2 * tmask(:,:,:) 440 trc3d(:,:,:,jp_pcs0_3d + 8) = zprorcad(:,:,:) * zrfact2 * tmask(:,:,:) * zysopt(:,:,:) 441 trc3d(:,:,:,jp_pcs0_3d + 9) = zprofed (:,:,:) * zrfact2 * tmask(:,:,:) 414 415 ! Total primary production per year 416 IF( iom_use( "tintpp" ) .OR. ( ln_check_mass .AND. kt == nitend .AND. knt == nrdttrc ) ) & 417 & tpp = glob_sum( ( zprorca(:,:,:) + zprorcad(:,:,:) ) * cvol(:,:,:) ) 418 419 IF( lk_iomput ) THEN 420 IF( knt == nrdttrc ) THEN 421 CALL wrk_alloc( jpi, jpj, zw2d ) 422 CALL wrk_alloc( jpi, jpj, jpk, zw3d ) 423 zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s 424 ! 425 IF( iom_use( "PPPHY" ) .OR. iom_use( "PPPHY2" ) ) THEN 426 zw3d(:,:,:) = zprorca (:,:,:) * zfact * tmask(:,:,:) ! primary production by nanophyto 427 CALL iom_put( "PPPHY" , zw3d ) 428 ! 429 zw3d(:,:,:) = zprorcad(:,:,:) * zfact * tmask(:,:,:) ! primary production by diatomes 430 CALL iom_put( "PPPHY2" , zw3d ) 431 ENDIF 432 IF( iom_use( "PPNEWN" ) .OR. iom_use( "PPNEWD" ) ) THEN 433 zw3d(:,:,:) = zpronew (:,:,:) * zfact * tmask(:,:,:) ! new primary production by nanophyto 434 CALL iom_put( "PPNEWN" , zw3d ) 435 ! 436 zw3d(:,:,:) = zpronewd(:,:,:) * zfact * tmask(:,:,:) ! new primary production by diatomes 437 CALL iom_put( "PPNEWD" , zw3d ) 438 ENDIF 439 IF( iom_use( "PBSi" ) ) THEN 440 zw3d(:,:,:) = zprorcad(:,:,:) * zfact * tmask(:,:,:) * zysopt(:,:,:) ! biogenic silica production 441 CALL iom_put( "PBSi" , zw3d ) 442 ENDIF 443 IF( iom_use( "PFeN" ) .OR. iom_use( "PFeD" ) ) THEN 444 zw3d(:,:,:) = zprofen(:,:,:) * zfact * tmask(:,:,:) ! biogenic iron production by nanophyto 445 CALL iom_put( "PFeN" , zw3d ) 446 ! 447 zw3d(:,:,:) = zprofed(:,:,:) * zfact * tmask(:,:,:) ! biogenic iron production by diatomes 448 CALL iom_put( "PFeD" , zw3d ) 449 ENDIF 450 IF( iom_use( "Mumax" ) ) THEN 451 zw3d(:,:,:) = prmax(:,:,:) * tmask(:,:,:) ! Maximum growth rate 452 CALL iom_put( "Mumax" , zw3d ) 453 ENDIF 454 IF( iom_use( "MuN" ) .OR. iom_use( "MuD" ) ) THEN 455 zw3d(:,:,:) = zprbio(:,:,:) * xlimphy(:,:,:) * tmask(:,:,:) ! Realized growth rate for nanophyto 456 CALL iom_put( "MuN" , zw3d ) 457 ! 458 zw3d(:,:,:) = zprdia(:,:,:) * xlimdia(:,:,:) * tmask(:,:,:) ! Realized growth rate for diatoms 459 CALL iom_put( "MuD" , zw3d ) 460 ENDIF 461 IF( iom_use( "LNlight" ) .OR. iom_use( "LDlight" ) ) THEN 462 zw3d(:,:,:) = zprbio (:,:,:) / (prmax(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 463 CALL iom_put( "LNlight" , zw3d ) 464 ! 465 zw3d(:,:,:) = zprdia (:,:,:) / (prmax(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 466 CALL iom_put( "LDlight" , zw3d ) 467 ENDIF 468 IF( iom_use( "TPP" ) ) THEN 469 zw3d(:,:,:) = ( zprorca(:,:,:) + zprorcad(:,:,:) ) * zfact * tmask(:,:,:) ! total primary production 470 CALL iom_put( "TPP" , zw3d ) 471 ENDIF 472 IF( iom_use( "TPNEW" ) ) THEN 473 zw3d(:,:,:) = ( zpronew(:,:,:) + zpronewd(:,:,:) ) * zfact * tmask(:,:,:) ! total new production 474 CALL iom_put( "TPNEW" , zw3d ) 475 ENDIF 476 IF( iom_use( "TPBFE" ) ) THEN 477 zw3d(:,:,:) = ( zprofen(:,:,:) + zprofed(:,:,:) ) * zfact * tmask(:,:,:) ! total biogenic iron production 478 CALL iom_put( "TPBFE" , zw3d ) 479 ENDIF 480 IF( iom_use( "INTPPPHY" ) .OR. iom_use( "INTPPPHY2" ) ) THEN 481 zw2d(:,:) = 0. 482 DO jk = 1, jpkm1 483 zw2d(:,:) = zw2d(:,:) + zprorca (:,:,jk) * fse3t(:,:,jk) * zfact * tmask(:,:,jk) ! vert. integrated primary produc. by nano 484 ENDDO 485 CALL iom_put( "INTPPPHY" , zw2d ) 486 ! 487 zw2d(:,:) = 0. 488 DO jk = 1, jpkm1 489 zw2d(:,:) = zw2d(:,:) + zprorcad(:,:,jk) * fse3t(:,:,jk) * zfact * tmask(:,:,jk) ! vert. integrated primary produc. by diatom 490 ENDDO 491 CALL iom_put( "INTPPPHY2" , zw2d ) 492 ENDIF 493 IF( iom_use( "INTPP" ) ) THEN 494 zw2d(:,:) = 0. 495 DO jk = 1, jpkm1 496 zw2d(:,:) = zw2d(:,:) + ( zprorca(:,:,jk) + zprorcad(:,:,jk) ) * fse3t(:,:,jk) * zfact * tmask(:,:,jk) ! vert. integrated pp 497 ENDDO 498 CALL iom_put( "INTPP" , zw2d ) 499 ENDIF 500 IF( iom_use( "INTPNEW" ) ) THEN 501 zw2d(:,:) = 0. 502 DO jk = 1, jpkm1 503 zw2d(:,:) = zw2d(:,:) + ( zpronew(:,:,jk) + zpronewd(:,:,jk) ) * fse3t(:,:,jk) * zfact * tmask(:,:,jk) ! vert. integrated new prod 504 ENDDO 505 CALL iom_put( "INTPNEW" , zw2d ) 506 ENDIF 507 IF( iom_use( "INTPBFE" ) ) THEN ! total biogenic iron production ( vertically integrated ) 508 zw2d(:,:) = 0. 509 DO jk = 1, jpkm1 510 zw2d(:,:) = zw2d(:,:) + ( zprofen(:,:,jk) + zprofed(:,:,jk) ) * fse3t(:,:,jk) * zfact * tmask(:,:,jk) ! vert integr. bfe prod 511 ENDDO 512 CALL iom_put( "INTPBFE" , zw2d ) 513 ENDIF 514 IF( iom_use( "INTPBSI" ) ) THEN ! total biogenic silica production ( vertically integrated ) 515 zw2d(:,:) = 0. 516 DO jk = 1, jpkm1 517 zw2d(:,:) = zw2d(:,:) + zprorcad(:,:,jk) * zysopt(:,:,jk) * fse3t(:,:,jk) * zfact * tmask(:,:,jk) ! vert integr. bsi prod 518 ENDDO 519 CALL iom_put( "INTPBSI" , zw2d ) 520 ENDIF 521 IF( iom_use( "tintpp" ) ) CALL iom_put( "tintpp" , tpp * zfact ) ! global total integrated primary production molC/s 522 ! 523 CALL wrk_dealloc( jpi, jpj, zw2d ) 524 CALL wrk_dealloc( jpi, jpj, jpk, zw3d ) 525 ENDIF 526 ELSE 527 IF( ln_diatrc ) THEN 528 zfact = 1.e+3 * rfact2r 529 trc3d(:,:,:,jp_pcs0_3d + 4) = zprorca (:,:,:) * zfact * tmask(:,:,:) 530 trc3d(:,:,:,jp_pcs0_3d + 5) = zprorcad(:,:,:) * zfact * tmask(:,:,:) 531 trc3d(:,:,:,jp_pcs0_3d + 6) = zpronew (:,:,:) * zfact * tmask(:,:,:) 532 trc3d(:,:,:,jp_pcs0_3d + 7) = zpronewd(:,:,:) * zfact * tmask(:,:,:) 533 trc3d(:,:,:,jp_pcs0_3d + 8) = zprorcad(:,:,:) * zfact * tmask(:,:,:) * zysopt(:,:,:) 534 trc3d(:,:,:,jp_pcs0_3d + 9) = zprofed (:,:,:) * zfact * tmask(:,:,:) 442 535 # if ! defined key_kriest 443 trc3d(:,:,:,jp_pcs0_3d + 10) = zprofen (:,:,:) * zrfact2* tmask(:,:,:)536 trc3d(:,:,:,jp_pcs0_3d + 10) = zprofen (:,:,:) * zfact * tmask(:,:,:) 444 537 # endif 445 ENDIF 446 ! 447 ENDIF 448 449 IF(ln_ctl) THEN ! print mean trends (used for debugging) 538 ENDIF 539 ENDIF 540 541 IF(ln_ctl) THEN ! print mean trends (used for debugging) 450 542 WRITE(charout, FMT="('prod')") 451 543 CALL prt_ctl_trc_info(charout) 452 544 CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm) 453 454 455 456 457 458 459 460 545 ENDIF 546 ! 547 CALL wrk_dealloc( jpi, jpj, zmixnano, zmixdiat, zstrn ) 548 CALL wrk_dealloc( jpi, jpj, jpk, zpislopead, zpislopead2, zprdia, zprbio, zprdch, zprnch, zysopt ) 549 CALL wrk_dealloc( jpi, jpj, jpk, zprorca, zprorcad, zprofed, zprofen, zprochln, zprochld, zpronew, zpronewd ) 550 ! 551 IF( nn_timing == 1 ) CALL timing_stop('p4z_prod') 552 ! 461 553 END SUBROUTINE p4z_prod 462 554 -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zrem.F90
r4624 r5581 59 59 CONTAINS 60 60 61 SUBROUTINE p4z_rem( kt, jnt )61 SUBROUTINE p4z_rem( kt, knt ) 62 62 !!--------------------------------------------------------------------- 63 63 !! *** ROUTINE p4z_rem *** … … 68 68 !!--------------------------------------------------------------------- 69 69 ! 70 INTEGER, INTENT(in) :: kt, jnt ! ocean time step70 INTEGER, INTENT(in) :: kt, knt ! ocean time step 71 71 ! 72 72 INTEGER :: ji, jj, jk … … 78 78 REAL(wp) :: zofer2 79 79 #endif 80 REAL(wp) :: zonitr, zstep, z rfact280 REAL(wp) :: zonitr, zstep, zfact 81 81 CHARACTER (len=25) :: charout 82 REAL(wp), POINTER, DIMENSION(:,: ) :: ztempbac 83 REAL(wp), POINTER, DIMENSION(:,:,:) :: zdepbac, zolimi, zdepprod 82 REAL(wp), POINTER, DIMENSION(:,: ) :: ztempbac 83 REAL(wp), POINTER, DIMENSION(:,:,:) :: zdepbac, zolimi, zdepprod, zw3d 84 84 !!--------------------------------------------------------------------- 85 85 ! … … 104 104 zdep = MAX( hmld(ji,jj), heup(ji,jj) ) 105 105 IF( fsdept(ji,jj,jk) < zdep ) THEN 106 zdepbac(ji,jj,jk) = MIN( 0.7 * ( tr n(ji,jj,jk,jpzoo) + 2.* trn(ji,jj,jk,jpmes) ), 4.e-6 )106 zdepbac(ji,jj,jk) = MIN( 0.7 * ( trb(ji,jj,jk,jpzoo) + 2.* trb(ji,jj,jk,jpmes) ), 4.e-6 ) 107 107 ztempbac(ji,jj) = zdepbac(ji,jj,jk) 108 108 ELSE … … 119 119 DO ji = 1, jpi 120 120 ! denitrification factor computed from O2 levels 121 nitrfac(ji,jj,jk) = MAX( 0.e0, 0.4 * ( 6.e-6 - tr n(ji,jj,jk,jpoxy) ) &122 & / ( oxymin + tr n(ji,jj,jk,jpoxy) ) )121 nitrfac(ji,jj,jk) = MAX( 0.e0, 0.4 * ( 6.e-6 - trb(ji,jj,jk,jpoxy) ) & 122 & / ( oxymin + trb(ji,jj,jk,jpoxy) ) ) 123 123 nitrfac(ji,jj,jk) = MIN( 1., nitrfac(ji,jj,jk) ) 124 124 END DO … … 140 140 ! Ammonification in oxic waters with oxygen consumption 141 141 ! ----------------------------------------------------- 142 zolimit = zremik * ( 1.- nitrfac(ji,jj,jk) ) * tr n(ji,jj,jk,jpdoc)143 zolimi(ji,jj,jk) = MIN( ( tr n(ji,jj,jk,jpoxy) - rtrn ) / o2ut, zolimit )142 zolimit = zremik * ( 1.- nitrfac(ji,jj,jk) ) * trb(ji,jj,jk,jpdoc) 143 zolimi(ji,jj,jk) = MIN( ( trb(ji,jj,jk,jpoxy) - rtrn ) / o2ut, zolimit ) 144 144 ! Ammonification in suboxic waters with denitrification 145 145 ! ------------------------------------------------------- 146 denitr(ji,jj,jk) = MIN( ( tr n(ji,jj,jk,jpno3) - rtrn ) / rdenit, &147 & zremik * nitrfac(ji,jj,jk) * tr n(ji,jj,jk,jpdoc) )146 denitr(ji,jj,jk) = MIN( ( trb(ji,jj,jk,jpno3) - rtrn ) / rdenit, & 147 & zremik * nitrfac(ji,jj,jk) * trb(ji,jj,jk,jpdoc) ) 148 148 ! 149 149 zolimi (ji,jj,jk) = MAX( 0.e0, zolimi (ji,jj,jk) ) … … 165 165 ! below 2 umol/L. Inhibited at strong light 166 166 ! ---------------------------------------------------------- 167 zonitr =nitrif * zstep * tr n(ji,jj,jk,jpnh4) / ( 1.+ emoy(ji,jj,jk) ) * ( 1.- nitrfac(ji,jj,jk) )168 denitnh4(ji,jj,jk) = nitrif * zstep * tr n(ji,jj,jk,jpnh4) * nitrfac(ji,jj,jk)167 zonitr =nitrif * zstep * trb(ji,jj,jk,jpnh4) / ( 1.+ emoy(ji,jj,jk) ) * ( 1.- nitrfac(ji,jj,jk) ) 168 denitnh4(ji,jj,jk) = nitrif * zstep * trb(ji,jj,jk,jpnh4) * nitrfac(ji,jj,jk) 169 169 ! Update of the tracers trends 170 170 ! ---------------------------- … … 192 192 ! ---------------------------------------------------------- 193 193 zbactfer = 10.e-6 * rfact2 * prmax(ji,jj,jk) * xlimbacl(ji,jj,jk) & 194 & * tr n(ji,jj,jk,jpfer) / ( 2.5E-10 + trn(ji,jj,jk,jpfer) ) &194 & * trb(ji,jj,jk,jpfer) / ( 2.5E-10 + trb(ji,jj,jk,jpfer) ) & 195 195 & * zdepprod(ji,jj,jk) * zdepbac(ji,jj,jk) 196 196 #if defined key_kriest … … 228 228 ! means a disaggregation constant about 0.5 the value in oxic zones 229 229 ! ----------------------------------------------------------------- 230 zorem = zremip * tr n(ji,jj,jk,jppoc)231 zofer = zremip * tr n(ji,jj,jk,jpsfe)230 zorem = zremip * trb(ji,jj,jk,jppoc) 231 zofer = zremip * trb(ji,jj,jk,jpsfe) 232 232 #if ! defined key_kriest 233 zorem2 = zremip * tr n(ji,jj,jk,jpgoc)234 zofer2 = zremip * tr n(ji,jj,jk,jpbfe)233 zorem2 = zremip * trb(ji,jj,jk,jpgoc) 234 zofer2 = zremip * trb(ji,jj,jk,jpbfe) 235 235 #else 236 zorem2 = zremip * tr n(ji,jj,jk,jpnum)236 zorem2 = zremip * trb(ji,jj,jk,jpnum) 237 237 #endif 238 238 … … 272 272 ! Remineralization rate of BSi depedant on T and saturation 273 273 ! --------------------------------------------------------- 274 zsatur = ( sio3eq(ji,jj,jk) - tr n(ji,jj,jk,jpsil) ) / ( sio3eq(ji,jj,jk) + rtrn )274 zsatur = ( sio3eq(ji,jj,jk) - trb(ji,jj,jk,jpsil) ) / ( sio3eq(ji,jj,jk) + rtrn ) 275 275 zsatur = MAX( rtrn, zsatur ) 276 276 zsatur2 = ( 1. + tsn(ji,jj,jk,jp_tem) / 400.)**37 … … 287 287 zfactdep = xsilab * EXP(-( xsiremlab - xsirem ) * znusil2 * zdep / wsbio2 ) * ztem / ( ztem + 10. ) 288 288 zsiremin = ( xsiremlab * zfactdep + xsirem * ( 1. - zfactdep ) ) * zstep * znusil 289 zosil = zsiremin * tr n(ji,jj,jk,jpgsi)289 zosil = zsiremin * trb(ji,jj,jk,jpgsi) 290 290 ! 291 291 tra(ji,jj,jk,jpgsi) = tra(ji,jj,jk,jpgsi) - zosil … … 315 315 END DO 316 316 317 IF( ln_diatrc .AND. lk_iomput .AND. jnt == nrdttrc ) THEN 318 zrfact2 = 1.e3 * rfact2r 319 CALL iom_put( "REMIN" , zolimi(:,:,:) * tmask(:,:,:) * zrfact2 ) ! Remineralisation rate 320 CALL iom_put( "DENIT" , denitr(:,:,:) * rdenit * rno3 * tmask(:,:,:) * zrfact2 ) ! Denitrification 321 ENDIF 317 IF( knt == nrdttrc ) THEN 318 CALL wrk_alloc( jpi, jpj, jpk, zw3d ) 319 zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s 320 ! 321 IF( iom_use( "REMIN" ) ) THEN 322 zw3d(:,:,:) = zolimi(:,:,:) * tmask(:,:,:) * zfact ! Remineralisation rate 323 CALL iom_put( "REMIN" , zw3d ) 324 ENDIF 325 IF( iom_use( "DENIT" ) ) THEN 326 zw3d(:,:,:) = denitr(:,:,:) * rdenit * rno3 * tmask(:,:,:) * zfact ! Denitrification 327 CALL iom_put( "DENIT" , zw3d ) 328 ENDIF 329 ! 330 CALL wrk_dealloc( jpi, jpj, jpk, zw3d ) 331 ENDIF 322 332 323 333 IF(ln_ctl) THEN ! print mean trends (used for debugging) -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsbc.F90
- Property svn:keywords set to Id
r4793 r5581 80 80 REAL(wp), PUBLIC :: rivdininput, rivdipinput, rivdsiinput 81 81 82 REAL(wp) :: ryyss !: number of seconds per year83 82 84 83 !!* Substitution … … 86 85 !!---------------------------------------------------------------------- 87 86 !! NEMO/TOP 3.3 , NEMO Consortium (2010) 88 !! $ Header:$87 !! $Id$ 89 88 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 90 89 !!---------------------------------------------------------------------- … … 118 117 IF( kt == nit000 .OR. ( kt /= nit000 .AND. ntimes_dust > 1 ) ) THEN 119 118 CALL fld_read( kt, 1, sf_dust ) 120 dust(:,:) = sf_dust(1)%fnow(:,:,1) 119 IF( nn_ice_tr == -1 .AND. .NOT. ln_ironice ) THEN 120 dust(:,:) = sf_dust(1)%fnow(:,:,1) 121 ELSE 122 dust(:,:) = sf_dust(1)%fnow(:,:,1) * ( 1.0 - fr_i(:,:) ) 123 ENDIF 121 124 ENDIF 122 125 ENDIF … … 137 140 DO jj = 1, jpj 138 141 DO ji = 1, jpi 139 zcoef = ryyss * cvol(ji,jj,1)142 zcoef = ryyss * e1e2t(ji,jj) * h_rnf(ji,jj) 140 143 rivalk(ji,jj) = sf_river(jr_dic)%fnow(ji,jj,1) & 141 144 & * 1.E3 / ( 12. * zcoef + rtrn ) … … 188 191 INTEGER :: ierr, ierr1, ierr2, ierr3 189 192 INTEGER :: ios ! Local integer output status for namelist read 193 INTEGER :: ik50 ! last level where depth less than 50 m 194 INTEGER :: isrow ! index for ORCA1 starting row 190 195 REAL(wp) :: zexpide, zdenitide, zmaskt 191 196 REAL(wp) :: ztimes_dust, ztimes_riv, ztimes_ndep … … 208 213 IF( nn_timing == 1 ) CALL timing_start('p4z_sbc_init') 209 214 ! 210 ryyss = nyear_len(1) * rday ! number of seconds per year and per month211 !212 215 ! !* set file information 213 216 REWIND( numnatp_ref ) ! Namelist nampissbc in reference namelist : Pisces external sources of nutrients … … 219 222 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampissbc in configuration namelist', lwp ) 220 223 IF(lwm) WRITE ( numonp, nampissbc ) 224 225 IF ( ( nn_ice_tr >= 0 ) .AND. ln_ironice ) THEN 226 IF(lwp) THEN 227 WRITE(numout,*) ' ln_ironice incompatible with nn_ice_tr = ', nn_ice_tr 228 WRITE(numout,*) ' Specify your sea ice iron concentration in nampisice instead ' 229 WRITE(numout,*) ' ln_ironice is forced to .FALSE. ' 230 ln_ironice = .FALSE. 231 ENDIF 232 ENDIF 221 233 222 234 IF(lwp) THEN … … 250 262 ENDIF 251 263 264 ! set the number of level over which river runoffs are applied 265 ! online configuration : computed in sbcrnf 266 IF( lk_offline ) THEN 267 nk_rnf(:,:) = 1 268 h_rnf (:,:) = fsdept(:,:,1) 269 ENDIF 270 252 271 ! dust input from the atmosphere 253 272 ! ------------------------------ … … 361 380 rivalkinput = 0._wp 362 381 END IF 363 364 382 ! nutrient input from dust 365 383 ! ------------------------ … … 413 431 CALL iom_close( numiron ) 414 432 ! 415 DO jk = 1, 5 433 ik50 = 5 ! last level where depth less than 50 m 434 DO jk = jpkm1, 1, -1 435 IF( gdept_1d(jk) > 50. ) ik50 = jk - 1 436 END DO 437 IF (lwp) WRITE(numout,*) 438 IF (lwp) WRITE(numout,*) ' Level corresponding to 50m depth ', ik50,' ', gdept_1d(ik50+1) 439 IF (lwp) WRITE(numout,*) 440 DO jk = 1, ik50 416 441 DO jj = 2, jpjm1 417 442 DO ji = fs_2, fs_jpim1 … … 424 449 END DO 425 450 END DO 426 IF( cp_cfg == 'orca' .AND. jp_cfg == 2 ) THEN 427 ii0 = 176 ; ii1 = 176 ! Southern Island : Kerguelen 428 ij0 = 37 ; ij1 = 37 ; zcmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 0.3_wp 429 ! 430 ii0 = 119 ; ii1 = 119 ! South Georgia 431 ij0 = 29 ; ij1 = 29 ; zcmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 0.3_wp 432 ! 433 ii0 = 111 ; ii1 = 111 ! Falklands 434 ij0 = 35 ; ij1 = 35 ; zcmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 0.3_wp 435 ! 436 ii0 = 168 ; ii1 = 168 ! Crozet 437 ij0 = 40 ; ij1 = 40 ; zcmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 0.3_wp 438 ! 439 ii0 = 119 ; ii1 = 119 ! South Orkney 440 ij0 = 28 ; ij1 = 28 ; zcmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 0.3_wp 441 ! 442 ii0 = 140 ; ii1 = 140 ! Bouvet Island 443 ij0 = 33 ; ij1 = 33 ; zcmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 0.3_wp 444 ! 445 ii0 = 178 ; ii1 = 178 ! Prince edwards 446 ij0 = 34 ; ij1 = 34 ; zcmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 0.3_wp 447 ! 448 ii0 = 43 ; ii1 = 43 ! Balleny islands 449 ij0 = 21 ; ij1 = 21 ; zcmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 0.3_wp 450 ENDIF 451 ! 451 452 CALL lbc_lnk( zcmask , 'T', 1. ) ! lateral boundary conditions on cmask (sign unchanged) 453 ! 452 454 DO jk = 1, jpk 453 455 DO jj = 1, jpj -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsed.F90
- Property svn:keywords set to Id
r4793 r5581 21 21 USE p4zopt ! optical model 22 22 USE p4zlim ! Co-limitations of differents nutrients 23 USE p4zrem ! Remineralisation of organic matter24 23 USE p4zsbc ! External source of nutrients 25 24 USE p4zint ! interpolation and computation of various fields … … 30 29 PRIVATE 31 30 32 PUBLIC p4z_sed 31 PUBLIC p4z_sed 32 PUBLIC p4z_sed_alloc 33 33 34 34 35 !! * Module variables 35 REAL(wp) :: ryyss !: number of seconds per year 36 REAL(wp) :: r1_ryyss !: inverse of ryyss 37 REAL(wp) :: rmtss !: number of seconds per month 36 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: nitrpot !: Nitrogen fixation 37 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,: ) :: sdenit !: Nitrate reduction in the sediments 38 38 REAL(wp) :: r1_rday !: inverse of rday 39 40 INTEGER :: numnit41 42 39 43 40 !!* Substitution … … 45 42 !!---------------------------------------------------------------------- 46 43 !! NEMO/TOP 3.3 , NEMO Consortium (2010) 47 !! $ Header:$44 !! $Id$ 48 45 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 49 46 !!---------------------------------------------------------------------- 50 47 CONTAINS 51 48 52 SUBROUTINE p4z_sed( kt, jnt )49 SUBROUTINE p4z_sed( kt, knt ) 53 50 !!--------------------------------------------------------------------- 54 51 !! *** ROUTINE p4z_sed *** … … 61 58 !!--------------------------------------------------------------------- 62 59 ! 63 INTEGER, INTENT(in) :: kt, jnt ! ocean time step60 INTEGER, INTENT(in) :: kt, knt ! ocean time step 64 61 INTEGER :: ji, jj, jk, ikt 65 62 #if ! defined key_sed … … 72 69 REAL(wp) :: zsiloss, zcaloss, zws3, zws4, zwsc, zdep, zwstpoc 73 70 REAL(wp) :: ztrfer, ztrpo4, zwdust, zlight 74 REAL(wp) :: zrdenittot, zsdenittot, znitrpottot75 71 ! 76 72 CHARACTER (len=25) :: charout 77 REAL(wp), POINTER, DIMENSION(:,: ) :: zpdep, zsidep, zwork1, zwork2, zwork3 , zwork473 REAL(wp), POINTER, DIMENSION(:,: ) :: zpdep, zsidep, zwork1, zwork2, zwork3 78 74 REAL(wp), POINTER, DIMENSION(:,: ) :: zdenit2d, zironice, zbureff 79 75 REAL(wp), POINTER, DIMENSION(:,: ) :: zwsbio3, zwsbio4, zwscal 80 REAL(wp), POINTER, DIMENSION(:,:,:) :: z nitrpot, zirondep, zsoufer76 REAL(wp), POINTER, DIMENSION(:,:,:) :: zirondep, zsoufer 81 77 !!--------------------------------------------------------------------- 82 78 ! 83 79 IF( nn_timing == 1 ) CALL timing_start('p4z_sed') 84 80 ! 85 IF( kt == nittrc000 .AND. jnt == 1 ) THEN 86 ryyss = nyear_len(1) * rday ! number of seconds per year and per month 87 rmtss = ryyss / raamo 88 r1_rday = 1. / rday 89 r1_ryyss = 1. / ryyss 90 IF( ln_check_mass .AND. lwp) & 91 & CALL ctl_opn( numnit, 'nitrogen.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea ) 92 ENDIF 81 IF( kt == nittrc000 .AND. knt == 1 ) r1_rday = 1. / rday 93 82 ! 94 83 ! Allocate temporary workspace 95 CALL wrk_alloc( jpi, jpj, zdenit2d, zwork1, zwork2, zwork3, z work4, zbureff )84 CALL wrk_alloc( jpi, jpj, zdenit2d, zwork1, zwork2, zwork3, zbureff ) 96 85 CALL wrk_alloc( jpi, jpj, zwsbio3, zwsbio4, zwscal ) 97 CALL wrk_alloc( jpi, jpj, jpk, z nitrpot, zsoufer )86 CALL wrk_alloc( jpi, jpj, jpk, zsoufer ) 98 87 99 88 zdenit2d(:,:) = 0.e0 100 89 zbureff (:,:) = 0.e0 90 zwork1 (:,:) = 0.e0 91 zwork2 (:,:) = 0.e0 92 zwork3 (:,:) = 0.e0 101 93 102 94 ! Iron input/uptake due to sea ice : Crude parameterization based on Lancelot et al. … … 110 102 zdep = rfact2 / fse3t(ji,jj,1) 111 103 zwflux = fmmflx(ji,jj) / 1000._wp 112 zfminus = MIN( 0._wp, -zwflux ) * tr n(ji,jj,1,jpfer) * zdep104 zfminus = MIN( 0._wp, -zwflux ) * trb(ji,jj,1,jpfer) * zdep 113 105 zfplus = MAX( 0._wp, -zwflux ) * icefeinput * zdep 114 106 zironice(ji,jj) = zfplus + zfminus … … 116 108 END DO 117 109 ! 118 tr n(:,:,1,jpfer) = trn(:,:,1,jpfer) + zironice(:,:)119 ! 120 IF( l n_diatrc .AND. lk_iomput .AND. jnt == nrdttrc) &110 tra(:,:,1,jpfer) = tra(:,:,1,jpfer) + zironice(:,:) 111 ! 112 IF( lk_iomput .AND. knt == nrdttrc .AND. iom_use( "Ironice" ) ) & 121 113 & CALL iom_put( "Ironice", zironice(:,:) * 1.e+3 * rfact2r * fse3t(:,:,1) * tmask(:,:,1) ) ! iron flux from ice 114 ! 122 115 CALL wrk_dealloc( jpi, jpj, zironice ) 123 116 ! … … 132 125 ! ! Iron and Si deposition at the surface 133 126 IF( ln_solub ) THEN 134 zirondep(:,:,1) = solub(:,:) * dust(:,:) * mfrac * rfact2 / fse3t(:,:,1) / ( 55.85 * rmtss )+ 3.e-10 * r1_ryyss127 zirondep(:,:,1) = solub(:,:) * dust(:,:) * mfrac * rfact2 / fse3t(:,:,1) / 55.85 + 3.e-10 * r1_ryyss 135 128 ELSE 136 zirondep(:,:,1) = dustsolub * dust(:,:) * mfrac * rfact2 / fse3t(:,:,1) / ( 55.85 * rmtss )+ 3.e-10 * r1_ryyss129 zirondep(:,:,1) = dustsolub * dust(:,:) * mfrac * rfact2 / fse3t(:,:,1) / 55.85 + 3.e-10 * r1_ryyss 137 130 ENDIF 138 zsidep(:,:) = 8.8 * 0.075 * dust(:,:) * mfrac * rfact2 / fse3t(:,:,1) / ( 28.1 * rmtss )139 zpdep (:,:) = 0.1 * 0.021 * dust(:,:) * mfrac * rfact2 / fse3t(:,:,1) / ( 31. * rmtss )/ po4r131 zsidep(:,:) = 8.8 * 0.075 * dust(:,:) * mfrac * rfact2 / fse3t(:,:,1) / 28.1 132 zpdep (:,:) = 0.1 * 0.021 * dust(:,:) * mfrac * rfact2 / fse3t(:,:,1) / 31. / po4r 140 133 ! ! Iron solubilization of particles in the water column 141 134 ! ! dust in kg/m2/s ---> 1/55.85 to put in mol/Fe ; wdust in m/j … … 145 138 END DO 146 139 ! ! Iron solubilization of particles in the water column 147 trn(:,:,1,jppo4) = trn(:,:,1,jppo4) + zpdep (:,:) 148 trn(:,:,1,jpsil) = trn(:,:,1,jpsil) + zsidep (:,:) 149 trn(:,:,:,jpfer) = trn(:,:,:,jpfer) + zirondep(:,:,:) 150 ! 151 IF( ln_diatrc ) THEN 152 zfact = 1.e+3 * rfact2r 153 IF( lk_iomput ) THEN 154 IF( jnt == nrdttrc ) THEN 155 CALL iom_put( "Irondep", zirondep(:,:,1) * zfact * fse3t(:,:,1) * tmask(:,:,1) ) ! surface downward dust depo of iron 156 CALL iom_put( "pdust" , dust(:,:) / ( wdust * rday ) * tmask(:,:,1) ) ! dust concentration at surface 157 ENDIF 158 ELSE 159 trc2d(:,:,jp_pcs0_2d + 11) = zirondep(:,:,1) * zfact * fse3t(:,:,1) * tmask(:,:,1) 140 tra(:,:,1,jppo4) = tra(:,:,1,jppo4) + zpdep (:,:) 141 tra(:,:,1,jpsil) = tra(:,:,1,jpsil) + zsidep (:,:) 142 tra(:,:,:,jpfer) = tra(:,:,:,jpfer) + zirondep(:,:,:) 143 ! 144 IF( lk_iomput ) THEN 145 IF( knt == nrdttrc ) THEN 146 IF( iom_use( "Irondep" ) ) & 147 & CALL iom_put( "Irondep", zirondep(:,:,1) * 1.e+3 * rfact2r * fse3t(:,:,1) * tmask(:,:,1) ) ! surface downward dust depo of iron 148 IF( iom_use( "pdust" ) ) & 149 & CALL iom_put( "pdust" , dust(:,:) / ( wdust * rday ) * tmask(:,:,1) ) ! dust concentration at surface 160 150 ENDIF 151 ELSE 152 IF( ln_diatrc ) & 153 & trc2d(:,:,jp_pcs0_2d + 11) = zirondep(:,:,1) * 1.e+3 * rfact2r * fse3t(:,:,1) * tmask(:,:,1) 161 154 ENDIF 162 155 CALL wrk_dealloc( jpi, jpj, zpdep, zsidep ) … … 168 161 ! ---------------------------------------------------------- 169 162 IF( ln_river ) THEN 170 trn(:,:,1,jppo4) = trn(:,:,1,jppo4) + rivdip(:,:) * rfact2 171 trn(:,:,1,jpno3) = trn(:,:,1,jpno3) + rivdin(:,:) * rfact2 172 trn(:,:,1,jpfer) = trn(:,:,1,jpfer) + rivdic(:,:) * 5.e-5 * rfact2 173 trn(:,:,1,jpsil) = trn(:,:,1,jpsil) + rivdsi(:,:) * rfact2 174 trn(:,:,1,jpdic) = trn(:,:,1,jpdic) + rivdic(:,:) * rfact2 175 trn(:,:,1,jptal) = trn(:,:,1,jptal) + ( rivalk(:,:) - rno3 * rivdin(:,:) ) * rfact2 163 DO jj = 1, jpj 164 DO ji = 1, jpi 165 DO jk = 1, nk_rnf(ji,jj) 166 tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + rivdip(ji,jj) * rfact2 167 tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) + rivdin(ji,jj) * rfact2 168 tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + rivdic(ji,jj) * 5.e-5 * rfact2 169 tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) + rivdsi(ji,jj) * rfact2 170 tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + rivdic(ji,jj) * rfact2 171 tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + ( rivalk(ji,jj) - rno3 * rivdin(ji,jj) ) * rfact2 172 ENDDO 173 ENDDO 174 ENDDO 176 175 ENDIF 177 176 … … 179 178 ! ---------------------------------------------------------- 180 179 IF( ln_ndepo ) THEN 181 tr n(:,:,1,jpno3) = trn(:,:,1,jpno3) + nitdep(:,:) * rfact2182 tr n(:,:,1,jptal) = trn(:,:,1,jptal) - rno3 * nitdep(:,:) * rfact2180 tra(:,:,1,jpno3) = tra(:,:,1,jpno3) + nitdep(:,:) * rfact2 181 tra(:,:,1,jptal) = tra(:,:,1,jptal) - rno3 * nitdep(:,:) * rfact2 183 182 ENDIF 184 183 … … 186 185 ! ------------------------------------------------------ 187 186 IF( ln_ironsed ) THEN 188 tr n(:,:,:,jpfer) = trn(:,:,:,jpfer) + ironsed(:,:,:) * rfact2187 tra(:,:,:,jpfer) = tra(:,:,:,jpfer) + ironsed(:,:,:) * rfact2 189 188 ! 190 IF( l n_diatrc .AND. lk_iomput .AND. jnt == nrdttrc) &189 IF( lk_iomput .AND. knt == nrdttrc .AND. iom_use( "Ironsed" ) ) & 191 190 & CALL iom_put( "Ironsed", ironsed(:,:,:) * 1.e+3 * tmask(:,:,:) ) ! iron inputs from sediments 192 191 ENDIF … … 195 194 ! ------------------------------------------------------ 196 195 IF( ln_hydrofe ) THEN 197 tr n(:,:,:,jpfer) = trn(:,:,:,jpfer) + hydrofe(:,:,:) * rfact2196 tra(:,:,:,jpfer) = tra(:,:,:,jpfer) + hydrofe(:,:,:) * rfact2 198 197 ! 199 IF( l n_diatrc .AND. lk_iomput .AND. jnt == nrdttrc) &198 IF( lk_iomput .AND. knt == nrdttrc .AND. iom_use( "HYDR" ) ) & 200 199 & CALL iom_put( "HYDR", hydrofe(:,:,:) * 1.e+3 * tmask(:,:,:) ) ! hydrothermal iron input 201 200 ENDIF 202 203 201 204 202 ! OA: Warning, the following part is necessary, especially with Kriest … … 224 222 ikt = mbkt(ji,jj) 225 223 # if defined key_kriest 226 zflx = tr n(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) * 1E3 * 1E6 / 1E4224 zflx = trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) * 1E3 * 1E6 / 1E4 227 225 # else 228 zflx = ( tr n(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj) &229 & + tr n(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) ) * 1E3 * 1E6 / 1E4226 zflx = ( trb(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj) & 227 & + trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) ) * 1E3 * 1E6 / 1E4 230 228 #endif 231 229 zflx = LOG10( MAX( 1E-3, zflx ) ) 232 zo2 = LOG10( MAX( 10. , tr n(ji,jj,ikt,jpoxy) * 1E6 ) )233 zno3 = LOG10( MAX( 1. , tr n(ji,jj,ikt,jpno3) * 1E6 * rno3 ) )230 zo2 = LOG10( MAX( 10. , trb(ji,jj,ikt,jpoxy) * 1E6 ) ) 231 zno3 = LOG10( MAX( 1. , trb(ji,jj,ikt,jpno3) * 1E6 * rno3 ) ) 234 232 zdep = LOG10( fsdepw(ji,jj,ikt+1) ) 235 233 zdenit2d(ji,jj) = -2.2567 - 1.185 * zflx - 0.221 * zflx**2 - 0.3995 * zno3 * zo2 + 1.25 * zno3 & … … 237 235 zdenit2d(ji,jj) = 10.0**( zdenit2d(ji,jj) ) 238 236 ! 239 zflx = ( tr n(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj) &240 & + tr n(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) ) * 1E6237 zflx = ( trb(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj) & 238 & + trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) ) * 1E6 241 239 zbureff(ji,jj) = 0.013 + 0.53 * zflx**2 / ( 7.0 + zflx )**2 242 240 ENDIF … … 250 248 DO jj = 1, jpj 251 249 DO ji = 1, jpi 252 ikt = mbkt(ji,jj) 250 IF( tmask(ji,jj,1) == 1 ) THEN 251 ikt = mbkt(ji,jj) 253 252 # if defined key_kriest 254 zwork1(ji,jj) = trn(ji,jj,ikt,jpgsi) * zwscal (ji,jj)255 zwork2(ji,jj) = trn(ji,jj,ikt,jppoc) * zwsbio3(ji,jj)253 zwork1(ji,jj) = trb(ji,jj,ikt,jpgsi) * zwscal (ji,jj) 254 zwork2(ji,jj) = trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) 256 255 # else 257 zwork1(ji,jj) = trn(ji,jj,ikt,jpgsi) * zwsbio4(ji,jj)258 zwork2(ji,jj) = trn(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj) + trn(ji,jj,ikt,jppoc) * zwsbio3(ji,jj)256 zwork1(ji,jj) = trb(ji,jj,ikt,jpgsi) * zwsbio4(ji,jj) 257 zwork2(ji,jj) = trb(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj) + trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) 259 258 # endif 260 ! For calcite, burial efficiency is made a function of saturation 261 zfactcal = MIN( excess(ji,jj,ikt), 0.2 ) 262 zfactcal = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) ) 263 zwork3(ji,jj) = trn(ji,jj,ikt,jpcal) * zwscal(ji,jj) * 2.e0 * zfactcal 259 ! For calcite, burial efficiency is made a function of saturation 260 zfactcal = MIN( excess(ji,jj,ikt), 0.2 ) 261 zfactcal = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) ) 262 zwork3(ji,jj) = trb(ji,jj,ikt,jpcal) * zwscal(ji,jj) * 2.e0 * zfactcal 263 ENDIF 264 264 END DO 265 265 END DO … … 279 279 DO ji = 1, jpi 280 280 ikt = mbkt(ji,jj) 281 zdep = xstep / fse3t(ji,jj,ikt) 281 zdep = xstep / fse3t(ji,jj,ikt) 282 282 zws4 = zwsbio4(ji,jj) * zdep 283 283 zwsc = zwscal (ji,jj) * zdep 284 284 # if defined key_kriest 285 zsiloss = tr n(ji,jj,ikt,jpgsi) * zws4285 zsiloss = trb(ji,jj,ikt,jpgsi) * zws4 286 286 # else 287 zsiloss = tr n(ji,jj,ikt,jpgsi) * zwsc287 zsiloss = trb(ji,jj,ikt,jpgsi) * zwsc 288 288 # endif 289 zcaloss = tr n(ji,jj,ikt,jpcal) * zwsc289 zcaloss = trb(ji,jj,ikt,jpcal) * zwsc 290 290 ! 291 tr n(ji,jj,ikt,jpgsi) = trn(ji,jj,ikt,jpgsi) - zsiloss292 tr n(ji,jj,ikt,jpcal) = trn(ji,jj,ikt,jpcal) - zcaloss291 tra(ji,jj,ikt,jpgsi) = tra(ji,jj,ikt,jpgsi) - zsiloss 292 tra(ji,jj,ikt,jpcal) = tra(ji,jj,ikt,jpcal) - zcaloss 293 293 #if ! defined key_sed 294 tr n(ji,jj,ikt,jpsil) = trn(ji,jj,ikt,jpsil) + zsiloss * zrivsil294 tra(ji,jj,ikt,jpsil) = tra(ji,jj,ikt,jpsil) + zsiloss * zrivsil 295 295 zfactcal = MIN( excess(ji,jj,ikt), 0.2 ) 296 296 zfactcal = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) ) 297 297 zrivalk = 1._wp - ( rivalkinput * r1_ryyss ) * zfactcal / ( zsumsedcal + rtrn ) 298 tr n(ji,jj,ikt,jptal) = trn(ji,jj,ikt,jptal) + zcaloss * zrivalk * 2.0299 tr n(ji,jj,ikt,jpdic) = trn(ji,jj,ikt,jpdic) + zcaloss * zrivalk298 tra(ji,jj,ikt,jptal) = tra(ji,jj,ikt,jptal) + zcaloss * zrivalk * 2.0 299 tra(ji,jj,ikt,jpdic) = tra(ji,jj,ikt,jpdic) + zcaloss * zrivalk 300 300 #endif 301 301 END DO … … 304 304 DO jj = 1, jpj 305 305 DO ji = 1, jpi 306 ikt 307 zdep = xstep / fse3t(ji,jj,ikt)306 ikt = mbkt(ji,jj) 307 zdep = xstep / fse3t(ji,jj,ikt) 308 308 zws4 = zwsbio4(ji,jj) * zdep 309 309 zws3 = zwsbio3(ji,jj) * zdep 310 310 zrivno3 = 1. - zbureff(ji,jj) 311 311 # if ! defined key_kriest 312 tr n(ji,jj,ikt,jpgoc) = trn(ji,jj,ikt,jpgoc) - trn(ji,jj,ikt,jpgoc) * zws4313 tr n(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) - trn(ji,jj,ikt,jppoc) * zws3314 tr n(ji,jj,ikt,jpbfe) = trn(ji,jj,ikt,jpbfe) - trn(ji,jj,ikt,jpbfe) * zws4315 tr n(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) - trn(ji,jj,ikt,jpsfe) * zws3316 zwstpoc = trn(ji,jj,ikt,jpgoc) * zws4 + trn(ji,jj,ikt,jppoc) * zws3312 tra(ji,jj,ikt,jpgoc) = tra(ji,jj,ikt,jpgoc) - trb(ji,jj,ikt,jpgoc) * zws4 313 tra(ji,jj,ikt,jppoc) = tra(ji,jj,ikt,jppoc) - trb(ji,jj,ikt,jppoc) * zws3 314 tra(ji,jj,ikt,jpbfe) = tra(ji,jj,ikt,jpbfe) - trb(ji,jj,ikt,jpbfe) * zws4 315 tra(ji,jj,ikt,jpsfe) = tra(ji,jj,ikt,jpsfe) - trb(ji,jj,ikt,jpsfe) * zws3 316 zwstpoc = trb(ji,jj,ikt,jpgoc) * zws4 + trb(ji,jj,ikt,jppoc) * zws3 317 317 # else 318 tr n(ji,jj,ikt,jpnum) = trn(ji,jj,ikt,jpnum) - trn(ji,jj,ikt,jpnum) * zws4319 tr n(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) - trn(ji,jj,ikt,jppoc) * zws3320 tr n(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) - trn(ji,jj,ikt,jpsfe) * zws3321 zwstpoc = tr n(ji,jj,ikt,jppoc) * zws3318 tra(ji,jj,ikt,jpnum) = tra(ji,jj,ikt,jpnum) - trb(ji,jj,ikt,jpnum) * zws4 319 tra(ji,jj,ikt,jppoc) = tra(ji,jj,ikt,jppoc) - trb(ji,jj,ikt,jppoc) * zws3 320 tra(ji,jj,ikt,jpsfe) = tra(ji,jj,ikt,jpsfe) - trb(ji,jj,ikt,jpsfe) * zws3 321 zwstpoc = trb(ji,jj,ikt,jppoc) * zws3 322 322 # endif 323 323 … … 325 325 ! The 0.5 factor in zpdenit and zdenitt is to avoid negative NO3 concentration after both denitrification 326 326 ! in the sediments and just above the sediments. Not very clever, but simpliest option. 327 zpdenit = MIN( 0.5 * ( tr n(ji,jj,ikt,jpno3) - rtrn ) / rdenit, zdenit2d(ji,jj) * zwstpoc * zrivno3 )327 zpdenit = MIN( 0.5 * ( trb(ji,jj,ikt,jpno3) - rtrn ) / rdenit, zdenit2d(ji,jj) * zwstpoc * zrivno3 ) 328 328 z1pdenit = zwstpoc * zrivno3 - zpdenit 329 zolimit = MIN( ( tr n(ji,jj,ikt,jpoxy) - rtrn ) / o2ut, z1pdenit * ( 1.- nitrfac(ji,jj,ikt) ) )330 zdenitt = MIN( 0.5 * ( tr n(ji,jj,ikt,jpno3) - rtrn ) / rdenit, z1pdenit * nitrfac(ji,jj,ikt) )331 tr n(ji,jj,ikt,jpdoc) = trn(ji,jj,ikt,jpdoc) + z1pdenit - zolimit - zdenitt332 tr n(ji,jj,ikt,jppo4) = trn(ji,jj,ikt,jppo4) + zpdenit + zolimit + zdenitt333 tr n(ji,jj,ikt,jpnh4) = trn(ji,jj,ikt,jpnh4) + zpdenit + zolimit + zdenitt334 tr n(ji,jj,ikt,jpno3) = trn(ji,jj,ikt,jpno3) - rdenit * (zpdenit + zdenitt)335 tr n(ji,jj,ikt,jpoxy) = trn(ji,jj,ikt,jpoxy) - zolimit * o2ut336 tr n(ji,jj,ikt,jptal) = trn(ji,jj,ikt,jptal) + rno3 * (zolimit + (1.+rdenit) * (zpdenit + zdenitt) )337 tr n(ji,jj,ikt,jpdic) = trn(ji,jj,ikt,jpdic) + zpdenit + zolimit + zdenitt338 zwork4(ji,jj) = rdenit * zpdenit * fse3t(ji,jj,ikt)329 zolimit = MIN( ( trb(ji,jj,ikt,jpoxy) - rtrn ) / o2ut, z1pdenit * ( 1.- nitrfac(ji,jj,ikt) ) ) 330 zdenitt = MIN( 0.5 * ( trb(ji,jj,ikt,jpno3) - rtrn ) / rdenit, z1pdenit * nitrfac(ji,jj,ikt) ) 331 tra(ji,jj,ikt,jpdoc) = tra(ji,jj,ikt,jpdoc) + z1pdenit - zolimit - zdenitt 332 tra(ji,jj,ikt,jppo4) = tra(ji,jj,ikt,jppo4) + zpdenit + zolimit + zdenitt 333 tra(ji,jj,ikt,jpnh4) = tra(ji,jj,ikt,jpnh4) + zpdenit + zolimit + zdenitt 334 tra(ji,jj,ikt,jpno3) = tra(ji,jj,ikt,jpno3) - rdenit * (zpdenit + zdenitt) 335 tra(ji,jj,ikt,jpoxy) = tra(ji,jj,ikt,jpoxy) - zolimit * o2ut 336 tra(ji,jj,ikt,jptal) = tra(ji,jj,ikt,jptal) + rno3 * (zolimit + (1.+rdenit) * (zpdenit + zdenitt) ) 337 tra(ji,jj,ikt,jpdic) = tra(ji,jj,ikt,jpdic) + zpdenit + zolimit + zdenitt 338 sdenit(ji,jj) = rdenit * zpdenit * fse3t(ji,jj,ikt) 339 339 #endif 340 340 END DO … … 356 356 #endif 357 357 ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) ) 358 ztrpo4 = tr n (ji,jj,jk,jppo4) / ( concnnh4 + trn(ji,jj,jk,jppo4) )359 zlight = ( 1.- EXP( -etot (ji,jj,jk) / diazolight ) )360 znitrpot(ji,jj,jk) = MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) * r1_rday ) &358 ztrpo4 = trb (ji,jj,jk,jppo4) / ( concnnh4 + trb (ji,jj,jk,jppo4) ) 359 zlight = ( 1.- EXP( -etot_ndcy(ji,jj,jk) / diazolight ) ) 360 nitrpot(ji,jj,jk) = MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) * r1_rday ) & 361 361 & * zfact * MIN( ztrfer, ztrpo4 ) * zlight 362 362 zsoufer(ji,jj,jk) = zlight * 2E-11 / (2E-11 + biron(ji,jj,jk)) … … 370 370 DO jj = 1, jpj 371 371 DO ji = 1, jpi 372 zfact = znitrpot(ji,jj,jk) * nitrfix373 tr n(ji,jj,jk,jpnh4) = trn(ji,jj,jk,jpnh4) + zfact374 tr n(ji,jj,jk,jptal) = trn(ji,jj,jk,jptal) + rno3 * zfact375 tr n(ji,jj,jk,jpoxy) = trn(ji,jj,jk,jpoxy) + o2nit * zfact376 tr n(ji,jj,jk,jppo4) = trn(ji,jj,jk,jppo4) + concdnh4 / ( concdnh4 + trn(ji,jj,jk,jppo4) ) &377 & * 0.002 * tr n(ji,jj,jk,jpdoc) * rfact2 / rday378 tr n(ji,jj,jk,jpfer) = trn(ji,jj,jk,jpfer) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * rfact2 / rday372 zfact = nitrpot(ji,jj,jk) * nitrfix 373 tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zfact 374 tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * zfact 375 tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + o2nit * zfact 376 tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + concdnh4 / ( concdnh4 + trb(ji,jj,jk,jppo4) ) & 377 & * 0.002 * trb(ji,jj,jk,jpdoc) * xstep 378 tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * xstep 379 379 END DO 380 380 END DO 381 381 END DO 382 382 383 384 IF( ln_check_mass) THEN385 ! Global budget of N SMS : denitrification in the water column and in the sediment386 ! nitrogen fixation by the diazotrophs387 ! --------------------------------------------------------------------------------388 zrdenittot = glob_sum ( denitr(:,:,:) * rdenit * xnegtr(:,:,:) * cvol(:,:,:) )389 zsdenittot = glob_sum ( zwork4(:,:) * e1e2t(:,:) )390 znitrpottot = glob_sum ( znitrpot(:,:,:) * nitrfix * cvol(:,:,:))391 IF( kt == nitend .AND. jnt == nrdttrc ) THEN392 zfact = 1.e+3 * rfact2r * rno3 * ryyss * 14. / 1e12393 IF(lwp) WRITE(numnit,9100) ndastp, znitrpottot * nitrfix * zfact, zrdenittot * zfact , zsdenittot * zfact383 IF( lk_iomput ) THEN 384 IF( knt == nrdttrc ) THEN 385 zfact = 1.e+3 * rfact2r * rno3 ! conversion from molC/l/kt to molN/m3/s 386 IF( iom_use("Nfix" ) ) CALL iom_put( "Nfix", nitrpot(:,:,:) * nitrfix * zfact * tmask(:,:,:) ) ! nitrogen fixation 387 IF( iom_use("INTNFIX") ) THEN ! nitrogen fixation rate in ocean ( vertically integrated ) 388 zwork1(:,:) = 0. 389 DO jk = 1, jpkm1 390 zwork1(:,:) = zwork1(:,:) + nitrpot(:,:,jk) * nitrfix * zfact * fse3t(:,:,jk) * tmask(:,:,jk) 391 ENDDO 392 CALL iom_put( "INTNFIX" , zwork1 ) 393 ENDIF 394 394 ENDIF 395 ENDIF 396 ! 397 IF( ln_diatrc ) THEN 398 zfact = 1.e+3 * rfact2r 399 IF( lk_iomput ) THEN 400 IF( jnt == nrdttrc ) THEN 401 CALL iom_put( "Nfix" , znitrpot(:,:,:) * nitrfix * rno3 * zfact * tmask(:,:,:) ) ! nitrogen fixation 402 CALL iom_put( "Sdenit", zwork4(:,:) * rno3 * zfact * tmask(:,:,1) ) ! Nitrate reduction in the sediments 403 ENDIF 404 ELSE 405 trc2d(:,:,jp_pcs0_2d + 12) = znitrpot(:,:,1) * nitrfix * zfact * fse3t(:,:,1) * tmask(:,:,1) 406 ENDIF 395 ELSE 396 IF( ln_diatrc ) & 397 & trc2d(:,:,jp_pcs0_2d + 12) = nitrpot(:,:,1) * nitrfix * rno3 * 1.e+3 * rfact2r * fse3t(:,:,1) * tmask(:,:,1) 407 398 ENDIF 408 399 ! … … 410 401 WRITE(charout, fmt="('sed ')") 411 402 CALL prt_ctl_trc_info(charout) 412 CALL prt_ctl_trc(tab4d=tr n, mask=tmask, clinfo=ctrcnm)413 ENDIF 414 ! 415 CALL wrk_dealloc( jpi, jpj, zdenit2d, zwork1, zwork2, zwork3, z work4, zbureff )403 CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm) 404 ENDIF 405 ! 406 CALL wrk_dealloc( jpi, jpj, zdenit2d, zwork1, zwork2, zwork3, zbureff ) 416 407 CALL wrk_dealloc( jpi, jpj, zwsbio3, zwsbio4, zwscal ) 417 CALL wrk_dealloc( jpi, jpj, jpk, z nitrpot, zsoufer )408 CALL wrk_dealloc( jpi, jpj, jpk, zsoufer ) 418 409 ! 419 410 IF( nn_timing == 1 ) CALL timing_stop('p4z_sed') … … 422 413 ! 423 414 END SUBROUTINE p4z_sed 415 416 417 INTEGER FUNCTION p4z_sed_alloc() 418 !!---------------------------------------------------------------------- 419 !! *** ROUTINE p4z_sed_alloc *** 420 !!---------------------------------------------------------------------- 421 ALLOCATE( nitrpot(jpi,jpj,jpk), sdenit(jpi,jpj), STAT=p4z_sed_alloc ) 422 ! 423 IF( p4z_sed_alloc /= 0 ) CALL ctl_warn('p4z_sed_alloc: failed to allocate arrays') 424 ! 425 END FUNCTION p4z_sed_alloc 426 424 427 425 428 #else -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsink.F90
r4793 r5581 41 41 #endif 42 42 43 INTEGER :: ik sed = 1043 INTEGER :: ik100 44 44 45 45 #if defined key_kriest … … 79 79 !!---------------------------------------------------------------------- 80 80 81 SUBROUTINE p4z_sink ( kt, jnt )81 SUBROUTINE p4z_sink ( kt, knt ) 82 82 !!--------------------------------------------------------------------- 83 83 !! *** ROUTINE p4z_sink *** … … 88 88 !! ** Method : - ??? 89 89 !!--------------------------------------------------------------------- 90 INTEGER, INTENT(in) :: kt, jnt90 INTEGER, INTENT(in) :: kt, knt 91 91 INTEGER :: ji, jj, jk, jit 92 92 INTEGER :: iiter1, iiter2 … … 94 94 REAL(wp) :: zagg , zaggfe, zaggdoc, zaggdoc2, zaggdoc3 95 95 REAL(wp) :: zfact, zwsmax, zmax, zstep 96 REAL(wp) :: zrfact297 INTEGER :: ik198 96 CHARACTER (len=25) :: charout 97 REAL(wp), POINTER, DIMENSION(:,:,:) :: zw3d 98 REAL(wp), POINTER, DIMENSION(:,: ) :: zw2d 99 99 !!--------------------------------------------------------------------- 100 100 ! … … 199 199 zfact = zstep * xdiss(ji,jj,jk) 200 200 ! Part I : Coagulation dependent on turbulence 201 zagg1 = 25.9 * zfact * tr n(ji,jj,jk,jppoc) * trn(ji,jj,jk,jppoc)202 zagg2 = 4452. * zfact * tr n(ji,jj,jk,jppoc) * trn(ji,jj,jk,jpgoc)201 zagg1 = 25.9 * zfact * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jppoc) 202 zagg2 = 4452. * zfact * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpgoc) 203 203 204 204 ! Part II : Differential settling 205 205 206 206 ! Aggregation of small into large particles 207 zagg3 = 47.1 * zstep * tr n(ji,jj,jk,jppoc) * trn(ji,jj,jk,jpgoc)208 zagg4 = 3.3 * zstep * tr n(ji,jj,jk,jppoc) * trn(ji,jj,jk,jppoc)207 zagg3 = 47.1 * zstep * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpgoc) 208 zagg4 = 3.3 * zstep * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jppoc) 209 209 210 210 zagg = zagg1 + zagg2 + zagg3 + zagg4 211 zaggfe = zagg * tr n(ji,jj,jk,jpsfe) / ( trn(ji,jj,jk,jppoc) + rtrn )211 zaggfe = zagg * trb(ji,jj,jk,jpsfe) / ( trb(ji,jj,jk,jppoc) + rtrn ) 212 212 213 213 ! Aggregation of DOC to POC : … … 215 215 ! 2nd term is shear aggregation of DOC-POC 216 216 ! 3rd term is differential settling of DOC-POC 217 zaggdoc = ( ( 0.369 * 0.3 * tr n(ji,jj,jk,jpdoc) + 102.4 * trn(ji,jj,jk,jppoc) ) * zfact &218 & + 2.4 * zstep * tr n(ji,jj,jk,jppoc) ) * 0.3 * trn(ji,jj,jk,jpdoc)217 zaggdoc = ( ( 0.369 * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4 * trb(ji,jj,jk,jppoc) ) * zfact & 218 & + 2.4 * zstep * trb(ji,jj,jk,jppoc) ) * 0.3 * trb(ji,jj,jk,jpdoc) 219 219 ! transfer of DOC to GOC : 220 220 ! 1st term is shear aggregation 221 221 ! 2nd term is differential settling 222 zaggdoc2 = ( 3.53E3 * zfact + 0.1 * zstep ) * tr n(ji,jj,jk,jpgoc) * 0.3 * trn(ji,jj,jk,jpdoc)222 zaggdoc2 = ( 3.53E3 * zfact + 0.1 * zstep ) * trb(ji,jj,jk,jpgoc) * 0.3 * trb(ji,jj,jk,jpdoc) 223 223 ! tranfer of DOC to POC due to brownian motion 224 zaggdoc3 = ( 5095. * tr n(ji,jj,jk,jppoc) + 114. * 0.3 * trn(ji,jj,jk,jpdoc) ) *zstep * 0.3 * trn(ji,jj,jk,jpdoc)224 zaggdoc3 = ( 5095. * trb(ji,jj,jk,jppoc) + 114. * 0.3 * trb(ji,jj,jk,jpdoc) ) *zstep * 0.3 * trb(ji,jj,jk,jpdoc) 225 225 226 226 ! Update the trends … … 235 235 END DO 236 236 237 ! Total primary production per year 238 t_oce_co2_exp = t_oce_co2_exp + glob_sum( ( sinking(:,:,iksed+1) + sinking2(:,:,iksed+1) ) * e1e2t(:,:) * tmask(:,:,1) ) 237 238 ! Total carbon export per year 239 IF( iom_use( "tcexp" ) .OR. ( ln_check_mass .AND. kt == nitend .AND. knt == nrdttrc ) ) & 240 & t_oce_co2_exp = glob_sum( ( sinking(:,:,ik100) + sinking2(:,:,ik100) ) * e1e2t(:,:) * tmask(:,:,1) ) 239 241 ! 240 IF( ln_diatrc ) THEN 241 zrfact2 = 1.e3 * rfact2r 242 ik1 = iksed + 1 243 IF( lk_iomput ) THEN 244 IF( jnt == nrdttrc ) THEN 245 CALL iom_put( "EPC100" , ( sinking(:,:,ik1) + sinking2(:,:,ik1) ) * zrfact2 * tmask(:,:,1) ) ! Export of carbon at 100m 246 CALL iom_put( "EPFE100" , ( sinkfer(:,:,ik1) + sinkfer2(:,:,ik1) ) * zrfact2 * tmask(:,:,1) ) ! Export of iron at 100m 247 CALL iom_put( "EPCAL100", sinkcal(:,:,ik1) * zrfact2 * tmask(:,:,1) ) ! Export of calcite at 100m 248 CALL iom_put( "EPSI100" , sinksil(:,:,ik1) * zrfact2 * tmask(:,:,1) ) ! Export of biogenic silica at 100m 249 ENDIF 250 ELSE 251 trc2d(:,:,jp_pcs0_2d + 4) = sinking (:,:,ik1) * zrfact2 * tmask(:,:,1) 252 trc2d(:,:,jp_pcs0_2d + 5) = sinking2(:,:,ik1) * zrfact2 * tmask(:,:,1) 253 trc2d(:,:,jp_pcs0_2d + 6) = sinkfer (:,:,ik1) * zrfact2 * tmask(:,:,1) 254 trc2d(:,:,jp_pcs0_2d + 7) = sinkfer2(:,:,ik1) * zrfact2 * tmask(:,:,1) 255 trc2d(:,:,jp_pcs0_2d + 8) = sinksil (:,:,ik1) * zrfact2 * tmask(:,:,1) 256 trc2d(:,:,jp_pcs0_2d + 9) = sinkcal (:,:,ik1) * zrfact2 * tmask(:,:,1) 242 IF( lk_iomput ) THEN 243 IF( knt == nrdttrc ) THEN 244 CALL wrk_alloc( jpi, jpj, zw2d ) 245 CALL wrk_alloc( jpi, jpj, jpk, zw3d ) 246 zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s 247 ! 248 IF( iom_use( "EPC100" ) ) THEN 249 zw2d(:,:) = ( sinking(:,:,ik100) + sinking2(:,:,ik100) ) * zfact * tmask(:,:,1) ! Export of carbon at 100m 250 CALL iom_put( "EPC100" , zw2d ) 251 ENDIF 252 IF( iom_use( "EPFE100" ) ) THEN 253 zw2d(:,:) = ( sinkfer(:,:,ik100) + sinkfer2(:,:,ik100) ) * zfact * tmask(:,:,1) ! Export of iron at 100m 254 CALL iom_put( "EPFE100" , zw2d ) 255 ENDIF 256 IF( iom_use( "EPCAL100" ) ) THEN 257 zw2d(:,:) = sinkcal(:,:,ik100) * zfact * tmask(:,:,1) ! Export of calcite at 100m 258 CALL iom_put( "EPCAL100" , zw2d ) 259 ENDIF 260 IF( iom_use( "EPSI100" ) ) THEN 261 zw2d(:,:) = sinksil(:,:,ik100) * zfact * tmask(:,:,1) ! Export of bigenic silica at 100m 262 CALL iom_put( "EPSI100" , zw2d ) 263 ENDIF 264 IF( iom_use( "EXPC" ) ) THEN 265 zw3d(:,:,:) = ( sinking(:,:,:) + sinking2(:,:,:) ) * zfact * tmask(:,:,:) ! Export of carbon in the water column 266 CALL iom_put( "EXPC" , zw3d ) 267 ENDIF 268 IF( iom_use( "EXPFE" ) ) THEN 269 zw3d(:,:,:) = ( sinkfer(:,:,:) + sinkfer2(:,:,:) ) * zfact * tmask(:,:,:) ! Export of iron 270 CALL iom_put( "EXPFE" , zw3d ) 271 ENDIF 272 IF( iom_use( "EXPCAL" ) ) THEN 273 zw3d(:,:,:) = sinkcal(:,:,:) * zfact * tmask(:,:,:) ! Export of calcite 274 CALL iom_put( "EXPCAL" , zw3d ) 275 ENDIF 276 IF( iom_use( "EXPSI" ) ) THEN 277 zw3d(:,:,:) = sinksil(:,:,:) * zfact * tmask(:,:,:) ! Export of bigenic silica 278 CALL iom_put( "EXPSI" , zw3d ) 279 ENDIF 280 IF( iom_use( "tcexp" ) ) CALL iom_put( "tcexp" , t_oce_co2_exp * zfact ) ! molC/s 281 ! 282 CALL wrk_dealloc( jpi, jpj, zw2d ) 283 CALL wrk_dealloc( jpi, jpj, jpk, zw3d ) 284 ENDIF 285 ELSE 286 IF( ln_diatrc ) THEN 287 zfact = 1.e3 * rfact2r 288 trc2d(:,:,jp_pcs0_2d + 4) = sinking (:,:,ik100) * zfact * tmask(:,:,1) 289 trc2d(:,:,jp_pcs0_2d + 5) = sinking2(:,:,ik100) * zfact * tmask(:,:,1) 290 trc2d(:,:,jp_pcs0_2d + 6) = sinkfer (:,:,ik100) * zfact * tmask(:,:,1) 291 trc2d(:,:,jp_pcs0_2d + 7) = sinkfer2(:,:,ik100) * zfact * tmask(:,:,1) 292 trc2d(:,:,jp_pcs0_2d + 8) = sinksil (:,:,ik100) * zfact * tmask(:,:,1) 293 trc2d(:,:,jp_pcs0_2d + 9) = sinkcal (:,:,ik100) * zfact * tmask(:,:,1) 257 294 ENDIF 258 295 ENDIF … … 272 309 !! *** ROUTINE p4z_sink_init *** 273 310 !!---------------------------------------------------------------------- 274 311 INTEGER :: jk 312 313 ik100 = 10 ! last level where depth less than 100 m 314 DO jk = jpkm1, 1, -1 315 IF( gdept_1d(jk) > 100. ) ik100 = jk - 1 316 END DO 317 IF (lwp) WRITE(numout,*) 318 IF (lwp) WRITE(numout,*) ' Level corresponding to 100m depth ', ik100 + 1 319 IF (lwp) WRITE(numout,*) 320 ! 275 321 t_oce_co2_exp = 0._wp 276 322 ! … … 282 328 !!---------------------------------------------------------------------- 283 329 284 SUBROUTINE p4z_sink ( kt, jnt )330 SUBROUTINE p4z_sink ( kt, knt ) 285 331 !!--------------------------------------------------------------------- 286 332 !! *** ROUTINE p4z_sink *** … … 292 338 !!--------------------------------------------------------------------- 293 339 ! 294 INTEGER, INTENT(in) :: kt, jnt340 INTEGER, INTENT(in) :: kt, knt 295 341 ! 296 342 INTEGER :: ji, jj, jk, jit, niter1, niter2 … … 300 346 REAL(wp) :: zdiv , zdiv1, zdiv2, zdiv3, zdiv4, zdiv5 301 347 REAL(wp) :: zval1, zval2, zval3, zval4 302 REAL(wp) :: z rfact2348 REAL(wp) :: zfact 303 349 INTEGER :: ik1 304 350 CHARACTER (len=25) :: charout 305 351 REAL(wp), POINTER, DIMENSION(:,:,:) :: znum3d 352 REAL(wp), POINTER, DIMENSION(:,:,:) :: zw3d 353 REAL(wp), POINTER, DIMENSION(:,: ) :: zw2d 306 354 !!--------------------------------------------------------------------- 307 355 ! … … 325 373 DO ji = 1, jpi 326 374 IF( tmask(ji,jj,jk) /= 0.e0 ) THEN 327 znum = tr n(ji,jj,jk,jppoc) / ( trn(ji,jj,jk,jpnum) + rtrn ) / xkr_massp375 znum = trb(ji,jj,jk,jppoc) / ( trb(ji,jj,jk,jpnum) + rtrn ) / xkr_massp 328 376 ! -------------- To avoid sinking speed over 50 m/day ------- 329 377 znum = MIN( xnumm(jk), znum ) … … 387 435 IF( tmask(ji,jj,jk) /= 0.e0 ) THEN 388 436 389 znum = tr n(ji,jj,jk,jppoc)/(trn(ji,jj,jk,jpnum)+rtrn) / xkr_massp437 znum = trb(ji,jj,jk,jppoc)/(trb(ji,jj,jk,jpnum)+rtrn) / xkr_massp 390 438 !-------------- To avoid sinking speed over 50 m/day ------- 391 439 znum = min(xnumm(jk),znum) … … 405 453 ! ---------------------------------------------- 406 454 407 zagg1 = 0.163 * tr n(ji,jj,jk,jpnum)**2 &455 zagg1 = 0.163 * trb(ji,jj,jk,jpnum)**2 & 408 456 & * 2.*( (zfm-1.)*(zfm*xkr_mass_max**3-xkr_mass_min**3) & 409 457 & * (zeps-1)/zdiv1 + 3.*(zfm*xkr_mass_max-xkr_mass_min) & 410 458 & * (zfm*xkr_mass_max**2-xkr_mass_min**2) & 411 459 & * (zeps-1.)**2/(zdiv2*zdiv3)) 412 zagg2 = 2*0.163*tr n(ji,jj,jk,jpnum)**2*zfm* &460 zagg2 = 2*0.163*trb(ji,jj,jk,jpnum)**2*zfm* & 413 461 & ((xkr_mass_max**3+3.*(xkr_mass_max**2 & 414 462 & *xkr_mass_min*(zeps-1.)/zdiv2 & … … 418 466 & (zeps-2.)+(zeps-1.)/zdiv3)+(zeps-1.)/zdiv1)) 419 467 420 zagg3 = 0.163*tr n(ji,jj,jk,jpnum)**2*zfm**2*8. * xkr_mass_max**3468 zagg3 = 0.163*trb(ji,jj,jk,jpnum)**2*zfm**2*8. * xkr_mass_max**3 421 469 422 470 ! Aggregation of small into large particles … … 424 472 ! ---------------------------------------------- 425 473 426 zagg4 = 2.*3.141*0.125*tr n(ji,jj,jk,jpnum)**2* &474 zagg4 = 2.*3.141*0.125*trb(ji,jj,jk,jpnum)**2* & 427 475 & xkr_wsbio_min*(zeps-1.)**2 & 428 476 & *(xkr_mass_min**2*((1.-zsm*zfm)/(zdiv3*zdiv4) & … … 431 479 & *xkr_eta)/(zdiv*zdiv3*zdiv5) ) 432 480 433 zagg5 = 2.*3.141*0.125*tr n(ji,jj,jk,jpnum)**2 &481 zagg5 = 2.*3.141*0.125*trb(ji,jj,jk,jpnum)**2 & 434 482 & *(zeps-1.)*zfm*xkr_wsbio_min & 435 483 & *(zsm*(xkr_mass_min**2-zfm*xkr_mass_max**2) & … … 441 489 ! ------------------------------------ 442 490 443 zfract = 2.*3.141*0.125*tr n(ji,jj,jk,jpmes)*12./0.12/0.06**3*trn(ji,jj,jk,jpnum) &491 zfract = 2.*3.141*0.125*trb(ji,jj,jk,jpmes)*12./0.12/0.06**3*trb(ji,jj,jk,jpnum) & 444 492 & * (0.01/xkr_mass_min)**(1.-zeps)*0.1**2 & 445 493 & * 10000.*xstep … … 448 496 ! -------------------------------------- 449 497 450 zaggdoc = 0.83 * tr n(ji,jj,jk,jpdoc) * xstep * xdiss(ji,jj,jk) * trn(ji,jj,jk,jpdoc) &451 & + 0.005 * 231. * tr n(ji,jj,jk,jpdoc) * xstep * trn(ji,jj,jk,jpdoc)452 zaggdoc1 = 271. * tr n(ji,jj,jk,jppoc) * xstep * xdiss(ji,jj,jk) * trn(ji,jj,jk,jpdoc) &453 & + 0.02 * 16706. * tr n(ji,jj,jk,jppoc) * xstep * trn(ji,jj,jk,jpdoc)498 zaggdoc = 0.83 * trb(ji,jj,jk,jpdoc) * xstep * xdiss(ji,jj,jk) * trb(ji,jj,jk,jpdoc) & 499 & + 0.005 * 231. * trb(ji,jj,jk,jpdoc) * xstep * trb(ji,jj,jk,jpdoc) 500 zaggdoc1 = 271. * trb(ji,jj,jk,jppoc) * xstep * xdiss(ji,jj,jk) * trb(ji,jj,jk,jpdoc) & 501 & + 0.02 * 16706. * trb(ji,jj,jk,jppoc) * xstep * trb(ji,jj,jk,jpdoc) 454 502 455 503 # if defined key_degrad … … 466 514 zagg = 0.5 * xkr_stick * ( zaggsh + zaggsi ) 467 515 ! 468 znumdoc = tr n(ji,jj,jk,jpnum) / ( trn(ji,jj,jk,jppoc) + rtrn )516 znumdoc = trb(ji,jj,jk,jpnum) / ( trb(ji,jj,jk,jppoc) + rtrn ) 469 517 tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zaggdoc + zaggdoc1 470 518 tra(ji,jj,jk,jpnum) = tra(ji,jj,jk,jpnum) + zfract + zaggdoc / xkr_massp - zagg … … 477 525 478 526 ! Total primary production per year 479 t_oce_co2_exp = t_oce_co2_exp + glob_sum( ( sinking(:,:, :) ) * cvol(:,:,:) )527 t_oce_co2_exp = t_oce_co2_exp + glob_sum( ( sinking(:,:,ik100) * e1e2t(:,:) * tmask(:,:,1) ) 480 528 ! 481 IF( ln_diatrc ) THEN 482 ! 483 ik1 = iksed + 1 484 zrfact2 = 1.e3 * rfact2r 485 IF( jnt == nrdttrc ) THEN 486 CALL iom_put( "POCFlx" , sinking (:,:,:) * zrfact2 * tmask(:,:,:) ) ! POC export 487 CALL iom_put( "NumFlx" , sinking2 (:,:,:) * zrfact2 * tmask(:,:,:) ) ! Num export 488 CALL iom_put( "SiFlx" , sinksil (:,:,:) * zrfact2 * tmask(:,:,:) ) ! Silica export 489 CALL iom_put( "CaCO3Flx", sinkcal (:,:,:) * zrfact2 * tmask(:,:,:) ) ! Calcite export 490 CALL iom_put( "xnum" , znum3d (:,:,:) * tmask(:,:,:) ) ! Number of particles in aggregats 491 CALL iom_put( "W1" , wsbio3 (:,:,:) * tmask(:,:,:) ) ! sinking speed of POC 492 CALL iom_put( "W2" , wsbio4 (:,:,:) * tmask(:,:,:) ) ! sinking speed of aggregats 529 IF( lk_iomput ) THEN 530 IF( knt == nrdttrc ) THEN 531 CALL wrk_alloc( jpi, jpj, zw2d ) 532 CALL wrk_alloc( jpi, jpj, jpk, zw3d ) 533 zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s 534 ! 535 IF( iom_use( "EPC100" ) ) THEN 536 zw2d(:,:) = sinking(:,:,ik100) * zfact * tmask(:,:,1) ! Export of carbon at 100m 537 CALL iom_put( "EPC100" , zw2d ) 538 ENDIF 539 IF( iom_use( "EPN100" ) ) THEN 540 zw2d(:,:) = sinking2(:,:,ik100) * zfact * tmask(:,:,1) ! Export of number of aggregates ? 541 CALL iom_put( "EPN100" , zw2d ) 542 ENDIF 543 IF( iom_use( "EPCAL100" ) ) THEN 544 zw2d(:,:) = sinkcal(:,:,ik100) * zfact * tmask(:,:,1) ! Export of calcite at 100m 545 CALL iom_put( "EPCAL100" , zw2d ) 546 ENDIF 547 IF( iom_use( "EPSI100" ) ) THEN 548 zw2d(:,:) = sinksil(:,:,ik100) * zfact * tmask(:,:,1) ! Export of bigenic silica at 100m 549 CALL iom_put( "EPSI100" , zw2d ) 550 ENDIF 551 IF( iom_use( "EXPC" ) ) THEN 552 zw3d(:,:,:) = sinking(:,:,:) * zfact * tmask(:,:,:) ! Export of carbon in the water column 553 CALL iom_put( "EXPC" , zw3d ) 554 ENDIF 555 IF( iom_use( "EXPN" ) ) THEN 556 zw3d(:,:,:) = sinking(:,:,:) * zfact * tmask(:,:,:) ! Export of carbon in the water column 557 CALL iom_put( "EXPN" , zw3d ) 558 ENDIF 559 IF( iom_use( "EXPCAL" ) ) THEN 560 zw3d(:,:,:) = sinkcal(:,:,:) * zfact * tmask(:,:,:) ! Export of calcite 561 CALL iom_put( "EXPCAL" , zw3d ) 562 ENDIF 563 IF( iom_use( "EXPSI" ) ) THEN 564 zw3d(:,:,:) = sinksil(:,:,:) * zfact * tmask(:,:,:) ! Export of bigenic silica 565 CALL iom_put( "EXPSI" , zw3d ) 566 ENDIF 567 IF( iom_use( "XNUM" ) ) THEN 568 zw3d(:,:,:) = znum3d(:,:,:) * tmask(:,:,:) ! Number of particles on aggregats 569 CALL iom_put( "XNUM" , zw3d ) 570 ENDIF 571 IF( iom_use( "WSC" ) ) THEN 572 zw3d(:,:,:) = wsbio3(:,:,:) * tmask(:,:,:) ! Sinking speed of carbon particles 573 CALL iom_put( "WSC" , zw3d ) 574 ENDIF 575 IF( iom_use( "WSN" ) ) THEN 576 zw3d(:,:,:) = wsbio4(:,:,:) * tmask(:,:,:) ! Sinking speed of particles number 577 CALL iom_put( "WSN" , zw3d ) 578 ENDIF 579 ! 580 CALL wrk_dealloc( jpi, jpj, zw2d ) 581 CALL wrk_dealloc( jpi, jpj, jpk, zw3d ) 582 ELSE 583 IF( ln_diatrc ) THEN 584 zfact = 1.e3 * rfact2r 585 trc2d(:,: ,jp_pcs0_2d + 4) = sinking (:,:,ik100) * zfact * tmask(:,:,1) 586 trc2d(:,: ,jp_pcs0_2d + 5) = sinking2(:,:,ik100) * zfact * tmask(:,:,1) 587 trc2d(:,: ,jp_pcs0_2d + 6) = sinkfer (:,:,ik100) * zfact * tmask(:,:,1) 588 trc2d(:,: ,jp_pcs0_2d + 7) = sinksil (:,:,ik100) * zfact * tmask(:,:,1) 589 trc2d(:,: ,jp_pcs0_2d + 8) = sinkcal (:,:,ik100) * zfact * tmask(:,:,1) 590 trc3d(:,:,:,jp_pcs0_3d + 11) = sinking (:,:,:) * zfact * tmask(:,:,:) 591 trc3d(:,:,:,jp_pcs0_3d + 12) = sinking2(:,:,:) * zfact * tmask(:,:,:) 592 trc3d(:,:,:,jp_pcs0_3d + 13) = sinksil (:,:,:) * zfact * tmask(:,:,:) 593 trc3d(:,:,:,jp_pcs0_3d + 14) = sinkcal (:,:,:) * zfact * tmask(:,:,:) 594 trc3d(:,:,:,jp_pcs0_3d + 15) = znum3d (:,:,:) * tmask(:,:,:) 595 trc3d(:,:,:,jp_pcs0_3d + 16) = wsbio3 (:,:,:) * tmask(:,:,:) 596 trc3d(:,:,:,jp_pcs0_3d + 17) = wsbio4 (:,:,:) * tmask(:,:,:) 493 597 ENDIF 494 # if ! defined key_iomput495 trc2d(:,: ,jp_pcs0_2d + 4) = sinking (:,:,ik1) * zrfact2 * tmask(:,:,1)496 trc2d(:,: ,jp_pcs0_2d + 5) = sinking2(:,:,ik1) * zrfact2 * tmask(:,:,1)497 trc2d(:,: ,jp_pcs0_2d + 6) = sinkfer (:,:,ik1) * zrfact2 * tmask(:,:,1)498 trc2d(:,: ,jp_pcs0_2d + 7) = sinksil (:,:,ik1) * zrfact2 * tmask(:,:,1)499 trc2d(:,: ,jp_pcs0_2d + 8) = sinkcal (:,:,ik1) * zrfact2 * tmask(:,:,1)500 trc3d(:,:,:,jp_pcs0_3d + 11) = sinking (:,:,:) * zrfact2 * tmask(:,:,:)501 trc3d(:,:,:,jp_pcs0_3d + 12) = sinking2(:,:,:) * zrfact2 * tmask(:,:,:)502 trc3d(:,:,:,jp_pcs0_3d + 13) = sinksil (:,:,:) * zrfact2 * tmask(:,:,:)503 trc3d(:,:,:,jp_pcs0_3d + 14) = sinkcal (:,:,:) * zrfact2 * tmask(:,:,:)504 trc3d(:,:,:,jp_pcs0_3d + 15) = znum3d (:,:,:) * tmask(:,:,:)505 trc3d(:,:,:,jp_pcs0_3d + 16) = wsbio3 (:,:,:) * tmask(:,:,:)506 trc3d(:,:,:,jp_pcs0_3d + 17) = wsbio4 (:,:,:) * tmask(:,:,:)507 # endif508 !509 598 ENDIF 599 510 600 ! 511 601 IF(ln_ctl) THEN ! print mean trends (used for debugging) … … 663 753 END DO 664 754 ! 755 ik100 = 10 ! last level where depth less than 100 m 756 DO jk = jpkm1, 1, -1 757 IF( gdept_1d(jk) > 100. ) iksed = jk - 1 758 END DO 759 IF (lwp) WRITE(numout,*) 760 IF (lwp) WRITE(numout,*) ' Level corresponding to 100m depth ', ik100 + 1 761 IF (lwp) WRITE(numout,*) 762 ! 665 763 t_oce_co2_exp = 0._wp 666 764 ! … … 702 800 ztraz(:,:,:) = 0.e0 703 801 zakz (:,:,:) = 0.e0 704 ztrb (:,:,:) = tr n(:,:,:,jp_tra)802 ztrb (:,:,:) = trb(:,:,:,jp_tra) 705 803 706 804 DO jk = 1, jpkm1 … … 717 815 ! first guess of the slopes interior values 718 816 DO jk = 2, jpkm1 719 ztraz(:,:,jk) = ( tr n(:,:,jk-1,jp_tra) - trn(:,:,jk,jp_tra) ) * tmask(:,:,jk)817 ztraz(:,:,jk) = ( trb(:,:,jk-1,jp_tra) - trb(:,:,jk,jp_tra) ) * tmask(:,:,jk) 720 818 END DO 721 819 ztraz(:,:,1 ) = 0.0 … … 748 846 zigma = zwsink2(ji,jj,jk+1) * zstep / fse3w(ji,jj,jk+1) 749 847 zew = zwsink2(ji,jj,jk+1) 750 psinkflx(ji,jj,jk+1) = -zew * ( tr n(ji,jj,jk,jp_tra) - 0.5 * ( 1 + zigma ) * zakz(ji,jj,jk) ) * zstep848 psinkflx(ji,jj,jk+1) = -zew * ( trb(ji,jj,jk,jp_tra) - 0.5 * ( 1 + zigma ) * zakz(ji,jj,jk) ) * zstep 751 849 END DO 752 850 END DO … … 761 859 DO ji = 1, jpi 762 860 zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / fse3t(ji,jj,jk) 763 tr n(ji,jj,jk,jp_tra) = trn(ji,jj,jk,jp_tra) + zflx861 trb(ji,jj,jk,jp_tra) = trb(ji,jj,jk,jp_tra) + zflx 764 862 END DO 765 863 END DO … … 777 875 END DO 778 876 779 tr n(:,:,:,jp_tra) = ztrb(:,:,:)877 trb(:,:,:,jp_tra) = ztrb(:,:,:) 780 878 psinkflx(:,:,:) = 2. * psinkflx(:,:,:) 781 879 ! -
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90
r4624 r5581 11 11 !! 'key_pisces' PISCES bio-model 12 12 !!---------------------------------------------------------------------- 13 !! p4zsms : Time loop of passive tracers sms13 !! p4zsms : Time loop of passive tracers sms 14 14 !!---------------------------------------------------------------------- 15 15 USE oce_trc ! shared variables between ocean and passive tracers … … 24 24 USE p4zsed ! Sedimentation 25 25 USE p4zint ! time interpolation 26 USE p4zrem ! remineralisation 26 27 USE iom ! I/O manager 27 USE trd mod_oce! Ocean trends variables28 USE trd mod_trc! TOP trends variables28 USE trd_oce ! Ocean trends variables 29 USE trdtrc ! TOP trends variables 29 30 USE sedmodel ! Sediment model 30 31 USE prtctl_trc ! print control for debugging … … 33 34 PRIVATE 34 35 35 PUBLIC p4z_sms_init ! called in p4zsms.F90 36 PUBLIC p4z_sms ! called in p4zsms.F90 37 38 REAL(wp) :: alkbudget, no3budget, silbudget, ferbudget 39 INTEGER :: numco2, numnut !: logical unit for co2 budget 40 36 PUBLIC p4z_sms_init ! called in p4zsms.F90 37 PUBLIC p4z_sms ! called in p4zsms.F90 38 39 REAL(wp) :: alkbudget, no3budget, silbudget, ferbudget, po4budget 40 REAL(wp) :: xfact1, xfact2, xfact3 41 INTEGER :: numco2, numnut, numnit !: logical unit for co2 budget 42 43 !!* Array used to indicate negative tracer values 44 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: xnegtr !: ??? 45 46 47 !! * Substitutions 48 # include "top_substitute.h90" 41 49 !!---------------------------------------------------------------------- 42 50 !! NEMO/TOP 3.3 , NEMO Consortium (2010) … … 61 69 INTEGER, INTENT( in ) :: kt ! ocean time-step index 62 70 !! 63 INTEGER :: jnt, jn, jl 71 INTEGER :: ji, jj, jk, jnt, jn, jl 72 REAL(wp) :: ztra 73 #if defined key_kriest 74 REAL(wp) :: zcoef1, zcoef2 75 #endif 64 76 CHARACTER (len=25) :: charout 65 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: ztrdpis66 77 !!--------------------------------------------------------------------- 67 78 ! 68 79 IF( nn_timing == 1 ) CALL timing_start('p4z_sms') 69 80 ! 70 IF( l_trdtrc ) THEN71 CALL wrk_alloc( jpi, jpj, jpk, jp_pisces, ztrdpis )72 DO jn = 1, jp_pisces73 jl = jn + jp_pcs0 - 174 ztrdpis(:,:,:,jn) = trn(:,:,:,jl)75 ENDDO76 ENDIF77 !78 81 IF( kt == nittrc000 ) THEN 82 ! 83 ALLOCATE( xnegtr(jpi,jpj,jpk) ) 79 84 ! 80 85 CALL p4z_che ! initialize the chemical constants … … 88 93 IF( ln_pisdmp .AND. MOD( kt - nn_dttrc, nn_pisdmp ) == 0 ) CALL p4z_dmp( kt ) ! Relaxation of some tracers 89 94 ! 95 ! ! set time step size (Euler/Leapfrog) 96 IF( ( neuler == 0 .AND. kt == nittrc000 ) .OR. ln_top_euler ) THEN ; rfact = rdttrc(1) ! at nittrc000 97 ELSEIF( kt <= nittrc000 + nn_dttrc ) THEN ; rfact = 2. * rdttrc(1) ! at nittrc000 or nittrc000+nn_dttrc (Leapfrog) 98 ENDIF 99 ! 100 IF( ( ln_top_euler .AND. kt == nittrc000 ) .OR. ( .NOT.ln_top_euler .AND. kt <= nittrc000 + nn_dttrc ) ) THEN 101 rfactr = 1. / rfact 102 rfact2 = rfact / FLOAT( nrdttrc ) 103 rfact2r = 1. / rfact2 104 xstep = rfact2 / rday ! Time step duration for biology 105 IF(lwp) WRITE(numout,*) 106 IF(lwp) WRITE(numout,*) ' Passive Tracer time step rfact = ', rfact, ' rdt = ', rdttra(1) 107 IF(lwp) write(numout,*) ' PISCES Biology time step rfact2 = ', rfact2 108 IF(lwp) WRITE(numout,*) 109 ENDIF 110 111 IF( ( neuler == 0 .AND. kt == nittrc000 ) .OR. ln_top_euler ) THEN 112 DO jn = jp_pcs0, jp_pcs1 ! SMS on tracer without Asselin time-filter 113 trb(:,:,:,jn) = trn(:,:,:,jn) 114 END DO 115 ENDIF 116 ! 90 117 IF( ndayflxtr /= nday_year ) THEN ! New days 91 118 ! … … 105 132 DO jnt = 1, nrdttrc ! Potential time splitting if requested 106 133 ! 107 CALL p4z_bio (kt, jnt) ! Biology 108 CALL p4z_sed (kt, jnt) ! Sedimentation 109 ! 134 CALL p4z_bio( kt, jnt ) ! Biology 135 CALL p4z_sed( kt, jnt ) ! Sedimentation 136 CALL p4z_lys( kt, jnt ) ! Compute CaCO3 saturation 137 CALL p4z_flx( kt, jnt ) ! Compute surface fluxes 138 ! 139 xnegtr(:,:,:) = 1.e0 110 140 DO jn = jp_pcs0, jp_pcs1 111 trb(:,:,:,jn) = trn(:,:,:,jn) 112 ENDDO 113 ! 141 DO jk = 1, jpk 142 DO jj = 1, jpj 143 DO ji = 1, jpi 144 IF( ( trb(ji,jj,jk,jn) + tra(ji,jj,jk,jn) ) < 0.e0 ) THEN 145 ztra = ABS( trb(ji,jj,jk,jn) ) / ( ABS( tra(ji,jj,jk,jn) ) + rtrn ) 146 xnegtr(ji,jj,jk) = MIN( xnegtr(ji,jj,jk), ztra ) 147 ENDIF 148 END DO 149 END DO 150 END DO 151 END DO 152 ! ! where at least 1 tracer concentration becomes negative 153 ! ! 154 DO jn = jp_pcs0, jp_pcs1 155 trb(:,:,:,jn) = trb(:,:,:,jn) + xnegtr(:,:,:) * tra(:,:,:,jn) 156 END DO 157 ! 158 DO jn = jp_pcs0, jp_pcs1 159 tra(:,:,:,jn) = 0._wp 160 END DO 161 ! 162 IF( ln_top_euler ) THEN 163 DO jn = jp_pcs0, jp_pcs1 164 trn(:,:,:,jn) = trb(:,:,:,jn) 165 END DO 166 ENDIF 114 167 END DO 115 168 116 IF( l_trdtrc ) THEN 117 DO jn = 1, jp_pisces 118 jl = jn + jp_pcs0 - 1 119 ztrdpis(:,:,:,jn) = ( ztrdpis(:,:,:,jn) - trn(:,:,:,jl) ) * rfact2r 120 ENDDO 121 ENDIF 122 123 CALL p4z_lys( kt ) ! Compute CaCO3 saturation 124 CALL p4z_flx( kt ) ! Compute surface fluxes 125 126 DO jn = jp_pcs0, jp_pcs1 127 CALL lbc_lnk( trn(:,:,:,jn), 'T', 1. ) 128 CALL lbc_lnk( trb(:,:,:,jn), 'T', 1. ) 129 CALL lbc_lnk( tra(:,:,:,jn), 'T', 1. ) 169 #if defined key_kriest 170 ! 171 zcoef1 = 1.e0 / xkr_massp 172 zcoef2 = 1.e0 / xkr_massp / 1.1 173 DO jk = 1,jpkm1 174 trb(:,:,jk,jpnum) = MAX( trb(:,:,jk,jpnum), trb(:,:,jk,jppoc) * zcoef1 / xnumm(jk) ) 175 trb(:,:,jk,jpnum) = MIN( trb(:,:,jk,jpnum), trb(:,:,jk,jppoc) * zcoef2 ) 130 176 END DO 131 177 ! 178 #endif 179 ! 180 ! 181 IF( l_trdtrc ) THEN 182 DO jn = jp_pcs0, jp_pcs1 183 CALL trd_trc( tra(:,:,:,jn), jn, jptra_sms, kt ) ! save trends 184 END DO 185 END IF 186 ! 132 187 IF( lk_sed ) THEN 133 188 ! … … 135 190 ! 136 191 DO jn = jp_pcs0, jp_pcs1 137 CALL lbc_lnk( tr n(:,:,:,jn), 'T', 1. )192 CALL lbc_lnk( trb(:,:,:,jn), 'T', 1. ) 138 193 END DO 139 194 ! … … 142 197 IF( lrst_trc ) CALL p4z_rst( kt, 'WRITE' ) !* Write PISCES informations in restart file 143 198 ! 144 IF( l_trdtrc ) THEN 145 DO jn = 1, jp_pisces 146 jl = jn + jp_pcs0 - 1 147 ztrdpis(:,:,:,jn) = ztrdpis(:,:,:,jn) + tra(:,:,:,jl) 148 CALL trd_mod_trc( ztrdpis(:,:,:,jn), jn, jptra_trd_sms, kt ) ! save trends 149 END DO 150 CALL wrk_dealloc( jpi, jpj, jpk, jp_pisces, ztrdpis ) 151 END IF 152 ! 153 CALL p4z_chk_mass( kt ) ! Mass conservation checking 199 200 IF( lk_iomput .OR. ln_check_mass ) CALL p4z_chk_mass( kt ) ! Mass conservation checking 154 201 155 202 IF ( lwm .AND. kt == nittrc000 ) CALL FLUSH ( numonp ) ! flush output namelist PISCES … … 281 328 ztmas = tmask(ji,jj,jk) 282 329 ztmas1 = 1. - tmask(ji,jj,jk) 283 zcaralk = tr n(ji,jj,jk,jptal) - borat(ji,jj,jk) / ( 1. + 1.E-8 / ( rtrn + akb3(ji,jj,jk) ) )284 zco3 = ( zcaralk - tr n(ji,jj,jk,jpdic) ) * ztmas + 0.5e-3 * ztmas1285 zbicarb = ( 2. * tr n(ji,jj,jk,jpdic) - zcaralk )330 zcaralk = trb(ji,jj,jk,jptal) - borat(ji,jj,jk) / ( 1. + 1.E-8 / ( rtrn + akb3(ji,jj,jk) ) ) 331 zco3 = ( zcaralk - trb(ji,jj,jk,jpdic) ) * ztmas + 0.5e-3 * ztmas1 332 zbicarb = ( 2. * trb(ji,jj,jk,jpdic) - zcaralk ) 286 333 hi(ji,jj,jk) = ( ak23(ji,jj,jk) * zbicarb / zco3 ) * ztmas + 1.e-9 * ztmas1 287 334 END DO … … 328 375 ENDIF 329 376 ! 377 IF( iom_varid( numrtr, 'tcflxcum', ldstop = .FALSE. ) > 0 ) THEN ! cumulative total flux of carbon 378 CALL iom_get( numrtr, 'tcflxcum' , t_oce_co2_flx_cum ) 379 ELSE 380 t_oce_co2_flx_cum = 0._wp 381 ENDIF 382 ! 330 383 ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN 331 384 IF( kt == nitrst ) THEN … … 337 390 CALL iom_rstput( kt, nitrst, numrtw, 'Silicalim', xksi(:,:) ) 338 391 CALL iom_rstput( kt, nitrst, numrtw, 'Silicamax', xksimax(:,:) ) 392 CALL iom_rstput( kt, nitrst, numrtw, 'tcflxcum', t_oce_co2_flx_cum ) 339 393 ENDIF 340 394 ! … … 355 409 REAL(wp) :: silmean = 91.51 ! mean value of silicate 356 410 ! 357 REAL(wp) :: zarea, zalksum, zpo4sum, zno3sum, zsilsum 411 REAL(wp) :: zarea, zalksumn, zpo4sumn, zno3sumn, zsilsumn 412 REAL(wp) :: zalksumb, zpo4sumb, zno3sumb, zsilsumb 358 413 !!--------------------------------------------------------------------- 359 414 … … 368 423 zarea = 1._wp / glob_sum( cvol(:,:,:) ) * 1e6 369 424 370 zalksum = glob_sum( trn(:,:,:,jptal) * cvol(:,:,:) ) * zarea371 zpo4sum = glob_sum( trn(:,:,:,jppo4) * cvol(:,:,:) ) * zarea * po4r372 zno3sum = glob_sum( trn(:,:,:,jpno3) * cvol(:,:,:) ) * zarea * rno3373 zsilsum = glob_sum( trn(:,:,:,jpsil) * cvol(:,:,:) ) * zarea425 zalksumn = glob_sum( trn(:,:,:,jptal) * cvol(:,:,:) ) * zarea 426 zpo4sumn = glob_sum( trn(:,:,:,jppo4) * cvol(:,:,:) ) * zarea * po4r 427 zno3sumn = glob_sum( trn(:,:,:,jpno3) * cvol(:,:,:) ) * zarea * rno3 428 zsilsumn = glob_sum( trn(:,:,:,jpsil) * cvol(:,:,:) ) * zarea 374 429 375 IF(lwp) WRITE(numout,*) ' TALK mean : ', zalksum 376 trn(:,:,:,jptal) = trn(:,:,:,jptal) * alkmean / zalksum 377 378 IF(lwp) WRITE(numout,*) ' PO4 mean : ', zpo4sum 379 trn(:,:,:,jppo4) = trn(:,:,:,jppo4) * po4mean / zpo4sum 380 381 IF(lwp) WRITE(numout,*) ' NO3 mean : ', zno3sum 382 trn(:,:,:,jpno3) = trn(:,:,:,jpno3) * no3mean / zno3sum 383 384 IF(lwp) WRITE(numout,*) ' SiO3 mean : ', zsilsum 385 trn(:,:,:,jpsil) = MIN( 400.e-6,trn(:,:,:,jpsil) * silmean / zsilsum ) 386 ! 387 ENDIF 388 430 IF(lwp) WRITE(numout,*) ' TALKN mean : ', zalksumn 431 trn(:,:,:,jptal) = trn(:,:,:,jptal) * alkmean / zalksumn 432 433 IF(lwp) WRITE(numout,*) ' PO4N mean : ', zpo4sumn 434 trn(:,:,:,jppo4) = trn(:,:,:,jppo4) * po4mean / zpo4sumn 435 436 IF(lwp) WRITE(numout,*) ' NO3N mean : ', zno3sumn 437 trn(:,:,:,jpno3) = trn(:,:,:,jpno3) * no3mean / zno3sumn 438 439 IF(lwp) WRITE(numout,*) ' SiO3N mean : ', zsilsumn 440 trn(:,:,:,jpsil) = MIN( 400.e-6,trn(:,:,:,jpsil) * silmean / zsilsumn ) 441 ! 442 ! 443 IF( .NOT. ln_top_euler ) THEN 444 zalksumb = glob_sum( trb(:,:,:,jptal) * cvol(:,:,:) ) * zarea 445 zpo4sumb = glob_sum( trb(:,:,:,jppo4) * cvol(:,:,:) ) * zarea * po4r 446 zno3sumb = glob_sum( trb(:,:,:,jpno3) * cvol(:,:,:) ) * zarea * rno3 447 zsilsumb = glob_sum( trb(:,:,:,jpsil) * cvol(:,:,:) ) * zarea 448 449 IF(lwp) WRITE(numout,*) ' ' 450 IF(lwp) WRITE(numout,*) ' TALKB mean : ', zalksumb 451 trb(:,:,:,jptal) = trb(:,:,:,jptal) * alkmean / zalksumb 452 453 IF(lwp) WRITE(numout,*) ' PO4B mean : ', zpo4sumb 454 trb(:,:,:,jppo4) = trb(:,:,:,jppo4) * po4mean / zpo4sumb 455 456 IF(lwp) WRITE(numout,*) ' NO3B mean : ', zno3sumb 457 trb(:,:,:,jpno3) = trb(:,:,:,jpno3) * no3mean / zno3sumb 458 459 IF(lwp) WRITE(numout,*) ' SiO3B mean : ', zsilsumb 460 trb(:,:,:,jpsil) = MIN( 400.e-6,trb(:,:,:,jpsil) * silmean / zsilsumb ) 461 ENDIF 462 ! 463 ENDIF 464 ! 389 465 END SUBROUTINE p4z_dmp 390 466 … … 399 475 ! 400 476 INTEGER, INTENT( in ) :: kt ! ocean time-step index 401 !! 477 REAL(wp) :: zrdenittot, zsdenittot, znitrpottot 478 CHARACTER(LEN=100) :: cltxt 479 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zvol 480 INTEGER :: jk 481 !!---------------------------------------------------------------------- 482 483 ! 402 484 !!--------------------------------------------------------------------- 403 485 … … 406 488 CALL ctl_opn( numco2, 'carbon.budget' , 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea ) 407 489 CALL ctl_opn( numnut, 'nutrient.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea ) 490 CALL ctl_opn( numnit, 'nitrogen.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea ) 491 xfact1 = rfact2r * 12. / 1.e15 * ryyss ! conversion molC/kt --> PgC/yr 492 xfact2 = 1.e+3 * rno3 * 14. / 1.e12 * ryyss ! conversion molC/l/s ----> TgN/m3/yr 493 xfact3 = 1.e+3 * rfact2r * rno3 ! conversion molC/l/kt ----> molN/m3/s 494 cltxt='time-step Alkalinity Nitrate Phosphorus Silicate Iron' 495 IF( lwp ) WRITE(numnut,*) TRIM(cltxt) 496 IF( lwp ) WRITE(numnut,*) 408 497 ENDIF 409 498 ENDIF 410 499 411 IF( ln_check_mass .AND. kt == nitend ) THEN ! Compute the budget of NO3, ALK, Si, Fer 500 ! 501 IF( iom_use( "pno3tot" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN 502 ! Compute the budget of NO3, ALK, Si, Fer 412 503 no3budget = glob_sum( ( trn(:,:,:,jpno3) + trn(:,:,:,jpnh4) & 413 504 & + trn(:,:,:,jpphy) + trn(:,:,:,jpdia) & … … 417 508 & + trn(:,:,:,jpgoc) & 418 509 #endif 419 & + trn(:,:,:,jpdoc) ) * cvol(:,:,:) ) 420 ! 510 & + trn(:,:,:,jpdoc) ) * cvol(:,:,:) ) 511 ! 512 no3budget = no3budget / areatot 513 CALL iom_put( "pno3tot", no3budget ) 514 ENDIF 515 ! 516 IF( iom_use( "ppo4tot" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN 517 po4budget = glob_sum( ( trn(:,:,:,jppo4) & 518 & + trn(:,:,:,jpphy) + trn(:,:,:,jpdia) & 519 & + trn(:,:,:,jpzoo) + trn(:,:,:,jpmes) & 520 & + trn(:,:,:,jppoc) & 521 #if ! defined key_kriest 522 & + trn(:,:,:,jpgoc) & 523 #endif 524 & + trn(:,:,:,jpdoc) ) * cvol(:,:,:) ) 525 po4budget = po4budget / areatot 526 CALL iom_put( "ppo4tot", po4budget ) 527 ENDIF 528 ! 529 IF( iom_use( "psiltot" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN 421 530 silbudget = glob_sum( ( trn(:,:,:,jpsil) + trn(:,:,:,jpgsi) & 422 531 & + trn(:,:,:,jpdsi) ) * cvol(:,:,:) ) 423 ! 532 ! 533 silbudget = silbudget / areatot 534 CALL iom_put( "psiltot", silbudget ) 535 ENDIF 536 ! 537 IF( iom_use( "palktot" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN 424 538 alkbudget = glob_sum( ( trn(:,:,:,jpno3) * rno3 & 425 539 & + trn(:,:,:,jptal) & 426 540 & + trn(:,:,:,jpcal) * 2. ) * cvol(:,:,:) ) 427 ! 541 ! 542 alkbudget = alkbudget / areatot 543 CALL iom_put( "palktot", alkbudget ) 544 ENDIF 545 ! 546 IF( iom_use( "pfertot" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN 428 547 ferbudget = glob_sum( ( trn(:,:,:,jpfer) + trn(:,:,:,jpnfe) & 429 548 & + trn(:,:,:,jpdfe) & … … 434 553 & + trn(:,:,:,jpzoo) * ferat3 & 435 554 & + trn(:,:,:,jpmes) * ferat3 ) * cvol(:,:,:) ) 436 437 ! 555 ! 556 ferbudget = ferbudget / areatot 557 CALL iom_put( "pfertot", ferbudget ) 558 ENDIF 559 ! 560 561 ! Global budget of N SMS : denitrification in the water column and in the sediment 562 ! nitrogen fixation by the diazotrophs 563 ! -------------------------------------------------------------------------------- 564 IF( iom_use( "tnfix" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN 565 znitrpottot = glob_sum ( nitrpot(:,:,:) * nitrfix * cvol(:,:,:) ) 566 CALL iom_put( "tnfix" , znitrpottot * 1.e+3 * rno3 ) ! Global nitrogen fixation molC/l to molN/m3 567 ENDIF 568 ! 569 IF( iom_use( "tdenit" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN 570 zrdenittot = glob_sum ( denitr(:,:,:) * rdenit * xnegtr(:,:,:) * cvol(:,:,:) ) 571 CALL iom_put( "tdenit" , zrdenittot * 1.e+3 * rno3 ) ! Total denitrification molC/l to molN/m3 572 ENDIF 573 ! 574 IF( iom_use( "Sdenit" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN 575 zsdenittot = glob_sum ( sdenit(:,:) * e1e2t(:,:) ) 576 CALL iom_put( "Sdenit", sdenit(:,:) * xfact3 * tmask(:,:,1) ) ! Nitrate reduction in the sediments 577 ENDIF 578 579 IF( ln_check_mass .AND. kt == nitend ) THEN ! Compute the budget of NO3, ALK, Si, Fer 438 580 t_atm_co2_flx = t_atm_co2_flx / glob_sum( e1e2t(:,:) ) 439 t_oce_co2_flx = t_oce_co2_flx * 12. / 1.e15* (-1 )440 tpp = tpp * 1000. * 12. / 1.E15441 t_oce_co2_exp = t_oce_co2_exp * 1000. * 12. / 1.E15442 !443 no3budget = no3budget / areatot444 silbudget = silbudget / areatot445 alkbudget = alkbudget / areatot446 ferbudget = ferbudget / areatot447 !448 IF(lwp) THEN449 WRITE(numco2,9000) ndastp, t_atm_co2_flx, t_oce_co2_flx, tpp, t_oce_co2_exp450 WRITE(numnut,9500) ndastp, alkbudget, no3budget, silbudget, ferbudget451 ENDIF452 ! 453 ENDIF 454 581 t_oce_co2_flx = t_oce_co2_flx * xfact1 * (-1 ) 582 tpp = tpp * 1000. * xfact1 583 t_oce_co2_exp = t_oce_co2_exp * 1000. * xfact1 584 IF( lwp ) WRITE(numco2,9000) ndastp, t_atm_co2_flx, t_oce_co2_flx, tpp, t_oce_co2_exp 585 IF( lwp ) WRITE(numnut,9100) ndastp, alkbudget * 1.e+06, & 586 & no3budget * rno3 * 1.e+06, & 587 & po4budget * po4r * 1.e+06, & 588 & silbudget * 1.e+06, & 589 & ferbudget * 1.e+09 590 ! 591 IF( lwp ) WRITE(numnit,9200) ndastp, znitrpottot * xfact2 , & 592 & zrdenittot * xfact2 , & 593 & zsdenittot * xfact2 594 595 ENDIF 596 ! 455 597 9000 FORMAT(i8,f10.5,e18.10,f10.5,f10.5) 456 9500 FORMAT(i8,4e18.10) 598 9100 FORMAT(i8,5e18.10) 599 9200 FORMAT(i8,3f10.5) 600 457 601 ! 458 602 END SUBROUTINE p4z_chk_mass
Note: See TracChangeset
for help on using the changeset viewer.
