Feel++ 0.91.0

Manipulating Function Spaces

Defining function spaces and functions

typedef

    //# marker1 #
    typedef FunctionSpace<mesh_type, bases<Lagrange<0,Scalar, Discontinuous> > > p0_space_type;
    typedef typename p0_space_type::element_type p0_element_type;
    //# endmarker1 #

P0 FunctionSpace typedef

    //# marker2 #
    typedef bases<Lagrange<Order,Scalar,Continuous> > basis_type;

    typedef FunctionSpace<mesh_type, basis_type> space_type;
    typedef boost::shared_ptr<space_type> space_ptrtype;
    typedef typename space_type::element_type element_type;
    //# endmarker2 #

Now we turn to the instantiation of the function space $ X_h $ of type functionspace_type. We have already instantiated a mesh (see this Section) of type mesh_ptrtype. The code looks like this

    //# marker31 #
    mesh_ptrtype mesh = createGMSHMesh( _mesh=new mesh_type,
                                        _update=MESH_CHECK|MESH_UPDATE_FACES|MESH_UPDATE_EDGES|MESH_RENUMBER,
                                        _desc=domain( _name= (boost::format( "%1%-%2%-%3%" ) % shape % Dim % Order).str() ,
                                                      _shape=shape,
                                                      _dim=Dim,
                                                      _order=Order,
                                                      _h=X[0] ) );

    //# endmarker31 #

Then we instantiate $ X_h$ using space_type::New() static member function and obtain a functionspace_ptrtype.

    space_ptrtype Xh = space_type::New( mesh );

We are now able to instantiate elements of the FunctionSpace $ X_h$. Two ways are presented, one with the auto keyword allowing to infer automatically the type of $ X_h$ elements

    auto u = Xh->element( "u" );

and one knowing the actual type of the element

    element_type v( Xh, "v" );

Note:
The strings "u" and "v" are the name we give to these elements of $ X_h$.

Using functions spaces and functions

Projection

First, we define some mathematical expression/functions

\begin{eqnarray} g(x,y,z) &=& \sin(\pi x/2)\ \cos(\pi y/2)\ \cos(\pi*z/2) \\ f(x,y,z) &=& (1-x^2)\ (1-y^2)\ (1-z^2)\ (x^2+y^2+z^2)^{\alpha/2.0},\quad \alpha=3 \end{eqnarray}

They are implemented using the new auto C++ keyword to infer the type of expression automatically.

    //# marker4 #
    auto g = sin(2*pi*Px())*cos(2*pi*Py())*cos(2*pi*Pz());
    auto f = (1-Px()*Px())*(1-Py()*Py())*(1-Pz()*Pz())*pow(trans(vf::P())*vf::P(),(alpha/2.0));
    //# endmarker4 #

Then we build the interpolant (Lagrange interpolant in this case since we chose Lagrange basis function), by calling the vf::project() function which can be applied on all or parts of the mesh thanks to the mesh iterators such as elements(mesh) or markedelements(mesh,marker). The return object is the interpolant of the function in the space $ X_h $ given as an argument to vf::project().

    //# marker5 #
    u = vf::project( Xh, elements(mesh), g );
    v = vf::project( Xh, elements(mesh), f );
    w = vf::project( Xh, elements(mesh), idv(u)-g );
    //# endmarker5 #

Norm computation

It is easy to compute norms

    //# marker6 #
    double L2g2 = integrate( elements(mesh), g*g ).evaluate()(0,0);
    double L2uerror2 = integrate( elements(mesh), (idv(u)-g)*(idv(u)-g) ).evaluate()(0,0);
    Log() << "||u-g||_0=" << math::sqrt( L2uerror2/L2g2 ) << "\n";
    double L2f2 = integrate( elements(mesh), f*f ).evaluate()(0,0);
    double L2verror2 = integrate( elements(mesh), (idv(v)-f)*(idv(v)-f) ).evaluate()(0,0);
    Log() << "||v-f||_0=" << math::sqrt( L2verror2/L2f2 ) << "\n";
    //# endmarker6 #

Interpolation

Exporting to Paraview or Gmsh

projection

    //# marker7 #
    exporter = export_ptrtype( export_type::New( this->vm(), (boost::format( "%1%-%2%-%3%-%4%" ) % this->about().appName() % shape % Dim % Order).str() ) );
    exporter->step(0)->setMesh( mesh );
    exporter->step(0)->add( "g", u );
    exporter->step(0)->add( "u-g", w );
    exporter->step(0)->add( "f", v );
    exporter->save();
    //# endmarker7 #

Execution

Execution on a Simplex

We execute this example on a simplex domain an export to the Gmsh format:

feel_doc_myfunctionspace --shape="simplex"--nochdir --exporter-format=gmsh

the output log of the execution of this example gives

Here are the graphical outputs on the d-simplex, d=1,2,3:

myfunctionspace-simplex-1-f.png

f on the Line

myfunctionspace-simplex-2-f.png

f on the Triangle

myfunctionspace-simplex-3-f.png

f on the Tetrahedron

myfunctionspace-simplex-1-g.png

g on the Line

myfunctionspace-simplex-2-g.png

g on the Triangle

myfunctionspace-simplex-3-g.png

g on the Tetrahedron

Execution on a Hypercube

We execute this example on a hypercube domain an export to the Gmsh format:

feel_doc_myfunctionspace --shape="hypercube"--nochdir --exporter-format=gmsh

the output log of the execution of this example gives

Here are the graphical outputs on the d-hypercube, d=1,2,3:

myfunctionspace-hypercube-1-f.png

f plot on Line

myfunctionspace-hypercube-2-f.png

f plot on Unit Square

myfunctionspace-hypercube-3-f.png

f plot on Unit Cube

myfunctionspace-hypercube-1-g.png

g plot on Line

myfunctionspace-hypercube-2-g.png

g plot on Unit Square

myfunctionspace-hypercube-3-g.png

g plot on Unit Cube