17template<
typename ValueType>
27 const auto exec = phi.
exec();
31 const auto [owner, neighbour, surfFaceCells] =
34 auto [refGradient, value, valueFraction, refValue] =
spans(
41 const auto [result, faceArea, fnGrad, vol] =
48 if (std::holds_alternative<SerialExecutor>(exec))
50 for (
size_t i = 0; i < nInternalFaces; i++)
52 ValueType flux = faceArea[i] * fnGrad[i];
53 result[
static_cast<size_t>(owner[i])] += flux;
54 result[
static_cast<size_t>(neighbour[i])] -= flux;
57 for (
size_t i = nInternalFaces; i < fnGrad.size(); i++)
59 auto own =
static_cast<size_t>(surfFaceCells[i - nInternalFaces]);
60 ValueType valueOwn = faceArea[i] * fnGrad[i];
61 result[own] += valueOwn;
64 for (
size_t celli = 0; celli < mesh.
nCells(); celli++)
66 result[celli] *= operatorScaling[celli] / vol[celli];
74 KOKKOS_LAMBDA(
const size_t i) {
75 ValueType flux = faceArea[i] * fnGrad[i];
76 Kokkos::atomic_add(&result[
static_cast<size_t>(owner[i])], flux);
77 Kokkos::atomic_sub(&result[
static_cast<size_t>(neighbour[i])], flux);
83 {nInternalFaces, fnGrad.size()},
84 KOKKOS_LAMBDA(
const size_t i) {
85 auto own =
static_cast<size_t>(surfFaceCells[i - nInternalFaces]);
86 ValueType valueOwn = faceArea[i] * fnGrad[i];
87 Kokkos::atomic_add(&result[own], valueOwn);
94 KOKKOS_LAMBDA(
const size_t celli) {
95 result[celli] *= operatorScaling[celli] / vol[celli];
101template<
typename ValueType>
110 static std::string
name() {
return "Gauss"; }
112 static std::string
doc() {
return "Gauss-Green Divergence"; }
114 static std::string
schema() {
return "none"; }
117 : Base(exec, mesh), surfaceInterpolation_(exec, mesh, inputs),
118 faceNormalGradient_(exec, mesh, inputs),
124 auto [A, b, sp] = ls.
view();
125 const auto& exec = A.exec();
129 Field<localIdx> mRowPtrs(exec, A.rowPtrs().data(), A.rowPtrs().size());
144 computeLaplacian<ValueType>(
145 faceNormalGradient_, gamma, phi, lapPhi.
internalField(), operatorScaling
156 computeLaplacian<ValueType>(faceNormalGradient_, gamma, phi, lapPhi, operatorScaling);
168 const auto exec = phi.
exec();
169 const auto [owner, neighbour, surfFaceCells, diagOffs, ownOffs, neiOffs] =
spans(
173 sparsityPattern_->diagOffset(),
174 sparsityPattern_->ownerOffset(),
175 sparsityPattern_->neighbourOffset()
178 const auto [sGamma, deltaCoeffs, magFaceArea] =
spans(
180 faceNormalGradient_.deltaCoeffs().internalField(),
184 auto [refGradient, value, valueFraction, refValue] =
spans(
191 const auto rowPtrs = ls.
matrix().rowPtrs();
192 const auto colIdxs = ls.
matrix().colIdxs();
193 auto values = ls.
matrix().values();
194 auto rhs = ls.
rhs().span();
199 KOKKOS_LAMBDA(
const size_t facei) {
200 scalar flux = deltaCoeffs[facei] * sGamma[facei] * magFaceArea[facei];
202 std::size_t own =
static_cast<std::size_t
>(owner[facei]);
203 std::size_t nei =
static_cast<std::size_t
>(neighbour[facei]);
206 std::size_t rowNeiStart = rowPtrs[nei];
207 std::size_t rowOwnStart = rowPtrs[own];
210 values[rowNeiStart + neiOffs[facei]] += flux * one<ValueType>();
211 Kokkos::atomic_sub(&values[rowOwnStart + diagOffs[own]], flux * one<ValueType>());
216 values[rowOwnStart + ownOffs[facei]] += flux * one<ValueType>();
217 Kokkos::atomic_sub(&values[rowNeiStart + diagOffs[nei]], flux * one<ValueType>());
223 {nInternalFaces, sGamma.size()},
224 KOKKOS_LAMBDA(
const size_t facei) {
225 std::size_t bcfacei = facei - nInternalFaces;
226 scalar flux = sGamma[facei] * magFaceArea[facei];
228 std::size_t own =
static_cast<std::size_t
>(surfFaceCells[bcfacei]);
229 std::size_t rowOwnStart = rowPtrs[own];
231 values[rowOwnStart + diagOffs[own]] -=
232 flux * valueFraction[bcfacei] * deltaCoeffs[facei] * one<ValueType>();
235 * (valueFraction[bcfacei] * deltaCoeffs[facei] * refValue[bcfacei]
236 + (1.0 - valueFraction[bcfacei]) * refGradient[bcfacei]));
243 KOKKOS_LAMBDA(
const size_t celli) {
244 rhs[celli] *= operatorScaling[celli];
245 for (
size_t i = rowPtrs[celli]; i < rowPtrs[celli + 1]; i++)
247 values[i] *= operatorScaling[celli];
253 std::unique_ptr<LaplacianOperatorFactory<ValueType>>
clone()
const override
255 return std::make_unique<GaussGreenLaplacian<ValueType>>(*this);
262 const std::shared_ptr<SparsityPattern> sparsityPattern_;
const NeoFOAM::Field< T > & refGrad() const
Get the view storing the Neumann boundary values.
const NeoFOAM::Field< T > & value() const
Get the view storing the computed values from the boundary condition.
const NeoFOAM::Field< T > & refValue() const
Get the view storing the Dirichlet boundary values.
const NeoFOAM::Field< scalar > & valueFraction() const
Get the view storing the fraction of the boundary value.
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.
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 scalarField & magFaceAreas() const
Get the field of magnitudes of face areas.
size_t nCells() const
Get the number of cells 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.
void faceNormalGrad(const VolumeField< ValueType > &volField, SurfaceField< ValueType > &surfaceField) const
static std::string schema()
la::LinearSystem< ValueType, localIdx > createEmptyLinearSystem() const override
virtual void laplacian(VolumeField< ValueType > &lapPhi, const SurfaceField< scalar > &gamma, VolumeField< ValueType > &phi, const dsl::Coeff operatorScaling) override
std::unique_ptr< LaplacianOperatorFactory< ValueType > > clone() const override
virtual void laplacian(la::LinearSystem< ValueType, localIdx > &ls, const SurfaceField< scalar > &gamma, VolumeField< ValueType > &phi, const dsl::Coeff operatorScaling) override
GaussGreenLaplacian(const Executor &exec, const UnstructuredMesh &mesh, const Input &inputs)
virtual void laplacian(Field< ValueType > &lapPhi, const SurfaceField< scalar > &gamma, VolumeField< ValueType > &phi, const dsl::Coeff operatorScaling) override
static std::string name()
const UnstructuredMesh & mesh() const
Returns a const reference to the unstructured mesh object.
const Executor & exec() const
Returns a const reference to the executor object.
const Field< ValueType > & internalField() const
Returns a const reference to the internal field.
const BoundaryFields< ValueType > & boundaryField() const
Returns a const reference to the boundary field.
Represents a surface field in a finite volume method.
Represents a volume field in a finite volume method.
A class representing a linear system of equations.
LinearSystemView< ValueType, IndexType > view()
Field< ValueType > & rhs()
std::string sparsityPattern() const
CSRMatrix< ValueType, IndexType > & matrix()
A template class for registering derived classes with a base class.
void computeLaplacian(const FaceNormalGradient< ValueType > &faceNormalGradient, const SurfaceField< scalar > &gamma, VolumeField< ValueType > &phi, Field< ValueType > &lapPhi, const dsl::Coeff operatorScaling)
std::variant< Dictionary, TokenList > Input
void parallelFor(const Executor &exec, std::pair< size_t, size_t > range, Kernel kernel)
std::variant< SerialExecutor, CPUExecutor, GPUExecutor > Executor
auto spans(Args &... fields)