add debug plots; step one dt

Author: juliasloan25

experiments/AMIP/coupler_driver.jl

@@ -140,9 +140,9 @@ parsed_args = parse_commandline(argparse_settings())
 
 ## setup coupler and model configurations
 # modify parsed args for fast testing from REPL #hide
+include("user_io/debug_plots.jl")
 pkg_dir = pkgdir(ClimaCoupler)
 if isinteractive()
-    include("user_io/debug_plots.jl")
     parsed_args["config_file"] =
         isnothing(parsed_args["config_file"]) ? joinpath(pkg_dir, "config/model_configs/interactive_debug.yml") :
         parsed_args["config_file"]
@@ -165,10 +165,11 @@ run_name = config_dict["run_name"]
 energy_check = config_dict["energy_check"]
 const FT = config_dict["FLOAT_TYPE"] == "Float64" ? Float64 : Float32
 land_sim_name = "bucket"
-t_end = Float64(time_to_seconds(config_dict["t_end"]))
+# t_end = Float64(time_to_seconds(config_dict["t_end"]))
 t_start = 0.0
-tspan = (t_start, t_end)
 Δt_cpl = Float64(config_dict["dt_cpl"])
+t_end = Δt_cpl
+tspan = (t_start, t_end)
 saveat = Float64(time_to_seconds(config_dict["dt_save_to_sol"]))
 date0 = date = DateTime(config_dict["start_date"], dateformat"yyyymmdd")
 mono_surface = config_dict["mono_surface"]
@@ -194,6 +195,9 @@ isdir(COUPLER_ARTIFACTS_DIR) ? nothing : mkpath(COUPLER_ARTIFACTS_DIR)
 @info COUPLER_OUTPUT_DIR
 config_dict["print_config_dict"] ? @info(config_dict) : nothing
 
+debug_dir = joinpath(COUPLER_ARTIFACTS_DIR, "debug_output/")
+isdir(debug_dir) ? nothing : mkpath(debug_dir)
+
 # get the paths to the necessary data files: land-sea mask, sst map, sea ice concentration
 include(joinpath(pkgdir(ClimaCoupler), "artifacts", "artifact_funcs.jl"))
 sst_data = joinpath(sst_dataset_path(), "sst.nc")
@@ -534,6 +538,8 @@ cs = CoupledSimulation{FT}(
     dir_paths,
 );
 
+debug(cs, debug_dir * "0_initialized_cs")
+
 #=
 ## Restart component model states if specified
 =#
@@ -563,12 +569,16 @@ import_combined_surface_fields!(cs.fields, cs.model_sims, cs.boundary_space, tur
 import_atmos_fields!(cs.fields, cs.model_sims, cs.boundary_space, turbulent_fluxes)
 update_model_sims!(cs.model_sims, cs.fields, turbulent_fluxes)
 
+debug(cs, "debug_output/1_init_exchange_cs")
+
 # 2) each surface component model calculates its own vapor specific humidity (q_sfc)
 # TODO: the q_sfc calculation follows the design of the bucket q_sfc, but it would be neater to abstract this from step!
 step!(land_sim, Δt_cpl)
 step!(ocean_sim, Δt_cpl)
 step!(ice_sim, Δt_cpl)
 
+debug(cs, "debug_output/2_step_sfc_models")
+
 # 3) coupler re-imports updated surface fields and calculates turbulent fluxes, while updating atmos sfc_conditions
 if turbulent_fluxes isa CombinedStateFluxes
     # calculate fluxes using combined surface states on the atmos grid
@@ -585,16 +595,22 @@ elseif turbulent_fluxes isa PartitionedStateFluxes
     atmos_sim.integrator.p.precomputed.sfc_conditions .= new_p.precomputed.sfc_conditions
 end
 
+debug(cs, "debug_output/3_init_flux_calc")
+
 # 4) given the new sfc_conditions, atmos calls the radiative flux callback
 reinit_model_sims!(cs.model_sims) # NB: for atmos this sets a nonzero radiation flux
 
+debug(cs, "debug_output/4_reinit")
+
 # 5) coupler re-imports updated atmos fluxes (radiative fluxes for both `turbulent_fluxes` types
 # and also turbulent fluxes if `turbulent_fluxes isa CombinedStateFluxes`,
 # and sends them to the surface component models. If `turbulent_fluxes isa PartitionedStateFluxes`
 # atmos receives the turbulent fluxes from the coupler.
 import_atmos_fields!(cs.fields, cs.model_sims, cs.boundary_space, turbulent_fluxes)
 update_model_sims!(cs.model_sims, cs.fields, turbulent_fluxes)
 
+debug(cs, "debug_output/5_import_atmos_fluxes")
+
 #=
 ## Coupling Loop
 =#
@@ -720,6 +736,7 @@ end #hide
 ## run the coupled simulation
 solve_coupler!(cs);
 
+debug(cs, "debug_output/6_after_solve")
 #=
 ## Postprocessing
 Currently all postprocessing is performed using the root process only.

experiments/AMIP/user_io/debug_plots.jl

@@ -24,11 +24,21 @@ end
 Plot useful coupler fields (in `field_names`) and save plots to a directory.
 """
 function debug(cs_fields::NamedTuple, dir)
-    field_names = (:F_turb_energy, :F_turb_moisture, :P_liq, :T_S, :ρ_sfc, :q_sfc)
+    field_names = (:albedo, :F_radiative, :F_turb_energy, :F_turb_moisture, :P_liq, :T_S, :ρ_sfc, :q_sfc, :beta, :z0b_S, :z0m_S)
     all_plots = []
     for field_name in field_names
         field = getproperty(cs_fields, field_name)
         push!(all_plots, Plots.plot(field, title = string(field_name) * print_extrema(field)))
+        if (field_name == :T_S) && (@isdefined debug_csf0)
+            push!(
+                all_plots,
+                Plots.plot(
+                    field .- debug_csf0.T_S,
+                    title = string(field_name) * "_anom" * print_extrema(field),
+                    color = :viridis,
+                ),
+            )
+        end
     end
     fig = Plots.plot(all_plots..., size = (1500, 800))
     Plots.png(joinpath(dir, "debug_coupler"))