- Timestamp:
- Dec 2, 2024, 5:02:15 PM (39 hours ago)
- Location:
- trunk/LMDZ.VENUS/libf/phyvenus
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.VENUS/libf/phyvenus/photolysis_mod.F90
r2925 r3530 40 40 real, dimension(nw), save :: xshocl ! hocl absorption cross-section (cm2) 41 41 real, dimension(nw), save :: xsso2_200, xsso2_298, xsso2_360 ! so2 absorption cross-section (cm2) 42 real, dimension(nw), save :: xsso 42 real, dimension(nw), save :: xsso_150, xsso_250 ! so absorption cross-section (cm2) 43 43 real, dimension(nw), save :: xsso3 ! so3 absorption cross-section (cm2) 44 44 real, dimension(nw), save :: xsclo ! clo absorption cross-section (cm2) … … 58 58 real, dimension(nw), save :: albedo ! surface albedo 59 59 60 !$OMP THREADPRIVATE(f, xsco2_195, xsco2_295, xsco2_370,yieldco2,xso2_150, xso2_200, xso2_250, xso2_300,yieldo2,xso3_218, xso3_298,xsh2o,xsh2o2,xsho2,xsh2,yieldh2,xsno2, xsno2_220, xsno2_294,yldno2_248,yldno2_298)61 !$OMP THREADPRIVATE(xsno, yieldno,xsn2,yieldn2,xshdo,albedo)60 !$OMP THREADPRIVATE(f, xsco2_195, xsco2_295, xsco2_370, yieldco2, xso2_150, xso2_200, xso2_250, xso2_300, yieldo2, xso3_218, xso3_298, xsh2o, xsh2o2, xsho2, xsh2, yieldh2, xsno2, xsno2_220, xsno2_294, yldno2_248, yldno2_298) 61 !$OMP THREADPRIVATE(xsno, yieldno, xsn2, yieldn2, xshdo, albedo) 62 62 63 63 contains … … 123 123 ! read and grid so cross-sections 124 124 125 call rdxsso(nw,wl,xsso )125 call rdxsso(nw,wl,xsso_150,xsso_250) 126 126 127 127 ! read and grid so3 cross-sections … … 2327 2327 !============================================================================== 2328 2328 2329 subroutine rdxsso(nw, wl, yg)2329 subroutine rdxsso(nw,wl,xsso_150,xsso_250) 2330 2330 2331 2331 !-----------------------------------------------------------------------------* 2332 2332 != PURPOSE: =* 2333 != Read SO cross-sections =*2334 != JPL 2006 recommendation=*2333 != Read SO cross-sections =* 2334 != Phillips (1981) or Heays et al. (2023) =* 2335 2335 !-----------------------------------------------------------------------------* 2336 2336 != PARAMETERS: =* … … 2342 2342 USE mod_phys_lmdz_transfert_para, ONLY: bcast 2343 2343 2344 IMPLICIT NONE2344 implicit none 2345 2345 2346 2346 ! input 2347 2347 2348 integer :: nw ! number of wavelength grid points 2348 real, dimension(nw) :: wl ! lower wavelength for each interval2349 real, dimension(nw) :: wl ! lower wavelength for each interval 2349 2350 2350 2351 ! output 2351 2352 2352 real, dimension(nw) :: yg ! SOcross-sections (cm2)2353 real, dimension(nw) :: xsso_150, xsso_250 ! so cross-sections (cm2) 2353 2354 2354 2355 ! local 2356 2357 integer, parameter :: kdata = 1000 2358 integer :: i, l, n, n1, n2, ierr, iopt 2359 integer :: kin, kout ! input/output logical units 2360 real, dimension(nw) :: yg, yg1, yg2 2361 real, dimension(kdata) :: x1, y1, x2, y2 2355 2362 real :: dummy 2356 2363 real, parameter :: deltax = 1.e-4 2357 integer, parameter :: kdata = 10002358 real, dimension(kdata) :: x1, y12359 integer :: i, n, ierr2360 2364 character*100 fil 2361 integer :: kin, kout ! input/output logical units2362 2365 2363 2366 kin = 10 2364 2367 2365 !*** cross sections from JPL [2006] 2366 2367 fil = 'cross_sections/so_marcq.txt' 2368 print*, 'section efficace SO: ', fil 2369 2370 if(is_master) then 2368 ! iopt = 1 : Phillips (1981) 2369 ! iopt = 2 : Heays et al., Molecular Physics, 2023 2370 2371 iopt = 2 2372 2373 if (iopt == 1) then 2374 2375 ! cross sections from Phillips (1981) 2376 2377 fil = 'cross_sections/so_marcq.txt' 2378 2379 if (is_master) then 2380 2381 print*, 'section efficace SO: ', fil 2371 2382 2372 2383 n = 851 2373 OPEN(kin,FILE=fil,STATUS='OLD') 2374 DO i = 1,n 2375 READ(kin,*) x1(i), dummy, dummy, dummy, dummy, y1(i) 2376 y1(i) = y1(i) * 0.88 2377 ENDDO 2378 CLOSE(kin) 2379 2380 CALL addpnt(x1,y1,kdata,n,x1(1)*(1.-deltax),0.) 2381 CALL addpnt(x1,y1,kdata,n, 0.,0.) 2382 CALL addpnt(x1,y1,kdata,n,x1(n)*(1.+deltax),0.) 2383 CALL addpnt(x1,y1,kdata,n, 1E38,0.) 2384 2385 CALL inter2(nw,wl,yg,n,x1,y1,ierr) 2384 OPEN(kin,FILE=fil,STATUS='OLD') 2385 DO i = 1,n 2386 READ(kin,*) x1(i), dummy, dummy, dummy, dummy, y1(i) 2387 y1(i) = y1(i)*0.88 ! from Mills, PhD, 1998 2388 ENDDO 2389 CLOSE(kin) 2390 2391 CALL addpnt(x1,y1,kdata,n,x1(1)*(1.-deltax),0.) 2392 CALL addpnt(x1,y1,kdata,n, 0.,0.) 2393 CALL addpnt(x1,y1,kdata,n,x1(n)*(1.+deltax),0.) 2394 CALL addpnt(x1,y1,kdata,n, 1E38,0.) 2395 2396 do i = 1,n 2397 print*, x1(i), y1(i) 2398 end do 2399 2400 CALL inter2(nw,wl,yg,n,x1,y1,ierr) 2386 2401 2387 IF (ierr .NE. 0) THEN 2388 WRITE(*,*) ierr, fil 2389 STOP 2390 ENDIF 2391 2392 endif !is_master 2402 IF (ierr .NE. 0) THEN 2403 WRITE(*,*) ierr, fil 2404 STOP 2405 ENDIF 2406 2407 DO l = 1, nw-1 2408 xsso_150(l) = yg(l) 2409 xsso_250(l) = yg(l) 2410 END DO 2411 2412 endif !is_master 2413 2414 else if (iopt == 2) then 2415 2416 ! cross sections from Heays et al. (2023) 2417 ! values given at 150 K and 250 K 2418 2419 fil = 'cross_sections/SO_heays_2023.txt' 2420 2421 if (is_master) then 2422 2423 print*, 'section efficace SO: ', fil 2424 2425 n1 = 525 2426 n2 = n1 2427 OPEN(kin,FILE=fil,STATUS='OLD') 2428 2429 DO i = 1, 3 2430 READ(kin,*) 2431 ENDDO 2432 2433 DO i = 1,n1 2434 READ(kin,*) x1(i), dummy, dummy, dummy, dummy, dummy, & 2435 dummy, dummy, dummy, y1(i), y2(i) 2436 x2(i) = x1(i) 2437 ENDDO 2438 CLOSE(kin) 2439 2440 CALL addpnt(x1,y1,kdata,n1,x1(1)*(1.-deltax),0.) 2441 CALL addpnt(x1,y1,kdata,n1, 0.,0.) 2442 CALL addpnt(x1,y1,kdata,n1,x1(n1)*(1.+deltax),0.) 2443 CALL addpnt(x1,y1,kdata,n1, 1.e+38,0.) 2444 2445 CALL inter2(nw,wl,yg1,n1,x1,y1,ierr) 2446 2447 IF (ierr .NE. 0) THEN 2448 WRITE(*,*) ierr, fil 2449 STOP 2450 ENDIF 2451 2452 CALL addpnt(x2,y2,kdata,n2,x2(1)*(1.-deltax),0.) 2453 CALL addpnt(x2,y2,kdata,n2, 0.,0.) 2454 CALL addpnt(x2,y2,kdata,n2,x2(n2)*(1.+deltax),0.) 2455 CALL addpnt(x2,y2,kdata,n2, 1.e+38,0.) 2456 2457 CALL inter2(nw,wl,yg2,n2,x2,y2,ierr) 2458 2459 IF (ierr .NE. 0) THEN 2460 WRITE(*,*) ierr, fil 2461 STOP 2462 ENDIF 2463 2464 DO l = 1, nw-1 2465 xsso_150(l) = yg1(l) 2466 xsso_250(l) = yg2(l) 2467 END DO 2468 2469 endif !is_master 2470 2471 end if ! iopt 2393 2472 2394 2473 call bcast(yg) -
trunk/LMDZ.VENUS/libf/phyvenus/photolysis_online.F
r2974 r3530 183 183 $ dtgas(1:nlayer,:,a_h2o2), sj) 184 184 185 ! so: 186 187 call setso(nphot, nlayer, nw, wc, temp, xsso_150, xsso_250, 188 $ j_so, colinc(1:nlayer), rm(:,i_so), 189 $ dtgas(1:nlayer,:,a_so), sj) 190 185 191 ! so2: 186 192 187 call setso2(nphot, nlayer, nw, wc, temp, xsso2_200, xsso2_298,188 $ xsso2_360, j_so2, colinc(1:nlayer), rm(:,i_so2),193 call setso2(nphot, nlayer, nw, wc, temp, xsso2_200, xsso2_298, 194 $ xsso2_360, j_so2, colinc(1:nlayer), rm(:,i_so2), 189 195 $ dtgas(1:nlayer,:,a_so2), sj) 190 196 … … 209 215 dtgas(ilay,iw,a_hocl) = colinc(ilay)*rm(ilay,i_hocl) 210 216 $ *xshocl(iw) 211 dtgas(ilay,iw,a_so) = colinc(ilay)*rm(ilay,i_so)*xsso(iw)212 217 dtgas(ilay,iw,a_so3) = colinc(ilay)*rm(ilay,i_so3)*xsso3(iw) 213 218 dtgas(ilay,iw,a_s2) = colinc(ilay)*rm(ilay,i_s2)*xss2(iw) … … 245 250 sj(ilay,iw,j_cl2) = xscl2(iw) ! cl2 246 251 sj(ilay,iw,j_hocl) = xshocl(iw) ! hocl 247 sj(ilay,iw,j_so) = xsso(iw) ! so248 252 sj(ilay,iw,j_s2) = xss2(iw) ! s2 249 253 sj(ilay,iw,j_so3) = xsso3(iw) ! so3 … … 647 651 648 652 end subroutine seto3 653 !============================================================================== 654 655 subroutine setso(nd, nlayer, nw, wc, tlay, xsso_150, xsso_250, 656 $ j_so, colinc, rm, dt, sj) 657 658 !-----------------------------------------------------------------------------* 659 != PURPOSE: =* 660 != Set up the so temperature dependent cross-sections and optical depth =* 661 !-----------------------------------------------------------------------------* 662 663 implicit none 664 665 ! input: 666 667 integer :: nd ! number of photolysis rates 668 integer :: nlayer ! number of vertical layers 669 integer :: nw ! number of wavelength grid points 670 integer :: j_so ! photolysis indexe 671 real, dimension(nw) :: wc ! central wavelength for each interval 672 real, dimension(nw) :: xsso_150, xsso_250 ! so cross-sections (cm2) 673 real, dimension(nlayer) :: tlay ! temperature (k) 674 real, dimension(nlayer) :: rm ! so mixing ratio 675 real, dimension(nlayer) :: colinc ! air column increment (molecule.cm-2) 676 677 ! output: 678 679 real, dimension(nlayer,nw) :: dt ! optical depth 680 real, dimension(nlayer,nw,nd) :: sj ! cross-section array (cm2) 681 682 ! local: 683 ! 684 integer :: ilev, iw 685 real :: temp, xsso 686 687 ! calculate temperature dependance 688 689 do ilev = 1,nlayer 690 temp = max(tlay(ilev),150.) 691 temp = min(temp, 250.) 692 do iw = 1, nw-1 693 if (xsso_150(iw) /= 0. .and. xsso_250(iw) /= 0.) then 694 xsso = log(xsso_150(iw)) 695 $ + (log(xsso_250(iw)) - log(xsso_150(iw))) 696 $ /(250. - 150.)*(temp - 150.) 697 698 sj(ilev,iw,j_so) = exp(xsso) 699 else 700 sj(ilev,iw,j_so) = 0. 701 end if 702 703 ! optical depth 704 705 dt(ilev,iw) = colinc(ilev)*rm(ilev)*sj(ilev,iw,j_so) 706 707 end do 708 end do 709 710 end subroutine setso 711 649 712 !============================================================================== 650 713
Note: See TracChangeset
for help on using the changeset viewer.