Initial conditions from file (rough draft, work-in-progress)

CliMA · Nov 26, 2024 · a021e40 · a021e40
1 parent 002b4f3
commit a021e40
Show file tree

Hide file tree

Showing 8 changed files with 177 additions and 2 deletions.
diff --git a/Artifacts.toml b/Artifacts.toml
@@ -43,3 +43,11 @@ git-tree-sha1 = "10742e0a2e343d13bb04df379e300a83402d4106"
     [[era5_cloud.download]]
     sha256 = "bb51e2f2d315b487e05a8d38944d4ad937ee4a40c43b68541220c5d54425e24a"
     url = "https://caltech.box.com/shared/static/b6ur4ap4vo04j09vdulem96z9fxqlgyn.gz"
+
+[DYAMOND_SUMMER_ICS_p98deg]
+git-tree-sha1 = "5ef4bb950eceb86430464c8637ede7f075a1449f"
+lazy = true
+
+    [[DYAMOND_SUMMER_ICS_p98deg.download]]
+    sha256 = "da552c32c4e66a053992f3b2c96ebaeb78c1aded1818b4e7c2deea77d521dda1"
+    url = "https://caltech.box.com/shared/static/ad8d774iwaezctwud1g1lqbp9irbjj38.gz"
diff --git a/artifacts/artifact_funcs.jl b/artifacts/artifact_funcs.jl
@@ -61,3 +61,16 @@ function gfdl_ogw_data_path()
     )
     return AW.get_data_folder(gfdl_data)
 end
+
+# DYAMOND summer initial conditions
+function DYAMOND_SUMMER_ICS_p98deg_dataset_path()
+    DYAMOND_SUMMER_ICS_p98deg_dataset = AW.ArtifactWrapper(
+        @__DIR__,
+        "DYAMOND_SUMMER_ICS_p98deg",
+        AW.ArtifactFile[AW.ArtifactFile(
+            url = "https://caltech.box.com/shared/static/1nqueb39ro33b2ps7zm2df9hpbkps6wy.gz",
+            filename = "DYAMOND_SUMMER_ICS_p98deg.nc",
+        ),],
+    )
+    return AW.get_data_folder(DYAMOND_SUMMER_ICS_p98deg_dataset)
+end
diff --git a/artifacts/download_artifacts.jl b/artifacts/download_artifacts.jl
@@ -7,6 +7,7 @@ function trigger_download(lazy_download = true)
     @info "topo dataset path:`$(topo_res_path())`"
     @info "MiMA convective gravity wave path:`$(mima_gwf_path())`"
     @info "GFDL OGWD test data:`$(gfdl_ogw_data_path())`"
+    @info "DYAMOND summer initial conditions path:`$(DYAMOND_SUMMER_ICS_p98deg_dataset_path())`"
     return nothing
 end
 trigger_download()
diff --git a/config/default_configs/default_config.yml b/config/default_configs/default_config.yml
@@ -179,7 +179,7 @@ surface_temperature:
   help: "Prescribed surface temperature functional form ['ZonallySymmetric' (default), 'ZonallyAsymmetric', 'RCEMIPII']"
   value: "ZonallySymmetric"
 initial_condition:
-  help: "Initial condition [`DryBaroclinicWave`, `MoistBaroclinicWave`, `DecayingProfile`, `IsothermalProfile`, `Bomex`, `DryDensityCurrentProfile`, `AgnesiHProfile`, `ScharProfile`, `RisingThermalBubbleProfile`, `ISDAC`]"
+  help: "Initial condition [`DryBaroclinicWave`, `MoistBaroclinicWave`, `DecayingProfile`, `IsothermalProfile`, `Bomex`, `DryDensityCurrentProfile`, `AgnesiHProfile`, `ScharProfile`, `RisingThermalBubbleProfile`, `ISDAC`, `DYAMONDSummer`]"
   value: "DecayingProfile"
 perturb_initstate:
   help: "Add a perturbation to the initial condition [`false`, `true` (default)]"

diff --git a/config/gpu_configs/gpu_aquaplanet_dyamond_summer.yml b/config/gpu_configs/gpu_aquaplanet_dyamond_summer.yml
@@ -0,0 +1,26 @@
+dt_save_state_to_disk: "Inf"
+dt_save_to_sol: "Inf"
+output_default_diagnostics: false
+h_elem: 30
+z_max: 60000.0
+z_elem: 63
+dz_bottom: 30.0
+rayleigh_sponge: true
+viscous_sponge: true
+moist: "equil"
+precip_model: "1M"
+rad: "allskywithclear"
+insolation: "timevarying"
+dt_rad: "1hours"
+dt_cloud_fraction: "1hours"
+vert_diff: "FriersonDiffusion"
+implicit_diffusion: true
+approximate_linear_solve_iters: 2
+surface_setup: "DefaultMoninObukhov"
+dt: "90secs"
+t_end: "1days"
+toml: [toml/longrun_aquaplanet.toml]
+prescribe_ozone: true
+aerosol_radiation: true
+prescribed_aerosols: ["CB1", "CB2", "DST01", "OC1", "OC2", "SO4", "SSLT01"]
+initial_condition: "DYAMONDSummer"
diff --git a/src/initial_conditions/initial_conditions.jl b/src/initial_conditions/initial_conditions.jl
@@ -160,6 +160,45 @@ function (initial_condition::DecayingProfile)(params)
     return local_state
 end
 
