Electrowetting on dielectric
Classes | Public Member Functions | Public Attributes | Protected Types | Protected Member Functions | Protected Attributes
CahnHilliard< dim > Class Template Reference
Inheritance diagram for CahnHilliard< dim >:
Inheritance graph
[legend]
Collaboration diagram for CahnHilliard< dim >:
Collaboration graph
[legend]

List of all members.

Classes

struct  ScratchData

Public Member Functions

 CahnHilliard (const Triangulation< dim > &tria, Material_Parameters &params, const unsigned deg, const TestCase _tc)
virtual ~CahnHilliard ()
void advance ()
double it_error ()

Public Attributes

AsFunction< dim > get_phase
 The phase field.
AsFunction< dim > get_potential
 The chemical potential.
AsFunction< dim > get_old_phase
 The old phase field.

Protected Types

typedef IteratorGroup
< 3, typename hp::DoFHandler
< dim >::active_cell_iterator > 
Iterator

Protected Member Functions

virtual void InitLADataSuffix ()
virtual void AssembleSystem (std::vector< AsFunction< dim > * > &data, const double time, const bool m_theaded=true)
void AssembleCell (const Iterator Its, ScratchData &scratch, typename Problem< dim >::PerTaskData &data)
virtual void ReinitPrec ()
virtual void SolveSystemPreffix ()
virtual void DoSolve (SolverControl &control)
virtual void SetInitialData ()
virtual void PreRefinementPreffix ()
virtual void PostRefinementSuffix (const std::vector< TrilinosWrappers::Vector > &sol_tmp)

Protected Attributes

const material_functionmobility
 Mobility of the phases.
const material_functionpermitivity
 Electric permittivity of the phases.
const material_functiondensity
 Density of the phases.
const GLPotentialW
 The derivative of the Ginzburg Landau potential.
const InterfaceEnergygamma_fs
 The derivative of the interface fee energy.
const double gamma
 Suface tension.
const double delta
 Interface thickness.
const double alpha
 Surface transport parameter.
const TestCase tc
 The test case.
TrilinosWrappers::Vector oldsol
 The old solution.
TrilinosWrappers::Vector olditsol
 Solution at the previous iterative step in the inner loop.

Detailed Description

template<int dim>
class CahnHilliard< dim >

A Cahn Hilliard equation with surface transport terms:

\[ \phi_t + \mathbf{u} \cdot \nabla \phi = \nabla\cdot( M(\phi) \nabla \mu ), \quad \text{ in } \Omega \]

