Logo
Finite Element Embedded Library and Language in C++
Feel++ Feel++ on Github Feel++ on Travis-CI Feel++ on Twitter Feel++ on YouTube Feel++ community
 All Classes Namespaces Files Functions Variables Typedefs Pages
Feel::DAR< Dim > Class Template Reference

Detailed Description

template<int Dim>
class Feel::DAR< Dim >

DAR Solver using continuous approximation spaces solve \( -\Delta u = f\) on \(\Omega\) and \(u= g\) on \(\Gamma\)

Template Parameters
Dimthe geometric dimension of the problem (e.g. Dim=1, 2 or 3)

Inherits Simget.

Public Types

typedef bases< Lagrange< Order,
Scalar > > 
basis_type
 the basis type of our approximation space
 
typedef Simplex< Dim > convex_type
 geometry entities type composing the mesh, here Simplex in Dimension Dim of Order 1
 
typedef space_type::element_type element_type
 an element type of the approximation function space
 
typedef boost::shared_ptr
< export_type
export_ptrtype
 the exporter factory (shared_ptr<> type)
 
typedef Exporter< mesh_typeexport_type
 the exporter factory type
 
typedef boost::shared_ptr
< mesh_type
mesh_ptrtype
 mesh shared_ptr<> type
 
typedef Mesh< convex_typemesh_type
 mesh type
 
typedef boost::shared_ptr
< space_type
space_ptrtype
 the approximation function space type (shared_ptr<> type)
 
typedef FunctionSpace
< mesh_type, basis_type
space_type
 the approximation function space type
 
typedef double value_type
 numerical type is double
 

Public Member Functions

 DAR ()
 
void run ()
 

Static Public Attributes

static const uint16_type Order = 2
 Polynomial order \(P_2\).
 

Constructor & Destructor Documentation

template<int Dim>
Feel::DAR< Dim >::DAR ( )
inline

Constructor

Member Function Documentation

template<int Dim>
void Feel::DAR< Dim >::run ( )

The function space and some associated elements(functions) are then defined

*/
space_ptrtype Xh = space_type::New( mesh );
element_type u( Xh, "u" );
element_type v( Xh, "v" );
element_type gproj( Xh, "v" );

define \(g\) the expression of the exact solution and \(f\) the expression of the right hand side such that \(g\) is the exact solution

*/
//# marker1 #
auto g = sin( pi*Px() )*cos( pi*Py() )*cos( pi*Pz() );
gproj = vf::project( Xh, elements( mesh ), g );
auto grad_g = vec( +pi*cos( pi*Px() )*cos( pi*Py() )*cos( pi*Pz() ),
-pi*sin( pi*Px() )*sin( pi*Py() )*cos( pi*Pz() ) );
auto beta = vec( cst( bx ),cst( by ) );
auto f = ( pi*pi*Dim*epsilon*g+trans( grad_g )*beta + mu*g );
//# endmarker1 #

Construction of the right hand side. F is the vector that holds the algebraic representation of the right habd side of the problem

*/
//# marker2 #
auto F = backend( _vm=this->vm() )->newVector( Xh );
form1( _test=Xh, _vector=F, _init=true ) =
integrate( _range=elements( mesh ), _expr=f*id( v ) );
//# endmarker2 #
if ( weak_dirichlet )
{
//# marker41 #
#if 1
form1( _test=Xh, _vector=F ) +=
integrate( _range=boundaryfaces( mesh ),
_expr=g*( -epsilon*grad( v )*vf::N()+penaldir*id( v )/hFace() ) );
#endif
//# endmarker41 #
}

create the matrix that will hold the algebraic representation of the left hand side

*/
size_type pattern=Pattern::COUPLED;
if ( stab )
pattern = Pattern::COUPLED|Pattern::EXTENDED;
auto D = backend()->newMatrix( _test=Xh, _trial=Xh, _pattern=pattern );

assemble $ u v$

*/
form2( _test=Xh, _trial=Xh, _matrix=D ) =
integrate( _range=elements( mesh ),
_expr=( epsilon*gradt( u )*trans( grad( v ) )+( gradt( u )*beta )*id( v )+mu*idt( u )*id( v ) ) );

weak dirichlet conditions treatment for the boundaries marked 1 and 3

  1. assemble \(\int_{\partial \Omega} -\nabla u \cdot \mathbf{n} v\)
  2. assemble \(\int_{\partial \Omega} -\nabla v \cdot \mathbf{n} u\)
  3. assemble \(\int_{\partial \Omega} \frac{\gamma}{h} u v\)
*/
//# marker10 #
#if 1
form2( _test=Xh, _trial=Xh, _matrix=D ) +=
integrate( _range=boundaryfaces( mesh ),
_expr= ( -( epsilon*gradt( u )*vf::N() )*id( v )
-( epsilon*grad( v )*vf::N() )*idt( u )
+penaldir*id( v )*idt( u )/hFace() ) );
#else
form2( _test=Xh, _trial=Xh, _matrix=D ) +=
integrate( _range=boundaryfaces( mesh ),
_expr= ( -( epsilon*gradt( u )*vf::N() )*id( v ) ) );
#endif
//# endmarker10 #

strong(algebraic) dirichlet conditions treatment for the boundaries marked 1 and 3

  1. first close the matrix (the matrix must be closed first before any manipulation )
  2. modify the matrix by cancelling out the rows and columns of D that are associated with the Dirichlet dof
*/
//# marker5 #
#if 1
form2( _test=Xh, _trial=Xh, _matrix=D ) +=
on( _range=boundaryfaces( mesh ),_element=u, _rhs=F, _expr=g );
#endif
//# endmarker5 #

solve the system

*/
//# marker6 #
backend( _rebuild=true,_vm=this->vm() )->solve( _matrix=D, _solution=u, _rhs=F );
//# endmarker6 #

compute the

*/
//# marker7 #
LOG(INFO) << "||error||_L2=" << normL2( _range=elements(mesh), _expr=idv(u)-g ) << "\n";
//# endmarker7 #

save the results

*/
element_type e( Xh, "e" );
e = vf::project( Xh, elements( mesh ), g );
export_ptrtype exporter( export_type::New( this->vm(),
( boost::format( "%1%-%2%-%3%" )
% this->about().appName()
% shape
% Dim ).str() ) );
if ( exporter->doExport() )
{
LOG(INFO) << "exportResults starts\n";
exporter->step( 0 )->setMesh( mesh );
exporter->step( 0 )->add( "u", u );
exporter->step( 0 )->add( "g", e );
exporter->save();
LOG(INFO) << "exportResults done\n";
}

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