diff --git a/.github/workflows/build_doc.yaml b/.github/workflows/build_doc.yaml index b01189cfc..59f2bfbee 100644 --- a/.github/workflows/build_doc.yaml +++ b/.github/workflows/build_doc.yaml @@ -5,8 +5,6 @@ env: on: pull_request: - branches: - - main types: [closed, synchronize, opened] push: tags: @@ -21,6 +19,7 @@ jobs: && pip3 install furo && pip3 install breathe && pip3 install sphinx-sitemap + && pip3 install sphinxcontrib-mermaid - name: Checkout repo uses: actions/checkout@v4 - name: Build docs @@ -46,6 +45,13 @@ jobs: with: folder: docs_build target-folder: latest + - name: Comment PR + uses: thollander/actions-comment-pull-request@v2 + with: + message: | + Deployed test documentation to https://exasim-project.com/NeoFOAM/Build_PR_${{ env.PR_NUMBER }} + comment_tag: build_url + - name: Echo Build PR URL run: | echo "Deploy to: https://exasim-project.com/NeoFOAM/Build_PR_${{ env.PR_NUMBER }}" diff --git a/doc/DSL/equation.rst b/doc/DSL/equation.rst new file mode 100644 index 000000000..1a059d01b --- /dev/null +++ b/doc/DSL/equation.rst @@ -0,0 +1,50 @@ +Equation +--------- + + +The `Equation` class in NeoFOAM holds, manages, builds and solves the expressed/programmed equation. +Its core responsibility lie in the answering of the following questions: + + - How to discretize the spatial terms? + - In OpenFOAM this information is provided in **fvSchemes** + - How to integrate the system in time? + - In OpenFOAM this information is provided in **fvSchemes** + - How to solve the system? + - In OpenFOAM this information is provided in **fvSolution** + +The main difference between NeoFOAM and OpenFOAM is that the DSL is evaluated lazily, i.e. evaluation is not performed on construction by default. +Since, the evaluation is not tied to the construction but rather delayed, other numerical integration strategies (e.g. RK methods or even sub-stepping within an the equation) are possible. + +To implement lazy evaluation, the `Equation` stores 3 vectors: + +.. mermaid:: + + classDiagram + class Operator { + +explicitOperation(...) + +implicitOperation(...) + } + class DivOperator { + +explicitOperation(...) + +implicitOperation(...) + } + class TemporalOperator { + +explicitOperation(...) + +implicitOperation(...) + } + class Others["..."] { + +explicitOperation(...) + +implicitOperation(...) + } + class Equation { + +temporalTerms_: vector~Operator~ + +implicitTerms_: vector~Operator~ + +explicitTerms_: vector~Operator~ + } + Operator <|-- DivOperator + Operator <|-- TemporalOperator + Operator <|-- Others + Equation <|-- Operator + +Thus, an `Equation` consists of multiple `Operators` which are either explicit, implicit, or temporal. +Consequently, addition, subtraction, and scaling with a field needs to be handled by the `Operator`. diff --git a/doc/DSL/index.rst b/doc/DSL/index.rst new file mode 100644 index 000000000..840f7712f --- /dev/null +++ b/doc/DSL/index.rst @@ -0,0 +1,71 @@ +.. _fvcc: + +Domain Specific Language (DSL) +============================== + +The concept of a Domain Specific Language (DSL) allows to dramatically simplify the process of implementing and solving equations in a given programming language like C++. +Engineers can express equations in a concise and readable form, closely resembling their mathematical representation and no or little knowledge of the numerical schemes and implementation is required. +This approach allows engineers to focus on the physics of the problem rather than the numerical implementation and helps in reducing the time and effort required to develop and maintain complex simulations. + +The Navier-Stokes equations can be expressed in the DSL in OpenFOAM as follows: + +.. code-block:: cpp + + fvVectorMatrix UEqn + ( + fvm::ddt(U) + + fvm::div(phi, U) + - fvm::laplacian(nu, U) + ); + + solve(UEqn == -fvc::grad(p)); + +To solve the continuity equation in OpenFOAM with the PISO or SIMPLE algorithm, the VectorMatrix, UEqn, needs to provide the diagonal and off-diagonal terms of the matrix. + +.. code-block:: cpp + + volScalarField rAU(1.0/UEqn.A()); + volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p)); + +This approach is readable and easy to understand for engineers familiar with OpenFOAM. However, it has several limitations due to its implementation in OpenFOAM: + + - the solution system is always a sparse matrix + - the sparse matrix is always a LDU matrix that is not supported by external linear solvers + - Only cell-centred discretisation is supported + +NeoFOAM DSL tries to address these issues by providing: + + - lazy evaluation of the equations terms. This allows for better optimisation of the resulting equation and can reduce the number of required temporaries. + - a more modular design + - Support for standard matrix formats like COO and CSR, to simplify the use of external linear solvers + +The use of standard matrix format combined with lazy evaluation allows for the use of external libraries to integrate PDEs in time and space. +The equation system can be passed to **sundials** and be integrated by **RK methods** and **BDF methods** on heterogeneous architectures. + +The NeoFOAM DSL is designed as drop in replacement for OpenFOAM DSL and the adoption should be possible with minimal effort and the same equation from above should read: + +.. code-block:: cpp + + dsl::Equation UEqn + ( + fvcc::impOp::ddt(U) + + fvcc::impOp::div(phi, U) + - fvcc::impOp::laplacian(nu, U) + ) + + solve(UEqn == -fvcc::expOp::grad(p)); + + +In contrast to OpenFOAM, here the majority of the work is done in the solve step. +That is 1. assemble the system and 2. solve the system. +After the system is assembled or solved, it provides access to the linear system for the SIMPLE and PISO algorithms. + + + + +.. toctree:: + :maxdepth: 2 + :glob: + + equation.rst + operator.rst diff --git a/doc/DSL/operator.rst b/doc/DSL/operator.rst new file mode 100644 index 000000000..30c0d3c15 --- /dev/null +++ b/doc/DSL/operator.rst @@ -0,0 +1,118 @@ +Operator +======= + + +The `Operator` class represents a term in an equation and can be instantiated with different value types. +An `Operator` is either explicit, implicit or temporal, and can be scalable by a an additional coefficient, for example a scalar value or a further field. +The `Operator` implementation uses Type Erasure (more details `[1] `_ `[2] `_ `[3] `_) to achieve polymorphism without inheritance. Consequently, the class needs only to implement the interface which is used in the DSL and which is shown in the below example: + +Example: + .. code-block:: cpp + + NeoFOAM::DSL::Operator divTerm = + Divergence(NeoFOAM::DSL::Operator::Type::Explicit, exec, ...); + + NeoFOAM::DSL::Operator ddtTerm = + TimeTerm(NeoFOAM::DSL::Operator::Type::Temporal, exec, ..); + + +To fit the specification of the EqnSystem (storage in a vector), the Operator needs to be able to be scaled: + +.. code-block:: cpp + + NeoFOAM::Field scalingField(exec, nCells, 2.0); + auto sF = scalingField.span(); + + dsl::Operator customTerm = + CustomTerm(dsl::Operator::Type::Explicit, exec, nCells, 1.0); + + auto constantScaledTerm = 2.0 * customTerm; // A constant scaling factor of 2 for the term. + auto fieldScaledTerm = scalingField * customTerm; // scalingField is used to scale the term. + + // Operator also supports a similar syntax as OpenFOAM + auto multiScaledTerm = (scale + scale + scale + scale) * customTerm; + + // Operator also supports the use of a lambda as scaling function to reduce the number of temporaries generated + auto lambdaScaledTerm = + (KOKKOS_LAMBDA(const NeoFOAM::size_t i) { return sF[i] + sF[i] + sF[i] + sF[i]; }) * customTerm; + +To add a user-defined `Operator`, a new derived class must be created, inheriting from `OperatorMixin`, + and provide the definitions of the below virtual functions that are required for the `Operator` interface: + + - build: build the term + - explicitOperation: perform the explicit operation + - implicitOperation: perform the implicit operation + - display: display the term + - getType: get the type of the term + - exec: get the executor + - nCells: get the number of cells + - volumeField: get the volume field + +An example is given below: + +.. code-block:: cpp + + class CustomOperator : public dsl::OperatorMixin + { + + public: + + // constructors .. + NeoFOAM::scalar read(const NeoFOAM::Input& input) + { + // .. + } + + void build(const NeoFOAM::Input& input) + { + value = read(input); + termEvaluated = true; + } + + std::string display() const { return "Laplacian"; } + + void explicitOperation(NeoFOAM::Field& source) + { + NeoFOAM::scalar setValue = value; + // scaleField is defined in OperatorMixin + // and accounts for the scaling of the terms + // and considers scaling by fields and scalars + auto scale = scaleField(); + auto sourceField = source.span(); + NeoFOAM::parallelFor( + source.exec(), + {0, source.size()}, + KOKKOS_LAMBDA(const size_t i) { sourceField[i] += scale[i] * setValue; } + ); + } + + // other helper functions + dsl::Operator::Type getType() const { return termType_; } + + const NeoFOAM::Executor& exec() const { return exec_; } + + std::size_t nCells() const { return nCells_; } + + fvcc::VolumeField* volumeField() { return nullptr; } + + dsl::Operator::Type termType_; + + + const NeoFOAM::Executor exec_; + std::size_t nCells_; + NeoFOAM::scalar value = 1.0; + }; + +The required scaling of the term is handle by the `scaleField` function, provided by `OperatorMixin`. The `scaleField` function returns the 'ScalingField' class that is used to scale by fields and scalars. + +.. code-block:: cpp + + template + class ScalingField + { + + // the span is only used if it is defined + KOKKOS_INLINE_FUNCTION + ValueType operator[](const size_t i) const { return useSpan ? values[i] * value : value; } + + } diff --git a/doc/conf.py b/doc/conf.py index 731d51a8b..56e66dd30 100644 --- a/doc/conf.py +++ b/doc/conf.py @@ -23,6 +23,7 @@ extensions = [ 'sphinx.ext.autodoc', + 'sphinxcontrib.mermaid', 'sphinx.ext.intersphinx', 'sphinx.ext.autosectionlabel', 'sphinx.ext.todo', diff --git a/doc/finiteVolume/DSL.rst b/doc/finiteVolume/DSL.rst deleted file mode 100644 index 514441de9..000000000 --- a/doc/finiteVolume/DSL.rst +++ /dev/null @@ -1,4 +0,0 @@ -.. _fvcc_DSL: - -Domain Specific Language -======================== diff --git a/doc/finiteVolume/cellCentred/index.rst b/doc/finiteVolume/cellCentred/index.rst index 21ce5a6d2..16b862c1e 100644 --- a/doc/finiteVolume/cellCentred/index.rst +++ b/doc/finiteVolume/cellCentred/index.rst @@ -9,7 +9,6 @@ cellCenteredFiniteVolume :glob: DSL.rst - fields.rst boundaryConditions.rst operators.rst stencil.rst diff --git a/doc/index.rst b/doc/index.rst index 0d7e3fef6..af18c2ae6 100644 --- a/doc/index.rst +++ b/doc/index.rst @@ -26,6 +26,7 @@ Table of Contents contributing basics/index finiteVolume/cellCentred/index + DSL/index api/index Compatibility with OpenFOAM