General concepts » First Steps with Rodin

Understanding the basics and creating your first Rodin program.

Now that you have Rodin installed, let's explore the basic concepts and create your first program.

Project Structure

A typical Rodin project consists of:

  • CMakeLists.txt - Build configuration
  • Source files (.cpp) - Your C++ code using Rodin
  • Mesh files - Geometry data (optional, can generate programmatically)

Basic CMakeLists.txt

Here's a minimal CMakeLists.txt for a Rodin project:

cmake_minimum_required(VERSION 3.16)
project(MyRodinProject CXX)

# Set C++20 standard
set(CMAKE_CXX_STANDARD 20)
set(CMAKE_CXX_STANDARD_REQUIRED ON)

# Find Rodin package
find_package(Rodin REQUIRED)

# Create executable
add_executable(my_program main.cpp)

# Link against Rodin libraries
target_link_libraries(my_program PRIVATE
  Rodin::Geometry
  Rodin::Variational
  Rodin::Solver
)

Rodin Modules

Rodin is organized into several modules. Each module has a well-defined responsibility, and you link only the modules you need:

ModuleHeaderPurposeKey classes
Geometry<Rodin/Geometry.h>Meshes, polytopes, connectivityMesh, Polytope, Connectivity
Variational<Rodin/Variational.h>Finite element spaces, forms, problemsP1, P0, H1, TrialFunction, TestFunction, GridFunction
Solver<Rodin/Solver.h>Linear system solversCG, SparseLU, UMFPack
IO<Rodin/IO/XDMF.h>Input/output and visualizationXDMF
Assembly(used internally)Matrix and vector assemblyAssemblers (not used directly)

Hello Rodin

Let's create a simple program that creates a mesh and prints information about it.

#include <Rodin/Geometry.h>
#include <iostream>

using namespace Rodin;
using namespace Rodin::Geometry;

int main() {
  // Create a 2D uniform triangular mesh
  Mesh mesh;
  mesh = mesh.UniformGrid(Polytope::Type::Triangle, {8, 8});

  // Print mesh information
  std::cout << "Created a mesh with:" << std::endl;
  std::cout << "  Vertices: " << mesh.getVertexCount() << std::endl;
  std::cout << "  Cells: " << mesh.getCellCount() << std::endl;
  std::cout << "  Dimension: " << mesh.getDimension() << std::endl;

  return 0;
}

Save this as main.cpp and build it:

cmake .
make
./my_program

You should see output like:

Created a mesh with:
  Vertices: 81
  Cells: 128
  Dimension: 2

Core Concepts

Working with Meshes

A Mesh represents the discretization of your domain $ \Omega $ into a collection of polytopes. You can create meshes in several ways:

Generate a uniform grid on the unit hypercube:**

Mesh mesh;
// 2D triangular mesh: 16×16 subdivisions on [0, 15] × [0, 15]
mesh = mesh.UniformGrid(Polytope::Type::Triangle, {16, 16});

// 3D tetrahedral mesh: 8×8×8 subdivisions
mesh = mesh.UniformGrid(Polytope::Type::Tetrahedron, {8, 8, 8});

Generate a boundary (surface) mesh:**

// 2D box boundary (segments): a rectangle boundary in 2D
mesh = mesh.Box(Polytope::Type::Segment, {4, 4});

// 3D box boundary (triangles): a surface mesh in 3D
mesh = mesh.Box(Polytope::Type::Triangle, {2, 2, 2});

Load from file:**

Mesh mesh;
mesh.load("domain.mesh", IO::FileFormat::MEDIT);   // MEDIT format
mesh.load("domain.mesh", IO::FileFormat::MFEM);     // MFEM format

Key properties:**

mesh.getDimension();          // Topological dimension (2 or 3)
mesh.getSpaceDimension();     // Ambient space dimension
mesh.getVertexCount();        // Total number of vertices
mesh.getCellCount();          // Total number of cells

Computing Connectivity

Before performing certain operations, you may need to compute connectivity information:

Mesh mesh;
mesh = mesh.UniformGrid(Polytope::Type::Triangle, {16, 16});

// Compute face-to-cell connectivity (needed for boundary operations)
mesh.getConnectivity().compute(1, 2);

The compute(d1, d2) method computes connectivity from polytopes of dimension d1 to polytopes of dimension d2. For a 2D mesh:

  • 0 = vertices
  • 1 = edges (faces)
  • 2 = triangles (cells)

Finite Element Spaces

A finite element space $ V_h $ defines the type of basis functions used for approximation. It is built on a mesh and determines how many degrees of freedom (DOFs) the discrete problem will have.

#include <Rodin/Variational.h>

using namespace Rodin::Variational;

Mesh mesh;
mesh = mesh.UniformGrid(Polytope::Type::Triangle, {16, 16});

// P1 (piecewise linear) finite element space
// One DOF per vertex → continuous across cells
P1 Vh(mesh);
std::cout << "P1 DOFs: " << Vh.getSize() << std::endl;

// P0 (piecewise constant) finite element space
// One DOF per cell → discontinuous across cells
P0 Wh(mesh);
std::cout << "P0 DOFs: " << Wh.getSize() << std::endl;

// Vector-valued P1 space (d components)
P1 Uh(mesh, mesh.getSpaceDimension());
std::cout << "Vector P1 DOFs: " << Uh.getSize() << std::endl;

For more details, see Finite Element Spaces.

Trial and Test Functions

To formulate variational problems, we create trial and test functions from a finite element space:

P1 Vh(mesh);

// Trial function: represents the unknown solution u
TrialFunction u(Vh);

// Test function: represents the weighting function v
TestFunction v(Vh);
  • The trial function u is the unknown we solve for. Its solution is a GridFunction accessible via u.getSolution().
  • The test function v is used to construct the equations: for each basis function of $ V_h $ , one equation of the linear system is generated.

These symbolic objects are combined with differential operators (Grad, Div, Jacobian) and integrals (Integral, BoundaryIntegral) to express the variational formulation of a PDE. See Variational Formulations for details.

Using Namespaces

To avoid typing Rodin::Geometry:: and Rodin::Variational:: repeatedly, most Rodin programs start with:

#include <Rodin/Geometry.h>
#include <Rodin/Variational.h>
#include <Rodin/Solver.h>

using namespace Rodin;
using namespace Rodin::Geometry;
using namespace Rodin::Variational;
using namespace Rodin::Solver;

Complete Example

Here's a complete program that demonstrates the basic workflow:

#include <Rodin/Geometry.h>
#include <Rodin/Variational.h>
#include <iostream>

using namespace Rodin;
using namespace Rodin::Geometry;
using namespace Rodin::Variational;

int main() {
  // 1. Create a mesh
  Mesh mesh;
  mesh = mesh.UniformGrid(Polytope::Type::Triangle, {8, 8});

  // 2. Compute connectivity (needed for boundary operations)
  mesh.getConnectivity().compute(1, 2);

  // 3. Create finite element space
  P1 Vh(mesh);

  // 4. Create trial and test functions
  TrialFunction u(Vh);
  TestFunction v(Vh);

  // 5. Print information
  std::cout << "Setup complete!" << std::endl;
  std::cout << "  Mesh cells: " << mesh.getCellCount() << std::endl;
  std::cout << "  DOFs: " << Vh.getSize() << std::endl;

  return 0;
}

What's Next?

You now understand the basic building blocks of a Rodin program. In the next section, Your First Problem: Solving the Poisson Equation, we'll put these concepts together to solve a complete PDE problem: the Poisson equation.