+"""
+    DYAMONDSummer(; perturb = true)
+
+An `InitialCondition` with a decaying temperature profile, and with an optional
+perturbation to the temperature.
+"""
+Base.@kwdef struct DYAMONDSummer <: InitialCondition
+    perturb::Bool = true
+end
+
+function (initial_condition::DYAMONDSummer)(params)
+    (; perturb) = initial_condition
+    function local_state(local_geometry)
+        FT = eltype(params)
+        grav = CAP.grav(params)
+        thermo_params = CAP.thermodynamics_params(params)
+        temp_profile = DecayingTemperatureProfile{FT}(
+            thermo_params,
+            FT(290),
+            FT(220),
+            FT(8e3),
+        )
+
+        (; z) = local_geometry.coordinates
+        coeff = perturb_coeff(local_geometry.coordinates)
+        T, p = temp_profile(thermo_params, z)
+        if perturb
+            T += coeff * FT(0.1) * (z < 5000)
+        end
+
+        return LocalState(;
+            params,
+            geometry = local_geometry,
+            thermo_state = TD.PhaseDry_pT(thermo_params, p, T),
+        )
+    end
+    return local_state
+end
+
 """
     AgnesiHProfile(; perturb = false)
 

diff --git a/src/solver/type_getters.jl b/src/solver/type_getters.jl
@@ -1,13 +1,17 @@
 using Dates: DateTime, @dateformat_str
 import Interpolations
 import NCDatasets
+import ClimaUtilities
 import ClimaUtilities.OutputPathGenerator
 import ClimaCore: InputOutput, Meshes, Spaces, Quadratures
 import ClimaAtmos.RRTMGPInterface as RRTMGPI
 import ClimaAtmos as CA
 import LinearAlgebra
 import ClimaCore.Fields
 import ClimaTimeSteppers as CTS
+import ClimaAtmos.AtmosArtifacts as AA
+#using ClimaParams
+import Thermodynamics as TD
 
 import ClimaDiagnostics
 
@@ -267,6 +271,10 @@ function get_initial_condition(parsed_args)
         return getproperty(ICs, Symbol(parsed_args["initial_condition"]))(
             parsed_args["perturb_initstate"],
         )
+    elseif parsed_args["initial_condition"] in ["DYAMONDSummer"]
+        return getproperty(ICs, Symbol(parsed_args["initial_condition"]))(
+            parsed_args["perturb_initstate"],
+        )
     elseif parsed_args["initial_condition"] in [
         "Nieuwstadt",
         "GABLS",
@@ -616,7 +624,6 @@ end
 function get_simulation(config::AtmosConfig)
     params = create_parameter_set(config)
     atmos = get_atmos(config, params)
-
     sim_info = get_sim_info(config)
     job_id = sim_info.job_id
     output_dir = sim_info.output_dir
@@ -661,6 +668,79 @@ function get_simulation(config::AtmosConfig)
 
     tracers = get_tracers(config.parsed_args)
 
+    if initial_condition isa ClimaAtmos.InitialConditions.DYAMONDSummer
+        @show "DYAMONDSummer, initialize from file"
+        file_path = AA.dyamond_summer_artifact_path(;
+            context = get_comms_context(config.parsed_args),
+        )
+        Pressure = ClimaUtilities.SpaceVaryingInputs.SpaceVaryingInput(
+            file_path,
+            "p",
+            spaces.center_space,
+        )
+        Temperature = ClimaUtilities.SpaceVaryingInputs.SpaceVaryingInput(
+            file_path,
+            "t",
+            spaces.center_space,
+        )
+        q =
+            TD.PhasePartition.(
+                ClimaUtilities.SpaceVaryingInputs.SpaceVaryingInput(
+                    file_path,
+                    "q",
+                    spaces.center_space,
+                )
+            )
+        Y.c.ρ .=
+            TD.air_density.(
+                params.thermodynamics_params,
+                Temperature,
+                Pressure,
+                q,
+            )
+        vel =
+            ClimaCore.Geometry.UVWVector.(
+                ClimaUtilities.SpaceVaryingInputs.SpaceVaryingInput(
+                    file_path,
+                    "u",
+                    spaces.center_space,
+                ),
+                ClimaUtilities.SpaceVaryingInputs.SpaceVaryingInput(
+                    file_path,
+                    "v",
+                    spaces.center_space,
+                ),
+                ClimaUtilities.SpaceVaryingInputs.SpaceVaryingInput(
+                    file_path,
+                    "w",
+                    spaces.center_space,
+                ),
+            )
+        Y.c.uₕ .=
+            ClimaCore.Geometry.Covariant12Vector.(
+                ClimaCore.Geometry.UVVector.(vel)
+            )
+        Y.f.u₃ .=
+            ᶠinterp.(
+                ClimaCore.Geometry.Covariant3Vector.(
+                    ClimaCore.Geometry.WVector.(vel)
+                )
+            )
+        e_kin = similar(Temperature)
+        compute_kinetic!(e_kin, Y.c.uₕ, Y.f.u₃)
+        e_pot =
+            ClimaCore.Fields.coordinate_field(Y.c).z .*
+            params.thermodynamics_params.grav
+        Y.c.ρe_tot .=
+            TD.total_energy.(
+                params.thermodynamics_params,
+                Temperature,
+                e_kin,
+                e_pot,
+            ) .* Y.c.ρ
+        @show "****************************"
+    end
+
     s = @timed_str begin
         p = build_cache(
             Y,

diff --git a/src/utils/AtmosArtifacts.jl b/src/utils/AtmosArtifacts.jl
@@ -82,4 +82,12 @@ function earth_orography_file_path(; context = nothing)
     )
 end
 
+function dyamond_summer_artifact_path(; context = nothing)
+    filename = "DYAMOND_SUMMER_ICS_p98deg.nc"
+    return joinpath(
+        @clima_artifact("DYAMOND_SUMMER_ICS_p98deg", context),
+        filename,
+    )
+end
+
 end