Changeset 14789 for NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU/src/OCE/TRA/traqsr.F90
- Timestamp:
- 2021-05-05T13:18:04+02:00 (3 years ago)
- Location:
- NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
src/OCE/TRA/traqsr.F90 (modified) (21 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU
- Property svn:externals
-
old new 3 3 ^/utils/build/mk@HEAD mk 4 4 ^/utils/tools@HEAD tools 5 ^/vendors/AGRIF/dev _r12970_AGRIF_CMEMSext/AGRIF5 ^/vendors/AGRIF/dev@HEAD ext/AGRIF 6 6 ^/vendors/FCM@HEAD ext/FCM 7 7 ^/vendors/IOIPSL@HEAD ext/IOIPSL 8 ^/vendors/PPR@HEAD ext/PPR 8 9 9 10 # SETTE 10 ^/utils/CI/sette@1 3559sette11 ^/utils/CI/sette@14244 sette
-
- Property svn:externals
-
NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU/src/OCE/TRA/traqsr.F90
r13497 r14789 9 9 !! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module 10 10 !! - ! 2005-11 (G. Madec) zco, zps, sco coordinate 11 !! 3.2 ! 2009-04 (G. Madec & NEMO team) 12 !! 3.6 ! 2012-05 (C. Rousset) store attenuation coef for use in ice model 11 !! 3.2 ! 2009-04 (G. Madec & NEMO team) 12 !! 3.6 ! 2012-05 (C. Rousset) store attenuation coef for use in ice model 13 13 !! 3.6 ! 2015-12 (O. Aumont, J. Jouanno, C. Ethe) use vertical profile of chlorophyll 14 !! 3.7 ! 2015-11 (G. Madec, A. Coward) remove optimisation for fix volume 14 !! 3.7 ! 2015-11 (G. Madec, A. Coward) remove optimisation for fix volume 15 15 !!---------------------------------------------------------------------- 16 16 17 17 !!---------------------------------------------------------------------- 18 !! tra_qsr : temperature trend due to the penetration of solar radiation 19 !! tra_qsr_init : initialization of the qsr penetration 18 !! tra_qsr : temperature trend due to the penetration of solar radiation 19 !! tra_qsr_init : initialization of the qsr penetration 20 20 !!---------------------------------------------------------------------- 21 21 USE oce ! ocean dynamics and active tracers 22 22 USE phycst ! physical constants 23 23 USE dom_oce ! ocean space and time domain 24 USE domtile 24 25 USE sbc_oce ! surface boundary condition: ocean 25 26 USE trc_oce ! share SMS/Ocean variables … … 44 45 ! !!* Namelist namtra_qsr: penetrative solar radiation 45 46 LOGICAL , PUBLIC :: ln_traqsr !: light absorption (qsr) flag 46 LOGICAL , PUBLIC :: ln_qsr_rgb !: Red-Green-Blue light absorption flag 47 LOGICAL , PUBLIC :: ln_qsr_rgb !: Red-Green-Blue light absorption flag 47 48 LOGICAL , PUBLIC :: ln_qsr_2bd !: 2 band light absorption flag 48 49 LOGICAL , PUBLIC :: ln_qsr_bio !: bio-model light absorption flag … … 53 54 ! 54 55 INTEGER , PUBLIC :: nksr !: levels below which the light cannot penetrate (depth larger than 391 m) 55 56 56 57 INTEGER, PARAMETER :: np_RGB = 1 ! R-G-B light penetration with constant Chlorophyll 57 58 INTEGER, PARAMETER :: np_RGBc = 2 ! R-G-B light penetration with Chlorophyll data … … 87 88 !! Considering the 2 wavebands case: 88 89 !! I(k) = Qsr*( rn_abs*EXP(z(k)/rn_si0) + (1.-rn_abs)*EXP(z(k)/rn_si1) ) 89 !! The temperature trend associated with the solar radiation penetration 90 !! The temperature trend associated with the solar radiation penetration 90 91 !! is given by : zta = 1/e3t dk[ I ] / (rho0*Cp) 91 92 !! At the bottom, boudary condition for the radiation is no flux : 92 93 !! all heat which has not been absorbed in the above levels is put 93 94 !! in the last ocean level. 94 !! The computation is only done down to the level where 95 !! I(k) < 1.e-15 W/m2 (i.e. over the top nksr levels) . 95 !! The computation is only done down to the level where 96 !! I(k) < 1.e-15 W/m2 (i.e. over the top nksr levels) . 96 97 !! 97 98 !! ** Action : - update ta with the penetrative solar radiation trend … … 107 108 ! 108 109 INTEGER :: ji, jj, jk ! dummy loop indices 109 INTEGER :: irgb 110 INTEGER :: irgb, isi, iei, isj, iej ! local integers 110 111 REAL(wp) :: zchl, zcoef, z1_2 ! local scalars 111 112 REAL(wp) :: zc0 , zc1 , zc2 , zc3 ! - - … … 120 121 IF( ln_timing ) CALL timing_start('tra_qsr') 121 122 ! 122 IF( kt == nit000 ) THEN 123 IF(lwp) WRITE(numout,*) 124 IF(lwp) WRITE(numout,*) 'tra_qsr : penetration of the surface solar radiation' 125 IF(lwp) WRITE(numout,*) '~~~~~~~' 123 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 124 IF( kt == nit000 ) THEN 125 IF(lwp) WRITE(numout,*) 126 IF(lwp) WRITE(numout,*) 'tra_qsr : penetration of the surface solar radiation' 127 IF(lwp) WRITE(numout,*) '~~~~~~~' 128 ENDIF 126 129 ENDIF 127 130 ! 128 131 IF( l_trdtra ) THEN ! trends diagnostic: save the input temperature trend 129 ALLOCATE( ztrdt(jpi,jpj,jpk) ) 132 ALLOCATE( ztrdt(jpi,jpj,jpk) ) 130 133 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) 131 134 ENDIF … … 134 137 ! ! before qsr induced heat content ! 135 138 ! !-----------------------------------! 139 IF( ntsi == Nis0 ) THEN ; isi = nn_hls ; ELSE ; isi = 0 ; ENDIF ! Avoid double-counting when using tiling 140 IF( ntsj == Njs0 ) THEN ; isj = nn_hls ; ELSE ; isj = 0 ; ENDIF 141 IF( ntei == Nie0 ) THEN ; iei = nn_hls ; ELSE ; iei = 0 ; ENDIF 142 IF( ntej == Nje0 ) THEN ; iej = nn_hls ; ELSE ; iej = 0 ; ENDIF 143 136 144 IF( kt == nit000 ) THEN !== 1st time step ==! 137 IF( ln_rstart .AND. iom_varid( numror, 'qsr_hc_b', ldstop = .FALSE. ) > 0 .AND. .NOT.l_1st_euler ) THEN ! read in restart 138 IF(lwp) WRITE(numout,*) ' nit000-1 qsr tracer content forcing field read in the restart file' 145 IF( ln_rstart .AND. .NOT.l_1st_euler ) THEN ! read in restart 139 146 z1_2 = 0.5_wp 140 CALL iom_get( numror, jpdom_auto, 'qsr_hc_b', qsr_hc_b, ldxios = lrxios ) ! before heat content trend due to Qsr flux 141 ELSE ! No restart or restart not found: Euler forward time stepping 147 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 148 IF(lwp) WRITE(numout,*) ' nit000-1 qsr tracer content forcing field read in the restart file' 149 CALL iom_get( numror, jpdom_auto, 'qsr_hc_b', qsr_hc_b ) ! before heat content trend due to Qsr flux 150 ENDIF 151 ELSE ! No restart or Euler forward at 1st time step 142 152 z1_2 = 1._wp 143 qsr_hc_b(:,:,:) = 0._wp 153 DO_3D( isi, iei, isj, iej, 1, jpk ) 154 qsr_hc_b(ji,jj,jk) = 0._wp 155 END_3D 144 156 ENDIF 145 157 ELSE !== Swap of qsr heat content ==! 146 158 z1_2 = 0.5_wp 147 qsr_hc_b(:,:,:) = qsr_hc(:,:,:) 159 DO_3D( isi, iei, isj, iej, 1, jpk ) 160 qsr_hc_b(ji,jj,jk) = qsr_hc(ji,jj,jk) 161 END_3D 148 162 ENDIF 149 163 ! … … 154 168 CASE( np_BIO ) !== bio-model fluxes ==! 155 169 ! 156 DO jk = 1, nksr157 qsr_hc( :,:,jk) = r1_rho0_rcp * ( etot3(:,:,jk) - etot3(:,:,jk+1) )158 END DO170 DO_3D( isi, iei, isj, iej, 1, nksr ) 171 qsr_hc(ji,jj,jk) = r1_rho0_rcp * ( etot3(ji,jj,jk) - etot3(ji,jj,jk+1) ) 172 END_3D 159 173 ! 160 174 CASE( np_RGB , np_RGBc ) !== R-G-B fluxes ==! 161 175 ! 162 ALLOCATE( ze0 ( jpi,jpj) , ze1 (jpi,jpj) , &163 & ze2 ( jpi,jpj) , ze3 (jpi,jpj) , &164 & ztmp3d( jpi,jpj,nksr + 1) )176 ALLOCATE( ze0 (A2D(nn_hls)) , ze1 (A2D(nn_hls)) , & 177 & ze2 (A2D(nn_hls)) , ze3 (A2D(nn_hls)) , & 178 & ztmp3d(A2D(nn_hls),nksr + 1) ) 165 179 ! 166 180 IF( nqsr == np_RGBc ) THEN !* Variable Chlorophyll 167 CALL fld_read( kt, 1, sf_chl ) ! Read Chl data and provides it at the current time step 181 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only for the full domain 182 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 183 CALL fld_read( kt, 1, sf_chl ) ! Read Chl data and provides it at the current time step 184 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 1 ) ! Revert to tile domain 185 ENDIF 168 186 ! 169 187 ! Separation in R-G-B depending on the surface Chl … … 172 190 ! most expensive calculations) 173 191 ! 174 DO_2D( 0, 0, 0, 0)192 DO_2D( isi, iei, isj, iej ) 175 193 ! zlogc = log(zchl) 176 zlogc = LOG ( MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) ) ) 194 zlogc = LOG ( MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) ) ) 177 195 ! zc1 : log(zCze) = log (1.12 * zchl**0.803) 178 zc1 = 0.113328685307 + 0.803 * zlogc 196 zc1 = 0.113328685307 + 0.803 * zlogc 179 197 ! zc2 : log(zCtot) = log(40.6 * zchl**0.459) 180 zc2 = 3.703768066608 + 0.459 * zlogc 198 zc2 = 3.703768066608 + 0.459 * zlogc 181 199 ! zc3 : log(zze) = log(568.2 * zCtot**(-0.746)) 182 zc3 = 6.34247346942 - 0.746 * zc2 200 zc3 = 6.34247346942 - 0.746 * zc2 183 201 ! IF( log(zze) > log(102.) ) log(zze) = log(200.0 * zCtot**(-0.293)) 184 IF( zc3 > 4.62497281328 ) zc3 = 5.298317366548 - 0.293 * zc2 185 ! 202 IF( zc3 > 4.62497281328 ) zc3 = 5.298317366548 - 0.293 * zc2 203 ! 186 204 ze0(ji,jj) = zlogc ! ze0 = log(zchl) 187 205 ze1(ji,jj) = EXP( zc1 ) ! ze1 = zCze … … 189 207 ze3(ji,jj) = EXP( - zc3 ) ! ze3 = 1/zze 190 208 END_2D 191 209 192 210 ! 193 DO_3D( 0, 0, 0, 0, 1, nksr + 1 )211 DO_3D( isi, iei, isj, iej, 1, nksr + 1 ) 194 212 ! zchl = ALOG( ze0(ji,jj) ) 195 213 zlogc = ze0(ji,jj) … … 211 229 ELSE !* constant chlorophyll 212 230 zchl = 0.05 213 ! NB. make sure constant value is such that: 231 ! NB. make sure constant value is such that: 214 232 zchl = MIN( 10. , MAX( 0.03, zchl ) ) 215 233 ! Convert chlorophyll value to attenuation coefficient look-up table index 216 234 zlui = 41 + 20.*LOG10(zchl) + 1.e-15 217 235 DO jk = 1, nksr + 1 218 ztmp3d(:,:,jk) = zlui 236 ztmp3d(:,:,jk) = zlui 219 237 END DO 220 238 ENDIF 221 239 ! 222 240 zcoef = ( 1. - rn_abs ) / 3._wp !* surface equi-partition in R-G-B 223 DO_2D( 0, 0, 0, 0)241 DO_2D( isi, iei, isj, iej ) 224 242 ze0(ji,jj) = rn_abs * qsr(ji,jj) 225 243 ze1(ji,jj) = zcoef * qsr(ji,jj) 226 244 ze2(ji,jj) = zcoef * qsr(ji,jj) 227 245 ze3(ji,jj) = zcoef * qsr(ji,jj) 228 ! store the surface SW radiation; re-use the surface ztmp3d array 246 ! store the surface SW radiation; re-use the surface ztmp3d array 229 247 ! since the surface attenuation coefficient is not used 230 248 ztmp3d(ji,jj,1) = qsr(ji,jj) … … 232 250 ! 233 251 ! !* interior equi-partition in R-G-B depending on vertical profile of Chl 234 DO_3D( 0, 0, 0, 0, 2, nksr + 1 )252 DO_3D( isi, iei, isj, iej, 2, nksr + 1 ) 235 253 ze3t = e3t(ji,jj,jk-1,Kmm) 236 254 irgb = NINT( ztmp3d(ji,jj,jk) ) … … 246 264 END_3D 247 265 ! 248 DO_3D( 0, 0, 0, 0, 1, nksr ) !* now qsr induced heat content266 DO_3D( isi, iei, isj, iej, 1, nksr ) !* now qsr induced heat content 249 267 qsr_hc(ji,jj,jk) = r1_rho0_rcp * ( ztmp3d(ji,jj,jk) - ztmp3d(ji,jj,jk+1) ) 250 268 END_3D 251 269 ! 252 DEALLOCATE( ze0 , ze1 , ze2 , ze3 , ztmp3d ) 270 DEALLOCATE( ze0 , ze1 , ze2 , ze3 , ztmp3d ) 253 271 ! 254 272 CASE( np_2BD ) !== 2-bands fluxes ==! … … 256 274 zz0 = rn_abs * r1_rho0_rcp ! surface equi-partition in 2-bands 257 275 zz1 = ( 1. - rn_abs ) * r1_rho0_rcp 258 DO_3D( 0, 0, 0, 0, 1, nksr ) ! solar heat absorbed at T-point in the top 400m276 DO_3D( isi, iei, isj, iej, 1, nksr ) !* now qsr induced heat content 259 277 zc0 = zz0 * EXP( -gdepw(ji,jj,jk ,Kmm)*xsi0r ) + zz1 * EXP( -gdepw(ji,jj,jk ,Kmm)*xsi1r ) 260 278 zc1 = zz0 * EXP( -gdepw(ji,jj,jk+1,Kmm)*xsi0r ) + zz1 * EXP( -gdepw(ji,jj,jk+1,Kmm)*xsi1r ) 261 qsr_hc(ji,jj,jk) = qsr(ji,jj) * ( zc0 * wmask(ji,jj,jk) - zc1 * wmask(ji,jj,jk+1) ) 279 qsr_hc(ji,jj,jk) = qsr(ji,jj) * ( zc0 * wmask(ji,jj,jk) - zc1 * wmask(ji,jj,jk+1) ) 262 280 END_3D 263 281 ! … … 274 292 ! 275 293 ! sea-ice: store the 1st ocean level attenuation coefficient 276 DO_2D( 0, 0, 0, 0)294 DO_2D( isi, iei, isj, iej ) 277 295 IF( qsr(ji,jj) /= 0._wp ) THEN ; fraqsr_1lev(ji,jj) = qsr_hc(ji,jj,1) / ( r1_rho0_rcp * qsr(ji,jj) ) 278 296 ELSE ; fraqsr_1lev(ji,jj) = 1._wp 279 297 ENDIF 280 298 END_2D 281 CALL lbc_lnk( 'traqsr', fraqsr_1lev(:,:), 'T', 1._wp ) 282 ! 283 IF( iom_use('qsr3d') ) THEN ! output the shortwave Radiation distribution 284 ALLOCATE( zetot(jpi,jpj,jpk) ) 285 zetot(:,:,nksr+1:jpk) = 0._wp ! below ~400m set to zero 286 DO jk = nksr, 1, -1 287 zetot(:,:,jk) = zetot(:,:,jk+1) + qsr_hc(:,:,jk) * rho0_rcp 288 END DO 289 CALL iom_put( 'qsr3d', zetot ) ! 3D distribution of shortwave Radiation 290 DEALLOCATE( zetot ) 291 ENDIF 292 ! 293 IF( lrst_oce ) THEN ! write in the ocean restart file 294 IF( lwxios ) CALL iom_swap( cwxios_context ) 295 CALL iom_rstput( kt, nitrst, numrow, 'qsr_hc_b' , qsr_hc , ldxios = lwxios ) 296 CALL iom_rstput( kt, nitrst, numrow, 'fraqsr_1lev', fraqsr_1lev, ldxios = lwxios ) 297 IF( lwxios ) CALL iom_swap( cxios_context ) 299 ! 300 ! TEMP: [tiling] This change not necessary and working array can use A2D(nn_hls) if using XIOS (subdomain support) 301 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 302 IF( iom_use('qsr3d') ) THEN ! output the shortwave Radiation distribution 303 ALLOCATE( zetot(jpi,jpj,jpk) ) 304 zetot(:,:,nksr+1:jpk) = 0._wp ! below ~400m set to zero 305 DO jk = nksr, 1, -1 306 zetot(:,:,jk) = zetot(:,:,jk+1) + qsr_hc(:,:,jk) * rho0_rcp 307 END DO 308 CALL iom_put( 'qsr3d', zetot ) ! 3D distribution of shortwave Radiation 309 DEALLOCATE( zetot ) 310 ENDIF 311 ENDIF 312 ! 313 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 314 IF( lrst_oce ) THEN ! write in the ocean restart file 315 CALL iom_rstput( kt, nitrst, numrow, 'qsr_hc_b' , qsr_hc ) 316 CALL iom_rstput( kt, nitrst, numrow, 'fraqsr_1lev', fraqsr_1lev ) 317 ENDIF 298 318 ENDIF 299 319 ! … … 301 321 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) - ztrdt(:,:,:) 302 322 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_qsr, ztrdt ) 303 DEALLOCATE( ztrdt ) 323 DEALLOCATE( ztrdt ) 304 324 ENDIF 305 325 ! ! print mean trends (used for debugging) … … 320 340 !! from two length scale of penetration (rn_si0,rn_si1) and a ratio 321 341 !! (rn_abs). These parameters are read in the namtra_qsr namelist. The 322 !! default values correspond to clear water (type I in Jerlov' 342 !! default values correspond to clear water (type I in Jerlov' 323 343 !! (1968) classification. 324 344 !! called by tra_qsr at the first timestep (nit000) … … 370 390 & ' 2 bands, 3 RGB bands or bio-model light penetration' ) 371 391 ! 372 IF( ln_qsr_rgb .AND. nn_chldta == 0 ) nqsr = np_RGB 392 IF( ln_qsr_rgb .AND. nn_chldta == 0 ) nqsr = np_RGB 373 393 IF( ln_qsr_rgb .AND. nn_chldta == 1 ) nqsr = np_RGBc 374 394 IF( ln_qsr_2bd ) nqsr = np_2BD … … 380 400 ! 381 401 SELECT CASE( nqsr ) 382 ! 402 ! 383 403 CASE( np_RGB , np_RGBc ) !== Red-Green-Blue light penetration ==! 384 ! 404 ! 385 405 IF(lwp) WRITE(numout,*) ' ==>>> R-G-B light penetration ' 386 406 ! 387 407 CALL trc_oce_rgb( rkrgb ) ! tabulated attenuation coef. 388 ! 408 ! 389 409 nksr = trc_oce_ext_lev( r_si2, 33._wp ) ! level of light extinction 390 410 ! … … 420 440 ! 421 441 CALL trc_oce_rgb( rkrgb ) ! tabulated attenuation coef. 422 ! 442 ! 423 443 nksr = trc_oce_ext_lev( r_si2, 33._wp ) ! level of light extinction 424 444 ! … … 431 451 ! 1st ocean level attenuation coefficient (used in sbcssm) 432 452 IF( iom_varid( numror, 'fraqsr_1lev', ldstop = .FALSE. ) > 0 ) THEN 433 CALL iom_get( numror, jpdom_auto, 'fraqsr_1lev' , fraqsr_1lev , ldxios = lrxios)453 CALL iom_get( numror, jpdom_auto, 'fraqsr_1lev' , fraqsr_1lev ) 434 454 ELSE 435 455 fraqsr_1lev(:,:) = 1._wp ! default : no penetration 436 456 ENDIF 437 457 ! 438 IF( lwxios ) THEN439 CALL iom_set_rstw_var_active('qsr_hc_b')440 CALL iom_set_rstw_var_active('fraqsr_1lev')441 ENDIF442 !443 458 END SUBROUTINE tra_qsr_init 444 459
Note: See TracChangeset
for help on using the changeset viewer.
