NeoFOAM
WIP Prototype of a modern OpenFOAM core
Loading...
Searching...
No Matches
surfaceIntegrate.hpp
Go to the documentation of this file.
1// SPDX-License-Identifier: MIT
2// SPDX-FileCopyrightText: 2025 NeoFOAM authors
3
4#pragma once
5
13
15{
16
17template<typename ValueType>
19 const Executor& exec,
20 size_t nInternalFaces,
21 std::span<const int> neighbour,
22 std::span<const int> owner,
23 std::span<const int> faceCells,
24 std::span<const ValueType> flux,
25 std::span<const scalar> V,
26 std::span<ValueType> res,
27 const dsl::Coeff operatorScaling
28)
29{
30 size_t nCells {V.size()};
31 const size_t nBoundaryFaces = faceCells.size();
32 // check if the executor is GPU
33 if (std::holds_alternative<SerialExecutor>(exec))
34 {
35 for (size_t i = 0; i < nInternalFaces; i++)
36 {
37 res[static_cast<size_t>(owner[i])] += flux[i];
38 res[static_cast<size_t>(neighbour[i])] -= flux[i];
39 }
40
41 for (size_t i = nInternalFaces; i < nInternalFaces + nBoundaryFaces; i++)
42 {
43 auto own = static_cast<size_t>(faceCells[i - nInternalFaces]);
44 res[own] += flux[i];
45 }
46
47 // TODO does it make sense to store invVol and multiply?
48 for (size_t celli = 0; celli < nCells; celli++)
49 {
50 res[celli] *= operatorScaling[celli] / V[celli];
51 }
52 }
53 else
54 {
56 exec,
57 {0, nInternalFaces},
58 KOKKOS_LAMBDA(const size_t i) {
59 Kokkos::atomic_add(&res[static_cast<size_t>(owner[i])], flux[i]);
60 Kokkos::atomic_sub(&res[static_cast<size_t>(neighbour[i])], flux[i]);
61 }
62 );
63
65 exec,
66 {nInternalFaces, nInternalFaces + nBoundaryFaces},
67 KOKKOS_LAMBDA(const size_t i) {
68 auto own = static_cast<size_t>(faceCells[i - nInternalFaces]);
69 Kokkos::atomic_add(&res[own], flux[i]);
70 }
71 );
72
74 exec,
75 {0, nCells},
76 KOKKOS_LAMBDA(const size_t celli) { res[celli] *= operatorScaling[celli] / V[celli]; }
77 );
78 }
79}
80
81
82template<typename ValueType>
84{
85
86public:
87
88 using FieldValueType = ValueType;
89
91 : flux_(flux), type_(dsl::Operator::Type::Explicit), coeffs_(1.0) {};
92
94 : flux_(surfaceIntegrate.flux_), type_(surfaceIntegrate.type_),
95 coeffs_(surfaceIntegrate.coeffs_) {};
96
97
98 void build(const Input& input)
99 {
100 // do nothing
101 }
102
104 {
105 NeoFOAM::Field<ValueType> tmpsource(source.exec(), source.size(), zero<ValueType>());
106 const auto operatorScaling = this->getCoefficient();
107
108 const UnstructuredMesh& mesh = flux_.mesh();
109 const auto exec = flux_.exec();
110
111 size_t nInternalFaces = mesh.nInternalFaces();
112 surfaceIntegrate<ValueType>(
113 exec,
114 nInternalFaces,
115 mesh.faceNeighbour().span(),
116 mesh.faceOwner().span(),
117 mesh.boundaryMesh().faceCells().span(),
118 this->flux_.internalField().span(),
119 mesh.cellVolumes().span(),
120 tmpsource.span(),
121 operatorScaling
122
123 );
124 source += tmpsource;
125 }
126
127 const Executor& exec() const { return flux_.exec(); }
128
129 dsl::Coeff& getCoefficient() { return coeffs_; }
130
131 const dsl::Coeff& getCoefficient() const { return coeffs_; }
132
133 dsl::Operator::Type getType() const { return type_; }
134
135 std::string getName() const { return "SurfaceIntegrate"; }
136
137private:
138
139 const SurfaceField<ValueType>& flux_;
141 dsl::Coeff coeffs_;
142};
143
144
145} // namespace NeoFOAM
const labelField & faceCells() const
Get the field of face cells.
A class to contain the data and executors for a field and define some basic operations.
Definition field.hpp:49
size_t size() const
Gets the size of the field.
Definition field.hpp:357
const Executor & exec() const
Gets the executor associated with the field.
Definition field.hpp:351
std::span< ValueType > span() &&=delete
Represents an unstructured mesh in NeoFOAM.
const BoundaryMesh & boundaryMesh() const
Get the boundary mesh.
const scalarField & cellVolumes() const
Get the field of cell volumes in the mesh.
size_t nInternalFaces() const
Get the number of internal faces in the mesh.
const labelField & faceOwner() const
Get the field of face owner cells.
const labelField & faceNeighbour() const
Get the field of face neighbour cells.
A class that represents a coefficient for the NeoFOAM dsl.
Definition coeff.hpp:22
Represents a surface field in a finite volume method.
SurfaceIntegrate(const SurfaceIntegrate &surfaceIntegrate)
SurfaceIntegrate(const SurfaceField< ValueType > &flux)
void surfaceIntegrate(const Executor &exec, size_t nInternalFaces, std::span< const int > neighbour, std::span< const int > owner, std::span< const int > faceCells, std::span< const ValueType > flux, std::span< const scalar > V, std::span< ValueType > res, const dsl::Coeff operatorScaling)
std::variant< Dictionary, TokenList > Input
Definition input.hpp:13
void parallelFor(const Executor &exec, std::pair< size_t, size_t > range, Kernel kernel)
std::variant< SerialExecutor, CPUExecutor, GPUExecutor > Executor
Definition executor.hpp:16