Changeset 2322 for trunk/LMDZ.MARS/libf
- Timestamp:
- May 13, 2020, 5:14:02 PM (5 years ago)
- Location:
- trunk/LMDZ.MARS/libf
- Files:
-
- 6 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/iniphysiq_mod.F90
r1682 r2322 11 11 12 12 use infotrac, only : nqtot, & ! number of advected tracers 13 tname ! tracer names 13 tname, & ! tracer names 14 nqdesc,iqfils,nqperes,nqfils! MVals: variables isotopes 14 15 use comgeomfi_h, only: ini_fillgeom 15 16 use temps_mod, only: day_ini, hour_ini … … 71 72 call phys_state_var_init(klon_omp,nlayer,nqtot,tname, & 72 73 day_ini,hour_ini,punjours,ptimestep, & 73 prad,pg,pr,pcpp) 74 prad,pg,pr,pcpp, & 75 nqdesc,iqfils,nqperes,nqfils) ! MVals: variables isotopes 74 76 call ini_fillgeom(klon_omp,latitude,longitude,cell_area) 75 77 ! work is needed to put what is in comgeomfi_h in geometry_mod? -
trunk/LMDZ.MARS/libf/phymars/callsedim_mod.F
r2312 r2322 22 22 & igcm_co2_ice, igcm_stormdust_mass, 23 23 & igcm_stormdust_number,igcm_topdust_mass, 24 & igcm_topdust_number 24 & igcm_topdust_number, 25 & iqfils,nqfils,qperemin,masseqmin ! MVals: variables isotopes 25 26 USE newsedim_mod, ONLY: newsedim 26 27 USE comcstfi_h, ONLY: g … … 91 92 real masse (ngrid,nlay) ! Layer mass (kg.m-2) 92 93 real epaisseur (ngrid,nlay) ! Layer thickness (m) 93 real wq(ngrid,nlay+1) ! displaced tracer mass(kg.m-2)94 real wq(ngrid,nlay+1),w(ngrid,nlay) ! displaced tracer mass wq (kg.m-2), MVals: displaced "pere" tracer mass w (kg.m-2) 94 95 real r0(ngrid,nlay) ! geometric mean radius used for 95 96 ! sedimentation (m) … … 105 106 REAl ccntyp 106 107 character(len=20),parameter :: modname="callsedim" 108 c MVals: transport of the isotopic ratio 109 REAL Ratio0(ngrid,nlay),Ratio(ngrid,nlay) 110 REAL masseq(ngrid,nlay) 111 REAL newmasse 112 REAL zq0(ngrid,nlay,nq) 113 INTEGER ifils,iq2 107 114 108 115 … … 388 395 if(radius(iq).gt.1.e-9 .and.(iq.ne.ico2_ice) .and. 389 396 & (iq .ne. iccnco2_mass) .and. (iq .ne. 390 & iccnco2_number) ) then! no sedim for gaz or CO2 clouds (done in microtimestep)391 397 & iccnco2_number) .and. ! no sedim for gaz or CO2 clouds (done in microtimestep) 398 & iq .ne. igcm_hdo_ice) then !MVals: hdo is transported by h2o 392 399 c ----------------------------------------------------------------- 393 400 c DOUBLEQ CASE … … 500 507 c ----------------------------------------------------------------- 501 508 else if ((iq .eq. iccn_mass) .or. (iq .eq. iccn_number) 502 & .or. (iq .eq. igcm_h2o_ice) 503 & .or. (iq .eq. igcm_hdo_ice)) then 509 & .or. (iq .eq. igcm_h2o_ice)) then 504 510 if (microphys) then 505 ! water ice sedimentation 511 c water ice sedimentation 512 c ~~~~~~~~~~ 506 513 call newsedim(ngrid,nlay,ngrid*nlay,ngrid*nlay, 507 514 & ptimestep,pplev,masse,epaisseur,zt,rsedcloud,rhocloud, 508 515 & zqi(1,1,iq),wq,beta) 516 c Surface Tendencies 517 c ~~~~~~~~~~ 518 do ig=1,ngrid 519 pdqs_sed(ig,iq)=wq(ig,1)/ptimestep 520 end do 509 521 else 510 ! water ice sedimentation 522 c water ice sedimentation 523 c ~~~~~~~~~~ 511 524 call newsedim(ngrid,nlay,ngrid*nlay,1, 512 525 & ptimestep,pplev,masse,epaisseur,zt,rsedcloud,rho_q(iq), 513 526 & zqi(1,1,iq),wq,beta) 527 c Surface tendencies 528 c ~~~~~~~~~~ 529 do ig=1,ngrid 530 pdqs_sed(ig,iq)=wq(ig,1)/ptimestep 531 end do 532 c Special case of isotopes 533 c ~~~~~~~~~~ 534 !MVals: Loop over the sons ("fils") 535 if (nqfils(iq).gt.0) then 536 iq2=iqfils(nqfils(iq),iq) ! for now always nqfils(iq)=1 (special case of HDO only son of H2O) 537 !MVals: input paramters in vlz_fi for hdo 538 do l=1,nlay 539 do ig=1,ngrid 540 if (zq0(ig,l,iq).gt.qperemin) then 541 Ratio0(ig,l)=zq0(ig,l,iq2)/zq0(ig,l,iq) 542 else 543 Ratio0(ig,l)=0. 544 endif 545 Ratio(ig,l)=Ratio0(ig,l) 546 masseq(ig,l)=max(masse(ig,l)*zq0(ig,l,iq),masseqmin) 547 w(ig,l)=wq(ig,l) !MVals: very important: hdo is transported by h2o (see vlsplt_p.F: correction bugg 15 mai 2015) 548 enddo !ig=1,ngrid 549 enddo !l=1,nlay 550 !MVals: no need to enter newsedim as the transporting mass w has been already calculated 551 call vlz_fi(ngrid,nlay,Ratio,2.,masseq,w,wq) 552 zqi(:,nlay,iq2)=zqi(:,nlay,iq)*Ratio0(:,nlay) 553 do l=1,nlay-1 554 do ig=1,ngrid 555 newmasse=max(masseq(ig,l)+w(ig,l+1)-w(ig,l),masseqmin) 556 Ratio(ig,l)=(Ratio0(ig,l)*masseq(ig,l) 557 & +wq(ig,l+1)-wq(ig,l))/newmasse 558 zqi(ig,l,iq2)=zqi(ig,l,iq)*Ratio(ig,l) 559 enddo 560 enddo !l=1,nlay-1 561 !MVals: hdo surface tendency 562 do ig=1,ngrid 563 if (w(ig,1).gt.masseqmin) then 564 pdqs_sed(ig,iq2)=pdqs_sed(ig,iq)*(wq(ig,1)/w(ig,1)) 565 else 566 pdqs_sed(ig,iq2)=pdqs_sed(ig,iq)*Ratio0(ig,1) 567 endif 568 end do 569 endif !(nqfils(iq).gt.0) 514 570 endif ! of if (microphys) 515 c Tendencies 516 c ~~~~~~~~~~ 517 do ig=1,ngrid 518 pdqs_sed(ig,iq)=wq(ig,1)/ptimestep 519 end do 571 520 572 c ----------------------------------------------------------------- 521 573 c GENERAL CASE … … 539 591 pdqsed(ig,l,iq)=(zqi(ig,l,iq)- 540 592 $ (pq(ig,l,iq) + pdqfi(ig,l,iq)*ptimestep))/ptimestep 593 !MVals: Special case of isotopes: for now only HDO 594 if (nqfils(iq).gt.0) then 595 iq2=iqfils(nqfils(iq),iq) 596 pdqsed(ig,l,iq2)=(zqi(ig,l,iq2)- 597 $ (pq(ig,l,iq2) + pdqfi(ig,l,iq2)*ptimestep))/ptimestep 598 endif 541 599 ENDDO 542 600 ENDDO -
trunk/LMDZ.MARS/libf/phymars/co2condens_mod.F
r2189 r2322 18 18 use tracer_mod, only: noms, igcm_h2o_ice, 19 19 & igcm_dust_mass, igcm_dust_number, 20 & igcm_ccn_mass, igcm_ccn_number 20 & igcm_ccn_mass, igcm_ccn_number, 21 & iqfils,nqperes,nqfils, ! MVals: variables isotopes 22 & qperemin,masseqmin 21 23 use surfdat_h, only: emissiv, phisfi 22 24 use geometry_mod, only: latitude, ! grid point latitudes (rad) … … 155 157 REAL ztm3D(ngrid,nlayer) 156 158 REAL zmflux3D(ngrid,nlayer) 159 160 c MVals: variables isotopes 161 REAL Ratio1(nlayer),Ratiom1(nlayer+1) 162 REAL masseq(nlayer),wq(nlayer+1) 163 INTEGER ifils,iq2 157 164 c---------------------------------------------------------------------- 158 165 … … 537 544 call vl1d(nlayer,zu ,2.,masse,w,zum) 538 545 call vl1d(nlayer,zv ,2.,masse,w,zvm) 539 do iq=1,nq 546 ! MVals: loop over the fathers ("peres") 547 do iq=1,nqperes 540 548 do l=1,nlayer 541 549 zq1(l)=zqc(l,iq) … … 547 555 zqm(l,iq)=zqm1(l) 548 556 enddo 549 enddo 557 ! MVals: loop over the sons ("fils") 558 if (nqfils(iq).gt.0) then 559 iq2=iqfils(nqfils(iq),iq) ! for now it is always nqfils(iq)=1 (special case of HDO only son of H2O) 560 do l=1,nlayer 561 if (zqc(l,iq).gt.qperemin) then 562 Ratio1(l)=zqc(l,iq2)/zqc(l,iq) 563 else 564 Ratio1(l)=0. 565 endif 566 masseq(l)=max(masse(l)*zqc(l,iq),masseqmin) 567 wq(l)=w(l)*zqm(l,iq) 568 enddo 569 Ratiom1(1)=zqm(1,iq2) 570 call vl1d(nlayer,Ratio1,2.,masseq,wq,Ratiom1) 571 zqm(1,iq2)=Ratiom1(1)*zqc(1,iq) 572 do l=2,nlayer 573 zqm(l,iq2)=Ratiom1(l)*zqm(l,iq) 574 enddo 575 endif !if (nqfils(iq).gt.0) then 576 enddo !iq=1,nqperes 550 577 551 578 c Surface condensation affects low winds -
trunk/LMDZ.MARS/libf/phymars/dyn1d/testphys1d.F
r2312 r2322 5 5 use dimphy, only : init_dimphy 6 6 use mod_grid_phy_lmdz, only : regular_lonlat 7 use infotrac, only: nqtot, tname 7 use infotrac, only: nqtot, tname, nqperes,nqdesc,iqfils,nqfils, 8 & iqpere, nqdesc_tot 8 9 use comsoil_h, only: volcapa, layer, mlayer, inertiedat, nsoilmx 9 10 use comgeomfi_h, only: sinlat, ini_fillgeom … … 127 128 REAL halfaxe, excentric, Lsperi 128 129 Logical paleomars 130 131 c MVals: isotopes as in the dynamics (CRisi) 132 INTEGER :: ifils,ipere,generation 133 CHARACTER(len=30), ALLOCATABLE, DIMENSION(:) :: tnom_transp ! transporting fluid short name 134 CHARACTER(len=80) :: line ! to store a line of text 135 INTEGER ierr0 136 LOGICAL :: continu 129 137 130 138 c======================================================================= … … 240 248 allocate(dqdyn(nlayer,nq)) 241 249 allocate(mqtot(nq)) 250 allocate(tnom_transp(nq)) 242 251 243 252 ! read tracer names from file traceur.def 244 253 do iq=1,nq 245 read(90, *,iostat=ierr) tname(iq)254 read(90,'(80a)',iostat=ierr) line ! store the line from traceur.def 246 255 if (ierr.ne.0) then 247 256 write(*,*) 'testphys1d: error reading tracer names...' 248 257 stop 249 258 endif 259 ! if format is tnom_0, tnom_transp (isotopes) 260 read(line,*,iostat=ierr0) tname(iq),tnom_transp(iq) 261 if (ierr0.ne.0) then 262 read(line,*) tname(iq) 263 tnom_transp(iq)='air' 264 endif 265 250 266 enddo 251 267 close(90) 268 269 ! Isotopes: as in the 3D case we have to determine father/son relations for isotopes and carrying fluid 270 ALLOCATE(nqfils(nqtot),nqdesc(nqtot)) 271 ALLOCATE(iqfils(nqtot,nqtot)) 272 ALLOCATE(iqpere(nqtot)) 273 nqperes=0 274 nqfils(:)=0 275 nqdesc(:)=0 276 iqfils(:,:)=0 277 iqpere(:)=0 278 nqdesc_tot=0 279 DO iq=1,nqtot 280 if (tnom_transp(iq) == 'air') then 281 ! ceci est un traceur père 282 WRITE(*,*) 'Le traceur',iq,', appele ', 283 & trim(tname(iq)),', est un pere' 284 nqperes=nqperes+1 285 iqpere(iq)=0 286 else !if (tnom_transp(iq) == 'air') then 287 ! ceci est un fils. Qui est son père? 288 WRITE(*,*) 'Le traceur',iq,', appele ', 289 & trim(tname(iq)),', est un fils' 290 continu=.true. 291 ipere=1 292 do while (continu) 293 if (tnom_transp(iq) .eq. tname(ipere)) then 294 ! Son père est ipere 295 WRITE(*,*) 'Le traceur',iq,'appele ', 296 & trim(tname(iq)),' est le fils de ', 297 & ipere,'appele ',trim(tname(ipere)) 298 nqfils(ipere)=nqfils(ipere)+1 299 iqfils(nqfils(ipere),ipere)=iq 300 iqpere(iq)=ipere 301 continu=.false. 302 else !if (tnom_transp(iq) == tnom_0(ipere)) then 303 ipere=ipere+1 304 if (ipere.gt.nqtot) then 305 WRITE(*,*) 'Le traceur',iq,'appele ', 306 & trim(tname(iq)),', est orpelin.' 307 CALL abort_gcm('infotrac_init', 308 & 'Un traceur est orphelin',1) 309 endif !if (ipere.gt.nqtot) then 310 endif !if (tnom_transp(iq) == tnom_0(ipere)) then 311 enddo !do while (continu) 312 endif !if (tnom_transp(iq) == 'air') then 313 enddo !DO iq=1,nqtot 314 WRITE(*,*) 'nqperes=',nqperes 315 WRITE(*,*) 'nqfils=',nqfils 316 WRITE(*,*) 'iqpere=',iqpere 317 WRITE(*,*) 'iqfils=',iqfils 318 ! Calculer le nombre de descendants à partir de iqfils et de nbfils 319 DO iq=1,nqtot 320 generation=0 321 continu=.true. 322 ifils=iq 323 do while (continu) 324 ipere=iqpere(ifils) 325 if (ipere.gt.0) then 326 nqdesc(ipere)=nqdesc(ipere)+1 327 nqdesc_tot=nqdesc_tot+1 328 iqfils(nqdesc(ipere),ipere)=iq 329 ifils=ipere 330 generation=generation+1 331 else !if (ipere.gt.0) then 332 continu=.false. 333 endif !if (ipere.gt.0) then 334 enddo !do while (continu) 335 WRITE(*,*) 'Le traceur ',iq,', appele ',trim(tname(iq)), 336 & ' est un traceur de generation: ',generation 337 ENDDO !DO iq=1,nqtot 338 WRITE(*,*) 'infotrac: nqdesc=',nqdesc 339 WRITE(*,*) 'iqfils=',iqfils 340 WRITE(*,*) 'nqdesc_tot=',nqdesc_tot 252 341 253 342 ! initialize tracers here: … … 583 672 call init_dimphy(1,nlayer) ! Initialize dimphy module 584 673 call phys_state_var_init(1,llm,nq,tname, 585 . day0,time,daysec,dtphys,rad,g,r,cpp) 674 . day0,time,daysec,dtphys,rad,g,r,cpp, 675 . nqdesc,iqfils,nqperes,nqfils)! MVals: variables isotopes 586 676 call ini_fillgeom(1,latitude,longitude,(/1.0/)) 587 677 call conf_phys(1,llm,nq) -
trunk/LMDZ.MARS/libf/phymars/phys_state_var_init_mod.F90
r2265 r2322 5 5 SUBROUTINE phys_state_var_init(ngrid,nlayer,nq,tname, & 6 6 day_ini,hour_ini,pdaysec,ptimestep, & 7 prad,pg,pr,pcpp) 7 prad,pg,pr,pcpp, & 8 dyn_nqdesc,dyn_iqfils,dyn_nqperes,dyn_nqfils)! MVals: variables isotopes 8 9 9 10 !======================================================================= … … 67 68 REAL,INTENT(IN) :: hour_ini 68 69 REAL,INTENT(IN) :: pdaysec,ptimestep,prad,pg,pr,pcpp 70 !MVals isotopes 71 INTEGER,INTENT(in) :: dyn_nqperes 72 INTEGER,INTENT(in) :: dyn_nqfils(nq) 73 INTEGER,INTENT(in) :: dyn_nqdesc(nq) 74 INTEGER,INTENT(in) :: dyn_iqfils(nq,nq) 69 75 70 76 ! set dimension and allocate arrays in tracer_mod 71 77 call end_tracer_mod 72 call ini_tracer_mod(nq,tname )78 call ini_tracer_mod(nq,tname,dyn_nqdesc,dyn_iqfils,dyn_nqperes,dyn_nqfils)! MVals: variables isotopes 73 79 74 80 -
trunk/LMDZ.MARS/libf/phymars/tracer_mod.F90
r2321 r2322 95 95 ! other tracers 96 96 integer,save :: igcm_ar_n2 ! for simulations using co2 +neutral gas 97 97 ! MVals: isotopes 98 integer, save :: nqperes ! numbers of tracers defined as "peres" 99 integer, allocatable, save :: nqfils(:) ! numbers of sons ("fils") of the considered tracer 100 integer, allocatable, save :: iqfils(:,:) ! indice of a son, ex: iqfils(nqfils(ipere),ipere) 101 real, parameter :: qperemin=1.e-16 ! threschold for the "pere" mixing ratio qpere to calculate Ratio=qfils/qpere 102 real, parameter :: masseqmin=1.e-16 ! threschold for the "pere" transporting masse 103 !integer, allocatable, save :: nqdesc(:) ! number of sons + all gran-sons over all generations: not useful for now in the martian case as there are no gran-sons 98 104 99 105 !----------------------------------------------------------------------- … … 101 107 contains 102 108 103 subroutine ini_tracer_mod(nq,tname )109 subroutine ini_tracer_mod(nq,tname,dyn_nqdesc,dyn_iqfils,dyn_nqperes,dyn_nqfils)! MVals: variables isotopes 104 110 implicit none 105 111 106 112 integer,intent(in) :: nq ! number of tracers 107 113 character(len=*),intent(in) :: tname(nq) ! tracer names 114 !MVals: variables isotopes 115 integer,intent(in) :: dyn_nqperes 116 integer,intent(in) :: dyn_nqfils(nq) 117 integer,intent(in) :: dyn_nqdesc(nq) 118 integer,intent(in) :: dyn_iqfils(nq,nq) 108 119 109 120 integer :: iq, count … … 117 128 write(*,*) "tracer_mod names : ", trim(noms(iq)) 118 129 enddo 130 131 !MVals: isotopes variables initialisation 132 do iq=1,nq 133 if (dyn_nqfils(iq).ne.dyn_nqdesc(iq)) then 134 write(*,*) ' for now all descendants must be sons: check the', & 135 ' relations between tracers in traceur.def !' 136 call abort_physic("ini_tracer_mod","relatives pattern between tracers not accepted",1) 137 endif 138 enddo 139 allocate(nqfils(nq))!,nqdesc(nq)) 140 allocate(iqfils(nq,nq)) 141 nqperes=dyn_nqperes 142 nqfils(:)=dyn_nqfils(:) 143 iqfils(:,:)=dyn_iqfils(:,:) 119 144 120 145 #ifndef MESOSCALE
Note: See TracChangeset
for help on using the changeset viewer.