dune-pdelab 2.7-git
Loading...
Searching...
No Matches
recipe-integrating-grid-functions.cc

See explanation at Integrating grid functions

#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <iostream>
#include <dune/common/filledarray.hh>
#include <dune/common/parallel/mpihelper.hh>
#include <dune/grid/yaspgrid.hh>
#include <dune/pdelab.hh>
int main(int argc, char** argv)
{
try{
// Initialize Mpi
Dune::MPIHelper::instance(argc, argv);
// need a grid in order to test grid functions
constexpr unsigned int dim = 2;
Dune::FieldVector<double,dim> L(1.0);
std::array<int,dim> N(Dune::filledArray<dim,int>(64));
typedef Dune::YaspGrid<dim> Grid;
Grid grid(L,N);
// [Defining an analytic grid function]
auto analyticFunction = Dune::PDELab::makeGridFunctionFromCallable (grid.leafGridView(), [&](const auto& i, const auto& x){
return exp(-(x*x));
});
// [Compute integral]
auto integral = Dune::PDELab::integrateGridFunction(analyticFunction,10);
// [Sum for parallel case]
integral = grid.leafGridView().comm().sum(integral);
std::cout << "Integral: " << integral << std::endl;
return 0;
}
catch (Dune::Exception &e){
std::cerr << "Dune reported error: " << e << std::endl;
return 1;
}
catch (...){
std::cerr << "Unknown exception thrown!" << std::endl;
return 1;
}
}
int main(int argc, char **argv)
Definition recipe-blocking.cc:42
static const int dim
Definition adaptivity.hh:84
GF::Traits::RangeType integrateGridFunction(const GF &gf, unsigned qorder=1)
Integrate a GridFunction.
Definition functionutilities.hh:51
WrapperConformingToGridFunctionInterface makeGridFunctionFromCallable(const GV &gv, const F &f)
Create a GridFunction adapter from a callable.
Definition callableadapter.hh:113