solverbase.h 6.83 KB
Newer Older
Jonas Kusch's avatar
Jonas Kusch committed
1
2
3
4
#ifndef SOLVER_H
#define SOLVER_H

// include Matrix, Vector definitions
5
6
#include "common/globalconstants.h"
#include "common/typedef.h"
7

8
9
10
// externals
#include "spdlog/spdlog.h"

steffen.schotthoefer's avatar
steffen.schotthoefer committed
11
// Forward Declarations
12
13
14
15
16
class NumericalFlux;
class Mesh;
class Config;
class ProblemBase;
class QuadratureBase;
17
class Reconstructor;
steffen.schotthoefer's avatar
steffen.schotthoefer committed
18

Jonas Kusch's avatar
Jonas Kusch committed
19
20
class Solver
{
21
  protected:
22
23
24
25
    Mesh* _mesh;           /*! @brief mesh object for writing out information */
    NumericalFlux* _g;     /*! @brief class for numerical flux */
    Config* _settings;     /*! @brief config class for global information */
    ProblemBase* _problem; /*! @brief problem class for initial conditions */
steffen.schotthoefer's avatar
steffen.schotthoefer committed
26

27
28
    // --------- Often used variables of member classes for faster access ----

29
30
31
32
33
34
    unsigned _nEnergies;          /*! @brief number of energy/time steps, number of nodal energy values for CSD */
    double _dE;                   /*! @brief energy/time step size */
    Vector _energies;             /*! @brief energy groups used in the simulation [keV] */
    std::vector<double> _density; /*! @brief patient density, dim(_density) = _nCells */
    Vector _s;                    /*! @brief stopping power, dim(_s) = _nTimeSteps    */
    std::vector<VectorVector> _Q; /*! @brief external source term */
35
36
37
38
39

    VectorVector _sigmaS; /*!  @brief scattering cross section for all energies */
    VectorVector _sigmaT; /*!  @brief total cross section for all energies */

