Using data from GEOS-FP
This example demonstrates how to use the GEOS-FP data loader in the EarthSciML ecosystem. The GEOS-FP data loader is used to load data from the GEOS-FP dataset.
First, let's initialize some packages and set up the GEOS-FP equation system.
using EarthSciData, EarthSciMLBase
using DomainSets, ModelingToolkit, MethodOfLines, DifferentialEquations
using ModelingToolkit: t, D
using Dates, Plots, DataFrames
using DynamicQuantities
using DynamicQuantities: dimension
# Set up system
@parameters(
lev,
lon, [unit=u"rad"],
lat, [unit=u"rad"],
)
geosfp, geosfp_updater = GEOSFP("4x5")
geosfp\[ \begin{align} A3dyn_{+}DTRAIN\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( DTRAIN \right), t, lon, lat, lev \right) \\ A3dyn_{+}OMEGA\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( OMEGA \right), t, lon, lat, lev \right) \\ A3dyn_{+}RH\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( RH \right), t, lon, lat, lev \right) \\ A3dyn_{+}U\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( U \right), t, lon, lat, lev \right) \\ A3dyn_{+}V\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( V \right), t, lon, lat, lev \right) \\ A3cld_{+}CLOUD\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( CLOUD \right), t, lon, lat, lev \right) \\ A3cld_{+}OPTDEPTH\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( OPTDEPTH \right), t, lon, lat, lev \right) \\ A3cld_{+}QCCU\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( QCCU \right), t, lon, lat, lev \right) \\ A3cld_{+}QI\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( QI \right), t, lon, lat, lev \right) \\ A3cld_{+}QL\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( QL \right), t, lon, lat, lev \right) \\ A3cld_{+}TAUCLI\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( TAUCLI \right), t, lon, lat, lev \right) \\ A3cld_{+}TAUCLW\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( TAUCLW \right), t, lon, lat, lev \right) \\ I3_{+}PS\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PS \right), t, lon, lat \right) \\ I3_{+}PV\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PV \right), t, lon, lat, lev \right) \\ I3_{+}QV\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( QV \right), t, lon, lat, lev \right) \\ I3_{+}T\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( T \right), t, lon, lat, lev \right) \\ A3mstE_{+}CMFMC\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( CMFMC \right), t, lon, lat, lev \right) \\ A3mstE_{+}PFICU\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PFICU \right), t, lon, lat, lev \right) \\ A3mstE_{+}PFILSAN\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PFILSAN \right), t, lon, lat, lev \right) \\ A3mstE_{+}PFLCU\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PFLCU \right), t, lon, lat, lev \right) \\ A3mstE_{+}PFLLSAN\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PFLLSAN \right), t, lon, lat, lev \right) \\ A1_{+}ALBEDO\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( ALBEDO \right), t, lon, lat \right) \\ A1_{+}CLDTOT\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( CLDTOT \right), t, lon, lat \right) \\ A1_{+}EFLUX\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( EFLUX \right), t, lon, lat \right) \\ A1_{+}EVAP\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( EVAP \right), t, lon, lat \right) \\ A1_{+}FRSEAICE\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( FRSEAICE \right), t, lon, lat \right) \\ A1_{+}FRSNO\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( FRSNO \right), t, lon, lat \right) \\ A1_{+}GRN\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( GRN \right), t, lon, lat \right) \\ A1_{+}GWETROOT\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( GWETROOT \right), t, lon, lat \right) \\ A1_{+}GWETTOP\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( GWETTOP \right), t, lon, lat \right) \\ A1_{+}HFLUX\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( HFLUX \right), t, lon, lat \right) \\ A1_{+}LAI\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( LAI \right), t, lon, lat \right) \\ A1_{+}LWI\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( LWI \right), t, lon, lat \right) \\ A1_{+}LWGNT\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( LWGNT \right), t, lon, lat \right) \\ A1_{+}LWTUP\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( LWTUP \right), t, lon, lat \right) \\ A1_{+}PARDF\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PARDF \right), t, lon, lat \right) \\ A1_{+}PARDR\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PARDR \right), t, lon, lat \right) \\ A1_{+}PBLH\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PBLH \right), t, lon, lat \right) \\ A1_{+}PRECANV\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PRECANV \right), t, lon, lat \right) \\ A1_{+}PRECCON\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PRECCON \right), t, lon, lat \right) \\ A1_{+}PRECLSC\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PRECLSC \right), t, lon, lat \right) \\ A1_{+}PRECSNO\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PRECSNO \right), t, lon, lat \right) \\ A1_{+}PRECTOT\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PRECTOT \right), t, lon, lat \right) \\ A1_{+}QV2M\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( QV2M \right), t, lon, lat \right) \\ A1_{+}SEAICE00\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE00 \right), t, lon, lat \right) \\ A1_{+}SEAICE10\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE10 \right), t, lon, lat \right) \\ A1_{+}SEAICE20\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE20 \right), t, lon, lat \right) \\ A1_{+}SEAICE30\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE30 \right), t, lon, lat \right) \\ A1_{+}SEAICE40\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE40 \right), t, lon, lat \right) \\ A1_{+}SEAICE50\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE50 \right), t, lon, lat \right) \\ A1_{+}SEAICE60\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE60 \right), t, lon, lat \right) \\ A1_{+}SEAICE70\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE70 \right), t, lon, lat \right) \\ A1_{+}SEAICE80\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE80 \right), t, lon, lat \right) \\ A1_{+}SEAICE90\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE90 \right), t, lon, lat \right) \\ A1_{+}SLP\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SLP \right), t, lon, lat \right) \\ A1_{+}SNODP\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SNODP \right), t, lon, lat \right) \\ A1_{+}SNOMAS\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SNOMAS \right), t, lon, lat \right) \\ A1_{+}SWGDN\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SWGDN \right), t, lon, lat \right) \\ A1_{+}TO3\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( TO3 \right), t, lon, lat \right) \\ A1_{+}TROPPT\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( TROPPT \right), t, lon, lat \right) \\ A1_{+}TS\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( TS \right), t, lon, lat \right) \\ A1_{+}T2M\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( T2M \right), t, lon, lat \right) \\ A1_{+}U10M\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( U10M \right), t, lon, lat \right) \\ A1_{+}USTAR\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( USTAR \right), t, lon, lat \right) \\ A1_{+}V10M\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( V10M \right), t, lon, lat \right) \\ A1_{+}Z0M\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( Z0M \right), t, lon, lat \right) \\ A1_{+}T10M\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( T10M \right), t, lon, lat \right) \\ A1_{+}Q850\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( Q850 \right), t, lon, lat \right) \\ A3mstC_{+}DQRCU\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( DQRCU \right), t, lon, lat, lev \right) \\ A3mstC_{+}DQRLSAN\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( DQRLSAN \right), t, lon, lat, lev \right) \\ A3mstC_{+}REEVAPCN\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( REEVAPCN \right), t, lon, lat, lev \right) \\ A3mstC_{+}REEVAPLS\left( t \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( REEVAPLS \right), t, lon, lat, lev \right) \\ P\left( t \right) &= P_{unit} \mathrm{DataInterpolations.LinearInterpolation{Vector{Float64}, UnitRange{Int64}, Vector{Float64}, Vector{Float64}, Float64}([0.0, 4.804826, 659.3752000000001, 1313.48, 1961.311, 2609.201, 3257.081, 3898.201, 4533.901, 5169.611, 5805.321, 6436.264, 7062.197999999999, 7883.4220000000005, 8909.992, 9936.521, 10918.17, 11895.86, 12869.59, 14291.0, 15626.0, 16960.9, 18161.9, 19309.7, 20325.899999999998, 21215.0, 21877.600000000002, 22389.8, 22436.3, 21686.5, 20119.2, 17693.0, 15039.3, 12783.7, 10866.3, 9236.572, 7851.231, 6660.340999999999, 5638.791, 4764.391, 4017.541, 3381.0009999999997, 2836.781, 2373.0409999999997, 1979.1599999999999, 1645.71, 1364.34, 1127.69, 929.2942, 761.9842, 621.6801, 504.68010000000004, 407.6571, 327.6431, 262.0211, 208.497, 165.079, 130.05100000000002, 101.94399999999999, 79.51341, 61.677910000000004, 47.58061, 36.50411, 27.85261, 21.134900000000002, 15.9495, 11.9703, 8.934502, 6.600001, 4.758501, 3.27, 2.0, 1.0], 1:73, Float64[], DataInterpolations.LinearParameterCache{Vector{Float64}}(Float64[]), false, FindFirstFunctions.Guesser{UnitRange{Int64}}(1:73, Base.RefValue{Int64}(1), true), false, true)}\left( lev \right) + \mathrm{DataInterpolations.LinearInterpolation{Vector{Float64}, UnitRange{Int64}, Vector{Float64}, Vector{Float64}, Float64}([1.0, 0.984952, 0.963406, 0.941865, 0.920387, 0.898908, 0.877429, 0.856018, 0.8346609, 0.8133039, 0.7919469, 0.7706375, 0.7493782, 0.721166, 0.6858999, 0.6506349, 0.6158184, 0.5810415, 0.5463042, 0.4945902, 0.4437402, 0.3928911, 0.3433811, 0.2944031, 0.2467411, 0.2003501, 0.1562241, 0.1136021, 0.06372006, 0.02801004, 0.006960025, 8.175413e-9, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], 1:73, Float64[], DataInterpolations.LinearParameterCache{Vector{Float64}}(Float64[]), false, FindFirstFunctions.Guesser{UnitRange{Int64}}(1:73, Base.RefValue{Int64}(1), true), false, true)}\left( lev \right) I3_{+}PS\left( t \right) \end{align} \]
Note that the GEOSFP function returns to things, an equation system and an object that can used to update the time in the underlying data loaders. We can see above the different variables that are available in the GEOS-FP dataset. But also, here they are in table form:
vars = unknowns(geosfp)
DataFrame(
:Name => [string(Symbolics.tosymbol(v, escape=false)) for v ∈ vars],
:Units => [dimension(ModelingToolkit.get_unit(v)) for v ∈ vars],
:Description => [ModelingToolkit.getdescription(v) for v ∈ vars],
)| Row | Name | Units | Description |
|---|---|---|---|
| String | Dimensio… | String | |
| 1 | A3dyn₊DTRAIN | m⁻² kg s⁻² | Detrainment cloud mass flux |
| 2 | A3dyn₊OMEGA | m⁻¹ kg s⁻³ | Vertical pressure velocity |
| 3 | A3dyn₊RH | Relative humidity | |
| 4 | A3dyn₊U | m s⁻¹ | Eastward component of wind |
| 5 | A3dyn₊V | m s⁻¹ | Northward component of wind |
| 6 | A3cld₊CLOUD | Total cloud fraction in grid box | |
| 7 | A3cld₊OPTDEPTH | Total in-cloud optical thickness (visible band) | |
| 8 | A3cld₊QCCU | Cloud condensate mixing ratio, convective updraft | |
| 9 | A3cld₊QI | Cloud ice water mixing ratio | |
| 10 | A3cld₊QL | Cloud liquid water mixing ratio | |
| 11 | A3cld₊TAUCLI | In-cloud ice optical thickness (visible band) | |
| 12 | A3cld₊TAUCLW | In-cloud water optical thickness (visible band) | |
| 13 | I3₊PS | m⁻¹ kg s⁻² | Pressure at the surface |
| 14 | I3₊PV | m⁻² kg⁻¹ s⁻¹ K | Ertel potential vorticity |
| 15 | I3₊QV | Specific humidity | |
| 16 | I3₊T | K | Temperature |
| 17 | A3mstE₊CMFMC | m⁻² kg s⁻² | Upward moist convective mass flux |
| 18 | A3mstE₊PFICU | m⁻² kg s⁻¹ | Downward flux of ice precipitation (convective) |
| 19 | A3mstE₊PFILSAN | m⁻² kg s⁻¹ | Downward flux of ice precipitation (large scale + anvil) |
| 20 | A3mstE₊PFLCU | m⁻² kg s⁻¹ | Downward flux of liquid precipitation (convective) |
| 21 | A3mstE₊PFLLSAN | m⁻² kg s⁻¹ | Downward flux of liquid precipitation (large scale + anvil) |
| 22 | A1₊ALBEDO | Surface albedo | |
| 23 | A1₊CLDTOT | Total cloud fraction | |
| 24 | A1₊EFLUX | kg s⁻³ | Latent heat flux (positive upward) |
| 25 | A1₊EVAP | m⁻² kg s⁻² | Surface evaporation |
| 26 | A1₊FRSEAICE | Fraction of sea ice on surface | |
| 27 | A1₊FRSNO | Fractional snow-covered area | |
| 28 | A1₊GRN | Vegetation greenness fraction | |
| 29 | A1₊GWETROOT | Root zone soil wetness | |
| 30 | A1₊GWETTOP | Top soil wetness | |
| 31 | A1₊HFLUX | kg s⁻³ | Sensible heat flux (positive upward) |
| 32 | A1₊LAI | Leaf area index | |
| 33 | A1₊LWI | Land/water/ice flags | |
| 34 | A1₊LWGNT | kg s⁻³ | Net longwave flux at the ground |
| 35 | A1₊LWTUP | kg s⁻³ | Upward longwave flux at top of atmosphere (TOA) |
| 36 | A1₊PARDF | kg s⁻³ | Surface downward PAR diffuse flux |
| 37 | A1₊PARDR | kg s⁻³ | Surface downward PAR beam flux |
| 38 | A1₊PBLH | m | Planetary boundary layer height above surface |
| 39 | A1₊PRECANV | m⁻² kg s⁻¹ | Surface precipitation flux from anvils |
| 40 | A1₊PRECCON | m⁻² kg s⁻¹ | Surface precipitation flux from convection |
| 41 | A1₊PRECLSC | m⁻² kg s⁻¹ | Surface precipitation flux from large-scale |
| 42 | A1₊PRECSNO | m⁻² kg s⁻¹ | Surface precipitation flux from snow |
| 43 | A1₊PRECTOT | m⁻² kg s⁻¹ | Total surface precipitation flux |
| 44 | A1₊QV2M | Specific humidity at 2m above the displacement height | |
| 45 | A1₊SEAICE00 | Fraction of grid box that has 0-10% sea ice coverage | |
| 46 | A1₊SEAICE10 | Fraction of grid box that has 10-20% sea ice coverage | |
| 47 | A1₊SEAICE20 | Fraction of grid box that has 20-30% sea ice coverage | |
| 48 | A1₊SEAICE30 | Fraction of grid box that has 30-40% sea ice coverage | |
| 49 | A1₊SEAICE40 | Fraction of grid box that has 40-50% sea ice coverage | |
| 50 | A1₊SEAICE50 | Fraction of grid box that has 50-60% sea ice coverage | |
| 51 | A1₊SEAICE60 | Fraction of grid box that has 60-70% sea ice coverage | |
| 52 | A1₊SEAICE70 | Fraction of grid box that has 70-80% sea ice coverage | |
| 53 | A1₊SEAICE80 | Fraction of grid box that has 80-90% sea ice coverage | |
| 54 | A1₊SEAICE90 | Fraction of grid box that has 90-100% sea ice coverage | |
| 55 | A1₊SLP | m⁻¹ kg s⁻² | Sea level pressure |
| 56 | A1₊SNODP | m | Snow depth |
| 57 | A1₊SNOMAS | m⁻² kg | Snow mass |
| 58 | A1₊SWGDN | kg s⁻³ | Surface incident shortwave flux |
| 59 | A1₊TO3 | m⁻² mol | Total column ozone |
| 60 | A1₊TROPPT | m⁻¹ kg s⁻² | Temperature-based tropopause pressure |
| 61 | A1₊TS | K | Surface skin temperature |
| 62 | A1₊T2M | K | Temperature 2m above displacement height |
| 63 | A1₊U10M | m s⁻¹ | Eastward wind 10m above displacement height |
| 64 | A1₊USTAR | m s⁻¹ | Friction velocity |
| 65 | A1₊V10M | m s⁻¹ | Northward wind 10m above displacement height |
| 66 | A1₊Z0M | m | Roughness length, momentum |
| 67 | A1₊T10M | K | Temperature at 10 m above the displacement height |
| 68 | A1₊Q850 | Specific humidity at 850 hPa | |
| 69 | A3mstC₊DQRCU | s⁻¹ | Precipitation production rate -- convective |
| 70 | A3mstC₊DQRLSAN | s⁻¹ | Precipitation production rate -- large scale + anvil |
| 71 | A3mstC₊REEVAPCN | s⁻¹ | Evaporation of precipitating convective condensate |
| 72 | A3mstC₊REEVAPLS | Evaporation of precipitating large-scale & anvil condensate | |
| 73 | P | m⁻¹ kg s⁻² | Pressure |
The GEOS-FP equation system isn't an ordinary differential equation (ODE) system, so we can't run it by itself. To fix this, we create another equation system that is an ODE. (We don't actually end up using this system for anything, it's just necessary to get the system to compile.)
function Example()
@variables c(t) = 5.0 [unit=u"s"]
ODESystem([D(c) ~ sin(lat * 6) + sin(lon * 6)], t, name=:Docs₊Example)
end
examplesys = Example()\[ \begin{align} \frac{\mathrm{d} c\left( t \right)}{\mathrm{d}t} &= \sin\left( 6 lon \right) + \sin\left( 6 lat \right) \end{align} \]
Now, let's couple these two systems together, and also add in advection and some information about the domain:
domain = DomainInfo(
partialderivatives_δxyδlonlat,
constIC(0.0, t ∈ Interval(Dates.datetime2unix(DateTime(2022, 1, 1)), Dates.datetime2unix(DateTime(2022, 1, 3)))),
zerogradBC(lat ∈ Interval(deg2rad(-80.0f0), deg2rad(80.0f0))),
periodicBC(lon ∈ Interval(deg2rad(-180.0f0), deg2rad(180.0f0))),
zerogradBC(lev ∈ Interval(1.0f0, 11.0f0)),
)
composed_sys = couple(examplesys, domain, geosfp, geosfp_updater)
pde_sys = convert(PDESystem, composed_sys)\[ \begin{align} \frac{\mathrm{d}}{\mathrm{d}t} Docs_{+}Example_{+}c\left( t, lat, lon, lev \right) &= \sin\left( 6 lon \right) + \sin\left( 6 lat \right) \\ GEOSFP_{+}A3dyn_{+}DTRAIN\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( DTRAIN \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3dyn_{+}OMEGA\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( OMEGA \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3dyn_{+}RH\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( RH \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3dyn_{+}U\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( U \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3dyn_{+}V\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( V \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3cld_{+}CLOUD\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( CLOUD \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3cld_{+}OPTDEPTH\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( OPTDEPTH \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3cld_{+}QCCU\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( QCCU \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3cld_{+}QI\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( QI \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3cld_{+}QL\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( QL \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3cld_{+}TAUCLI\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( TAUCLI \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3cld_{+}TAUCLW\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( TAUCLW \right), t, lon, lat, lev \right) \\ GEOSFP_{+}I3_{+}PS\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PS \right), t, lon, lat \right) \\ GEOSFP_{+}I3_{+}PV\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PV \right), t, lon, lat, lev \right) \\ GEOSFP_{+}I3_{+}QV\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( QV \right), t, lon, lat, lev \right) \\ GEOSFP_{+}I3_{+}T\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( T \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3mstE_{+}CMFMC\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( CMFMC \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3mstE_{+}PFICU\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PFICU \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3mstE_{+}PFILSAN\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PFILSAN \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3mstE_{+}PFLCU\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PFLCU \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3mstE_{+}PFLLSAN\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PFLLSAN \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A1_{+}ALBEDO\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( ALBEDO \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}CLDTOT\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( CLDTOT \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}EFLUX\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( EFLUX \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}EVAP\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( EVAP \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}FRSEAICE\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( FRSEAICE \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}FRSNO\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( FRSNO \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}GRN\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( GRN \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}GWETROOT\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( GWETROOT \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}GWETTOP\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( GWETTOP \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}HFLUX\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( HFLUX \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}LAI\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( LAI \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}LWI\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( LWI \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}LWGNT\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( LWGNT \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}LWTUP\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( LWTUP \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}PARDF\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PARDF \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}PARDR\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PARDR \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}PBLH\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PBLH \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}PRECANV\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PRECANV \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}PRECCON\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PRECCON \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}PRECLSC\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PRECLSC \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}PRECSNO\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PRECSNO \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}PRECTOT\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( PRECTOT \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}QV2M\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( QV2M \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SEAICE00\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE00 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SEAICE10\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE10 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SEAICE20\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE20 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SEAICE30\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE30 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SEAICE40\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE40 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SEAICE50\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE50 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SEAICE60\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE60 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SEAICE70\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE70 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SEAICE80\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE80 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SEAICE90\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SEAICE90 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SLP\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SLP \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SNODP\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SNODP \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SNOMAS\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SNOMAS \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}SWGDN\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( SWGDN \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}TO3\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( TO3 \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}TROPPT\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( TROPPT \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}TS\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( TS \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}T2M\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( T2M \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}U10M\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( U10M \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}USTAR\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( USTAR \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}V10M\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( V10M \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}Z0M\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( Z0M \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}T10M\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( T10M \right), t, lon, lat \right) \\ GEOSFP_{+}A1_{+}Q850\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( Q850 \right), t, lon, lat \right) \\ GEOSFP_{+}A3mstC_{+}DQRCU\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( DQRCU \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3mstC_{+}DQRLSAN\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( DQRLSAN \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3mstC_{+}REEVAPCN\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( REEVAPCN \right), t, lon, lat, lev \right) \\ GEOSFP_{+}A3mstC_{+}REEVAPLS\left( t, lat, lon, lev \right) &= \mathrm{interp}_{unsafe}\left( \mathrm{GEOSFPFileSet}\left( REEVAPLS \right), t, lon, lat, lev \right) \\ GEOSFP_{+}P\left( t, lat, lon, lev \right) &= GEOSFP_{+}P_{unit} \mathrm{DataInterpolations.LinearInterpolation{Vector{Float64}, UnitRange{Int64}, Vector{Float64}, Vector{Float64}, Float64}([0.0, 4.804826, 659.3752000000001, 1313.48, 1961.311, 2609.201, 3257.081, 3898.201, 4533.901, 5169.611, 5805.321, 6436.264, 7062.197999999999, 7883.4220000000005, 8909.992, 9936.521, 10918.17, 11895.86, 12869.59, 14291.0, 15626.0, 16960.9, 18161.9, 19309.7, 20325.899999999998, 21215.0, 21877.600000000002, 22389.8, 22436.3, 21686.5, 20119.2, 17693.0, 15039.3, 12783.7, 10866.3, 9236.572, 7851.231, 6660.340999999999, 5638.791, 4764.391, 4017.541, 3381.0009999999997, 2836.781, 2373.0409999999997, 1979.1599999999999, 1645.71, 1364.34, 1127.69, 929.2942, 761.9842, 621.6801, 504.68010000000004, 407.6571, 327.6431, 262.0211, 208.497, 165.079, 130.05100000000002, 101.94399999999999, 79.51341, 61.677910000000004, 47.58061, 36.50411, 27.85261, 21.134900000000002, 15.9495, 11.9703, 8.934502, 6.600001, 4.758501, 3.27, 2.0, 1.0], 1:73, Float64[], DataInterpolations.LinearParameterCache{Vector{Float64}}(Float64[]), false, FindFirstFunctions.Guesser{UnitRange{Int64}}(1:73, Base.RefValue{Int64}(1), true), false, true)}\left( lev \right) + \mathrm{DataInterpolations.LinearInterpolation{Vector{Float64}, UnitRange{Int64}, Vector{Float64}, Vector{Float64}, Float64}([1.0, 0.984952, 0.963406, 0.941865, 0.920387, 0.898908, 0.877429, 0.856018, 0.8346609, 0.8133039, 0.7919469, 0.7706375, 0.7493782, 0.721166, 0.6858999, 0.6506349, 0.6158184, 0.5810415, 0.5463042, 0.4945902, 0.4437402, 0.3928911, 0.3433811, 0.2944031, 0.2467411, 0.2003501, 0.1562241, 0.1136021, 0.06372006, 0.02801004, 0.006960025, 8.175413e-9, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], 1:73, Float64[], DataInterpolations.LinearParameterCache{Vector{Float64}}(Float64[]), false, FindFirstFunctions.Guesser{UnitRange{Int64}}(1:73, Base.RefValue{Int64}(1), true), false, true)}\left( lev \right) GEOSFP_{+}I3_{+}PS\left( t, lat, lon, lev \right) \end{align} \]
You can see above that we add both geosfp and geosfp_updater to our coupled system. If we didn't include the updater, the resulting model would not give the correct results.
Now, finally, we can run the simulation and plot the GEOS-FP wind fields in the result:
(The code below is commented out because it is very slow right now. A faster solution is coming soon!)
# discretization = MOLFiniteDifference([lat => 10, lon => 10, lev => 10], t, approx_order=2)
# @time pdeprob = discretize(pde_sys, discretization)
# pdesol = solve(pdeprob, Tsit5(), saveat=3600.0)
# discrete_lon = pdesol[lon]
# discrete_lat = pdesol[lat]
# discrete_lev = pdesol[lev]
# discrete_t = pdesol[t]
# @variables meanwind₊u(..) meanwind₊v(..) examplesys₊c(..)
# sol_u = pdesol[meanwind₊u(t, lat, lon, lev)]
# sol_v = pdesol[meanwind₊v(t, lat, lon, lev)]
# sol_c = pdesol[examplesys₊c(t, lat, lon, lev)]
# anim = @animate for k in 1:length(discrete_t)
# p1 = heatmap(discrete_lon, discrete_lat, sol_c[k, 1:end, 1:end, 2], clim=(minimum(sol_c[:, :, :, 2]), maximum(sol_c[:, :, :, 2])),
# xlabel="Longitude", ylabel="Latitude", title="examplesys.c: $(Dates.unix2datetime(discrete_t[k]))")
# p2 = heatmap(discrete_lon, discrete_lat, sol_u[k, 1:end, 1:end, 2], clim=(minimum(sol_u[:, :, :, 2]), maximum(sol_u[:, :, :, 2])),
# title="U")
# p3 = heatmap(discrete_lon, discrete_lat, sol_v[k, 1:end, 1:end, 2], clim=(minimum(sol_v[:, :, :, 2]), maximum(sol_v[:, :, :, 2])),
# title="V")
# plot(p1, p2, p3, size=(800, 500))
# end
# gif(anim, "animation.gif", fps = 8)