#include <euler_solver.h>
typedef UnsteadySolver Parent
typedef DifferentiableSystem sys_type
EulerSolver (sys_type &s)
virtual ~EulerSolver ()
virtual Real error_order () const
virtual bool element_residual (bool request_jacobian, DiffContext &)
virtual bool side_residual (bool request_jacobian, DiffContext &)
virtual void init ()
virtual void solve ()
virtual void advance_timestep ()
Number old_nonlinear_solution (const unsigned int global_dof_number) const
virtual Real du (const SystemNorm &norm) const
virtual void reinit ()
virtual void adjoint_solve ()
virtual void adjoint_recede_timestep ()
virtual void before_timestep ()
const sys_type & system () const
virtual AutoPtr< DiffSolver > & diff_solver ()
virtual AutoPtr< DiffSolver > & adjoint_solver ()
static std::string get_info ()
static void print_info ()
static unsigned int n_objects ()
Real theta
AutoPtr< NumericVector< Number > > old_local_nonlinear_solution
bool quiet
unsigned int reduce_deltat_on_diffsolver_failure
typedef std::map< std::string, std::pair< unsigned int, unsigned int > > Counts
sys_type & system ()
void increment_constructor_count (const std::string &name)
void increment_destructor_count (const std::string &name)
bool first_solve
AutoPtr< DiffSolver > _diff_solver
AutoPtr< DiffSolver > _adjoint_solver
sys_type & _system
static Counts _counts
static Threads::atomic< unsigned int > _n_objects
static Threads::spin_mutex _mutex
This class defines a theta-method Euler (defaulting to Backward Euler with theta = 1.0) solver to handle time integration of DifferentiableSystems.
This class is part of the new DifferentiableSystem framework, which is still experimental. Users of this framework should beware of bugs and future API changes.
Author:
Definition at line 44 of file euler_solver.h.
Definition at line 105 of file reference_counter.h.
Definition at line 50 of file euler_solver.h.
Reimplemented in EigenTimeSolver, and SteadySolver.
Definition at line 62 of file time_solver.h.
Definition at line 27 of file euler_solver.C.
: UnsteadySolver(s), theta(1.)
{
}
Definition at line 34 of file euler_solver.C.
{
}
Definition at line 90 of file time_solver.C.
{
}
Definition at line 77 of file time_solver.C.
References TimeSolver::_adjoint_solver.
{
_adjoint_solver->adjoint_solve();
}
Definition at line 157 of file time_solver.h.
References TimeSolver::_adjoint_solver.
{ return _adjoint_solver; }
Reimplemented from TimeSolver.
Reimplemented in AdaptiveTimeSolver.
Definition at line 124 of file unsteady_solver.C.
References TimeSolver::_system, DifferentiableSystem::deltat, UnsteadySolver::first_solve, System::get_vector(), UnsteadySolver::old_nonlinear_solution(), System::solution, and DifferentiableSystem::time.
Referenced by UnsteadySolver::solve().
{
NumericVector<Number> &old_nonlinear_solution =
_system.get_vector('_old_nonlinear_solution');
NumericVector<Number> &nonlinear_solution =
*(_system.solution);
old_nonlinear_solution = nonlinear_solution;
if (!first_solve)
_system.time += _system.deltat;
}
Definition at line 142 of file time_solver.h.
{}
Reimplemented in AdaptiveTimeSolver.
Definition at line 152 of file time_solver.h.
References TimeSolver::_diff_solver.
{ return _diff_solver; }
Note that, while you can always call this function, its result may or may not be very meaningful. For example, if you call this function right after calling advance_timestep() then you'll get a result of zero since old_nonlinear_solution is set equal to nonlinear_solution in this function.
Implements TimeSolver.
Definition at line 150 of file unsteady_solver.C.
References TimeSolver::_system, System::calculate_norm(), System::get_vector(), and System::solution.
{
AutoPtr<NumericVector<Number> > solution_copy =
_system.solution->clone();
solution_copy->add(-1., _system.get_vector('_old_nonlinear_solution'));
solution_copy->close();
return _system.calculate_norm(*solution_copy, norm);
}
Implements TimeSolver.
Definition at line 50 of file euler_solver.C.
References TimeSolver::_system, DenseVector< T >::add(), DifferentiableSystem::deltat, DiffContext::dof_indices, DiffContext::elem_fixed_solution, DiffContext::elem_jacobian, DiffContext::elem_reinit(), DiffContext::elem_residual, DiffContext::elem_solution, DiffContext::elem_solution_derivative, DifferentiableSystem::element_constraint(), DifferentiableSystem::element_time_derivative(), DifferentiableSystem::eulerian_residual(), DiffContext::fixed_solution_derivative, DifferentiableSystem::mass_residual(), UnsteadySolver::old_nonlinear_solution(), DenseVector< T >::size(), DenseMatrix< T >::swap(), DenseVector< T >::swap(), theta, DifferentiableSystem::use_fixed_solution, DenseVector< T >::zero(), and DenseMatrix< T >::zero().
{
unsigned int n_dofs = context.elem_solution.size();
// Local nonlinear solution at old timestep
DenseVector<Number> old_elem_solution(n_dofs);
for (unsigned int i=0; i != n_dofs; ++i)
old_elem_solution(i) =
old_nonlinear_solution(context.dof_indices[i]);
// Local nonlinear solution at time t_theta
DenseVector<Number> theta_solution(context.elem_solution);
theta_solution *= theta;
theta_solution.add(1. - theta, old_elem_solution);
// Technically the elem_solution_derivative is either theta
// or -1.0 in this implementation, but we scale the former part
// ourselves
context.elem_solution_derivative = 1.0;
// Technically the fixed_solution_derivative is always theta,
// but we're scaling a whole jacobian by theta after these first
// evaluations
context.fixed_solution_derivative = 1.0;
// If a fixed solution is requested, we'll use theta_solution
if (_system.use_fixed_solution)
context.elem_fixed_solution = theta_solution;
// Move theta_->elem_, elem_->theta_
context.elem_solution.swap(theta_solution);
// Move the mesh into place first if necessary
context.elem_reinit(theta);
// We're going to compute just the change in elem_residual
// (and possibly elem_jacobian), then add back the old values
DenseVector<Number> old_elem_residual(context.elem_residual);
DenseMatrix<Number> old_elem_jacobian;
if (request_jacobian)
{
old_elem_jacobian = context.elem_jacobian;
context.elem_jacobian.zero();
}
context.elem_residual.zero();
// Get the time derivative at t_theta
bool jacobian_computed =
_system.element_time_derivative(request_jacobian, context);
// For a moving mesh problem we may need the pseudoconvection term too
jacobian_computed =
_system.eulerian_residual(jacobian_computed, context) && jacobian_computed;
// Scale the time-dependent residual and jacobian correctly
context.elem_residual *= _system.deltat;
if (jacobian_computed)
context.elem_jacobian *= (theta * _system.deltat);
// The fixed_solution_derivative is always theta,
// and now we're done scaling jacobians
context.fixed_solution_derivative = theta;
// We evaluate mass_residual with the change in solution
// to get the mass matrix, reusing old_elem_solution to hold that
// delta_solution. We're solving dt*F(u) - du = 0, so our
// delta_solution is old_solution - new_solution.
// We're still keeping elem_solution in theta_solution for now
old_elem_solution -= theta_solution;
// Move old_->elem_, theta_->old_
context.elem_solution.swap(old_elem_solution);
// We do a trick here to avoid using a non-1
// elem_solution_derivative:
context.elem_jacobian *= -1.0;
jacobian_computed = _system.mass_residual(jacobian_computed, context) &&
jacobian_computed;
context.elem_jacobian *= -1.0;
// Move elem_->elem_, old_->theta_
context.elem_solution.swap(theta_solution);
// Restore the elem position if necessary
context.elem_reinit(1.);
// Add the constraint term
jacobian_computed = _system.element_constraint(jacobian_computed, context) &&
jacobian_computed;
// Add back the old residual and jacobian
context.elem_residual += old_elem_residual;
if (request_jacobian)
{
if (jacobian_computed)
context.elem_jacobian += old_elem_jacobian;
else
context.elem_jacobian.swap(old_elem_jacobian);
}
return jacobian_computed;
}
Implements UnsteadySolver.
Definition at line 40 of file euler_solver.C.
References theta.
{
if (theta == 0.5)
return 2.;
return 1.;
}
Definition at line 45 of file reference_counter.C.
References ReferenceCounter::_counts, and Quality::name().
Referenced by ReferenceCounter::print_info().
{
#if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
std::ostringstream out;
out << '
<< ' ----------------------------------------------------------------------------
<< '| Reference count information |
<< ' ---------------------------------------------------------------------------- ;
for (Counts::iterator it = _counts.begin();
it != _counts.end(); ++it)
{
const std::string name(it->first);
const unsigned int creations = it->second.first;
const unsigned int destructions = it->second.second;
out << '| ' << name << ' reference count information:
<< '| Creations: ' << creations << '
<< '| Destructions: ' << destructions << ';
}
out << ' ---------------------------------------------------------------------------- ;
return out.str();
#else
return '';
#endif
}
Definition at line 149 of file reference_counter.h.
References ReferenceCounter::_counts, Quality::name(), and Threads::spin_mtx.
Referenced by ReferenceCountedObject< Value >::ReferenceCountedObject().
{
Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
std::pair<unsigned int, unsigned int>& p = _counts[name];
p.first++;
}
Definition at line 167 of file reference_counter.h.
References ReferenceCounter::_counts, Quality::name(), and Threads::spin_mtx.
Referenced by ReferenceCountedObject< Value >::~ReferenceCountedObject().
{
Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
std::pair<unsigned int, unsigned int>& p = _counts[name];
p.second++;
}
Reimplemented from TimeSolver.
Reimplemented in AdaptiveTimeSolver.
Definition at line 44 of file unsteady_solver.C.
References TimeSolver::_system, and System::add_vector().
{
TimeSolver::init();
_system.add_vector('_old_nonlinear_solution');
}
Definition at line 76 of file reference_counter.h.
References ReferenceCounter::_n_objects.
Referenced by System::read_serialized_blocked_dof_objects(), and System::write_serialized_blocked_dof_objects().
{ return _n_objects; }
Definition at line 139 of file unsteady_solver.C.
References TimeSolver::_system, System::get_dof_map(), DofMap::n_dofs(), and UnsteadySolver::old_local_nonlinear_solution.
Referenced by UnsteadySolver::advance_timestep(), AdaptiveTimeSolver::advance_timestep(), element_residual(), Euler2Solver::element_residual(), FEMSystem::eulerian_residual(), side_residual(), and Euler2Solver::side_residual().
{
libmesh_assert (global_dof_number < _system.get_dof_map().n_dofs());
libmesh_assert (global_dof_number < old_local_nonlinear_solution->size());
return (*old_local_nonlinear_solution)(global_dof_number);
}
Definition at line 83 of file reference_counter.C.
References ReferenceCounter::get_info().
{
#if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
std::cout << ReferenceCounter::get_info();
#endif
}
Reimplemented in AdaptiveTimeSolver, and EigenTimeSolver.
Definition at line 46 of file time_solver.C.
References TimeSolver::_adjoint_solver, and TimeSolver::_diff_solver.
{
_diff_solver->reinit();
_adjoint_solver->reinit();
}
Implements TimeSolver.
Definition at line 156 of file euler_solver.C.
References TimeSolver::_system, DenseVector< T >::add(), DifferentiableSystem::deltat, DiffContext::dof_indices, DiffContext::elem_fixed_solution, DiffContext::elem_jacobian, DiffContext::elem_residual, DiffContext::elem_side_reinit(), DiffContext::elem_solution, DiffContext::elem_solution_derivative, DiffContext::fixed_solution_derivative, UnsteadySolver::old_nonlinear_solution(), DifferentiableSystem::side_constraint(), DifferentiableSystem::side_mass_residual(), DifferentiableSystem::side_time_derivative(), DenseVector< T >::size(), DenseMatrix< T >::swap(), DenseVector< T >::swap(), theta, DifferentiableSystem::use_fixed_solution, DenseVector< T >::zero(), and DenseMatrix< T >::zero().
{
unsigned int n_dofs = context.elem_solution.size();
// Local nonlinear solution at old timestep
DenseVector<Number> old_elem_solution(n_dofs);
for (unsigned int i=0; i != n_dofs; ++i)
old_elem_solution(i) =
old_nonlinear_solution(context.dof_indices[i]);
// Local nonlinear solution at time t_theta
DenseVector<Number> theta_solution(context.elem_solution);
theta_solution *= theta;
theta_solution.add(1. - theta, old_elem_solution);
// Technically the elem_solution_derivative is either theta
// or 1.0 in this implementation, but we scale the former part
// ourselves
context.elem_solution_derivative = 1.0;
// Technically the fixed_solution_derivative is always theta,
// but we're scaling a whole jacobian by theta after these first
// evaluations
context.fixed_solution_derivative = 1.0;
// If a fixed solution is requested, we'll use theta_solution
if (_system.use_fixed_solution)
context.elem_fixed_solution = theta_solution;
// Move theta_->elem_, elem_->theta_
context.elem_solution.swap(theta_solution);
// Move the mesh into place first if necessary
context.elem_side_reinit(theta);
// We're going to compute just the change in elem_residual
// (and possibly elem_jacobian), then add back the old values
DenseVector<Number> old_elem_residual(context.elem_residual);
DenseMatrix<Number> old_elem_jacobian;
if (request_jacobian)
{
old_elem_jacobian = context.elem_jacobian;
context.elem_jacobian.zero();
}
context.elem_residual.zero();
// Get the time derivative at t_theta
bool jacobian_computed =
_system.side_time_derivative(request_jacobian, context);
// Scale the time-dependent residual and jacobian correctly
context.elem_residual *= _system.deltat;
if (jacobian_computed)
context.elem_jacobian *= (theta * _system.deltat);
// The fixed_solution_derivative is always theta,
// and now we're done scaling jacobians
context.fixed_solution_derivative = theta;
// We evaluate side_mass_residual with the change in solution
// to get the mass matrix, reusing old_elem_solution to hold that
// delta_solution. We're solving dt*F(u) - du = 0, so our
// delta_solution is old_solution - new_solution.
// We're still keeping elem_solution in theta_solution for now
old_elem_solution -= theta_solution;
// Move old_->elem_, theta_->old_
context.elem_solution.swap(old_elem_solution);
// We do a trick here to avoid using a non-1
// elem_solution_derivative:
context.elem_jacobian *= -1.0;
jacobian_computed = _system.side_mass_residual(jacobian_computed, context) &&
jacobian_computed;
context.elem_jacobian *= -1.0;
// Move elem_->elem_, old_->theta_
context.elem_solution.swap(theta_solution);
// Restore the elem position if necessary
context.elem_side_reinit(1.);
// Add the constraint term
jacobian_computed = _system.side_constraint(jacobian_computed, context) &&
jacobian_computed;
// Add back the old residual and jacobian
context.elem_residual += old_elem_residual;
if (request_jacobian)
{
if (jacobian_computed)
context.elem_jacobian += old_elem_jacobian;
else
context.elem_jacobian.swap(old_elem_jacobian);
}
return jacobian_computed;
}
Reimplemented from TimeSolver.
Reimplemented in AdaptiveTimeSolver, and TwostepTimeSolver.
Definition at line 53 of file unsteady_solver.C.
References TimeSolver::_diff_solver, TimeSolver::_system, UnsteadySolver::advance_timestep(), DifferentiableSystem::deltat, DiffSolver::DIVERGED_BACKTRACKING_FAILURE, DiffSolver::DIVERGED_MAX_NONLINEAR_ITERATIONS, UnsteadySolver::first_solve, System::get_dof_map(), DofMap::get_send_list(), System::get_vector(), libMeshEnums::GHOSTED, NumericVector< T >::localize(), System::n_dofs(), System::n_local_dofs(), UnsteadySolver::old_local_nonlinear_solution, TimeSolver::quiet, TimeSolver::reduce_deltat_on_diffsolver_failure, and libMeshEnums::SERIAL.
{
if (first_solve)
{
advance_timestep();
first_solve = false;
}
#ifdef LIBMESH_ENABLE_GHOSTED
old_local_nonlinear_solution->init (_system.n_dofs(), _system.n_local_dofs(),
_system.get_dof_map().get_send_list(), false,
GHOSTED);
#else
old_local_nonlinear_solution->init (_system.n_dofs(), false, SERIAL);
#endif
_system.get_vector('_old_nonlinear_solution').localize
(*old_local_nonlinear_solution,
_system.get_dof_map().get_send_list());
unsigned int solve_result = _diff_solver->solve();
// If we requested the UnsteadySolver to attempt reducing dt after a
// failed DiffSolver solve, check the results of the solve now.
if (reduce_deltat_on_diffsolver_failure)
{
bool backtracking_failed =
solve_result & DiffSolver::DIVERGED_BACKTRACKING_FAILURE;
bool max_iterations =
solve_result & DiffSolver::DIVERGED_MAX_NONLINEAR_ITERATIONS;
if (backtracking_failed || max_iterations)
{
// Cut timestep in half
for (unsigned int nr=0; nr<reduce_deltat_on_diffsolver_failure; ++nr)
{
_system.deltat *= 0.5;
std::cout << 'Newton backtracking failed. Trying with smaller timestep, dt='
<< _system.deltat << std::endl;
solve_result = _diff_solver->solve();
// Check solve results with reduced timestep
bool backtracking_failed =
solve_result & DiffSolver::DIVERGED_BACKTRACKING_FAILURE;
bool max_iterations =
solve_result & DiffSolver::DIVERGED_MAX_NONLINEAR_ITERATIONS;
if (!backtracking_failed && !max_iterations)
{
if (!quiet)
std::cout << 'Reduced dt solve succeeded.' << std::endl;
return;
}
}
// If we made it here, we still couldn't converge the solve after
// reducing deltat
std::cout << 'DiffSolver::solve() did not succeed after '
<< reduce_deltat_on_diffsolver_failure
<< ' attempts.' << std::endl;
libmesh_convergence_failure();
} // end if (backtracking_failed || max_iterations)
} // end if (reduce_deltat_on_diffsolver_failure)
}
Definition at line 202 of file time_solver.h.
References TimeSolver::_system.
{ return _system; }
Definition at line 147 of file time_solver.h.
References TimeSolver::_system.
{ return _system; }
Definition at line 197 of file time_solver.h.
Referenced by TimeSolver::adjoint_solve(), TimeSolver::adjoint_solver(), TimeSolver::init(), and TimeSolver::reinit().
Definition at line 110 of file reference_counter.h.
Referenced by ReferenceCounter::get_info(), ReferenceCounter::increment_constructor_count(), and ReferenceCounter::increment_destructor_count().
Definition at line 192 of file time_solver.h.
Referenced by TimeSolver::diff_solver(), TimeSolver::init(), TimeSolver::reinit(), UnsteadySolver::solve(), and TimeSolver::solve().
Definition at line 123 of file reference_counter.h.
Definition at line 118 of file reference_counter.h.
Referenced by ReferenceCounter::n_objects(), ReferenceCounter::ReferenceCounter(), and ReferenceCounter::~ReferenceCounter().
Definition at line 207 of file time_solver.h.
Referenced by UnsteadySolver::advance_timestep(), AdaptiveTimeSolver::advance_timestep(), UnsteadySolver::du(), SteadySolver::element_residual(), element_residual(), Euler2Solver::element_residual(), EigenTimeSolver::element_residual(), UnsteadySolver::init(), TimeSolver::init(), EigenTimeSolver::init(), UnsteadySolver::old_nonlinear_solution(), SteadySolver::side_residual(), side_residual(), Euler2Solver::side_residual(), EigenTimeSolver::side_residual(), UnsteadySolver::solve(), TwostepTimeSolver::solve(), EigenTimeSolver::solve(), and TimeSolver::system().
Reimplemented from TimeSolver.
Definition at line 124 of file unsteady_solver.h.
Referenced by UnsteadySolver::advance_timestep(), AdaptiveTimeSolver::advance_timestep(), UnsteadySolver::solve(), and TwostepTimeSolver::solve().
Reimplemented from TimeSolver.
Definition at line 105 of file unsteady_solver.h.
Referenced by AdaptiveTimeSolver::AdaptiveTimeSolver(), AdaptiveTimeSolver::init(), UnsteadySolver::old_nonlinear_solution(), UnsteadySolver::solve(), and AdaptiveTimeSolver::~AdaptiveTimeSolver().
Definition at line 162 of file time_solver.h.
Referenced by UnsteadySolver::solve(), TwostepTimeSolver::solve(), and EigenTimeSolver::solve().
Definition at line 185 of file time_solver.h.
Referenced by UnsteadySolver::solve(), and TwostepTimeSolver::solve().
Definition at line 91 of file euler_solver.h.
Referenced by element_residual(), error_order(), and side_residual().
Generated automatically by Doxygen for libMesh from the source code.