Interpolation between different meshes (#2245)

* add neighbor version

* minor fix

* small simplifications

* fix type of points

* add timers

* update demo

* rewrite cell collision

* some small updates

* fix function call

* fix tests

* minor adjustments

* Add test of 2D->3D and 3D->2D interpolation + one std::size_t compiler fix

* create function

* fix request

* fix funciton call

* remove unused variable

* Update cpp/dolfinx/fem/interpolate.h

Co-authored-by: Jørgen Schartum Dokken <dokken92@gmail.com>

* Apply suggestions from code review

Co-authored-by: Jørgen Schartum Dokken <dokken92@gmail.com>

* Loop optimisation

* minor tweaks in demo

* Add new function that distributes points

* Separate functions into sub-routines

* std::size_t->std::int32_t

* Remove extra input parameter

* Improve docs

* Regenerate demos and various span fixes

* Remove comment

* Remove const span reference and fix in,out

* Debug

* Works in serial, crashes in parallel

* Refactoring code. Remove xtensor. Various documation fixes.

* Doc update

* Applying suggestions by Garth

* Remove copy_N

* Tentatively updating python wrapper

Co-authored-by: IgorBaratta <igorbaratta@gmail.com>
Co-authored-by: Jørgen Dokken <jsd55@cam.ac.uk>
Co-authored-by: Jørgen Schartum Dokken <dokken92@gmail.com>
Co-authored-by: Massimiliano Leoni <mleoni@bcamath.org>
Co-authored-by: Chris Richardson <chris@bpi.cam.ac.uk>
Co-authored-by: Garth N. Wells <gnw20@cam.ac.uk>
Co-authored-by: Matthew Scroggs <matthew.w.scroggs@gmail.com>

cpp/demo/interpolation-io/CMakeLists.txt

@@ -16,6 +16,8 @@ if (NOT TARGET dolfinx)
   find_package(DOLFINX REQUIRED)
 endif()
 
+set(CMAKE_INCLUDE_CURRENT_DIR ON)
+
 add_executable(${PROJECT_NAME} main.cpp)
 target_link_libraries(${PROJECT_NAME} dolfinx)
 

cpp/demo/interpolation_different_meshes/CMakeLists.txt

@@ -0,0 +1,35 @@
+# This file was generated by running
+#
+#     python cmake/scripts/generate-cmakefiles from dolfinx/cpp
+#
+cmake_minimum_required(VERSION 3.16)
+
+set(PROJECT_NAME demo_interpolation_different_meshes)
+project(${PROJECT_NAME} LANGUAGES C CXX)
+
+# Set C++20 standard
+set(CMAKE_CXX_STANDARD 20)
+set(CMAKE_CXX_STANDARD_REQUIRED ON)
+set(CMAKE_CXX_EXTENSIONS OFF)
+
+if (NOT TARGET dolfinx)
+  find_package(DOLFINX REQUIRED)
+endif()
+
+set(CMAKE_INCLUDE_CURRENT_DIR ON)
+
+add_executable(${PROJECT_NAME} main.cpp)
+target_link_libraries(${PROJECT_NAME} dolfinx)
+
+# Do not throw error for 'multi-line comments' (these are typical in
+# rst which includes LaTeX)
+include(CheckCXXCompilerFlag)
+CHECK_CXX_COMPILER_FLAG("-Wno-comment" HAVE_NO_MULTLINE)
+set_source_files_properties(main.cpp PROPERTIES COMPILE_FLAGS "$<$<BOOL:${HAVE_NO_MULTLINE}>:-Wno-comment -Wall -Wextra -pedantic -Werror>")
+
+# Test targets (used by DOLFINx testing system)
+set(TEST_PARAMETERS2 -np 2 ${MPIEXEC_PARAMS} "./${PROJECT_NAME}")
+set(TEST_PARAMETERS3 -np 3 ${MPIEXEC_PARAMS} "./${PROJECT_NAME}")
+add_test(NAME ${PROJECT_NAME}_mpi_2 COMMAND "mpirun" ${TEST_PARAMETERS2})
+add_test(NAME ${PROJECT_NAME}_mpi_3 COMMAND "mpirun" ${TEST_PARAMETERS3})
+add_test(NAME ${PROJECT_NAME}_serial COMMAND ${PROJECT_NAME})

cpp/demo/interpolation_different_meshes/main.cpp

@@ -0,0 +1,72 @@
+// Copyright (C) 2022 Igor A. Baratta and Massimiliano Leoni
+//
+// This file is part of DOLFINx (https://www.fenicsproject.org)
+//
+// SPDX-License-Identifier:    LGPL-3.0-or-later
+
+#include <basix/e-lagrange.h>
+#include <dolfinx/common/timing.h>
+#include <dolfinx/fem/utils.h>
+#include <dolfinx/io/ADIOS2Writers.h>
+#include <dolfinx/mesh/generation.h>
+#include <memory>
+
+using namespace dolfinx;
+
+using T = double;
+
+int main(int argc, char* argv[])
+{
+  dolfinx::init_logging(argc, argv);
+  MPI_Init(&argc, &argv);
+  {
+    MPI_Comm comm{MPI_COMM_WORLD};
+
+    // Create a tetrahedral mesh
+    auto mesh_tet = std::make_shared<mesh::Mesh>(
+        mesh::create_box(comm, {{{0, 0, 0}, {1, 1, 1}}}, {10, 10, 10},
+                         mesh::CellType::tetrahedron, mesh::GhostMode::none));
+    // Create a hexahedral mesh
+    auto mesh_hex = std::make_shared<mesh::Mesh>(
+        mesh::create_box(comm, {{{0, 0, 0}, {1, 1, 1}}}, {9, 8, 7},
+                         mesh::CellType::hexahedron, mesh::GhostMode::none));
+
+    basix::FiniteElement eL = basix::element::create_lagrange(
+        mesh::cell_type_to_basix_type(mesh_tet->topology().cell_type()), 1,
+        basix::element::lagrange_variant::equispaced, false);
+    auto V_tet = std::make_shared<fem::FunctionSpace>(
+        fem::create_functionspace(mesh_tet, eL, 3));
+
+    basix::FiniteElement eR = basix::element::create_lagrange(
+        mesh::cell_type_to_basix_type(mesh_hex->topology().cell_type()), 2,
+        basix::element::lagrange_variant::equispaced, false);
+    auto V_hex = std::make_shared<fem::FunctionSpace>(
+        fem::create_functionspace(mesh_hex, eR, 3));
+
+    auto u_tet = std::make_shared<fem::Function<T>>(V_tet);
+    auto u_hex = std::make_shared<fem::Function<T>>(V_hex);
+
+    auto fun = [](auto& x)
+    {
+      auto r = xt::zeros_like(x);
+      xt::row(r, 0) = xt::cos(10 * xt::row(x, 0)) * xt::sin(10 * xt::row(x, 2));
+      xt::row(r, 1) = xt::sin(10 * xt::row(x, 0)) * xt::sin(10 * xt::row(x, 2));
+      xt::row(r, 2) = xt::cos(10 * xt::row(x, 0)) * xt::cos(10 * xt::row(x, 2));
+      return r;
+    };
+
+    u_tet->interpolate(fun);
+    u_hex->interpolate(*u_tet);
+
+#ifdef HAS_ADIOS2
+    io::VTXWriter write_tet(mesh_tet->comm(), "u_tet.vtx", {u_tet});
+    write_tet.write(0.0);
+
+    io::VTXWriter write_hex(mesh_hex->comm(), "u_hex.vtx", {u_hex});
+    write_hex.write(0.0);
+#endif
+  }
+  MPI_Finalize();
+
+  return 0;
+}

cpp/demo/poisson/CMakeLists.txt

@@ -7,7 +7,7 @@ cmake_minimum_required(VERSION 3.16)
 set(PROJECT_NAME demo_poisson)
 project(${PROJECT_NAME} LANGUAGES C CXX)
 
-# Set standard
+# Set C++20 standard
 set(CMAKE_CXX_STANDARD 20)
 set(CMAKE_CXX_STANDARD_REQUIRED ON)
 set(CMAKE_CXX_EXTENSIONS OFF)

cpp/demo/poisson_matrix_free/CMakeLists.txt

@@ -25,6 +25,8 @@ add_custom_command(
   COMMAND ffcx ${CMAKE_CURRENT_SOURCE_DIR}/poisson.py ${SCALAR_TYPE}
   VERBATIM DEPENDS poisson.py COMMENT "Compile poisson.py using FFCx")
 
+set(CMAKE_INCLUDE_CURRENT_DIR ON)
+
 add_executable(${PROJECT_NAME} main.cpp ${CMAKE_CURRENT_BINARY_DIR}/poisson.c)
 target_link_libraries(${PROJECT_NAME} dolfinx)