Page MenuHomec4science

transport_stagger.cpp
No OneTemporary

File Metadata

Created
Fri, Jul 19, 10:49

transport_stagger.cpp

#include "transport_stagger.hpp"
#include "variables.hpp"
#include "reactmicp/solver/staggers_base/stagger_structs.hpp"
namespace specmicp {
namespace reactmicp {
namespace systems {
namespace satdiff {
// class SaturatedTransportStagger::
namespace internal {
//! \brief Translate a ParabolicDriverReturnCode to a StaggerReturnCode
solver::StaggerReturnCode Parabolic2StaggerReturnCode(dfpmsolver::ParabolicDriverReturnCode retcode)
{
using ParRC = dfpmsolver::ParabolicDriverReturnCode;
using StaRC = solver::StaggerReturnCode;
switch (retcode)
{
case ParRC::ResidualMinimized:
return StaRC::ResidualMinimized;
break;
case ParRC::ErrorMinimized:
return StaRC::ErrorMinimized;
break;
default: // error
switch(retcode)
{
case ParRC::MaxIterations:
return StaRC::MaximumIterationsReached;
break;
case ParRC::StationaryPoint:
return StaRC::StationaryPoint;
break;
case ParRC::NotConvergedYet:
return StaRC::NotConvergedYet;
break;
default:
return StaRC::UnknownError;
}
}
}
} // end namespace internal
SaturatedTransportStagger::SaturatedTransportStagger(SaturatedVariablesPtr variables,
std::vector<index_t> list_fixed_nodes):
m_program(variables, list_fixed_nodes),
m_solver(m_program)
{
}
SaturatedTransportStagger::SaturatedTransportStagger(
SaturatedVariablesPtr variables,
std::vector<index_t> list_fixed_nodes,
std::map<index_t, index_t> list_slave_nodes,
std::vector<index_t> list_immobile_components):
m_program(variables, list_fixed_nodes, list_slave_nodes, list_immobile_components),
m_solver(m_program)
{
}
//! \brief Initialize the stagger at the beginning of an iteration
void SaturatedTransportStagger::initialize_timestep(scalar_t dt, VariablesBasePtr var)
{
m_dt = dt;
SaturatedVariablesPtr true_var = cast_var_from_base(var);
true_var->predictor() = true_var->displacement();
true_var->velocity().setZero();
true_var->transport_rate().setZero();
m_solver.initialize_timestep(dt, true_var->displacement());
Vector tmp = true_var->chemistry_rate();
true_var->chemistry_rate().setZero();
Eigen::VectorXd residuals;
m_program.compute_residuals(true_var->displacement(), true_var->velocity(), residuals);
m_residual_0 = residuals.norm();
true_var->chemistry_rate() = tmp;
m_program.set_variables(true_var);
}
//! \brief Solve the equation for the timestep
solver::StaggerReturnCode SaturatedTransportStagger::restart_timestep(VariablesBasePtr var)
{
SaturatedVariablesPtr true_var = cast_var_from_base(var);
m_solver.velocity() = true_var->velocity();
dfpmsolver::ParabolicDriverReturnCode retcode = m_solver.restart_timestep(true_var->displacement());
// copy variables if successful
if (retcode > dfpmsolver::ParabolicDriverReturnCode::NotConvergedYet)
{
true_var->velocity() = m_solver.velocity();
}
return internal::Parabolic2StaggerReturnCode(retcode);
}
scalar_t SaturatedTransportStagger::get_update(VariablesBasePtr var)
{
return cast_var_from_base(var)->velocity().norm();
}
//! \brief Compute the residuals norm
scalar_t SaturatedTransportStagger::get_residual(VariablesBasePtr var)
{
SaturatedVariablesPtr true_var = cast_var_from_base(var);
Eigen::VectorXd residuals;
m_program.compute_residuals(true_var->displacement(), true_var->velocity(), residuals);
return residuals.norm();
}
} // end namespace satdiff
} // end namespace systems
} // end namespace reactmicp
} // end namespace specmicp

Event Timeline