本文整理汇总了C++中DiffContext::get_elem_solution方法的典型用法代码示例。如果您正苦于以下问题:C++ DiffContext::get_elem_solution方法的具体用法?C++ DiffContext::get_elem_solution怎么用?C++ DiffContext::get_elem_solution使用的例子?那么, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在类DiffContext
的用法示例。
在下文中一共展示了DiffContext::get_elem_solution方法的4个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的C++代码示例。
示例1:
bool Euler2Solver::_general_residual (bool request_jacobian,
DiffContext & context,
ResFuncType mass,
ResFuncType damping,
ResFuncType time_deriv,
ResFuncType constraint,
ReinitFuncType reinit_func,
bool compute_second_order_eqns)
{
unsigned int n_dofs = context.get_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.get_dof_indices()[i]);
// Local time derivative of solution
context.get_elem_solution_rate() = context.get_elem_solution();
context.get_elem_solution_rate() -= old_elem_solution;
context.elem_solution_rate_derivative = 1 / _system.deltat;
context.get_elem_solution_rate() *=
context.elem_solution_rate_derivative;
// Our first evaluations are at the final elem_solution
context.elem_solution_derivative = 1.0;
// If a fixed solution is requested, we'll use the elem_solution
// at the new timestep
// FIXME - should this be the theta solution instead?
if (_system.use_fixed_solution)
context.get_elem_fixed_solution() = context.get_elem_solution();
context.fixed_solution_derivative = 1.0;
// We need to save the old jacobian and old residual since we'll be
// multiplying some of the new contributions by theta or 1-theta
DenseMatrix<Number> old_elem_jacobian(n_dofs, n_dofs);
DenseVector<Number> old_elem_residual(n_dofs);
old_elem_residual.swap(context.get_elem_residual());
if (request_jacobian)
old_elem_jacobian.swap(context.get_elem_jacobian());
// Local time derivative of solution
context.get_elem_solution_rate() = context.get_elem_solution();
context.get_elem_solution_rate() -= old_elem_solution;
context.elem_solution_rate_derivative = 1 / _system.deltat;
context.get_elem_solution_rate() *=
context.elem_solution_rate_derivative;
// If we are asked to compute residuals for second order variables,
// we also populate the acceleration part so the user can use that.
if (compute_second_order_eqns)
this->prepare_accel(context);
// Move the mesh into place first if necessary, set t = t_{n+1}
(context.*reinit_func)(1.);
// First, evaluate time derivative at the new timestep.
// The element should already be in the proper place
// even for a moving mesh problem.
bool jacobian_computed =
(_system.get_physics()->*time_deriv)(request_jacobian, context);
// Next, evaluate the mass residual at the new timestep
jacobian_computed = (_system.get_physics()->*mass)(jacobian_computed, context) &&
jacobian_computed;
// If we have second-order variables, we need to get damping terms
// and the velocity equations
if (compute_second_order_eqns)
{
jacobian_computed = (_system.get_physics()->*damping)(jacobian_computed, context) &&
jacobian_computed;
jacobian_computed = this->compute_second_order_eqns(jacobian_computed, context) &&
jacobian_computed;
}
// Add the constraint term
jacobian_computed = (_system.get_physics()->*constraint)(jacobian_computed, context) &&
jacobian_computed;
// The new solution's contribution is scaled by theta
context.get_elem_residual() *= theta;
context.get_elem_jacobian() *= theta;
// Save the new solution's term
DenseMatrix<Number> elem_jacobian_newterm(n_dofs, n_dofs);
DenseVector<Number> elem_residual_newterm(n_dofs);
elem_residual_newterm.swap(context.get_elem_residual());
if (request_jacobian)
elem_jacobian_newterm.swap(context.get_elem_jacobian());
// Add the time-dependent term for the old solution
// Make sure elem_solution is set up for elem_reinit to use
// Move elem_->old_, old_->elem_
context.get_elem_solution().swap(old_elem_solution);
//.........这里部分代码省略.........
示例2: element_residual
bool EulerSolver::element_residual (bool request_jacobian,
DiffContext &context)
{
unsigned int n_dofs = context.get_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.get_dof_indices()[i]);
// Local nonlinear solution at time t_theta
DenseVector<Number> theta_solution(context.get_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.get_elem_fixed_solution() = theta_solution;
// Move theta_->elem_, elem_->theta_
context.get_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.get_elem_residual());
DenseMatrix<Number> old_elem_jacobian;
if (request_jacobian)
{
old_elem_jacobian = context.get_elem_jacobian();
context.get_elem_jacobian().zero();
}
context.get_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.get_elem_residual() *= _system.deltat;
if (jacobian_computed)
context.get_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.get_elem_solution().swap(old_elem_solution);
// We do a trick here to avoid using a non-1
// elem_solution_derivative:
context.get_elem_jacobian() *= -1.0;
context.fixed_solution_derivative *= -1.0;
jacobian_computed = _system.mass_residual(jacobian_computed, context) &&
jacobian_computed;
context.get_elem_jacobian() *= -1.0;
context.fixed_solution_derivative *= -1.0;
// Move elem_->elem_, old_->theta_
context.get_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.get_elem_residual() += old_elem_residual;
if (request_jacobian)
{
if (jacobian_computed)
context.get_elem_jacobian() += old_elem_jacobian;
else
//.........这里部分代码省略.........
示例3: _general_residual
bool NewmarkSolver::_general_residual (bool request_jacobian,
DiffContext & context,
ResFuncType mass,
ResFuncType damping,
ResFuncType time_deriv,
ResFuncType constraint)
{
unsigned int n_dofs = context.get_elem_solution().size();
// We might need to save the old jacobian in case one of our physics
// terms later is unable to update it analytically.
DenseMatrix<Number> old_elem_jacobian(n_dofs, n_dofs);
// Local velocity at old time step
DenseVector<Number> old_elem_solution_rate(n_dofs);
for (unsigned int i=0; i != n_dofs; ++i)
old_elem_solution_rate(i) =
old_solution_rate(context.get_dof_indices()[i]);
// The user is computing the initial acceleration
// So upstream we've swapped _system.solution and _old_local_solution_accel
// So we need to give the context the correct entries since we're solving for
// acceleration here.
if( _is_accel_solve )
{
// System._solution is actually the acceleration right now so we need
// to reset the elem_solution to the right thing, which in this case
// is the initial guess for displacement, which got swapped into
// _old_solution_accel vector
DenseVector<Number> old_elem_solution(n_dofs);
for (unsigned int i=0; i != n_dofs; ++i)
old_elem_solution(i) =
old_solution_accel(context.get_dof_indices()[i]);
context.elem_solution_derivative = 0.0;
context.elem_solution_rate_derivative = 0.0;
context.elem_solution_accel_derivative = 1.0;
// Acceleration is currently the unknown so it's already sitting
// in elem_solution() thanks to FEMContext::pre_fe_reinit
context.get_elem_solution_accel() = context.get_elem_solution();
// Now reset elem_solution() to what the user is expecting
context.get_elem_solution() = old_elem_solution;
context.get_elem_solution_rate() = old_elem_solution_rate;
// The user's Jacobians will be targeting derivatives w.r.t. u_{n+1}.
// Although the vast majority of cases will have the correct analytic
// Jacobians in this iteration, since we reset elem_solution_derivative*,
// if there are coupled/overlapping problems, there could be
// mismatches in the Jacobian. So we force finite differencing for
// the first iteration.
request_jacobian = false;
}
// Otherwise, the unknowns are the displacements and everything is straight
// foward and is what you think it is
else
{
if (request_jacobian)
old_elem_jacobian.swap(context.get_elem_jacobian());
// Local displacement 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.get_dof_indices()[i]);
// Local acceleration at old time step
DenseVector<Number> old_elem_solution_accel(n_dofs);
for (unsigned int i=0; i != n_dofs; ++i)
old_elem_solution_accel(i) =
old_solution_accel(context.get_dof_indices()[i]);
// Convenience
libMesh::Real dt = _system.deltat;
context.elem_solution_derivative = 1.0;
// Local velocity at current time step
// v_{n+1} = gamma/(beta*Delta t)*(x_{n+1}-x_n)
// + (1-(gamma/beta))*v_n
// + (1-gamma/(2*beta))*(Delta t)*a_n
context.elem_solution_rate_derivative = (_gamma/(_beta*dt));
context.get_elem_solution_rate() = context.get_elem_solution();
context.get_elem_solution_rate() -= old_elem_solution;
context.get_elem_solution_rate() *= context.elem_solution_rate_derivative;
context.get_elem_solution_rate().add( (1.0-_gamma/_beta), old_elem_solution_rate);
context.get_elem_solution_rate().add( (1.0-_gamma/(2.0*_beta))*dt, old_elem_solution_accel);
// Local acceleration at current time step
// a_{n+1} = (1/(beta*(Delta t)^2))*(x_{n+1}-x_n)
// - 1/(beta*Delta t)*v_n
// - (1/(2*beta)-1)*a_n
context.elem_solution_accel_derivative = 1.0/(_beta*dt*dt);
context.get_elem_solution_accel() = context.get_elem_solution();
//.........这里部分代码省略.........
示例4: _general_residual
bool EulerSolver::_general_residual (bool request_jacobian,
DiffContext &context,
ResFuncType mass,
ResFuncType time_deriv,
ResFuncType constraint,
ReinitFuncType reinit_func)
{
unsigned int n_dofs = context.get_elem_solution().size();
// We might need to save the old jacobian in case one of our physics
// terms later is unable to update it analytically.
DenseMatrix<Number> old_elem_jacobian(n_dofs, n_dofs);
if (request_jacobian)
old_elem_jacobian.swap(context.get_elem_jacobian());
// 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.get_dof_indices()[i]);
// Local time derivative of solution
context.get_elem_solution_rate() = context.get_elem_solution();
context.get_elem_solution_rate() -= old_elem_solution;
context.elem_solution_rate_derivative = 1 / _system.deltat;
context.get_elem_solution_rate() *=
context.elem_solution_rate_derivative;
// Local nonlinear solution at time t_theta
DenseVector<Number> theta_solution(context.get_elem_solution());
theta_solution *= theta;
theta_solution.add(1. - theta, old_elem_solution);
context.elem_solution_derivative = theta;
context.fixed_solution_derivative = theta;
// If a fixed solution is requested, we'll use theta_solution
if (_system.use_fixed_solution)
context.get_elem_fixed_solution() = theta_solution;
// Move theta_->elem_, elem_->theta_
context.get_elem_solution().swap(theta_solution);
// Move the mesh into place first if necessary
(context.*reinit_func)(theta);
// Get the time derivative at t_theta
bool jacobian_computed =
(_system.*time_deriv)(request_jacobian, context);
jacobian_computed = (_system.*mass)(jacobian_computed, context) &&
jacobian_computed;
// Restore the elem position if necessary
(context.*reinit_func)(1);
// Move elem_->elem_, theta_->theta_
context.get_elem_solution().swap(theta_solution);
context.elem_solution_derivative = 1;
// Add the constraint term
jacobian_computed = (_system.*constraint)(jacobian_computed, context) &&
jacobian_computed;
// Add back (or restore) the old jacobian
if (request_jacobian)
{
if (jacobian_computed)
context.get_elem_jacobian() += old_elem_jacobian;
else
context.get_elem_jacobian().swap(old_elem_jacobian);
}
return jacobian_computed;
}