Changeset 1958 for trunk/LMDZ.TITAN/libf
- Timestamp:
- Jul 3, 2018, 12:50:23 PM (6 years ago)
- Location:
- trunk/LMDZ.TITAN/libf
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.TITAN/libf/chimtitan/gptitan.c
r1956 r1958 12 12 char CORPS[][10], double Y[][NLEV], 13 13 double *FIN, double *LAT, double *MASS, double MD[][NLEV], 14 double *KEDD, double *botCH4, double *fluxCH4, doubleKRATE[][NLEV],14 double *KEDD, double KRATE[][NLEV], 15 15 int reactif[][5], int *nom_prod, int *nom_perte, 16 16 int prod[][200], int perte[][200][2], int *aerprod, int *utilaer, … … 22 22 int i,j,k,l; 23 23 int ireac,ncom1,ncom2; 24 int in2, i ch4, ih, ih2, in4s;24 int in2, ih, ih2, in4s; 25 25 double ***a,***b,**c; 26 26 double *fl,*fp,*mu,**jac,**ym1,**f; 27 double dyCH4surf;28 27 double conv,delta,deltamax; 29 28 double cm,cp,dim,dip,dm,dp,dym,dyp,km,kp,r,dra,dram,drap; … … 100 99 for( i = 0; i <= ST-1; i++ ) 101 100 { 102 if( strcmp(corps[i], "CH4") == 0 ) ich4 = i;103 101 if( strcmp(corps[i], "H" ) == 0 ) ih = i; 104 102 if( strcmp(corps[i], "H2" ) == 0 ) ih2 = i; … … 107 105 } 108 106 109 /* initialisation mu , CH4 au sol*/107 /* initialisation mu */ 110 108 111 109 for( j = 0; j <= NLEV-1; j++ ) … … 114 112 for( i = 0; i <= ST-1; i++ ) 115 113 { 116 if( ( i == ich4 ) && ( Y[i][j] <= *botCH4 ) && ( j == 0 ) )117 {118 dyCH4surf = (*botCH4 - Y[i][j]);119 Y[i][j] = *botCH4;120 }121 114 mu[j] += ( MASS[i] * Y[i][j] ); 122 115 } … … 789 782 /* end inversion by blocs: */ 790 783 791 /* CH4 au sol */792 /* ---------- */793 if( Y[ich4][0] < *botCH4 )794 {795 dyCH4surf += (*botCH4 - Y[ich4][0]*10.0e0);796 Y[ich4][0] = *botCH4;797 }798 799 784 /* ------------------ */ 800 785 /* Tests et evolution */ … … 888 873 /* ------------ */ 889 874 890 *fluxCH4=dyCH4surf;891 892 875 for( j = 0; j <= NLEV-1; j++ ) 893 876 { -
trunk/LMDZ.TITAN/libf/phytitan/calchim.F90
r1956 r1958 1 1 SUBROUTINE calchim(ngrid,qy_c,declin,dtchim, & 2 ctemp,cpphi,cplay,cplev,czlay,czlev,dqyc ,zdyevapCH4)2 ctemp,cpphi,cplay,cplev,czlay,czlev,dqyc) 3 3 4 4 !--------------------------------------------------------------------------------------------------------- … … 86 86 87 87 REAL*8, DIMENSION(ngrid,klev,nkim), INTENT(OUT) :: dqyc ! Chemical species tendencies on GCM layers (mol/mol/s). 88 REAL*8, DIMENSION(ngrid), INTENT(OUT) :: zdyevapCH4 ! Diagnostic surface methane pseudo-evaporation flux (mol/mol/s).89 88 90 89 ! Local variables : … … 117 116 REAL*8, DIMENSION(nlaykim_tot) :: rinter ! Inter-layer distance (km) to planetographic center (RA grid in chem. module). 118 117 ! NB : rinter is on nlaykim_tot too, we don't care of the uppermost layer upper boundary altitude. 119 120 REAL(c_double) :: fluxCH4 ! Surface "evaporation" flux (mol/mol)121 118 122 119 ! Saved variables initialized at firstcall … … 469 466 nomqy_c,cqy, & 470 467 dtchim,latitude(ig)*180./pi,mass,md, & 471 kedd, botCH4,fluxCH4,krate,reactif,&468 kedd,krate,reactif, & 472 469 nom_prod,nom_perte,prod,perte, & 473 470 aerprod,utilaer,cmaer,cprodaer,ccsn,ccsh, & 474 471 htoh2,surfhaze) 475 476 zdyevapCH4(ig) = fluxCH4 / dtchim ! Diagnostic pseudo-evaporation ( due to readjustement to botCH4 value ) (mol/mol/s)477 472 478 473 ! 5. Calculates tendencies on composition for advected tracers … … 495 490 496 491 ELSE ! In 2D chemistry, if following grid point at same latitude, same zonal mean so don't do calculations again ! 497 zdyevapCH4(ig) = zdyevapCH4(igm1)498 492 dqyc(ig,:,:) = dqyc(igm1,:,:) ! will be put back in 3D with longitudinal variations assuming same relative tendencies within a lat band 499 493 ykim_up(:,ig,:) = ykim_up(:,igm1,:) ! no horizontal mixing in upper layers -> no longitudinal variations -
trunk/LMDZ.TITAN/libf/phytitan/physiq_mod.F90
r1957 r1958 17 17 use radcommon_h, only: sigma, gzlat, gzlat_ig, Cmk, grav, BWNV 18 18 use surfdat_h, only: phisfi, zmea, zstd, zsig, zgam, zthe 19 use comchem_h, only: nkim, cnames, nlaykim_up, ykim_up, ykim_tot 19 use comchem_h, only: nkim, cnames, nlaykim_up, ykim_up, ykim_tot, botCH4 20 20 use comdiurn_h, only: coslat, sinlat, coslon, sinlon 21 21 use comsaison_h, only: mu0, fract, dist_star, declin, right_ascen … … 287 287 288 288 real zdqchi(ngrid,nlayer,nq) ! Chemical tendency ( chemistry routine ). 289 real zdqcond(ngrid,nlayer,nq) ! Condensation tendency ( chemistry routine ).290 289 291 290 real zdqmufi(ngrid,nlayer,nq) ! Microphysical tendency. … … 376 375 real, dimension(ngrid,nlayer,nkim) :: ychimbar ! For 2D chemistry 377 376 378 ! Molar fraction tendencies ( chemistry and condensation ) for tracers (mol/mol/s) 379 real, dimension(ngrid,nlayer,nq) :: dyccond ! All tracers, we want to use indx on it. 380 real, dimension(ngrid,nlayer,nq) :: dyccondbar ! For 2D chemistry 381 real, dimension(:,:,:), allocatable, save :: dycchi ! Only for chem tracers. Saved since chemistry is not called every step. 377 ! Molar fraction tendencies ( chemistry, condensation and evaporation ) for tracers (mol/mol/s) 378 real, dimension(:,:,:), allocatable, save :: dycchi ! NB : Only for chem tracers. Saved since chemistry is not called every step. 382 379 !$OMP THREADPRIVATE(dycchi) 380 real, dimension(ngrid,nlayer,nq) :: dyccond ! Condensation flux. NB : for all tracers, as we want to use indx on it. 381 real, dimension(ngrid,nlayer,nq) :: dyccondbar ! For 2D chemistry 382 real, dimension(ngrid) :: dycevapCH4 ! Surface "pseudo-evaporation" flux (forcing constant surface humidity). 383 383 384 384 ! Saturation profiles … … 389 389 ! Surface methane 390 390 real, dimension(:), allocatable, save :: tankCH4 ! Depth of surface methane tank (m) 391 real, dimension(:), allocatable, save :: zdyevapCH4 ! Surface pseudo-evaporation flux (chemistry keeping constant surface humidity) (mol/mol/s). 392 !$OMP THREADPRIVATE(tankCH4,zdyevapCH4) 391 !$OMP THREADPRIVATE(tankCH4) 392 393 393 394 394 ! -----******----- FOR MUPHYS OPTICS -----******----- … … 523 523 allocate(dycchi(ngrid,nlayer,nkim)) ! only for chemical tracers 524 524 allocate(qysat(nlayer,nkim)) 525 allocate(zdyevapCH4(ngrid))526 525 527 526 ! Chemistry timestep … … 540 539 541 540 zdqchi(:,:,:) = 0.0 542 zdqcond(:,:,:) = 0.0543 541 544 542 endif … … 1146 1144 endif 1147 1145 1148 ! i. Condensation after the transport1149 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 1146 ! i. Condensation of the 3D tracers after the transport 1147 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 1150 1148 do iq=1,nkim 1151 1149 do l=1,nlayer … … 1164 1162 do iq=1,nkim 1165 1163 ychim(:,:,iq) = ychim(:,:,iq) + dyccond(:,:,iq+nmicro) ! update molar ychim for following calchim 1166 zdqcond(:,:,iq+nmicro) = dyccond(:,:,iq+nmicro)*rat_mmol(iq+nmicro) ! convert tendencies to mass mixing ratio 1164 1165 pdq(:,:,iq+nmicro) = pdq(:,:,iq+nmicro) + dyccond(:,:,iq+nmicro)*rat_mmol(iq+nmicro) ! convert tendencies to mass mixing ratio 1167 1166 enddo 1168 1167 1169 pdq(:,:,:) = pdq(:,:,:) + zdqcond(:,:,:) 1170 1171 ! 2D zonally averaged fields needed to condense before photochem1168 1169 ! ii. 2D zonally averaged fields need to condense and evap before photochem 1170 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 1172 1171 if ( moyzon_ch .and. mod(icount-1,ichim).eq. 0 ) then 1172 1173 1173 do iq = 1,nkim 1174 1174 do l=1,nlayer … … 1189 1189 enddo 1190 1190 1191 ! Pseudo-evap ( forcing constant surface humidity ) 1192 do ig=1,ngrid 1193 if ( ychimbar(ig,1,7+nmicro) .lt. botCH4 ) ychimbar(ig,1,7+nmicro) = botCH4 1194 enddo 1195 1191 1196 endif ! if ( moyzon_ch .and. mod(icount-1,ichim).eq. 0 ) 1192 1197 1193 ! ii . Photochemistry ( must be call after condensation)1194 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 1198 ! iii. Photochemistry ( must be call after condensation (and evap of 2D) ) 1199 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 1195 1200 if( mod(icount-1,ichim).eq.0. ) then 1196 1201 … … 1201 1206 ! Here we send zonal average fields ( corrected with cond ) from dynamics to chem. module 1202 1207 call calchim(ngrid,ychimbar,declin,ctimestep,ztfibar,zphibar, & 1203 zplaybar,zplevbar,zzlaybar,zzlevbar,dycchi ,zdyevapCH4)1208 zplaybar,zplevbar,zzlaybar,zzlevbar,dycchi) 1204 1209 else ! 3D chemistry (or 1D run) 1205 1210 call calchim(ngrid,ychim,declin,ctimestep,pt,pphi, & 1206 pplay,pplev,zzlay,zzlev,dycchi ,zdyevapCH4)1211 pplay,pplev,zzlay,zzlev,dycchi) 1207 1212 endif ! if moyzon 1208 1213 1209 1214 endif 1210 1215 1211 ! Add diagnostic-only surface pseudo-evapoaration in condensation tendency for bottom layer 1212 zdqcond(:,1,chimi_indx(7)) = zdyevapCH4(:)*rat_mmol(chimi_indx(7)) 1213 1216 ! iv. Surface pseudo-evaporation 1217 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 1218 do ig=1,ngrid 1219 if ( (ychim(ig,1,7+nmicro)+dycchi(ig,1,7+nmicro)*ptimestep) .lt. botCH4 ) then ! +dycchi because ychim not yet updated 1220 dycevapCH4(ig) = ( -ychim(ig,1,7+nmicro)+botCH4 ) / ptimestep - dycchi(ig,1,7+nmicro) 1221 else 1222 dycevapCH4(ig) = 0.0 1223 endif 1224 enddo 1225 1226 pdq(:,1,7+nmicro) = pdq(:,1,7+nmicro) + dycevapCH4(:)*rat_mmol(7+nmicro) 1227 1228 ! v. Updates and positivity check 1229 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 1214 1230 do iq=1,nkim 1215 1231 zdqchi(:,:,iq+nmicro) = dycchi(:,:,iq)*rat_mmol(iq+nmicro) ! convert tendencies to mass mixing ratio … … 1564 1580 call writediagfi(ngrid,cnames(iq),cnames(iq),'mol/mol',3,zq(:,:,iq+nmicro)/rat_mmol(iq+nmicro)) 1565 1581 enddo 1566 call writediagfi(ngrid,"fluxCH4","Surface CH4 pseudo-evaporation",'mol/mol/s',2, zdyevapCH4)1582 call writediagfi(ngrid,"fluxCH4","Surface CH4 pseudo-evaporation",'mol/mol/s',2,dycevapCH4) 1567 1583 endif 1568 1584 … … 1644 1660 1645 1661 ! Condensation tendencies ( mol/mol/s ) 1646 CALL send_xios_field("dqcond_CH4",zdqcond(:,:,7+nmicro)/rat_mmol(7+nmicro)) 1647 CALL send_xios_field("dqcond_C2H2",zdqcond(:,:,10+nmicro)/rat_mmol(10+nmicro)) 1648 CALL send_xios_field("dqcond_C2H4",zdqcond(:,:,12+nmicro)/rat_mmol(12+nmicro)) 1649 CALL send_xios_field("dqcond_C2H6",zdqcond(:,:,14+nmicro)/rat_mmol(14+nmicro)) 1650 CALL send_xios_field("dqcond_C3H6",zdqcond(:,:,17+nmicro)/rat_mmol(17+nmicro)) 1651 CALL send_xios_field("dqcond_C4H4",zdqcond(:,:,21+nmicro)/rat_mmol(21+nmicro)) 1652 CALL send_xios_field("dqcond_CH3CCH",zdqcond(:,:,24+nmicro)/rat_mmol(24+nmicro)) 1653 CALL send_xios_field("dqcond_C3H8",zdqcond(:,:,25+nmicro)/rat_mmol(25+nmicro)) 1654 CALL send_xios_field("dqcond_C4H2",zdqcond(:,:,26+nmicro)/rat_mmol(26+nmicro)) 1655 CALL send_xios_field("dqcond_C4H6",zdqcond(:,:,27+nmicro)/rat_mmol(27+nmicro)) 1656 CALL send_xios_field("dqcond_C4H10",zdqcond(:,:,28+nmicro)/rat_mmol(28+nmicro)) 1657 CALL send_xios_field("dqcond_AC6H6",zdqcond(:,:,29+nmicro)/rat_mmol(29+nmicro)) 1658 CALL send_xios_field("dqcond_HCN",zdqcond(:,:,36+nmicro)/rat_mmol(36+nmicro)) 1659 CALL send_xios_field("dqcond_CH3CN",zdqcond(:,:,40+nmicro)/rat_mmol(40+nmicro)) 1660 CALL send_xios_field("dqcond_HC3N",zdqcond(:,:,42+nmicro)/rat_mmol(42+nmicro)) 1661 CALL send_xios_field("dqcond_NCCN",zdqcond(:,:,43+nmicro)/rat_mmol(43+nmicro)) 1662 CALL send_xios_field("dqcond_C4N2",zdqcond(:,:,44+nmicro)/rat_mmol(44+nmicro)) 1662 CALL send_xios_field("dqcond_CH4",dyccond(:,:,7+nmicro)) 1663 CALL send_xios_field("dqcond_C2H2",dyccond(:,:,10+nmicro)) 1664 CALL send_xios_field("dqcond_C2H4",dyccond(:,:,12+nmicro)) 1665 CALL send_xios_field("dqcond_C2H6",dyccond(:,:,14+nmicro)) 1666 CALL send_xios_field("dqcond_C3H6",dyccond(:,:,17+nmicro)) 1667 CALL send_xios_field("dqcond_C4H4",dyccond(:,:,21+nmicro)) 1668 CALL send_xios_field("dqcond_CH3CCH",dyccond(:,:,24+nmicro)) 1669 CALL send_xios_field("dqcond_C3H8",dyccond(:,:,25+nmicro)) 1670 CALL send_xios_field("dqcond_C4H2",dyccond(:,:,26+nmicro)) 1671 CALL send_xios_field("dqcond_C4H6",dyccond(:,:,27+nmicro)) 1672 CALL send_xios_field("dqcond_C4H10",dyccond(:,:,28+nmicro)) 1673 CALL send_xios_field("dqcond_AC6H6",dyccond(:,:,29+nmicro)) 1674 CALL send_xios_field("dqcond_HCN",dyccond(:,:,36+nmicro)) 1675 CALL send_xios_field("dqcond_CH3CN",dyccond(:,:,40+nmicro)) 1676 CALL send_xios_field("dqcond_HC3N",dyccond(:,:,42+nmicro)) 1677 CALL send_xios_field("dqcond_NCCN",dyccond(:,:,43+nmicro)) 1678 CALL send_xios_field("dqcond_C4N2",dyccond(:,:,44+nmicro)) 1679 1680 ! Pseudo-evaporation flux (mol/mol/s) 1681 CALL send_xios_field("dqevapCH4",dycevapCH4(:)) 1663 1682 1664 1683 endif ! of 'if callchim'
Note: See TracChangeset
for help on using the changeset viewer.