\[ \mu = \gamma \left( \frac1\delta \mathcal{W}'(\phi) - \delta \Delta \phi \right) - \frac12\varepsilon'(\phi)|\nabla V|^2 + \frac12 \rho'(\phi)|\mathbf{u}|^2, \quad \text{ in } \Omega \]

with boundary conditions

\[ \alpha\left( \phi_t + \mathbf{u}_{\boldsymbol{\tau}} \partial_{\boldsymbol{\tau}} \phi \right) + \gamma \left( \Theta_{fs}'(\phi) + \delta \partial_{\mathbf{n}} \phi \right) = 0, \quad M(\phi) \partial_n \mu = 0, \quad \text{on }\Gamma, \]

Definition at line 58 of file CahnHilliard.h.


Member Typedef Documentation

template<int dim>
typedef IteratorGroup<3,typename hp::DoFHandler<dim>::active_cell_iterator> CahnHilliard< dim >::Iterator [protected]

We group the iterators of this problem and the problems this one takes information from as follows:

Definition at line 107 of file CahnHilliard.h.


Constructor & Destructor Documentation

template<int dim>
CahnHilliard< dim >::CahnHilliard ( const Triangulation< dim > &  tria,
Material_Parameters params,
const unsigned  deg,
const TestCase  _tc 
)

Constructor. Copy the provided parameters to the members:

Parameters:
tria: A reference to the mesh where the problem will live.
params: The collection of material parameters.
deg: The polynomial degree for the finite element spaces where the phase field and chemical potential are defined.
_tc: The testcase. Needed to know which initial condition to invoke.
template<int dim>
virtual CahnHilliard< dim >::~CahnHilliard ( ) [virtual]

Clear the data.


Member Function Documentation

template<int dim>
void CahnHilliard< dim >::advance ( )

Advance in time.

template<int dim>
double CahnHilliard< dim >::it_error ( )

compute the iterative error.

template<int dim>
virtual void CahnHilliard< dim >::InitLADataSuffix ( ) [protected, virtual]

What to do after the LA data is initialized

Reimplemented from Problem< dim >.

template<int dim>
virtual void CahnHilliard< dim >::AssembleSystem ( std::vector< AsFunction< dim > * > &  data,
const double  time,
const bool  m_theaded = true 
) [protected, virtual]

Assemble the system matrix and right hand side for the current status of the problem

Parameters:
data: The views of the other problems this one needs information from to be able to assemble.
  1. Voltage.
  2. Velocity.
time: The current time.
m_theaded: flag that indicates if the assembly should be done using multithreading.

Implements Problem< dim >.

template<int dim>
void CahnHilliard< dim >::AssembleCell ( const Iterator  Its,
ScratchData scratch,
typename Problem< dim >::PerTaskData data 
) [protected]

Assemble the local problems.

Parameters:
Its: The cell iterators for the current problem and the data.
scratch: The scratch data.
data: The per task data.
template<int dim>
virtual void CahnHilliard< dim >::ReinitPrec ( ) [protected, virtual]

Does nothing. Needed for compatibility

Implements Problem< dim >.

template<int dim>
virtual void CahnHilliard< dim >::SolveSystemPreffix ( ) [protected, virtual]

What to do before we solve the system

Reimplemented from Problem< dim >.

template<int dim>
virtual void CahnHilliard< dim >::DoSolve ( SolverControl control) [protected, virtual]

Solve the system to given accuracy

Implements Problem< dim >.

template<int dim>
virtual void CahnHilliard< dim >::SetInitialData ( ) [protected, virtual]

Set the initial data

Implements Problem< dim >.

template<int dim>
virtual void CahnHilliard< dim >::PreRefinementPreffix ( ) [protected, virtual]

Used to handle changes in the triangulation: Copy to the solution transfer

Implements Problem< dim >.

template<int dim>
virtual void CahnHilliard< dim >::PostRefinementSuffix ( const std::vector< TrilinosWrappers::Vector > &  sol_tmp) [protected, virtual]

Used to handle changes in the triangulation: Copy from the solution transfer

Implements Problem< dim >.


Member Data Documentation

template<int dim>
const material_function& CahnHilliard< dim >::mobility [protected]

Mobility of the phases.

Definition at line 83 of file CahnHilliard.h.

template<int dim>
const material_function & CahnHilliard< dim >::permitivity [protected]

Electric permittivity of the phases.

Definition at line 83 of file CahnHilliard.h.

template<int dim>
const material_function & CahnHilliard< dim >::density [protected]

Density of the phases.

Definition at line 83 of file CahnHilliard.h.

template<int dim>
const GLPotential& CahnHilliard< dim >::W [protected]

The derivative of the Ginzburg Landau potential.

Definition at line 86 of file CahnHilliard.h.

template<int dim>
const InterfaceEnergy& CahnHilliard< dim >::gamma_fs [protected]

The derivative of the interface fee energy.

Definition at line 87 of file CahnHilliard.h.

template<int dim>
const double CahnHilliard< dim >::gamma [protected]

Suface tension.

Definition at line 88 of file CahnHilliard.h.

template<int dim>
const double CahnHilliard< dim >::delta [protected]

Interface thickness.

Definition at line 88 of file CahnHilliard.h.

template<int dim>
const double CahnHilliard< dim >::alpha [protected]

Surface transport parameter.

Definition at line 88 of file CahnHilliard.h.

template<int dim>
const TestCase CahnHilliard< dim >::tc [protected]

The test case.

Definition at line 91 of file CahnHilliard.h.

template<int dim>
TrilinosWrappers::Vector CahnHilliard< dim >::oldsol [protected]

The old solution.

Definition at line 92 of file CahnHilliard.h.

template<int dim>
TrilinosWrappers::Vector CahnHilliard< dim >::olditsol [protected]

Solution at the previous iterative step in the inner loop.

Definition at line 92 of file CahnHilliard.h.

template<int dim>
AsFunction<dim> CahnHilliard< dim >::get_phase

The phase field.

Definition at line 202 of file CahnHilliard.h.

template<int dim>
AsFunction<dim> CahnHilliard< dim >::get_potential

The chemical potential.

Definition at line 202 of file CahnHilliard.h.

template<int dim>
AsFunction<dim> CahnHilliard< dim >::get_old_phase

The old phase field.

Definition at line 202 of file CahnHilliard.h.


The documentation for this class was generated from the following file: