NeoN
A framework for CFD software
Loading...
Searching...
No Matches
solver.hpp
Go to the documentation of this file.
1// SPDX-FileCopyrightText: 2023 - 2025 NeoN authors
2//
3// SPDX-License-Identifier: MIT
4
5#pragma once
6
7#include <iostream>
8#include <memory>
9#include <type_traits>
10#include <utility>
11#include <concepts>
12
13#include "NeoN/fields/field.hpp"
15#include "NeoN/core/input.hpp"
19
23
24
25namespace NeoN::dsl
26{
27
28/* @brief solve an expression
29 *
30 * @param exp - Expression which is to be solved/updated.
31 * @param solution - Solution field, where the solution will be 'written to'.
32 * @param t - the time at the start of the time step.
33 * @param dt - time step for the temporal integration
34 * @param fvSchemes - Dictionary containing spatial operator and time integration properties
35 * @param fvSolution - Dictionary containing linear solver properties
36 */
37template<typename VectorType>
38void solve(
40 VectorType& solution,
41 scalar t,
42 scalar dt,
43 const Dictionary& fvSchemes,
44 const Dictionary& fvSolution
45)
46{
47 // TODO:
48 if (exp.temporalOperators().size() == 0 && exp.spatialOperators().size() == 0)
49 {
50 NF_ERROR_EXIT("No temporal or implicit terms to solve.");
51 }
52 exp.read(fvSchemes);
53 if (exp.temporalOperators().size() > 0)
54 {
55 // integrate equations in time
57 fvSchemes.subDict("ddtSchemes"), fvSolution
58 );
59 timeIntegrator.solve(exp, solution, t, dt);
60 }
61 else
62 {
63 // solve sparse matrix system
64 using ValueType = typename VectorType::ElementType;
65
66 auto sparsity = la::SparsityPattern(solution.mesh());
67 auto ls = la::createEmptyLinearSystem<ValueType, localIdx>(solution.mesh(), sparsity);
68
69 exp.implicitOperation(ls);
70 auto expTmp = exp.explicitOperation(solution.mesh().nCells());
71
72 auto [vol, expSource, rhs] = views(solution.mesh().cellVolumes(), expTmp, ls.rhs());
73
74 // subtract the explicit source term from the rhs
76 solution.exec(),
77 {0, rhs.size()},
78 KOKKOS_LAMBDA(const localIdx i) { rhs[i] -= expSource[i] * vol[i]; }
79 );
80
81 auto solver = la::Solver(solution.exec(), fvSolution);
82 solver.solve(ls, solution.internalVector());
83 }
84}
85
86} // namespace dsl
A class representing a dictionary that stores key-value pairs.
Dictionary & subDict(const std::string &key)
Retrieves a sub-dictionary associated with the given key.
const std::vector< TemporalOperator< ValueType > > & temporalOperators() const
void implicitOperation(la::LinearSystem< ValueType, localIdx > &ls)
Vector< ValueType > explicitOperation(localIdx nCells) const
const std::vector< SpatialOperator< ValueType > > & spatialOperators() const
void read(const Dictionary &input)
void solve(Expression &eqn, SolutionVectorType &sol, scalar t, scalar dt)
#define NF_ERROR_EXIT(message)
Macro for printing an error message and aborting the program.
Definition error.hpp:110
void solve(Expression< typename VectorType::ElementType > &exp, VectorType &solution, scalar t, scalar dt, const Dictionary &fvSchemes, const Dictionary &fvSolution)
Definition solver.hpp:38
void parallelFor(const Executor &exec, std::pair< localIdx, localIdx > range, Kernel kernel, std::string name="parallelFor")
int32_t localIdx
Definition label.hpp:32
float scalar
Definition scalar.hpp:16
auto views(Types &... args)
Unpacks all views of the passed classes.
Definition view.hpp:110