    // quadrature related numbers
40
41
42
    QuadratureBase* _quadrature; /*! @brief quadrature to create members below */
    unsigned _nq;                /*! @brief number of quadrature points */

43
44
45
46
47
48
49
50
51
    // Mesh related members
    unsigned _nCells;                          /*! @brief number of spatial cells */
    std::vector<BOUNDARY_TYPE> _boundaryCells; /*! boundary type for all cells, dim(_boundary) = (_NCells) */
    std::vector<double> _areas;                /*! @brief surface area of all spatial cells, dim(_areas) = _NCells */
    /*! @brief edge normals multiplied by edge length, dim(_normals) = (_NCells,nEdgesPerCell,spatialDim) */
    std::vector<std::vector<Vector>> _normals;
    /*! @brief edge neighbor cell ids, dim(_neighbors) = (_NCells,nEdgesPerCell) */
    std::vector<std::vector<unsigned>> _neighbors;

Tianbai Xiao's avatar
Tianbai Xiao committed
52
    // slope related params
53
54
55
56
57
58
    Reconstructor* _reconstructor; /*! @brief reconstructor object for high-order scheme */
    unsigned _reconsOrder; /*! @brief reconstruction order (current: 1 & 2) */
    VectorVector _psiDx; /*! @brief slope of solutions in X direction */
    VectorVector _psiDy; /*! @brief slope of solutions in Y direction */
    VectorVector _cellMidPoints; /*! @brief middle point locations of elements */
    std::vector<std::vector<Vector>> _interfaceMidPoints; /*! @brief middle point locations of edges */
59

60
61
    // Solution related members
    VectorVector _sol;                 /*! @brief solution of the PDE, e.g. angular flux or moments */
62
    std::vector<double> _solverOutput; /*! @brief LEGACY: Outputfield for solver ==> Will be replaced by _outputFields in the near future */
63

64
    // Output related members
65
    std::vector<std::vector<std::vector<double>>> _outputFields; /*! @brief: Solver Output: dimensions (GroupID,FieldID,CellID).*/
66
    std::vector<std::vector<std::string>> _outputFieldNames;     /*! @brief: Names of the outputFields: dimensions (GroupID,FieldID) */
67
    // we will have to add a further dimension for quadPoints and weights once we start with multilevel SN
68

69
70
71
72
    // Output related members
    std::vector<double> _screenOutputFields;          /*! @brief: Solver Output: dimensions (FieldID). */
    std::vector<std::string> _screenOutputFieldNames; /*! @brief: Names of the outputFields: dimensions (FieldID) */

73
74
75
76
    // Output related members
    std::vector<double> _historyOutputFields;          /*! @brief: Solver Output: dimensions (FieldID). */
    std::vector<std::string> _historyOutputFieldNames; /*! @brief: Names of the outputFields: dimensions (FieldID) */

77
78
79
80
81
    // Internal Members
    VectorVector _solNew; /*! @brief: VectorVector to store the new flux and later the new solution per iteration */    // REPLACES psiNEW
    Vector _fluxNew; /*! @brief: Vector to store the new Flux */
    Vector _flux;    /*! @brief: Vector to store the old Flux */

82
83
    // ---- Member functions ----

84
    // Solver
85
86
    /*! @brief Performs preprocessing for the current solver iteration */
    virtual void IterPreprocessing() = 0;
87
88
    /*! @brief Performs postprocessing for the current solver iteration */
    virtual void IterPostprocessing() = 0;
89
    /*! @brief Constructs  the flux update for the current iteration and stores it in psiNew*/
90
    virtual void FluxUpdate() = 0;
91
    /*! @brief Computes the finite Volume update step for the current iteration */
92
    virtual void FVMUpdate( unsigned idx_energy ) = 0;
93

94
    // Helper
95
    /*! @brief ComputeTimeStep calculates the maximal stable time step */
96
    double ComputeTimeStep( double cfl ) const;
97
98
99
    /*! @brief: Computes the flux of the solution to check conservation properties */
    virtual void ComputeRadFlux() = 0;

100
101
102
103
104
    // IO
    /*! @brief Initializes the output groups and fields of this solver and names the fields */
    virtual void PrepareVolumeOutput() = 0;
    /*! @brief Function that prepares VTK export and csv export of the current solver iteration  */
    virtual void WriteVolumeOutput( unsigned iteration ) = 0;
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
105
106
    /*! @brief Save Output solution at given energy (pseudo time) to VTK file. Write frequency is given by
               option VOLUME_OUTPUT_FREQUENCY. Always prints last iteration without iteration affix.*/
107
    void PrintVolumeOutput( int currEnergy ) const;
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
108
    /*! @brief: Initialized the output fields and their Names for the screenoutput */
109
    void PrepareScreenOutput();
110
111
    /*! @brief Function that writes screen and history output fields */
    void WriteScalarOutput( unsigned iteration );
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
112
113
    /*! @brief Prints ScreenOutputFields to Screen and to logger. Write frequency is given by
               option SCREEN_OUTPUT_FREQUENCY. Always prints last iteration. */
114
    void PrintScreenOutput( unsigned iteration );
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
115
116
    /*! @brief: Initialized the historyOutputFields and their Names for history output. Write frequency is given by
               option HISTORY_OUTPUT_FREQUENCY. Always prints last iteration. */
117
    void PrepareHistoryOutput();
118
119
    /*! @brief Prints HistoryOutputFields to logger */
    void PrintHistoryOutput( unsigned iteration );
120

121
  public:
122
123
    /*! @brief Solver constructor
     *  @param settings stores all needed information */
124
    Solver( Config* settings );
Jonas Kusch's avatar
Jonas Kusch committed
125

126
127
    ~Solver();

128
129
130
    /*! @brief Create constructor
     *  @param settings stores all needed information
     *  @return pointer to Solver */
131
    static Solver* Create( Config* settings );
Jonas Kusch's avatar
Jonas Kusch committed
132

133
    /*! @brief Solve functions runs main time loop */
134
    virtual void Solve();
135

136
    /*! @brief Save Output solution to VTK file */
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
137
    void PrintVolumeOutput() const;    // Only for debugging purposes.
Jonas Kusch's avatar
Jonas Kusch committed
138
};
139
#endif    // SOLVER_H