Report a Bug · Request a Feature · Ask a Question

The use of mesh adaptation methods in numerical simulation allows to drastically reduce the memory footprint and the computational costs. There are different kinds of methods: AMR patch-based, AMR cell-based, multiresolution cell-based or point-based, ...

Different open source software is available to the community to manage mesh adaptation: AMReX for patch-based AMR, p4est and pablo for cell-based adaptation.

The strength of samurai is that it allows to implement all the above mentioned mesh adaptation methods from the same data structure. The mesh is represented as intervals and a set algebra allows to efficiently search for subsets among these intervals. Samurai also offers a flexible and pleasant interface to easily implement numerical methods.

Table of Contents

• Get started
  • The advection equation
  • The projection operator
  • There's more
• Features
• Installation
  • From conda
  • From Conan Center
  • From source
• Get help
• Project assistance
• Contributing
• License

Get started

In this section, we propose two examples: the first one solves a 2D advection equation with mesh adaptation using multiresolution, the second one shows the use of set algebra on intervals.

The advection equation

We want to solve the 2D advection equation given by

$$ \partial_t u + a \cdot \nabla u = 0 \; \text{in} \; [0, 1]\times [0, 1] $$

with homogeneous Dirichlet boundary conditions and $a = (1, 1)$. The initial solution is given by

