Electrowetting on dielectric
|
Classes | |
struct | PerTaskData |
Public Member Functions | |
Problem (const Triangulation< dim > &tria, const unsigned _n_transfers, const double _eps, bool _prec_needs_reinit=false, const unsigned _update_prec=10) | |
virtual | ~Problem () |
void | set_dt (const double _dt) |
void | init () |
virtual void | solve (std::vector< AsFunction< dim > * > &data, const double time) |
const hp::FECollection< dim > & | get_fe () const |
unsigned | size () const |
Protected Member Functions | |
void | SetupDoFs () |
virtual void | SetupDoFsSuffix () |
void | InitLAData () |
virtual void | InitLADataSuffix () |
virtual void | AssembleSystem (std::vector< AsFunction< dim > * > &data, const double time, const bool m_threaded=true)=0 |
void | CopyToGlob (const PerTaskData &data) |
void | SolveSystem () |
virtual void | SolveSystemPreffix () |
virtual void | ReinitPrec ()=0 |
virtual void | DoSolve (SolverControl &control)=0 |
virtual void | SolveSystemSuffix () |
virtual void | SetInitialData ()=0 |
void | PreRefinement () |
virtual void | PreRefinementPreffix ()=0 |
void | PostRefinement () |
virtual void | PostRefinementSuffix (const std::vector< TrilinosWrappers::Vector > &sol_tmp)=0 |
bool | isNeeded (const typename hp::DoFHandler< dim >::active_cell_iterator &c, const unsigned face) const |
Protected Attributes | |
const Triangulation< dim > & | tri |
Reference to the triangulation. | |
bool | mesh_has_changed |
Flag to indicate that the mesh changed so that we can reinit the prec. | |
hp::DoFHandler< dim > | dh |
The DoF handler object. | |
ConstraintMatrix | constraints |
The hanging node (and other) constraints. | |
hp::FECollection< dim > | fe |
the finite element | |
unsigned | n_dofs |
number of degrees of freedom | |
double | dt |
the time-step | |
unsigned | steps |
the number of steps | |
TrilinosWrappers::Vector | sol |
The solution to the problem. | |
TrilinosWrappers::Vector | rhs |
The RHS. | |
TrilinosWrappers::SparseMatrix | K |
The system matrix. | |
const double | eps |
Tolerance for the iterative solver. | |
const unsigned | update_prec |
How often to update the preconditoner. | |
const bool | prec_needs_reinit |
Flag to indicate if the preconditioner needs reinitialization. | |
const unsigned | n_of_transfers |
Indicates how many vectors we want to transfer on each mesh. | |
std::vector < TrilinosWrappers::Vector > | x_sol |
The vectors that are going to be transferred. | |
SolutionTransfer< dim, TrilinosWrappers::Vector, hp::DoFHandler< dim > > | transfer |
The transferring mechanism. |
Once a problem is discretized in time and space, its solution basically amounts to defining a finite element space, assembling a matrix and solving a linear system. This class is an attempt at maximizing the reuse of the common parts of this process.
The way we try to maximize this is by wrapping around all the common features of this procedures in this class. There is, however, some caveats. It is not possible -- or at least we do not know how -- to fully abstract some parts of this process. For instace, different solvers might be needed depending on whether the problem at hand is symmetric or not. For this reason, we chose to set these as pure virtual functions and the child classes should implement them.
Problem< dim >::Problem | ( | const Triangulation< dim > & | tria, |
const unsigned | _n_transfers, | ||
const double | _eps, | ||
bool | _prec_needs_reinit = false , |
||
const unsigned | _update_prec = 10 |
||
) |
Copy the data to members and/or initialize them
tria | : A reference to the triangulation where this problem is defined |
_n_transfers | : We want to be able to use mesh adaptation. Since most likely our problem is time dependent we want to be able to keep the solution when doing this. In addition, there might be some auxiliary quantities that we want to keep as well. This will indicate how many of these extra quantities we have. |
_eps | : The tolerance for the iterative solver. |
_prec_needs_reinit | : A flag that indicates if the matrix of this problem will change during the existence of the object. If it does then we will need to reinitialize the preconditioner every once in a while. |
_update_prec | : If the preconditioner needs reinitialization, then this indicates what exactly we mean by "once in a while". |
Initialize the problem. This requires:
virtual void Problem< dim >::solve | ( | std::vector< AsFunction< dim > * > & | data, |
const double | time | ||
) | [virtual] |
We assemble and solve the system.
data | : The views from other problems this one uses as data or coefficients. |
time | : The current time |
const hp::FECollection<dim>& Problem< dim >::get_fe | ( | ) | const |
Returns a constant reference to the used FiniteElement object
We distribute the degrees of freedom, set the number of degrees of freedom and make the hanging node constraints.
virtual void Problem< dim >::SetupDoFsSuffix | ( | ) | [protected, virtual] |
It is possible that some constraints must be added to the problem at hand. For instance, if the problem has Dirichlet or Neumann boundary conditions. It is also possible that something else needs to be done. We will call this method just from SetupDoFs() just before closing the ConstraintMatrix so that the child classes can do what they want there.
Reimplemented in Velocity< dim >, Voltage< dim >, and Pressure< dim >.
void Problem< dim >::InitLAData | ( | ) | [protected] |
We initialize the members that are related to the linear algebra part of the solution process. That is, we create the sparsity pattern, which in turn is needed to initialize the system matrix. Then we set the solution and right hand side vectors to the correct size.
virtual void Problem< dim >::InitLADataSuffix | ( | ) | [protected, virtual] |
In a particular problem there might be more linear algebra related objects that need initialization. For instace if we are keeping extrapolations of the solution or the solution at previous time-steps. This cannot be known in advance and so, we add a suffix that will be called at the end InitLAData() and that the child classes will overload.
Reimplemented in Velocity< dim >, CahnHilliard< dim >, and Pressure< dim >.
virtual void Problem< dim >::AssembleSystem | ( | std::vector< AsFunction< dim > * > & | data, |
const double | time, | ||
const bool | m_threaded = true |
||
) | [protected, pure virtual] |
Assemble the system. There is a lot of common structure in the assembly of the system, since the only thing that is really problem-specific is the expression to assemble the local matrices and vectors. However, due to the structure used by the deal.II parallel routines, i.e. division in Scratch and PerTask data, this method needs to be made pure virtual, because there is now way of knowing beforehand what kind of scratch data might be needed for the assembly of the local problems.
data | : The AsFunction objects that the problem uses as coefficients and or right hand side. |
time | : The current time |
m_threaded | : Indicates if the assembly routine must be performed in parallel |
Implemented in Velocity< dim >, CahnHilliard< dim >, Voltage< dim >, Pressure< dim >, and Charge< dim >.
void Problem< dim >::CopyToGlob | ( | const PerTaskData & | data | ) | [protected] |
This IS the assembly procedure. That is, it adds the local matrices and right hand sides to the global ones, taking into account the constraints in the system.
data | : This bundles the local matrix, vector and the local to global map. |
void Problem< dim >::SolveSystem | ( | ) | [protected] |
Do all the bookkeeping needed before solving the system, actually solve it and do all the bookkeeping that is needed after the system is solved. In other words, this method does the following:
virtual void Problem< dim >::SolveSystemPreffix | ( | ) | [protected, virtual] |
What to do before solving the system. For instance we might need to save the old solution.
Reimplemented in Velocity< dim >, and CahnHilliard< dim >.
virtual void Problem< dim >::ReinitPrec | ( | ) | [protected, pure virtual] |
Reinit the preconditioner. The type of preconditioner to be used must be problem specific, and so it should be its initialization. This justifies this function being pure virtual.
Implemented in Velocity< dim >, CahnHilliard< dim >, Voltage< dim >, Pressure< dim >, and Charge< dim >.
virtual void Problem< dim >::DoSolve | ( | SolverControl & | control | ) | [protected, pure virtual] |
Call the actual solver. The solution procedure must be problem specific, so this is pure virtual
Implemented in Velocity< dim >, CahnHilliard< dim >, Voltage< dim >, Pressure< dim >, and Charge< dim >.
virtual void Problem< dim >::SolveSystemSuffix | ( | ) | [protected, virtual] |
What to do after the system is solved. For instance, we might need to update extrapolations.
Reimplemented in Pressure< dim >.
virtual void Problem< dim >::SetInitialData | ( | ) | [protected, pure virtual] |
Set the initial data
Implemented in Velocity< dim >, CahnHilliard< dim >, Voltage< dim >, Pressure< dim >, and Charge< dim >.
void Problem< dim >::PreRefinement | ( | ) | [protected] |
Executed before the refinement
virtual void Problem< dim >::PreRefinementPreffix | ( | ) | [protected, pure virtual] |
Since we do not know how to fill x_sol, we deferr it to the child classes
Implemented in Velocity< dim >, CahnHilliard< dim >, Voltage< dim >, Pressure< dim >, and Charge< dim >.
void Problem< dim >::PostRefinement | ( | ) | [protected] |
Executed after the refinement. We call SetupDofs() and InitLAData(). Then we check if we are at the beginning of time. If we are we call SetInitialData(), otherwise we transfer the vectors to the new mesh.
virtual void Problem< dim >::PostRefinementSuffix | ( | const std::vector< TrilinosWrappers::Vector > & | sol_tmp | ) | [protected, pure virtual] |
We do not know where x_sol came from, so its transfer is deferred to the child classes
Implemented in Velocity< dim >, CahnHilliard< dim >, Voltage< dim >, Pressure< dim >, and Charge< dim >.
bool Problem< dim >::isNeeded | ( | const typename hp::DoFHandler< dim >::active_cell_iterator & | c, |
const unsigned | face | ||
) | const [protected] |
In our application, special care must be taken at the interface between the fluid and the plates. This method determines wether a face is there.
c | : The cell |
face | : The face |
const Triangulation<dim>& Problem< dim >::tri [protected] |
bool Problem< dim >::mesh_has_changed [protected] |
hp::DoFHandler<dim> Problem< dim >::dh [protected] |
ConstraintMatrix Problem< dim >::constraints [protected] |
hp::FECollection<dim> Problem< dim >::fe [protected] |
TrilinosWrappers::Vector Problem< dim >::sol [protected] |
TrilinosWrappers::Vector Problem< dim >::rhs [protected] |
TrilinosWrappers::SparseMatrix Problem< dim >::K [protected] |
const unsigned Problem< dim >::update_prec [protected] |
const bool Problem< dim >::prec_needs_reinit [protected] |
const unsigned Problem< dim >::n_of_transfers [protected] |
std::vector<TrilinosWrappers::Vector> Problem< dim >::x_sol [protected] |
SolutionTransfer<dim, TrilinosWrappers::Vector, hp::DoFHandler<dim> > Problem< dim >::transfer [protected] |