21 " Testing program for the finite element implementation of the SSA.\n"
22 " Does a time-independent calculation. Does not run IceModel or a derived\n"
23 " class thereof. Uses verification test J. Also may be used in a PISM\n"
24 " software (regression) test.\n\n";
28#include "pism/stressbalance/ssa/SSATestCase.hh"
29#include "pism/util/Mask.hh"
30#include "pism/util/Context.hh"
31#include "pism/util/error_handling.hh"
32#include "pism/util/petscwrappers/PetscInitializer.hh"
33#include "pism/util/pism_options.hh"
34#include "pism/verification/tests/exactTestsIJ.h"
37namespace stressbalance {
39std::shared_ptr<Grid>
ssa_test_j_grid(std::shared_ptr<Context> ctx,
int Mx,
int My) {
57 double x,
double y,
double *u,
double *v);
66 double ocean_rho =
m_config->get_number(
"constants.sea_water.density"),
67 ice_rho =
m_config->get_number(
"constants.ice.density");
69 const double H0 = 500.;
74 m_ssa->strength_extension->set_notional_strength(nu0 *
H0);
75 m_ssa->strength_extension->set_min_thickness(800);
79 for (
auto p =
m_grid->points(); p; p.next()) {
80 const int i = p.i(), j = p.j();
92 if ((i == ((
int)
m_grid->Mx()) / 2) and
93 (j == ((
int)
m_grid->My()) / 2)) {
109 double *u,
double *v) {
118int main(
int argc,
char *argv[]) {
120 using namespace pism;
123 MPI_Comm com = MPI_COMM_WORLD;
126 com = PETSC_COMM_WORLD;
130 std::shared_ptr<Context> ctx = context_from_options(com,
"ssa_testj");
133 std::string usage =
"\n"
134 "usage of SSA_TESTJ:\n"
135 " run ssafe_test -Mx <number> -My <number> -ssa_method <fd|fem>\n"
138 bool stop = show_usage_check_req_opts(*ctx->log(),
"ssa_testj", {}, usage);
146 unsigned int Mx = config->get_number(
"grid.Mx");
147 unsigned int My = config->get_number(
"grid.My");
149 auto method = config->get_string(
"stress_balance.ssa.method");
150 auto output_file = config->get_string(
"output.file");
152 bool write_output = config->get_string(
"output.size") !=
"none";
154 config->set_flag(
"basal_resistance.pseudo_plastic.enabled",
false);
155 config->set_string(
"stress_balance.ssa.flow_law",
"isothermal_glen");
157 auto grid = ssa_test_j_grid(ctx, Mx, My);
158 SSATestCaseJ testcase(SSATestCase::solver(grid, method));
163 testcase.
write(output_file);
167 handle_fatal_errors(com);
std::shared_ptr< Config > Ptr
double enthalpy(double T, double omega, double P) const
Compute enthalpy from absolute temperature, liquid water fraction, and pressure.
Converts between specific enthalpy and temperature or liquid content.
array::Scalar2 ice_surface_elevation
array::CellType2 cell_type
array::Scalar2 ice_thickness
array::Scalar2 bed_elevation
Makes sure that we call begin_access() and end_access() for all accessed array::Arrays.
void set(double c)
Result: v[j] <- c for all j.
void update_ghosts()
Updates ghost points.
virtual void initializeSSACoefficients()
Set up the coefficient variables as appropriate for the test case.
SSATestCaseJ(std::shared_ptr< SSA > ssa)
virtual void exactSolution(int i, int j, double x, double y, double *u, double *v)
std::shared_ptr< const pism::Grid > m_grid
virtual void init()
Initialize the test case at the start of a run.
const units::System::Ptr m_sys
array::Array3D m_ice_enthalpy
static std::shared_ptr< Grid > grid(std::shared_ptr< Context > ctx, int Mx, int My, double Lx, double Ly, grid::Registration registration, grid::Periodicity periodicity)
array::Vector2 m_bc_values
virtual void report(const std::string &testname)
Report on the generated solution.
virtual void run()
Solve the SSA.
const Config::ConstPtr m_config
std::shared_ptr< SSA > m_ssa
virtual void write(const std::string &filename)
Save the computation and data to a file.
struct TestJParameters exactJ(const double x, const double y)
std::shared_ptr< Grid > ssa_test_j_grid(std::shared_ptr< Context > ctx, int Mx, int My)
Stress balance models and related diagnostics.
double convert(System::Ptr system, double input, const std::string &spec1, const std::string &spec2)
Convert a quantity from unit1 to unit2.
int main(int argc, char *argv[])