$$ u_0(x, y) = \left\{ \begin{align*} 1 & \; \text{in} \; [0.4, 0.6]\times [0.4, 0.6], \\ 0 & \; \text{elsewhere}. \end{align*} \right. $$

To solve this equation, we use the well known upwind scheme.

The following steps describe how to solve this problem with samurai. It is important to note that these steps are generally the same whatever the equations we want to solve.

• Define the configuration of the problem

  constexpr size_t dim = 2;
  using Config = samurai::MRConfig<dim>;
  std::size_t min_level = 2, max_level = 8;

• Create the Cartesian mesh

  const samurai::Box<double, dim> box({0., 0.}, {1., 1.});
  samurai::MRMesh<Config> mesh(box, min_level, max_level);

• Create the field on this mesh

  auto u = samurai::make_field<double, 1>("u", mesh);
  samurai::make_bc<samurai::Dirichlet<1>>(u, 0.);

• Initialization of this field

  samurai::for_each_cell(mesh, [&](const auto& cell)
  {
      double length = 0.2;
      if (xt::all(xt::abs(cell.center() - 0.5) <= 0.5*length))
      {
          u[cell] = 1;
      }
  });

• Create the adaptation method

  auto MRadaptation = samurai::make_MRAdapt(u);

• Time loop

  double dx = mesh.cell_length(max_level);
  double dt = 0.5*dx;
  auto unp1 = samurai::make_field<double, 1>("u", mesh);
  
  // Time loop
  for (std::size_t nite = 0; nite < 50; ++nite)
  {
      // adapt u
      MRadaptation(1e-4, 2);
  
      // update the ghosts used by the upwind scheme
      samurai::update_ghost_mr(u);
  
      // upwind scheme
      samurai::for_each_interval(mesh, [&](std::size_t level, const auto& i, const auto& index)
      {
          double dx = mesh.cell_length(level);
          auto j = index[0];
  
          unp1(level, i, j) = u(level, i, j) - dt / dx * (u(level, i, j) - u(level, i - 1, j)
                                                        + u(level, i, j) - u(level, i, j - 1));
      });
  
      std::swap(unp1.array(), u.array());
  }

The whole example can be found here.

The projection operator

When manipulating grids of different resolution levels, it is often necessary to transmit the solution of a level $l$ to a level $l+1$ and vice versa. We are interested here in the projection operator defined by

$$ u(l, i, j) = \frac{1}{4}\sum_{k_i=0}^1\sum_{k_j=0}^1 u(l+1, 2i + k_i, 2j + k_j) $$

This operator allows to compute the cell-average value of the solution at a grid node at level $l$ from cell-average values of the solution known on children-nodes at grid level $l + 1$ for a 2D problem.

We assume that we already have a samurai mesh with several level defined in the variable mesh. To access to a level, we use the operator mesh[level]. We also assume that we created a field on this mesh using the make_field and initialized it.

The following steps describe how to implement the projection operator with samurai.

• Create a subset of the mesh using set algebra

auto set = samurai::intersection(mesh[level], mesh[level+1]).on(level);

• Apply an operator on this subset

set([&](const auto& i, const auto index)
{
    auto j = index[0];
    u(level, i, j) = 0.25*(u(level+1, 2*i, 2*j)
                         + u(level+1, 2*i+1, 2*j)
                         + u(level+1, 2*i, 2*j+1)
                         + u(level+1, 2*i+1, 2*j+1));
});

The multi dimensional projection operator can be found here.

There's more

If you want to learn more about samurai skills by looking at examples, we encourage you to browse the demos directory.

The tutorial directory is a good first step followed by the FiniteVolume directory.

Features

• Facilitate data manipulation by using the formalism on a uniform Cartesian grid
• Facilitate the implementation of complex operators between grid levels
• High memory compression of an adapted mesh
• Complex mesh creation using a set of meshes
• Finite volume methods using flux construction
• Lattice Boltzmann methods examples
• Finite difference methods
• Discontinuous Galerkin methods
• Matrix assembling of the discrete operators using PETSc
• AMR cell-based methods
• AMR patch-based and block-based methods
• MRA cell-based methods
• MRA point-based methods
• HDF5 ouput format support
• MPI implementation

Installation

From conda

mamba install samurai

For compiling purposes, you have to install a C++ compiler, cmake, and (optionaly) make:

mamba install cxx-compiler cmake [make]

If you have to use PETSc to assemble the matrix of your problem, you need to install it:

mamba install petsc pkg-config

For parallel computation,

mamba install libboost-mpi libboost-devel libboost-headers 'hdf5=*=mpi*'

From Conan Center

If you want to install samurai from Conan, you can use the following command:

conan install --requires=samurai/0.13.0

From source

Run the cmake configuration

• With mamba or conda

  First, you need to create the environment with all the dependencies installed, run

  mamba env create --file conda/environment.yml

  for sequential computation, or

  mamba env create --file conda/mpi-environment.yml

  for parallel computation. Then

  mamba activate samurai-env

  cmake . -B build -DCMAKE_BUILD_TYPE=Release -DBUILD_DEMOS=ON

• With vcpkg

  cmake . -B ./build -DENABLE_VCPKG=ON -DBUILD_DEMOS=ON

• With conan

  cmake . -B ./build -DCMAKE_BUILD_TYPE=Release -DENABLE_CONAN_OPTION=ON -DBUILD_DEMOS=ON

Build the demos

cmake --build ./build --config Release

CMake configuration

Here is a minimal example of CMakeLists.txt:

cmake_minimum_required(VERSION 3.15)
set(CMAKE_CXX_STANDARD 17)

project(my_samurai_project CXX)

find_package(samurai CONFIG REQUIRED)

add_executable(my_samurai_project main.cpp)
target_link_libraries(my_samurai_project PRIVATE samurai::samurai samurai::libdeps)

Get help

For a better understanding of all the components of samurai, you can consult the documentation https://hpc-math-samurai.readthedocs.io.

If you have any question or remark, you can write a message on github discussions and we will be happy do help you or to discuss with you.

Project assistance

If you want to say thank you or/and support active development of samurai:

• Add a GitHub Star to the project.
• Tweet about samurai.
• Write interesting articles about the project on Dev.to, Medium or your personal blog.

Together, we can make samurai better!

Contributing

First off, thanks for taking the time to contribute! Contributions are what make the open-source community such an amazing place to learn, inspire, and create. Any contributions you make will benefit everybody else and are greatly appreciated.

Please read our contribution guidelines, and thank you for being involved!

License

This project is licensed under the BSD license.

See LICENSE for more information.