sourceTerm.hpp

Go to the documentation of this file.

// SPDX-License-Identifier: MIT
// SPDX-FileCopyrightText: 2023 NeoFOAM authors
 
#pragma once
 
#include "NeoFOAM/fields/field.hpp"
#include "NeoFOAM/core/executor/executor.hpp"
#include "NeoFOAM/core/input.hpp"
#include "NeoFOAM/dsl/spatialOperator.hpp"
#include "NeoFOAM/mesh/unstructured.hpp"
#include "NeoFOAM/finiteVolume/cellCentred/fields/volumeField.hpp"
#include "NeoFOAM/finiteVolume/cellCentred/operators/sparsityPattern.hpp"
 
namespace NeoFOAM::finiteVolume::cellCentred
{
 
 
template<typename ValueType>
class SourceTerm : public dsl::OperatorMixin<VolumeField<ValueType>>
{
 
public:
 
    using FieldValueType = ValueType;
 
    SourceTerm(
        dsl::Operator::Type termType,
        VolumeField<scalar>& coefficients,
        VolumeField<ValueType>& field
    )
        : dsl::OperatorMixin<VolumeField<ValueType>>(
            field.exec(), dsl::Coeff(1.0), field, termType
        ),
          coefficients_(coefficients),
          sparsityPattern_(SparsityPattern::readOrCreate(field.mesh())) {};
    SourceTerm( {…}
 
    void explicitOperation(Field<ValueType>& source)
    {
        auto operatorScaling = this->getCoefficient();
        const auto vol = coefficients_.mesh().cellVolumes().span();
        auto [sourceSpan, fieldSpan, coeff] =
            spans(source, this->field_.internalField(), coefficients_.internalField());
        NeoFOAM::parallelFor(
            source.exec(),
            source.range(),
            KOKKOS_LAMBDA(const size_t celli) {
                sourceSpan[celli] += operatorScaling[celli] * coeff[celli] * fieldSpan[celli];
            }
        );
    }
    void explicitOperation(Field<ValueType>& source) {…}
 
    la::LinearSystem<ValueType, localIdx> createEmptyLinearSystem() const
    {
        la::LinearSystem<scalar, localIdx> ls(sparsityPattern_->linearSystem());
        auto [A, b, sp] = ls.view();
        const auto& exec = A.exec();
 
        Field<ValueType> values(exec, A.nNonZeros(), zero<ValueType>());
        Field<localIdx> mColIdxs(exec, A.colIdxs().data(), A.nNonZeros());
        Field<localIdx> mRowPtrs(exec, A.rowPtrs().data(), A.rowPtrs().size());
 
        la::CSRMatrix<ValueType, localIdx> matrix(values, mColIdxs, mRowPtrs);
        Field<ValueType> rhs(exec, b.size(), zero<ValueType>());
 
        return {matrix, rhs, ls.sparsityPattern()};
    }
    la::LinearSystem<ValueType, localIdx> createEmptyLinearSystem() const {…}
 
    void implicitOperation(la::LinearSystem<ValueType, localIdx>& ls)
    {
        const auto operatorScaling = this->getCoefficient();
        const auto vol = coefficients_.mesh().cellVolumes().span();
        const auto [diagOffs, coeff] =
            spans(sparsityPattern_->diagOffset(), coefficients_.internalField());
        auto [values, cols, rows] = ls.matrix().span().span();
 
        NeoFOAM::parallelFor(
            ls.exec(),
            {0, coeff.size()},
            KOKKOS_LAMBDA(const size_t celli) {
                std::size_t idx = rows[celli] + diagOffs[celli];
                values[idx] +=
                    operatorScaling[celli] * coeff[celli] * vol[celli] * one<ValueType>();
            }
        );
    }
    void implicitOperation(la::LinearSystem<ValueType, localIdx>& ls) {…}
 
 
    void build(const Input& input)
    {
        // do nothing
    }
    void build(const Input& input) {…}
 
    std::string getName() const { return "DivOperator"; }
 
private:
 
    const VolumeField<scalar>& coefficients_;
    const std::shared_ptr<SparsityPattern> sparsityPattern_;
};
class SourceTerm : public dsl::OperatorMixin<VolumeField<ValueType>> {…};
 
 
} // namespace NeoFOAM