|
Reference documentation for deal.II version 9.2.0
|
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\)
\(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\)
\(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\)
\(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Go to the documentation of this file.
16 #ifndef dealii_eigen_h
17 #define dealii_eigen_h
54 template <
typename VectorType = Vector<
double>>
86 const AdditionalData & data = AdditionalData());
95 template <
typename MatrixType>
129 template <
typename VectorType = Vector<
double>>
173 const AdditionalData & data = AdditionalData());
182 template <
typename MatrixType>
197 template <
class VectorType>
202 , additional_data(data)
207 template <
class VectorType>
208 template <
typename MatrixType>
223 double length = x.l2_norm();
224 double old_length = 0.;
233 y.add(additional_data.shift, x);
237 length = y.l2_norm();
243 double thresh = length / x.size();
254 value = (entry * x(i) < 0.) ? -length : length;
255 value -= additional_data.shift;
258 x.equ(1 / length, y);
265 conv = this->iteration_status(iter,
266 std::fabs(1. / length - 1. / old_length),
273 iter,
std::fabs(1. / length - 1. / old_length)));
280 template <
class VectorType>
285 , additional_data(data)
290 template <
class VectorType>
291 template <
typename MatrixType>
303 double current_shift = -
value;
312 unsigned int goal = additional_data.start_adaption;
322 double length = x.l2_norm();
323 double old_value =
value;
332 solver.
solve(A_s, y, x, prec);
335 length = y.l2_norm();
341 double thresh = length / x.size();
352 value = (entry * x(i) < 0. ? -1. : 1.) / length - current_shift;
356 const auto & relaxation = additional_data.relaxation;
357 const double new_shift =
358 relaxation * (-
value) + (1. - relaxation) * current_shift;
361 current_shift = new_shift;
367 x.equ(1. / length, y);
369 if (additional_data.use_residual)
373 r.sadd(-1.,
value, x);
376 conv = this->iteration_status(iter, res, x);
381 conv = this->iteration_status(iter, res, x);
LinearOperator< Range, Range, Payload > identity_operator(const std::function< void(Range &, bool)> &reinit_vector)
EigenPower(SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
AdditionalData(double relaxation=1., unsigned int start_adaption=6, bool use_residual=true)
unsigned int start_adaption
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
AdditionalData additional_data
Expression fabs(const Expression &x)
EigenInverse(SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
@ iterate
Continue iteration.
LinearOperator< Range, Domain, Payload > linear_operator(const Matrix &matrix)
void solve(const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &preconditioner)
unsigned int global_dof_index
#define DEAL_II_NAMESPACE_OPEN
AdditionalData(const double shift=0.)
void solve(double &value, const MatrixType &A, VectorType &x)
void solve(double &value, const MatrixType &A, VectorType &x)
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
AdditionalData additional_data
@ success
Stop iteration, goal reached.
#define DEAL_II_NAMESPACE_CLOSE
#define AssertThrow(cond, exc)
T max(const T &t, const MPI_Comm &mpi_communicator)