artificial-tsunami.py

Go to the documentation of this file.

# This file is part of the ExaHyPE2 project. For conditions of distribution and
# use, please see the copyright notice at www.peano-framework.org
import os
import sys
 
import peano4
import exahype2
 
sys.path.insert(0, os.path.abspath("../../../../applications/exahype2/shallow-water"))
from PDE import *
 
end_time = 10.0
size = [10000, 10000]  # [m]
offset = [-5000, -5000]  # [m]
min_depth = 4
amr_levels = 0
dg_order = 5
number_of_snapshots = 0
build_mode = "Release"
output = "solutions"
 
initial_conditions = """
  for (int i = 0; i < NumberOfUnknowns + NumberOfAuxiliaryVariables; i++) {
    Q[i] = 0.0;
  }
 
  constexpr double INITIAL_WATER_HEIGHT = 100.0;
  constexpr double INITIAL_BATHYMETRY_BEFORE_EARTHQUAKE = -100.0;
 
  Q[Shortcuts::h] = INITIAL_WATER_HEIGHT;
 
  auto displacementX = [](double x) { return std::sin((x / 500.0 + 1.0) * tarch::la::PI); };
  auto displacementY = [](double y) { return -std::pow(y / 500.0, 2) + 1.0; };
  auto displacement = [&](double x, double y) { return 5.0 * displacementX(x) * displacementY(y); };
 
  const double bathymetryAfterEarthquake = INITIAL_BATHYMETRY_BEFORE_EARTHQUAKE + displacement(x(0), x(1));
 
  if (std::abs(x(0)) <= 500 and std::abs(x(1)) <= 500) { // Center of domain is [0, 0]
    Q[Shortcuts::z] = bathymetryAfterEarthquake;
  } else {
    Q[Shortcuts::z] = INITIAL_BATHYMETRY_BEFORE_EARTHQUAKE;
  }
"""
 
boundary_conditions = """
  Qoutside[Shortcuts::h]  = Qinside[Shortcuts::h];
  Qoutside[Shortcuts::hu] = -Qinside[Shortcuts::hu];
  Qoutside[Shortcuts::hv] = -Qinside[Shortcuts::hv];
  Qoutside[Shortcuts::z]  = Qinside[Shortcuts::z];
"""
 
refinement_criterion = """
  auto result = ::exahype2::RefinementCommand::Keep;
  if (std::abs(x(0)) <= 500 and std::abs(x(1)) <= 500) { // Center of domain is [0, 0]
    result = ::exahype2::RefinementCommand::Refine;
  }
  return result;
"""
 
constants = {
    "g": [9.81, "double"],
    "hThreshold": [1e-5, "double"],
}
 
max_h = 1.1 * min(size) / (3.0**min_depth)
min_h = max_h * 3.0 ** (-amr_levels)
 
aderdg_solver = exahype2.solvers.aderdg.GlobalAdaptiveTimeStep(
    name="ADERDGSolver",
    order=dg_order,
    unknowns={"h": 1, "hu": 1, "hv": 1, "z": 1},
    auxiliary_variables=0,
    min_cell_h=min_h,
    max_cell_h=max_h,
    time_step_relaxation=0.5,
)
 
aderdg_solver.set_implementation(
    initial_conditions=initial_conditions,
    boundary_conditions=boundary_conditions,
    refinement_criterion=refinement_criterion,
    flux=flux,
    ncp=nonconservative_product,
    max_eigenvalue=eigenvalue + """
  return sFlux;
""",
)
 
aderdg_solver.add_kernel_optimisations(is_linear=False, polynomials=exahype2.solvers.aderdg.Polynomials.Gauss_Legendre)
 
project = exahype2.Project(
    namespace=["benchmarks", "exahype2", "swe"],
    project_name="ArtificialTsunami",
    directory=".",
    executable=f"ArtificialTsunami",
)
project.add_solver(aderdg_solver)
 
if number_of_snapshots <= 0:
    time_in_between_plots = 0.0
else:
    time_in_between_plots = end_time / number_of_snapshots
    project.set_output_path(output)
 
project.set_global_simulation_parameters(
    dimensions=2,
    size=size,
    offset=offset,
    min_end_time=end_time,
    max_end_time=end_time,
    first_plot_time_stamp=0.0,
    time_in_between_plots=time_in_between_plots,
    periodic_BC=[False, False],
)
 
project.set_load_balancer("new ::exahype2::LoadBalancingConfiguration")
project.set_Peano4_installation("../../../../", mode=peano4.output.string_to_mode(build_mode))
project = project.generate_Peano4_project(verbose=False)
for const_name, const_info in constants.items():
    const_val, const_type = const_info
    project.constants.export_constexpr_with_type(
        const_name, str(const_val), const_type
    )
project.set_fenv_handler(True)
project.build(make=True, make_clean_first=True, throw_away_data_after_build=True)