Changeset 2566
- Timestamp:
- Sep 24, 2021, 1:46:05 PM (3 years ago)
- Location:
- trunk/LMDZ.MARS
- Files:
-
- 3 edited
-
README (modified) (1 diff)
-
libf/phymars/co2condens_mod.F (modified) (14 diffs)
-
libf/phymars/physiq_mod.F (modified) (3 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.MARS/README
r2563 r2566 3477 3477 wstats will be in the stats.nc file (which matches the behaviour prior to 3478 3478 this improvement). 3479 3480 == 24/09/2021 == CM 3481 merge co2condens_mod4micro.F with co2condens_mod.F 3482 variable used: 3483 if (co2clouds) 3484 condens_layer, condens_column 3485 else 3486 zcondicea, zfallice -
trunk/LMDZ.MARS/libf/phymars/co2condens_mod.F
r2524 r2566 21 21 & igcm_hdo_ice, igcm_hdo_vap, 22 22 & nqperes,nqfils, ! MVals: variables isotopes 23 & qperemin,masseqmin 23 & qperemin,masseqmin, 24 & igcm_co2 24 25 use surfdat_h, only: emissiv, phisfi 25 26 use geometry_mod, only: latitude, ! grid point latitudes (rad) … … 114 115 REAL zcondices(ngrid) ! condensation rate on the ground (kg/m2/s) 115 116 REAL zfallice(ngrid,nlayer+1) ! amount of ice falling from layer l (kg/m2/s) 117 REAL condens_layer(ngrid,nlayer) ! co2clouds: condensation rate in layer l (kg/m2/s) 118 REAL condens_column(ngrid) ! co2clouds: sum(condens_layer(ig,:)) (kg/m2/s) 116 119 REAL zfallheat 117 120 REAL zmflux(nlayer+1) … … 326 329 327 330 ELSE ! if co2 clouds 328 DO ig=1,ngrid 329 zfallice(ig,1) = zdqssed_co2(ig) 330 ENDDO 331 condens_layer(:,:) = 0. 332 condens_column(:) = 0. 331 333 DO l=nlayer , 1, -1 332 334 DO ig=1,ngrid 333 zcondicea(ig,l) = pcondicea_co2microp(ig,l)*334 & (pplev(ig,l) - pplev(ig,l+1))/g335 condens_layer(ig,l) = pcondicea_co2microp(ig,l)* 336 & (pplev(ig,l) - pplev(ig,l+1))/g 335 337 ENDDO 336 338 ENDDO 337 339 DO ig=1,ngrid 340 condens_column(ig) = sum(condens_layer(ig,:)) 341 piceco2(ig) = piceco2(ig) + zdqssed_co2(ig)*ptimestep 342 ENDDO 338 343 ENDIF ! end of if co2clouds 339 344 … … 368 373 c (compute zcondices and pdtsrfc) 369 374 c 375 c No microphysics of CO2 clouds 376 if (.not.co2clouds)then 370 377 DO ig=1,ngrid 371 378 IF(latitude(ig).lt.0) THEN … … 392 399 zfallheat=0. 393 400 ENDIF 394 395 401 c condensation or partial sublimation of CO2 ice 396 402 c """"""""""""""""""""""""""""""""""""""""""""""" … … 399 405 pdtsrfc(ig) = (ztcondsol(ig) - ztsrf(ig))/ptimestep 400 406 zdiceco2(ig) = zcondices(ig) + zfallice(ig,1) 401 402 407 #ifdef MESOSCALE 403 408 print*, "not enough CO2 tracer in 1st layer to condense" … … 409 414 IF(ico2.ne.0) then 410 415 c Available CO2 tracer in layer 1 at end of timestep (kg/m2) 416 #ifndef MESOSCALE 411 417 availco2= pq(ig,1,ico2)*((ap(1)-ap(2))+ 412 418 & (bp(1)-bp(2))*(pplev(ig,1)/g-zdiceco2(ig)*ptimestep)) 413 419 #else 420 availco2 = pq(ig,1,igcm_co2) 421 PRINT*, "MESOSCALE: CO2 tracer AT FIRST LEVEL IS NOT' 422 & ' CORRECTED FROM SIGMA LEVELS" 423 #endif 414 424 IF ((zcondices(ig) + zcondicea(ig,1))*ptimestep 415 425 & .gt.availco2) then … … 418 428 pdtsrfc(ig)=(latcond/pcapcal(ig))* 419 429 & (zcondices(ig)+zfallheat) 420 ENDIF 421 ENDIF 430 ENDIF 431 ENDIF 422 432 #endif 423 433 … … 425 435 c """""""""""""""""""""""""""""""""""""""""""""""""""" 426 436 c (including what has just condensed in the atmosphere) 427 428 437 IF((piceco2(ig)/ptimestep+zfallice(ig,1)).LE. 429 438 & -zcondices(ig))THEN … … 436 445 c Changing CO2 ice amount and pressure : 437 446 c """""""""""""""""""""""""""""""""""" 438 439 447 piceco2(ig) = piceco2(ig) + zdiceco2(ig)*ptimestep 440 448 pdpsrf(ig) = -zdiceco2(ig)*g … … 453 461 END IF ! if there is condensation/sublimation 454 462 ENDDO ! of DO ig=1,ngrid 455 463 else 464 do ig=1,ngrid 465 if ((ztsrf(ig)<ztcondsol(ig)) .OR. (zdqssed_co2(ig)/=0.) .OR. 466 & ((ztsrf(ig)>ztcondsol(ig)) .AND. (piceco2(ig)/=0.))) then 467 condsub(ig) = .true. 468 zfallheat = 0. 469 !----------------------------------------------------------------------------------------------------------------------! 470 ! Compute direct condensation/sublimation of CO2 ice 471 !----------------------------------------------------------------------------------------------------------------------! 472 zcondices(ig) = pcapcal(ig) * (ztcondsol(ig)-ztsrf(ig)) / 473 & (latcond*ptimestep) - zfallheat 474 pdtsrfc(ig) = (ztcondsol(ig) - ztsrf(ig)) / ptimestep 475 !----------------------------------------------------------------------------------------------------------------------! 476 ! If there is not enough CO2 tracer in 1st layer to condense 477 !----------------------------------------------------------------------------------------------------------------------! 478 ! Available CO2 tracer in layer 1 at end of timestep (kg/m2) 479 #ifndef MESOSCALE 480 availco2 = pq(ig,1,igcm_co2) * ( (ap(1)-ap(2)) + 481 & (bp(1)-bp(2)) * (pplev(ig,1)/g - zcondices(ig) 482 & *ptimestep) ) 483 #else 484 availco2 = pq(ig,1,igcm_co2) 485 PRINT*, "MESOSCALE: CO2 tracer AT FIRST LEVEL IS NOT CORRECTED FROM SIGMA LEVELS" 486 #endif 487 if ( zcondices(ig) * ptimestep>availco2 ) then 488 zcondices(ig) = availco2/ptimestep 489 zdiceco2(ig) = zcondices(ig) + zdqssed_co2(ig) 490 pdtsrfc(ig) = (latcond/pcapcal(ig)) * 491 & (zcondices(ig)+zfallheat) 492 end if 493 !----------------------------------------------------------------------------------------------------------------------! 494 ! If the entire CO2 ice layer sublimes (including what has just condensed in the atmosphere) 495 !----------------------------------------------------------------------------------------------------------------------! 496 if ( (piceco2(ig)/ptimestep) <= -zcondices(ig) ) then 497 zcondices(ig) = -piceco2(ig)/ptimestep 498 pdtsrfc(ig) = (latcond/pcapcal(ig)) * 499 & (zcondices(ig)+zfallheat) 500 end if 501 !----------------------------------------------------------------------------------------------------------------------! 502 ! Changing CO2 ice amount and pressure 503 !----------------------------------------------------------------------------------------------------------------------! 504 zdiceco2(ig) = zcondices(ig) + zdqssed_co2(ig) + 505 & condens_column(ig) 506 piceco2(ig) = piceco2(ig) + zcondices(ig)*ptimestep 507 pdpsrf(ig) = -zdiceco2(ig) * g 508 509 if (abs(pdpsrf(ig)*ptimestep)>pplev(ig,1)) then 510 print *, 'STOP in condens' 511 print *, 'condensing more than total mass' 512 print *, 'Grid point ', ig 513 print *, 'Longitude(degrees): ', longitude_deg(ig) 514 print *, 'Latitude(degrees): ', latitude_deg(ig) 515 print *, 'Ps = ', pplev(ig,1) 516 print *, 'd Ps = ', pdpsrf(ig) 517 call abort_physic('co2condens', 518 & 'condensing more than total mass', 1) 519 end if 520 end if 521 end do 522 end if 456 523 c ******************************************************************** 457 524 c Surface albedo and emissivity of the surface below the snow (emisref) … … 498 565 c Mass fluxes through the sigma levels (kg.m-2.s-1) (>0 when up) 499 566 c """""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" 500 501 zmflux(1) = -zcondices(ig)502 DO l=1,nlayer503 zmflux(l+1) = zmflux(l) -zcondicea(ig,l)567 if (.not. co2clouds) then 568 zmflux(1) = -zcondices(ig) 569 DO l=1,nlayer 570 zmflux(l+1) = zmflux(l) -zcondicea(ig,l) 504 571 #ifndef MESOSCALE 505 & + (bp(l)-bp(l+1))*(zfallice(ig,1)-zmflux(1)) 506 c zmflux set to 0 if very low to avoid: top layer is disappearing in v1ld 507 if (abs(zmflux(l+1)).lt.1E-13.OR.bp(l+1).eq.0.) zmflux(l+1)=0. 572 & + (bp(l)-bp(l+1))*(zfallice(ig,1)-zmflux(1)) 573 c zmflux set to 0 if very low to avoid: top layer is disappearing in v1ld 574 if (abs(zmflux(l+1)).lt.1E-13.OR.bp(l+1).eq.0.) then 575 zmflux(l+1)=0. 576 end if 508 577 #else 509 if (abs(zmflux(l+1)).lt.1E-13) zmflux(l+1)=0. 578 zmflux(l+1) = zmflux(l) - zcondicea(ig,l) 579 if (abs(zmflux(l+1)).lt.1E-13) zmflux(l+1)=0. 580 PRINT*, "MESOSCALE: FLUX THROUGH SIGMA LEVELS from"// 581 & "dPS HAVE TO BE IMPLEMENTED" 510 582 #endif 511 END DO 512 583 END DO 584 else 585 zmflux(1) = - zcondices(ig) - zdqssed_co2(ig) 586 do l = 1, nlayer 587 #ifndef MESOSCALE 588 zmflux(l+1) = zmflux(l) - condens_layer(ig,l) 589 & + (bp(l)-bp(l+1)) * (-pdpsrf(ig)/g) 590 ! zmflux set to 0 if very low to avoid: top layer is disappearing in v1ld 591 if (abs(zmflux(l+1))<1E-13.OR.bp(l+1)==0.) then 592 zmflux(l+1) = 0. 593 end if 594 #else 595 zmflux(l+1) = zmflux(l) - condens_layer(ig,l) 596 if (abs(zmflux(l+1))<1E-13) then 597 zmflux(l+1) = 0. 598 end if 599 PRINT*, "MESOSCALE: FLUX THROUGH SIGMA LEVELS from"// 600 & "dPS HAVE TO BE IMPLEMENTED" 601 #endif 602 end do 603 end if 513 604 #ifdef MESOSCALE 514 605 print*,"absurd mass set because bp not available" … … 642 733 #endif 643 734 735 if (.not. co2clouds) then 736 do iq=1,nq 737 ! if (noms(iq).eq.'co2') then 738 if (iq.eq.ico2) then 739 c SPECIAL Case when the tracer IS CO2 : 740 DO l=1,nlayer 741 pdqc(ig,l,iq)= (1/masse(l)) * 742 & ( zmflux(l)*(zqm(l,iq) - zqc(l,iq)) 743 & - zmflux(l+1)*(zqm(l+1,iq) - zqc(l,iq)) 744 & + zcondicea(ig,l)*(zqc(l,iq)-1.) ) 745 END DO 746 else 747 DO l=1,nlayer 748 pdqc(ig,l,iq)= (1/masse(l)) * 749 & ( zmflux(l)*(zqm(l,iq) - zqc(l,iq)) 750 & - zmflux(l+1)*(zqm(l+1,iq) - zqc(l,iq)) 751 & + zcondicea(ig,l)*zqc(l,iq) ) 752 753 pdqc(ig,l,iq)=pdqc(ig,l,iq)+zdq_scav(ig,l,iq) ! zdq_scav=0 if watercloud=false 754 END DO 755 end if 756 enddo 757 else 644 758 c Tendencies on Q 645 do iq=1,nq646 ! if (noms(iq).eq.'co2') then647 if (iq.eq.ico2) then648 c SPECIAL Case when the tracer IS CO2 :649 DO l=1,nlayer650 pdqc(ig,l,iq)= (1/masse(l)) *651 & ( zmflux(l)*(zqm(l,iq) - zqc(l,iq))652 & - zmflux(l+1)*(zqm(l+1,iq) - zqc(l,iq))653 & + zcondicea(ig,l)*(zqc(l,iq)-1.) )654 END DO655 else656 DO l=1,nlayer657 pdqc(ig,l,iq)= (1/masse(l)) *658 & ( zmflux(l)*(zqm(l,iq) - zqc(l,iq))659 & - zmflux(l+1)*(zqm(l+1,iq) - zqc(l,iq))660 & + zcondicea(ig,l)*zqc(l,iq) )661 662 pdqc(ig,l,iq)=pdqc(ig,l,iq)+zdq_scav(ig,l,iq) ! zdq_scav=0 if watercloud=false663 END DO664 665 enddo 759 do iq=1,nq 760 ! if (noms(iq).eq.'co2') then 761 if (iq.eq.ico2) then 762 c SPECIAL Case when the tracer IS CO2 : 763 DO l=1,nlayer 764 pdqc(ig,l,iq)= (1/masse(l)) * 765 & ( zmflux(l)*(zqm(l,iq) - zqc(l,iq)) 766 & - zmflux(l+1)*(zqm(l+1,iq) - zqc(l,iq)) 767 & + condens_layer(ig,l)*(zqc(l,iq)-1.) ) 768 END DO 769 else 770 DO l=1,nlayer 771 pdqc(ig,l,iq)= (1/masse(l)) * 772 & ( zmflux(l)*(zqm(l,iq) - zqc(l,iq)) 773 & - zmflux(l+1)*(zqm(l+1,iq) - zqc(l,iq)) 774 & + condens_layer(ig,l)*zqc(l,iq) ) 775 END DO 776 end if 777 enddo 778 end if 779 666 780 667 781 end if ! if (condsub) … … 671 785 c CO2 snow / clouds scheme 672 786 c *************************************************************** 673 674 call co2snow(ngrid,nlayer,ptimestep,emisref,condsub,pplev, 675 & zcondicea,zcondices,zfallice,pemisurf) 676 787 if (.not. co2clouds) then 788 call co2snow(ngrid,nlayer,ptimestep,emisref,condsub,pplev, 789 & zcondicea,zcondices,zfallice,pemisurf) 790 else 791 call co2snow(ngrid,nlayer,ptimestep,emisref,condsub,pplev, 792 & zcondicea,zcondices,zdqssed_co2,pemisurf) 793 end if 677 794 c *************************************************************** 678 795 c Ecriture des diagnostiques -
trunk/LMDZ.MARS/libf/phymars/physiq_mod.F
r2563 r2566 1999 1999 zdvc(:,:) = 0. 2000 2000 zdqc(:,:,:) = 0. 2001 2002 if (co2clouds) then 2003 CALL co2condens4micro(ngrid,nlayer,nq,ptimestep, 2004 $ capcal,zplay,zplev,tsurf,pt, 2005 $ pphi,pdt,pdu,pdv,zdtsurf,pu,pv,pq,pdq, 2006 $ co2ice,albedo,emis, 2007 $ zdtc,zdtsurfc,pdpsrf,zduc,zdvc,zdqc, 2008 $ fluxsurf_sw,zls, 2009 $ zdqssed_co2,zcondicea_co2microp 2010 & ) 2011 ! no scavenging yet 2012 zdqsc(:,:) = 0. 2013 else 2014 CALL co2condens(ngrid,nlayer,nq,ptimestep, 2001 CALL co2condens(ngrid,nlayer,nq,ptimestep, 2015 2002 $ capcal,zplay,zplev,tsurf,pt, 2016 2003 $ pphi,pdt,pdu,pdv,zdtsurf,pu,pv,pq,pdq, … … 2020 2007 $ zdqssed_co2,zcondicea_co2microp, 2021 2008 & zdqsc) 2022 2009 DO iq=1, nq 2023 2010 DO ig=1,ngrid 2024 2011 dqsurf(ig,iq)=dqsurf(ig,iq)+zdqsc(ig,iq) 2025 2012 ENDDO ! (ig) 2026 ENDDO ! (iq) 2027 end if 2028 DO l=1,nlayer 2013 ENDDO ! (iq) 2014 DO l=1,nlayer 2029 2015 DO ig=1,ngrid 2030 2016 pdt(ig,l)=pdt(ig,l)+zdtc(ig,l) … … 2032 2018 pdu(ig,l)=pdu(ig,l)+zduc(ig,l) 2033 2019 ENDDO 2034 ENDDO2035 DO ig=1,ngrid2020 ENDDO 2021 DO ig=1,ngrid 2036 2022 zdtsurf(ig) = zdtsurf(ig) + zdtsurfc(ig) 2037 ENDDO2023 ENDDO 2038 2024 2039 2025 IF (tracer) THEN
Note: See TracChangeset
for help on using the changeset viewer.
