The next step is to construct a function space over the mesh. The source code is available in myexpression.cpp
.
Step by step explanations
We start by loading a Mesh in 2D
auto mesh = loadMesh( _mesh=new Mesh<Simplex<2>> );
then we define some expression through the command line or config file: g
is a scalar field and f
is a vector field
auto g = expr(soption(_name="functions.g"));
std::cout << "g=" << g << std::endl;
auto f = expr<2,1>(soption(_name="functions.f"));
std::cout << "f=" << f << std::endl;
here is an example how to enter them, more are available below feelpp_doc_myexpression
–functions.g="x*y:x:y" –functions.f="{sin(pi*x),cos(pi*y)}:x:y"
then we compute the gradient of g
and f
auto grad_g=grad<2>(g);
auto grad_f=grad(f);
std::cout << "grad(g)=" << grad_g << std::endl;
std::cout << "grad(f)=" << grad_f << std::endl;
Notice that template argument are given to grad
to specify the shape of the gradient: in the case of \(\nabla g\) it is \(1\times2\) and \(\nabla f\) \(2\times 2\) since we are in 2D.
then we compute the laplacian of g
and f
auto laplacian_g=laplacian(g);
std::cout << "laplacian(g)=" << laplacian_g << std::endl;
auto laplacian_f=laplacian(f);
std::cout << "laplacian(f)=" << laplacian_f << std::endl;
then we compute the divergence of f
auto div_f=div(f);
std::cout << "div(f)=" << div_f << std::endl;
and the curl of f
auto curl_f=curl(f);
std::cout << "curl(f)=" << curl_f << std::endl;
Finally we evaluate these expression at one point given by the option x
and y
std::cout << "Evaluation at (" << doption("x") << "," << doption("y") << "):" << std::endl;
std::cout << " g(x,y)=" << g.evaluate() << std::endl;
std::cout << " f(x,y)=" << f.evaluate() << std::endl;
std::cout << "Gradient:\n";
std::cout << " grad(g)(x,y)=" << grad_g.evaluate() << std::endl;
std::cout << " grad(f)(x,y)=" << grad_f.evaluate() << std::endl;
std::cout << "Divergence:\n";
std::cout << " div(f)(x,y)=" << div_f.evaluate() << std::endl;
std::cout << "Curl:\n";
std::cout << " curl(f)(x,y)=" << curl_f.evaluate() << std::endl;
std::cout << "Laplacian:\n";
std::cout << "laplacian(g)(x,y)=" << laplacian_g.evaluate() << std::endl;
std::cout << "laplacian(f)(x,y)=" << laplacian_f.evaluate() << std::endl;
Some results
We start with the following function \(g=1\) and \(f=(1,1)\).
shell> ./feelpp_doc_myexpression --functions.g=1:x:y --functions.f="{1,1}:x:y
g=1
f={x,-y}
grad(g)=[[0]]
grad(f)=[[1,0],[0,-1]]
laplacian(g)=[[0]]
laplacian(f)=[[0],[0]]
div(f)=[[0]]
curl(f)=[[0]]
Evaluation at (0,0):
g(x,y)=1
f(x,y)= 0
-0
Gradient:
grad(g)(x,y)= 0 -0
grad(f)(x,y)= 1 0
0 -1
Divergence:
div(f)(x,y)=0
Curl:
curl(f)(x,y)=-3.14159
Laplacian:
laplacian(g)(x,y)=0
laplacian(f)(x,y)=0
0
The symbolic calculus system worked as expected.
Complete code
The complete code reads as follows
#include <feel/feelcore/environment.hpp>
#include <feel/feelfilters/loadmesh.hpp>
#include <feel/feelvf/ginac.hpp>
using namespace Feel;
int main(int argc, char**argv )
{
Environment env( _argc=argc, _argv=argv,
_about=about(_name="myexpression",
_author="Feel++ Consortium",
_email="feelpp-devel@feelpp.org"));
auto mesh = loadMesh(_mesh=new Mesh<Simplex<2>>);
auto g = expr(soption(_name="functions.g"));
std::cout << "g=" << g << std::endl;
auto f = expr<2,1>(soption(_name="functions.f"));
std::cout << "f=" << f << std::endl;
auto grad_g=grad<2>(g);
auto grad_f=grad(f);
std::cout << "grad(g)=" << grad_g << std::endl;
std::cout << "grad(f)=" << grad_f << std::endl;
auto laplacian_g=laplacian(g);
std::cout << "laplacian(g)=" << laplacian_g << std::endl;
auto laplacian_f=laplacian(f);
std::cout << "laplacian(f)=" << laplacian_f << std::endl;
auto div_f=div(f);
std::cout << "div(f)=" << div_f << std::endl;
auto curl_f=curl(f);
std::cout << "curl(f)=" << curl_f << std::endl;
std::cout << "Evaluation at (" << doption("x") << "," << doption("y") << "):" << std::endl;
std::cout << " g(x,y)=" << g.evaluate() << std::endl;
std::cout << " f(x,y)=" << f.evaluate() << std::endl;
std::cout << "Gradient:\n";
std::cout << " grad(g)(x,y)=" << grad_g.evaluate() << std::endl;
std::cout << " grad(f)(x,y)=" << grad_f.evaluate() << std::endl;
std::cout << "Divergence:\n";
std::cout << " div(f)(x,y)=" << div_f.evaluate() << std::endl;
std::cout << "Curl:\n";
std::cout << " curl(f)(x,y)=" << curl_f.evaluate() << std::endl;
std::cout << "Laplacian:\n";
std::cout << "laplacian(g)(x,y)=" << laplacian_g.evaluate() << std::endl;
std::cout << "laplacian(f)(x,y)=" << laplacian_f.evaluate() << std::endl;
}
to compile just type
make feelpp_doc_myexpression