- Timestamp:
- Jan 4, 2022, 9:04:10 AM (3 years ago)
- Location:
- trunk/LMDZ.MARS
- Files:
-
- 3 edited
-
README (modified) (1 diff)
-
libf/phymars/co2condens_mod.F (modified) (17 diffs)
-
libf/phymars/physiq_mod.F (modified) (3 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.MARS/README
r2597 r2599 3556 3556 Also cleaned up and commented comsaison_h in the process. 3557 3557 3558 == 04/01/2021 == CM 3559 - Clean co2condens_mod.F 3560 - remove dqsurf duplication after call co2condens -
trunk/LMDZ.MARS/libf/phymars/co2condens_mod.F
r2566 r2599 323 323 ENDDO 324 324 325 condens_layer(:,:) = zcondicea(:,:) 326 condens_column(:) = 0. 327 325 328 if (scavco2cond) then 326 329 call scavenging_by_co2(ngrid,nlayer,nq,ptimestep,pplev,zq, … … 339 342 DO ig=1,ngrid 340 343 condens_column(ig) = sum(condens_layer(ig,:)) 344 zfallice(ig,1) = zdqssed_co2(ig) 341 345 piceco2(ig) = piceco2(ig) + zdqssed_co2(ig)*ptimestep 342 346 ENDDO … … 347 351 & " ",3,pdtc) 348 352 349 call WRITEdiagfi(ngrid," zcondicea",353 call WRITEdiagfi(ngrid,"condens_layer", 350 354 & "", 351 & " ",3, zcondicea)355 & " ",3,condens_layer) 352 356 353 357 call WRITEdiagfi(ngrid,"zfallice", … … 374 378 c 375 379 c No microphysics of CO2 clouds 376 if (.not.co2clouds)then377 380 DO ig=1,ngrid 378 381 IF(latitude(ig).lt.0) THEN … … 383 386 icap=1 384 387 ENDIF 385 386 388 c 387 389 c Loop on where we have condensation/ sublimation … … 390 392 $ ((ztsrf(ig) .GT. ztcondsol(ig)) .AND. ! ground sublim. 391 393 $ ((piceco2(ig)+zfallice(ig,1)*ptimestep) .NE. 0.))) THEN 392 condsub(ig) = .true. 393 394 IF (zfallice(ig,1).gt.0) then 395 zfallheat=zfallice(ig,1)* 396 & (pphi(ig,1)- phisfi(ig) + 397 & cpice*(ztcond(ig,1)-ztcondsol(ig)))/latcond 394 condsub(ig) = .true. 395 396 IF (zfallice(ig,1).gt.0 .and. .not. co2clouds) then 397 zfallheat = zfallice(ig,1) * (pphi(ig,1)- phisfi(ig) + 398 & cpice*(ztcond(ig,1)-ztcondsol(ig)))/latcond 398 399 ELSE 399 zfallheat=0.400 zfallheat = 0. 400 401 ENDIF 401 402 c condensation or partial sublimation of CO2 ice 402 403 c """"""""""""""""""""""""""""""""""""""""""""""" 403 zcondices(ig) =pcapcal(ig)*(ztcondsol(ig)-ztsrf(ig))404 & /(latcond*ptimestep)- zfallheat404 zcondices(ig) = pcapcal(ig)*(ztcondsol(ig)-ztsrf(ig)) 405 & /(latcond*ptimestep) - zfallheat 405 406 pdtsrfc(ig) = (ztcondsol(ig) - ztsrf(ig))/ptimestep 406 zdiceco2(ig) = zcondices(ig) + zfallice(ig,1)407 407 #ifdef MESOSCALE 408 408 print*, "not enough CO2 tracer in 1st layer to condense" … … 413 413 c """""""""""""""""""""""""""""""""""""""""""""""""""""" 414 414 IF(ico2.ne.0) then 415 c 415 c Available CO2 tracer in layer 1 at end of timestep (kg/m2) 416 416 #ifndef MESOSCALE 417 availco2= pq(ig,1,ico2)*((ap(1)-ap(2))+ 418 & (bp(1)-bp(2))*(pplev(ig,1)/g-zdiceco2(ig)*ptimestep)) 417 availco2 = pq(ig,1,ico2)*((ap(1)-ap(2))+ 418 & (bp(1)-bp(2))*(pplev(ig,1)/g - 419 & (zcondices(ig) + zfallice(ig,1))*ptimestep)) 420 if ((zcondices(ig) + condens_layer(ig,1))*ptimestep 421 & .gt.availco2) then 422 zcondices(ig) = availco2/ptimestep - 423 & condens_layer(ig,1) 424 pdtsrfc(ig) = (latcond/pcapcal(ig))* 425 & (zcondices(ig)+zfallheat) 426 end if 419 427 #else 420 428 availco2 = pq(ig,1,igcm_co2) … … 422 430 & ' CORRECTED FROM SIGMA LEVELS" 423 431 #endif 424 IF ((zcondices(ig) + zcondicea(ig,1))*ptimestep425 & .gt.availco2) then426 zcondices(ig) = availco2/ptimestep -zcondicea(ig,1)427 zdiceco2(ig) = zcondices(ig) + zfallice(ig,1)428 pdtsrfc(ig)=(latcond/pcapcal(ig))*429 & (zcondices(ig)+zfallheat)430 ENDIF431 432 ENDIF 432 433 #endif … … 437 438 IF((piceco2(ig)/ptimestep+zfallice(ig,1)).LE. 438 439 & -zcondices(ig))THEN 439 zcondices(ig) = -piceco2(ig)/ptimestep - zfallice(ig,1) 440 zcondices(ig) = -piceco2(ig)/ptimestep - zfallice(ig,1) 440 441 pdtsrfc(ig)=(latcond/pcapcal(ig))* 441 442 & (zcondices(ig)+zfallheat) 442 zdiceco2(ig) = zcondices(ig) + zfallice(ig,1)443 443 END IF 444 444 445 445 c Changing CO2 ice amount and pressure : 446 446 c """""""""""""""""""""""""""""""""""" 447 zdiceco2(ig) = zcondices(ig) + zfallice(ig,1) 448 & + condens_column(ig) 449 if (co2clouds) then 450 ! add here only direct condensation/sublimation 451 piceco2(ig) = piceco2(ig) + zcondices(ig)*ptimestep 452 else 453 ! add here also CO2 ice in the atmosphere 447 454 piceco2(ig) = piceco2(ig) + zdiceco2(ig)*ptimestep 455 end if 448 456 pdpsrf(ig) = -zdiceco2(ig)*g 449 457 … … 461 469 END IF ! if there is condensation/sublimation 462 470 ENDDO ! of DO ig=1,ngrid 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 471 523 472 c ******************************************************************** 524 473 c Surface albedo and emissivity of the surface below the snow (emisref) … … 565 514 c Mass fluxes through the sigma levels (kg.m-2.s-1) (>0 when up) 566 515 c """""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" 567 if (.not. co2clouds) then 568 zmflux(1) = -zcondices(ig) 516 zmflux(1) = -zcondices(ig) - zdqssed_co2(ig) 569 517 DO l=1,nlayer 570 zmflux(l+1) = zmflux(l) - zcondicea(ig,l)518 zmflux(l+1) = zmflux(l) - condens_layer(ig,l) 571 519 #ifndef MESOSCALE 572 & + (bp(l)-bp(l+1))*( zfallice(ig,1)-zmflux(1))520 & + (bp(l)-bp(l+1))*(-pdpsrf(ig)/g) 573 521 c zmflux set to 0 if very low to avoid: top layer is disappearing in v1ld 574 522 if (abs(zmflux(l+1)).lt.1E-13.OR.bp(l+1).eq.0.) then … … 576 524 end if 577 525 #else 578 zmflux(l+1) = zmflux(l) - zcondicea(ig,l)526 zmflux(l+1) = zmflux(l) - condens_layer(ig,l) 579 527 if (abs(zmflux(l+1)).lt.1E-13) zmflux(l+1)=0. 580 528 PRINT*, "MESOSCALE: FLUX THROUGH SIGMA LEVELS from"// … … 582 530 #endif 583 531 END DO 584 else585 zmflux(1) = - zcondices(ig) - zdqssed_co2(ig)586 do l = 1, nlayer587 #ifndef MESOSCALE588 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 v1ld591 if (abs(zmflux(l+1))<1E-13.OR.bp(l+1)==0.) then592 zmflux(l+1) = 0.593 end if594 #else595 zmflux(l+1) = zmflux(l) - condens_layer(ig,l)596 if (abs(zmflux(l+1))<1E-13) then597 zmflux(l+1) = 0.598 end if599 PRINT*, "MESOSCALE: FLUX THROUGH SIGMA LEVELS from"//600 & "dPS HAVE TO BE IMPLEMENTED"601 #endif602 end do603 end if604 532 #ifdef MESOSCALE 605 533 print*,"absurd mass set because bp not available" … … 705 633 & ( zmflux(l)*(ztm(l) - ztc(l)) 706 634 & - zmflux(l+1)*(ztm(l+1) - ztc(l)) 707 & + zcondicea(ig,l)*(ztcond(ig,l)-ztc(l)) )635 & + condens_layer(ig,l)*(ztcond(ig,l)-ztc(l)) ) 708 636 ELSE 709 637 zdtsig(ig,l) = (1/masse(l)) * … … 733 661 #endif 734 662 735 if (.not. co2clouds) then736 663 do iq=1,nq 737 664 ! if (noms(iq).eq.'co2') then 738 if (iq.eq.ico2) then739 c SPECIAL Case when the tracer IS CO2 :740 DO l=1,nlayer741 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 DO746 else747 DO l=1,nlayer748 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=false754 END DO755 end if756 enddo757 else758 c Tendencies on Q759 do iq=1,nq760 ! if (noms(iq).eq.'co2') then761 665 if (iq.eq.ico2) then 762 666 c SPECIAL Case when the tracer IS CO2 : … … 770 674 DO l=1,nlayer 771 675 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) ) 676 & ( zmflux(l)*(zqm(l,iq) - zqc(l,iq)) 677 & - zmflux(l+1)*(zqm(l+1,iq) - zqc(l,iq)) 678 & + condens_layer(ig,l)*zqc(l,iq) ) 679 680 pdqc(ig,l,iq)=pdqc(ig,l,iq)+zdq_scav(ig,l,iq) ! zdq_scav=0 if watercloud=false 775 681 END DO 776 682 end if 777 683 enddo 778 end if779 780 684 781 685 end if ! if (condsub) … … 785 689 c CO2 snow / clouds scheme 786 690 c *************************************************************** 787 if (.not. co2clouds) then788 691 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 692 & condens_layer,zcondices,zfallice,pemisurf) 794 693 c *************************************************************** 795 694 c Ecriture des diagnostiques -
trunk/LMDZ.MARS/libf/phymars/physiq_mod.F
r2597 r2599 874 874 co2totA = co2totA + co2ice(ig) 875 875 end do 876 else 877 co2totA = 0. 878 do ig=1,ngrid 879 do l=1,nlayer 880 co2totA = co2totA + (zplev(ig,l)-zplev(ig,l+1))/g* 881 & (pq(ig,l,igcm_co2) 882 & +pdq(ig,l,igcm_co2)*ptimestep) 883 end do 884 co2totA = co2totA + co2ice(ig) 885 end do 876 886 endif ! of if (igcm_co2_ice.ne.0) 877 887 c----------------------------------------------------------------------- … … 2059 2069 ENDDO 2060 2070 ENDDO 2061 2062 DO iq=1, nq2063 DO ig=1,ngrid2064 dqsurf(ig,iq)=dqsurf(ig,iq)+zdqsc(ig,iq)2065 ENDDO ! (ig)2066 ENDDO ! (iq)2067 2068 2071 ENDIF ! of IF (tracer) 2069 2072 … … 3953 3956 co2totB = co2totB + co2ice(ig) 3954 3957 enddo 3955 3958 else 3959 co2totB = 0. ! added by C.M. 3960 do ig=1,ngrid 3961 do l=1,nlayer 3962 co2totB = co2totB + (zplev(ig,l)-zplev(ig,l+1))/g* 3963 & (pq(ig,l,igcm_co2)+pdq(ig,l,igcm_co2)*ptimestep) 3964 enddo 3965 co2totB = co2totB + co2ice(ig) 3966 enddo 3967 endif ! of if (igcm_co2_ice.ne.0) 3956 3968 call WRITEDIAGFI(ngrid,'co2conservation', 3957 3969 & 'Total CO2 mass conservation in physic', 3958 3970 & '%',0,[(co2totA-co2totB)/co2totA]) 3959 endif ! of if (igcm_co2_ice.ne.0)3960 3971 ! XIOS outputs 3961 3972 #ifdef CPP_XIOS
Note: See TracChangeset
for help on using the changeset viewer.
