- Timestamp:
- 2020-07-02T11:27:33+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/TRA/traqsr.F90
r12724 r13217 111 111 REAL(wp) :: zc0 , zc1 , zc2 , zc3 ! - - 112 112 REAL(wp) :: zzc0, zzc1, zzc2, zzc3 ! - - 113 REAL(wp) :: zz0 , zz1 ! - - 114 REAL(wp) :: zCb, zCmax, zze, zpsi, zpsimax, zdelpsi, zCtot, zCze 115 REAL(wp) :: zlogc, zlogc2, zlogc3 116 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zekb, zekg, zekr 117 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ze0, ze1, ze2, ze3, zea, ztrdt 118 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zetot, zchl3d 113 REAL(wp) :: zz0 , zz1 , ze3t, zlui ! - - 114 REAL(wp) :: zCb, zCmax, zpsi, zpsimax, zrdpsi, zCze 115 REAL(wp) :: zlogc, zlogze, zlogCtot, zlogCze 116 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: ze0, ze1, ze2, ze3 117 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdt, zetot, ztmp3d 119 118 !!---------------------------------------------------------------------- 120 119 ! … … 161 160 CASE( np_RGB , np_RGBc ) !== R-G-B fluxes ==! 162 161 ! 163 ALLOCATE( ze kb(jpi,jpj) , zekg(jpi,jpj) , zekr (jpi,jpj) ,&164 & ze 0 (jpi,jpj,jpk) , ze1 (jpi,jpj,jpk) , ze2 (jpi,jpj,jpk) ,&165 & z e3 (jpi,jpj,jpk) , zea (jpi,jpj,jpk) , zchl3d(jpi,jpj,jpk))162 ALLOCATE( ze0 (jpi,jpj) , ze1 (jpi,jpj) , & 163 & ze2 (jpi,jpj) , ze3 (jpi,jpj) , & 164 & ztmp3d(jpi,jpj,nksr + 1) ) 166 165 ! 167 166 IF( nqsr == np_RGBc ) THEN !* Variable Chlorophyll 168 167 CALL fld_read( kt, 1, sf_chl ) ! Read Chl data and provides it at the current time step 168 ! 169 ! Separation in R-G-B depending on the surface Chl 170 ! perform and store as many of the 2D calculations as possible 171 ! before the 3D loop (use the temporary 2D arrays to replace the 172 ! most expensive calculations) 173 ! 174 DO_2D_00_00 175 ! zlogc = log(zchl) 176 zlogc = LOG ( MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) ) ) 177 ! zc1 : log(zCze) = log (1.12 * zchl**0.803) 178 zc1 = 0.113328685307 + 0.803 * zlogc 179 ! zc2 : log(zCtot) = log(40.6 * zchl**0.459) 180 zc2 = 3.703768066608 + 0.459 * zlogc 181 ! zc3 : log(zze) = log(568.2 * zCtot**(-0.746)) 182 zc3 = 6.34247346942 - 0.746 * zc2 183 ! 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 ! 186 ze0(ji,jj) = zlogc ! ze0 = log(zchl) 187 ze1(ji,jj) = EXP( zc1 ) ! ze1 = zCze 188 ze2(ji,jj) = 1._wp / ( 0.710 + zlogc * ( 0.159 + zlogc * 0.021 ) ) ! ze2 = 1/zdelpsi 189 ze3(ji,jj) = EXP( - zc3 ) ! ze3 = 1/zze 190 END_2D 191 192 ! 193 DO_3D_00_00 ( 1, nksr + 1 ) 194 ! zchl = ALOG( ze0(ji,jj) ) 195 zlogc = ze0(ji,jj) 196 ! 197 zCb = 0.768 + zlogc * ( 0.087 - zlogc * ( 0.179 + zlogc * 0.025 ) ) 198 zCmax = 0.299 - zlogc * ( 0.289 - zlogc * 0.579 ) 199 zpsimax = 0.6 - zlogc * ( 0.640 - zlogc * ( 0.021 + zlogc * 0.115 ) ) 200 ! zdelpsi = 0.710 + zlogc * ( 0.159 + zlogc * 0.021 ) 201 ! 202 zCze = ze1(ji,jj) 203 zrdpsi = ze2(ji,jj) ! 1/zdelpsi 204 zpsi = ze3(ji,jj) * gdepw(ji,jj,jk,Kmm) ! gdepw/zze 205 ! 206 ! NB. make sure zchl value is such that: zchl = MIN( 10. , MAX( 0.03, zchl ) ) 207 zchl = MIN( 10. , MAX( 0.03, zCze * ( zCb + zCmax * EXP( -( (zpsi - zpsimax) * zrdpsi )**2 ) ) ) ) 208 ! Convert chlorophyll value to attenuation coefficient look-up table index 209 ztmp3d(ji,jj,jk) = 41 + 20.*LOG10(zchl) + 1.e-15 210 END_3D 211 ELSE !* constant chlorophyll 212 zchl = 0.05 213 ! NB. make sure constant value is such that: 214 zchl = MIN( 10. , MAX( 0.03, zchl ) ) 215 ! Convert chlorophyll value to attenuation coefficient look-up table index 216 zlui = 41 + 20.*LOG10(zchl) + 1.e-15 169 217 DO jk = 1, nksr + 1 170 DO jj = 2, jpjm1 ! Separation in R-G-B depending of the surface Chl 171 DO ji = 2, jpim1 172 zchl = MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) ) 173 zCtot = 40.6 * zchl**0.459 174 zze = 568.2 * zCtot**(-0.746) 175 IF( zze > 102. ) zze = 200.0 * zCtot**(-0.293) 176 zpsi = gdepw(ji,jj,jk,Kmm) / zze 177 ! 178 zlogc = LOG( zchl ) 179 zlogc2 = zlogc * zlogc 180 zlogc3 = zlogc * zlogc * zlogc 181 zCb = 0.768 + 0.087 * zlogc - 0.179 * zlogc2 - 0.025 * zlogc3 182 zCmax = 0.299 - 0.289 * zlogc + 0.579 * zlogc2 183 zpsimax = 0.6 - 0.640 * zlogc + 0.021 * zlogc2 + 0.115 * zlogc3 184 zdelpsi = 0.710 + 0.159 * zlogc + 0.021 * zlogc2 185 zCze = 1.12 * (zchl)**0.803 186 ! 187 zchl3d(ji,jj,jk) = zCze * ( zCb + zCmax * EXP( -( (zpsi - zpsimax) / zdelpsi )**2 ) ) 188 END DO 189 ! 190 END DO 218 ztmp3d(:,:,jk) = zlui 191 219 END DO 192 ELSE !* constant chrlorophyll193 DO jk = 1, nksr + 1194 zchl3d(:,:,jk) = 0.05195 ENDDO196 220 ENDIF 197 221 ! 198 222 zcoef = ( 1. - rn_abs ) / 3._wp !* surface equi-partition in R-G-B 199 223 DO_2D_00_00 200 ze0(ji,jj,1) = rn_abs * qsr(ji,jj) 201 ze1(ji,jj,1) = zcoef * qsr(ji,jj) 202 ze2(ji,jj,1) = zcoef * qsr(ji,jj) 203 ze3(ji,jj,1) = zcoef * qsr(ji,jj) 204 zea(ji,jj,1) = qsr(ji,jj) 224 ze0(ji,jj) = rn_abs * qsr(ji,jj) 225 ze1(ji,jj) = zcoef * qsr(ji,jj) 226 ze2(ji,jj) = zcoef * qsr(ji,jj) 227 ze3(ji,jj) = zcoef * qsr(ji,jj) 228 ! store the surface SW radiation; re-use the surface ztmp3d array 229 ! since the surface attenuation coefficient is not used 230 ztmp3d(ji,jj,1) = qsr(ji,jj) 205 231 END_2D 206 232 ! 207 DO jk = 2, nksr+1 !* interior equi-partition in R-G-B depending of vertical profile of Chl 208 DO_2D_00_00 209 zchl = MIN( 10. , MAX( 0.03, zchl3d(ji,jj,jk) ) ) 210 irgb = NINT( 41 + 20.*LOG10(zchl) + 1.e-15 ) 211 zekb(ji,jj) = rkrgb(1,irgb) 212 zekg(ji,jj) = rkrgb(2,irgb) 213 zekr(ji,jj) = rkrgb(3,irgb) 214 END_2D 215 216 DO_2D_00_00 217 zc0 = ze0(ji,jj,jk-1) * EXP( - e3t(ji,jj,jk-1,Kmm) * xsi0r ) 218 zc1 = ze1(ji,jj,jk-1) * EXP( - e3t(ji,jj,jk-1,Kmm) * zekb(ji,jj) ) 219 zc2 = ze2(ji,jj,jk-1) * EXP( - e3t(ji,jj,jk-1,Kmm) * zekg(ji,jj) ) 220 zc3 = ze3(ji,jj,jk-1) * EXP( - e3t(ji,jj,jk-1,Kmm) * zekr(ji,jj) ) 221 ze0(ji,jj,jk) = zc0 222 ze1(ji,jj,jk) = zc1 223 ze2(ji,jj,jk) = zc2 224 ze3(ji,jj,jk) = zc3 225 zea(ji,jj,jk) = ( zc0 + zc1 + zc2 + zc3 ) * wmask(ji,jj,jk) 226 END_2D 227 END DO 233 !* interior equi-partition in R-G-B depending on vertical profile of Chl 234 DO_3D_00_00 ( 2, nksr + 1 ) 235 ze3t = e3t(ji,jj,jk-1,Kmm) 236 irgb = NINT( ztmp3d(ji,jj,jk) ) 237 zc0 = ze0(ji,jj) * EXP( - ze3t * xsi0r ) 238 zc1 = ze1(ji,jj) * EXP( - ze3t * rkrgb(1,irgb) ) 239 zc2 = ze2(ji,jj) * EXP( - ze3t * rkrgb(2,irgb) ) 240 zc3 = ze3(ji,jj) * EXP( - ze3t * rkrgb(3,irgb) ) 241 ze0(ji,jj) = zc0 242 ze1(ji,jj) = zc1 243 ze2(ji,jj) = zc2 244 ze3(ji,jj) = zc3 245 ztmp3d(ji,jj,jk) = ( zc0 + zc1 + zc2 + zc3 ) * wmask(ji,jj,jk) 246 END_3D 228 247 ! 229 248 DO_3D_00_00( 1, nksr ) 230 qsr_hc(ji,jj,jk) = r1_rho0_rcp * ( z ea(ji,jj,jk) - zea(ji,jj,jk+1) )249 qsr_hc(ji,jj,jk) = r1_rho0_rcp * ( ztmp3d(ji,jj,jk) - ztmp3d(ji,jj,jk+1) ) 231 250 END_3D 232 251 ! 233 DEALLOCATE( ze kb , zekg , zekr , ze0 , ze1 , ze2 , ze3 , zea , zchl3d )252 DEALLOCATE( ze0 , ze1 , ze2 , ze3 , ztmp3d ) 234 253 ! 235 254 CASE( np_2BD ) !== 2-bands fluxes ==!
Note: See TracChangeset
for help on using the changeset viewer.