Field Architecture

[greole]

I am currently thinking about what the right level of abstraction and inheritance for expressing fields could be. For OpenFOAM we have several field classes:

1. Field basically a container of data allowing simple arithmetic functions like +,- and scalar *,/ build on top of List and UList
2. DimensionedField add runtime dimension checks
3. [primitiveField]() e.g. Field<scalar>
4. GeometricField<T, Patch, GeoMesh>, adds boundaries and reference to geometric data.

Additionally various typedefs exists for convenience like:

// volFieldsFwd.H
typedef GeometricField<scalar, fvPatchField, volMesh> volScalarField;
typedef GeometricField<vector, fvPatchField, volMesh> volVectorField;

The following code snipped shows how boundary and internalField are handled

// GeometricFieldFunctions.C
template                                                                       \
<class Form, class Type, template<class> class PatchField, class GeoMesh>      \
void opFunc                                                                    \
(                                                                              \
    GeometricField                                                             \
    <typename product<Type, Form>::type, PatchField, GeoMesh>& gf,             \
    const GeometricField<Type, PatchField, GeoMesh>& gf1,                      \
    const dimensioned<Form>& dvs                                               \
)                                                                              \
{                                                                              \
    Foam::opFunc(gf.primitiveFieldRef(), gf1.primitiveField(), dvs.value());   \
    Foam::opFunc(gf.boundaryFieldRef(), gf1.boundaryField(), dvs.value());     \
    gf.oriented() = gf1.oriented();                                            \
}                  

Required Components

1. Memory handling and an underlying container class (gko::array which has an executor)

Proposed changes

1. Use composition over inheritance.

Open questions:

1. Naming: Use ScalarField instead of volScalarField to make naming more concise; (GeometricScalarField)
2. What to expose to public API, for example low level container and memory related things might not required to be exposed; Most likely users should only need access to fields like scalarField and volScalarField.
3. Do we need dimension checks at runTime;
   • probably easier to implement compared to compile time unit checking especially since we want to read field data and variables at run time.
   • Use mixin for that
   • (for implementation of POC not needed yet, )
4. How to handle scalarFields vs volScalarFields, both have access to arithmetic functions but volScalarFields contain boundaries;
5. Should we implement a factory approach for field construction as in Ginkgo -> look at
6. How to incorporate distributed functionality
   • does it know it is one of many

Notes

• look at domain field

[HenningScheufler]

One of the requirements should be that you can switch between different parallel implementation: mpi, openmp+mpi and GPU+mpi by changing a keyword. Consequently, the field needs to be allocated on the different corresponding device.

A solution to achive could be the executor model of ginkgo.

Ginkgo implements multiple executors:

// define the omp, cuda, hip and reference kernels which will be bound to the
// operation
namespace kernels {
namespace omp {
void my_kernel(int x) {
     // omp code
}
}
namespace cuda {
void my_kernel(int x) {
     // cuda code
}
}
namespace hip {
void my_kernel(int x) {
     // hip code
}
}
namespace dpcpp {
void my_kernel(int x) {
     // dpcpp code
}
}
namespace reference {
void my_kernel(int x) {
    // reference code
}
}

// Bind the kernels to the operation
[GKO_REGISTER_OPERATION](https://ginkgo-project.github.io/ginkgo-generated-documentation/doc/develop/group__Executor.html#ga7f2c119ff9f4f51c7f8e2d5dbbbbd044)(my_op, my_kernel);

int main() {
    // create executors
    auto omp = OmpExecutor::create();
    auto cuda = CudaExecutor::create(0, omp);
    auto hip = HipExecutor::create(0, omp);
    auto dpcpp = DpcppExecutor::create(0, omp);
    auto ref = ReferenceExecutor::create();

However, this would allow to implement the kernel differently for different GPU architectures. But it would also requires the implementation and maintenance of the same kernel for multiple architectures.

I could be sufficient to have a GPU, OpenMP, Serial exectutor. The type of GPU would be decided at compile time e.g. CUDA, HIP, SYCL.

Another question is how would you distribute the binaries (parameters: architecture, implementation)?

[HenningScheufler]

@bevanwsjones

[HenningScheufler]

regarding 3 (units): This library could be a candidate https://github.com/mpusz/mp-units

[greole]

Also potentially interesting talk https://www.youtube.com/watch?v=aF3samjRzD4

[greole]

Use thrust function calls to get std::algorithm kind functionality. For executor model ->

[greole]

Discussing it we came up with the names Domain for DomainField and MeshedDomain for GeometricFields