*
* @author Guido Kanschat, 2006, 2010
*/
- template <class VECTOR>
- class Newton : public Operator<VECTOR>
+ template <typename VectorType>
+ class Newton : public Operator<VectorType>
{
public:
/**
* Constructor, receiving the applications computing the residual and
* solving the linear problem, respectively.
*/
- Newton (Operator<VECTOR> &residual, Operator<VECTOR> &inverse_derivative);
+ Newton (Operator<VectorType> &residual, Operator<VectorType> &inverse_derivative);
/**
* Declare the parameters applicable to Newton's method.
/**
* Initialize the pointer data_out for debugging.
*/
- void initialize (OutputOperator<VECTOR> &output);
+ void initialize (OutputOperator<VectorType> &output);
/**
* The actual Newton iteration. The initial value is in <tt>out(0)</tt>,
/**
* The operator computing the residual.
*/
- SmartPointer<Operator<VECTOR>, Newton<VECTOR> > residual;
+ SmartPointer<Operator<VectorType>, Newton<VectorType> > residual;
/**
* The operator applying the inverse derivative to the residual.
*/
- SmartPointer<Operator<VECTOR>, Newton<VECTOR> > inverse_derivative;
+ SmartPointer<Operator<VectorType>, Newton<VectorType> > inverse_derivative;
/**
* The operator handling the output in case the debug_vectors is true.
* Call the initialize function first.
*/
- SmartPointer<OutputOperator<VECTOR>, Newton<VECTOR> > data_out;
+ SmartPointer<OutputOperator<VectorType>, Newton<VectorType> > data_out;
/**
* This flag is set by the function assemble(), indicating that the matrix
namespace Algorithms
{
- template <class VECTOR>
- Newton<VECTOR>::Newton(Operator<VECTOR> &residual, Operator<VECTOR> &inverse_derivative)
+ template <typename VectorType>
+ Newton<VectorType>::Newton(Operator<VectorType> &residual, Operator<VectorType> &inverse_derivative)
:
residual(&residual), inverse_derivative(&inverse_derivative),
assemble_now(false),
{}
- template <class VECTOR>
+ template <typename VectorType>
void
- Newton<VECTOR>::declare_parameters(ParameterHandler ¶m)
+ Newton<VectorType>::declare_parameters(ParameterHandler ¶m)
{
param.enter_subsection("Newton");
ReductionControl::declare_parameters (param);
param.leave_subsection();
}
- template <class VECTOR>
+ template <typename VectorType>
void
- Newton<VECTOR>::parse_parameters (ParameterHandler ¶m)
+ Newton<VectorType>::parse_parameters (ParameterHandler ¶m)
{
param.enter_subsection("Newton");
control.parse_parameters (param);
param.leave_subsection ();
}
- template <class VECTOR>
+ template <typename VectorType>
void
- Newton<VECTOR>::initialize (OutputOperator<VECTOR> &output)
+ Newton<VectorType>::initialize (OutputOperator<VectorType> &output)
{
data_out = &output;
}
- template <class VECTOR>
+ template <typename VectorType>
void
- Newton<VECTOR>::notify(const Event &e)
+ Newton<VectorType>::notify(const Event &e)
{
residual->notify(e);
inverse_derivative->notify(e);
}
- template <class VECTOR>
+ template <typename VectorType>
double
- Newton<VECTOR>::threshold(const double thr)
+ Newton<VectorType>::threshold(const double thr)
{
const double t = assemble_threshold;
assemble_threshold = thr;
}
- template <class VECTOR>
+ template <typename VectorType>
void
- Newton<VECTOR>::operator() (AnyData &out, const AnyData &in)
+ Newton<VectorType>::operator() (AnyData &out, const AnyData &in)
{
Assert (out.size() == 1, ExcNotImplemented());
deallog.push ("Newton");
- VECTOR &u = *out.entry<VECTOR *>(0);
+ VectorType &u = *out.entry<VectorType *>(0);
if (debug>2)
deallog << "u: " << u.l2_norm() << std::endl;
- GrowingVectorMemory<VECTOR> mem;
- typename VectorMemory<VECTOR>::Pointer Du(mem);
- typename VectorMemory<VECTOR>::Pointer res(mem);
+ GrowingVectorMemory<VectorType> mem;
+ typename VectorMemory<VectorType>::Pointer Du(mem);
+ typename VectorMemory<VectorType>::Pointer res(mem);
res->reinit(u);
AnyData src1;
AnyData src2;
- src1.add<const VECTOR *>(&u, "Newton iterate");
+ src1.add<const VectorType *>(&u, "Newton iterate");
src1.merge(in);
- src2.add<const VECTOR *>(res, "Newton residual");
+ src2.add<const VectorType *>(res, "Newton residual");
src2.merge(src1);
AnyData out1;
- out1.add<VECTOR *>(res, "Residual");
+ out1.add<VectorType *>(res, "Residual");
AnyData out2;
- out2.add<VECTOR *>(Du, "Update");
+ out2.add<VectorType *>(Du, "Update");
unsigned int step = 0;
// fill res with (f(u), v)
if (debug_vectors)
{
AnyData out;
- VECTOR *p = &u;
- out.add<const VECTOR *>(p, "solution");
+ VectorType *p = &u;
+ out.add<const VectorType *>(p, "solution");
p = Du;
- out.add<const VECTOR *>(p, "update");
+ out.add<const VectorType *>(p, "update");
p = res;
- out.add<const VECTOR *>(p, "residual");
+ out.add<const VectorType *>(p, "residual");
*data_out << step;
*data_out << out;
}
if (debug_vectors)
{
AnyData out;
- VECTOR *p = &u;
- out.add<const VECTOR *>(p, "solution");
+ VectorType *p = &u;
+ out.add<const VectorType *>(p, "solution");
p = Du;
- out.add<const VECTOR *>(p, "update");
+ out.add<const VectorType *>(p, "update");
p = res;
- out.add<const VECTOR *>(p, "residual");
+ out.add<const VectorType *>(p, "residual");
*data_out << step;
*data_out << out;
}
*
* @author Guido Kanschat, 2010
*/
- template <class VECTOR>
+ template <typename VectorType>
class Operator : public OperatorBase
{
public:
*
* @author Guido Kanschat, 2010
*/
- template <class VECTOR>
+ template <typename VectorType>
class OutputOperator : public Subscriptor
{
- OutputOperator(const OutputOperator<VECTOR> &);
+ OutputOperator(const OutputOperator<VectorType> &);
public:
OutputOperator ();
/**
/**
* Output all the vectors in AnyData.
*/
- virtual OutputOperator<VECTOR> &operator<< (const AnyData &vectors);
+ virtual OutputOperator<VectorType> &operator<< (const AnyData &vectors);
protected:
unsigned int step;
std::ostream *os;
};
- template <class VECTOR>
+ template <typename VectorType>
inline
void
- OutputOperator<VECTOR>::set_step (const unsigned int s)
+ OutputOperator<VectorType>::set_step (const unsigned int s)
{
step = s;
}
*
* @relates OutputOperator
*/
- template <class VECTOR>
+ template <typename VectorType>
inline
- OutputOperator<VECTOR> &
- operator<< (OutputOperator<VECTOR> &out, unsigned int step)
+ OutputOperator<VectorType> &
+ operator<< (OutputOperator<VectorType> &out, unsigned int step)
{
out.set_step(step);
return out;
namespace Algorithms
{
- template <class VECTOR>
- Operator<VECTOR>::Operator()
+ template <typename VectorType>
+ Operator<VectorType>::Operator()
{}
- template <class VECTOR>
- OutputOperator<VECTOR>::~OutputOperator()
+ template <typename VectorType>
+ OutputOperator<VectorType>::~OutputOperator()
{}
- template <class VECTOR>
- OutputOperator<VECTOR>::OutputOperator()
+ template <typename VectorType>
+ OutputOperator<VectorType>::OutputOperator()
:
os(0)
{}
- template <class VECTOR>
- void OutputOperator<VECTOR>::initialize_stream(std::ostream &stream)
+ template <typename VectorType>
+ void OutputOperator<VectorType>::initialize_stream(std::ostream &stream)
{
os =&stream;
}
- template <class VECTOR>
- OutputOperator<VECTOR> &
- OutputOperator<VECTOR>::operator<< (const AnyData &vectors)
+ template <typename VectorType>
+ OutputOperator<VectorType> &
+ OutputOperator<VectorType>::operator<< (const AnyData &vectors)
{
if (os == 0)
{
deallog << "Step " << step << std::endl;
for (unsigned int i=0; i<vectors.size(); ++i)
{
- const VECTOR *v = vectors.try_read_ptr<VECTOR>(i);
+ const VectorType *v = vectors.try_read_ptr<VectorType>(i);
if (v == 0) continue;
deallog << vectors.name(i);
for (unsigned int j=0; j<v->size(); ++j)
(*os) << ' ' << step;
for (unsigned int i=0; i<vectors.size(); ++i)
{
- const VECTOR *v = vectors.try_read_ptr<VECTOR>(i);
+ const VectorType *v = vectors.try_read_ptr<VectorType>(i);
if (v == 0) continue;
for (unsigned int j=0; j<v->size(); ++j)
(*os) << ' ' << (*v)(j);
* @author Guido Kanschat
* @date 2010
*/
- template <class VECTOR>
- class ThetaTimestepping : public Operator<VECTOR>
+ template <typename VectorType>
+ class ThetaTimestepping : public Operator<VectorType>
{
public:
/**
* #op_implicit. For their meaning, see the description of those
* variables.
*/
- ThetaTimestepping (Operator<VECTOR> &op_explicit,
- Operator<VECTOR> &op_implicit);
+ ThetaTimestepping (Operator<VectorType> &op_explicit,
+ Operator<VectorType> &op_implicit);
/**
* The timestepping scheme.
* AnyData objects used as input for the operators #op_explicit and
* #op_implicit.
*
- * @param out in its first argument must contain a pointer to a `VECTOR`,
- * which contains the initial value when the operator is called. It
- * contains the final value when the operator returns.
+ * @param out in its first argument must contain a pointer to a VectorType
+ * instance, which contains the initial value when the operator is
+ * called. It contains the final value when the operator returns.
*/
virtual void operator() (AnyData &out, const AnyData &in);
* Define an operator which will output the result in each step. Note that
* no output will be generated without this.
*/
- void set_output(OutputOperator<VECTOR> &output);
+ void set_output(OutputOperator<VectorType> &output);
/**
* Declare parameters in a parameter handler.
* vector, $M$ the mass matrix, $F$ the operator in space and $c$ is the
* adjusted time step size $(1-\theta) \Delta t$.
*/
- SmartPointer<Operator<VECTOR>, ThetaTimestepping<VECTOR> > op_explicit;
+ SmartPointer<Operator<VectorType>, ThetaTimestepping<VectorType> > op_explicit;
/**
* The operator solving the implicit part of the scheme. It will receive
* the input data, <i>M</i> the mass matrix, <i>F</i> the operator in
* space and <i>c</i> is the adjusted time step size $ \theta \Delta t$
*/
- SmartPointer<Operator<VECTOR>, ThetaTimestepping<VECTOR> > op_implicit;
+ SmartPointer<Operator<VectorType>, ThetaTimestepping<VectorType> > op_implicit;
/**
* The operator writing the output in each time step
*/
- SmartPointer<OutputOperator<VECTOR>, ThetaTimestepping<VECTOR> > output;
+ SmartPointer<OutputOperator<VectorType>, ThetaTimestepping<VectorType> > output;
};
- template <class VECTOR>
+ template <typename VectorType>
inline
const TimestepData &
- ThetaTimestepping<VECTOR>::explicit_data () const
+ ThetaTimestepping<VectorType>::explicit_data () const
{
return d_explicit;
}
- template <class VECTOR>
+ template <typename VectorType>
inline
const TimestepData &
- ThetaTimestepping<VECTOR>::implicit_data () const
+ ThetaTimestepping<VectorType>::implicit_data () const
{
return d_implicit;
}
- template <class VECTOR>
+ template <typename VectorType>
inline
TimestepControl &
- ThetaTimestepping<VECTOR>::timestep_control ()
+ ThetaTimestepping<VectorType>::timestep_control ()
{
return control;
}
- template <class VECTOR>
+ template <typename VectorType>
inline
- void ThetaTimestepping<VECTOR>::set_output (OutputOperator<VECTOR> &out)
+ void ThetaTimestepping<VectorType>::set_output (OutputOperator<VectorType> &out)
{
output = &out;
}
- template <class VECTOR>
+ template <typename VectorType>
inline
- double ThetaTimestepping<VECTOR>::theta () const
+ double ThetaTimestepping<VectorType>::theta () const
{
return vtheta;
}
- template <class VECTOR>
+ template <typename VectorType>
inline
- double ThetaTimestepping<VECTOR>::theta (double new_theta)
+ double ThetaTimestepping<VectorType>::theta (double new_theta)
{
const double tmp = vtheta;
vtheta = new_theta;
}
- template <class VECTOR>
+ template <typename VectorType>
inline
- double ThetaTimestepping<VECTOR>::current_time () const
+ double ThetaTimestepping<VectorType>::current_time () const
{
return control.now();
}
namespace Algorithms
{
- template <class VECTOR>
- ThetaTimestepping<VECTOR>::ThetaTimestepping (Operator<VECTOR> &e, Operator<VECTOR> &i)
+ template <typename VectorType>
+ ThetaTimestepping<VectorType>::ThetaTimestepping (Operator<VectorType> &e, Operator<VectorType> &i)
: vtheta(0.5), adaptive(false), op_explicit(&e), op_implicit(&i)
{}
- template <class VECTOR>
+ template <typename VectorType>
void
- ThetaTimestepping<VECTOR>::notify(const Event &e)
+ ThetaTimestepping<VectorType>::notify(const Event &e)
{
op_explicit->notify(e);
op_implicit->notify(e);
}
- template <class VECTOR>
+ template <typename VectorType>
void
- ThetaTimestepping<VECTOR>::declare_parameters(ParameterHandler ¶m)
+ ThetaTimestepping<VectorType>::declare_parameters(ParameterHandler ¶m)
{
param.enter_subsection("ThetaTimestepping");
TimestepControl::declare_parameters (param);
param.leave_subsection();
}
- template <class VECTOR>
+ template <typename VectorType>
void
- ThetaTimestepping<VECTOR>::parse_parameters (ParameterHandler ¶m)
+ ThetaTimestepping<VectorType>::parse_parameters (ParameterHandler ¶m)
{
param.enter_subsection("ThetaTimestepping");
control.parse_parameters (param);
}
- template <class VECTOR>
+ template <typename VectorType>
void
- ThetaTimestepping<VECTOR>::operator() (AnyData &out, const AnyData &in)
+ ThetaTimestepping<VectorType>::operator() (AnyData &out, const AnyData &in)
{
Assert(!adaptive, ExcNotImplemented());
deallog.push ("Theta");
- VECTOR &solution = *out.entry<VECTOR *>(0);
- GrowingVectorMemory<VECTOR> mem;
- typename VectorMemory<VECTOR>::Pointer aux(mem);
+ VectorType &solution = *out.entry<VectorType *>(0);
+ GrowingVectorMemory<VectorType> mem;
+ typename VectorMemory<VectorType>::Pointer aux(mem);
aux->reinit(solution);
control.restart();
// vector associated with the old
// timestep
AnyData src1;
- src1.add<const VECTOR *>(&solution, "Previous iterate");
+ src1.add<const VectorType *>(&solution, "Previous iterate");
src1.add<const double *>(&d_explicit.time, "Time");
src1.add<const double *>(&d_explicit.step, "Timestep");
src1.add<const double *>(&vtheta, "Theta");
AnyData src2;
AnyData out1;
- out1.add<VECTOR *>(aux, "Solution");
+ out1.add<VectorType *>(aux, "Solution");
// The data provided to the inner solver
- src2.add<const VECTOR *>(aux, "Previous time");
- src2.add<const VECTOR *>(&solution, "Previous iterate");
+ src2.add<const VectorType *>(aux, "Previous time");
+ src2.add<const VectorType *>(&solution, "Previous iterate");
src2.add<const double *>(&d_implicit.time, "Time");
src2.add<const double *>(&d_implicit.step, "Timestep");
src2.add<const double *>(&vtheta, "Theta");
* Abstract class for time stepping methods. These methods assume that the
* equation has the form: $ \frac{\partial y}{\partial t} = f(t,y) $.
*/
- template <typename VECTOR>
+ template <typename VectorType>
class TimeStepping
{
public:
* and a vector. The output is the value of function at this point. This
* function returns the time at the end of the time step.
*/
- virtual double evolve_one_time_step(
- std::vector<std_cxx11::function<VECTOR (const double, const VECTOR &)> > &F,
- std::vector<std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> > &J_inverse,
- double t,
- double delta_t,
- VECTOR &y) = 0;
+ virtual double evolve_one_time_step
+ (std::vector<std_cxx11::function<VectorType (const double, const VectorType &)> > &F,
+ std::vector<std_cxx11::function<VectorType (const double, const double, const VectorType &)> > &J_inverse,
+ double t,
+ double delta_t,
+ VectorType &y) = 0;
/**
* Empty structure used to store information.
* @author Damien Lebrun-Grandie, Bruno Turcksin
* @date 2014
*/
- template <typename VECTOR>
- class RungeKutta : public TimeStepping<VECTOR>
+ template <typename VectorType>
+ class RungeKutta : public TimeStepping<VectorType>
{
public:
/**
* returns the time at the end of the time step. When using Runge-Kutta
* methods, @p F and @ J_inverse can only contain one element.
*/
- double evolve_one_time_step(
- std::vector<std_cxx11::function<VECTOR (const double, const VECTOR &)> > &F,
- std::vector<std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> > &J_inverse,
- double t,
- double delta_t,
- VECTOR &y);
+ double evolve_one_time_step
+ (std::vector<std_cxx11::function<VectorType (const double, const VectorType &)> > &F,
+ std::vector<std_cxx11::function<VectorType (const double, const double, const VectorType &)> > &J_inverse,
+ double t,
+ double delta_t,
+ VectorType &y);
/**
* Purely virtual function. This function is used to advance from time @p
* vector. The output is the value of function at this point.
* evolve_one_time_step returns the time at the end of the time step.
*/
- virtual double evolve_one_time_step(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> id_minus_tau_J_inverse,
- double t,
- double delta_t,
- VECTOR &y) = 0;
+ virtual double evolve_one_time_step
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ std_cxx11::function<VectorType (const double, const double, const VectorType &)> id_minus_tau_J_inverse,
+ double t,
+ double delta_t,
+ VectorType &y) = 0;
protected:
/**
* ExplicitRungeKutta is derived from RungeKutta and implement the explicit
* methods.
*/
- template <typename VECTOR>
- class ExplicitRungeKutta : public RungeKutta<VECTOR>
+ template <typename VectorType>
+ class ExplicitRungeKutta : public RungeKutta<VectorType>
{
public:
- using RungeKutta<VECTOR>::evolve_one_time_step;
+ using RungeKutta<VectorType>::evolve_one_time_step;
/**
* Default constructor. initialize(runge_kutta_method) needs to be called
* of function at this point. evolve_one_time_step returns the time at the
* end of the time step.
*/
- double evolve_one_time_step(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> id_minus_tau_J_inverse,
- double t,
- double delta_t,
- VECTOR &y);
+ double evolve_one_time_step
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ std_cxx11::function<VectorType (const double, const double, const VectorType &)> id_minus_tau_J_inverse,
+ double t,
+ double delta_t,
+ VectorType &y);
/**
* This function is used to advance from time @p t to t+ @p delta_t. This
* methods. evolve_one_time_step returns the time at the end of the time
* step.
*/
- double evolve_one_time_step(std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- double t,
- double delta_t,
- VECTOR &y);
+ double evolve_one_time_step
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ double t,
+ double delta_t,
+ VectorType &y);
/**
* This structure stores the name of the method used.
*/
- struct Status : public TimeStepping<VECTOR>::Status
+ struct Status : public TimeStepping<VectorType>::Status
{
runge_kutta_method method;
};
/**
* Compute the different stages needed.
*/
- void compute_stages(std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- const double t,
- const double delta_t,
- const VECTOR &y,
- std::vector<VECTOR> &f_stages) const;
+ void compute_stages
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ const double t,
+ const double delta_t,
+ const VectorType &y,
+ std::vector<VectorType> &f_stages) const;
/**
* Status structure of the object.
* This class is derived from RungeKutta and implement the implicit methods.
* This class works only for Diagonal Implicit Runge-Kutta (DIRK) methods.
*/
- template <typename VECTOR>
- class ImplicitRungeKutta : public RungeKutta<VECTOR>
+ template <typename VectorType>
+ class ImplicitRungeKutta : public RungeKutta<VectorType>
{
public:
- using RungeKutta<VECTOR>::evolve_one_time_step;
+ using RungeKutta<VectorType>::evolve_one_time_step;
/**
* Default constructor. initialize(runge_kutta_method) and
* The output is the value of function at this point. evolve_one_time_step
* returns the time at the end of the time step.
*/
- double evolve_one_time_step(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> id_minus_tau_J_inverse,
- double t,
- double delta_t,
- VECTOR &y);
+ double evolve_one_time_step
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ std_cxx11::function<VectorType (const double, const double, const VectorType &)> id_minus_tau_J_inverse,
+ double t,
+ double delta_t,
+ VectorType &y);
/**
* Set the maximum number of iterations and the tolerance used by the
* Structure that stores the name of the method, the number of Newton
* iterations and the norm of the residual when exiting the Newton solver.
*/
- struct Status : public TimeStepping<VECTOR>::Status
+ struct Status : public TimeStepping<VectorType>::Status
{
runge_kutta_method method;
unsigned int n_iterations;
/**
* Compute the different stages needed.
*/
- void compute_stages(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> id_minus_tau_J_inverse,
- double t,
- double delta_t,
- VECTOR &y,
- std::vector<VECTOR> &f_stages);
+ void compute_stages
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ std_cxx11::function<VectorType (const double, const double, const VectorType &)> id_minus_tau_J_inverse,
+ double t,
+ double delta_t,
+ VectorType &y,
+ std::vector<VectorType> &f_stages);
/**
* Newton solver used for the implicit stages.
*/
- void newton_solve(std_cxx11::function<void (const VECTOR &,VECTOR &)> get_residual,
- std_cxx11::function<VECTOR (const VECTOR &)> id_minus_tau_J_inverse,
- VECTOR &y);
+ void newton_solve(std_cxx11::function<void (const VectorType &,VectorType &)> get_residual,
+ std_cxx11::function<VectorType (const VectorType &)> id_minus_tau_J_inverse,
+ VectorType &y);
/**
* Compute the residual needed by the Newton solver.
*/
- void compute_residual(std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- double t,
- double delta_t,
- const VECTOR &old_y,
- const VECTOR &y,
- VECTOR &tendency,
- VECTOR &residual) const;
+ void compute_residual(std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ double t,
+ double delta_t,
+ const VectorType &old_y,
+ const VectorType &y,
+ VectorType &tendency,
+ VectorType &residual) const;
/**
* When using SDIRK, there is no need to compute the linear combination of
* This is class is derived from RungeKutta and implement embedded explicit
* methods.
*/
- template <typename VECTOR>
- class EmbeddedExplicitRungeKutta : public RungeKutta<VECTOR>
+ template <typename VectorType>
+ class EmbeddedExplicitRungeKutta : public RungeKutta<VectorType>
{
public:
- using RungeKutta<VECTOR>::evolve_one_time_step;
+ using RungeKutta<VectorType>::evolve_one_time_step;
/**
* Default constructor. initialize(runge_kutta_method) and
* value of function at this point. evolve_one_time_step returns the time
* at the end of the time step.
*/
- double evolve_one_time_step(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> id_minus_tau_J_inverse,
- double t,
- double delta_t,
- VECTOR &y);
+ double evolve_one_time_step
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ std_cxx11::function<VectorType (const double, const double, const VectorType &)> id_minus_tau_J_inverse,
+ double t,
+ double delta_t,
+ VectorType &y);
/**
* This function is used to advance from time @p t to t+ @p delta_t. This
* methods. evolve_one_time_step returns the time at the end of the time
* step.
*/
- double evolve_one_time_step(std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- double t,
- double delta_t,
- VECTOR &y);
+ double evolve_one_time_step
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ double t,
+ double delta_t,
+ VectorType &y);
/**
* Set the parameters necessary for the time adaptation.
* guess of what the next time step should be, and an estimate of the norm
* of the error.
*/
- struct Status : public TimeStepping<VECTOR>::Status
+ struct Status : public TimeStepping<VectorType>::Status
{
runge_kutta_method method;
embedded_runge_kutta_time_step exit_delta_t;
/**
* Compute the different stages needed.
*/
- void compute_stages(std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- const double t,
- const double delta_t,
- const VECTOR &y,
- std::vector<VECTOR> &f_stages);
+ void compute_stages(std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ const double t,
+ const double delta_t,
+ const VectorType &y,
+ std::vector<VectorType> &f_stages);
/**
* This parameter is the factor (>1) by which the time step is multiplied
* If the last_same_as_first flag is set to true, the last stage is saved
* and reused as the first stage of the next time step.
*/
- VECTOR *last_stage;
+ VectorType *last_stage;
/**
* Status structure of the object.
// RungeKutta
// ----------------------------------------------------------------------
- template <typename VECTOR>
- double RungeKutta<VECTOR>::evolve_one_time_step(
- std::vector<std_cxx11::function<VECTOR (const double, const VECTOR &)> > &F,
- std::vector<std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> > &J_inverse,
+ template <typename VectorType>
+ double RungeKutta<VectorType>::evolve_one_time_step(
+ std::vector<std_cxx11::function<VectorType (const double, const VectorType &)> > &F,
+ std::vector<std_cxx11::function<VectorType (const double, const double, const VectorType &)> > &J_inverse,
double t,
double delta_t,
- VECTOR &y)
+ VectorType &y)
{
AssertThrow(F.size()==0,
ExcMessage("RungeKutta methods cannot handle more that one function to integate."));
// ExplicitRungeKutta
// ----------------------------------------------------------------------
- template <typename VECTOR>
- ExplicitRungeKutta<VECTOR>::ExplicitRungeKutta(runge_kutta_method method)
+ template <typename VectorType>
+ ExplicitRungeKutta<VectorType>::ExplicitRungeKutta(runge_kutta_method method)
{
initialize(method);
}
- template <typename VECTOR>
- void ExplicitRungeKutta<VECTOR>::initialize(runge_kutta_method method)
+ template <typename VectorType>
+ void ExplicitRungeKutta<VectorType>::initialize(runge_kutta_method method)
{
status.method = method;
- template <typename VECTOR>
- double ExplicitRungeKutta<VECTOR>::evolve_one_time_step(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> /*id_minus_tau_J_inverse*/,
- double t,
- double delta_t,
- VECTOR &y)
+ template <typename VectorType>
+ double ExplicitRungeKutta<VectorType>::evolve_one_time_step
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ std_cxx11::function<VectorType (const double, const double, const VectorType &)> /*id_minus_tau_J_inverse*/,
+ double t,
+ double delta_t,
+ VectorType &y)
{
return evolve_one_time_step(f,t,delta_t,y);
}
- template <typename VECTOR>
- double ExplicitRungeKutta<VECTOR>::evolve_one_time_step(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- double t,
- double delta_t,
- VECTOR &y)
+ template <typename VectorType>
+ double ExplicitRungeKutta<VectorType>::evolve_one_time_step
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ double t,
+ double delta_t,
+ VectorType &y)
{
- std::vector<VECTOR> f_stages(this->n_stages,y);
+ std::vector<VectorType> f_stages(this->n_stages,y);
// Compute the different stages needed.
compute_stages(f,t,delta_t,y,f_stages);
- template <typename VECTOR>
- const typename ExplicitRungeKutta<VECTOR>::Status &ExplicitRungeKutta<VECTOR>::get_status() const
+ template <typename VectorType>
+ const typename ExplicitRungeKutta<VectorType>::Status &ExplicitRungeKutta<VectorType>::get_status() const
{
return status;
}
- template <typename VECTOR>
- void ExplicitRungeKutta<VECTOR>::compute_stages(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- const double t,
- const double delta_t,
- const VECTOR &y,
- std::vector<VECTOR> &f_stages) const
+ template <typename VectorType>
+ void ExplicitRungeKutta<VectorType>::compute_stages
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ const double t,
+ const double delta_t,
+ const VectorType &y,
+ std::vector<VectorType> &f_stages) const
{
for (unsigned int i=0; i<this->n_stages; ++i)
{
- VECTOR Y(y);
+ VectorType Y(y);
for (unsigned int j=0; j<i; ++j)
Y.sadd(1.,delta_t *this->a[i][j],f_stages[j]);
// Evaluate the function f at the point (t+c[i]*delta_t,Y).
// ImplicitRungeKutta
// ----------------------------------------------------------------------
- template <typename VECTOR>
- ImplicitRungeKutta<VECTOR>::ImplicitRungeKutta(runge_kutta_method method,
+ template <typename VectorType>
+ ImplicitRungeKutta<VectorType>::ImplicitRungeKutta(runge_kutta_method method,
unsigned int max_it,
double tolerance)
:
- RungeKutta<VECTOR> (),
+ RungeKutta<VectorType> (),
skip_linear_combi(false),
max_it(max_it),
tolerance(tolerance)
- template <typename VECTOR>
- void ImplicitRungeKutta<VECTOR>::initialize(runge_kutta_method method)
+ template <typename VectorType>
+ void ImplicitRungeKutta<VectorType>::initialize(runge_kutta_method method)
{
status.method = method;
- template <typename VECTOR>
- double ImplicitRungeKutta<VECTOR>::evolve_one_time_step(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> id_minus_tau_J_inverse,
- double t,
- double delta_t,
- VECTOR &y)
+ template <typename VectorType>
+ double ImplicitRungeKutta<VectorType>::evolve_one_time_step
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ std_cxx11::function<VectorType (const double, const double, const VectorType &)> id_minus_tau_J_inverse,
+ double t,
+ double delta_t,
+ VectorType &y)
{
- VECTOR old_y(y);
- std::vector<VECTOR> f_stages(this->n_stages,y);
+ VectorType old_y(y);
+ std::vector<VectorType> f_stages(this->n_stages,y);
// Compute the different stages needed.
compute_stages(f,id_minus_tau_J_inverse,t,delta_t,y,f_stages);
- template <typename VECTOR>
- void ImplicitRungeKutta<VECTOR>::set_newton_solver_parameters(unsigned int max_it_, double tolerance_)
+ template <typename VectorType>
+ void ImplicitRungeKutta<VectorType>::set_newton_solver_parameters(unsigned int max_it_, double tolerance_)
{
max_it = max_it_;
tolerance = tolerance_;
- template <typename VECTOR>
- const typename ImplicitRungeKutta<VECTOR>::Status &ImplicitRungeKutta<VECTOR>::get_status() const
+ template <typename VectorType>
+ const typename ImplicitRungeKutta<VectorType>::Status &ImplicitRungeKutta<VectorType>::get_status() const
{
return status;
}
- template <typename VECTOR>
- void ImplicitRungeKutta<VECTOR>::compute_stages(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> id_minus_tau_J_inverse,
+ template <typename VectorType>
+ void ImplicitRungeKutta<VectorType>::compute_stages(
+ std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ std_cxx11::function<VectorType (const double, const double, const VectorType &)> id_minus_tau_J_inverse,
double t,
double delta_t,
- VECTOR &y,
- std::vector<VECTOR> &f_stages)
+ VectorType &y,
+ std::vector<VectorType> &f_stages)
{
- VECTOR z(y);
+ VectorType z(y);
for (unsigned int i=0; i<this->n_stages; ++i)
{
- VECTOR old_y(z);
+ VectorType old_y(z);
for (unsigned int j=0; j<i; ++j)
old_y.sadd(1.,delta_t *this->a[i][j],f_stages[j]);
// Solve the nonlinear system using Newton's method
const double new_t = t+this->c[i]*delta_t;
const double new_delta_t = this->a[i][i]*delta_t;
- newton_solve(std_cxx11::bind(&ImplicitRungeKutta<VECTOR>::compute_residual,this,f,new_t,new_delta_t,
+ newton_solve(std_cxx11::bind(&ImplicitRungeKutta<VectorType>::compute_residual,this,f,new_t,new_delta_t,
std_cxx11::cref(old_y),std_cxx11::_1,std_cxx11::ref(f_stages[i]),std_cxx11::_2),
std_cxx11::bind(id_minus_tau_J_inverse,new_t,new_delta_t,std_cxx11::_1),y);
}
- template <typename VECTOR>
- void ImplicitRungeKutta<VECTOR>::newton_solve(
- std_cxx11::function<void (const VECTOR &,VECTOR &)> get_residual,
- std_cxx11::function<VECTOR (const VECTOR &)> id_minus_tau_J_inverse,
- VECTOR &y)
+ template <typename VectorType>
+ void ImplicitRungeKutta<VectorType>::newton_solve(
+ std_cxx11::function<void (const VectorType &,VectorType &)> get_residual,
+ std_cxx11::function<VectorType (const VectorType &)> id_minus_tau_J_inverse,
+ VectorType &y)
{
- VECTOR residual(y);
+ VectorType residual(y);
get_residual(y,residual);
unsigned int i=0;
const double initial_residual_norm = residual.l2_norm();
- template <typename VECTOR>
- void ImplicitRungeKutta<VECTOR>::compute_residual(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- double t,
- double delta_t,
- const VECTOR &old_y,
- const VECTOR &y,
- VECTOR &tendency,
- VECTOR &residual) const
+ template <typename VectorType>
+ void ImplicitRungeKutta<VectorType>::compute_residual
+ (std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ double t,
+ double delta_t,
+ const VectorType &old_y,
+ const VectorType &y,
+ VectorType &tendency,
+ VectorType &residual) const
{
// The tendency is stored to save one evaluation of f.
tendency = f(t,y);
// EmbeddedExplicitRungeKutta
// ----------------------------------------------------------------------
- template <typename VECTOR>
- EmbeddedExplicitRungeKutta<VECTOR>::EmbeddedExplicitRungeKutta(runge_kutta_method method,
- double coarsen_param,
- double refine_param,
- double min_delta,
- double max_delta,
- double refine_tol,
- double coarsen_tol)
+ template <typename VectorType>
+ EmbeddedExplicitRungeKutta<VectorType>::EmbeddedExplicitRungeKutta
+ (runge_kutta_method method,
+ double coarsen_param,
+ double refine_param,
+ double min_delta,
+ double max_delta,
+ double refine_tol,
+ double coarsen_tol)
:
coarsen_param(coarsen_param),
refine_param(refine_param),
- template <typename VECTOR>
- void EmbeddedExplicitRungeKutta<VECTOR>::initialize(runge_kutta_method method)
+ template <typename VectorType>
+ void EmbeddedExplicitRungeKutta<VectorType>::initialize(runge_kutta_method method)
{
status.method = method;
- template <typename VECTOR>
- void EmbeddedExplicitRungeKutta<VECTOR>::free_memory()
+ template <typename VectorType>
+ void EmbeddedExplicitRungeKutta<VectorType>::free_memory()
{
if (last_stage!=NULL)
delete last_stage;
- template <typename VECTOR>
- double EmbeddedExplicitRungeKutta<VECTOR>::evolve_one_time_step(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- std_cxx11::function<VECTOR (const double, const double, const VECTOR &)> /*id_minus_tau_J_inverse*/,
+ template <typename VectorType>
+ double EmbeddedExplicitRungeKutta<VectorType>::evolve_one_time_step(
+ std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ std_cxx11::function<VectorType (const double, const double, const VectorType &)> /*id_minus_tau_J_inverse*/,
double t,
double delta_t,
- VECTOR &y)
+ VectorType &y)
{
return evolve_one_time_step(f,t,delta_t,y);
}
- template <typename VECTOR>
- double EmbeddedExplicitRungeKutta<VECTOR>::evolve_one_time_step(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
- double t, double delta_t, VECTOR &y)
+ template <typename VectorType>
+ double EmbeddedExplicitRungeKutta<VectorType>::evolve_one_time_step(
+ std_cxx11::function<VectorType (const double, const VectorType &)> f,
+ double t, double delta_t, VectorType &y)
{
bool done = false;
unsigned int count = 0;
double error_norm = 0.;
- VECTOR old_y(y);
- VECTOR error(y);
- std::vector<VECTOR> f_stages(this->n_stages,y);
+ VectorType old_y(y);
+ VectorType error(y);
+ std::vector<VectorType> f_stages(this->n_stages,y);
while (!done)
{
if (last_same_as_first==true)
{
if (last_stage==NULL)
- last_stage = new VECTOR(f_stages.back());
+ last_stage = new VectorType(f_stages.back());
else
*last_stage = f_stages.back();
}
- template <typename VECTOR>
- void EmbeddedExplicitRungeKutta<VECTOR>::set_time_adaptation_parameters(double coarsen_param_,
+ template <typename VectorType>
+ void EmbeddedExplicitRungeKutta<VectorType>::set_time_adaptation_parameters(double coarsen_param_,
double refine_param_,
double min_delta_,
double max_delta_,
- template <typename VECTOR>
- const typename EmbeddedExplicitRungeKutta<VECTOR>::Status &EmbeddedExplicitRungeKutta<VECTOR>::get_status() const
+ template <typename VectorType>
+ const typename EmbeddedExplicitRungeKutta<VectorType>::Status &EmbeddedExplicitRungeKutta<VectorType>::get_status() const
{
return status;
}
- template <typename VECTOR>
- void EmbeddedExplicitRungeKutta<VECTOR>::compute_stages(
- std_cxx11::function<VECTOR (const double, const VECTOR &)> f,
+ template <typename VectorType>
+ void EmbeddedExplicitRungeKutta<VectorType>::compute_stages(
+ std_cxx11::function<VectorType (const double, const VectorType &)> f,
const double t,
const double delta_t,
- const VECTOR &y,
- std::vector<VECTOR> &f_stages)
+ const VectorType &y,
+ std::vector<VectorType> &f_stages)
{
- VECTOR Y(y);
+ VectorType Y(y);
unsigned int i = 0;
// If the last stage is the same as the first, we can skip the evaluation
* (unsigned int)</tt> and a function <tt>size() const</tt>.
*
* The use of this object is straightforward. It duplicates the random access
- * operator of the <tt>VECTOR</tt> and adds an offset to every index.
+ * operator of the <tt>VectorType</tt> and adds an offset to every index.
*
* Some precautions have to be taken if it is used for a constant vector: the
* VectorSlice object has to be constant, too. The appropriate initialization
* @ingroup data
* @author Guido Kanschat, 2004
*/
-template <class VECTOR>
+template <typename VectorType>
class VectorSlice
{
public:
* just put in the vector itself as argument and let this constructor make a
* slice for you.
*/
- VectorSlice(VECTOR &v);
+ VectorSlice(VectorType &v);
/**
* The real constructor for a vector slice, allowing you to specify the
* start index and the length of the slice.
*/
- VectorSlice(VECTOR &v,
+ VectorSlice(VectorType &v,
unsigned int start,
unsigned int length);
* Access an element of the slice using the same interface as
* <tt>std::vector</tt>.
*/
- typename VECTOR::reference operator[] (unsigned int i);
+ typename VectorType::reference operator[] (unsigned int i);
/**
* Access an element of a constant slice using the same interface as
* <tt>std::vector</tt>.
*/
- typename VECTOR::const_reference operator[] (unsigned int i) const;
+ typename VectorType::const_reference operator[] (unsigned int i) const;
/**
* Standard-conforming iterator function.
*/
- typename VECTOR::iterator begin();
+ typename VectorType::iterator begin();
/**
* Standard-conforming iterator function.
*/
- typename VECTOR::const_iterator begin() const;
+ typename VectorType::const_iterator begin() const;
/**
* Standard-conforming iterator function.
*/
- typename VECTOR::iterator end();
+ typename VectorType::iterator end();
/**
* Standard-conforming iterator function.
*/
- typename VECTOR::const_iterator end() const;
+ typename VectorType::const_iterator end() const;
private:
/**
* The vector we extract from.
*/
- VECTOR &v;
+ VectorType &v;
/**
* The start index of the slice.
*/
* @relates VectorSlice
* @author Guido Kanschat, 2004
*/
-template <class VECTOR>
+template <typename VectorType>
inline
-const VectorSlice<const VECTOR>
-make_slice (VECTOR &v)
+const VectorSlice<const VectorType>
+make_slice (VectorType &v)
{
- const VectorSlice<const VECTOR> r(v);
+ const VectorSlice<const VectorType> r(v);
return r;
}
* @relates VectorSlice
* @author Guido Kanschat, 2004
*/
-template <class VECTOR>
+template <typename VectorType>
inline
-const VectorSlice<const VECTOR>
-make_slice (VECTOR &v,
+const VectorSlice<const VectorType>
+make_slice (VectorType &v,
const unsigned int start,
const unsigned int length)
{
- const VectorSlice<const VECTOR> r(v, start, length);
+ const VectorSlice<const VectorType> r(v, start, length);
return r;
}
//---------------------------------------------------------------------------
-template <class VECTOR>
+template <typename VectorType>
inline
-VectorSlice<VECTOR>::VectorSlice(VECTOR &v)
+VectorSlice<VectorType>::VectorSlice(VectorType &v)
:
v(v), start(0), length(v.size())
{}
-template <class VECTOR>
+template <typename VectorType>
inline
-VectorSlice<VECTOR>::VectorSlice(VECTOR &v,
- unsigned int start,
- unsigned int length)
+VectorSlice<VectorType>::VectorSlice(VectorType &v,
+ unsigned int start,
+ unsigned int length)
:
v(v), start(start), length(length)
{
}
-template <class VECTOR>
+template <typename VectorType>
inline
unsigned int
-VectorSlice<VECTOR>::size() const
+VectorSlice<VectorType>::size() const
{
return length;
}
-template <class VECTOR>
+template <typename VectorType>
inline
-typename VECTOR::reference
-VectorSlice<VECTOR>::operator[](unsigned int i)
+typename VectorType::reference
+VectorSlice<VectorType>::operator[](unsigned int i)
{
Assert ((i<length), ExcIndexRange(i, 0, length));
}
-template <class VECTOR>
+template <typename VectorType>
inline
-typename VECTOR::const_reference
-VectorSlice<VECTOR>::operator[](unsigned int i) const
+typename VectorType::const_reference
+VectorSlice<VectorType>::operator[](unsigned int i) const
{
Assert ((i<length), ExcIndexRange(i, 0, length));
}
-template <class VECTOR>
+template <typename VectorType>
inline
-typename VECTOR::const_iterator
-VectorSlice<VECTOR>::begin() const
+typename VectorType::const_iterator
+VectorSlice<VectorType>::begin() const
{
return v.begin()+start;
}
-template <class VECTOR>
+template <typename VectorType>
inline
-typename VECTOR::iterator
-VectorSlice<VECTOR>::begin()
+typename VectorType::iterator
+VectorSlice<VectorType>::begin()
{
return v.begin()+start;
}
-template <class VECTOR>
+template <typename VectorType>
inline
-typename VECTOR::const_iterator
-VectorSlice<VECTOR>::end() const
+typename VectorType::const_iterator
+VectorSlice<VectorType>::end() const
{
return v.begin()+start+length;
}
-template <class VECTOR>
+template <typename VectorType>
inline
-typename VECTOR::iterator
-VectorSlice<VECTOR>::end()
+typename VectorType::iterator
+VectorSlice<VectorType>::end()
{
return v.begin()+start+length;
}
* interpolate() or deserialize() you need to supply distributed vectors
* without ghost elements.
*
- * <h3>Transferring a solution</h3> Here VECTOR is your favorite vector
+ * <h3>Transferring a solution</h3> Here VectorType is your favorite vector
* type, e.g. PETScWrappers::MPI::Vector, TrilinosWrappers::MPI::Vector,
* or corresponding blockvectors.
* @code
- * SolutionTransfer<dim, VECTOR> soltrans(dof_handler);
+ * SolutionTransfer<dim, VectorType> soltrans(dof_handler);
* // flag some cells for refinement
* // and coarsening, e.g.
* GridRefinement::refine_and_coarsen_fixed_fraction(
* // redistribute dofs,
* dof_handler.distribute_dofs (fe);
* // and interpolate the solution
- * VECTOR interpolated_solution;
- * //create VECTOR in the right size here
+ * VectorType interpolated_solution;
+ * //create VectorType in the right size here
* soltrans.interpolate(interpolated_solution);
* @endcode
*
* follows:
* *@code
*
- * parallel::distributed::SolutionTransfer<dim,VECTOR> sol_trans(dof_handler);
+ * parallel::distributed::SolutionTransfer<dim,VectorType> sol_trans(dof_handler);
* sol_trans.prepare_serialization (vector);
*
* triangulation.save(filename);
* //[create coarse mesh...]
* triangulation.load(filename);
*
- * parallel::distributed::SolutionTransfer<dim,VECTOR> sol_trans(dof_handler);
+ * parallel::distributed::SolutionTransfer<dim,VectorType> sol_trans(dof_handler);
* sol_trans.deserialize (distributed_vector);
* @endcode
*
* @ingroup distributed
* @author Timo Heister, 2009-2011
*/
- template<int dim, typename VECTOR, class DH=DoFHandler<dim> >
+ template<int dim, typename VectorType, class DH=DoFHandler<dim> >
class SolutionTransfer
{
public:
* includes all vectors that are to be interpolated onto the new
* (refined and/or coarsened) grid.
*/
- void prepare_for_coarsening_and_refinement (const std::vector<const VECTOR *> &all_in);
+ void prepare_for_coarsening_and_refinement (const std::vector<const VectorType *> &all_in);
/**
* Same as previous function but for only one discrete function to be
* interpolated.
*/
- void prepare_for_coarsening_and_refinement (const VECTOR &in);
+ void prepare_for_coarsening_and_refinement (const VectorType &in);
/**
* Interpolate the data previously stored in this object before the mesh
* prepare_for_coarsening_and_refinement() and write the result into the
* given set of vectors.
*/
- void interpolate (std::vector<VECTOR *> &all_out);
+ void interpolate (std::vector<VectorType *> &all_out);
/**
* Same as the previous function. It interpolates only one function. It
* several functions can be performed in one step by using
* <tt>interpolate (all_in, all_out)</tt>
*/
- void interpolate (VECTOR &out);
+ void interpolate (VectorType &out);
/**
* on the locally active DoFs (it must be ghosted). See documentation of
* this class for more information.
*/
- void prepare_serialization(const VECTOR &in);
+ void prepare_serialization(const VectorType &in);
/**
* Same as the function above, only for a list of vectors.
*/
- void prepare_serialization(const std::vector<const VECTOR *> &all_in);
+ void prepare_serialization(const std::vector<const VectorType *> &all_in);
/**
* fully distributed vector without ghost elements. See documentation of
* this class for more information.
*/
- void deserialize(VECTOR &in);
+ void deserialize(VectorType &in);
/**
* Same as the function above, only for a list of vectors.
*/
- void deserialize(std::vector<VECTOR *> &all_in);
+ void deserialize(std::vector<VectorType *> &all_in);
private:
/**
* Pointer to the degree of freedom handler to work with.
*/
- SmartPointer<const DH,SolutionTransfer<dim,VECTOR,DH> > dof_handler;
+ SmartPointer<const DH,SolutionTransfer<dim,VectorType,DH> > dof_handler;
/**
* A vector that stores pointers to all the vectors we are supposed to
* copy over from the old to the new mesh.
*/
- std::vector<const VECTOR *> input_vectors;
+ std::vector<const VectorType *> input_vectors;
/**
* The offset that the Triangulation has assigned to this object
void unpack_callback(const typename Triangulation<dim,dim>::cell_iterator &cell,
const typename Triangulation<dim,dim>::CellStatus status,
const void *data,
- std::vector<VECTOR *> &all_out);
+ std::vector<VectorType *> &all_out);
/**
* @author Luca Heltai, Marco Tezzele 2013, 2015
*/
template <int dim, int spacedim=dim,
- class VECTOR=Vector<double>,
+ typename VectorType=Vector<double>,
class DH=DoFHandler<dim,spacedim> >
class MappingFEField : public Mapping<dim,spacedim>
{
public:
/**
- * Constructor. The first argument is a VECTOR that specifies the
+ * Constructor. The first argument is a VectorType that specifies the
* transformation of the domain from the reference to the current
* configuration.
*
*
* If an incompatible mask is passed, an exception is thrown.
*/
- MappingFEField (const DH &euler_dof_handler,
- const VECTOR &euler_vector,
+ MappingFEField (const DH &euler_dof_handler,
+ const VectorType &euler_vector,
const ComponentMask mask=ComponentMask());
/**
* Copy constructor.
*/
- MappingFEField (const MappingFEField<dim,spacedim,VECTOR,DH> &mapping);
+ MappingFEField (const MappingFEField<dim,spacedim,VectorType,DH> &mapping);
/**
* Return a pointer to a copy of the present object. The caller of this copy
/**
* Reference to the vector of shifts.
*/
- SmartPointer<const VECTOR, MappingFEField<dim,spacedim,VECTOR,DH> > euler_vector;
+ SmartPointer<const VectorType, MappingFEField<dim,spacedim,VectorType,DH> > euler_vector;
/**
* A FiniteElement object which is only needed in 3D, since it knows how to
* prevent construction in 1D and 2D, but since memory and time requirements
* are not particularly high this seems unnecessary at the moment.
*/
- SmartPointer<const FiniteElement<dim,spacedim>, MappingFEField<dim,spacedim,VECTOR,DH> > fe;
+ SmartPointer<const FiniteElement<dim,spacedim>, MappingFEField<dim,spacedim,VectorType,DH> > fe;
/**
* Pointer to the DoFHandler to which the mapping vector is associated.
*/
- SmartPointer<const DH,MappingFEField<dim,spacedim,VECTOR,DH> > euler_dof_handler;
+ SmartPointer<const DH,MappingFEField<dim,spacedim,VectorType,DH> > euler_dof_handler;
private:
/**
*
* @author Michael Stadler, 2001
*/
-template <int dim, class VECTOR = Vector<double>, int spacedim=dim >
+template <int dim, typename VectorType = Vector<double>, int spacedim=dim >
class MappingQ1Eulerian : public MappingQGeneric<dim,spacedim>
{
public:
* problem. Alternatively, the @p Vector can be initialized by
* <tt>DoFAccessor::set_dof_values()</tt>.
*/
- MappingQ1Eulerian (const VECTOR &euler_transform_vectors,
+ MappingQ1Eulerian (const VectorType &euler_transform_vectors,
const DoFHandler<dim,spacedim> &shiftmap_dof_handler);
/**
* then assumes ownership of it.
*/
virtual
- MappingQ1Eulerian<dim,VECTOR,spacedim> *clone () const;
+ MappingQ1Eulerian<dim,VectorType,spacedim> *clone () const;
/**
* Always returns @p false because MappingQ1Eulerian does not in general
/**
* Reference to the vector of shifts.
*/
- SmartPointer<const VECTOR, MappingQ1Eulerian<dim,VECTOR,spacedim> > euler_transform_vectors;
+ SmartPointer<const VectorType, MappingQ1Eulerian<dim,VectorType,spacedim> > euler_transform_vectors;
/**
* Pointer to the DoFHandler to which the mapping vector is associated.
*/
- SmartPointer<const DoFHandler<dim,spacedim>,MappingQ1Eulerian<dim,VECTOR,spacedim> > shiftmap_dof_handler;
+ SmartPointer<const DoFHandler<dim,spacedim>,MappingQ1Eulerian<dim,VectorType,spacedim> > shiftmap_dof_handler;
};
/*@}*/
#ifndef DOXYGEN
-template <int dim, class VECTOR, int spacedim>
+template <int dim, typename VectorType, int spacedim>
inline
bool
-MappingQ1Eulerian<dim,VECTOR,spacedim>::preserves_vertex_locations () const
+MappingQ1Eulerian<dim,VectorType,spacedim>::preserves_vertex_locations () const
{
return false;
}
*
* @author Joshua White, 2008
*/
-template <int dim, class VECTOR = Vector<double>, int spacedim=dim >
+template <int dim, typename VectorType = Vector<double>, int spacedim=dim >
class MappingQEulerian : public MappingQ<dim, spacedim>
{
public:
* interpreted as the displacement we use in defining the mapping,
* relative to the location of cells of the underlying triangulation.
*/
- MappingQEulerian (const unsigned int degree,
+ MappingQEulerian (const unsigned int degree,
const DoFHandler<dim,spacedim> &euler_dof_handler,
- const VECTOR &euler_vector);
+ const VectorType &euler_vector);
/**
* @deprecated Use the constructor with the reverse order of second and
* third argument.
*/
- MappingQEulerian (const unsigned int degree,
- const VECTOR &euler_vector,
+ MappingQEulerian (const unsigned int degree,
+ const VectorType &euler_vector,
const DoFHandler<dim,spacedim> &euler_dof_handler) DEAL_II_DEPRECATED;
/**
/**
* Reference to the vector of shifts.
*/
- SmartPointer<const VECTOR, MappingQEulerian<dim,VECTOR,spacedim> > euler_vector;
+ SmartPointer<const VectorType, MappingQEulerian<dim,VectorType,spacedim> > euler_vector;
/**
* Pointer to the DoFHandler to which the mapping vector is associated.
*/
- SmartPointer<const DoFHandler<dim,spacedim>,MappingQEulerian<dim,VECTOR,spacedim> > euler_dof_handler;
+ SmartPointer<const DoFHandler<dim,spacedim>,MappingQEulerian<dim,VectorType,spacedim> > euler_dof_handler;
private:
#ifndef DOXYGEN
-template <int dim, class VECTOR, int spacedim>
+template <int dim, typename VectorType, int spacedim>
inline
bool
-MappingQEulerian<dim,VECTOR,spacedim>::preserves_vertex_locations () const
+MappingQEulerian<dim,VectorType,spacedim>::preserves_vertex_locations () const
{
return false;
}
#endif // dealii__mapping_q_eulerian_h
-
* Leave it at its default zero, which will be reset to the size of
* <code>eigenvalues</code> internally.
*/
- template <typename VECTOR, typename MATRIX1,
+ template <typename VectorType, typename MATRIX1,
typename MATRIX2, typename INVERSE>
- void solve(
- const MATRIX1 &A,
- const MATRIX2 &B,
- const INVERSE &inverse,
- std::vector<std::complex<double> > &eigenvalues,
- std::vector<VECTOR> &eigenvectors,
- const unsigned int n_eigenvalues = 0);
+ void solve (const MATRIX1 &A,
+ const MATRIX2 &B,
+ const INVERSE &inverse,
+ std::vector<std::complex<double> > &eigenvalues,
+ std::vector<VectorType> &eigenvectors,
+ const unsigned int n_eigenvalues = 0);
protected:
{}
-template <typename VECTOR, typename MATRIX1,
+template <typename VectorType, typename MATRIX1,
typename MATRIX2, typename INVERSE>
inline
-void ArpackSolver::solve (
- const MATRIX1 &system_matrix,
- const MATRIX2 &mass_matrix,
- const INVERSE &inverse,
- std::vector<std::complex<double> > &eigenvalues,
- std::vector<VECTOR> &eigenvectors,
- const unsigned int n_eigenvalues)
+void ArpackSolver::solve (const MATRIX1 &system_matrix,
+ const MATRIX2 &mass_matrix,
+ const INVERSE &inverse,
+ std::vector<std::complex<double> > &eigenvalues,
+ std::vector<VectorType> &eigenvectors,
+ const unsigned int n_eigenvalues)
{
//inside the routines of ARPACK the
//values change magically, so store
case -1:
{
- VECTOR src,dst,tmp;
+ VectorType src,dst,tmp;
src.reinit(eigenvectors[0]);
dst.reinit(src);
tmp.reinit(src);
case 1:
{
- VECTOR src,dst,tmp, tmp2;
+ VectorType src,dst,tmp, tmp2;
src.reinit(eigenvectors[0]);
dst.reinit(src);
tmp.reinit(src);
case 2:
{
- VECTOR src,dst;
+ VectorType src,dst;
src.reinit(eigenvectors[0]);
dst.reinit(src);
DEAL_II_NAMESPACE_OPEN
template <typename number> class Vector;
-template <class VECTOR> class FilteredMatrixBlock;
+template <class VectorType> class FilteredMatrixBlock;
/*! @addtogroup Matrix2
*
* @author Wolfgang Bangerth 2001, Luca Heltai 2006, Guido Kanschat 2007, 2008
*/
-template <class VECTOR>
+template <typename VectorType>
class FilteredMatrix : public Subscriptor
{
public:
* Constructor. Since we use accessors only for read access, a const
* matrix pointer is sufficient.
*/
- Accessor (const FilteredMatrix<VECTOR> *matrix,
- const size_type index);
+ Accessor (const FilteredMatrix<VectorType> *matrix,
+ const size_type index);
public:
/**
/**
* The matrix accessed.
*/
- const FilteredMatrix<VECTOR> *matrix;
+ const FilteredMatrix<VectorType> *matrix;
/**
* Current row number.
/**
* Constructor.
*/
- const_iterator(const FilteredMatrix<VECTOR> *matrix,
- const size_type index);
+ const_iterator(const FilteredMatrix<VectorType> *matrix,
+ const size_type index);
/**
* Prefix increment.
* the second parameter to @p true to use a faster algorithm. Note: This
* method is deprecated as matrix_is_symmetric parameter is no longer used.
*/
- void apply_constraints (VECTOR &v,
- const bool matrix_is_symmetric) const DEAL_II_DEPRECATED;
+ void apply_constraints (VectorType &v,
+ const bool matrix_is_symmetric) const DEAL_II_DEPRECATED;
/**
* Apply the constraints to a right hand side vector. This needs to be done
* before starting to solve with the filtered matrix.
*/
- void apply_constraints (VECTOR &v) const;
+ void apply_constraints (VectorType &v) const;
/**
* Matrix-vector multiplication: this operation performs pre_filter(),
* multiplication with the stored matrix and post_filter() in that order.
*/
- void vmult (VECTOR &dst,
- const VECTOR &src) const;
+ void vmult (VectorType &dst,
+ const VectorType &src) const;
/**
* Matrix-vector multiplication: this operation performs pre_filter(),
* transposed multiplication with the stored matrix and post_filter() in
* that order.
*/
- void Tvmult (VECTOR &dst,
- const VECTOR &src) const;
+ void Tvmult (VectorType &dst,
+ const VectorType &src) const;
/**
* Adding matrix-vector multiplication.
* entries set to zero, independent of the previous value of <tt>dst</tt>.
* We excpect that in most cases this is the required behavior.
*/
- void vmult_add (VECTOR &dst,
- const VECTOR &src) const;
+ void vmult_add (VectorType &dst,
+ const VectorType &src) const;
/**
* Adding transpose matrix-vector multiplication:
* entries set to zero, independent of the previous value of <tt>dst</tt>.
* We excpect that in most cases this is the required behavior.
*/
- void Tvmult_add (VECTOR &dst,
- const VECTOR &src) const;
+ void Tvmult_add (VectorType &dst,
+ const VectorType &src) const;
//@}
/**
/**
* Pointer to the sparsity pattern used for this matrix.
*/
- std_cxx11::shared_ptr<PointerMatrixBase<VECTOR> > matrix;
+ std_cxx11::shared_ptr<PointerMatrixBase<VectorType> > matrix;
/**
* Sorted list of pairs denoting the index of the variable and the value to
* Do the pre-filtering step, i.e. zero out those components that belong to
* constrained degrees of freedom.
*/
- void pre_filter (VECTOR &v) const;
+ void pre_filter (VectorType &v) const;
/**
* Do the postfiltering step, i.e. set constrained degrees of freedom to the
* value of the input vector, as the matrix contains only ones on the
* diagonal for these degrees of freedom.
*/
- void post_filter (const VECTOR &in,
- VECTOR &out) const;
+ void post_filter (const VectorType &in,
+ VectorType &out) const;
friend class Accessor;
/**
* FilteredMatrixBlock accesses pre_filter() and post_filter().
*/
- friend class FilteredMatrixBlock<VECTOR>;
+ friend class FilteredMatrixBlock<VectorType>;
};
/*@}*/
//--------------------------------Iterators--------------------------------------//
-template<class VECTOR>
+template<typename VectorType>
inline
-FilteredMatrix<VECTOR>::Accessor::Accessor(
- const FilteredMatrix<VECTOR> *matrix,
- const size_type index)
+FilteredMatrix<VectorType>::Accessor::Accessor
+(const FilteredMatrix<VectorType> *matrix,
+ const size_type index)
:
matrix(matrix),
index(index)
-template<class VECTOR>
+template<typename VectorType>
inline
types::global_dof_index
-FilteredMatrix<VECTOR>::Accessor::row() const
+FilteredMatrix<VectorType>::Accessor::row() const
{
return matrix->constraints[index].first;
}
-template<class VECTOR>
+template<typename VectorType>
inline
types::global_dof_index
-FilteredMatrix<VECTOR>::Accessor::column() const
+FilteredMatrix<VectorType>::Accessor::column() const
{
return matrix->constraints[index].first;
}
-template<class VECTOR>
+template<typename VectorType>
inline
double
-FilteredMatrix<VECTOR>::Accessor::value() const
+FilteredMatrix<VectorType>::Accessor::value() const
{
return matrix->constraints[index].second;
}
-template<class VECTOR>
+template<typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::Accessor::advance()
+FilteredMatrix<VectorType>::Accessor::advance()
{
Assert (index < matrix->constraints.size(), ExcIteratorPastEnd());
++index;
-template<class VECTOR>
+template<typename VectorType>
inline
-FilteredMatrix<VECTOR>::const_iterator::const_iterator(
- const FilteredMatrix<VECTOR> *matrix,
- const size_type index)
+FilteredMatrix<VectorType>::const_iterator::const_iterator
+(const FilteredMatrix<VectorType> *matrix,
+ const size_type index)
:
accessor(matrix, index)
{}
-template<class VECTOR>
+template<typename VectorType>
inline
-typename FilteredMatrix<VECTOR>::const_iterator &
-FilteredMatrix<VECTOR>::const_iterator::operator++ ()
+typename FilteredMatrix<VectorType>::const_iterator &
+FilteredMatrix<VectorType>::const_iterator::operator++ ()
{
accessor.advance();
return *this;
}
-template <class VECTOR>
+template <typename VectorType>
inline
bool
-FilteredMatrix<VECTOR>::PairComparison::
+FilteredMatrix<VectorType>::PairComparison::
operator () (const IndexValuePair &i1,
const IndexValuePair &i2) const
{
-template <class VECTOR>
+template <typename VectorType>
template <class MATRIX>
inline
void
-FilteredMatrix<VECTOR>::initialize (const MATRIX &m, bool ecs)
+FilteredMatrix<VectorType>::initialize (const MATRIX &m, bool ecs)
{
- matrix.reset (new_pointer_matrix_base(m, VECTOR()));
+ matrix.reset (new_pointer_matrix_base(m, VectorType()));
expect_constrained_source = ecs;
}
-template <class VECTOR>
+template <typename VectorType>
inline
-FilteredMatrix<VECTOR>::FilteredMatrix ()
+FilteredMatrix<VectorType>::FilteredMatrix ()
{}
-template <class VECTOR>
+template <typename VectorType>
inline
-FilteredMatrix<VECTOR>::FilteredMatrix (const FilteredMatrix &fm)
+FilteredMatrix<VectorType>::FilteredMatrix (const FilteredMatrix &fm)
:
Subscriptor(),
expect_constrained_source(fm.expect_constrained_source),
-template <class VECTOR>
+template <typename VectorType>
template <class MATRIX>
inline
-FilteredMatrix<VECTOR>::
+FilteredMatrix<VectorType>::
FilteredMatrix (const MATRIX &m, bool ecs)
{
initialize (m, ecs);
-template <class VECTOR>
+template <typename VectorType>
inline
-FilteredMatrix<VECTOR> &
-FilteredMatrix<VECTOR>::operator = (const FilteredMatrix &fm)
+FilteredMatrix<VectorType> &
+FilteredMatrix<VectorType>::operator = (const FilteredMatrix &fm)
{
matrix = fm.matrix;
expect_constrained_source = fm.expect_constrained_source;
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::add_constraint (const size_type index, const double value)
+FilteredMatrix<VectorType>::add_constraint (const size_type index, const double value)
{
// add new constraint to end
constraints.push_back(IndexValuePair(index, value));
-template <class VECTOR>
+template <typename VectorType>
template <class ConstraintList>
inline
void
-FilteredMatrix<VECTOR>::add_constraints (const ConstraintList &new_constraints)
+FilteredMatrix<VectorType>::add_constraints (const ConstraintList &new_constraints)
{
// add new constraints to end
const size_type old_size = constraints.size();
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::clear_constraints ()
+FilteredMatrix<VectorType>::clear_constraints ()
{
// swap vectors to release memory
std::vector<IndexValuePair> empty;
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::clear ()
+FilteredMatrix<VectorType>::clear ()
{
clear_constraints();
matrix.reset();
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::apply_constraints (
- VECTOR &v,
- const bool /* matrix_is_symmetric */) const
+FilteredMatrix<VectorType>::apply_constraints
+(VectorType &v,
+ const bool /* matrix_is_symmetric */) const
{
apply_constraints(v);
}
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::apply_constraints (
- VECTOR &v) const
+FilteredMatrix<VectorType>::apply_constraints (VectorType &v) const
{
- GrowingVectorMemory<VECTOR> mem;
- typename VectorMemory<VECTOR>::Pointer tmp_vector(mem);
+ GrowingVectorMemory<VectorType> mem;
+ typename VectorMemory<VectorType>::Pointer tmp_vector(mem);
tmp_vector->reinit(v);
const_index_value_iterator i = constraints.begin();
const const_index_value_iterator e = constraints.end();
}
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::pre_filter (VECTOR &v) const
+FilteredMatrix<VectorType>::pre_filter (VectorType &v) const
{
// iterate over all constraints and
// zero out value
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::post_filter (const VECTOR &in,
- VECTOR &out) const
+FilteredMatrix<VectorType>::post_filter (const VectorType &in,
+ VectorType &out) const
{
// iterate over all constraints and
// set value correctly
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::vmult (VECTOR &dst, const VECTOR &src) const
+FilteredMatrix<VectorType>::vmult (VectorType &dst, const VectorType &src) const
{
if (!expect_constrained_source)
{
- GrowingVectorMemory<VECTOR> mem;
- VECTOR *tmp_vector = mem.alloc();
+ GrowingVectorMemory<VectorType> mem;
+ VectorType *tmp_vector = mem.alloc();
// first copy over src vector and
// pre-filter
tmp_vector->reinit(src, true);
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::Tvmult (VECTOR &dst, const VECTOR &src) const
+FilteredMatrix<VectorType>::Tvmult (VectorType &dst, const VectorType &src) const
{
if (!expect_constrained_source)
{
- GrowingVectorMemory<VECTOR> mem;
- VECTOR *tmp_vector = mem.alloc();
+ GrowingVectorMemory<VectorType> mem;
+ VectorType *tmp_vector = mem.alloc();
// first copy over src vector and
// pre-filter
tmp_vector->reinit(src, true);
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::vmult_add (VECTOR &dst, const VECTOR &src) const
+FilteredMatrix<VectorType>::vmult_add (VectorType &dst, const VectorType &src) const
{
if (!expect_constrained_source)
{
- GrowingVectorMemory<VECTOR> mem;
- VECTOR *tmp_vector = mem.alloc();
+ GrowingVectorMemory<VectorType> mem;
+ VectorType *tmp_vector = mem.alloc();
// first copy over src vector and
// pre-filter
tmp_vector->reinit(src, true);
-template <class VECTOR>
+template <typename VectorType>
inline
void
-FilteredMatrix<VECTOR>::Tvmult_add (VECTOR &dst, const VECTOR &src) const
+FilteredMatrix<VectorType>::Tvmult_add (VectorType &dst, const VectorType &src) const
{
if (!expect_constrained_source)
{
- GrowingVectorMemory<VECTOR> mem;
- VECTOR *tmp_vector = mem.alloc();
+ GrowingVectorMemory<VectorType> mem;
+ VectorType *tmp_vector = mem.alloc();
// first copy over src vector and
// pre-filter
tmp_vector->reinit(src, true);
-template <class VECTOR>
+template <typename VectorType>
inline
std::size_t
-FilteredMatrix<VECTOR>::memory_consumption () const
+FilteredMatrix<VectorType>::memory_consumption () const
{
return (MemoryConsumption::memory_consumption (matrix) +
MemoryConsumption::memory_consumption (constraints));
/**
* A wrapper to least_squares(), implementing the standard MATRIX interface.
*/
- template<class VECTOR>
- void vmult (VECTOR &dst, const VECTOR &src) const;
+ template<class VectorType>
+ void vmult (VectorType &dst, const VectorType &src) const;
- template<class VECTOR>
- void Tvmult (VECTOR &dst, const VECTOR &src) const;
+ template<class VectorType>
+ void Tvmult (VectorType &dst, const VectorType &src) const;
private:
template <typename number>
-template <class VECTOR>
+template <class VectorType>
void
-Householder<number>::vmult (VECTOR &dst, const VECTOR &src) const
+Householder<number>::vmult (VectorType &dst, const VectorType &src) const
{
least_squares (dst, src);
}
template <typename number>
-template <class VECTOR>
+template <class VectorType>
void
-Householder<number>::Tvmult (VECTOR &, const VECTOR &) const
+Householder<number>::Tvmult (VectorType &, const VectorType &) const
{
Assert(false, ExcNotImplemented());
}
DEAL_II_NAMESPACE_CLOSE
#endif
-
* @author Guido Kanschat
* @date 2010
*/
-template <class VECTOR>
+template <typename VectorType>
class IterativeInverse : public Subscriptor
{
public:
/**
* Solve for right hand side <tt>src</tt>.
*/
- void vmult (VECTOR &dst, const VECTOR &src) const;
+ void vmult (VectorType &dst, const VectorType &src) const;
/**
* Solve for right hand side <tt>src</tt>, but allow for the fact that the
* The solver, which allows selection of the actual solver as well as
* adjustment of parameters.
*/
- SolverSelector<VECTOR> solver;
+ SolverSelector<VectorType> solver;
private:
/**
* The matrix in use.
*/
- std_cxx11::shared_ptr<PointerMatrixBase<VECTOR> > matrix;
+ std_cxx11::shared_ptr<PointerMatrixBase<VectorType> > matrix;
/**
* The preconditioner to use.
*/
- std_cxx11::shared_ptr<PointerMatrixBase<VECTOR> > preconditioner;
+ std_cxx11::shared_ptr<PointerMatrixBase<VectorType> > preconditioner;
};
-template <class VECTOR>
+template <typename VectorType>
template <class MATRIX, class PRECONDITION>
inline
void
-IterativeInverse<VECTOR>::initialize(const MATRIX &m, const PRECONDITION &p)
+IterativeInverse<VectorType>::initialize(const MATRIX &m, const PRECONDITION &p)
{
// dummy variable
- VECTOR *v = 0;
- matrix = std_cxx11::shared_ptr<PointerMatrixBase<VECTOR> > (new_pointer_matrix_base(m, *v));
- preconditioner = std_cxx11::shared_ptr<PointerMatrixBase<VECTOR> > (new_pointer_matrix_base(p, *v));
+ VectorType *v = 0;
+ matrix = std_cxx11::shared_ptr<PointerMatrixBase<VectorType> > (new_pointer_matrix_base(m, *v));
+ preconditioner = std_cxx11::shared_ptr<PointerMatrixBase<VectorType> > (new_pointer_matrix_base(p, *v));
}
-template <class VECTOR>
+template <typename VectorType>
inline
void
-IterativeInverse<VECTOR>::clear()
+IterativeInverse<VectorType>::clear()
{
matrix = 0;
preconditioner = 0;
}
-template <class VECTOR>
+template <typename VectorType>
inline void
-IterativeInverse<VECTOR>::vmult (VECTOR &dst, const VECTOR &src) const
+IterativeInverse<VectorType>::vmult (VectorType &dst, const VectorType &src) const
{
Assert(matrix.get() != 0, ExcNotInitialized());
Assert(preconditioner.get() != 0, ExcNotInitialized());
}
-template <class VECTOR>
+template <typename VectorType>
template <class VECTOR2>
inline void
-IterativeInverse<VECTOR>::vmult (VECTOR2 &dst, const VECTOR2 &src) const
+IterativeInverse<VectorType>::vmult (VECTOR2 &dst, const VECTOR2 &src) const
{
Assert(matrix.get() != 0, ExcNotInitialized());
Assert(preconditioner.get() != 0, ExcNotInitialized());
- GrowingVectorMemory<VECTOR> mem;
- typename VectorMemory<VECTOR>::Pointer sol(mem);
- typename VectorMemory<VECTOR>::Pointer rhs(mem);
+ GrowingVectorMemory<VectorType> mem;
+ typename VectorMemory<VectorType>::Pointer sol(mem);
+ typename VectorMemory<VectorType>::Pointer rhs(mem);
sol->reinit(dst);
*rhs = src;
solver.solve(*matrix, *sol, *rhs, *preconditioner);
DEAL_II_NAMESPACE_CLOSE
#endif
-
-
* No index computations are done, thus, the vectors need to have sizes
* matching #matrix.
*/
- template<class VECTOR>
- void vmult (VECTOR &w, const VECTOR &v) const;
+ template<class VectorType>
+ void vmult (VectorType &w, const VectorType &v) const;
/**
* Matrix-vector-multiplication, forwarding to the same function in MATRIX.
* No index computations are done, thus, the vectors need to have sizes
* matching #matrix.
*/
- template<class VECTOR>
- void vmult_add (VECTOR &w, const VECTOR &v) const;
+ template<class VectorType>
+ void vmult_add (VectorType &w, const VectorType &v) const;
/**
* Matrix-vector-multiplication, forwarding to the same function in MATRIX.
* No index computations are done, thus, the vectors need to have sizes
* matching #matrix.
*/
- template<class VECTOR>
- void Tvmult (VECTOR &w, const VECTOR &v) const;
+ template<class VectorType>
+ void Tvmult (VectorType &w, const VectorType &v) const;
/**
* Matrix-vector-multiplication, forwarding to the same function in MATRIX.
* No index computations are done, thus, the vectors need to have sizes
* matching #matrix.
*/
- template<class VECTOR>
- void Tvmult_add (VECTOR &w, const VECTOR &v) const;
+ template<class VectorType>
+ void Tvmult_add (VectorType &w, const VectorType &v) const;
/**
* The memory used by this object.
template <class MATRIX>
-template <class VECTOR>
+template <class VectorType>
inline
void
-MatrixBlock<MATRIX>::vmult (VECTOR &w, const VECTOR &v) const
+MatrixBlock<MATRIX>::vmult (VectorType &w, const VectorType &v) const
{
matrix.vmult(w,v);
}
template <class MATRIX>
-template <class VECTOR>
+template <class VectorType>
inline
void
-MatrixBlock<MATRIX>::vmult_add (VECTOR &w, const VECTOR &v) const
+MatrixBlock<MATRIX>::vmult_add (VectorType &w, const VectorType &v) const
{
matrix.vmult_add(w,v);
}
template <class MATRIX>
-template <class VECTOR>
+template <class VectorType>
inline
void
-MatrixBlock<MATRIX>::Tvmult (VECTOR &w, const VECTOR &v) const
+MatrixBlock<MATRIX>::Tvmult (VectorType &w, const VectorType &v) const
{
matrix.Tvmult(w,v);
}
template <class MATRIX>
-template <class VECTOR>
+template <class VectorType>
inline
void
-MatrixBlock<MATRIX>::Tvmult_add (VECTOR &w, const VECTOR &v) const
+MatrixBlock<MATRIX>::Tvmult_add (VectorType &w, const VectorType &v) const
{
matrix.Tvmult_add(w,v);
}
*
* @author Denis Davydov, 2015.
*/
-template <typename VECTOR>
+template <typename VectorType>
class PArpackSolver : public Subscriptor
{
public:
/**
* Apply <code>A-sigma * B</code>
*/
- void vmult (VECTOR &dst, const VECTOR &src) const
+ void vmult (VectorType &dst, const VectorType &src) const
{
B.vmult(dst,src);
dst *= (-sigma);
/**
* Apply <code>A^T-sigma * B^T</code>
*/
- void Tvmult (VECTOR &dst, const VECTOR &src) const
+ void Tvmult (VectorType &dst, const VectorType &src) const
{
B.Tvmult(dst,src);
dst *= (-sigma);
*/
template <typename MATRIX1,
typename MATRIX2, typename INVERSE>
- void solve(
- const MATRIX1 &A,
- const MATRIX2 &B,
- const INVERSE &inverse,
- std::vector<std::complex<double> > &eigenvalues,
- std::vector<VECTOR> &eigenvectors,
- const unsigned int n_eigenvalues);
+ void solve
+ (const MATRIX1 &A,
+ const MATRIX2 &B,
+ const INVERSE &inverse,
+ std::vector<std::complex<double> > &eigenvalues,
+ std::vector<VectorType> &eigenvectors,
+ const unsigned int n_eigenvalues);
std::size_t memory_consumption() const;
/**
* Temporary vectors used between Arpack and deal.II
*/
- VECTOR src,dst,tmp;
+ VectorType src,dst,tmp;
/**
* Indices of local degrees of freedom.
<< " Arnoldi vectors.");
};
-template <typename VECTOR>
+template <typename VectorType>
std::size_t
-PArpackSolver<VECTOR>::memory_consumption() const
+PArpackSolver<VectorType>::memory_consumption() const
{
return MemoryConsumption::memory_consumption (double()) *
(workl.size() +
MemoryConsumption::memory_consumption (types::global_dof_index()) * local_indices.size();
}
-template <typename VECTOR>
-PArpackSolver<VECTOR>::AdditionalData::
-AdditionalData (const unsigned int number_of_arnoldi_vectors,
+template <typename VectorType>
+PArpackSolver<VectorType>::AdditionalData::
+AdditionalData (const unsigned int number_of_arnoldi_vectors,
const WhichEigenvalues eigenvalue_of_interest,
- const bool symmetric)
+ const bool symmetric)
:
number_of_arnoldi_vectors(number_of_arnoldi_vectors),
eigenvalue_of_interest(eigenvalue_of_interest),
symmetric(symmetric)
{}
-template <typename VECTOR>
-PArpackSolver<VECTOR>::PArpackSolver (SolverControl &control,
- const MPI_Comm &mpi_communicator,
- const AdditionalData &data)
+template <typename VectorType>
+PArpackSolver<VectorType>::PArpackSolver (SolverControl &control,
+ const MPI_Comm &mpi_communicator,
+ const AdditionalData &data)
:
solver_control (control),
additional_data (data),
{}
-template <typename VECTOR>
-void PArpackSolver<VECTOR>::set_shift(const double s )
+template <typename VectorType>
+void PArpackSolver<VectorType>::set_shift(const double s )
{
shift_value = s;
}
-template <typename VECTOR>
-void PArpackSolver<VECTOR>::reinit(const dealii::IndexSet &locally_owned_dofs)
+template <typename VectorType>
+void PArpackSolver<VectorType>::reinit(const dealii::IndexSet &locally_owned_dofs)
{
// store local indices to write to vectors
locally_owned_dofs.fill_index_vector(local_indices);
}
-template <typename VECTOR>
+template <typename VectorType>
template <typename MATRIX1,typename MATRIX2, typename INVERSE>
-void PArpackSolver<VECTOR>::solve (
- const MATRIX1 &/*system_matrix*/,
- const MATRIX2 &mass_matrix,
- const INVERSE &inverse,
- std::vector<std::complex<double> > &eigenvalues,
- std::vector<VECTOR> &eigenvectors,
- const unsigned int n_eigenvalues)
+void PArpackSolver<VectorType>::solve
+(const MATRIX1 &/*system_matrix*/,
+ const MATRIX2 &mass_matrix,
+ const INVERSE &inverse,
+ std::vector<std::complex<double> > &eigenvalues,
+ std::vector<VectorType> &eigenvectors,
+ const unsigned int n_eigenvalues)
{
Assert (n_eigenvalues <= eigenvectors.size(),
}
-template <typename VECTOR>
-SolverControl &PArpackSolver<VECTOR>::control () const
+template <typename VectorType>
+SolverControl &PArpackSolver<VectorType>::control () const
{
return solver_control;
}
DEAL_II_NAMESPACE_OPEN
-template<class VECTOR> class VectorMemory;
+template<typename VectorType> class VectorMemory;
class IdentityMatrix;
template <typename number> class FullMatrix;
*
* @author Guido Kanschat, 2000, 2001, 2002
*/
-template<class VECTOR>
+template<typename VectorType>
class PointerMatrixBase : public Subscriptor
{
public:
* This was defined to make this matrix a possible template argument to
* BlockMatrixArray.
*/
- typedef typename VECTOR::value_type value_type;
+ typedef typename VectorType::value_type value_type;
/**
* Virtual destructor. Does nothing except making sure that the destructor
/**
* Matrix-vector product.
*/
- virtual void vmult (VECTOR &dst,
- const VECTOR &src) const = 0;
+ virtual void vmult (VectorType &dst,
+ const VectorType &src) const = 0;
/**
* Transposed matrix-vector product.
*/
- virtual void Tvmult (VECTOR &dst,
- const VECTOR &src) const = 0;
+ virtual void Tvmult (VectorType &dst,
+ const VectorType &src) const = 0;
/**
* Matrix-vector product, adding to <tt>dst</tt>.
*/
- virtual void vmult_add (VECTOR &dst,
- const VECTOR &src) const = 0;
+ virtual void vmult_add (VectorType &dst,
+ const VectorType &src) const = 0;
/**
* Transposed matrix-vector product, adding to <tt>dst</tt>.
*/
- virtual void Tvmult_add (VECTOR &dst,
- const VECTOR &src) const = 0;
+ virtual void Tvmult_add (VectorType &dst,
+ const VectorType &src) const = 0;
};
*
* @author Guido Kanschat 2000, 2001, 2002
*/
-template<class MATRIX, class VECTOR>
-class PointerMatrix : public PointerMatrixBase<VECTOR>
+template<class MATRIX, typename VectorType>
+class PointerMatrix : public PointerMatrixBase<VectorType>
{
public:
/**
/**
* Matrix-vector product.
*/
- virtual void vmult (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void vmult (VectorType &dst,
+ const VectorType &src) const;
/**
* Transposed matrix-vector product.
*/
- virtual void Tvmult (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void Tvmult (VectorType &dst,
+ const VectorType &src) const;
/**
* Matrix-vector product, adding to <tt>dst</tt>.
*/
- virtual void vmult_add (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void vmult_add (VectorType &dst,
+ const VectorType &src) const;
/**
* Transposed matrix-vector product, adding to <tt>dst</tt>.
*/
- virtual void Tvmult_add (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void Tvmult_add (VectorType &dst,
+ const VectorType &src) const;
private:
/**
* The pointer to the actual matrix.
*/
- SmartPointer<const MATRIX,PointerMatrix<MATRIX,VECTOR> > m;
+ SmartPointer<const MATRIX,PointerMatrix<MATRIX,VectorType> > m;
};
*
* @author Guido Kanschat 2006
*/
-template<class MATRIX, class VECTOR>
-class PointerMatrixAux : public PointerMatrixBase<VECTOR>
+template<class MATRIX, typename VectorType>
+class PointerMatrixAux : public PointerMatrixBase<VectorType>
{
public:
/**
*
* If <tt>mem</tt> is zero, then GrowingVectorMemory is used.
*/
- PointerMatrixAux (VectorMemory<VECTOR> *mem = 0,
- const MATRIX *M=0);
+ PointerMatrixAux (VectorMemory<VectorType> *mem = 0,
+ const MATRIX *M=0);
/**
* Constructor not using a matrix.
* argument to this function is used to this end, i.e., you can in essence
* assign a name to the current PointerMatrix object.
*/
- PointerMatrixAux(VectorMemory<VECTOR> *mem,
- const char *name);
+ PointerMatrixAux(VectorMemory<VectorType> *mem,
+ const char *name);
/**
* Constructor. <tt>M</tt> points to a matrix which must live longer than
* argument to this function is used to this end, i.e., you can in essence
* assign a name to the current PointerMatrix object.
*/
- PointerMatrixAux(VectorMemory<VECTOR> *mem,
- const MATRIX *M,
- const char *name);
+ PointerMatrixAux(VectorMemory<VectorType> *mem,
+ const MATRIX *M,
+ const char *name);
// Use doc from base class
virtual void clear();
/**
* Assign a new VectorMemory object for getting auxiliary vectors.
*/
- void set_memory(VectorMemory<VECTOR> *mem);
+ void set_memory(VectorMemory<VectorType> *mem);
/**
* Assign a new matrix pointer. Deletes the old pointer and releases its
/**
* Matrix-vector product.
*/
- virtual void vmult (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void vmult (VectorType &dst,
+ const VectorType &src) const;
/**
* Transposed matrix-vector product.
*/
- virtual void Tvmult (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void Tvmult (VectorType &dst,
+ const VectorType &src) const;
/**
* Matrix-vector product, adding to <tt>dst</tt>.
*/
- virtual void vmult_add (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void vmult_add (VectorType &dst,
+ const VectorType &src) const;
/**
* Transposed matrix-vector product, adding to <tt>dst</tt>.
*/
- virtual void Tvmult_add (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void Tvmult_add (VectorType &dst,
+ const VectorType &src) const;
private:
/**
* The backup memory if none was provided.
*/
- mutable GrowingVectorMemory<VECTOR> my_memory;
+ mutable GrowingVectorMemory<VectorType> my_memory;
/**
* Object for getting the auxiliary vector.
*/
- mutable SmartPointer<VectorMemory<VECTOR>,PointerMatrixAux<MATRIX,VECTOR> > mem;
+ mutable SmartPointer<VectorMemory<VectorType>,PointerMatrixAux<MATRIX,VectorType> > mem;
/**
* The pointer to the actual matrix.
*/
- SmartPointer<const MATRIX,PointerMatrixAux<MATRIX,VECTOR> > m;
+ SmartPointer<const MATRIX,PointerMatrixAux<MATRIX,VectorType> > m;
};
* should overload the function in order to save memory and time.
*
* The result is a PointerMatrixBase* pointing to <tt>matrix</tt>. The
- * <tt>VECTOR</tt> argument is a dummy just used to determine the template
+ * <tt>VectorType</tt> argument is a dummy just used to determine the template
* arguments.
*
* @relates PointerMatrixBase @relates PointerMatrixAux
*/
-template <class VECTOR, class MATRIX>
+template <typename VectorType, class MATRIX>
inline
-PointerMatrixBase<VECTOR> *
-new_pointer_matrix_base(MATRIX &matrix, const VECTOR &, const char *name = "PointerMatrixAux")
+PointerMatrixBase<VectorType> *
+new_pointer_matrix_base(MATRIX &matrix, const VectorType &, const char *name = "PointerMatrixAux")
{
- return new PointerMatrixAux<MATRIX, VECTOR>(0, &matrix, name);
+ return new PointerMatrixAux<MATRIX, VectorType>(0, &matrix, name);
}
/**
*
* @relates PointerMatrixBase @relates PointerMatrix
*/
-template <class VECTOR, typename numberm>
-PointerMatrixBase<VECTOR> *
-new_pointer_matrix_base(const BlockSparseMatrix<numberm> &matrix, const VECTOR &, const char *name = "PointerMatrix")
+template <typename VectorType, typename numberm>
+PointerMatrixBase<VectorType> *
+new_pointer_matrix_base(const BlockSparseMatrix<numberm> &matrix, const VectorType &, const char *name = "PointerMatrix")
{
- return new PointerMatrix<BlockSparseMatrix<numberm>, VECTOR>(&matrix, name);
+ return new PointerMatrix<BlockSparseMatrix<numberm>, VectorType>(&matrix, name);
}
*
* @relates PointerMatrixBase @relates PointerMatrix
*/
-template <class VECTOR, typename numberm>
-PointerMatrixBase<VECTOR> *
-new_pointer_matrix_base(const BlockSparseMatrixEZ<numberm> &matrix, const VECTOR &, const char *name = "PointerMatrix")
+template <typename VectorType, typename numberm>
+PointerMatrixBase<VectorType> *
+new_pointer_matrix_base(const BlockSparseMatrixEZ<numberm> &matrix, const VectorType &, const char *name = "PointerMatrix")
{
- return new PointerMatrix<BlockSparseMatrixEZ<numberm>, VECTOR>(&matrix, name);
+ return new PointerMatrix<BlockSparseMatrixEZ<numberm>, VectorType>(&matrix, name);
}
*
* @relates PointerMatrixBase @relates PointerMatrix
*/
-template <typename numberv, typename numberm, typename BLOCK_VECTOR>
-PointerMatrixBase<BLOCK_VECTOR> *
-new_pointer_matrix_base(const BlockMatrixArray<numberm,BLOCK_VECTOR> &matrix, const BLOCK_VECTOR &, const char *name = "PointerMatrix")
+template <typename numberv, typename numberm, typename BLOCK_VectorType>
+PointerMatrixBase<BLOCK_VectorType> *
+new_pointer_matrix_base(const BlockMatrixArray<numberm,BLOCK_VectorType> &matrix, const BLOCK_VectorType &, const char *name = "PointerMatrix")
{
- return new PointerMatrix<BlockMatrixArray<numberm,BLOCK_VECTOR>, BlockVector<numberv> >(&matrix, name);
+ return new PointerMatrix<BlockMatrixArray<numberm,BLOCK_VectorType>, BlockVector<numberv> >(&matrix, name);
}
/*@}*/
//---------------------------------------------------------------------------
-template<class VECTOR>
+template<typename VectorType>
inline
-PointerMatrixBase<VECTOR>::~PointerMatrixBase ()
+PointerMatrixBase<VectorType>::~PointerMatrixBase ()
{}
//----------------------------------------------------------------------//
-template<class MATRIX, class VECTOR>
-PointerMatrix<MATRIX, VECTOR>::PointerMatrix (const MATRIX *M)
+template<class MATRIX, typename VectorType>
+PointerMatrix<MATRIX, VectorType>::PointerMatrix (const MATRIX *M)
: m(M, typeid(*this).name())
{}
-template<class MATRIX, class VECTOR>
-PointerMatrix<MATRIX, VECTOR>::PointerMatrix (const char *name)
+template<class MATRIX, typename VectorType>
+PointerMatrix<MATRIX, VectorType>::PointerMatrix (const char *name)
: m(0, name)
{}
-template<class MATRIX, class VECTOR>
-PointerMatrix<MATRIX, VECTOR>::PointerMatrix (
- const MATRIX *M,
- const char *name)
+template<class MATRIX, typename VectorType>
+PointerMatrix<MATRIX, VectorType>::PointerMatrix (const MATRIX *M,
+ const char *name)
: m(M, name)
{}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrix<MATRIX, VECTOR>::clear ()
+PointerMatrix<MATRIX, VectorType>::clear ()
{
m = 0;
}
-template<class MATRIX, class VECTOR>
-inline const PointerMatrix<MATRIX, VECTOR> &
-PointerMatrix<MATRIX, VECTOR>::operator= (const MATRIX *M)
+template<class MATRIX, typename VectorType>
+inline const PointerMatrix<MATRIX, VectorType> &
+PointerMatrix<MATRIX, VectorType>::operator= (const MATRIX *M)
{
m = M;
return *this;
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline bool
-PointerMatrix<MATRIX, VECTOR>::empty () const
+PointerMatrix<MATRIX, VectorType>::empty () const
{
if (m == 0)
return true;
return m->empty();
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrix<MATRIX, VECTOR>::vmult (VECTOR &dst,
- const VECTOR &src) const
+PointerMatrix<MATRIX, VectorType>::vmult (VectorType &dst,
+ const VectorType &src) const
{
Assert (m != 0, ExcNotInitialized());
m->vmult (dst, src);
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrix<MATRIX, VECTOR>::Tvmult (VECTOR &dst,
- const VECTOR &src) const
+PointerMatrix<MATRIX, VectorType>::Tvmult (VectorType &dst,
+ const VectorType &src) const
{
Assert (m != 0, ExcNotInitialized());
m->Tvmult (dst, src);
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrix<MATRIX, VECTOR>::vmult_add (VECTOR &dst,
- const VECTOR &src) const
+PointerMatrix<MATRIX, VectorType>::vmult_add (VectorType &dst,
+ const VectorType &src) const
{
Assert (m != 0, ExcNotInitialized());
m->vmult_add (dst, src);
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrix<MATRIX, VECTOR>::Tvmult_add (VECTOR &dst,
- const VECTOR &src) const
+PointerMatrix<MATRIX, VectorType>::Tvmult_add (VectorType &dst,
+ const VectorType &src) const
{
Assert (m != 0, ExcNotInitialized());
m->Tvmult_add (dst, src);
//----------------------------------------------------------------------//
-template<class MATRIX, class VECTOR>
-PointerMatrixAux<MATRIX, VECTOR>::PointerMatrixAux (
- VectorMemory<VECTOR> *mem,
+template<class MATRIX, typename VectorType>
+PointerMatrixAux<MATRIX, VectorType>::PointerMatrixAux (
+ VectorMemory<VectorType> *mem,
const MATRIX *M)
: mem(mem, typeid(*this).name()),
m(M, typeid(*this).name())
}
-template<class MATRIX, class VECTOR>
-PointerMatrixAux<MATRIX, VECTOR>::PointerMatrixAux (
- VectorMemory<VECTOR> *mem,
+template<class MATRIX, typename VectorType>
+PointerMatrixAux<MATRIX, VectorType>::PointerMatrixAux (
+ VectorMemory<VectorType> *mem,
const char *name)
: mem(mem, name),
m(0, name)
}
-template<class MATRIX, class VECTOR>
-PointerMatrixAux<MATRIX, VECTOR>::PointerMatrixAux (
- VectorMemory<VECTOR> *mem,
+template<class MATRIX, typename VectorType>
+PointerMatrixAux<MATRIX, VectorType>::PointerMatrixAux (
+ VectorMemory<VectorType> *mem,
const MATRIX *M,
const char *name)
: mem(mem, name),
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrixAux<MATRIX, VECTOR>::clear ()
+PointerMatrixAux<MATRIX, VectorType>::clear ()
{
m = 0;
}
-template<class MATRIX, class VECTOR>
-inline const PointerMatrixAux<MATRIX, VECTOR> &
-PointerMatrixAux<MATRIX, VECTOR>::operator= (const MATRIX *M)
+template<class MATRIX, typename VectorType>
+inline const PointerMatrixAux<MATRIX, VectorType> &
+PointerMatrixAux<MATRIX, VectorType>::operator= (const MATRIX *M)
{
m = M;
return *this;
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrixAux<MATRIX, VECTOR>::set_memory(VectorMemory<VECTOR> *M)
+PointerMatrixAux<MATRIX, VectorType>::set_memory(VectorMemory<VectorType> *M)
{
mem = M;
if (mem == 0)
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline bool
-PointerMatrixAux<MATRIX, VECTOR>::empty () const
+PointerMatrixAux<MATRIX, VectorType>::empty () const
{
if (m == 0)
return true;
return m->empty();
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrixAux<MATRIX, VECTOR>::vmult (VECTOR &dst,
- const VECTOR &src) const
+PointerMatrixAux<MATRIX, VectorType>::vmult (VectorType &dst,
+ const VectorType &src) const
{
if (mem == 0)
mem = &my_memory;
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrixAux<MATRIX, VECTOR>::Tvmult (VECTOR &dst,
- const VECTOR &src) const
+PointerMatrixAux<MATRIX, VectorType>::Tvmult (VectorType &dst,
+ const VectorType &src) const
{
if (mem == 0)
mem = &my_memory;
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrixAux<MATRIX, VECTOR>::vmult_add (VECTOR &dst,
- const VECTOR &src) const
+PointerMatrixAux<MATRIX, VectorType>::vmult_add (VectorType &dst,
+ const VectorType &src) const
{
if (mem == 0)
mem = &my_memory;
Assert (mem != 0, ExcNotInitialized());
Assert (m != 0, ExcNotInitialized());
- VECTOR *v = mem->alloc();
+ VectorType *v = mem->alloc();
v->reinit(dst);
m->vmult (*v, src);
dst += *v;
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-PointerMatrixAux<MATRIX, VECTOR>::Tvmult_add (VECTOR &dst,
- const VECTOR &src) const
+PointerMatrixAux<MATRIX, VectorType>::Tvmult_add (VectorType &dst,
+ const VectorType &src) const
{
if (mem == 0)
mem = &my_memory;
Assert (mem != 0, ExcNotInitialized());
Assert (m != 0, ExcNotInitialized());
- VECTOR *v = mem->alloc();
+ VectorType *v = mem->alloc();
v->reinit(dst);
m->Tvmult (*v, src);
dst += *v;
* LinearOperator class: Jean-Paul Pelteret, 2015
*/
template <class MATRIX = SparseMatrix<double>,
- class VECTOR = dealii::Vector<double> >
+ typename VectorType = dealii::Vector<double> >
class PreconditionSelector : public Subscriptor
{
public:
* Constructor. @p omega denotes the damping parameter of the
* preconditioning.
*/
- PreconditionSelector(const std::string &preconditioning,
- const typename VECTOR::value_type &omega=1.);
+ PreconditionSelector (const std::string &preconditioning,
+ const typename VectorType::value_type &omega=1.);
/**
* Destructor.
* Precondition procedure. Calls the preconditioning that was specified in
* the constructor.
*/
- virtual void vmult (VECTOR &dst, const VECTOR &src) const;
+ virtual void vmult (VectorType &dst, const VectorType &src) const;
/**
* Transpose precondition procedure. Calls the preconditioning that was
* specified in the constructor.
*/
- virtual void Tvmult (VECTOR &dst, const VECTOR &src) const;
+ virtual void Tvmult (VectorType &dst, const VectorType &src) const;
/**
* Get the names of all implemented preconditionings.
* Matrix that is used for the matrix-builtin preconditioning function. cf.
* also @p PreconditionUseMatrix.
*/
- SmartPointer<const MATRIX,PreconditionSelector<MATRIX,VECTOR> > A;
+ SmartPointer<const MATRIX,PreconditionSelector<MATRIX,VectorType> > A;
/**
* Stores the damping parameter of the preconditioner.
*/
- const typename VECTOR::value_type omega;
+ const typename VectorType::value_type omega;
};
/*@}*/
/* --------------------- Inline and template functions ------------------- */
-template <class MATRIX, class VECTOR>
-PreconditionSelector<MATRIX,VECTOR>
-::PreconditionSelector(const std::string &preconditioning,
- const typename VECTOR::value_type &omega) :
+template <class MATRIX, typename VectorType>
+PreconditionSelector<MATRIX,VectorType>
+::PreconditionSelector(const std::string &preconditioning,
+ const typename VectorType::value_type &omega) :
preconditioning(preconditioning),
omega(omega) {}
-template <class MATRIX, class VECTOR>
-PreconditionSelector<MATRIX,VECTOR>::~PreconditionSelector()
+template <class MATRIX, typename VectorType>
+PreconditionSelector<MATRIX,VectorType>::~PreconditionSelector()
{
// release the matrix A
A=0;
}
-template <class MATRIX, class VECTOR>
-void PreconditionSelector<MATRIX,VECTOR>::use_matrix(const MATRIX &M)
+template <class MATRIX, typename VectorType>
+void PreconditionSelector<MATRIX,VectorType>::use_matrix(const MATRIX &M)
{
A=&M;
}
+<<<<<<< fc87b7c22812a5fc7751b36c66b20e6fa54df72c
-template <class MATRIX, class VECTOR>
-inline typename PreconditionSelector<MATRIX,VECTOR>::size_type
-PreconditionSelector<MATRIX,VECTOR>::m () const
+template <class MATRIX, typename VectorType>
+inline typename PreconditionSelector<MATRIX,VectorType>::size_type
+PreconditionSelector<MATRIX,VectorType>::m () const
{
Assert(A!=0, ExcNoMatrixGivenToUse());
return A->m();
}
-template <class MATRIX, class VECTOR>
-inline typename PreconditionSelector<MATRIX,VECTOR>::size_type
-PreconditionSelector<MATRIX,VECTOR>::n () const
+template <class MATRIX, typename VectorType>
+inline typename PreconditionSelector<MATRIX,VectorType>::size_type
+PreconditionSelector<MATRIX,VectorType>::n () const
{
Assert(A!=0, ExcNoMatrixGivenToUse());
return A->n();
}
-template <class MATRIX, class VECTOR>
-void PreconditionSelector<MATRIX,VECTOR>::vmult (VECTOR &dst,
- const VECTOR &src) const
+
+template <class MATRIX, typename VectorType>
+void PreconditionSelector<MATRIX,VectorType>::vmult (VectorType &dst,
+ const VectorType &src) const
{
if (preconditioning=="none")
{
}
-template <class MATRIX, class VECTOR>
-void PreconditionSelector<MATRIX,VECTOR>::Tvmult (VECTOR &dst,
- const VECTOR &src) const
+template <class MATRIX, typename VectorType>
+void PreconditionSelector<MATRIX,VectorType>::Tvmult (VectorType &dst,
+ const VectorType &src) const
{
if (preconditioning=="none")
{
}
-template <class MATRIX, class VECTOR>
-std::string PreconditionSelector<MATRIX,VECTOR>::get_precondition_names()
+template <class MATRIX, typename VectorType>
+std::string PreconditionSelector<MATRIX,VectorType>::get_precondition_names()
{
return "none|jacobi|sor|ssor";
}
/**
* Matrix-vector-product.
*/
- template <class VECTOR>
- void vmult (VECTOR &dst, const VECTOR &src) const;
+ template <typename VectorType>
+ void vmult (VectorType &dst, const VectorType &src) const;
/**
* Residual.
*/
- template <class VECTOR>
- double residual (VECTOR &dst, const VECTOR &src, const VECTOR &rhs) const;
+ template <typename VectorType>
+ double residual (VectorType &dst, const VectorType &src, const VectorType &rhs) const;
private:
/**
/**
* Auxiliary vector.
*/
- // VECTOR aux;
+ // VectorType aux;
/**
* Shift parameter.
*/
*
* @author Guido Kanschat, 2001
*/
-template<class MATRIX, class MASSMATRIX, class VECTOR>
+template<class MATRIX, class MASSMATRIX, class VectorType>
class ShiftedMatrixGeneralized
{
public:
/**
* Matrix-vector-product.
*/
- void vmult (VECTOR &dst, const VECTOR &src) const;
+ void vmult (VectorType &dst, const VectorType &src) const;
/**
* Residual.
*/
- double residual (VECTOR &dst, const VECTOR &src, const VECTOR &rhs) const;
+ double residual (VectorType &dst, const VectorType &src, const VectorType &rhs) const;
private:
/**
* Storage for base matrix.
*/
- SmartPointer<const MATRIX,ShiftedMatrixGeneralized<MATRIX,MASSMATRIX,VECTOR> > A;
+ SmartPointer<const MATRIX,ShiftedMatrixGeneralized<MATRIX,MASSMATRIX,VectorType> > A;
/**
* Storage for mass matrix.
*/
- SmartPointer<const MASSMATRIX,ShiftedMatrixGeneralized<MATRIX,MASSMATRIX,VECTOR> > M;
+ SmartPointer<const MASSMATRIX,ShiftedMatrixGeneralized<MATRIX,MASSMATRIX,VectorType> > M;
/**
* Auxiliary vector.
*/
- VECTOR aux;
+ VectorType aux;
/**
* Shift parameter.
template <class MATRIX>
-template <class VECTOR>
+template <class VectorType>
inline void
-ShiftedMatrix<MATRIX>::vmult (VECTOR &dst, const VECTOR &src) const
+ShiftedMatrix<MATRIX>::vmult (VectorType &dst, const VectorType &src) const
{
A->vmult(dst, src);
if (sigma != 0.)
template <class MATRIX>
-template <class VECTOR>
+template <class VectorType>
inline double
-ShiftedMatrix<MATRIX>::residual (VECTOR &dst,
- const VECTOR &src,
- const VECTOR &rhs) const
+ShiftedMatrix<MATRIX>::residual (VectorType &dst,
+ const VectorType &src,
+ const VectorType &rhs) const
{
A->vmult(dst, src);
if (sigma != 0.)
//---------------------------------------------------------------------------
-template <class MATRIX, class MASSMATRIX, class VECTOR>
+template <class MATRIX, class MASSMATRIX, class VectorType>
inline
-ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VECTOR>
-::ShiftedMatrixGeneralized (const MATRIX &A,
+ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VectorType>
+::ShiftedMatrixGeneralized (const MATRIX &A,
const MASSMATRIX &M,
- const double sigma)
+ const double sigma)
:
A(&A), M(&M), sigma(sigma)
{}
-template <class MATRIX, class MASSMATRIX, class VECTOR>
+template <class MATRIX, class MASSMATRIX, class VectorType>
inline void
-ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VECTOR>::shift (const double s)
+ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VectorType>::shift (const double s)
{
sigma = s;
}
-template <class MATRIX, class MASSMATRIX, class VECTOR>
+template <class MATRIX, class MASSMATRIX, class VectorType>
inline double
-ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VECTOR>::shift () const
+ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VectorType>::shift () const
{
return sigma;
}
-template <class MATRIX, class MASSMATRIX, class VECTOR>
+template <class MATRIX, class MASSMATRIX, class VectorType>
inline void
-ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VECTOR>::vmult (VECTOR &dst,
- const VECTOR &src) const
+ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VectorType>::vmult (VectorType &dst,
+ const VectorType &src) const
{
A->vmult(dst, src);
if (sigma != 0.)
}
-template <class MATRIX, class MASSMATRIX, class VECTOR>
+template <class MATRIX, class MASSMATRIX, class VectorType>
inline double
-ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VECTOR>::residual (VECTOR &dst,
- const VECTOR &src,
- const VECTOR &rhs) const
+ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VectorType>::residual (VectorType &dst,
+ const VectorType &src,
+ const VectorType &rhs) const
{
A->vmult(dst, src);
if (sigma != 0.)
* public:
* // Application of matrix to vector src.
* // Write result into dst
- * void vmult (VECTOR &dst,
- * const VECTOR &src) const;
+ * void vmult (VectorType &dst,
+ * const VectorType &src) const;
*
* // Application of transpose to a vector.
* // Only used by some iterative methods.
- * void Tvmult (VECTOR &dst,
- * const VECTOR &src) const;
+ * void Tvmult (VectorType &dst,
+ * const VectorType &src) const;
* };
*
*
* @endcode
*
* In addition, for some solvers there has to be a global function
- * <tt>swap(VECTOR &a, VECTOR &b)</tt> that exchanges the values of the two
- * vectors.
+ * <tt>swap(VectorType &a, VectorType &b)</tt> that exchanges the values of
+ * the two vectors.
*
* Finally, the solvers also expect an instantiation of
- * GrowingVectorMemory@<VECTOR@>. These instantiations are provided by the
+ * GrowingVectorMemory@<VectorType@>. These instantiations are provided by the
* deal.II library for the built-in vector types, but must be explicitly added
* for user-provided vector classes. Otherwise, the linker will complain that
* it cannot find the constructors and destructors of GrowingVectorMemory that
* @author Wolfgang Bangerth, Guido Kanschat, Ralf Hartmann, 1997-2001, 2005,
* 2014
*/
-template <class VECTOR = Vector<double> >
+template <class VectorType = Vector<double> >
class Solver : public Subscriptor
{
public:
/**
* A typedef for the underlying vector type
*/
- typedef VECTOR vector_type;
+ typedef VectorType vector_type;
/**
* Constructor. Takes a control object which evaluates the conditions for
* responsibility to guarantee that the lifetime of the two arguments is at
* least as long as that of the solver object.
*/
- Solver (SolverControl &solver_control,
- VectorMemory<VECTOR> &vector_memory);
+ Solver (SolverControl &solver_control,
+ VectorMemory<VectorType> &vector_memory);
/**
* Constructor. Takes a control object which evaluates the conditions for
*/
boost::signals2::connection
connect (const std_cxx11::function<SolverControl::State (const unsigned int iteration,
- const double check_value,
- const VECTOR ¤t_iterate)> &slot);
+ const double check_value,
+ const VectorType ¤t_iterate)> &slot);
* A static vector memory object to be used whenever no such object has been
* given to the constructor.
*/
- mutable GrowingVectorMemory<VECTOR> static_vector_memory;
+ mutable GrowingVectorMemory<VectorType> static_vector_memory;
/**
* A reference to an object that provides memory for auxiliary vectors.
*/
- VectorMemory<VECTOR> &memory;
+ VectorMemory<VectorType> &memory;
private:
/**
* signal's return value indicates that the iteration should be terminated.
*/
boost::signals2::signal<SolverControl::State (const unsigned int iteration,
- const double check_value,
- const VECTOR ¤t_iterate),
+ const double check_value,
+ const VectorType ¤t_iterate),
StateCombiner> iteration_status;
};
/*-------------------------------- Inline functions ------------------------*/
-template <class VECTOR>
+template <class VectorType>
inline
SolverControl::State
-Solver<VECTOR>::StateCombiner::operator ()(const SolverControl::State state1,
+Solver<VectorType>::StateCombiner::operator ()(const SolverControl::State state1,
const SolverControl::State state2) const
{
if ((state1 == SolverControl::failure)
}
-template <class VECTOR>
+template <class VectorType>
template <typename Iterator>
inline
SolverControl::State
-Solver<VECTOR>::StateCombiner::operator ()(const Iterator begin,
+Solver<VectorType>::StateCombiner::operator ()(const Iterator begin,
const Iterator end) const
{
Assert (begin != end, ExcMessage ("You can't combine iterator states if no state is given."));
}
-template<class VECTOR>
+template<class VectorType>
inline
-Solver<VECTOR>::Solver (SolverControl &solver_control,
- VectorMemory<VECTOR> &vector_memory)
+Solver<VectorType>::Solver (SolverControl &solver_control,
+ VectorMemory<VectorType> &vector_memory)
:
memory(vector_memory)
{
-template<class VECTOR>
+template<class VectorType>
inline
-Solver<VECTOR>::Solver (SolverControl &solver_control)
+Solver<VectorType>::Solver (SolverControl &solver_control)
:
// use the static memory object this class owns
memory(static_vector_memory)
-template<class VECTOR>
+template<class VectorType>
inline
boost::signals2::connection
-Solver<VECTOR>::
+Solver<VectorType>::
connect (const std_cxx11::function<SolverControl::State (const unsigned int iteration,
- const double check_value,
- const VECTOR ¤t_iterate)> &slot)
+ const double check_value,
+ const VectorType ¤t_iterate)> &slot)
{
return iteration_status.connect (slot);
}
* to observe the progress of the iteration.
*
*/
-template <class VECTOR = Vector<double> >
-class SolverBicgstab : public Solver<VECTOR>
+template <typename VectorType = Vector<double> >
+class SolverBicgstab : public Solver<VectorType>
{
public:
/**
/**
* Constructor.
*/
- SolverBicgstab (SolverControl &cn,
- VectorMemory<VECTOR> &mem,
- const AdditionalData &data=AdditionalData());
+ SolverBicgstab (SolverControl &cn,
+ VectorMemory<VectorType> &mem,
+ const AdditionalData &data=AdditionalData());
/**
* Constructor. Use an object of type GrowingVectorMemory as a default to
*/
template<class MATRIX, class PRECONDITIONER>
void
- solve (const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
+ solve (const MATRIX &A,
+ VectorType &x,
+ const VectorType &b,
const PRECONDITIONER &precondition);
protected:
* Computation of the stopping criterion.
*/
template <class MATRIX>
- double criterion (const MATRIX &A, const VECTOR &x, const VECTOR &b);
+ double criterion (const MATRIX &A, const VectorType &x, const VectorType &b);
/**
* Interface for derived class. This function gets the current iteration
* for a graphical output of the convergence history.
*/
virtual void print_vectors(const unsigned int step,
- const VECTOR &x,
- const VECTOR &r,
- const VECTOR &d) const;
+ const VectorType &x,
+ const VectorType &r,
+ const VectorType &d) const;
/**
* Auxiliary vector.
*/
- VECTOR *Vx;
+ VectorType *Vx;
/**
* Auxiliary vector.
*/
- VECTOR *Vr;
+ VectorType *Vr;
/**
* Auxiliary vector.
*/
- VECTOR *Vrbar;
+ VectorType *Vrbar;
/**
* Auxiliary vector.
*/
- VECTOR *Vp;
+ VectorType *Vp;
/**
* Auxiliary vector.
*/
- VECTOR *Vy;
+ VectorType *Vy;
/**
* Auxiliary vector.
*/
- VECTOR *Vz;
+ VectorType *Vz;
/**
* Auxiliary vector.
*/
- VECTOR *Vt;
+ VectorType *Vt;
/**
* Auxiliary vector.
*/
- VECTOR *Vv;
+ VectorType *Vv;
/**
* Right hand side vector.
*/
- const VECTOR *Vb;
+ const VectorType *Vb;
/**
* Auxiliary value.
#ifndef DOXYGEN
-template<class VECTOR>
-SolverBicgstab<VECTOR>::IterationResult::IterationResult(const bool breakdown,
- const SolverControl::State state,
- const unsigned int last_step,
- const double last_residual)
+template<typename VectorType>
+SolverBicgstab<VectorType>::IterationResult::IterationResult
+(const bool breakdown,
+ const SolverControl::State state,
+ const unsigned int last_step,
+ const double last_residual)
:
breakdown (breakdown),
state (state),
{}
-template<class VECTOR>
-SolverBicgstab<VECTOR>::SolverBicgstab (SolverControl &cn,
- VectorMemory<VECTOR> &mem,
- const AdditionalData &data)
+template<typename VectorType>
+SolverBicgstab<VectorType>::SolverBicgstab (SolverControl &cn,
+ VectorMemory<VectorType> &mem,
+ const AdditionalData &data)
:
- Solver<VECTOR>(cn,mem),
+ Solver<VectorType>(cn,mem),
additional_data(data)
{}
-template<class VECTOR>
-SolverBicgstab<VECTOR>::SolverBicgstab (SolverControl &cn,
- const AdditionalData &data)
+template<typename VectorType>
+SolverBicgstab<VectorType>::SolverBicgstab (SolverControl &cn,
+ const AdditionalData &data)
:
- Solver<VECTOR>(cn),
+ Solver<VectorType>(cn),
additional_data(data)
{}
-template<class VECTOR>
-SolverBicgstab<VECTOR>::~SolverBicgstab ()
+template<typename VectorType>
+SolverBicgstab<VectorType>::~SolverBicgstab ()
{}
-template <class VECTOR>
+template <typename VectorType>
template <class MATRIX>
double
-SolverBicgstab<VECTOR>::criterion (const MATRIX &A, const VECTOR &x, const VECTOR &b)
+SolverBicgstab<VectorType>::criterion (const MATRIX &A, const VectorType &x, const VectorType &b)
{
A.vmult(*Vt, x);
Vt->add(-1.,b);
-template <class VECTOR >
+template <typename VectorType >
template <class MATRIX>
SolverControl::State
-SolverBicgstab<VECTOR>::start(const MATRIX &A)
+SolverBicgstab<VectorType>::start(const MATRIX &A)
{
A.vmult(*Vr, *Vx);
Vr->sadd(-1.,1.,*Vb);
-template<class VECTOR>
+template<typename VectorType>
void
-SolverBicgstab<VECTOR>::print_vectors(const unsigned int,
- const VECTOR &,
- const VECTOR &,
- const VECTOR &) const
+SolverBicgstab<VectorType>::print_vectors(const unsigned int,
+ const VectorType &,
+ const VectorType &,
+ const VectorType &) const
{}
-template<class VECTOR>
+template<typename VectorType>
template<class MATRIX, class PRECONDITIONER>
-typename SolverBicgstab<VECTOR>::IterationResult
-SolverBicgstab<VECTOR>::iterate(const MATRIX &A,
- const PRECONDITIONER &precondition)
+typename SolverBicgstab<VectorType>::IterationResult
+SolverBicgstab<VectorType>::iterate(const MATRIX &A,
+ const PRECONDITIONER &precondition)
{
//TODO:[GK] Implement "use the length of the computed orthogonal residual" in the BiCGStab method.
SolverControl::State state = SolverControl::iterate;
alpha = omega = rho = 1.;
- VECTOR &r = *Vr;
- VECTOR &rbar = *Vrbar;
- VECTOR &p = *Vp;
- VECTOR &y = *Vy;
- VECTOR &z = *Vz;
- VECTOR &t = *Vt;
- VECTOR &v = *Vv;
+ VectorType &r = *Vr;
+ VectorType &rbar = *Vrbar;
+ VectorType &p = *Vp;
+ VectorType &y = *Vy;
+ VectorType &z = *Vz;
+ VectorType &t = *Vt;
+ VectorType &v = *Vv;
rbar = r;
bool startup = true;
}
-template<class VECTOR>
+template<typename VectorType>
template<class MATRIX, class PRECONDITIONER>
void
-SolverBicgstab<VECTOR>::solve(const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
- const PRECONDITIONER &precondition)
+SolverBicgstab<VectorType>::solve(const MATRIX &A,
+ VectorType &x,
+ const VectorType &b,
+ const PRECONDITIONER &precondition)
{
deallog.push("Bicgstab");
Vr = this->memory.alloc();
* automatically, avoiding restart.
*/
- template <class VECTOR>
+ template <typename VectorType>
class TmpVectors
{
public:
/**
- * Constructor. Prepares an array of @p VECTOR of length @p max_size.
+ * Constructor. Prepares an array of @p VectorType of length @p max_size.
*/
- TmpVectors(const unsigned int max_size,
- VectorMemory<VECTOR> &vmem);
+ TmpVectors(const unsigned int max_size,
+ VectorMemory<VectorType> &vmem);
/**
* Delete all allocated vectors.
* Get vector number @p i. If this vector was unused before, an error
* occurs.
*/
- VECTOR &operator[] (const unsigned int i) const;
+ VectorType &operator[] (const unsigned int i) const;
/**
* Get vector number @p i. Allocate it if necessary.
* If a vector must be allocated, @p temp is used to reinit it to the
* proper dimensions.
*/
- VECTOR &operator() (const unsigned int i,
- const VECTOR &temp);
+ VectorType &operator() (const unsigned int i,
+ const VectorType &temp);
private:
/**
* Pool were vectors are obtained from.
*/
- VectorMemory<VECTOR> &mem;
+ VectorMemory<VectorType> &mem;
/**
* Field for storing the vectors.
*/
- std::vector<VECTOR *> data;
+ std::vector<VectorType *> data;
/**
* Offset of the first vector. This is for later when vector rotation
*
* @author Wolfgang Bangerth, Guido Kanschat, Ralf Hartmann.
*/
-template <class VECTOR = Vector<double> >
-class SolverGMRES : public Solver<VECTOR>
+template <class VectorType = Vector<double> >
+class SolverGMRES : public Solver<VectorType>
{
public:
/**
/**
* Constructor.
*/
- SolverGMRES (SolverControl &cn,
- VectorMemory<VECTOR> &mem,
- const AdditionalData &data=AdditionalData());
+ SolverGMRES (SolverControl &cn,
+ VectorMemory<VectorType> &mem,
+ const AdditionalData &data=AdditionalData());
/**
* Constructor. Use an object of type GrowingVectorMemory as a default to
template<class MATRIX, class PRECONDITIONER>
void
solve (const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
+ VectorType &x,
+ const VectorType &b,
const PRECONDITIONER &precondition);
/**
* All subsequent iterations use re-orthogonalization.
*/
static double
- modified_gram_schmidt (const internal::SolverGMRES::TmpVectors<VECTOR> &orthogonal_vectors,
- const unsigned int dim,
- const unsigned int accumulated_iterations,
- VECTOR &vv,
- Vector<double> &h,
- bool &re_orthogonalize);
+ modified_gram_schmidt (const internal::SolverGMRES::TmpVectors<VectorType> &orthogonal_vectors,
+ const unsigned int dim,
+ const unsigned int accumulated_iterations,
+ VectorType &vv,
+ Vector<double> &h,
+ bool &re_orthogonalize);
/**
* Estimates the eigenvalues from the Hessenberg matrix, H_orig, generated
/**
* No copy constructor.
*/
- SolverGMRES (const SolverGMRES<VECTOR> &);
+ SolverGMRES (const SolverGMRES<VectorType> &);
};
/**
*
* @author Guido Kanschat, 2003
*/
-template <class VECTOR = Vector<double> >
-class SolverFGMRES : public Solver<VECTOR>
+template <class VectorType = Vector<double> >
+class SolverFGMRES : public Solver<VectorType>
{
public:
/**
/**
* Constructor.
*/
- SolverFGMRES (SolverControl &cn,
- VectorMemory<VECTOR> &mem,
- const AdditionalData &data=AdditionalData());
+ SolverFGMRES (SolverControl &cn,
+ VectorMemory<VectorType> &mem,
+ const AdditionalData &data=AdditionalData());
/**
* Constructor. Use an object of type GrowingVectorMemory as a default to
template<class MATRIX, class PRECONDITIONER>
void
solve (const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
+ VectorType &x,
+ const VectorType &b,
const PRECONDITIONER &precondition);
private:
{
namespace SolverGMRES
{
- template <class VECTOR>
+ template <class VectorType>
inline
- TmpVectors<VECTOR>::
- TmpVectors (const unsigned int max_size,
- VectorMemory<VECTOR> &vmem)
+ TmpVectors<VectorType>::
+ TmpVectors (const unsigned int max_size,
+ VectorMemory<VectorType> &vmem)
:
mem(vmem),
data (max_size, 0),
{}
- template <class VECTOR>
+ template <class VectorType>
inline
- TmpVectors<VECTOR>::~TmpVectors ()
+ TmpVectors<VectorType>::~TmpVectors ()
{
- for (typename std::vector<VECTOR *>::iterator v = data.begin();
+ for (typename std::vector<VectorType *>::iterator v = data.begin();
v != data.end(); ++v)
if (*v != 0)
mem.free(*v);
}
- template <class VECTOR>
- inline VECTOR &
- TmpVectors<VECTOR>::operator[] (const unsigned int i) const
+ template <class VectorType>
+ inline VectorType &
+ TmpVectors<VectorType>::operator[] (const unsigned int i) const
{
Assert (i+offset<data.size(),
ExcIndexRange(i, -offset, data.size()-offset));
}
- template <class VECTOR>
- inline VECTOR &
- TmpVectors<VECTOR>::operator() (const unsigned int i,
- const VECTOR &temp)
+ template <class VectorType>
+ inline VectorType &
+ TmpVectors<VectorType>::operator() (const unsigned int i,
+ const VectorType &temp)
{
Assert (i+offset<data.size(),
ExcIndexRange(i,-offset, data.size()-offset));
-template <class VECTOR>
+template <class VectorType>
inline
-SolverGMRES<VECTOR>::AdditionalData::
+SolverGMRES<VectorType>::AdditionalData::
AdditionalData (const unsigned int max_n_tmp_vectors,
const bool right_preconditioning,
const bool use_default_residual,
-template <class VECTOR>
+template <class VectorType>
inline
-SolverGMRES<VECTOR>::AdditionalData::
+SolverGMRES<VectorType>::AdditionalData::
AdditionalData (const unsigned int max_n_tmp_vectors,
const bool right_preconditioning,
const bool use_default_residual,
-template <class VECTOR>
-SolverGMRES<VECTOR>::SolverGMRES (SolverControl &cn,
- VectorMemory<VECTOR> &mem,
- const AdditionalData &data)
+template <class VectorType>
+SolverGMRES<VectorType>::SolverGMRES (SolverControl &cn,
+ VectorMemory<VectorType> &mem,
+ const AdditionalData &data)
:
- Solver<VECTOR> (cn,mem),
+ Solver<VectorType> (cn,mem),
additional_data(data)
{}
-template <class VECTOR>
-SolverGMRES<VECTOR>::SolverGMRES (SolverControl &cn,
- const AdditionalData &data) :
- Solver<VECTOR> (cn),
+template <class VectorType>
+SolverGMRES<VectorType>::SolverGMRES (SolverControl &cn,
+ const AdditionalData &data) :
+ Solver<VectorType> (cn),
additional_data(data)
{}
-template <class VECTOR>
+template <class VectorType>
inline
void
-SolverGMRES<VECTOR>::givens_rotation (Vector<double> &h,
- Vector<double> &b,
- Vector<double> &ci,
- Vector<double> &si,
- int col) const
+SolverGMRES<VectorType>::givens_rotation (Vector<double> &h,
+ Vector<double> &b,
+ Vector<double> &ci,
+ Vector<double> &si,
+ int col) const
{
for (int i=0 ; i<col ; i++)
{
-template <class VECTOR>
+template <class VectorType>
inline
double
-SolverGMRES<VECTOR>::modified_gram_schmidt (const internal::SolverGMRES::TmpVectors<VECTOR> &orthogonal_vectors,
- const unsigned int dim,
- const unsigned int accumulated_iterations,
- VECTOR &vv,
- Vector<double> &h,
- bool &re_orthogonalize)
+SolverGMRES<VectorType>::modified_gram_schmidt
+(const internal::SolverGMRES::TmpVectors<VectorType> &orthogonal_vectors,
+ const unsigned int dim,
+ const unsigned int accumulated_iterations,
+ VectorType &vv,
+ Vector<double> &h,
+ bool &re_orthogonalize)
{
Assert(dim > 0, ExcInternalError());
const unsigned int inner_iteration = dim - 1;
if (re_orthogonalize == false && inner_iteration % 5 == 4)
{
if (norm_vv > 10. * norm_vv_start *
- std::sqrt(std::numeric_limits<typename VECTOR::value_type>::epsilon()))
+ std::sqrt(std::numeric_limits<typename VectorType::value_type>::epsilon()))
return norm_vv;
else
-template<class VECTOR>
+template<class VectorType>
inline void
-SolverGMRES<VECTOR>::compute_eigs_and_cond(
- const FullMatrix<double> &H_orig,
- const unsigned int dim,
- const boost::signals2::signal<void (const std::vector<std::complex<double> > &)> &eigenvalues_signal,
- const boost::signals2::signal<void (double)> &cond_signal,
- const bool log_eigenvalues)
+SolverGMRES<VectorType>::compute_eigs_and_cond
+(const FullMatrix<double> &H_orig,
+ const unsigned int dim,
+ const boost::signals2::signal<void (const std::vector<std::complex<double> > &)> &eigenvalues_signal,
+ const boost::signals2::signal<void (double)> &cond_signal,
+ const bool log_eigenvalues)
{
//Avoid copying the Hessenberg matrix if it isn't needed.
if (!eigenvalues_signal.empty() || !cond_signal.empty() || log_eigenvalues )
-template<class VECTOR>
+template<class VectorType>
template<class MATRIX, class PRECONDITIONER>
void
-SolverGMRES<VECTOR>::solve (const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
- const PRECONDITIONER &precondition)
+SolverGMRES<VectorType>::solve (const MATRIX &A,
+ VectorType &x,
+ const VectorType &b,
+ const PRECONDITIONER &precondition)
{
// this code was written a very long time ago by people not associated with
// deal.II. we don't make any guarantees to its optimality or that it even
const unsigned int n_tmp_vectors = additional_data.max_n_tmp_vectors;
// Generate an object where basis vectors are stored.
- internal::SolverGMRES::TmpVectors<VECTOR> tmp_vectors (n_tmp_vectors, this->memory);
+ internal::SolverGMRES::TmpVectors<VectorType> tmp_vectors (n_tmp_vectors, this->memory);
// number of the present iteration; this
// number is not reset to zero upon a
const bool use_default_residual = additional_data.use_default_residual;
// define two aliases
- VECTOR &v = tmp_vectors(0, x);
- VECTOR &p = tmp_vectors(n_tmp_vectors-1, x);
+ VectorType &v = tmp_vectors(0, x);
+ VectorType &p = tmp_vectors(n_tmp_vectors-1, x);
// Following vectors are needed
// when not the default residuals
// are used as stopping criterion
- VECTOR *r=0;
- VECTOR *x_=0;
+ VectorType *r=0;
+ VectorType *x_=0;
dealii::Vector<double> *gamma_=0;
if (!use_default_residual)
{
{
++accumulated_iterations;
// yet another alias
- VECTOR &vv = tmp_vectors(inner_iteration+1, x);
+ VectorType &vv = tmp_vectors(inner_iteration+1, x);
if (left_precondition)
{
-template<class VECTOR>
+template<class VectorType>
boost::signals2::connection
-SolverGMRES<VECTOR>::connect_condition_number_slot(
- const std_cxx11::function<void(double)> &slot,
- const bool every_iteration)
+SolverGMRES<VectorType>::connect_condition_number_slot
+(const std_cxx11::function<void(double)> &slot,
+ const bool every_iteration)
{
if (every_iteration)
{
-template<class VECTOR>
+template<class VectorType>
boost::signals2::connection
-SolverGMRES<VECTOR>::connect_eigenvalues_slot(
- const std_cxx11::function<void (const std::vector<std::complex<double> > &)> &slot,
- const bool every_iteration)
+SolverGMRES<VectorType>::connect_eigenvalues_slot
+(const std_cxx11::function<void (const std::vector<std::complex<double> > &)> &slot,
+ const bool every_iteration)
{
if (every_iteration)
{
-template<class VECTOR>
+template<class VectorType>
double
-SolverGMRES<VECTOR>::criterion ()
+SolverGMRES<VectorType>::criterion ()
{
// dummy implementation. this function is not needed for the present
// implementation of gmres
//----------------------------------------------------------------------//
-template <class VECTOR>
-SolverFGMRES<VECTOR>::SolverFGMRES (SolverControl &cn,
- VectorMemory<VECTOR> &mem,
- const AdditionalData &data)
+template <class VectorType>
+SolverFGMRES<VectorType>::SolverFGMRES (SolverControl &cn,
+ VectorMemory<VectorType> &mem,
+ const AdditionalData &data)
:
- Solver<VECTOR> (cn, mem),
+ Solver<VectorType> (cn, mem),
additional_data(data)
{}
-template <class VECTOR>
-SolverFGMRES<VECTOR>::SolverFGMRES (SolverControl &cn,
- const AdditionalData &data)
+template <class VectorType>
+SolverFGMRES<VectorType>::SolverFGMRES (SolverControl &cn,
+ const AdditionalData &data)
:
- Solver<VECTOR> (cn),
+ Solver<VectorType> (cn),
additional_data(data)
{}
-template<class VECTOR>
+template<class VectorType>
template<class MATRIX, class PRECONDITIONER>
void
-SolverFGMRES<VECTOR>::solve (
- const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
- const PRECONDITIONER &precondition)
+SolverFGMRES<VectorType>::solve (const MATRIX &A,
+ VectorType &x,
+ const VectorType &b,
+ const PRECONDITIONER &precondition)
{
deallog.push("FGMRES");
const unsigned int basis_size = additional_data.max_basis_size;
// Generate an object where basis vectors are stored.
- typename internal::SolverGMRES::TmpVectors<VECTOR> v (basis_size, this->memory);
- typename internal::SolverGMRES::TmpVectors<VECTOR> z (basis_size, this->memory);
+ typename internal::SolverGMRES::TmpVectors<VectorType> v (basis_size, this->memory);
+ typename internal::SolverGMRES::TmpVectors<VectorType> z (basis_size, this->memory);
// number of the present iteration; this number is not reset to zero upon a
// restart
// Iteration starts here
double res = -std::numeric_limits<double>::max();
- VECTOR *aux = this->memory.alloc();
+ VectorType *aux = this->memory.alloc();
aux->reinit(x);
do
{
*
* @author Thomas Richter, 2000, Luca Heltai, 2006
*/
-template <class VECTOR = Vector<double> >
-class SolverMinRes : public Solver<VECTOR>
+template <class VectorType = Vector<double> >
+class SolverMinRes : public Solver<VectorType>
{
public:
/**
/**
* Constructor.
*/
- SolverMinRes (SolverControl &cn,
- VectorMemory<VECTOR> &mem,
- const AdditionalData &data=AdditionalData());
+ SolverMinRes (SolverControl &cn,
+ VectorMemory<VectorType> &mem,
+ const AdditionalData &data=AdditionalData());
/**
* Constructor. Use an object of type GrowingVectorMemory as a default to
template<class MATRIX, class PRECONDITIONER>
void
solve (const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
+ VectorType &x,
+ const VectorType &b,
const PRECONDITIONER &precondition);
/**
* for a graphical output of the convergence history.
*/
virtual void print_vectors(const unsigned int step,
- const VECTOR &x,
- const VECTOR &r,
- const VECTOR &d) const;
+ const VectorType &x,
+ const VectorType &r,
+ const VectorType &d) const;
/**
* Temporary vectors, allocated through the @p VectorMemory object at the
* start of the actual solution process and deallocated at the end.
*/
- VECTOR *Vu0, *Vu1, *Vu2;
- VECTOR *Vm0, *Vm1, *Vm2;
- VECTOR *Vv;
+ VectorType *Vu0, *Vu1, *Vu2;
+ VectorType *Vm0, *Vm1, *Vm2;
+ VectorType *Vv;
/**
* Within the iteration loop, the square of the residual vector is stored in
#ifndef DOXYGEN
-template<class VECTOR>
-SolverMinRes<VECTOR>::SolverMinRes (SolverControl &cn,
- VectorMemory<VECTOR> &mem,
- const AdditionalData &)
+template<class VectorType>
+SolverMinRes<VectorType>::SolverMinRes (SolverControl &cn,
+ VectorMemory<VectorType> &mem,
+ const AdditionalData &)
:
- Solver<VECTOR>(cn,mem)
+ Solver<VectorType>(cn,mem)
{}
-template<class VECTOR>
-SolverMinRes<VECTOR>::SolverMinRes (SolverControl &cn,
- const AdditionalData &)
+template<class VectorType>
+SolverMinRes<VectorType>::SolverMinRes (SolverControl &cn,
+ const AdditionalData &)
:
- Solver<VECTOR>(cn)
+ Solver<VectorType>(cn)
{}
-template<class VECTOR>
-SolverMinRes<VECTOR>::~SolverMinRes ()
+template<class VectorType>
+SolverMinRes<VectorType>::~SolverMinRes ()
{}
-template<class VECTOR>
+template<class VectorType>
double
-SolverMinRes<VECTOR>::criterion()
+SolverMinRes<VectorType>::criterion()
{
return res2;
}
-template<class VECTOR>
+template<class VectorType>
void
-SolverMinRes<VECTOR>::print_vectors(const unsigned int,
- const VECTOR &,
- const VECTOR &,
- const VECTOR &) const
+SolverMinRes<VectorType>::print_vectors(const unsigned int,
+ const VectorType &,
+ const VectorType &,
+ const VectorType &) const
{}
-template<class VECTOR>
+template<class VectorType>
template<class MATRIX, class PRECONDITIONER>
void
-SolverMinRes<VECTOR>::solve (const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
- const PRECONDITIONER &precondition)
+SolverMinRes<VectorType>::solve (const MATRIX &A,
+ VectorType &x,
+ const VectorType &b,
+ const PRECONDITIONER &precondition)
{
SolverControl::State conv=SolverControl::iterate;
Vm1 = this->memory.alloc();
Vm2 = this->memory.alloc();
// define some aliases for simpler access
- typedef VECTOR *vecptr;
+ typedef VectorType *vecptr;
vecptr u[3] = {Vu0, Vu1, Vu2};
vecptr m[3] = {Vm0, Vm1, Vm2};
- VECTOR &v = *Vv;
+ VectorType &v = *Vv;
// resize the vectors, but do not set
// the values since they'd be overwritten
// soon anyway.
*
* @author Guido Kanschat, 1999
*/
-template <class VECTOR = Vector<double> >
-class SolverQMRS : public Solver<VECTOR>
+template <typename VectorType = Vector<double> >
+class SolverQMRS : public Solver<VectorType>
{
public:
/**
/**
* Constructor.
*/
- SolverQMRS (SolverControl &cn,
- VectorMemory<VECTOR> &mem,
- const AdditionalData &data=AdditionalData());
+ SolverQMRS (SolverControl &cn,
+ VectorMemory<VectorType> &mem,
+ const AdditionalData &data=AdditionalData());
/**
* Constructor. Use an object of type GrowingVectorMemory as a default to
template<class MATRIX, class PRECONDITIONER>
void
solve (const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
+ VectorType &x,
+ const VectorType &b,
const PRECONDITIONER &precondition);
/**
* vector, the residual and the update vector in each step. It can be used
* for a graphical output of the convergence history.
*/
- virtual void print_vectors(const unsigned int step,
- const VECTOR &x,
- const VECTOR &r,
- const VECTOR &d) const;
+ virtual void print_vectors (const unsigned int step,
+ const VectorType &x,
+ const VectorType &r,
+ const VectorType &d) const;
protected:
/**
* Implementation of the computation of the norm of the residual.
* Temporary vectors, allocated through the @p VectorMemory object at the
* start of the actual solution process and deallocated at the end.
*/
- VECTOR *Vv;
- VECTOR *Vp;
- VECTOR *Vq;
- VECTOR *Vt;
- VECTOR *Vd;
+ VectorType *Vv;
+ VectorType *Vp;
+ VectorType *Vq;
+ VectorType *Vt;
+ VectorType *Vd;
/**
* Iteration vector.
*/
- VECTOR *Vx;
+ VectorType *Vx;
/**
* RHS vector.
*/
- const VECTOR *Vb;
+ const VectorType *Vb;
/**
* Within the iteration loop, the square of the residual vector is stored in
#ifndef DOXYGEN
-template<class VECTOR>
-SolverQMRS<VECTOR>::IterationResult::IterationResult(const SolverControl::State state,
- const double last_residual)
+template<class VectorType>
+SolverQMRS<VectorType>::IterationResult::IterationResult (const SolverControl::State state,
+ const double last_residual)
:
state (state),
last_residual (last_residual)
{}
-template<class VECTOR>
-SolverQMRS<VECTOR>::SolverQMRS(SolverControl &cn,
- VectorMemory<VECTOR> &mem,
- const AdditionalData &data)
+template<class VectorType>
+SolverQMRS<VectorType>::SolverQMRS (SolverControl &cn,
+ VectorMemory<VectorType> &mem,
+ const AdditionalData &data)
:
- Solver<VECTOR>(cn,mem),
+ Solver<VectorType>(cn,mem),
additional_data(data)
{}
-template<class VECTOR>
-SolverQMRS<VECTOR>::SolverQMRS(SolverControl &cn,
- const AdditionalData &data)
+template<class VectorType>
+SolverQMRS<VectorType>::SolverQMRS(SolverControl &cn,
+ const AdditionalData &data)
:
- Solver<VECTOR>(cn),
+ Solver<VectorType>(cn),
additional_data(data)
{}
-template<class VECTOR>
+template<class VectorType>
double
-SolverQMRS<VECTOR>::criterion()
+SolverQMRS<VectorType>::criterion()
{
return std::sqrt(res2);
}
-template<class VECTOR>
+template<class VectorType>
void
-SolverQMRS<VECTOR>::print_vectors(const unsigned int,
- const VECTOR &,
- const VECTOR &,
- const VECTOR &) const
+SolverQMRS<VectorType>::print_vectors(const unsigned int,
+ const VectorType &,
+ const VectorType &,
+ const VectorType &) const
{}
-template<class VECTOR>
+template<class VectorType>
template<class MATRIX, class PRECONDITIONER>
void
-SolverQMRS<VECTOR>::solve (const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
- const PRECONDITIONER &precondition)
+SolverQMRS<VectorType>::solve (const MATRIX &A,
+ VectorType &x,
+ const VectorType &b,
+ const PRECONDITIONER &precondition)
{
deallog.push("QMRS");
-template<class VECTOR>
+template<class VectorType>
template<class MATRIX, class PRECONDITIONER>
-typename SolverQMRS<VECTOR>::IterationResult
-SolverQMRS<VECTOR>::iterate(const MATRIX &A,
- const PRECONDITIONER &precondition)
+typename SolverQMRS<VectorType>::IterationResult
+SolverQMRS<VectorType>::iterate(const MATRIX &A,
+ const PRECONDITIONER &precondition)
{
/* Remark: the matrix A in the article is the preconditioned matrix.
* Therefore, we have to precondition x before we compute the first residual.
SolverControl::State state = SolverControl::iterate;
// define some aliases for simpler access
- VECTOR &v = *Vv;
- VECTOR &p = *Vp;
- VECTOR &q = *Vq;
- VECTOR &t = *Vt;
- VECTOR &d = *Vd;
- VECTOR &x = *Vx;
- const VECTOR &b = *Vb;
+ VectorType &v = *Vv;
+ VectorType &p = *Vp;
+ VectorType &q = *Vq;
+ VectorType &t = *Vt;
+ VectorType &d = *Vd;
+ VectorType &x = *Vx;
+ const VectorType &b = *Vb;
int it=0;
* @author Guido Kanschat
* @date 2010
*/
-template <class VECTOR = Vector<double> >
-class SolverRelaxation : public Solver<VECTOR>
+template <typename VectorType = Vector<double> >
+class SolverRelaxation : public Solver<VectorType>
{
public:
/**
*/
template<class MATRIX, class RELAXATION>
void
- solve (const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
+ solve (const MATRIX &A,
+ VectorType &x,
+ const VectorType &b,
const RELAXATION &R);
};
//----------------------------------------------------------------------//
-template <class VECTOR>
-SolverRelaxation<VECTOR>::SolverRelaxation(SolverControl &cn,
- const AdditionalData &)
+template <class VectorType>
+SolverRelaxation<VectorType>::SolverRelaxation (SolverControl &cn,
+ const AdditionalData &)
:
- Solver<VECTOR> (cn)
+ Solver<VectorType> (cn)
{}
-template <class VECTOR>
-SolverRelaxation<VECTOR>::~SolverRelaxation()
+template <class VectorType>
+SolverRelaxation<VectorType>::~SolverRelaxation()
{}
-template <class VECTOR>
+template <class VectorType>
template <class MATRIX, class RELAXATION>
void
-SolverRelaxation<VECTOR>::solve (
- const MATRIX &A,
- VECTOR &x,
- const VECTOR &b,
- const RELAXATION &R)
+SolverRelaxation<VectorType>::solve (const MATRIX &A,
+ VectorType &x,
+ const VectorType &b,
+ const RELAXATION &R)
{
- GrowingVectorMemory<VECTOR> mem;
+ GrowingVectorMemory<VectorType> mem;
SolverControl::State conv=SolverControl::iterate;
// Memory allocation
- typename VectorMemory<VECTOR>::Pointer Vr(mem);
- VECTOR &r = *Vr;
+ typename VectorMemory<VectorType>::Pointer Vr(mem);
+ VectorType &r = *Vr;
r.reinit(x);
- typename VectorMemory<VECTOR>::Pointer Vd(mem);
- VECTOR &d = *Vd;
+ typename VectorMemory<VectorType>::Pointer Vd(mem);
+ VectorType &d = *Vd;
d.reinit(x);
deallog.push("Relaxation");
template <typename number> class SparseMatrix;
template <typename number> class SparseLUDecomposition;
template <typename number> class SparseILU;
-template <class VECTOR> class VectorSlice;
+template <typename VectorType> class VectorSlice;
namespace ChunkSparsityPatternIterators
{
* @ingroup Matrix2
* @author Guido Kanschat, 2006
*/
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
class
- TransposeMatrix : public PointerMatrixBase<VECTOR>
+ TransposeMatrix : public PointerMatrixBase<VectorType>
{
public:
/**
/**
* Matrix-vector product.
*/
- virtual void vmult (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void vmult (VectorType &dst,
+ const VectorType &src) const;
/**
* Transposed matrix-vector product.
*/
- virtual void Tvmult (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void Tvmult (VectorType &dst,
+ const VectorType &src) const;
/**
* Matrix-vector product, adding to <tt>dst</tt>.
*/
- virtual void vmult_add (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void vmult_add (VectorType &dst,
+ const VectorType &src) const;
/**
* Transposed matrix-vector product, adding to <tt>dst</tt>.
*/
- virtual void Tvmult_add (VECTOR &dst,
- const VECTOR &src) const;
+ virtual void Tvmult_add (VectorType &dst,
+ const VectorType &src) const;
private:
/**
* The pointer to the actual matrix.
*/
- SmartPointer<const MATRIX,TransposeMatrix<MATRIX,VECTOR> > m;
+ SmartPointer<const MATRIX,TransposeMatrix<MATRIX,VectorType> > m;
};
//----------------------------------------------------------------------//
-template<class MATRIX, class VECTOR>
-TransposeMatrix<MATRIX, VECTOR>::TransposeMatrix (const MATRIX *M)
+template<class MATRIX, typename VectorType>
+TransposeMatrix<MATRIX, VectorType>::TransposeMatrix (const MATRIX *M)
: m(M)
{}
-template<class MATRIX, class VECTOR>
-TransposeMatrix<MATRIX, VECTOR>::TransposeMatrix (const char *name)
+template<class MATRIX, typename VectorType>
+TransposeMatrix<MATRIX, VectorType>::TransposeMatrix (const char *name)
: m(0, name)
{}
-template<class MATRIX, class VECTOR>
-TransposeMatrix<MATRIX, VECTOR>::TransposeMatrix (
+template<class MATRIX, typename VectorType>
+TransposeMatrix<MATRIX, VectorType>::TransposeMatrix (
const MATRIX *M,
const char *name)
: m(M, name)
{}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-TransposeMatrix<MATRIX, VECTOR>::clear ()
+TransposeMatrix<MATRIX, VectorType>::clear ()
{
m = 0;
}
-template<class MATRIX, class VECTOR>
-inline const TransposeMatrix<MATRIX, VECTOR> &
-TransposeMatrix<MATRIX, VECTOR>::operator= (const MATRIX *M)
+template<class MATRIX, typename VectorType>
+inline const TransposeMatrix<MATRIX, VectorType> &
+TransposeMatrix<MATRIX, VectorType>::operator= (const MATRIX *M)
{
m = M;
return *this;
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline bool
-TransposeMatrix<MATRIX, VECTOR>::empty () const
+TransposeMatrix<MATRIX, VectorType>::empty () const
{
if (m == 0)
return true;
return m->empty();
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-TransposeMatrix<MATRIX, VECTOR>::vmult (VECTOR &dst,
- const VECTOR &src) const
+TransposeMatrix<MATRIX, VectorType>::vmult (VectorType &dst,
+ const VectorType &src) const
{
Assert (m != 0, ExcNotInitialized());
m->Tvmult (dst, src);
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-TransposeMatrix<MATRIX, VECTOR>::Tvmult (VECTOR &dst,
- const VECTOR &src) const
+TransposeMatrix<MATRIX, VectorType>::Tvmult (VectorType &dst,
+ const VectorType &src) const
{
Assert (m != 0, ExcNotInitialized());
m->vmult (dst, src);
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-TransposeMatrix<MATRIX, VECTOR>::vmult_add (VECTOR &dst,
- const VECTOR &src) const
+TransposeMatrix<MATRIX, VectorType>::vmult_add (VectorType &dst,
+ const VectorType &src) const
{
Assert (m != 0, ExcNotInitialized());
m->Tvmult_add (dst, src);
}
-template<class MATRIX, class VECTOR>
+template<class MATRIX, typename VectorType>
inline void
-TransposeMatrix<MATRIX, VECTOR>::Tvmult_add (VECTOR &dst,
- const VECTOR &src) const
+TransposeMatrix<MATRIX, VectorType>::Tvmult_add (VectorType &dst,
+ const VectorType &src) const
{
Assert (m != 0, ExcNotInitialized());
m->vmult_add (dst, src);
DEAL_II_NAMESPACE_OPEN
-template <typename VECTOR>
-typename GrowingVectorMemory<VECTOR>::Pool GrowingVectorMemory<VECTOR>::pool;
+template <typename VectorType>
+typename GrowingVectorMemory<VectorType>::Pool GrowingVectorMemory<VectorType>::pool;
-template <typename VECTOR>
-Threads::Mutex GrowingVectorMemory<VECTOR>::mutex;
+template <typename VectorType>
+Threads::Mutex GrowingVectorMemory<VectorType>::mutex;
-template <typename VECTOR>
+template <typename VectorType>
inline
-GrowingVectorMemory<VECTOR>::Pool::Pool()
+GrowingVectorMemory<VectorType>::Pool::Pool()
:
data(0)
{}
-template <typename VECTOR>
+template <typename VectorType>
inline
-GrowingVectorMemory<VECTOR>::Pool::~Pool()
+GrowingVectorMemory<VectorType>::Pool::~Pool()
{
// Nothing to do if memory was unused.
if (data == 0) return;
}
-template <typename VECTOR>
+template <typename VectorType>
inline
void
-GrowingVectorMemory<VECTOR>::Pool::initialize(const size_type size)
+GrowingVectorMemory<VectorType>::Pool::initialize(const size_type size)
{
if (data == 0)
{
++i)
{
i->first = false;
- i->second = new VECTOR;
+ i->second = new VectorType;
}
}
}
-template <typename VECTOR>
+template <typename VectorType>
inline
-GrowingVectorMemory<VECTOR>::GrowingVectorMemory (const size_type initial_size,
- const bool log_statistics)
+GrowingVectorMemory<VectorType>::GrowingVectorMemory (const size_type initial_size,
+ const bool log_statistics)
:
total_alloc(0),
}
-template<typename VECTOR>
+template<typename VectorType>
inline
-GrowingVectorMemory<VECTOR>::~GrowingVectorMemory()
+GrowingVectorMemory<VectorType>::~GrowingVectorMemory()
{
AssertNothrow(current_alloc == 0,
StandardExceptions::ExcMemoryLeak(current_alloc));
-template<typename VECTOR>
+template<typename VectorType>
inline
-VECTOR *
-GrowingVectorMemory<VECTOR>::alloc ()
+VectorType *
+GrowingVectorMemory<VectorType>::alloc ()
{
Threads::Mutex::ScopedLock lock(mutex);
++total_alloc;
// no free vector found, so let's
// just allocate a new one
- const entry_type t (true, new VECTOR);
+ const entry_type t (true, new VectorType);
pool.data->push_back(t);
return t.second;
-template<typename VECTOR>
+template<typename VectorType>
inline
void
-GrowingVectorMemory<VECTOR>::free(const VECTOR *const v)
+GrowingVectorMemory<VectorType>::free(const VectorType *const v)
{
Threads::Mutex::ScopedLock lock(mutex);
for (typename std::vector<entry_type>::iterator i=pool.data->begin();
return;
}
}
- Assert(false, typename VectorMemory<VECTOR>::ExcNotAllocatedHere());
+ Assert(false, typename VectorMemory<VectorType>::ExcNotAllocatedHere());
}
-template<typename VECTOR>
+template<typename VectorType>
inline
void
-GrowingVectorMemory<VECTOR>::release_unused_memory ()
+GrowingVectorMemory<VectorType>::release_unused_memory ()
{
Threads::Mutex::ScopedLock lock(mutex);
-template<typename VECTOR>
+template<typename VectorType>
inline
std::size_t
-GrowingVectorMemory<VECTOR>::memory_consumption () const
+GrowingVectorMemory<VectorType>::memory_consumption () const
{
Threads::Mutex::ScopedLock lock(mutex);
*
* @author Guido Kanschat, 2002
*/
-template <class VECTOR>
+template <typename VectorType>
class MGMatrixBase : public Subscriptor
{
public:
* Matrix-vector-multiplication on a certain level.
*/
virtual void vmult (const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const = 0;
+ VectorType &dst,
+ const VectorType &src) const = 0;
/**
* Adding matrix-vector-multiplication on a certain level.
*/
virtual void vmult_add (const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const = 0;
+ VectorType &dst,
+ const VectorType &src) const = 0;
/**
* Transpose matrix-vector-multiplication on a certain level.
*/
virtual void Tvmult (const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const = 0;
+ VectorType &dst,
+ const VectorType &src) const = 0;
/**
* Adding transpose matrix-vector-multiplication on a certain level.
*/
virtual void Tvmult_add (const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const = 0;
+ VectorType &dst,
+ const VectorType &src) const = 0;
};
*
* @author Guido Kanschat, 2002
*/
-template <class VECTOR>
+template <typename VectorType>
class MGCoarseGridBase : public Subscriptor
{
public:
/**
* Solution operator.
*/
- virtual void operator() (const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const = 0;
+ virtual void operator() (const unsigned int level,
+ VectorType &dst,
+ const VectorType &src) const = 0;
};
*
* @author Wolfgang Bangerth, Guido Kanschat, 1999, 2002, 2007
*/
-template <class VECTOR>
+template <typename VectorType>
class MGTransferBase : public Subscriptor
{
public:
* finer level.
*/
virtual void prolongate (const unsigned int to_level,
- VECTOR &dst,
- const VECTOR &src) const = 0;
+ VectorType &dst,
+ const VectorType &src) const = 0;
/**
* Restrict a vector from level <tt>from_level</tt> to level
*
*/
virtual void restrict_and_add (const unsigned int from_level,
- VECTOR &dst,
- const VECTOR &src) const = 0;
+ VectorType &dst,
+ const VectorType &src) const = 0;
};
*
* @author Guido Kanschat, 2002
*/
-template <class VECTOR>
+template <typename VectorType>
class MGSmootherBase : public Subscriptor
{
public:
* Smoothing function. This is the function used in multigrid methods.
*/
virtual void smooth (const unsigned int level,
- VECTOR &u,
- const VECTOR &rhs) const = 0;
+ VectorType &u,
+ const VectorType &rhs) const = 0;
};
/*@}*/
*
* @author Guido Kanschat, 1999, Ralf Hartmann, 2002.
*/
-template<class SOLVER, class VECTOR = Vector<double> >
-class MGCoarseGridLACIteration : public MGCoarseGridBase<VECTOR>
+template<class SOLVER, class VectorType = Vector<double> >
+class MGCoarseGridLACIteration : public MGCoarseGridBase<VectorType>
{
public:
/**
* Implementation of the abstract function. Calls the solver method with
* matrix, vectors and preconditioner.
*/
- void operator() (const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const;
+ void operator() (const unsigned int level,
+ VectorType &dst,
+ const VectorType &src) const;
/**
* Sets the matrix. This gives the possibility to replace the matrix that
/**
* Reference to the solver.
*/
- SmartPointer<SOLVER,MGCoarseGridLACIteration<SOLVER,VECTOR> > solver;
+ SmartPointer<SOLVER,MGCoarseGridLACIteration<SOLVER,VectorType> > solver;
/**
* Reference to the matrix.
*/
- PointerMatrixBase<VECTOR> *matrix;
+ PointerMatrixBase<VectorType> *matrix;
/**
* Reference to the preconditioner.
*/
- PointerMatrixBase<VECTOR> *precondition;
+ PointerMatrixBase<VectorType> *precondition;
};
*
* @author Guido Kanschat, 2003, 2012
*/
-template<typename number = double, class VECTOR = Vector<number> >
-class MGCoarseGridHouseholder : public MGCoarseGridBase<VECTOR>
+template<typename number = double, class VectorType = Vector<number> >
+class MGCoarseGridHouseholder : public MGCoarseGridBase<VectorType>
{
public:
/**
*/
void initialize (const FullMatrix<number> &A);
- void operator() (const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const;
+ void operator() (const unsigned int level,
+ VectorType &dst,
+ const VectorType &src) const;
private:
/**
*
* @author Guido Kanschat, 2003, 2012
*/
-template<typename number = double, class VECTOR = Vector<number> >
-class MGCoarseGridSVD : public MGCoarseGridBase<VECTOR>
+template<typename number = double, class VectorType = Vector<number> >
+class MGCoarseGridSVD : public MGCoarseGridBase<VectorType>
{
public:
/**
*/
void initialize (const FullMatrix<number> &A, const double threshold = 0);
- void operator() (const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const;
+ void operator() (const unsigned int level,
+ VectorType &dst,
+ const VectorType &src) const;
/**
* Write the singular values to @p deallog.
/* ------------------ Functions for MGCoarseGridLACIteration ------------ */
-template<class SOLVER, class VECTOR>
-MGCoarseGridLACIteration<SOLVER, VECTOR>
+template<class SOLVER, class VectorType>
+MGCoarseGridLACIteration<SOLVER, VectorType>
::MGCoarseGridLACIteration()
:
solver(0, typeid(*this).name()),
{}
-template<class SOLVER, class VECTOR>
+template<class SOLVER, class VectorType>
template<class MATRIX, class PRECOND>
-MGCoarseGridLACIteration<SOLVER, VECTOR>
+MGCoarseGridLACIteration<SOLVER, VectorType>
::MGCoarseGridLACIteration(SOLVER &s,
const MATRIX &m,
const PRECOND &p)
:
solver(&s, typeid(*this).name())
{
- matrix = new PointerMatrix<MATRIX, VECTOR>(&m);
- precondition = new PointerMatrix<PRECOND, VECTOR>(&p);
+ matrix = new PointerMatrix<MATRIX, VectorType>(&m);
+ precondition = new PointerMatrix<PRECOND, VectorType>(&p);
}
-template<class SOLVER, class VECTOR>
-MGCoarseGridLACIteration<SOLVER, VECTOR>
+template<class SOLVER, class VectorType>
+MGCoarseGridLACIteration<SOLVER, VectorType>
::~MGCoarseGridLACIteration()
{
clear();
}
-template<class SOLVER, class VECTOR>
+template<class SOLVER, class VectorType>
template<class MATRIX, class PRECOND>
void
-MGCoarseGridLACIteration<SOLVER, VECTOR>
+MGCoarseGridLACIteration<SOLVER, VectorType>
::initialize(SOLVER &s,
const MATRIX &m,
const PRECOND &p)
solver = &s;
if (matrix)
delete matrix;
- matrix = new PointerMatrix<MATRIX, VECTOR>(&m);
+ matrix = new PointerMatrix<MATRIX, VectorType>(&m);
if (precondition)
delete precondition;
- precondition = new PointerMatrix<PRECOND, VECTOR>(&p);
+ precondition = new PointerMatrix<PRECOND, VectorType>(&p);
}
-template<class SOLVER, class VECTOR>
+template<class SOLVER, class VectorType>
void
-MGCoarseGridLACIteration<SOLVER, VECTOR>
+MGCoarseGridLACIteration<SOLVER, VectorType>
::clear()
{
solver = 0;
}
-template<class SOLVER, class VECTOR>
+template<class SOLVER, class VectorType>
void
-MGCoarseGridLACIteration<SOLVER, VECTOR>
-::operator() (const unsigned int /* level */,
- VECTOR &dst,
- const VECTOR &src) const
+MGCoarseGridLACIteration<SOLVER, VectorType>
+::operator() (const unsigned int /* level */,
+ VectorType &dst,
+ const VectorType &src) const
{
Assert(solver!=0, ExcNotInitialized());
Assert(matrix!=0, ExcNotInitialized());
}
-template<class SOLVER, class VECTOR>
+template<class SOLVER, class VectorType>
template<class MATRIX>
void
-MGCoarseGridLACIteration<SOLVER, VECTOR>
+MGCoarseGridLACIteration<SOLVER, VectorType>
::set_matrix(const MATRIX &m)
{
if (matrix)
delete matrix;
- matrix = new PointerMatrix<MATRIX, VECTOR>(&m);
+ matrix = new PointerMatrix<MATRIX, VectorType>(&m);
}
//---------------------------------------------------------------------------
-template<typename number, class VECTOR>
-MGCoarseGridHouseholder<number, VECTOR>::MGCoarseGridHouseholder(
+template<typename number, class VectorType>
+MGCoarseGridHouseholder<number, VectorType>::MGCoarseGridHouseholder(
const FullMatrix<number> *A)
{
if (A != 0) householder.initialize(*A);
-template<typename number, class VECTOR>
+template<typename number, class VectorType>
void
-MGCoarseGridHouseholder<number, VECTOR>::initialize(
+MGCoarseGridHouseholder<number, VectorType>::initialize(
const FullMatrix<number> &A)
{
householder.initialize(A);
-template<typename number, class VECTOR>
+template<typename number, class VectorType>
void
-MGCoarseGridHouseholder<number, VECTOR>::operator() (
+MGCoarseGridHouseholder<number, VectorType>::operator() (
const unsigned int /*level*/,
- VECTOR &dst,
- const VECTOR &src) const
+ VectorType &dst,
+ const VectorType &src) const
{
householder.least_squares(dst, src);
}
//---------------------------------------------------------------------------
-template<typename number, class VECTOR>
+template<typename number, class VectorType>
inline
-MGCoarseGridSVD<number, VECTOR>::MGCoarseGridSVD()
+MGCoarseGridSVD<number, VectorType>::MGCoarseGridSVD()
{}
-template<typename number, class VECTOR>
+template<typename number, class VectorType>
void
-MGCoarseGridSVD<number, VECTOR>::initialize(
+MGCoarseGridSVD<number, VectorType>::initialize(
const FullMatrix<number> &A,
double threshold)
{
}
-template<typename number, class VECTOR>
+template<typename number, class VectorType>
void
-MGCoarseGridSVD<number, VECTOR>::operator() (
+MGCoarseGridSVD<number, VectorType>::operator() (
const unsigned int /*level*/,
- VECTOR &dst,
- const VECTOR &src) const
+ VectorType &dst,
+ const VectorType &src) const
{
matrix.vmult(dst, src);
}
-template<typename number, class VECTOR>
+template<typename number, class VectorType>
void
-MGCoarseGridSVD<number, VECTOR>::log() const
+MGCoarseGridSVD<number, VectorType>::log() const
{
const unsigned int n = std::min(matrix.n_rows(), matrix.n_cols());
* @author Guido Kanschat
* @date 2002, 2010
*/
- template <class VECTOR = Vector<double> >
+ template <typename VectorType = Vector<double> >
class Matrix
- : public MGMatrixBase<VECTOR>
+ : public MGMatrixBase<VectorType>
{
public:
/**
/**
* Access matrix on a level.
*/
- const PointerMatrixBase<VECTOR> &operator[] (unsigned int level) const;
+ const PointerMatrixBase<VectorType> &operator[] (unsigned int level) const;
- virtual void vmult (const unsigned int level, VECTOR &dst, const VECTOR &src) const;
- virtual void vmult_add (const unsigned int level, VECTOR &dst, const VECTOR &src) const;
- virtual void Tvmult (const unsigned int level, VECTOR &dst, const VECTOR &src) const;
- virtual void Tvmult_add (const unsigned int level, VECTOR &dst, const VECTOR &src) const;
+ virtual void vmult (const unsigned int level, VectorType &dst, const VectorType &src) const;
+ virtual void vmult_add (const unsigned int level, VectorType &dst, const VectorType &src) const;
+ virtual void Tvmult (const unsigned int level, VectorType &dst, const VectorType &src) const;
+ virtual void Tvmult_add (const unsigned int level, VectorType &dst, const VectorType &src) const;
/**
* Memory used by this object.
*/
std::size_t memory_consumption () const;
private:
- MGLevelObject<std_cxx11::shared_ptr<PointerMatrixBase<VECTOR> > > matrices;
+ MGLevelObject<std_cxx11::shared_ptr<PointerMatrixBase<VectorType> > > matrices;
};
}
namespace mg
{
- template <class VECTOR>
+ template <typename VectorType>
template <class MATRIX>
inline
void
- Matrix<VECTOR>::initialize (const MGLevelObject<MATRIX> &p)
+ Matrix<VectorType>::initialize (const MGLevelObject<MATRIX> &p)
{
matrices.resize(p.min_level(), p.max_level());
for (unsigned int level=p.min_level(); level <= p.max_level(); ++level)
- matrices[level] = std_cxx11::shared_ptr<PointerMatrixBase<VECTOR> > (new_pointer_matrix_base(p[level], VECTOR()));
+ matrices[level] = std_cxx11::shared_ptr<PointerMatrixBase<VectorType> >
+ (new_pointer_matrix_base(p[level], VectorType()));
}
- template <class VECTOR>
+ template <typename VectorType>
template <class MATRIX>
inline
- Matrix<VECTOR>::Matrix (const MGLevelObject<MATRIX> &p)
+ Matrix<VectorType>::Matrix (const MGLevelObject<MATRIX> &p)
{
initialize(p);
}
- template <class VECTOR>
+ template <typename VectorType>
inline
- Matrix<VECTOR>::Matrix ()
+ Matrix<VectorType>::Matrix ()
{}
- template <class VECTOR>
+ template <typename VectorType>
inline
- const PointerMatrixBase<VECTOR> &
- Matrix<VECTOR>::operator[] (unsigned int level) const
+ const PointerMatrixBase<VectorType> &
+ Matrix<VectorType>::operator[] (unsigned int level) const
{
return *matrices[level];
}
- template <class VECTOR>
+ template <typename VectorType>
void
- Matrix<VECTOR>::vmult (
- const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const
+ Matrix<VectorType>::vmult (const unsigned int level,
+ VectorType &dst,
+ const VectorType &src) const
{
matrices[level]->vmult(dst, src);
}
- template <class VECTOR>
+ template <typename VectorType>
void
- Matrix<VECTOR>::vmult_add (
- const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const
+ Matrix<VectorType>::vmult_add (const unsigned int level,
+ VectorType &dst,
+ const VectorType &src) const
{
matrices[level]->vmult_add(dst, src);
}
- template <class VECTOR>
+ template <typename VectorType>
void
- Matrix<VECTOR>::Tvmult (const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const
+ Matrix<VectorType>::Tvmult (const unsigned int level,
+ VectorType &dst,
+ const VectorType &src) const
{
matrices[level]->Tvmult(dst, src);
}
- template <class VECTOR>
+ template <typename VectorType>
void
- Matrix<VECTOR>::Tvmult_add (const unsigned int level,
- VECTOR &dst,
- const VECTOR &src) const
+ Matrix<VectorType>::Tvmult_add (const unsigned int level,
+ VectorType &dst,
+ const VectorType &src) const
{
matrices[level]->Tvmult_add(dst, src);
}
- template <class VECTOR>
+ template <typename VectorType>
inline
std::size_t
- Matrix<VECTOR>::memory_consumption () const
+ Matrix<VectorType>::memory_consumption () const
{
return sizeof(*this) + matrices->memory_consumption();
}
*
* @author Guido Kanschat 2009
*/
-template <class VECTOR>
-class MGSmoother : public MGSmootherBase<VECTOR>
+template <typename VectorType>
+class MGSmoother : public MGSmootherBase<VectorType>
{
public:
/**
* The object is marked as mutable since we will need to use it to allocate
* temporary vectors also in functions that are const.
*/
- mutable GrowingVectorMemory<VECTOR> vector_memory;
+ mutable GrowingVectorMemory<VectorType> vector_memory;
/**
* Number of smoothing steps on the finest level. If no #variable smoothing
*
* @author Guido Kanschat, 1999, 2002
*/
-template <class VECTOR>
-class MGSmootherIdentity : public MGSmootherBase<VECTOR>
+template <typename VectorType>
+class MGSmootherIdentity : public MGSmootherBase<VectorType>
{
public:
/**
* that the the smoothing operator equals the null operator.
*/
virtual void smooth (const unsigned int level,
- VECTOR &u,
- const VECTOR &rhs) const;
+ VectorType &u,
+ const VectorType &rhs) const;
virtual void clear ();
};
*
* A relaxation class is an object that has two member functions,
* @code
- * void step(VECTOR& x, const VECTOR& d) const;
- * void Tstep(VECTOR& x, const VECTOR& d) const;
+ * void step(VectorType& x, const VectorType& d) const;
+ * void Tstep(VectorType& x, const VectorType& d) const;
* @endcode
* performing one step of the smoothing scheme.
*
* @author Guido Kanschat,
* @date 2003, 2009, 2010
*/
- template<class RELAX, class VECTOR>
- class SmootherRelaxation : public MGLevelObject<RELAX>, public MGSmoother<VECTOR>
+ template<class RELAX, typename VectorType>
+ class SmootherRelaxation : public MGLevelObject<RELAX>, public MGSmoother<VectorType>
{
public:
/**
* The actual smoothing method.
*/
virtual void smooth (const unsigned int level,
- VECTOR &u,
- const VECTOR &rhs) const;
+ VectorType &u,
+ const VectorType &rhs) const;
/**
* Memory used by this object.
*
* A relaxation class is an object that has two member functions,
* @code
- * void step(VECTOR& x, const VECTOR& d) const;
- * void Tstep(VECTOR& x, const VECTOR& d) const;
+ * void step(VectorType& x, const VectorType& d) const;
+ * void Tstep(VectorType& x, const VectorType& d) const;
* @endcode
* performing one step of the smoothing scheme.
*
*
* @author Guido Kanschat, 2003
*/
-template<class MATRIX, class RELAX, class VECTOR>
-class MGSmootherRelaxation : public MGSmoother<VECTOR>
+template<class MATRIX, class RELAX, typename VectorType>
+class MGSmootherRelaxation : public MGSmoother<VectorType>
{
public:
/**
* The actual smoothing method.
*/
virtual void smooth (const unsigned int level,
- VECTOR &u,
- const VECTOR &rhs) const;
+ VectorType &u,
+ const VectorType &rhs) const;
/**
* Object containing relaxation methods.
/**
* Pointer to the matrices.
*/
- MGLevelObject<PointerMatrix<MATRIX, VECTOR> > matrices;
+ MGLevelObject<PointerMatrix<MATRIX, VectorType> > matrices;
};
*
* @author Guido Kanschat, 2009
*/
-template<class MATRIX, class PRECONDITIONER, class VECTOR>
-class MGSmootherPrecondition : public MGSmoother<VECTOR>
+template<class MATRIX, class PRECONDITIONER, typename VectorType>
+class MGSmootherPrecondition : public MGSmoother<VectorType>
{
public:
/**
* The actual smoothing method.
*/
virtual void smooth (const unsigned int level,
- VECTOR &u,
- const VECTOR &rhs) const;
+ VectorType &u,
+ const VectorType &rhs) const;
/**
* Object containing relaxation methods.
/**
* Pointer to the matrices.
*/
- MGLevelObject<PointerMatrix<MATRIX, VECTOR> > matrices;
+ MGLevelObject<PointerMatrix<MATRIX, VectorType> > matrices;
};
#ifndef DOXYGEN
-template <class VECTOR>
+template <typename VectorType>
inline void
-MGSmootherIdentity<VECTOR>::smooth (
- const unsigned int, VECTOR &,
- const VECTOR &) const
+MGSmootherIdentity<VectorType>::smooth (const unsigned int,
+ VectorType &,
+ const VectorType &) const
{}
-template <class VECTOR>
+template <typename VectorType>
inline void
-MGSmootherIdentity<VECTOR>::clear ()
+MGSmootherIdentity<VectorType>::clear ()
{}
//---------------------------------------------------------------------------
-template <class VECTOR>
+template <typename VectorType>
inline
-MGSmoother<VECTOR>::MGSmoother(
+MGSmoother<VectorType>::MGSmoother(
const unsigned int steps,
const bool variable,
const bool symmetric,
{}
-template <class VECTOR>
+template <typename VectorType>
inline void
-MGSmoother<VECTOR>::set_steps (const unsigned int s)
+MGSmoother<VectorType>::set_steps (const unsigned int s)
{
steps = s;
}
-template <class VECTOR>
+template <typename VectorType>
inline void
-MGSmoother<VECTOR>::set_debug (const unsigned int s)
+MGSmoother<VectorType>::set_debug (const unsigned int s)
{
debug = s;
}
-template <class VECTOR>
+template <typename VectorType>
inline void
-MGSmoother<VECTOR>::set_variable (const bool flag)
+MGSmoother<VectorType>::set_variable (const bool flag)
{
variable = flag;
}
-template <class VECTOR>
+template <typename VectorType>
inline void
-MGSmoother<VECTOR>::set_symmetric (const bool flag)
+MGSmoother<VectorType>::set_symmetric (const bool flag)
{
symmetric = flag;
}
-template <class VECTOR>
+template <typename VectorType>
inline void
-MGSmoother<VECTOR>::set_transpose (const bool flag)
+MGSmoother<VectorType>::set_transpose (const bool flag)
{
transpose = flag;
}
namespace mg
{
- template <class RELAX, class VECTOR>
+ template <class RELAX, typename VectorType>
inline
- SmootherRelaxation<RELAX, VECTOR>::SmootherRelaxation(
+ SmootherRelaxation<RELAX, VectorType>::SmootherRelaxation(
const unsigned int steps,
const bool variable,
const bool symmetric,
const bool transpose)
- : MGSmoother<VECTOR>(steps, variable, symmetric, transpose)
+ : MGSmoother<VectorType>(steps, variable, symmetric, transpose)
{}
- template <class RELAX, class VECTOR>
+ template <class RELAX, typename VectorType>
inline void
- SmootherRelaxation<RELAX, VECTOR>::clear ()
+ SmootherRelaxation<RELAX, VectorType>::clear ()
{
MGLevelObject<RELAX>::clear();
}
- template <class RELAX, class VECTOR>
+ template <class RELAX, typename VectorType>
template <class MATRIX2>
inline void
- SmootherRelaxation<RELAX, VECTOR>::initialize (
+ SmootherRelaxation<RELAX, VectorType>::initialize (
const MGLevelObject<MATRIX2> &m,
const typename RELAX::AdditionalData &data)
{
}
- template <class RELAX, class VECTOR>
+ template <class RELAX, typename VectorType>
template <class MATRIX2, class DATA>
inline void
- SmootherRelaxation<RELAX, VECTOR>::initialize (
+ SmootherRelaxation<RELAX, VectorType>::initialize (
const MGLevelObject<MATRIX2> &m,
const MGLevelObject<DATA> &data)
{
}
- template <class RELAX, class VECTOR>
+ template <class RELAX, typename VectorType>
inline void
- SmootherRelaxation<RELAX, VECTOR>::smooth(
- const unsigned int level,
- VECTOR &u,
- const VECTOR &rhs) const
+ SmootherRelaxation<RELAX, VectorType>::smooth (const unsigned int level,
+ VectorType &u,
+ const VectorType &rhs) const
{
unsigned int maxlevel = this->max_level();
unsigned int steps2 = this->steps;
}
- template <class RELAX, class VECTOR>
+ template <class RELAX, typename VectorType>
inline
std::size_t
- SmootherRelaxation<RELAX, VECTOR>::
+ SmootherRelaxation<RELAX, VectorType>::
memory_consumption () const
{
return sizeof(*this)
//----------------------------------------------------------------------//
-template <class MATRIX, class RELAX, class VECTOR>
+template <class MATRIX, class RELAX, typename VectorType>
inline
-MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::MGSmootherRelaxation(
+MGSmootherRelaxation<MATRIX, RELAX, VectorType>::MGSmootherRelaxation(
const unsigned int steps,
const bool variable,
const bool symmetric,
const bool transpose)
:
- MGSmoother<VECTOR>(steps, variable, symmetric, transpose)
+ MGSmoother<VectorType>(steps, variable, symmetric, transpose)
{}
-template <class MATRIX, class RELAX, class VECTOR>
+template <class MATRIX, class RELAX, typename VectorType>
inline void
-MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::clear ()
+MGSmootherRelaxation<MATRIX, RELAX, VectorType>::clear ()
{
smoothers.clear();
}
-template <class MATRIX, class RELAX, class VECTOR>
+template <class MATRIX, class RELAX, typename VectorType>
template <class MATRIX2>
inline void
-MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::initialize (
+MGSmootherRelaxation<MATRIX, RELAX, VectorType>::initialize (
const MGLevelObject<MATRIX2> &m,
const typename RELAX::AdditionalData &data)
{
}
}
-template <class MATRIX, class RELAX, class VECTOR>
+template <class MATRIX, class RELAX, typename VectorType>
template <class MATRIX2, class DATA>
inline void
-MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::initialize (
+MGSmootherRelaxation<MATRIX, RELAX, VectorType>::initialize (
const MGLevelObject<MATRIX2> &m,
const MGLevelObject<DATA> &data)
{
}
}
-template <class MATRIX, class RELAX, class VECTOR>
+template <class MATRIX, class RELAX, typename VectorType>
template <class MATRIX2, class DATA>
inline void
-MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::initialize (
+MGSmootherRelaxation<MATRIX, RELAX, VectorType>::initialize (
const MGLevelObject<MATRIX2> &m,
const DATA &data,
const unsigned int row,
}
}
-template <class MATRIX, class RELAX, class VECTOR>
+template <class MATRIX, class RELAX, typename VectorType>
template <class MATRIX2, class DATA>
inline void
-MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::initialize (
+MGSmootherRelaxation<MATRIX, RELAX, VectorType>::initialize (
const MGLevelObject<MATRIX2> &m,
const MGLevelObject<DATA> &data,
const unsigned int row,
}
-template <class MATRIX, class RELAX, class VECTOR>
+template <class MATRIX, class RELAX, typename VectorType>
inline void
-MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::smooth(
- const unsigned int level,
- VECTOR &u,
- const VECTOR &rhs) const
+MGSmootherRelaxation<MATRIX, RELAX, VectorType>::smooth (const unsigned int level,
+ VectorType &u,
+ const VectorType &rhs) const
{
unsigned int maxlevel = smoothers.max_level();
unsigned int steps2 = this->steps;
-template <class MATRIX, class RELAX, class VECTOR>
+template <class MATRIX, class RELAX, typename VectorType>
inline
std::size_t
-MGSmootherRelaxation<MATRIX, RELAX, VECTOR>::
+MGSmootherRelaxation<MATRIX, RELAX, VectorType>::
memory_consumption () const
{
return sizeof(*this)
//----------------------------------------------------------------------//
-template <class MATRIX, class PRECONDITIONER, class VECTOR>
+template <class MATRIX, class PRECONDITIONER, typename VectorType>
inline
-MGSmootherPrecondition<MATRIX, PRECONDITIONER, VECTOR>::MGSmootherPrecondition(
+MGSmootherPrecondition<MATRIX, PRECONDITIONER, VectorType>::MGSmootherPrecondition(
const unsigned int steps,
const bool variable,
const bool symmetric,
const bool transpose)
:
- MGSmoother<VECTOR>(steps, variable, symmetric, transpose)
+ MGSmoother<VectorType>(steps, variable, symmetric, transpose)
{}
-template <class MATRIX, class PRECONDITIONER, class VECTOR>
+template <class MATRIX, class PRECONDITIONER, typename VectorType>
inline void
-MGSmootherPrecondition<MATRIX, PRECONDITIONER, VECTOR>::clear ()
+MGSmootherPrecondition<MATRIX, PRECONDITIONER, VectorType>::clear ()
{
smoothers.clear();
-template <class MATRIX, class PRECONDITIONER, class VECTOR>
+template <class MATRIX, class PRECONDITIONER, typename VectorType>
template <class MATRIX2>
inline void
-MGSmootherPrecondition<MATRIX, PRECONDITIONER, VECTOR>::initialize (
+MGSmootherPrecondition<MATRIX, PRECONDITIONER, VectorType>::initialize (
const MGLevelObject<MATRIX2> &m,
const typename PRECONDITIONER::AdditionalData &data)
{
-template <class MATRIX, class PRECONDITIONER, class VECTOR>
+template <class MATRIX, class PRECONDITIONER, typename VectorType>
template <class MATRIX2, class DATA>
inline void
-MGSmootherPrecondition<MATRIX, PRECONDITIONER, VECTOR>::initialize (
+MGSmootherPrecondition<MATRIX, PRECONDITIONER, VectorType>::initialize (
const MGLevelObject<MATRIX2> &m,
const MGLevelObject<DATA> &data)
{
-template <class MATRIX, class PRECONDITIONER, class VECTOR>
+template <class MATRIX, class PRECONDITIONER, typename VectorType>
template <class MATRIX2, class DATA>
inline void
-MGSmootherPrecondition<MATRIX, PRECONDITIONER, VECTOR>::initialize (
+MGSmootherPrecondition<MATRIX, PRECONDITIONER, VectorType>::initialize (
const MGLevelObject<MATRIX2> &m,
const DATA &data,
const unsigned int row,
-template <class MATRIX, class PRECONDITIONER, class VECTOR>
+template <class MATRIX, class PRECONDITIONER, typename VectorType>
template <class MATRIX2, class DATA>
inline void
-MGSmootherPrecondition<MATRIX, PRECONDITIONER, VECTOR>::initialize (
+MGSmootherPrecondition<MATRIX, PRECONDITIONER, VectorType>::initialize (
const MGLevelObject<MATRIX2> &m,
const MGLevelObject<DATA> &data,
const unsigned int row,
-template <class MATRIX, class PRECONDITIONER, class VECTOR>
+template <class MATRIX, class PRECONDITIONER, typename VectorType>
inline void
-MGSmootherPrecondition<MATRIX, PRECONDITIONER, VECTOR>::smooth(
- const unsigned int level,
- VECTOR &u,
- const VECTOR &rhs) const
+MGSmootherPrecondition<MATRIX, PRECONDITIONER, VectorType>::smooth
+(const unsigned int level,
+ VectorType &u,
+ const VectorType &rhs) const
{
unsigned int maxlevel = matrices.max_level();
unsigned int steps2 = this->steps;
if (this->variable)
steps2 *= (1<<(maxlevel-level));
- typename VectorMemory<VECTOR>::Pointer r(this->vector_memory);
- typename VectorMemory<VECTOR>::Pointer d(this->vector_memory);
+ typename VectorMemory<VectorType>::Pointer r(this->vector_memory);
+ typename VectorMemory<VectorType>::Pointer d(this->vector_memory);
r->reinit(u,true);
d->reinit(u,true);
-template <class MATRIX, class PRECONDITIONER, class VECTOR>
+template <class MATRIX, class PRECONDITIONER, typename VectorType>
inline
std::size_t
-MGSmootherPrecondition<MATRIX, PRECONDITIONER, VECTOR>::
+MGSmootherPrecondition<MATRIX, PRECONDITIONER, VectorType>::
memory_consumption () const
{
return sizeof(*this)
namespace internal
{
- template <class VECTOR>
+ template <typename VectorType>
struct MatrixSelector
{
typedef ::dealii::SparsityPattern Sparsity;
- typedef ::dealii::SparseMatrix<typename VECTOR::value_type> Matrix;
+ typedef ::dealii::SparseMatrix<typename VectorType::value_type> Matrix;
template <class DSP, class DH>
static void reinit(Matrix &matrix, Sparsity &sparsity, int level, const DSP &dsp, const DH &)
* @author Wolfgang Bangerth, Guido Kanschat
* @date 1999, 2000, 2001, 2002, 2003, 2004, 2012
*/
-template <class VECTOR>
-class MGTransferPrebuilt : public MGTransferBase<VECTOR>
+template <typename VectorType>
+class MGTransferPrebuilt : public MGTransferBase<VectorType>
{
public:
/**
template <int dim, int spacedim>
void build_matrices (const DoFHandler<dim,spacedim> &mg_dof);
- virtual void prolongate (const unsigned int to_level,
- VECTOR &dst,
- const VECTOR &src) const;
+ virtual void prolongate (const unsigned int to_level,
+ VectorType &dst,
+ const VectorType &src) const;
- virtual void restrict_and_add (const unsigned int from_level,
- VECTOR &dst,
- const VECTOR &src) const;
+ virtual void restrict_and_add (const unsigned int from_level,
+ VectorType &dst,
+ const VectorType &src) const;
/**
* Transfer from a vector on the global grid to vectors defined on each of
template <int dim, class InVector, int spacedim>
void
copy_to_mg (const DoFHandler<dim,spacedim> &mg_dof,
- MGLevelObject<VECTOR> &dst,
- const InVector &src) const;
+ MGLevelObject<VectorType> &dst,
+ const InVector &src) const;
/**
* Transfer from multi-level vector to normal vector.
*/
template <int dim, class OutVector, int spacedim>
void
- copy_from_mg (const DoFHandler<dim,spacedim> &mg_dof,
- OutVector &dst,
- const MGLevelObject<VECTOR> &src) const;
+ copy_from_mg (const DoFHandler<dim,spacedim> &mg_dof,
+ OutVector &dst,
+ const MGLevelObject<VectorType> &src) const;
/**
* Add a multi-level vector to a normal vector.
*/
template <int dim, class OutVector, int spacedim>
void
- copy_from_mg_add (const DoFHandler<dim,spacedim> &mg_dof,
- OutVector &dst,
- const MGLevelObject<VECTOR> &src) const;
+ copy_from_mg_add (const DoFHandler<dim,spacedim> &mg_dof,
+ OutVector &dst,
+ const MGLevelObject<VectorType> &src) const;
/**
* If this object operates on BlockVector objects, we need to describe how
/**
* Sparsity patterns for transfer matrices.
*/
- std::vector<std_cxx11::shared_ptr<typename internal::MatrixSelector<VECTOR>::Sparsity> > prolongation_sparsities;
+ std::vector<std_cxx11::shared_ptr<typename internal::MatrixSelector<VectorType>::Sparsity> > prolongation_sparsities;
/**
* The actual prolongation matrix. column indices belong to the dof indices
* of the mother cell, i.e. the coarse level. while row indices belong to
* the child cell, i.e. the fine level.
*/
- std::vector<std_cxx11::shared_ptr<typename internal::MatrixSelector<VECTOR>::Matrix> > prolongation_matrices;
+ std::vector<std_cxx11::shared_ptr<typename internal::MatrixSelector<VectorType>::Matrix> > prolongation_matrices;
/**
* Mapping for the copy_to_mg() and copy_from_mg() functions. Here only
/**
* The constraints of the global system.
*/
- SmartPointer<const ConstraintMatrix, MGTransferPrebuilt<VECTOR> > constraints;
+ SmartPointer<const ConstraintMatrix, MGTransferPrebuilt<VectorType> > constraints;
/**
* The mg_constrained_dofs of the level systems.
*/
- SmartPointer<const MGConstrainedDoFs, MGTransferPrebuilt<VECTOR> > mg_constrained_dofs;
+ SmartPointer<const MGConstrainedDoFs, MGTransferPrebuilt<VectorType> > mg_constrained_dofs;
};
-template <class VECTOR>
+template <typename VectorType>
template <int dim, class InVector, int spacedim>
void
-MGTransferPrebuilt<VECTOR>::copy_to_mg (
- const DoFHandler<dim,spacedim> &mg_dof_handler,
- MGLevelObject<VECTOR> &dst,
- const InVector &src) const
+MGTransferPrebuilt<VectorType>::copy_to_mg
+(const DoFHandler<dim,spacedim> &mg_dof_handler,
+ MGLevelObject<VectorType> &dst,
+ const InVector &src) const
{
reinit_vector(mg_dof_handler, component_to_block_map, dst);
bool first = true;
for (unsigned int level=mg_dof_handler.get_tria().n_global_levels(); level != 0;)
{
--level;
- VECTOR &dst_level = dst[level];
+ VectorType &dst_level = dst[level];
#ifdef DEBUG_OUTPUT
MPI_Barrier(MPI_COMM_WORLD);
-template <class VECTOR>
+template <typename VectorType>
template <int dim, class OutVector, int spacedim>
void
-MGTransferPrebuilt<VECTOR>::copy_from_mg(
- const DoFHandler<dim,spacedim> &mg_dof_handler,
- OutVector &dst,
- const MGLevelObject<VECTOR> &src) const
+MGTransferPrebuilt<VectorType>::copy_from_mg
+(const DoFHandler<dim,spacedim> &mg_dof_handler,
+ OutVector &dst,
+ const MGLevelObject<VectorType> &src) const
{
// For non-DG: degrees of
// freedom in the refinement
-template <class VECTOR>
+template <typename VectorType>
template <int dim, class OutVector, int spacedim>
void
-MGTransferPrebuilt<VECTOR>::copy_from_mg_add (
- const DoFHandler<dim,spacedim> &mg_dof_handler,
- OutVector &dst,
- const MGLevelObject<VECTOR> &src) const
+MGTransferPrebuilt<VectorType>::copy_from_mg_add
+(const DoFHandler<dim,spacedim> &mg_dof_handler,
+ OutVector &dst,
+ const MGLevelObject<VectorType> &src) const
{
// For non-DG: degrees of
// freedom in the refinement
-template <class VECTOR>
+template <typename VectorType>
void
-MGTransferPrebuilt<VECTOR>::
+MGTransferPrebuilt<VectorType>::
set_component_to_block_map (const std::vector<unsigned int> &map)
{
component_to_block_map = map;
}
-template <class VECTOR>
+template <typename VectorType>
std::size_t
-MGTransferPrebuilt<VECTOR>::memory_consumption () const
+MGTransferPrebuilt<VectorType>::memory_consumption () const
{
std::size_t result = sizeof(*this);
result += sizeof(unsigned int) * sizes.size();
*
* @author Guido Kanschat, 1999 - 2005
*/
-template <class VECTOR>
+template <typename VectorType>
class Multigrid : public Subscriptor
{
public:
f_cycle
};
- typedef VECTOR vector_type;
- typedef const VECTOR const_vector_type;
+ typedef VectorType vector_type;
+ typedef const VectorType const_vector_type;
/**
* Constructor. The DoFHandler is used to determine the highest possible
* this type as late as possible.
*/
template <int dim>
- Multigrid(const DoFHandler<dim> &mg_dof_handler,
- const MGMatrixBase<VECTOR> &matrix,
- const MGCoarseGridBase<VECTOR> &coarse,
- const MGTransferBase<VECTOR> &transfer,
- const MGSmootherBase<VECTOR> &pre_smooth,
- const MGSmootherBase<VECTOR> &post_smooth,
- Cycle cycle = v_cycle);
+ Multigrid(const DoFHandler<dim> &mg_dof_handler,
+ const MGMatrixBase<VectorType> &matrix,
+ const MGCoarseGridBase<VectorType> &coarse,
+ const MGTransferBase<VectorType> &transfer,
+ const MGSmootherBase<VectorType> &pre_smooth,
+ const MGSmootherBase<VectorType> &post_smooth,
+ Cycle cycle = v_cycle);
/**
* Experimental constructor for cases in which no DoFHandler is available.
*
* @warning Not intended for general use.
*/
- Multigrid(const unsigned int minlevel,
- const unsigned int maxlevel,
- const MGMatrixBase<VECTOR> &matrix,
- const MGCoarseGridBase<VECTOR> &coarse,
- const MGTransferBase<VECTOR> &transfer,
- const MGSmootherBase<VECTOR> &pre_smooth,
- const MGSmootherBase<VECTOR> &post_smooth,
- Cycle cycle = v_cycle);
+ Multigrid(const unsigned int minlevel,
+ const unsigned int maxlevel,
+ const MGMatrixBase<VectorType> &matrix,
+ const MGCoarseGridBase<VectorType> &coarse,
+ const MGTransferBase<VectorType> &transfer,
+ const MGSmootherBase<VectorType> &pre_smooth,
+ const MGSmootherBase<VectorType> &post_smooth,
+ Cycle cycle = v_cycle);
/**
* Reinit this class according to #minlevel and #maxlevel.
* <tt>edge_in</tt>. In particular, for symmetric operators, both arguments
* can refer to the same matrix, saving assembling of one of them.
*/
- void set_edge_matrices (const MGMatrixBase<VECTOR> &edge_out,
- const MGMatrixBase<VECTOR> &edge_in);
+ void set_edge_matrices (const MGMatrixBase<VectorType> &edge_out,
+ const MGMatrixBase<VectorType> &edge_in);
/**
* Set additional matrices to correct residual computation at refinement
* <tt>edge_up</tt>. In particular, for symmetric operators, both arguments
* can refer to the same matrix, saving assembling of one of them.
*/
- void set_edge_flux_matrices (const MGMatrixBase<VECTOR> &edge_down,
- const MGMatrixBase<VECTOR> &edge_up);
+ void set_edge_flux_matrices (const MGMatrixBase<VectorType> &edge_down,
+ const MGMatrixBase<VectorType> &edge_up);
/**
* Return the finest level for multigrid.
* Input vector for the cycle. Contains the defect of the outer method
* projected to the multilevel vectors.
*/
- MGLevelObject<VECTOR> defect;
+ MGLevelObject<VectorType> defect;
/**
* The solution update after the multigrid step.
*/
- MGLevelObject<VECTOR> solution;
+ MGLevelObject<VectorType> solution;
private:
/**
* Auxiliary vector.
*/
- MGLevelObject<VECTOR> t;
+ MGLevelObject<VectorType> t;
/**
* Auxiliary vector for W- and F-cycles. Left uninitialized in V-cycle.
*/
- MGLevelObject<VECTOR> defect2;
+ MGLevelObject<VectorType> defect2;
/**
* The matrix for each level.
*/
- SmartPointer<const MGMatrixBase<VECTOR>,Multigrid<VECTOR> > matrix;
+ SmartPointer<const MGMatrixBase<VectorType>,Multigrid<VectorType> > matrix;
/**
* The matrix for each level.
*/
- SmartPointer<const MGCoarseGridBase<VECTOR>,Multigrid<VECTOR> > coarse;
+ SmartPointer<const MGCoarseGridBase<VectorType>,Multigrid<VectorType> > coarse;
/**
* Object for grid tranfer.
*/
- SmartPointer<const MGTransferBase<VECTOR>,Multigrid<VECTOR> > transfer;
+ SmartPointer<const MGTransferBase<VectorType>,Multigrid<VectorType> > transfer;
/**
* The pre-smoothing object.
*/
- SmartPointer<const MGSmootherBase<VECTOR>,Multigrid<VECTOR> > pre_smooth;
+ SmartPointer<const MGSmootherBase<VectorType>,Multigrid<VectorType> > pre_smooth;
/**
* The post-smoothing object.
*/
- SmartPointer<const MGSmootherBase<VECTOR>,Multigrid<VECTOR> > post_smooth;
+ SmartPointer<const MGSmootherBase<VectorType>,Multigrid<VectorType> > post_smooth;
/**
* Edge matrix from the interior of the refined part to the refinement edge.
*
* @note Only <tt>vmult</tt> is used for these matrices.
*/
- SmartPointer<const MGMatrixBase<VECTOR> > edge_out;
+ SmartPointer<const MGMatrixBase<VectorType> > edge_out;
/**
* Transpose edge matrix from the refinement edge to the interior of the
*
* @note Only <tt>Tvmult</tt> is used for these matrices.
*/
- SmartPointer<const MGMatrixBase<VECTOR> > edge_in;
+ SmartPointer<const MGMatrixBase<VectorType> > edge_in;
/**
* Edge matrix from fine to coarse.
*
* @note Only <tt>vmult</tt> is used for these matrices.
*/
- SmartPointer<const MGMatrixBase<VECTOR>,Multigrid<VECTOR> > edge_down;
+ SmartPointer<const MGMatrixBase<VectorType>,Multigrid<VectorType> > edge_down;
/**
* Transpose edge matrix from coarse to fine.
*
* @note Only <tt>Tvmult</tt> is used for these matrices.
*/
- SmartPointer<const MGMatrixBase<VECTOR>,Multigrid<VECTOR> > edge_up;
+ SmartPointer<const MGMatrixBase<VectorType>,Multigrid<VectorType> > edge_up;
/**
* Level for debug output. Defaults to zero and can be set by set_debug().
* multi-level preconditioning and provide the standard interface for LAC
* iterative methods.
*
- * Furthermore, it needs functions <tt>void copy_to_mg(const VECTOR&)</tt> to
+ * Furthermore, it needs functions <tt>void copy_to_mg(const VectorType&)</tt> to
* store @p src in the right hand side of the multi-level method and <tt>void
- * copy_from_mg(VECTOR&)</tt> to store the result of the v-cycle in @p dst.
+ * copy_from_mg(VectorType&)</tt> to store the result of the v-cycle in @p dst.
*
* @author Guido Kanschat, 1999, 2000, 2001, 2002
*/
-template<int dim, class VECTOR, class TRANSFER>
+template<int dim, typename VectorType, class TRANSFER>
class PreconditionMG : public Subscriptor
{
public:
* Constructor. Arguments are the multigrid object, pre-smoother, post-
* smoother and coarse grid solver.
*/
- PreconditionMG(const DoFHandler<dim> &dof_handler,
- Multigrid<VECTOR> &mg,
- const TRANSFER &transfer);
+ PreconditionMG(const DoFHandler<dim> &dof_handler,
+ Multigrid<VectorType> &mg,
+ const TRANSFER &transfer);
/**
* Dummy function needed by other classes.
/**
* Associated @p DoFHandler.
*/
- SmartPointer<const DoFHandler<dim>,PreconditionMG<dim,VECTOR,TRANSFER> > dof_handler;
+ SmartPointer<const DoFHandler<dim>,PreconditionMG<dim,VectorType,TRANSFER> > dof_handler;
/**
* The multigrid object.
*/
- SmartPointer<Multigrid<VECTOR>,PreconditionMG<dim,VECTOR,TRANSFER> > multigrid;
+ SmartPointer<Multigrid<VectorType>,PreconditionMG<dim,VectorType,TRANSFER> > multigrid;
/**
* Object for grid tranfer.
*/
- SmartPointer<const TRANSFER,PreconditionMG<dim,VECTOR,TRANSFER> > transfer;
+ SmartPointer<const TRANSFER,PreconditionMG<dim,VectorType,TRANSFER> > transfer;
};
/*@}*/
/* --------------------------- inline functions --------------------- */
-template <class VECTOR>
+template <typename VectorType>
template <int dim>
-Multigrid<VECTOR>::Multigrid (const DoFHandler<dim> &mg_dof_handler,
- const MGMatrixBase<VECTOR> &matrix,
- const MGCoarseGridBase<VECTOR> &coarse,
- const MGTransferBase<VECTOR> &transfer,
- const MGSmootherBase<VECTOR> &pre_smooth,
- const MGSmootherBase<VECTOR> &post_smooth,
- Cycle cycle)
+Multigrid<VectorType>::Multigrid (const DoFHandler<dim> &mg_dof_handler,
+ const MGMatrixBase<VectorType> &matrix,
+ const MGCoarseGridBase<VectorType> &coarse,
+ const MGTransferBase<VectorType> &transfer,
+ const MGSmootherBase<VectorType> &pre_smooth,
+ const MGSmootherBase<VectorType> &post_smooth,
+ Cycle cycle)
:
cycle_type(cycle),
minlevel(0),
-template <class VECTOR>
+template <typename VectorType>
inline
unsigned int
-Multigrid<VECTOR>::get_maxlevel () const
+Multigrid<VectorType>::get_maxlevel () const
{
return maxlevel;
}
-template <class VECTOR>
+template <typename VectorType>
inline
unsigned int
-Multigrid<VECTOR>::get_minlevel () const
+Multigrid<VectorType>::get_minlevel () const
{
return minlevel;
}
/* --------------------------- inline functions --------------------- */
-template<int dim, class VECTOR, class TRANSFER>
-PreconditionMG<dim, VECTOR, TRANSFER>
-::PreconditionMG(const DoFHandler<dim> &dof_handler,
- Multigrid<VECTOR> &mg,
- const TRANSFER &transfer)
+template<int dim, typename VectorType, class TRANSFER>
+PreconditionMG<dim, VectorType, TRANSFER>
+::PreconditionMG(const DoFHandler<dim> &dof_handler,
+ Multigrid<VectorType> &mg,
+ const TRANSFER &transfer)
:
dof_handler(&dof_handler),
multigrid(&mg),
transfer(&transfer)
{}
-template<int dim, class VECTOR, class TRANSFER>
+template<int dim, typename VectorType, class TRANSFER>
inline bool
-PreconditionMG<dim, VECTOR, TRANSFER>::empty () const
+PreconditionMG<dim, VectorType, TRANSFER>::empty () const
{
return false;
}
-template<int dim, class VECTOR, class TRANSFER>
+template<int dim, typename VectorType, class TRANSFER>
template<class VECTOR2>
void
-PreconditionMG<dim, VECTOR, TRANSFER>::vmult (
+PreconditionMG<dim, VectorType, TRANSFER>::vmult (
VECTOR2 &dst,
const VECTOR2 &src) const
{
}
-template<int dim, class VECTOR, class TRANSFER>
+template<int dim, typename VectorType, class TRANSFER>
template<class VECTOR2>
void
-PreconditionMG<dim, VECTOR, TRANSFER>::vmult_add (
+PreconditionMG<dim, VectorType, TRANSFER>::vmult_add (
VECTOR2 &dst,
const VECTOR2 &src) const
{
}
-template<int dim, class VECTOR, class TRANSFER>
+template<int dim, typename VectorType, class TRANSFER>
template<class VECTOR2>
void
-PreconditionMG<dim, VECTOR, TRANSFER>::Tvmult (
+PreconditionMG<dim, VectorType, TRANSFER>::Tvmult (
VECTOR2 &,
const VECTOR2 &) const
{
}
-template<int dim, class VECTOR, class TRANSFER>
+template<int dim, typename VectorType, class TRANSFER>
template<class VECTOR2>
void
-PreconditionMG<dim, VECTOR, TRANSFER>::Tvmult_add (
+PreconditionMG<dim, VectorType, TRANSFER>::Tvmult_add (
VECTOR2 &,
const VECTOR2 &) const
{
DEAL_II_NAMESPACE_OPEN
-template <class VECTOR>
-Multigrid<VECTOR>::Multigrid (const unsigned int minlevel,
- const unsigned int maxlevel,
- const MGMatrixBase<VECTOR> &matrix,
- const MGCoarseGridBase<VECTOR> &coarse,
- const MGTransferBase<VECTOR> &transfer,
- const MGSmootherBase<VECTOR> &pre_smooth,
- const MGSmootherBase<VECTOR> &post_smooth,
- typename Multigrid<VECTOR>::Cycle cycle)
+template <typename VectorType>
+Multigrid<VectorType>::Multigrid (const unsigned int minlevel,
+ const unsigned int maxlevel,
+ const MGMatrixBase<VectorType> &matrix,
+ const MGCoarseGridBase<VectorType> &coarse,
+ const MGTransferBase<VectorType> &transfer,
+ const MGSmootherBase<VectorType> &pre_smooth,
+ const MGSmootherBase<VectorType> &post_smooth,
+ typename Multigrid<VectorType>::Cycle cycle)
:
cycle_type(cycle),
minlevel(minlevel),
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::reinit(const unsigned int min_level,
- const unsigned int max_level)
+Multigrid<VectorType>::reinit (const unsigned int min_level,
+ const unsigned int max_level)
{
minlevel=min_level;
maxlevel=max_level;
}
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::set_maxlevel (const unsigned int l)
+Multigrid<VectorType>::set_maxlevel (const unsigned int l)
{
Assert (l <= maxlevel, ExcIndexRange(l,minlevel,maxlevel+1));
Assert (l >= minlevel, ExcIndexRange(l,minlevel,maxlevel+1));
}
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::set_minlevel (const unsigned int l,
- const bool relative)
+Multigrid<VectorType>::set_minlevel (const unsigned int l,
+ const bool relative)
{
Assert (l <= maxlevel, ExcIndexRange(l,minlevel,maxlevel+1));
minlevel = (relative)
}
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::set_cycle(typename Multigrid<VECTOR>::Cycle c)
+Multigrid<VectorType>::set_cycle(typename Multigrid<VectorType>::Cycle c)
{
cycle_type = c;
}
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::set_debug (const unsigned int d)
+Multigrid<VectorType>::set_debug (const unsigned int d)
{
debug = d;
}
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::set_edge_matrices (const MGMatrixBase<VECTOR> &down,
- const MGMatrixBase<VECTOR> &up)
+Multigrid<VectorType>::set_edge_matrices (const MGMatrixBase<VectorType> &down,
+ const MGMatrixBase<VectorType> &up)
{
edge_out = &down;
edge_in = &up;
}
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::set_edge_flux_matrices (const MGMatrixBase<VECTOR> &down,
- const MGMatrixBase<VECTOR> &up)
+Multigrid<VectorType>::set_edge_flux_matrices (const MGMatrixBase<VectorType> &down,
+ const MGMatrixBase<VectorType> &up)
{
edge_down = &down;
edge_up = &up;
}
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::level_v_step(const unsigned int level)
+Multigrid<VectorType>::level_v_step (const unsigned int level)
{
if (debug>0)
deallog << "V-cycle entering level " << level << std::endl;
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::level_step(const unsigned int level,
- Cycle cycle)
+Multigrid<VectorType>::level_step(const unsigned int level,
+ Cycle cycle)
{
char cychar = '?';
switch (cycle)
}
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::cycle()
+Multigrid<VectorType>::cycle()
{
// The defect vector has been
// initialized by copy_to_mg. Now
}
-template <class VECTOR>
+template <typename VectorType>
void
-Multigrid<VECTOR>::vcycle()
+Multigrid<VectorType>::vcycle()
{
// The defect vector has been
// initialized by copy_to_mg. Now
* includes all of the usual vector types, but also IndexSet (see step-41
* for a use of this).
*/
- template <class VECTOR>
- void add_data_vector (const VECTOR &data,
+ template <class VectorType>
+ void add_data_vector (const VectorType &data,
const std::vector<std::string> &names,
- const DataVectorType type = type_automatic,
+ const DataVectorType type = type_automatic,
const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation
= std::vector<DataComponentInterpretation::DataComponentInterpretation>());
* which FEValues can extract values on a cell using the
* FEValuesBase::get_function_values() function.
*/
- template <class VECTOR>
- void add_data_vector (const VECTOR &data,
+ template <class VectorType>
+ void add_data_vector (const VectorType &data,
const std::string &name,
- const DataVectorType type = type_automatic,
+ const DataVectorType type = type_automatic,
const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation
= std::vector<DataComponentInterpretation::DataComponentInterpretation>());
* represents dof data, the data vector type argument present in the other
* methods above is skipped.
*/
- template <class VECTOR>
+ template <class VectorType>
void add_data_vector (const DH &dof_handler,
- const VECTOR &data,
+ const VectorType &data,
const std::vector<std::string> &names,
const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation
= std::vector<DataComponentInterpretation::DataComponentInterpretation>());
* This function is an abbreviation of the function above with only a scalar
* @p dof_handler given and a single data name.
*/
- template <class VECTOR>
- void add_data_vector (const DH &dof_handler,
- const VECTOR &data,
- const std::string &name,
+ template <class VectorType>
+ void add_data_vector (const DH &dof_handler,
+ const VectorType &data,
+ const std::string &name,
const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation
= std::vector<DataComponentInterpretation::DataComponentInterpretation>());
* error by declaring the data postprocessor variable before the DataOut
* variable as objects are destroyed in reverse order of declaration.
*/
- template <class VECTOR>
- void add_data_vector (const VECTOR &data,
+ template <class VectorType>
+ void add_data_vector (const VectorType &data,
const DataPostprocessor<DH::space_dimension> &data_postprocessor);
/**
* postprocessor can only read data from the given DoFHandler and solution
* vector, not other solution vectors or DoFHandlers.
*/
- template <class VECTOR>
- void add_data_vector (const DH &dof_handler,
- const VECTOR &data,
+ template <class VectorType>
+ void add_data_vector (const DH &dof_handler,
+ const VectorType &data,
const DataPostprocessor<DH::space_dimension> &data_postprocessor);
/**
* An output operator writing a separate file in each step and writing the
* vectors as finite element functions with respect to a given DoFHandler.
*/
- template <class VECTOR, int dim, int spacedim=dim>
- class DoFOutputOperator : public OutputOperator<VECTOR>
+ template <typename VectorType, int dim, int spacedim=dim>
+ class DoFOutputOperator : public OutputOperator<VectorType>
{
public:
/*
void parse_parameters(ParameterHandler ¶m);
void initialize (const DoFHandler<dim, spacedim> &dof_handler);
- virtual OutputOperator<VECTOR> &
+ virtual OutputOperator<VectorType> &
operator << (const AnyData &vectors);
private:
SmartPointer<const DoFHandler<dim, spacedim>,
- DoFOutputOperator<VECTOR, dim, spacedim> > dof;
+ DoFOutputOperator<VectorType, dim, spacedim> > dof;
const std::string filename_base;
const unsigned int digits;
DataOut<dim> out;
};
- template <class VECTOR, int dim, int spacedim>
+ template <typename VectorType, int dim, int spacedim>
inline void
- DoFOutputOperator<VECTOR, dim, spacedim>::initialize(const DoFHandler<dim, spacedim> &dof_handler)
+ DoFOutputOperator<VectorType, dim, spacedim>::initialize(const DoFHandler<dim, spacedim> &dof_handler)
{
dof = &dof_handler;
}
namespace Algorithms
{
- template <class VECTOR, int dim, int spacedim>
- DoFOutputOperator<VECTOR, dim, spacedim>::DoFOutputOperator (
+ template <typename VectorType, int dim, int spacedim>
+ DoFOutputOperator<VectorType, dim, spacedim>::DoFOutputOperator (
const std::string filename_base,
const unsigned int digits)
:
}
- template <class VECTOR, int dim, int spacedim>
+ template <typename VectorType, int dim, int spacedim>
void
- DoFOutputOperator<VECTOR, dim, spacedim>::parse_parameters(ParameterHandler ¶m)
+ DoFOutputOperator<VectorType, dim, spacedim>::parse_parameters(ParameterHandler ¶m)
{
out.parse_parameters(param);
}
- template <class VECTOR, int dim, int spacedim>
- OutputOperator<VECTOR> &
- DoFOutputOperator<VECTOR, dim, spacedim>::operator<<(
+ template <typename VectorType, int dim, int spacedim>
+ OutputOperator<VectorType> &
+ DoFOutputOperator<VectorType, dim, spacedim>::operator<<(
const AnyData &data)
{
Assert ((dof!=0), ExcNotInitialized());
out.attach_dof_handler (*dof);
for (unsigned int i=0; i<data.size(); ++i)
{
- const VECTOR *p = data.try_read_ptr<VECTOR>(i);
+ const VectorType *p = data.try_read_ptr<VectorType>(i);
if (p!=0)
{
out.add_data_vector (*p, data.name(i));
* @ingroup output
* @author Guido Kanschat, 2000
*/
-template<int dim, class SOLVER, class VECTOR = Vector<double> >
+template<int dim, class SOLVER, class VectorType = Vector<double> >
class DoFPrintSolverStep : public SOLVER
{
public:
* produced for each iteration step.
*/
DoFPrintSolverStep (SolverControl &control,
- VectorMemory<VECTOR> &mem,
- DataOut<dim> &data_out,
- const std::string &basename);
+ VectorMemory<VectorType> &mem,
+ DataOut<dim> &data_out,
+ const std::string &basename);
/**
* Call-back function for the iterative method.
*/
virtual void print_vectors (const unsigned int step,
- const VECTOR &x,
- const VECTOR &r,
- const VECTOR &d) const;
+ const VectorType &x,
+ const VectorType &r,
+ const VectorType &d) const;
private:
/**
* Output object.
/* ----------------------- template functions --------------- */
-template<int dim, class SOLVER, class VECTOR>
-DoFPrintSolverStep<dim, SOLVER, VECTOR>::DoFPrintSolverStep (SolverControl &control,
- VectorMemory<VECTOR> &mem,
- DataOut<dim> &data_out,
- const std::string &basename)
+template<int dim, class SOLVER, class VectorType>
+DoFPrintSolverStep<dim, SOLVER, VectorType>::DoFPrintSolverStep
+(SolverControl &control,
+ VectorMemory<VectorType> &mem,
+ DataOut<dim> &data_out,
+ const std::string &basename)
: SOLVER (control, mem),
out (data_out),
basename (basename)
{}
-template<int dim, class SOLVER, class VECTOR>
+template<int dim, class SOLVER, class VectorType>
void
-DoFPrintSolverStep<dim, SOLVER, VECTOR>::print_vectors (const unsigned int step,
- const VECTOR &x,
- const VECTOR &r,
- const VECTOR &d) const
+DoFPrintSolverStep<dim, SOLVER, VectorType>::print_vectors (const unsigned int step,
+ const VectorType &x,
+ const VectorType &r,
+ const VectorType &d) const
{
out.clear_data_vectors();
out.add_data_vector(x, "solution");
*/
template <int dim,
typename DH=DoFHandler<dim>,
- typename VECTOR=Vector<double> >
+ typename VectorType=Vector<double> >
class FEFieldFunction : public Function<dim>
{
public:
* lay. Otherwise the standard Q1 mapping is used.
*/
FEFieldFunction (const DH &dh,
- const VECTOR &data_vector,
+ const VectorType &data_vector,
const Mapping<dim> &mapping = StaticMappingQ1<dim>::mapping);
/**
/**
* Pointer to the dof handler.
*/
- SmartPointer<const DH,FEFieldFunction<dim, DH, VECTOR> > dh;
+ SmartPointer<const DH,FEFieldFunction<dim, DH, VectorType> > dh;
/**
* A reference to the actual data vector.
*/
- const VECTOR &data_vector;
+ const VectorType &data_vector;
/**
* A reference to the mapping being used.
namespace Functions
{
- template <int dim, typename DH, typename VECTOR>
- FEFieldFunction<dim, DH, VECTOR>::FEFieldFunction (const DH &mydh,
- const VECTOR &myv,
- const Mapping<dim> &mymapping)
+ template <int dim, typename DH, typename VectorType>
+ FEFieldFunction<dim, DH, VectorType>::FEFieldFunction (const DH &mydh,
+ const VectorType &myv,
+ const Mapping<dim> &mymapping)
:
Function<dim>(mydh.get_fe().n_components()),
dh(&mydh, "FEFieldFunction"),
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
void
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
set_active_cell(const typename DH::active_cell_iterator &newcell)
{
cell_hint.get() = newcell;
- template <int dim, typename DH, typename VECTOR>
- void FEFieldFunction<dim, DH, VECTOR>::vector_value (const Point<dim> &p,
- Vector<double> &values) const
+ template <int dim, typename DH, typename VectorType>
+ void FEFieldFunction<dim, DH, VectorType>::vector_value (const Point<dim> &p,
+ Vector<double> &values) const
{
Assert (values.size() == n_components,
ExcDimensionMismatch(values.size(), n_components));
FEValues<dim> fe_v(mapping, cell->get_fe(), quad,
update_values);
fe_v.reinit(cell);
- std::vector< Vector<typename VECTOR::value_type> >
- vvalues (1, Vector<typename VECTOR::value_type>(values.size()));
+ std::vector< Vector<typename VectorType::value_type> >
+ vvalues (1, Vector<typename VectorType::value_type>(values.size()));
fe_v.get_function_values(data_vector, vvalues);
values = vvalues[0];
}
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
double
- FEFieldFunction<dim, DH, VECTOR>::value (const Point<dim> &p,
- const unsigned int comp) const
+ FEFieldFunction<dim, DH, VectorType>::value (const Point<dim> &p,
+ const unsigned int comp) const
{
Vector<double> values(n_components);
vector_value(p, values);
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
void
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
vector_gradient (const Point<dim> &p,
std::vector<Tensor<1,dim> > &gradients) const
{
// FIXME: we need a temp argument because get_function_values wants to put
// its data into an object storing the correct scalar type, but this
// function wants to return everything in a vector<double>
- std::vector< std::vector<Tensor<1,dim,typename VECTOR::value_type> > > vgrads
- (1, std::vector<Tensor<1,dim,typename VECTOR::value_type> >(n_components) );
+ std::vector< std::vector<Tensor<1,dim,typename VectorType::value_type> > > vgrads
+ (1, std::vector<Tensor<1,dim,typename VectorType::value_type> >(n_components) );
fe_v.get_function_gradients(data_vector, vgrads);
gradients = std::vector<Tensor<1,dim> >(vgrads[0].begin(), vgrads[0].end());
}
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
Tensor<1,dim>
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
gradient (const Point<dim> &p,
const unsigned int comp) const
{
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
void
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
vector_laplacian (const Point<dim> &p,
Vector<double> &values) const
{
FEValues<dim> fe_v(mapping, cell->get_fe(), quad,
update_hessians);
fe_v.reinit(cell);
- std::vector< Vector<typename VECTOR::value_type> >
- vvalues (1, Vector<typename VECTOR::value_type>(values.size()));
+ std::vector< Vector<typename VectorType::value_type> >
+ vvalues (1, Vector<typename VectorType::value_type>(values.size()));
fe_v.get_function_laplacians(data_vector, vvalues);
values = vvalues[0];
}
- template <int dim, typename DH, typename VECTOR>
- double FEFieldFunction<dim, DH, VECTOR>::laplacian
+ template <int dim, typename DH, typename VectorType>
+ double FEFieldFunction<dim, DH, VectorType>::laplacian
(const Point<dim> &p, const unsigned int comp) const
{
Vector<double> lap(n_components);
// Now the list versions
// ==============================
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
void
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
vector_value_list (const std::vector<Point< dim > > &points,
std::vector< Vector<double> > &values) const
{
{
fe_v.reinit(cells[i], i, 0);
const unsigned int nq = qpoints[i].size();
- std::vector< Vector<typename VECTOR::value_type> > vvalues (nq, Vector<typename VECTOR::value_type>(n_components));
+ std::vector< Vector<typename VectorType::value_type> > vvalues (nq, Vector<typename VectorType::value_type>(n_components));
fe_v.get_present_fe_values ().get_function_values(data_vector, vvalues);
for (unsigned int q=0; q<nq; ++q)
values[maps[i][q]] = vvalues[q];
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
void
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
value_list (const std::vector<Point< dim > > &points,
- std::vector< double > &values,
- const unsigned int component) const
+ std::vector< double > &values,
+ const unsigned int component) const
{
Assert(points.size() == values.size(),
ExcDimensionMismatch(points.size(), values.size()));
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
void
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
vector_gradient_list (const std::vector<Point< dim > > &points,
std::vector<
std::vector< Tensor<1,dim> > > &values) const
{
fe_v.reinit(cells[i], i, 0);
const unsigned int nq = qpoints[i].size();
- std::vector< std::vector<Tensor<1,dim,typename VECTOR::value_type> > >
- vgrads (nq, std::vector<Tensor<1,dim,typename VECTOR::value_type> >(n_components));
+ std::vector< std::vector<Tensor<1,dim,typename VectorType::value_type> > >
+ vgrads (nq, std::vector<Tensor<1,dim,typename VectorType::value_type> >(n_components));
fe_v.get_present_fe_values ().get_function_gradients(data_vector, vgrads);
for (unsigned int q=0; q<nq; ++q)
{
}
}
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
void
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
gradient_list (const std::vector<Point< dim > > &points,
std::vector< Tensor<1,dim> > &values,
const unsigned int component) const
}
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
void
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
vector_laplacian_list (const std::vector<Point< dim > > &points,
std::vector< Vector<double> > &values) const
{
{
fe_v.reinit(cells[i], i, 0);
const unsigned int nq = qpoints[i].size();
- std::vector< Vector<typename VECTOR::value_type> > vvalues (nq, Vector<typename VECTOR::value_type>(n_components));
+ std::vector< Vector<typename VectorType::value_type> > vvalues (nq, Vector<typename VectorType::value_type>(n_components));
fe_v.get_present_fe_values ().get_function_laplacians(data_vector, vvalues);
for (unsigned int q=0; q<nq; ++q)
values[maps[i][q]] = vvalues[q];
}
}
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
void
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
laplacian_list (const std::vector<Point< dim > > &points,
std::vector< double > &values,
const unsigned int component) const
- template <int dim, typename DH, typename VECTOR>
- unsigned int FEFieldFunction<dim, DH, VECTOR>::
+ template <int dim, typename DH, typename VectorType>
+ unsigned int FEFieldFunction<dim, DH, VectorType>::
compute_point_locations(const std::vector<Point<dim> > &points,
std::vector<typename DH::active_cell_iterator > &cells,
std::vector<std::vector<Point<dim> > > &qpoints,
}
- template <int dim, typename DH, typename VECTOR>
+ template <int dim, typename DH, typename VectorType>
boost::optional<Point<dim> >
- FEFieldFunction<dim, DH, VECTOR>::
+ FEFieldFunction<dim, DH, VectorType>::
get_reference_coordinates (const typename DH::active_cell_iterator &cell,
const Point<dim> &point) const
{
/**
- * Put another mnemonic string (and hence @p VECTOR) into the class. This
+ * Put another mnemonic string (and hence @p VectorType) into the class. This
* method adds storage space for variables equal to the number of true
* values in component_mask. This also adds extra entries for points that
* are already in the class, so @p add_field_name and @p add_points can be
const ComponentMask &component_mask = ComponentMask());
/**
- * Put another mnemonic string (and hence @p VECTOR) into the class. This
+ * Put another mnemonic string (and hence @p VectorType) into the class. This
* method adds storage space for n_components variables. This also adds
* extra entries for points that are already in the class, so @p
* add_field_name and @p add_points can be called in any order. This method
/**
- * Extract values at the stored points from the VECTOR supplied and add them
- * to the new dataset in vector_name. The component mask supplied when the
- * field was added is used to select components to extract. If a @p
+ * Extract values at the stored points from the VectorType supplied and add
+ * them to the new dataset in vector_name. The component mask supplied when
+ * the field was added is used to select components to extract. If a @p
* DoFHandler is used, one (and only one) evaluate_field method must be
* called for each dataset (time step, iteration, etc) for each vector_name,
* otherwise a @p ExcDataLostSync error can occur.
*/
- template <class VECTOR>
+ template <class VectorType>
void evaluate_field(const std::string &name,
- const VECTOR &solution);
+ const VectorType &solution);
/**
- * Compute values using a @p DataPostprocessor object with the @p VECTOR
+ * Compute values using a @p DataPostprocessor object with the @p VectorType
* supplied and add them to the new dataset in vector_name. The
* component_mask supplied when the field was added is used to select
* components to extract from the @p DataPostprocessor return vector. This
* method must be called for each dataset (time step, iteration, etc) for
* each vector_name, otherwise a @p ExcDataLostSync error can occur.
*/
- template <class VECTOR>
+ template <class VectorType>
void evaluate_field(const std::vector <std::string> &names,
- const VECTOR &solution,
- const DataPostprocessor<dim> &data_postprocessor,
- const Quadrature<dim> &quadrature);
+ const VectorType &solution,
+ const DataPostprocessor<dim> &data_postprocessor,
+ const Quadrature<dim> &quadrature);
/**
* Construct a std::vector <std::string> containing only vector_name and
* call the above function. The above function is more efficient if multiple
* fields use the same @p DataPostprocessor object.
*/
- template <class VECTOR>
- void evaluate_field(const std::string &name,
- const VECTOR &solution,
+ template <class VectorType>
+ void evaluate_field(const std::string &name,
+ const VectorType &solution,
const DataPostprocessor<dim> &data_postprocessor,
- const Quadrature<dim> &quadrature);
+ const Quadrature<dim> &quadrature);
/**
- * Extract values at the points actually requested from the VECTOR supplied
+ * Extract values at the points actually requested from the VectorType supplied
* and add them to the new dataset in vector_name. Unlike the other
* evaluate_field methods this method does not care if the dof_handler has
* been modified because it uses calls to @p VectorTools::point_value to
* called for each dataset (time step, iteration, etc) for each vector_name,
* otherwise a @p ExcDataLostSync error can occur.
*/
- template <class VECTOR>
+ template <class VectorType>
void evaluate_field_at_requested_location(const std::string &name,
- const VECTOR &solution);
+ const VectorType &solution);
/**
* @author Ralf Hartmann, 1999, Oliver Kayser-Herold and Wolfgang Bangerth,
* 2006, Wolfgang Bangerth 2014
*/
-template<int dim, typename VECTOR=Vector<double>, class DH=DoFHandler<dim> >
+template<int dim, typename VectorType=Vector<double>, class DH=DoFHandler<dim> >
class SolutionTransfer
{
public:
* vectors that are to be interpolated onto the new (refined and/or
* coarsenend) grid.
*/
- void prepare_for_coarsening_and_refinement (const std::vector<VECTOR> &all_in);
+ void prepare_for_coarsening_and_refinement (const std::vector<VectorType> &all_in);
/**
* Same as previous function but for only one discrete function to be
* interpolated.
*/
- void prepare_for_coarsening_and_refinement (const VECTOR &in);
+ void prepare_for_coarsening_and_refinement (const VectorType &in);
/**
* This function interpolates the discrete function @p in, which is a vector
* is called and the refinement is executed before. Multiple calling of this
* function is allowed. e.g. for interpolating several functions.
*/
- void refine_interpolate (const VECTOR &in,
- VECTOR &out) const;
+ void refine_interpolate (const VectorType &in,
+ VectorType &out) const;
/**
* This function interpolates the discrete functions that are stored in @p
* the right size (@p n_dofs_refined). Otherwise an assertion will be
* thrown.
*/
- void interpolate (const std::vector<VECTOR> &all_in,
- std::vector<VECTOR> &all_out) const;
+ void interpolate (const std::vector<VectorType> &all_in,
+ std::vector<VectorType> &all_out) const;
/**
* Same as the previous function. It interpolates only one function. It
* functions can be performed in one step by using <tt>interpolate (all_in,
* all_out)</tt>
*/
- void interpolate (const VECTOR &in,
- VECTOR &out) const;
+ void interpolate (const VectorType &in,
+ VectorType &out) const;
/**
* Determine an estimate for the memory consumption (in bytes) of this
/**
* Pointer to the degree of freedom handler to work with.
*/
- SmartPointer<const DH,SolutionTransfer<dim,VECTOR,DH> > dof_handler;
+ SmartPointer<const DH,SolutionTransfer<dim,VectorType,DH> > dof_handler;
/**
* Stores the number of DoFs before the refinement and/or coarsening.
indices_ptr(indices_ptr_in),
dof_values_ptr (0),
active_fe_index(active_fe_index_in) {};
- Pointerstruct(std::vector<Vector<typename VECTOR::value_type> > *dof_values_ptr_in,
+ Pointerstruct(std::vector<Vector<typename VectorType::value_type> > *dof_values_ptr_in,
const unsigned int active_fe_index_in = 0) :
indices_ptr (0),
dof_values_ptr(dof_values_ptr_in),
std::size_t memory_consumption () const;
std::vector<types::global_dof_index> *indices_ptr;
- std::vector<Vector<typename VECTOR::value_type> > *dof_values_ptr;
+ std::vector<Vector<typename VectorType::value_type> > *dof_values_ptr;
unsigned int active_fe_index;
};
* Is used for @p prepare_for_refining_and_coarsening The interpolated dof
* values of all cells that'll be coarsened will be stored in this vector.
*/
- std::vector<std::vector<Vector<typename VECTOR::value_type> > > dof_values_on_cell;
+ std::vector<std::vector<Vector<typename VectorType::value_type> > > dof_values_on_cell;
};
* @todo The @p mapping argument should be replaced by a
* hp::MappingCollection in case of a hp::DoFHandler.
*/
- template <class VECTOR, int dim, int spacedim, template <int,int> class DH>
- void interpolate (const Mapping<dim,spacedim> &mapping,
- const DH<dim,spacedim> &dof,
- const Function<spacedim,double> &function,
- VECTOR &vec);
+ template <typename VectorType, int dim, int spacedim, template <int,int> class DH>
+ void interpolate (const Mapping<dim,spacedim> &mapping,
+ const DH<dim,spacedim> &dof,
+ const Function<spacedim,double> &function,
+ VectorType &vec);
/**
* Calls the @p interpolate() function above with
* <tt>mapping=MappingQGeneric1@<dim>@()</tt>.
*/
- template <class VECTOR, class DH>
- void interpolate (const DH &dof,
- const Function<DH::space_dimension,double> &function,
- VECTOR &vec);
+ template <typename VectorType, class DH>
+ void interpolate (const DH &dof,
+ const Function<DH::space_dimension,double> &function,
+ VectorType &vec);
/**
* Interpolate different finite element spaces. The interpolation of vector
*
* @author Valentin Zingan, 2013
*/
- template<typename VECTOR, typename DH>
+ template<typename VectorType, typename DH>
void
- interpolate_based_on_material_id(const Mapping<DH::dimension, DH::space_dimension> &mapping,
- const DH &dof_handler,
- const std::map< types::material_id, const Function<DH::space_dimension,double>* > &function_map,
- VECTOR &dst,
- const ComponentMask &component_mask = ComponentMask());
+ interpolate_based_on_material_id
+ (const Mapping<DH::dimension, DH::space_dimension> &mapping,
+ const DH &dof_handler,
+ const std::map< types::material_id, const Function<DH::space_dimension,double>* > &function_map,
+ VectorType &dst,
+ const ComponentMask &component_mask = ComponentMask());
/**
* Gives the interpolation of a @p dof1-function @p u1 to a @p dof2-function
*/
template <int dim, int spacedim,
template <int,int> class DH,
- class VECTOR>
+ typename VectorType>
void
interpolate_to_different_mesh (const DH<dim, spacedim> &dof1,
- const VECTOR &u1,
+ const VectorType &u1,
const DH<dim, spacedim> &dof2,
- VECTOR &u2);
+ VectorType &u2);
/**
* Gives the interpolation of a @p dof1-function @p u1 to a @p dof2-function
*/
template <int dim, int spacedim,
template <int,int> class DH,
- class VECTOR>
+ typename VectorType>
void
interpolate_to_different_mesh (const DH<dim, spacedim> &dof1,
- const VECTOR &u1,
+ const VectorType &u1,
const DH<dim, spacedim> &dof2,
const ConstraintMatrix &constraints,
- VECTOR &u2);
+ VectorType &u2);
/**
*/
template <int dim, int spacedim,
template <int,int> class DH,
- class VECTOR>
+ typename VectorType>
void
interpolate_to_different_mesh (const InterGridMap<DH<dim, spacedim> > &intergridmap,
- const VECTOR &u1,
+ const VectorType &u1,
const ConstraintMatrix &constraints,
- VECTOR &u2);
+ VectorType &u2);
/**
* Compute the projection of @p function to the finite element space.
* In 1d, the default value of the boundary quadrature formula is an invalid
* object since integration on the boundary doesn't happen in 1d.
*/
- template <int dim, class VECTOR, int spacedim>
- void project (const Mapping<dim, spacedim> &mapping,
- const DoFHandler<dim,spacedim> &dof,
- const ConstraintMatrix &constraints,
- const Quadrature<dim> &quadrature,
- const Function<spacedim,double> &function,
- VECTOR &vec,
- const bool enforce_zero_boundary = false,
- const Quadrature<dim-1> &q_boundary = (dim > 1 ?
- QGauss<dim-1>(2) :
- Quadrature<dim-1>(0)),
- const bool project_to_boundary_first = false);
+ template <int dim, typename VectorType, int spacedim>
+ void project (const Mapping<dim, spacedim> &mapping,
+ const DoFHandler<dim,spacedim> &dof,
+ const ConstraintMatrix &constraints,
+ const Quadrature<dim> &quadrature,
+ const Function<spacedim,double> &function,
+ VectorType &vec,
+ const bool enforce_zero_boundary = false,
+ const Quadrature<dim-1> &q_boundary = (dim > 1 ?
+ QGauss<dim-1>(2) :
+ Quadrature<dim-1>(0)),
+ const bool project_to_boundary_first = false);
/**
* Calls the project() function above, with
* <tt>mapping=MappingQGeneric@<dim@>(1)</tt>.
*/
- template <int dim, class VECTOR, int spacedim>
- void project (const DoFHandler<dim,spacedim> &dof,
- const ConstraintMatrix &constraints,
- const Quadrature<dim> &quadrature,
- const Function<spacedim,double> &function,
- VECTOR &vec,
- const bool enforce_zero_boundary = false,
- const Quadrature<dim-1> &q_boundary = (dim > 1 ?
+ template <int dim, typename VectorType, int spacedim>
+ void project (const DoFHandler<dim,spacedim> &dof,
+ const ConstraintMatrix &constraints,
+ const Quadrature<dim> &quadrature,
+ const Function<spacedim,double> &function,
+ VectorType &vec,
+ const bool enforce_zero_boundary = false,
+ const Quadrature<dim-1> &q_boundary = (dim > 1 ?
QGauss<dim-1>(2) :
Quadrature<dim-1>(0)),
- const bool project_to_boundary_first = false);
+ const bool project_to_boundary_first = false);
/**
* Same as above, but for arguments of type hp::DoFHandler,
* hp::QuadratureCollection, hp::MappingCollection
*/
- template <int dim, class VECTOR, int spacedim>
- void project (const hp::MappingCollection<dim, spacedim> &mapping,
- const hp::DoFHandler<dim,spacedim> &dof,
- const ConstraintMatrix &constraints,
- const hp::QCollection<dim> &quadrature,
- const Function<spacedim,double> &function,
- VECTOR &vec,
- const bool enforce_zero_boundary = false,
- const hp::QCollection<dim-1> &q_boundary = hp::QCollection<dim-1>(dim > 1 ?
- QGauss<dim-1>(2) :
- Quadrature<dim-1>(0)),
- const bool project_to_boundary_first = false);
+ template <int dim, typename VectorType, int spacedim>
+ void project (const hp::MappingCollection<dim, spacedim> &mapping,
+ const hp::DoFHandler<dim,spacedim> &dof,
+ const ConstraintMatrix &constraints,
+ const hp::QCollection<dim> &quadrature,
+ const Function<spacedim,double> &function,
+ VectorType &vec,
+ const bool enforce_zero_boundary = false,
+ const hp::QCollection<dim-1> &q_boundary = hp::QCollection<dim-1>(dim > 1 ?
+ QGauss<dim-1>(2) :
+ Quadrature<dim-1>(0)),
+ const bool project_to_boundary_first = false);
/**
* Calls the project() function above, with a collection of
* $Q_1$ mapping objects, i.e., with hp::StaticMappingQ1::mapping_collection.
*/
- template <int dim, class VECTOR, int spacedim>
- void project (const hp::DoFHandler<dim,spacedim> &dof,
- const ConstraintMatrix &constraints,
- const hp::QCollection<dim> &quadrature,
- const Function<spacedim,double> &function,
- VECTOR &vec,
- const bool enforce_zero_boundary = false,
- const hp::QCollection<dim-1> &q_boundary = hp::QCollection<dim-1>(dim > 1 ?
- QGauss<dim-1>(2) :
- Quadrature<dim-1>(0)),
- const bool project_to_boundary_first = false);
+ template <int dim, typename VectorType, int spacedim>
+ void project (const hp::DoFHandler<dim,spacedim> &dof,
+ const ConstraintMatrix &constraints,
+ const hp::QCollection<dim> &quadrature,
+ const Function<spacedim,double> &function,
+ VectorType &vec,
+ const bool enforce_zero_boundary = false,
+ const hp::QCollection<dim-1> &q_boundary = hp::QCollection<dim-1>(dim > 1 ?
+ QGauss<dim-1>(2) :
+ Quadrature<dim-1>(0)),
+ const bool project_to_boundary_first = false);
/**
* Compute Dirichlet boundary conditions. This function makes up a map of
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void point_difference (const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
const Function<spacedim,double> &exact_solution,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void point_difference (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_value (const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_value (const hp::DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double
point_value (const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double
point_value (const hp::DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_value (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_value (const hp::MappingCollection<dim, spacedim> &mapping,
const hp::DoFHandler<dim,spacedim> &dof,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double
point_value (const Mapping<dim,spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double
point_value (const hp::MappingCollection<dim,spacedim> &mapping,
const hp::DoFHandler<dim,spacedim> &dof,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_gradient (const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_gradient (const hp::DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
Tensor<1, spacedim, typename VectorType::value_type>
point_gradient (const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
Tensor<1, spacedim, typename VectorType::value_type>
point_gradient (const hp::DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_gradient (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_gradient (const hp::MappingCollection<dim, spacedim> &mapping,
const hp::DoFHandler<dim,spacedim> &dof,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
Tensor<1, spacedim, typename VectorType::value_type>
point_gradient (const Mapping<dim,spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
* @note If the cell in which the point is found is not locally owned, an
* exception of type VectorTools::ExcPointNotAvailableHere is thrown.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
Tensor<1, spacedim, typename VectorType::value_type>
point_gradient (const hp::MappingCollection<dim,spacedim> &mapping,
const hp::DoFHandler<dim,spacedim> &dof,
* not equal to $(1,1,\ldots,1)^T$. For such elements, a different procedure
* has to be used when subtracting the mean value.
*/
- template <class VECTOR>
- void subtract_mean_value(VECTOR &v,
+ template <typename VectorType>
+ void subtract_mean_value(VectorType &v,
const std::vector<bool> &p_select = std::vector<bool>());
* Lagrangian elements. For all other elements, you will need to compute the
* mean value and subtract it right inside the evaluation routine.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double compute_mean_value (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
const Quadrature<dim> &quadrature,
* Calls the other compute_mean_value() function, see above, with
* <tt>mapping=MappingQGeneric@<dim@>(1)</tt>.
*/
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double compute_mean_value (const DoFHandler<dim,spacedim> &dof,
const Quadrature<dim> &quadrature,
const VectorType &v,
* points of a FE_Q() finite element of the same degree as the
* degree of the required components.
*
- * Curved manifold are respected, and the resulting VECTOR will be
+ * Curved manifold are respected, and the resulting VectorType will be
* geometrically consistent. The resulting map is guaranteed to be
* interpolatory at the support points of a FE_Q() finite element of
* the same degree as the degree of the required components.
*
* If the underlying finite element is an FE_Q(1)^spacedim, then the
- * resulting VECTOR is a finite element field representation of the vertices
- * of the Triangulation.
+ * resulting @p VectorType is a finite element field representation of the
+ * vertices of the Triangulation.
*
* The optional ComponentMask argument can be used to specify what
* components of the FiniteElement to use to describe the geometry. If no
*
* @author Luca Heltai, 2015
*/
- template<class DH, class VECTOR>
+ template<class DH, typename VectorType>
void get_position_vector(const DH &dh,
- VECTOR &vector,
+ VectorType &vector,
const ComponentMask &mask=ComponentMask());
//@}
namespace VectorTools
{
- template <class VECTOR, int dim, int spacedim, template <int,int> class DH>
- void interpolate (const Mapping<dim,spacedim> &mapping,
- const DH<dim,spacedim> &dof,
- const Function<spacedim> &function,
- VECTOR &vec)
+ template <typename VectorType, int dim, int spacedim, template <int,int> class DH>
+ void interpolate (const Mapping<dim,spacedim> &mapping,
+ const DH<dim,spacedim> &dof,
+ const Function<spacedim> &function,
+ VectorType &vec)
{
Assert (vec.size() == dof.n_dofs(),
ExcDimensionMismatch (vec.size(), dof.n_dofs()));
}
- template <class VECTOR, class DH>
- void interpolate (const DH &dof,
- const Function<DH::space_dimension> &function,
- VECTOR &vec)
+ template <typename VectorType, class DH>
+ void interpolate (const DH &dof,
+ const Function<DH::space_dimension> &function,
+ VectorType &vec)
{
interpolate(StaticMappingQ1<DH::dimension, DH::space_dimension>::mapping,
dof, function, vec);
}
- template<typename VECTOR, typename DH>
+ template<typename VectorType, typename DH>
void
interpolate_based_on_material_id(const Mapping<DH::dimension, DH::space_dimension> &mapping,
const DH &dof,
const std::map< types::material_id, const Function<DH::space_dimension>* > &function_map,
- VECTOR &dst,
+ VectorType &dst,
const ComponentMask &component_mask)
{
const unsigned int dim = DH::dimension;
- template <int dim, int spacedim,
- template <int,int> class DH,
- class VectorType>
+ template <int dim, int spacedim, template <int,int> class DH, typename VectorType>
void
interpolate_to_different_mesh (const DH<dim, spacedim> &dof1,
const VectorType &u1,
- template <int dim, int spacedim,
- template <int,int> class DH,
- class VectorType>
+ template <int dim, int spacedim, template <int,int> class DH, typename VectorType>
void
interpolate_to_different_mesh (const DH<dim, spacedim> &dof1,
const VectorType &u1,
}
- template <int dim, int spacedim,
- template <int,int> class DH,
- class VectorType>
+ template <int dim, int spacedim, template <int,int> class DH, typename VectorType>
void
interpolate_to_different_mesh (const InterGridMap<DH<dim, spacedim> > &intergridmap,
const VectorType &u1,
/**
* Compute the boundary values to be used in the project() functions.
*/
- template <int dim, int spacedim,
- template <int,int> class DH,
- template <int,int> class M_or_MC,
- template <int> class Q_or_QC>
+ template <int dim, int spacedim, template <int,int> class DH,
+ template <int,int> class M_or_MC, template <int> class Q_or_QC>
void project_compute_b_v (const M_or_MC<dim, spacedim> &mapping,
const DH<dim,spacedim> &dof,
const Function<spacedim> &function,
/**
* Generic implementation of the project() function
*/
- template <int dim, int spacedim,
- class VectorType,
- template <int,int> class DH,
- template <int,int> class M_or_MC,
- template <int> class Q_or_QC>
- void do_project (const M_or_MC<dim, spacedim> &mapping,
- const DH<dim,spacedim> &dof,
- const ConstraintMatrix &constraints,
- const Q_or_QC<dim> &quadrature,
- const Function<spacedim> &function,
- VectorType &vec_result,
- const bool enforce_zero_boundary,
- const Q_or_QC<dim-1> &q_boundary,
- const bool project_to_boundary_first)
+ template <int dim, int spacedim, typename VectorType, template <int,int> class DH,
+ template <int,int> class M_or_MC, template <int> class Q_or_QC>
+ void do_project (const M_or_MC<dim, spacedim> &mapping,
+ const DH<dim,spacedim> &dof,
+ const ConstraintMatrix &constraints,
+ const Q_or_QC<dim> &quadrature,
+ const Function<spacedim> &function,
+ VectorType &vec_result,
+ const bool enforce_zero_boundary,
+ const Q_or_QC<dim-1> &q_boundary,
+ const bool project_to_boundary_first)
{
Assert (dof.get_fe().n_components() == function.n_components,
ExcDimensionMismatch(dof.get_fe().n_components(),
}
- template <int dim, class VECTOR, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void project (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
- const ConstraintMatrix &constraints,
- const Quadrature<dim> &quadrature,
- const Function<spacedim> &function,
- VECTOR &vec_result,
- const bool enforce_zero_boundary,
- const Quadrature<dim-1> &q_boundary,
- const bool project_to_boundary_first)
+ const ConstraintMatrix &constraints,
+ const Quadrature<dim> &quadrature,
+ const Function<spacedim> &function,
+ VectorType &vec_result,
+ const bool enforce_zero_boundary,
+ const Quadrature<dim-1> &q_boundary,
+ const bool project_to_boundary_first)
{
do_project (mapping, dof, constraints, quadrature,
function, vec_result,
}
- template <int dim, class VECTOR, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void project (const DoFHandler<dim,spacedim> &dof,
- const ConstraintMatrix &constraints,
- const Quadrature<dim> &quadrature,
- const Function<spacedim> &function,
- VECTOR &vec,
- const bool enforce_zero_boundary,
- const Quadrature<dim-1> &q_boundary,
- const bool project_to_boundary_first)
+ const ConstraintMatrix &constraints,
+ const Quadrature<dim> &quadrature,
+ const Function<spacedim> &function,
+ VectorType &vec,
+ const bool enforce_zero_boundary,
+ const Quadrature<dim-1> &q_boundary,
+ const bool project_to_boundary_first)
{
project(StaticMappingQ1<dim,spacedim>::mapping, dof, constraints, quadrature, function, vec,
enforce_zero_boundary, q_boundary, project_to_boundary_first);
- template <int dim, class VECTOR, int spacedim>
- void project (const hp::MappingCollection<dim, spacedim> &mapping,
- const hp::DoFHandler<dim,spacedim> &dof,
- const ConstraintMatrix &constraints,
- const hp::QCollection<dim> &quadrature,
- const Function<spacedim> &function,
- VECTOR &vec_result,
- const bool enforce_zero_boundary,
- const hp::QCollection<dim-1> &q_boundary,
- const bool project_to_boundary_first)
+ template <int dim, typename VectorType, int spacedim>
+ void project (const hp::MappingCollection<dim, spacedim> &mapping,
+ const hp::DoFHandler<dim,spacedim> &dof,
+ const ConstraintMatrix &constraints,
+ const hp::QCollection<dim> &quadrature,
+ const Function<spacedim> &function,
+ VectorType &vec_result,
+ const bool enforce_zero_boundary,
+ const hp::QCollection<dim-1> &q_boundary,
+ const bool project_to_boundary_first)
{
do_project (mapping, dof, constraints, quadrature,
function, vec_result,
}
- template <int dim, class VECTOR, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void project (const hp::DoFHandler<dim,spacedim> &dof,
- const ConstraintMatrix &constraints,
- const hp::QCollection<dim> &quadrature,
- const Function<spacedim> &function,
- VECTOR &vec,
- const bool enforce_zero_boundary,
- const hp::QCollection<dim-1> &q_boundary,
- const bool project_to_boundary_first)
+ const ConstraintMatrix &constraints,
+ const hp::QCollection<dim> &quadrature,
+ const Function<spacedim> &function,
+ VectorType &vec,
+ const bool enforce_zero_boundary,
+ const hp::QCollection<dim-1> &q_boundary,
+ const bool project_to_boundary_first)
{
project(hp::StaticMappingQ1<dim,spacedim>::mapping_collection,
dof, constraints, quadrature, function, vec,
// faces are points and it is far
// easier to simply work on
// individual vertices
- template <class DH,
- template <int,int> class M_or_MC>
+ template <class DH, template <int,int> class M_or_MC>
static inline
void do_interpolate_boundary_values (const M_or_MC<DH::dimension, DH::space_dimension> &,
const DH &dof,
// dim_, it is clearly less specialized than the 1d function above and
// whenever possible (i.e., if dim==1), the function template above
// will be used
- template <class DH,
- template <int,int> class M_or_MC,
- int dim_>
+ template <class DH, template <int,int> class M_or_MC, int dim_>
static inline
void
do_interpolate_boundary_values (const M_or_MC<DH::dimension, DH::space_dimension> &mapping,
namespace
{
- template <int dim, int spacedim,
- template <int,int> class DH,
- template <int,int> class M_or_MC,
- template <int> class Q_or_QC>
+ template <int dim, int spacedim, template <int,int> class DH,
+ template <int,int> class M_or_MC, template <int> class Q_or_QC>
void
do_project_boundary_values (const M_or_MC<dim, spacedim> &mapping,
const DH<dim, spacedim> &dof,
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_difference (const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_difference (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_value (const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_value (const hp::DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double
point_value (const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double
point_value (const hp::DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_value (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_value (const hp::MappingCollection<dim, spacedim> &mapping,
const hp::DoFHandler<dim,spacedim> &dof,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double
point_value (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double
point_value (const hp::MappingCollection<dim, spacedim> &mapping,
const hp::DoFHandler<dim,spacedim> &dof,
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_gradient (const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_gradient (const hp::DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
Tensor<1, spacedim, typename VectorType::value_type>
point_gradient (const DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
Tensor<1, spacedim, typename VectorType::value_type>
point_gradient (const hp::DoFHandler<dim,spacedim> &dof,
const VectorType &fe_function,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_gradient (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
void
point_gradient (const hp::MappingCollection<dim, spacedim> &mapping,
const hp::DoFHandler<dim,spacedim> &dof,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
Tensor<1, spacedim, typename VectorType::value_type>
point_gradient (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
Tensor<1, spacedim, typename VectorType::value_type>
point_gradient (const hp::MappingCollection<dim, spacedim> &mapping,
const hp::DoFHandler<dim,spacedim> &dof,
- template <class VECTOR>
+ template <typename VectorType>
void
- subtract_mean_value(VECTOR &v,
+ subtract_mean_value(VectorType &v,
const std::vector<bool> &p_select)
{
if (p_select.size() == 0)
Assert(p_select.size() == n,
ExcDimensionMismatch(p_select.size(), n));
- typename VECTOR::value_type s = 0.;
+ typename VectorType::value_type s = 0.;
unsigned int counter = 0;
for (unsigned int i=0; i<n; ++i)
if (p_select[i])
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double
compute_mean_value (const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim,spacedim> &dof,
}
- template <int dim, class VectorType, int spacedim>
+ template <int dim, typename VectorType, int spacedim>
double
compute_mean_value (const DoFHandler<dim,spacedim> &dof,
const Quadrature<dim> &quadrature,
}
- template<class DH, class VECTOR>
+ template<class DH, typename VectorType>
void get_position_vector(const DH &dh,
- VECTOR &vector,
+ VectorType &vector,
const ComponentMask &mask)
{
AssertDimension(vector.size(), dh.n_dofs());