From e67d98810014098ca29c4ca383f92ca58ac1d4a2 Mon Sep 17 00:00:00 2001 From: Stuart Daines Date: Fri, 21 Jun 2024 10:56:08 +0100 Subject: [PATCH] add some 'romglb' tests that got missed --- Project.toml | 4 +- examples/Project.toml | 1 + examples/romglb/runtests.jl | 80 +++++++++++++++++++++++++++++++++++++ 3 files changed, 84 insertions(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index be4125e..05bb991 100644 --- a/Project.toml +++ b/Project.toml @@ -37,6 +37,7 @@ Plots = "1.38" Preferences = "1.3" SIMD = "3.4" SnoopPrecompile = "1.0" +SparseDiffTools = "2.0" SpecialFunctions = "1.0, 2.0" Sundials = "4.0" TestEnv = "1.0" @@ -52,9 +53,10 @@ Logging = "56ddb016-857b-54e1-b83d-db4d58db5568" PALEOcopse = "4a6ed817-0e58-48c6-8452-9e9afc8cb508" PALEOmodel = "bf7b4fbe-ccb1-42c5-83c2-e6e9378b660c" Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" +SparseDiffTools = "47a9eef4-7e08-11e9-0b38-333d64bd3804" Sundials = "c3572dad-4567-51f8-b174-8c6c989267f4" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" ZipFile = "a5390f91-8eb1-5f08-bee0-b1d1ffed6cea" [targets] -test = ["DataFrames", "DiffEqBase", "Documenter", "Downloads", "Logging", "PALEOcopse", "PALEOmodel", "Plots", "Sundials", "Test", "ZipFile"] +test = ["DataFrames", "DiffEqBase", "Documenter", "Downloads", "Logging", "PALEOcopse", "PALEOmodel", "Plots", "Sundials", "SparseDiffTools", "Test", "ZipFile"] diff --git a/examples/Project.toml b/examples/Project.toml index 380f933..2ffed7a 100644 --- a/examples/Project.toml +++ b/examples/Project.toml @@ -17,6 +17,7 @@ PALEOsediment = "e0a37952-6f01-4236-91ff-62fdc855f67b" PlotlyJS = "f0f68f2c-4968-5e81-91da-67840de0976a" Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" Revise = "295af30f-e4ad-537b-8983-00126c2a3abe" +SparseDiffTools = "47a9eef4-7e08-11e9-0b38-333d64bd3804" SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" Sundials = "c3572dad-4567-51f8-b174-8c6c989267f4" UnicodePlots = "b8865327-cd53-5732-bb35-84acbb429228" diff --git a/examples/romglb/runtests.jl b/examples/romglb/runtests.jl index 69bd82f..166c97a 100644 --- a/examples/romglb/runtests.jl +++ b/examples/romglb/runtests.jl @@ -3,6 +3,8 @@ using Logging using DiffEqBase using Sundials import DataFrames +import SparseArrays +import SparseDiffTools import PALEOboxes as PB @@ -14,6 +16,7 @@ import PALEOcopse @testset "romglb examples" begin skipped_testsets = [ + # "airsea_O2", # "O2_only", # "P_O2", # "P_O2_S_Carb_open", @@ -32,6 +35,81 @@ include("SedimentationRate_dev.jl") include("config_ocean_romglb_expts.jl") +!("airsea_O2" in skipped_testsets) && @testset "airsea_O2" begin + + model = PB.create_model_from_config( + joinpath(@__DIR__, "PALEO_examples_romglb_cfg.yaml"), "romglb_abiotic_O2"; + modelpars=Dict( + "matdir"=>matdir, + ) + ) + + # Test OceanBase domain configuration + @test PB.get_num_domains(model) == 6 + + global_domain = PB.get_domain(model, "global") + @test PB.get_length(global_domain) == 1 + + ocean_domain = PB.get_domain(model, "ocean") + @test PB.get_length(ocean_domain) == 79 + + + # test OceanBase variables + + initial_state, modeldata = PALEOmodel.initialize!(model) + + ocean_modelcreated_vars_dict = Dict([(var.name, var) for var in PB.get_variables(ocean_domain, hostdep=false)]) + + println("ocean model created variables after initialize!:") + for (name, var ) in ocean_modelcreated_vars_dict + println("\t", PB.fullname(var), " = ", PB.get_data(var, modeldata)) + end + + # bodge a test for ocean tracer with single non-zero cell + ocean_T_data = PB.get_data(PB.get_variable(ocean_domain,"T"), modeldata) + ocean_T_data .= 0.0 + ocean_T_data[1] = 1.0 + # update initial_state with our bodged values + initial_state = PALEOmodel.get_statevar(modeldata.solver_view_all) + + # Check model derivative + + PB.do_deriv(modeldata.dispatchlists_all) + + println("state, sms variables after check model derivative:") + for (state_var, sms_var) in PB.IteratorUtils.zipstrict(PB.get_vars(modeldata.solver_view_all.stateexplicit), + PB.get_vars(modeldata.solver_view_all.stateexplicit_deriv)) + println(PB.fullname(state_var), " ", PB.get_data(state_var, modeldata)) + println(PB.fullname(sms_var), " ", PB.get_data(sms_var, modeldata)) + end + + # check conservation + ocean_T_sms_data = PB.get_data(PB.get_variable(ocean_domain,"T_sms"), modeldata) + sum_T_sms = sum(ocean_T_sms_data) + println("check conservation: sum(Tracer_sms)=",sum_T_sms) + @test abs(sum_T_sms) < 1e-15 + + # check jacobian sparsity calculation + jac_prototype = PALEOmodel.JacobianAD.calcJacobianSparsitySparsityTracing!(model, modeldata, initial_state, 0.0) + @test SparseArrays.nnz(jac_prototype) == 819 + jac_proto_fill = PALEOmodel.SparseUtils.fill_sparse_jac(jac_prototype) + colors = SparseDiffTools.matrix_colors(copy(jac_proto_fill)) + @test maximum(colors) == 23 + + # integrate to approx steady state + paleorun = PALEOmodel.Run(model=model, output = PALEOmodel.OutputWriters.OutputMemory()) + PALEOmodel.ODE.integrateForwardDiff(paleorun, initial_state, modeldata, (0, 1e5), solvekwargs=(reltol=1e-5,)) # first run includes JIT time + + # check conservation + T_total = PB.get_data(paleorun.output, "ocean.T_total") + @test abs(T_total[1] - 1.0) < 1e-16 + @test abs(T_total[end] - T_total[1]) < 1e-4 + + total_O2 = PB.get_data(paleorun.output, "global.total_O2") + @test abs(total_O2[end] - total_O2[1]) < 1e-7*total_O2[1] + +end + !("O2_only" in skipped_testsets) && @testset "O2_only" begin model = PB.create_model_from_config( @@ -216,3 +294,5 @@ end end + +nothing # so no output printed when run from REPL \ No newline at end of file