1! ====================================================================================================== 2! This kernel represents a distillation of *part* of 3! the taumol03 calculation in the gas optics part of the PSRAD 4! atmospheric 5! radiation code. 6! 7! It is meant to show conceptually how one might "SIMD-ize" swaths of 8! the taumol03 code related to calculating the 9! taug term, so that the impact of the conditional expression on 10! specparm could be reduced and at least partial vectorization 11! across columns could be achieved. 12! 13! I consider it at this point to be "compiling pseudo-code". 14! 15! By this I mean that the code as written compiles under ifort, but has 16! not been tested 17! for correctness, nor I have written a driver routine for it. It does 18! not contain everything 19! that is going on in the taug parent taumol03 code, but I don't claim 20! to actually completely 21! understand the physical meaning of all or even most of the inputs 22! required to make it run. 23! 24! It has been written to vectorize, but apparently does not actually do 25! that 26! under the ifort V13 compiler with the -xHost -O3 level of 27! optimization, even with !dir$ assume_aligned directives. 28! I hypothesize that the compiler is baulking to do so for the indirect 29! addressed calls into the absa 30! look-up table, either that or 4 byte integers may be causing alignment 31! issues relative to 8 byte reals. Anyway, 32! it seems to complain about the key loop being too complex. 33! ====================================================================================================== 34MODULE mo_taumol03 35 USE mo_kind, only:wp 36 USE mo_lrtm_setup, ONLY: nspa 37 USE mo_lrtm_setup, ONLY: nspb 38 USE rrlw_planck, ONLY: chi_mls 39 USE rrlw_kg03, ONLY: selfref 40 USE rrlw_kg03, ONLY: forref 41 USE rrlw_kg03, ONLY: ka_mn2o 42 USE rrlw_kg03, ONLY: absa 43 USE rrlw_kg03, ONLY: fracrefa 44 USE rrlw_kg03, ONLY: kb_mn2o 45 USE rrlw_kg03, ONLY: absb 46 USE rrlw_kg03, ONLY: fracrefb 47 IMPLICIT NONE 48 PRIVATE 49 PUBLIC taumol03_lwr,taumol03_upr 50 CONTAINS 51 SUBROUTINE taumol03_lwr(ncol, startCol, laytrop, nlayers, & 52 rat_h2oco2, colco2, colh2o, coln2o, coldry, & 53 fac0, fac1, minorfrac, & 54 selffac,selffrac,forfac,forfrac, & 55 jp, jt, ig, indself, & 56 indfor, indminor, & 57 taug, fracs) 58 IMPLICIT NONE 59 60 real(kind=wp), PARAMETER :: oneminus = 1.0_wp - 1.0e-06_wp 61 integer, intent(in) :: ncol ! number of simd columns 62 integer, intent(in) :: startCol !starting index column 63 integer, intent(in) :: laytrop ! number of layers forwer atmosphere kernel 64 integer, intent(in) :: nlayers ! total number of layers 65 real(kind=wp), intent(in), dimension(startCol:ncol,0:1,nlayers) :: rat_h2oco2,fac0,fac1 ! not sure of ncol depend 66 real(kind=wp), intent(in), dimension(startCol:ncol,nlayers) :: colco2,colh2o,coln2o,coldry ! these appear to be gas concentrations 67 68 real(kind=wp), intent(in), dimension(startCol:ncol,nlayers) :: selffac,selffrac ! not sure of ncol depend 69 real(kind=wp), intent(in), dimension(startCol:ncol,nlayers) :: forfac,forfrac ! not sure of ncol depend 70 real(kind=wp), intent(in), dimension(startCol:ncol,nlayers) :: minorfrac ! not sure of ncol depend 71 72 ! Look up tables and related lookup indices 73 ! I assume all lookup indices depend on 3D position 74 ! ================================================= 75 76 integer, intent(in) :: jp(startCol:ncol,nlayers) ! I assume jp depends on ncol 77 integer, intent(in) :: jt(startCol:ncol,0:1,nlayers) ! likewise for jt 78 integer, intent(in) :: ig ! ig indexes into lookup tables 79 integer, intent(in) :: indself(startCol:ncol,nlayers) ! self index array 80 integer, intent(in) :: indfor(startCol:ncol,nlayers) ! for index array 81 integer, intent(in) :: indminor(startCol:ncol,nlayers) ! ka_mn2o index array 82 real(kind=wp), intent(out), dimension(startCol:ncol,nlayers) :: taug ! kernel result 83 real(kind=wp), intent(out), dimension(startCol:ncol,nlayers) :: fracs ! kernel result 84 85 ! Local variable 86 ! ============== 87 88 integer :: lay ! layer index 89 integer :: i ! specparm types index 90 integer :: icol ! column index 91 ! vector temporaries 92 ! ==================== 93 94 integer, dimension(1:3,1:3) :: caseTypeOperations 95 integer, dimension(startCol:ncol) :: caseType 96 real(kind=wp), dimension(startCol:ncol) :: p, p4, fs 97 real(kind=wp), dimension(startCol:ncol) :: fmn2o,fmn2omf 98 real(kind=wp), dimension(startCol:ncol) :: fpl 99 real(kind=wp), dimension(startCol:ncol) :: specmult, speccomb, specparm 100 real(kind=wp), dimension(startCol:ncol) :: specmult_mn2o, speccomb_mn2o,specparm_mn2o 101 real(kind=wp), dimension(startCol:ncol) :: specmult_planck, speccomb_planck,specparm_planck 102 real(kind=wp), dimension(startCol:ncol) :: n2om1,n2om2,absn2o,adjcoln2o,adjfac,chi_n2o,ratn2o 103 real(kind=wp), dimension(startCol:ncol,0:1) :: tau_major 104 real(kind=wp), dimension(startCol:ncol) :: taufor,tauself 105 integer, dimension(startCol:ncol) :: js, ind0, ind00, ind01, ind02, jmn2o, jpl 106 107 ! Register temporaries 108 ! ==================== 109 110 real(kind=wp) :: p2,fk0,fk1,fk2 111 real(kind=wp) :: refrat_planck_a, refrat_planck_b 112 real(kind=wp) :: refrat_m_a, refrat_m_b 113 integer :: rrpk_counter=0 114 115 !dir$ assume_aligned jp:64 116 !dir$ assume_aligned jt:64 117 !dir$ assume_aligned rat_h2oco2:64 118 !dir$ assume_aligned colco2:64 119 !dir$ assume_aligned colh2o:64 120 !dir$ assume_aligned fac0:64 121 !dir$ assume_aligned fac1:64 122 !dir$ assume_aligned taug:64 123 124 !dir$ assume_aligned p:64 125 !dir$ assume_aligned p4:64 126 !dir$ assume_aligned specmult:64 127 !dir$ assume_aligned speccomb:64 128 !dir$ assume_aligned specparm:64 129 !dir$ assume_aligned specmult_mn2o:64 130 !dir$ assume_aligned speccomb_mn2o:64 131 !dir$ assume_aligned specparm_mn2o:64 132 !dir$ assume_aligned specmult_planck:64 133 !dir$ assume_aligned speccomb_planck:64 134 !dir$ assume_aligned specparm_planck:64 135 !dir$ assume_aligned indself:64 136 !dir$ assume_aligned indfor:64 137 !dir$ assume_aligned indminor:64 138 !dir$ assume_aligned fs:64 139 !dir$ assume_aligned tau_major:64 140 141 !dir$ assume_aligned js:64 142 !dir$ assume_aligned ind0:64 143 !dir$ assume_aligned ind00:64 144 !dir$ assume_aligned ind01:64 145 !dir$ assume_aligned ind02:64 146 147 !dir$ assume_aligned caseTypeOperations:64 148 !dir$ assume_aligned caseType:64 149 150 ! Initialize Case type operations 151 !================================= 152 153 caseTypeOperations(1:3,1) = (/0, 1, 2/) 154 caseTypeOperations(1:3,2) = (/1, 0,-1/) 155 caseTypeOperations(1:3,3) = (/0, 1, 1/) 156 157 ! Minor gas mapping levels: 158 ! lower - n2o, p = 706.272 mbar, t = 278.94 k 159 ! upper - n2o, p = 95.58 mbar, t = 215.7 k 160 161 ! P = 212.725 mb 162 refrat_planck_a = chi_mls(1,9)/chi_mls(2,9) 163 164 ! P = 95.58 mb 165 refrat_planck_b = chi_mls(1,13)/chi_mls(2,13) 166 167 ! P = 706.270mb 168 refrat_m_a = chi_mls(1,3)/chi_mls(2,3) 169 170 ! P = 95.58 mb 171 refrat_m_b = chi_mls(1,13)/chi_mls(2,13) 172 173 ! Lower atmosphere loop 174 ! ===================== 175 176 DO lay = 1,laytrop ! loop over layers 177 178 ! Compute tau_major term 179 ! ====================== 180 181 DO i=0,1 ! loop over specparm types 182 183 ! This loop should vectorize 184 ! ============================= 185 !dir$ SIMD 186 DO icol=startCol,ncol ! Vectorizes with dir - 14.0.2 187 speccomb(icol) = colh2o(icol,lay) + rat_h2oco2(icol,i,lay)*colco2(icol,lay) 188 specparm(icol) = colh2o(icol,lay)/speccomb(icol) 189 IF (specparm(icol) .GE. oneminus) specparm(icol) = oneminus 190 specmult(icol) = 8._wp*(specparm(icol)) 191 js(icol) = 1 + INT(specmult(icol)) 192 fs(icol) = MOD(specmult(icol),1.0_wp) 193 END DO 194 195 ! The only conditional loop 196 ! ========================= 197 198 DO icol=startCol,ncol ! Vectorizes as is 14.0.2 199 IF (specparm(icol) .LT. 0.125_wp) THEN 200 caseType(icol)=1 201 p(icol) = fs(icol) - 1.0_wp 202 p2 = p(icol)*p(icol) 203 p4(icol) = p2*p2 204 ELSE IF (specparm(icol) .GT. 0.875_wp) THEN 205 caseType(icol)=2 206 p(icol) = -fs(icol) 207 p2 = p(icol)*p(icol) 208 p4(icol) = p2*p2 209 ELSE 210 caseType(icol)=3 211 ! SIMD way of writing fk0= 1-fs and fk1 = fs, fk2=zero 212 ! =========================================================== 213 214 p4(icol) = 1.0_wp - fs(icol) 215 p(icol) = -p4(icol) ! complicated way of getting fk2 = zero for ELSE case 216 ENDIF 217 END DO 218 219 ! Vector/SIMD index loop calculation 220 ! ================================== 221 222 !dir$ SIMD 223 DO icol=startCol,ncol ! Vectorizes with dir 14.0.2 224 ind0(icol) = ((jp(icol,lay)-(1*MOD(i+1,2)))*5+(jt(icol,i,lay)-1))*nspa(3) +js(icol) 225 ind00(icol) = ind0(icol) + caseTypeOperations(1,caseType(icol)) 226 ind01(icol) = ind0(icol) + caseTypeOperations(2,caseType(icol)) 227 ind02(icol) = ind0(icol) + caseTypeOperations(3,caseType(icol)) 228 END DO 229 230 ! What we've been aiming for a nice flop intensive 231 ! SIMD/vectorizable loop! 232 ! 17 flops 233 ! 234 ! Albeit at the cost of a couple extra flops for the fk2 term 235 ! and a repeated lookup table access for the fk2 term in the 236 ! the ELSE case when specparm or specparm1 is (> .125 && < .875) 237 ! =============================================================== 238 239 !dir$ SIMD 240 DO icol=startCol,ncol ! Vectorizes with dir 14.0.2 241 242 fk0 = p4(icol) 243 fk1 = 1.0_wp - p(icol) - 2.0_wp*p4(icol) 244 fk2 = p(icol) + p4(icol) 245 tau_major(icol,i) = speccomb(icol) * ( & 246 fac0(icol,i,lay)*(fk0*absa(ind00(icol),ig) + & 247 fk1*absa(ind01(icol),ig) + & 248 fk2*absa(ind02(icol),ig)) + & 249 fac1(icol,i,lay)*(fk0*absa(ind00(icol)+9,ig) + & 250 fk1*absa(ind01(icol)+9,ig) + & 251 fk2*absa(ind02(icol)+9,ig))) 252 END DO 253 254 END DO ! end loop over specparm types for tau_major calculation 255 256 ! Compute taufor and tauself terms: 257 ! Note the use of 1D bilinear interpolation of selfref and forref 258 ! lookup table values 259 ! =================================================================================== 260 261 !dir$ SIMD 262 DO icol=startCol,ncol ! Vectorizes with dir 14.0.2 263 tauself(icol) = selffac(icol,lay)*(selfref(indself(icol,lay),ig) +& 264 selffrac(icol,lay)*(selfref(indself(icol,lay)+1,ig)- selfref(indself(icol,lay),ig))) 265 taufor(icol) = forfac(icol,lay)*( forref(indfor(icol,lay),ig) +& 266 forfrac(icol,lay)*( forref(indfor(icol,lay)+1,ig) -forref(indfor(icol,lay),ig))) 267 END DO 268 269 ! Compute absn2o term: 270 ! Note the use of 2D bilinear interpolation ka_mn2o lookup table 271 ! values 272 ! ===================================================================== 273 274 !dir$ SIMD 275 DO icol=startCol,ncol !vectorizes with dir 14.0.2 276 speccomb_mn2o(icol) = colh2o(icol,lay) +refrat_m_a*colco2(icol,lay) 277 specparm_mn2o(icol) = colh2o(icol,lay)/speccomb_mn2o(icol) 278 END DO 279 280 do icol=startCol,ncol ! vectorizes as is 14.0.2 281 IF (specparm_mn2o(icol) .GE. oneminus) specparm_mn2o(icol) =oneminus 282 end do 283 284 !dir$ SIMD ! vectorizes with dir 14.0.2 285 DO icol=startCol,ncol 286 specmult_mn2o(icol) = 8.0_wp*specparm_mn2o(icol) 287 jmn2o(icol) = 1 + INT(specmult_mn2o(icol)) 288 fmn2o(icol) = MOD(specmult_mn2o(icol),1.0_wp) 289 fmn2omf(icol) = minorfrac(icol,lay)*fmn2o(icol) 290 END DO 291 292 ! 293 ! 2D bilinear interpolation 294 ! ========================= 295 296 !dir$ SIMD 297 do icol=startCol,ncol ! vectorizes with dir 14.0.2 298 n2om1(icol) = ka_mn2o(jmn2o(icol),indminor(icol,lay) ,ig) + & 299 fmn2o(icol)*(ka_mn2o(jmn2o(icol)+1,indminor(icol,lay),ig) - & 300 ka_mn2o(jmn2o(icol),indminor(icol,lay) ,ig)) 301 n2om2(icol) = ka_mn2o(jmn2o(icol),indminor(icol,lay)+1,ig) + & 302 fmn2o(icol)*(ka_mn2o(jmn2o(icol)+1,indminor(icol,lay)+1,ig)- & 303 ka_mn2o(jmn2o(icol),indminor(icol,lay)+1,ig)) 304 absn2o(icol) = n2om1(icol) + minorfrac(icol,lay)*(n2om2(icol) -n2om1(icol)) 305 end do 306 307 ! In atmospheres where the amount of N2O is too great to be 308 ! considered 309 ! a minor species, adjust the column amount of N2O by an empirical 310 ! factor 311 ! to obtain the proper contribution. 312 ! ======================================================================== 313 314 !dir$ SIMD 315 do icol=startCol,ncol ! vectorized with dir 14.0.2 316 chi_n2o(icol) = coln2o(icol,lay)/coldry(icol,lay) 317 ratn2o(icol) = 1.e20*chi_n2o(icol)/chi_mls(4,jp(icol,lay)+1) 318 end do 319 320 do icol=startCol,ncol ! vectorizes as is 14.0.2 321 IF (ratn2o(icol) .GT. 1.5_wp) THEN 322 adjfac(icol) = 0.5_wp+(ratn2o(icol)-0.5_wp)**0.65_wp 323 adjcoln2o(icol) =adjfac(icol)*chi_mls(4,jp(icol,lay)+1)*coldry(icol,lay)*1.e-20_wp 324 ELSE 325 adjcoln2o(icol) = coln2o(icol,lay) 326 ENDIF 327 end do 328 329 ! Compute taug, one of two terms returned by the lower atmosphere 330 ! kernel (the other is fracs) 331 ! This loop could be parallelized over specparm types (i) but might 332 ! produce 333 ! different results for different thread counts 334 ! =========================================================================================== 335 336 !dir$ SIMD ! DOES NOT VECTORIZE even with SIMD dir 14.0.2 337 DO icol=startCol,ncol 338 taug(icol,lay) = tau_major(icol,0) + tau_major(icol,1) +tauself(icol) + taufor(icol) + adjcoln2o(icol)*absn2o(icol) 339 END DO 340 341 !dir$ SIMD ! vectorizes with dir 14.0.2 342 DO icol=startCol,ncol 343 speccomb_planck(icol) = colh2o(icol,lay)+refrat_planck_a*colco2(icol,lay) 344 specparm_planck(icol) = colh2o(icol,lay)/speccomb_planck(icol) 345 END DO 346 347 DO icol=startCol,ncol ! vectorizes as is 14.0.2 348 IF (specparm_planck(icol) .GE. oneminus) specparm_planck(icol)=oneminus 349 END DO 350 351 !dir$ SIMD 352 DO icol=startCol,ncol !vectorizes with dir 14.0.2 353 specmult_planck(icol) = 8.0_wp*specparm_planck(icol) 354 jpl(icol)= 1 + INT(specmult_planck(icol)) 355 fpl(icol) = MOD(specmult_planck(icol),1.0_wp) 356 fracs(icol,lay) = fracrefa(ig,jpl(icol)) + fpl(icol)*(fracrefa(ig,jpl(icol)+1)-fracrefa(ig,jpl(icol))) 357 END DO 358 rrpk_counter=rrpk_counter+1 359 END DO ! end lower atmosphere loop 360 END SUBROUTINE taumol03_lwr 361 362 363 SUBROUTINE taumol03_upr(ncol, startCol, laytrop, nlayers, & 364 rat_h2oco2, colco2, colh2o, coln2o, coldry, & 365 fac0, fac1, minorfrac, & 366 forfac,forfrac, & 367 jp, jt, ig, & 368 indfor, indminor, & 369 taug, fracs) 370 371 use mo_kind, only : wp 372 USE mo_lrtm_setup, ONLY: nspa 373 USE mo_lrtm_setup, ONLY: nspb 374 USE rrlw_planck, ONLY: chi_mls 375 USE rrlw_kg03, ONLY: selfref 376 USE rrlw_kg03, ONLY: forref 377 USE rrlw_kg03, ONLY: ka_mn2o 378 USE rrlw_kg03, ONLY: absa 379 USE rrlw_kg03, ONLY: fracrefa 380 USE rrlw_kg03, ONLY: kb_mn2o 381 USE rrlw_kg03, ONLY: absb 382 USE rrlw_kg03, ONLY: fracrefb 383 384 IMPLICIT NONE 385 386 real(kind=wp), PARAMETER :: oneminus = 1.0_wp - 1.0e-06_wp 387 388 integer, intent(in) :: ncol ! number of simd columns 389 integer, intent(in) :: startCol ! starting index for iterations in order to support parallelization across architectures 390 integer, intent(in) :: laytrop ! number of layers for lower atmosphere kernel 391 integer, intent(in) :: nlayers ! total number of layers 392 real(kind=wp), intent(in), dimension(startCol:ncol,0:1,nlayers) :: rat_h2oco2,fac0,fac1 ! not sure of ncol depend 393 real(kind=wp), intent(in), dimension(startCol:ncol,nlayers) :: colco2,colh2o,coln2o,coldry ! these appear to be gas concentrations 394 real(kind=wp), intent(in), dimension(startCol:ncol,nlayers) :: forfac,forfrac ! not sure of ncol depend 395 real(kind=wp), intent(in), dimension(startCol:ncol,nlayers) :: minorfrac ! not sure of ncol depend 396 397 ! Look up tables and related lookup indices 398 ! I assume all lookup indices depend on 3D position 399 ! ================================================= 400 401 integer, intent(in) :: jp(startCol:ncol,nlayers) ! I assume jp depends on ncol 402 integer, intent(in) :: jt(startCol:ncol,0:1,nlayers) ! likewise for jt 403 integer, intent(in) :: ig ! ig indexes into lookup tables 404 integer, intent(in) :: indfor(startCol:ncol,nlayers) ! for index array 405 integer, intent(in) :: indminor(startCol:ncol,nlayers) ! ka_mn2o index array 406 real(kind=wp), intent(out), dimension(startCol:ncol,nlayers) :: taug ! kernel result 407 real(kind=wp), intent(out), dimension(startCol:ncol,nlayers) :: fracs ! kernel result 408 409 ! Local variable 410 ! ============== 411 412 integer :: lay ! layer index 413 integer :: i ! specparm types index 414 integer :: icol ! column index 415 416 ! vector temporaries 417 ! ==================== 418 419 real(kind=wp), dimension(startCol:ncol) :: fs 420 real(kind=wp), dimension(startCol:ncol) :: fmn2o,fmn2omf 421 real(kind=wp), dimension(startCol:ncol) :: fpl 422 real(kind=wp), dimension(startCol:ncol) :: specmult, speccomb, specparm 423 real(kind=wp), dimension(startCol:ncol) :: specmult_mn2o, speccomb_mn2o, specparm_mn2o 424 real(kind=wp), dimension(startCol:ncol) :: specmult_planck, speccomb_planck,specparm_planck 425 real(kind=wp), dimension(startCol:ncol) :: n2om1,n2om2,absn2o,adjcoln2o,adjfac,chi_n2o,ratn2o 426 real(kind=wp), dimension(startCol:ncol,0:1) :: tau_major 427 real(kind=wp), dimension(startCol:ncol) :: taufor,tauself 428 integer, dimension(startCol:ncol) :: js, ind0, jmn2o, jpl 429 430 ! Register temporaries 431 ! ==================== 432 433 real(kind=wp) :: p2,fk0,fk1,fk2 434 real(kind=wp) :: refrat_planck_a, refrat_planck_b 435 real(kind=wp) :: refrat_m_a, refrat_m_b 436 437 !dir$ assume_aligned jp:64 438 !dir$ assume_aligned jt:64 439 !dir$ assume_aligned rat_h2oco2:64 440 !dir$ assume_aligned colco2:64 441 !dir$ assume_aligned colh2o:64 442 !dir$ assume_aligned fac0:64 443 !dir$ assume_aligned fac1:64 444 !dir$ assume_aligned taug:64 445 446 !dir$ assume_aligned specmult:64 447 !dir$ assume_aligned speccomb:64 448 !dir$ assume_aligned specparm:64 449 !dir$ assume_aligned specmult_mn2o:64 450 !dir$ assume_aligned speccomb_mn2o:64 451 !dir$ assume_aligned specparm_mn2o:64 452 !dir$ assume_aligned specmult_planck:64 453 !dir$ assume_aligned speccomb_planck:64 454 !dir$ assume_aligned specparm_planck:64 455 !dir$ assume_aligned fs:64 456 !dir$ assume_aligned tau_major:64 457 !dir$ assume_aligned chi_n2o:64 458 459 !dir$ assume_aligned js:64 460 !dir$ assume_aligned ind0:64 461 !dir$ assume_aligned jpl:64 462 !dir$ assume_aligned fpl:64 463 464 !dir$ assume_aligned absn2o:64 465 !dir$ assume_aligned adjcoln2o:64 466 !dir$ assume_aligned adjfac:64 467 !dir$ assume_aligned ratn2o:64 468 469 ! Upper atmosphere loop 470 ! ======================== 471 refrat_planck_b = chi_mls(1,13)/chi_mls(2,13) 472 refrat_m_b = chi_mls(1,13)/chi_mls(2,13) 473 DO lay = laytrop+1, nlayers 474 475 DO i=0,1 ! loop over specparm types 476 477 ! This loop should vectorize 478 ! ============================= 479 480 !dir$ SIMD 481 DO icol=startCol,ncol ! Vectorizes with dir 14.0.2 482 speccomb(icol) = colh2o(icol,lay) + rat_h2oco2(icol,i,lay)*colco2(icol,lay) 483 specparm(icol) = colh2o(icol,lay)/speccomb(icol) 484 IF (specparm(icol) .ge. oneminus) specparm(icol) = oneminus 485 specmult(icol) = 4.0_wp*(specparm(icol)) 486 js(icol) = 1 + INT(specmult(icol)) 487 fs(icol) = MOD(specmult(icol),1.0_wp) 488 ind0(icol) = ((jp(icol,lay)-13+i)*5+(jt(icol,i,lay)-1))*nspb(3) +js(icol) 489 END DO 490 491 !dir$ SIMD 492 DO icol=startCol,ncol ! Vectorizes with dir 14.0.2 493 tau_major(icol,i) = speccomb(icol) * & 494 ((1.0_wp - fs(icol))*fac0(icol,i,lay)*absb(ind0(icol) ,ig) + & 495 fs(icol) *fac0(icol,i,lay)*absb(ind0(icol)+1,ig) + & 496 (1.0_wp - fs(icol))*fac1(icol,i,lay)*absb(ind0(icol)+5,ig) + & 497 fs(icol) *fac1(icol,i,lay)*absb(ind0(icol)+6,ig)) 498 END DO 499 500 END DO ! end loop over specparm types for tau_major calculation 501 502 ! Compute taufor terms 503 ! Note the use of 1D bilinear interpolation of selfref and forref lookup 504 ! table values 505 ! =================================================================================== 506 !dir$ SIMD 507 DO icol=startCol,ncol ! Vectorizes with dir 14.0.2 508 taufor(icol) = forfac(icol,lay)*( forref(indfor(icol,lay),ig) + & 509 forfrac(icol,lay)*( forref(indfor(icol,lay)+1,ig) - forref(indfor(icol,lay),ig))) 510 END DO 511 512 ! Compute absn2o term: 513 ! Note the use of 2D bilinear interpolation ka_mn2o lookup table values 514 ! ===================================================================== 515 !$DIR SIMD 516 DO icol=startCol,ncol ! Vectorizes with dir 14.0.2 517 speccomb_mn2o(icol) = colh2o(icol,lay) + refrat_m_b*colco2(icol,lay) 518 specparm_mn2o(icol) = colh2o(icol,lay)/speccomb_mn2o(icol) 519 IF (specparm_mn2o(icol) .GE. oneminus) specparm_mn2o(icol) = oneminus 520 specmult_mn2o(icol) = 4.0_wp*specparm_mn2o(icol) 521 jmn2o(icol) = 1 + INT(specmult_mn2o(icol)) 522 fmn2o(icol) = MOD(specmult_mn2o(icol),1.0_wp) 523 fmn2omf(icol) = minorfrac(icol,lay)*fmn2o(icol) 524 END DO 525 526 ! In atmospheres where the amount of N2O is too great to be considered 527 ! a minor species, adjust the column amount of N2O by an empirical factor 528 ! to obtain the proper contribution. 529 ! ======================================================================== 530 531 !dir$ SIMD 532 DO icol=startCol,ncol ! loop vectorized with directive 14.0.2 533 chi_n2o(icol) = coln2o(icol,lay)/coldry(icol,lay) 534 ratn2o(icol) = 1.e20*chi_n2o(icol)/chi_mls(4,jp(icol,lay)+1) 535 END DO 536 537 DO icol=startCol,ncol ! Loop vectorized as is 14.0.2 538 IF (ratn2o(icol) .GT. 1.5_wp) THEN 539 adjfac(icol) = 0.5_wp+(ratn2o(icol)-0.5_wp)**0.65_wp 540 adjcoln2o(icol) = adjfac(icol)*chi_mls(4,jp(icol,lay)+1)*coldry(icol,lay)*1.e-20_wp 541 ELSE 542 adjcoln2o(icol) = coln2o(icol,lay) 543 ENDIF 544 END DO 545 546 ! 547 ! 2D bilinear interpolation 548 ! ========================= 549 550 !dir$ SIMD 551 DO icol=startCol,ncol ! loop vectorizes with directive 14.0.2 552 n2om1(icol) = kb_mn2o(jmn2o(icol),indminor(icol,lay) ,ig) + & 553 fmn2o(icol)*(kb_mn2o(jmn2o(icol)+1,indminor(icol,lay),ig) - & 554 kb_mn2o(jmn2o(icol) ,indminor(icol,lay) ,ig)) 555 n2om2(icol) = kb_mn2o(jmn2o(icol),indminor(icol,lay)+1,ig) + & 556 fmn2o(icol)*(kb_mn2o(jmn2o(icol)+1,indminor(icol,lay)+1,ig)- & 557 kb_mn2o(jmn2o(icol),indminor(icol,lay)+1,ig)) 558 absn2o(icol) = n2om1(icol) + minorfrac(icol,lay)*(n2om2(icol) -n2om1(icol)) 559 END DO 560 561 !dir$ SIMD 562 DO icol=startCol,ncol ! loop vectorizes with directive 14.0.2 563 taug(icol,lay) = tau_major(icol,0) + tau_major(icol,1) + taufor(icol) + adjcoln2o(icol)*absn2o(icol) 564 END DO 565 566 !dir$ SIMD 567 DO icol=startCol,ncol ! loop vectorizes with directive 14.0.2 568 speccomb_planck(icol) = colh2o(icol,lay)+refrat_planck_b*colco2(icol,lay) 569 specparm_planck(icol) = colh2o(icol,lay)/speccomb_planck(icol) 570 END DO 571 572 !dir$ SIMD 573 DO icol=startCol,ncol ! loop vectorizes with directive 14.0.2 574 IF (specparm_planck(icol) .GE. oneminus) specparm_planck(icol)=oneminus 575 specmult_planck(icol) = 4.0_wp*specparm_planck(icol) 576 jpl(icol)= 1 + INT(specmult_planck(icol)) 577 fpl(icol) = MOD(specmult_planck(icol),1.0_wp) 578 fracs(icol,lay) = fracrefb(ig,jpl(icol)) + fpl(icol)*(fracrefb(ig,jpl(icol)+1)-fracrefb(ig,jpl(icol))) 579 END DO 580 END DO ! nlayers loop 581 582 END SUBROUTINE taumol03_upr 583 584END MODULE mo_taumol03 585