Homec4science
Diffusion SpecMiCP / ReactMiCP (specmicp_bin) Browse

specmicp/tests/reactmicp/systems/saturated_react/transport_program.cpp
a4ba90c69915specmicp_bin

Last Change
Raw File

transport_program.cpp
View Options

#include "catch.hpp"
#include "speciation_system.hpp"
#include "specmicp/adimensional/adimensional_system_solution.hpp"
#include "specmicp/adimensional/adimensional_system_solution_extractor.hpp"
#include "dfpm/mesh.hpp"
#include "reactmicp/systems/saturated_react/variables.hpp"
#include "reactmicp/systems/saturated_react/init_variables.hpp"
#include "reactmicp/systems/saturated_react/transport_program.hpp"
#include "dfpmsolver/parabolic_driver.hpp"
#include "reactmicp/systems/saturated_react/transport_stagger.hpp"
TEST_CASE("Saturated transport program", "[SaturatedReact, transport, program, stagger]") {
specmicp::units::UnitsSet the_units;
the_units.mass = specmicp::units::MassUnit::kilogram;
the_units.length = specmicp::units::LengthUnit::centimeter;
specmicp::mesh::Mesh1DPtr the_mesh = specmicp::mesh::uniform_mesh1d(5, 0.1, 5.0);
specmicp::RawDatabasePtr the_database = leaching_database();
std::vector<specmicp::AdimensionalSystemSolution> list_initial_composition;
list_initial_composition.push_back(initial_leaching_pb(the_database, 0.005, the_units));
list_initial_composition.push_back(initial_blank_leaching_pb(the_database, 0.005, the_units));
std::vector<specmicp::index_t> list_fixed_nodes = {0, };
std::vector<int> index_initial_state = {1, 0, 0, 0, 0};
specmicp::reactmicp::systems::satdiff::SaturatedVariablesPtr variables =
specmicp::reactmicp::systems::satdiff::init_variables(the_mesh, the_database, the_units,
list_fixed_nodes,
list_initial_composition,
index_initial_state);
for (specmicp::index_t node: the_mesh->range_nodes())
{
variables->porosity(node) = 0.2;
variables->diffusion_coefficient(node) = 1.0e-8;
}
SECTION("Initialisation") {
specmicp::reactmicp::systems::satdiff::SaturatedDiffusion program(variables, list_fixed_nodes);
specmicp::Vector residuals;
REQUIRE(program.get_ndf() == 2*4);
REQUIRE(program.get_tot_ndf() == 2*4*5);
REQUIRE(program.get_neq() == 4*4);
program.compute_residuals(variables->displacement(), variables->velocity(), residuals);
}
SECTION("Solving the problem manually") {
specmicp::reactmicp::systems::satdiff::SaturatedDiffusion program(variables, list_fixed_nodes);
specmicp::dfpmsolver::ParabolicDriver<specmicp::reactmicp::systems::satdiff::SaturatedDiffusion> solver(program);
variables->predictor() = variables->displacement();
int retcode = (int) solver.solve_timestep(100, variables->displacement());
REQUIRE(retcode > 0);
}
SECTION("Transport stagger") {
specmicp::reactmicp::systems::satdiff::SaturatedTransportStagger stagger(variables, list_fixed_nodes);
stagger.initialize_timestep(100, variables);
specmicp::reactmicp::solver::StaggerReturnCode retcode = stagger.restart_timestep(variables);
REQUIRE((int) retcode > 0);
}
}
c4science · Help