Changeset 719 for trunk/LMDZ.MARS/libf
- Timestamp:
- Jul 12, 2012, 1:42:13 PM (12 years ago)
- Location:
- trunk/LMDZ.MARS/libf/phymars
- Files:
-
- 2 edited
-
physiq.F (modified) (13 diffs)
-
wstats.F90 (modified) (2 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.MARS/libf/phymars/physiq.F
r708 r719 322 322 ! Qabs instead of Qext 323 323 ! (direct comparison with TES) 324 325 REAL dqdustsurf(ngridmx) ! surface q dust flux (kg/m2/s) 326 REAL dndustsurf(ngridmx) ! surface n dust flux (number/m2/s) 324 327 325 328 c Test 1d/3d scavenging 326 329 real h2otot(ngridmx) 327 REAL satu(ngridmx,nlayermx) ! satu ratio for output328 REAL zqsat(ngridmx,nlayermx) 330 REAL satu(ngridmx,nlayermx) ! satu ratio for output 331 REAL zqsat(ngridmx,nlayermx) ! saturation 329 332 330 333 REAL time_phys … … 1474 1477 1475 1478 if (tracer) then 1479 1480 if(doubleq) then 1481 do ig=1,ngrid 1482 dqdustsurf(ig) = 1483 & dqsurf(ig,igcm_dust_mass)*tauscaling(ig) 1484 dndustsurf(ig) = 1485 & dqsurf(ig,igcm_dust_number)*tauscaling(ig) 1486 enddo 1487 if (scavenging) then 1488 do ig=1,ngrid 1489 dqdustsurf(ig) = dqdustsurf(ig) + 1490 & dqsurf(ig,igcm_ccn_mass)*tauscaling(ig) 1491 dndustsurf(ig) = dndustsurf(ig) + 1492 & dqsurf(ig,igcm_ccn_number)*tauscaling(ig) 1493 enddo 1494 endif 1495 endif 1496 1476 1497 if (water) then 1477 1478 1498 mtot(:)=0 1479 1499 icetot(:)=0 … … 1588 1608 if (water) then 1589 1609 vmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_vap) 1590 & *mugaz/mmol(igcm_h2o_vap)1591 call wstats(ngrid,"vmr_h2ovap or",1610 & *mmean(1:ngridmx,1:nlayermx)/mmol(igcm_h2o_vap) 1611 call wstats(ngrid,"vmr_h2ovap", 1592 1612 & "H2O vapor volume mixing ratio","mol/mol", 1593 1613 & 3,vmr) 1594 1614 vmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_ice) 1595 & *mugaz/mmol(igcm_h2o_ice)1615 & *mmean(1:ngridmx,1:nlayermx)/mmol(igcm_h2o_ice) 1596 1616 call wstats(ngrid,"vmr_h2oice", 1597 1617 & "H2O ice volume mixing ratio","mol/mol", 1598 1618 & 3,vmr) 1599 1619 vmr=zqsat(1:ngridmx,1:nlayermx) 1600 & *mugaz/mmol(igcm_h2o_vap)1620 & *mmean(1:ngridmx,1:nlayermx)/mmol(igcm_h2o_vap) 1601 1621 call wstats(ngrid,"vmr_h2osat", 1602 1622 & "saturation volume mixing ratio","mol/mol", … … 1605 1625 & "surface h2o_ice","kg/m2", 1606 1626 & 2,qsurf(1,igcm_h2o_ice)) 1607 call wstats(ngrid,'albedo',1608 & 'albedo',1609 & '',2,albedo(1:ngridmx,1))1627 c call wstats(ngrid,'albedo', 1628 c & 'albedo', 1629 c & '',2,albedo(1:ngridmx,1)) 1610 1630 call wstats(ngrid,"mtot", 1611 1631 & "total mass of water vapor","kg/m2", … … 1637 1657 1638 1658 endif ! of if (water) 1659 1660 1661 if (dustbin.ne.0) then 1662 1663 call wstats(ngridmx,'tau','taudust','SI',2,tau(1,1)) 1664 1665 if (doubleq) then 1666 c call wstats(ngridmx,'qsurf','qsurf', 1667 c & 'kg.m-2',2,qsurf(1,igcm_dust_mass)) 1668 c call wstats(ngridmx,'Nsurf','N particles', 1669 c & 'N.m-2',2,qsurf(1,igcm_dust_number)) 1670 c call wstats(ngridmx,'dqsdev','ddevil lift', 1671 c & 'kg.m-2.s-1',2,zdqsdev(1,1)) 1672 c call wstats(ngridmx,'dqssed','sedimentation', 1673 c & 'kg.m-2.s-1',2,zdqssed(1,1)) 1674 c call wstats(ngridmx,'dqsdif','diffusion', 1675 c & 'kg.m-2.s-1',2,zdqsdif(1,1)) 1676 call wstats(ngridmx,'dqsdust', 1677 & 'deposited surface dust mass', 1678 & 'kg.m-2.s-1',2,dqdustsurf) 1679 call wstats(ngridmx,'dqndust', 1680 & 'deposited surface dust number', 1681 & 'number.m-2.s-1',2,dndustsurf) 1682 call wstats(ngridmx,'reffdust','reffdust', 1683 & 'm',3,rdust*ref_r0) 1684 call wstats(ngridmx,'dustq','Dust mass mr', 1685 & 'kg/kg',3,pq(1,1,igcm_dust_mass)) 1686 call wstats(ngridmx,'dustN','Dust number', 1687 & 'part/kg',3,pq(1,1,igcm_dust_number)) 1688 else 1689 do iq=1,dustbin 1690 write(str2(1:2),'(i2.2)') iq 1691 call wstats(ngridmx,'q'//str2,'mix. ratio', 1692 & 'kg/kg',3,zq(1,1,iq)) 1693 call wstats(ngridmx,'qsurf'//str2,'qsurf', 1694 & 'kg.m-2',2,qsurf(1,iq)) 1695 end do 1696 endif ! (doubleq) 1697 1698 if (scavenging) then 1699 call wstats(ngridmx,'ccnq','CCN mass mr', 1700 & 'kg/kg',3,pq(1,1,igcm_ccn_mass)) 1701 call wstats(ngridmx,'ccnN','CCN number', 1702 & 'part/kg',3,pq(1,1,igcm_ccn_number)) 1703 endif ! (scavenging) 1704 1705 endif ! (dustbin.ne.0) 1706 1639 1707 1640 1708 if (thermochem.or.photochem) then … … 1643 1711 $ noms(iq) .ne. "dust_number" .and. 1644 1712 $ noms(iq) .ne. "ccn_mass" .and. 1645 $ noms(iq) .ne. "ccn_number") then 1713 $ noms(iq) .ne. "ccn_number" .and. 1714 $ noms(iq) .ne. "h2o_vap" .and. 1715 $ noms(iq) .ne. "h2o_ice") then 1646 1716 vmr(1:ngrid,1:nlayer)=zq(1:ngrid,1:nlayer,iq) 1647 1717 & *mmean(1:ngrid,1:nlayer)/mmol(iq) … … 1710 1780 qsurfice(1:ngrid) = qsurf(1:ngrid,igcm_h2o_ice) 1711 1781 vmr=1.e6 * zq(1:ngrid,1:nlayer,igcm_h2o_ice) 1712 . * mugaz/ mmol(igcm_h2o_ice)1782 . *mmean(1:ngridmx,1:nlayermx) / mmol(igcm_h2o_ice) 1713 1783 ENDIF 1714 1784 !! Dust quantity integration along the vertical axe … … 1745 1815 & tsurf) 1746 1816 call WRITEDIAGFI(ngrid,"ps","surface pressure","Pa",2,ps) 1747 call WRITEDIAGFI(ngrid,"co2ice","co2 ice thickness","kg.m-2",2,1817 call WRITEDIAGFI(ngrid,"co2ice","co2 ice thickness","kg.m-2",2, 1748 1818 & co2ice) 1749 1819 1750 ccall WRITEDIAGFI(ngrid,"temp7","temperature in layer 7",1751 c& "K",2,zt(1,7))1752 ccall WRITEDIAGFI(ngrid,"fluxsurf_lw","fluxsurf_lw","W.m-2",2,1753 c& fluxsurf_lw)1754 ccall WRITEDIAGFI(ngrid,"fluxsurf_sw","fluxsurf_sw","W.m-2",2,1755 c& fluxsurf_sw_tot)1756 ccall WRITEDIAGFI(ngrid,"fluxtop_lw","fluxtop_lw","W.m-2",2,1757 c& fluxtop_lw)1758 ccall WRITEDIAGFI(ngrid,"fluxtop_sw","fluxtop_sw","W.m-2",2,1759 c& fluxtop_sw_tot)1760 ccall WRITEDIAGFI(ngrid,"temp","temperature","K",3,zt)1761 ccall WRITEDIAGFI(ngrid,"u","Zonal wind","m.s-1",3,zu)1762 ccall WRITEDIAGFI(ngrid,"v","Meridional wind","m.s-1",3,zv)1763 ccall WRITEDIAGFI(ngrid,"w","Vertical wind","m.s-1",3,pw)1764 ccall WRITEDIAGFI(ngrid,"rho","density","none",3,rho)1820 call WRITEDIAGFI(ngrid,"temp7","temperature in layer 7", 1821 & "K",2,zt(1,7)) 1822 call WRITEDIAGFI(ngrid,"fluxsurf_lw","fluxsurf_lw","W.m-2",2, 1823 & fluxsurf_lw) 1824 call WRITEDIAGFI(ngrid,"fluxsurf_sw","fluxsurf_sw","W.m-2",2, 1825 & fluxsurf_sw_tot) 1826 call WRITEDIAGFI(ngrid,"fluxtop_lw","fluxtop_lw","W.m-2",2, 1827 & fluxtop_lw) 1828 call WRITEDIAGFI(ngrid,"fluxtop_sw","fluxtop_sw","W.m-2",2, 1829 & fluxtop_sw_tot) 1830 call WRITEDIAGFI(ngrid,"temp","temperature","K",3,zt) 1831 call WRITEDIAGFI(ngrid,"u","Zonal wind","m.s-1",3,zu) 1832 call WRITEDIAGFI(ngrid,"v","Meridional wind","m.s-1",3,zv) 1833 call WRITEDIAGFI(ngrid,"w","Vertical wind","m.s-1",3,pw) 1834 call WRITEDIAGFI(ngrid,"rho","density","none",3,rho) 1765 1835 c call WRITEDIAGFI(ngrid,"q2","q2","kg.m-3",3,q2) 1766 1836 c call WRITEDIAGFI(ngrid,'Teta','T potentielle','K',3,zh) 1767 ccall WRITEDIAGFI(ngrid,"pressure","Pressure","Pa",3,pplay)1837 call WRITEDIAGFI(ngrid,"pressure","Pressure","Pa",3,pplay) 1768 1838 c call WRITEDIAGFI(ngrid,"ssurf","Surface stress","N.m-2",2, 1769 1839 c & zstress) … … 1808 1878 ! call WRITEDIAGFI(ngrid,"co2_l1","co2 mix. ratio in 1st layer", 1809 1879 ! & "kg/kg",2,zq(1,1,igcm_co2)) 1810 !call WRITEDIAGFI(ngrid,"co2","co2 mass mixing ratio",1811 !& "kg/kg",3,zq(1,1,igcm_co2))1880 call WRITEDIAGFI(ngrid,"co2","co2 mass mixing ratio", 1881 & "kg/kg",3,zq(1,1,igcm_co2)) 1812 1882 1813 1883 ! Compute co2 column … … 1844 1914 #endif 1845 1915 1846 1916 CALL WRITEDIAGFI(ngridmx,'mtot', 1847 1917 & 'total mass of water vapor', 1848 1918 & 'kg/m2',2,mtot) 1849 1919 CALL WRITEDIAGFI(ngridmx,'icetot', 1850 1920 & 'total mass of water ice', 1851 1921 & 'kg/m2',2,icetot) 1852 cvmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_ice)1853 c & *mugaz/mmol(igcm_h2o_ice)1854 ccall WRITEDIAGFI(ngridmx,'vmr_h2oice','h2o ice vmr',1855 c& 'mol/mol',3,vmr)1856 cvmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_vap)1857 c & *mugaz/mmol(igcm_h2o_vap)1858 ccall WRITEDIAGFI(ngridmx,'vmr_h2ovap','h2o vap vmr',1859 c& 'mol/mol',3,vmr)1860 1922 vmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_ice) 1923 & *mmean(1:ngridmx,1:nlayermx)/mmol(igcm_h2o_ice) 1924 call WRITEDIAGFI(ngridmx,'vmr_h2oice','h2o ice vmr', 1925 & 'mol/mol',3,vmr) 1926 vmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_vap) 1927 & *mmean(1:ngridmx,1:nlayermx)/mmol(igcm_h2o_vap) 1928 call WRITEDIAGFI(ngridmx,'vmr_h2ovap','h2o vap vmr', 1929 & 'mol/mol',3,vmr) 1930 CALL WRITEDIAGFI(ngridmx,'reffice', 1861 1931 & 'Mean reff', 1862 1932 & 'm',2,rave) 1863 1933 CALL WRITEDIAGFI(ngrid,"Nccntot", 1864 1934 & "condensation nuclei","Nbr/m2", 1865 1935 & 2,Nccntot) 1866 cCALL WRITEDIAGFI(ngrid,"Mccntot",1867 c& "mass condensation nuclei","kg/m2",1868 c& 2,Mccntot)1869 ccall WRITEDIAGFI(ngridmx,'rice','Ice particle size',1870 c& 'm',3,rice)1871 1936 CALL WRITEDIAGFI(ngrid,"Mccntot", 1937 & "mass condensation nuclei","kg/m2", 1938 & 2,Mccntot) 1939 call WRITEDIAGFI(ngridmx,'rice','Ice particle size', 1940 & 'm',3,rice) 1941 call WRITEDIAGFI(ngridmx,'h2o_ice_s', 1872 1942 & 'surface h2o_ice', 1873 1943 & 'kg.m-2',2,qsurf(1,igcm_h2o_ice)) … … 1902 1972 1903 1973 if (tracer.and.(dustbin.ne.0)) then 1904 ccall WRITEDIAGFI(ngridmx,'tau','taudust','SI',2,tau(1,1))1974 call WRITEDIAGFI(ngridmx,'tau','taudust','SI',2,tau(1,1)) 1905 1975 if (doubleq) then 1906 1976 c call WRITEDIAGFI(ngridmx,'qsurf','qsurf', … … 1914 1984 c call WRITEDIAGFI(ngridmx,'dqsdif','diffusion', 1915 1985 c & 'kg.m-2.s-1',2,zdqsdif(1,1)) 1916 c call WRITEDIAGFI(ngridmx,'reffdust','reffdust', 1917 c & 'm',3,rdust*ref_r0) 1918 c call WRITEDIAGFI(ngridmx,'dustq','Dust mass mr', 1919 c & 'kg/kg',3,pq(1,1,igcm_dust_mass)) 1920 c call WRITEDIAGFI(ngridmx,'dustN','Dust number', 1921 c & 'part/kg',3,pq(1,1,igcm_dust_number)) 1986 call WRITEDIAGFI(ngridmx,'dqsdust', 1987 & 'deposited surface dust mass', 1988 & 'kg.m-2.s-1',2,dqdustsurf) 1989 call WRITEDIAGFI(ngridmx,'dqndust', 1990 & 'deposited surface dust number', 1991 & 'number.m-2.s-1',2,dndustsurf) 1992 call WRITEDIAGFI(ngridmx,'reffdust','reffdust', 1993 & 'm',3,rdust*ref_r0) 1994 call WRITEDIAGFI(ngridmx,'dustq','Dust mass mr', 1995 & 'kg/kg',3,pq(1,1,igcm_dust_mass)) 1996 call WRITEDIAGFI(ngridmx,'dustN','Dust number', 1997 & 'part/kg',3,pq(1,1,igcm_dust_number)) 1922 1998 #ifdef MESOINI 1923 1999 call WRITEDIAGFI(ngridmx,'dustq','Dust mass mr', … … 1941 2017 1942 2018 if (scavenging) then 1943 ccall WRITEDIAGFI(ngridmx,'ccnq','CCN mass mr',1944 c& 'kg/kg',3,pq(1,1,igcm_ccn_mass))1945 ccall WRITEDIAGFI(ngridmx,'ccnN','CCN number',1946 c& 'part/kg',3,pq(1,1,igcm_ccn_number))2019 call WRITEDIAGFI(ngridmx,'ccnq','CCN mass mr', 2020 & 'kg/kg',3,pq(1,1,igcm_ccn_mass)) 2021 call WRITEDIAGFI(ngridmx,'ccnN','CCN number', 2022 & 'part/kg',3,pq(1,1,igcm_ccn_number)) 1947 2023 endif ! (scavenging) 1948 2024 -
trunk/LMDZ.MARS/libf/phymars/wstats.F90
r38 r719 18 18 character (len=50) :: namebis 19 19 character (len=50), save :: firstvar 20 integer :: ierr, varid,nbdim,nid20 integer :: ierr,ierr2,varid,nbdim,nid 21 21 integer :: meanid,sdid 22 22 integer, dimension(4) :: id,start,size … … 24 24 integer :: l,i,j,ig0 25 25 integer,save :: index 26 ! Added to use stats.def to select output variable 27 logical,save :: stats_def 28 logical :: getout 29 integer,save :: n_nom_def 30 integer :: n 31 integer,parameter :: n_nom_def_max=99 32 character(len=20),save :: nom_def(n_nom_def_max) 33 logical, save :: notfoundfirst=.TRUE. 26 34 27 35 integer, save :: step=0 28 36 37 29 38 30 39 if (firstcall) then 31 40 firstcall=.false. 32 firstvar=trim((nom))33 41 call inistats(ierr) 34 endif 42 ! at very first call, check if there is a "stats.def" to use and read it 43 open(99,file="stats.def",status='old',form='formatted',iostat=ierr2) 44 if(ierr2.eq.0) then 45 stats_def=.true. 46 write(*,*) "******************" 47 write(*,*) "Reading stats.def" 48 write(*,*) "******************" 49 do n=1,n_nom_def_max 50 read(99,fmt='(a)',end=88) nom_def(n) 51 write(*,*) 'Output in stats: ', trim(nom_def(n)) 52 end do 53 88 continue 54 if (n.ge.n_nom_def_max) then 55 write(*,*)"n_nom_def_max too small in wstats.F90:",n 56 stop 57 end if 58 n_nom_def=n-1 59 close(99) 60 else 61 firstvar=trim((nom)) 62 stats_def=.false. 63 notfoundfirst=.false. 64 endif 65 endif 66 67 ! find firstvar that is in stats.def 68 if (notfoundfirst) then 69 do n=1,n_nom_def_max 70 if(trim((nom_def(n))).eq.trim((nom))) then 71 firstvar=trim((nom)) 72 notfoundfirst=.false. 73 endif 74 end do 75 endif 76 77 ! get out of wstats if there is stats.def AND variable not listed 78 if (stats_def) then 79 getout=.true. 80 do n=1,n_nom_def 81 if(trim(nom_def(n)).eq.nom) getout=.false. 82 end do 83 if (getout) return 84 end if 35 85 36 86 if (firstvar==nom) then ! If we're back to the first variable
Note: See TracChangeset
for help on using the changeset viewer.
