This project implements a 1D overlapping Schwarz solver library for linear diffusion–reaction problems using a tridiagonal (Thomas) factorization on each subdomain. Multiple test cases and utilities are provided to try out the solver and visualize results.
The project is built with the vision of expansion in the future for the inclusion of 2D domains and cluster execution using MPI in conjunction with OpenMP.
- Objective of the Projects
- How to read the docs
- Project structure
- Class Structure
- Example of Usage
- Authors
The objective is to create a solver for a 1 dimensional Elliptic PDE with constant diffusion and reaction coefficients. We expect to expand it in the future to non-constant coefficients.
This page only contains overall documentation for the project, for in depth documentation see the appropriate pages:
The project structure is as follows:
01-DD/
├── bin/
├── build/
├── docs/
│ ├── Doxygen/
│ └── pdfs/
├── include/
├── lib/
├── outputs/
├── src/
└── test/
binwill contain the executables (for example tests) created by cmake.buildwill be used by cmake, please see Compilation and execution for more information.docs/Doxygenwill contain the doxygen generated documentation, please see Doxygen documentation for more information.docs/pdfscontains useful pds references pertaining DD methods.includecontains all the.hppfiles used by the library.libswill contain the compiled.alibrary files.outputis a commodity directory used by test to save.vtkfiles.srccontains the actual body of the library and a useful python executable used for visualization, see Visualization for more information.testcontains premade test to stress the library and ensure the correct behaviours of all its components.
This project uses inheritance as a way to avoid repeated code, not much room was left for polymorphism.
All the classes are templates based on the enum Dimension
this is because we intend to extend the library to 2D in the future
Every class in the project inherits from Types<dim> which contains the definition for essential types.
For example in 1D the type Domain is simply a pair of 2 real numbers.
-
The
PDESolver<dim>class is not meant to be used, it's contains useful functions that are used by bothSubdomainSolver<dim>andDiscreteSolver<dim>. It also contains attributes that must be known to both classes to correctly compute the solution. -
The
DiscreteSolver<dim>class construct the approximation of the PDE and handles parallelism by dividing work and different subdomains. It is the class to be called to actually obtain a solution, as it is in this class that we advance computation and check for the correctness of the solution (based on the difference between iterates). -
The
SubdomainSolver<dim>handles the construction of the algebraic system and finds the solution on the assigned subdomain. -
The
FactorizedTridiag<dim>class contains the tridiagonal matrix and solves the system$Au = f$ . Before solving the system the matrix is factorised ($A = L*U$ ) we then solve 2 system in sequence: $$\begin{cases} Ly=b \ Ux = y \end{cases}$$ We callftd.solve(b)to solve the system on the r.h.s.b, please see Thomas Algorithm for more general info on the algorith used to solve the system.
For more information please read the compiled documentation.
To solve your own 1D problem, you need to create the following three small structs and select a grid size:
-
PDEParams— information about the real PDE problemmu(diffusion coefficient),c(reaction coefficient)omega = {a, b}domain limitsdirichlet = {u_a, u_b}boundary values ataandbf(x)right-hand side (lambda or functor)
-
SchwarzParams— information about the domain decomposition parametersNnumber of subdomainsdeltaoverlap width (physical length)
-
SolverParams— iteration controlsmax_itermaximum number of Schwarz iterationsepsstopping tolerance on the sup-norm difference between iterates
-
Grid spacing
h = (omega.b - omega.a) / (N_nodes - 1);N_nodesis your chosen number of grid points.
Typical usage to compute the solution is (as in test/test_solver.cpp):
PDEParams pde{mu, c, f, omega, dirichlet};
SchwarzParams sp{N_subdomains, delta};
SolverParams sol{eps, max_iter};
Real h = (omega.b - omega.a) / (N_nodes - 1);
DiscreteSolver<Line> solver(pde, sp, sol, h);
solver.solve();
// check solver status (optional)
auto u = solver.get_solution();
// Optional: write VTK for visualization
solver.print_to_file(file);Notes:
N_subdomains * subdomain_widthshould match the domain length; choosedeltaso overlaps are a few grid cells wide.fis evaluated pointwise on the grid the solver builds fromomegaandh.- Check
solver.statusfor convergence;solver.iterreports iterations. - Use
OMP_NUM_THREADSenvironment variable to control OpenMP parallelism. - VTK files are written to
outputs/; view withsrc/visualization_pipeline.pyin ParaView.
| Name | PC | GitHub Name |
|---|---|---|
| Samuele Allegranza | 10766615 | samueleallegranza |
| Giulio Enzo Donninelli | 10823453 | gdonninelli |
| Valeriia Potrebina | 11114749 | ValeryPotrebina |
| Alessia Rigoni | 10859832 | alessiarigoni |
| Vale Turco | 10809855 | PurpleVale |