Revision 381d913240790774432fb3130f59687a10db4df9 authored by Vincent Chabannes on 20 August 2014, 10:14:00 UTC, committed by Vincent Chabannes on 20 August 2014, 10:14:00 UTC
1 parent e69c9ea
Raw File
README.org
#+OPTIONS: LaTeX:nil
#+STYLE: <script type="text/x-mathjax-config">  MathJax.Hub.Config({tex2jax: {inlineMath: [['$','$'], ['\\(','\\)']]}});</script><script type="text/javascript"  src="http://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML"></script>

* QuickStart

This directory contains a quickstart example that solves a Laplacian in 2D using
P1 Lagrange basis functions.

*Note*: you need to have [[http://www.paraview.org][Paraview]] installed to visualize the results of this example.

** Code

Denote $\Omega=[0,1]^{2}$, we consider the following problem, find $u$ such that

$$-\Delta u = f \mbox{ in } \Omega, \quad u=g \mbox{ on } \partial \Omega$$,

Here is the Feel++ code to solve this problem
#+BEGIN_SRC C++
#include <feel/feel.hpp>

int main(int argc, char**argv )
{
    using namespace Feel;
	Environment env( _argc=argc, _argv=argv,
                     _desc=feel_options(),
                     _about=about(_name="qs_laplacian",
                                  _author="Feel++ Consortium",
                                  _email="feelpp-devel@feelpp.org"));

    auto mesh = unitSquare();
    auto Vh = Pch<1>( mesh );
    auto u = Vh->element();
    auto v = Vh->element();

    auto l = form1( _test=Vh );
    l = integrate(_range=elements(mesh),
                  _expr=id(v));

    auto a = form2( _trial=Vh, _test=Vh );
    a = integrate(_range=elements(mesh),
                  _expr=gradt(u)*trans(grad(v)) );
    a+=on(_range=boundaryfaces(mesh), _rhs=l, _element=u,
          _expr=constant(0.) );
    a.solve(_rhs=l,_solution=u);

    auto e = exporter( _mesh=mesh, _name="qs_laplacian" );
    e->add( "u", u );
    e->save();
    return 0;
}
#+END_SRC

*Note:* =Pch<N>( mesh )= generates the function space

$$V_h = P_{c,h} = { v \in C^0(\Omega) such that v_{|K} \in P_N(K), \forall K \in
\mathcal{T}_h }.$$

where $\mathcal{T}_h=\{K\}$ is the mesh (a set of cells).

** Compilation

To compile the quickstart example, just type

#+BEGIN_SRC shell
make
#+END_SRC

It generates the executable =feelpp_qs_laplacian= (the Feel++ QuickStart Laplacian).

** Execution

To run the example, just type
#+BEGIN_SRC shell
./feelpp_qs_laplacian
#+END_SRC

The example above is run using 1 core. There won't be any output in the
console. Type =ls= (Unix) or =dir= (Windows) and see that several files have
been generated.

 - =square.geo= Gmsh geo file describing the geometry of the problem
 - =square.msh= Gmsh mesh file generated by Gmsh from the .geo file
 - =qs_laplacian-1.sos= (Ensight Server of Server) sos data file to be loaded by
   Paraview, describing the simulation for 1 core
 - =qs_laplacian-1_0.case= (Ensight Case) file loaded by Paraview
 - =qs_laplacian-1_0.geo001= Ensight mesh file loaded by Paraview
 - =qs_laplacian.pid-1_0.001= Ensight data file loaded by Paraview containing the
   process id of mesh element
 - =qs_laplacian.u-1_0.001= Ensight data file loaded by Paraview containing the
   solution of the example

You can also run the example in parallel, if you have, say 8 cores, then type

#+BEGIN_SRC shell
mpirun -np 8 ./feelpp_qs_laplacian
#+END_SRC

Load the file =qs_laplacian-8.sos= in Paraview, it describes the data set to be
loaded by Paraview for the simulation run on 8 cores.

** Playing further

You can easily change the dimension of the problem, for example replace
#+BEGIN_SRC c++
auto mesh = unitSquare();
#+END_SRC
by
#+BEGIN_SRC c++
auto mesh = unitCube();
#+END_SRC
and the example now solves the Laplacian in the unit cube [0,1]^3.


You can also change the polynomial order approximation,
for example replace
#+BEGIN_SRC C++
auto Vh = Pch<1>( mesh );
#+END_SRC
into
#+BEGIN_SRC C++
auto Vh = Pch<3>( mesh );
#+END_SRC
to get order 3 approximation.

back to top