Changeset 10302 for branches/UKMO/dev_r5518_GO6_package_FOAMv14/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_fin.F90
- Timestamp:
- 2018-11-13T18:21:16+01:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5518_GO6_package_FOAMv14/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_fin.F90
r8442 r10302 6 6 !! History : 7 7 !! - ! 2017-04 (M. Stringer) Code taken from trcbio_medusa.F90 8 !! - ! 2017-08 (A. Yool) Amend bethic reservoir updating 8 9 !!---------------------------------------------------------------------- 9 10 #if defined key_medusa … … 32 33 !!---------------------------------------------------------------------- 33 34 USE bio_medusa_mod 34 USE dom_oce, ONLY: atfp, atfp1, neuler, rdt, e3t_n,tmask35 USE dom_oce, ONLY: atfp, atfp1, neuler, rdt, tmask 35 36 USE in_out_manager, ONLY: lwp, numout 36 # if defined key_iomput 37 USE iom, ONLY: iom_put, lk_iomput 38 # endif 37 USE iom, ONLY: iom_put 39 38 USE lbclnk, ONLY: lbc_lnk 39 USE oce, ONLY: chloro_out_cpl 40 40 USE par_medusa, ONLY: jp_medusa_2d, jp_medusa_3d, & 41 jp_medusa_trd 41 jp_medusa_trd, jpchd, jpchn 42 42 USE par_oce, ONLY: jpi, jpim1, jpj, jpjm1, jpk 43 43 USE phycst, ONLY: rsmall 44 USE sbc_oce, ONLY: lk_oasis 44 45 USE sms_medusa, ONLY: jinorgben, jorgben, & 45 46 f3_co3, f3_h2co3, f3_hco3, & 46 47 f3_omarg, f3_omcal, f3_pH, & 48 ln_foam_medusa, mld_max, pgrow_avg, & 49 ploss_avg, phyt_avg, & 47 50 za_sed_c, za_sed_ca, za_sed_fe, & 48 51 za_sed_n, za_sed_si, & … … 50 53 zb_sed_n, zb_sed_si, & 51 54 zn_sed_c, zn_sed_ca, zn_sed_fe, & 52 zn_sed_n, zn_sed_si 53 USE trc, ONLY: med_diag, nittrc000 55 zn_sed_n, zn_sed_si, zn_chl_srf, & 56 scl_chl, chl_out 57 USE trc, ONLY: med_diag, nittrc000, trn 54 58 USE trcnam_trp, ONLY: ln_trcadv_cen2, ln_trcadv_tvd 59 USE zdfmxl, ONLY: hmld 55 60 56 61 !! time (integer timestep) … … 62 67 REAL(wp) :: fq0,fq1,fq2,fq3 63 68 69 # if defined key_roam 70 !!---------------------------------------------------------------------- 71 !! AXY (09/08/17): fix benthic submodel 64 72 !!---------------------------------------------------------------------- 65 73 !! Process benthic in/out fluxes 66 74 !! These can be handled outside of the 3D calculations since the 67 !! benthic pools (and fluxes) are 2D in nature; this code is68 !! (shamelessly) borrowed from corresponding code in the LOBSTER69 !! model75 !! benthic pools (and fluxes) are 2D in nature; this code was 76 !! developed with help from George Nurser (NOC); it cannot be run 77 !! in a configuration with variable time-stepping with depth 70 78 !!---------------------------------------------------------------------- 71 79 !! 72 !! IF(lwp) WRITE(numout,*) 'AXY: rdt = ', rdt80 !! time-step calculation 73 81 if (jorgben.eq.1) then 74 za_sed_n(:,:) = zn_sed_n(:,:) + & 75 ( f_sbenin_n(:,:) + f_fbenin_n(:,:) - & 76 f_benout_n(:,:) ) * (rdt / 86400.) 82 za_sed_n(:,:) = zb_sed_n(:,:) + ((2. * (rdt / 86400.)) * & 83 ( f_sbenin_n(:,:) + f_fbenin_n(:,:) - f_benout_n(:,:) )) 84 za_sed_fe(:,:) = zb_sed_fe(:,:) + ((2. * (rdt / 86400.)) * & 85 ( f_sbenin_fe(:,:) + f_fbenin_fe(:,:) - f_benout_fe(:,:) )) 86 za_sed_c(:,:) = zb_sed_c(:,:) + ((2. * (rdt / 86400.)) * & 87 ( f_sbenin_c(:,:) + f_fbenin_c(:,:) - f_benout_c(:,:) )) 88 endif 89 if (jinorgben.eq.1) then 90 za_sed_si(:,:) = zb_sed_si(:,:) + ((2. * (rdt / 86400.)) * & 91 ( f_fbenin_si(:,:) - f_benout_si(:,:) )) 92 za_sed_ca(:,:) = zb_sed_ca(:,:) + ((2. * (rdt / 86400.)) * & 93 ( f_fbenin_ca(:,:) - f_benout_ca(:,:) )) 94 endif 95 !! 96 !! time-level calculation 97 if (jorgben.eq.1) then 98 zb_sed_n(:,:) = zn_sed_n(:,:) + (atfp * & 99 ( za_sed_n(:,:) - (2. * zn_sed_n(:,:)) + zb_sed_n(:,:) )) 77 100 zn_sed_n(:,:) = za_sed_n(:,:) 78 !! 79 za_sed_fe(:,:) = zn_sed_fe(:,:) + & 80 ( f_sbenin_fe(:,:) + f_fbenin_fe(:,:) - & 81 f_benout_fe(:,:) ) * (rdt / 86400.) 101 zb_sed_fe(:,:) = zn_sed_fe(:,:) + (atfp * & 102 ( za_sed_fe(:,:) - (2. * zn_sed_fe(:,:)) + zb_sed_fe(:,:) )) 82 103 zn_sed_fe(:,:) = za_sed_fe(:,:) 83 !! 84 za_sed_c(:,:) = zn_sed_c(:,:) + & 85 ( f_sbenin_c(:,:) + f_fbenin_c(:,:) - & 86 f_benout_c(:,:) ) * (rdt / 86400.) 104 zb_sed_c(:,:) = zn_sed_c(:,:) + (atfp * & 105 ( za_sed_c(:,:) - (2. * zn_sed_c(:,:)) + zb_sed_c(:,:) )) 87 106 zn_sed_c(:,:) = za_sed_c(:,:) 88 107 endif 89 108 if (jinorgben.eq.1) then 90 za_sed_si(:,:) = zn_sed_si(:,:) + & 91 ( f_fbenin_si(:,:) - f_benout_si(:,:) ) * & 92 (rdt / 86400.) 109 zb_sed_si(:,:) = zn_sed_si(:,:) + (atfp * & 110 ( za_sed_si(:,:) - (2. * zn_sed_si(:,:)) + zb_sed_si(:,:) )) 93 111 zn_sed_si(:,:) = za_sed_si(:,:) 94 !! 95 za_sed_ca(:,:) = zn_sed_ca(:,:) + & 96 ( f_fbenin_ca(:,:) - f_benout_ca(:,:) ) * & 97 (rdt / 86400.) 112 zb_sed_ca(:,:) = zn_sed_ca(:,:) + (atfp * & 113 ( za_sed_ca(:,:) - (2. * zn_sed_ca(:,:)) + zb_sed_ca(:,:) )) 98 114 zn_sed_ca(:,:) = za_sed_ca(:,:) 99 115 endif 100 !! 101 if (ibenthic.eq.2) then 102 !! The code below (in this if ... then ... endif loop) is 103 !! effectively commented out because it does not work as 104 !! anticipated; it can be deleted at a later date 105 if (jorgben.eq.1) then 106 za_sed_n(:,:) = ( f_sbenin_n(:,:) + f_fbenin_n(:,:) - & 107 f_benout_n(:,:) ) * rdt 108 za_sed_fe(:,:) = ( f_sbenin_fe(:,:) + f_fbenin_fe(:,:) - & 109 f_benout_fe(:,:) ) * rdt 110 za_sed_c(:,:) = ( f_sbenin_c(:,:) + f_fbenin_c(:,:) - & 111 f_benout_c(:,:) ) * rdt 112 endif 113 if (jinorgben.eq.1) then 114 za_sed_si(:,:) = ( f_fbenin_si(:,:) - f_benout_si(:,:) ) * rdt 115 za_sed_ca(:,:) = ( f_fbenin_ca(:,:) - f_benout_ca(:,:) ) * rdt 116 endif 116 # endif 117 118 IF ( ln_foam_medusa ) THEN 119 !!---------------------------------------------------------------------- 120 !! Diagnostics required for ocean colour assimilation: 121 !! Mixed layer average phytoplankton growth, loss and concentration 122 !! Maximum mixed layer depth 123 !!---------------------------------------------------------------------- 117 124 !! 118 !! Leap-frog scheme - only in explicit case, otherwise the 119 !! time stepping is already being done in trczdf 120 !! IF( l_trczdf_exp .AND. (ln_trcadv_cen2 .OR. ln_trcadv_tvd) ) THEN 121 !! zfact = 2. * rdttra(jk) * FLOAT( ndttrc ) 122 !! IF( neuler == 0 .AND. kt == nittrc000 ) zfact = rdttra(jk) * 123 !! FLOAT(ndttrc) 124 !! if (jorgben.eq.1) then 125 !! za_sed_n(:,:) = zb_sed_n(:,:) + ( zfact * za_sed_n(:,:) ) 126 !! za_sed_fe(:,:) = zb_sed_fe(:,:) + ( zfact * za_sed_fe(:,:) ) 127 !! za_sed_c(:,:) = zb_sed_c(:,:) + ( zfact * za_sed_c(:,:) ) 128 !! endif 129 !! if (jinorgben.eq.1) then 130 !! za_sed_si(:,:) = zb_sed_si(:,:) + ( zfact * za_sed_si(:,:) ) 131 !! za_sed_ca(:,:) = zb_sed_ca(:,:) + ( zfact * za_sed_ca(:,:) ) 132 !! endif 133 !! ENDIF 134 !! 135 !! Time filter and swap of arrays 136 IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd ) THEN ! centred or tvd scheme 137 IF( neuler == 0 .AND. kt == nittrc000 ) THEN 138 if (jorgben.eq.1) then 139 zb_sed_n(:,:) = zn_sed_n(:,:) 140 zn_sed_n(:,:) = za_sed_n(:,:) 141 za_sed_n(:,:) = 0.0 142 !! 143 zb_sed_fe(:,:) = zn_sed_fe(:,:) 144 zn_sed_fe(:,:) = za_sed_fe(:,:) 145 za_sed_fe(:,:) = 0.0 146 !! 147 zb_sed_c(:,:) = zn_sed_c(:,:) 148 zn_sed_c(:,:) = za_sed_c(:,:) 149 za_sed_c(:,:) = 0.0 150 endif 151 if (jinorgben.eq.1) then 152 zb_sed_si(:,:) = zn_sed_si(:,:) 153 zn_sed_si(:,:) = za_sed_si(:,:) 154 za_sed_si(:,:) = 0.0 155 !! 156 zb_sed_ca(:,:) = zn_sed_ca(:,:) 157 zn_sed_ca(:,:) = za_sed_ca(:,:) 158 za_sed_ca(:,:) = 0.0 159 endif 160 ELSE 161 if (jorgben.eq.1) then 162 zb_sed_n(:,:) = (atfp * & 163 ( zb_sed_n(:,:) + za_sed_n(:,:) )) + & 164 (atfp1 * zn_sed_n(:,:) ) 165 zn_sed_n(:,:) = za_sed_n(:,:) 166 za_sed_n(:,:) = 0.0 167 !! 168 zb_sed_fe(:,:) = (atfp * & 169 ( zb_sed_fe(:,:) + za_sed_fe(:,:) )) + & 170 (atfp1 * zn_sed_fe(:,:)) 171 zn_sed_fe(:,:) = za_sed_fe(:,:) 172 za_sed_fe(:,:) = 0.0 173 !! 174 zb_sed_c(:,:) = (atfp * & 175 ( zb_sed_c(:,:) + za_sed_c(:,:) )) + & 176 (atfp1 * zn_sed_c(:,:) ) 177 zn_sed_c(:,:) = za_sed_c(:,:) 178 za_sed_c(:,:) = 0.0 179 endif 180 if (jinorgben.eq.1) then 181 zb_sed_si(:,:) = (atfp * & 182 ( zb_sed_si(:,:) + za_sed_si(:,:) )) + & 183 (atfp1 * zn_sed_si(:,:)) 184 zn_sed_si(:,:) = za_sed_si(:,:) 185 za_sed_si(:,:) = 0.0 186 !! 187 zb_sed_ca(:,:) = (atfp * & 188 ( zb_sed_ca(:,:) + za_sed_ca(:,:) )) + & 189 (atfp1 * zn_sed_ca(:,:)) 190 zn_sed_ca(:,:) = za_sed_ca(:,:) 191 za_sed_ca(:,:) = 0.0 192 endif 193 ENDIF 194 ELSE ! case of smolar scheme or muscl 195 if (jorgben.eq.1) then 196 zb_sed_n(:,:) = za_sed_n(:,:) 197 zn_sed_n(:,:) = za_sed_n(:,:) 198 za_sed_n(:,:) = 0.0 199 !! 200 zb_sed_fe(:,:) = za_sed_fe(:,:) 201 zn_sed_fe(:,:) = za_sed_fe(:,:) 202 za_sed_fe(:,:) = 0.0 203 !! 204 zb_sed_c(:,:) = za_sed_c(:,:) 205 zn_sed_c(:,:) = za_sed_c(:,:) 206 za_sed_c(:,:) = 0.0 207 endif 208 if (jinorgben.eq.1) then 209 zb_sed_si(:,:) = za_sed_si(:,:) 210 zn_sed_si(:,:) = za_sed_si(:,:) 211 za_sed_si(:,:) = 0.0 212 !! 213 zb_sed_ca(:,:) = za_sed_ca(:,:) 214 zn_sed_ca(:,:) = za_sed_ca(:,:) 215 za_sed_ca(:,:) = 0.0 216 endif 217 ENDIF 218 endif 125 DO jj = 2,jpjm1 126 DO ji = 2,jpim1 127 IF ( hmld(ji,jj) .GT. 0.0 ) THEN 128 pgrow_avg(ji,jj) = pgrow_avg(ji,jj) / hmld(ji,jj) 129 ploss_avg(ji,jj) = ploss_avg(ji,jj) / hmld(ji,jj) 130 phyt_avg(ji,jj) = phyt_avg(ji,jj) / hmld(ji,jj) 131 IF ( hmld(ji,jj) .GT. mld_max(ji,jj) ) THEN 132 mld_max(ji,jj) = hmld(ji,jj) 133 ENDIF 134 ENDIF 135 END DO 136 END DO 137 ENDIF 219 138 220 139 # if defined key_debug_medusa … … 253 172 fq1 = f_sbenin_n(ji,jj) + f_fbenin_n(ji,jj) 254 173 fq2 = fq0 + fq1 255 IF (lwp) write (numout,'(a,2i3,a,3f15.10)') & 256 'AXY N cons: (i,j)=',ji,jj,', (flx,ben,err)=', & 257 fq0,fq1,fq2 174 fq3 = f_benout_n(ji,jj) 175 if (lwp) write (numout,'a,2i3,a,4f15,5)') & 176 'AXY N cons: (i,j)=',ji,jj,', (flx,ben,err,out)=', & 177 fq0,fq1,fq2,fq3 258 178 ENDIF 259 179 ENDDO … … 266 186 fq1 = f_fbenin_si(ji,jj) 267 187 fq2 = fq0 + fq1 268 IF (lwp) write (numout,'(a,2i3,a,3f15.10)') & 269 'AXY Si cons: (i,j)=',ji,jj,', (flx,ben,err)=', & 270 fq0,fq1,fq2 188 fq3 = f_benout_si(ji,jj) 189 if (lwp) write (numout,'a,2i3,a,4f15,5)') & 190 'AXY Si cons: (i,j)=',ji,jj,', (flx,ben,err,out)=', & 191 fq0,fq1,fq2,fq3 271 192 ENDIF 272 193 ENDDO … … 278 199 fq0 = fflx_c(ji,jj) 279 200 fq1 = f_sbenin_c(ji,jj) + f_fbenin_c(ji,jj) + f_fbenin_ca(ji,jj) 280 fq2 = f_co2flux(ji,jj) * e3t_n(ji,jj,1)201 fq2 = f_co2flux(ji,jj) * fse3t(ji,jj,1) 281 202 fq3 = fq0 + fq1 282 IF (lwp) write (numout,'(a,2i3,a,4f15.10)') & 283 'AXY C cons: (i,j)=',ji,jj,', (flx,ben,asf,err)=', & 284 fq0,fq1,fq2,fq3 285 ENDIF 286 ENDDO 287 ENDDO 288 !! alkalinity 289 DO jj = 2,jpjm1 290 DO ji = 2,jpim1 291 if (tmask(ji,jj,1) == 1) then 292 fq0 = fflx_a(ji,jj) 293 fq1 = 2.0 * f_fbenin_ca(ji,jj) 294 fq2 = fq0 + fq1 295 IF (lwp) write (numout,'(a,2i3,a,3f15.10)') & 296 'AXY alk cons: (i,j)=',ji,jj,', (flx,ben,err)=', & 297 fq0,fq1,fq2 203 fq4 = f_benout_c(ji,jj) + f_benout_ca(ji,jj) 204 if (lwp) write (numout,'a,2i3,a,5f15,5)') & 205 'AXY C cons: (i,j)=',ji,jj,', (flx,ben,asf,err,out)=', & 206 fq0,fq1,fq2,fq3,fq4 207 ENDIF 208 ENDDO 209 ENDDO 210 !! alkalinity 211 DO jj = 2,jpjm1 212 DO ji = 2,jpim1 213 if (tmask(ji,jj,1) == 1) then 214 fq0 = fflx_a(ji,jj) 215 fq1 = 2.0 * f_fbenin_ca(ji,jj) 216 fq2 = fq0 + fq1 217 fq3 = 2.0 * f_benout_ca(ji,jj) 218 if (lwp) write (numout,'a,2i3,a,4f15,5)') & 219 'AXY alk cons: (i,j)=',ji,jj,', (flx,ben,err,out)=', & 220 fq0,fq1,fq2,fq3 298 221 ENDIF 299 222 ENDDO … … 333 256 ENDDO 334 257 ENDDO 258 259 !!!--------------------------------------------------------------- 260 !! Calculates Chl diag for UM coupling 261 !!!--------------------------------------------------------------- 262 !! JPALM -- 02-06-2017 -- 263 !! add Chl surf coupling 264 !! no need to output, just pass to cpl var 265 IF (lk_oasis) THEN 266 IF (chl_out.eq.1) THEN 267 !! export and scale surface chl 268 zn_chl_srf(:,:) = MAX( 0.0, (trn(:,:,1,jpchd) + trn(:,:,1,jpchn)) * 1.0E-6 ) 269 !! surf Chl in Kg-chl/m3 as needed for cpl 270 ELSEIF (chl_out.eq.2) THEN 271 !! export and scale mld chl 272 zn_chl_srf(:,:) = MAX( 0.0, fchl_ml(:,:) * 1.0E-6 ) 273 !! mld Chl in Kg-chl/m3 as needed for cpl 274 ENDIF 275 chloro_out_cpl(:,:) = zn_chl_srf(:,:) * scl_chl !! Coupling Chl 276 END IF 277 335 278 !!---------------------------------------------------------------- 336 279 !! Add in XML diagnostics stuff … … 360 303 CALL iom_put( "OCAL_LVL" , fccd ) 361 304 ENDIF 305 IF ( med_diag%CHL_MLD%dgsave ) THEN 306 CALL iom_put( "CHL_MLD" , fchl_ml ) 307 ENDIF 308 IF (lk_oasis) THEN 309 IF ( med_diag%CHL_CPL%dgsave ) THEN 310 CALL iom_put( "CHL_CPL" , chloro_out_cpl ) 311 ENDIF 312 ENDIF 362 313 IF ( med_diag%PN_JLIM%dgsave ) THEN 363 314 CALL iom_put( "PN_JLIM" , fjln2d )
Note: See TracChangeset
for help on using the changeset viewer.
