pnsolver.h 4.47 KB
Newer Older
1
2
3
#ifndef PNSOLVER_H
#define PNSOLVER_H

Tianbai Xiao's avatar
Tianbai Xiao committed
4
#include "solvers/solverbase.h"
5
6
7
8

class PNSolver : public Solver
{
  public:
9
10
    /*! @brief PNSolver constructor
     *  @param settings stores all needed information
11
12
13
     */
    PNSolver( Config* settings );

14
15
16
17
    /*! @brief PNSolver destructor */
    ~PNSolver();

  private:
18
19
    unsigned _nTotalEntries; /*! @brief: total number of equations in the system */
    unsigned _LMaxDegree;    /*! @brief: maximal degree of the spherical harmonics basis*/
20
21
22

    // System Matrix for x, y and z flux
    //    ==> not needed after computation of A+ and A- ==> maybe safe only temporarly and remove as member?
23
24
25
    SymMatrix _Ax; /*! @brief:  Flux Jacbioan in x direction */
    SymMatrix _Ay; /*! @brief:  Flux Jacbioan in x direction */
    SymMatrix _Az; /*! @brief:  Flux Jacbioan in x direction */
26
27

    // Upwinding Matrices
28
29
30
31
32
33
34
35
36
37
    Matrix _AxPlus;  /*! @brief:  Flux Jacbioan in x direction, positive part */
    Matrix _AxMinus; /*! @brief:  Flux Jacbioan in x direction, negative part */
    Matrix _AxAbs;   /*! @brief:  Flux Jacbioan in x direction, absolute part */
    Matrix _AyPlus;  /*! @brief:  Flux Jacbioan in y direction, positive part */
    Matrix _AyMinus; /*! @brief:  Flux Jacbioan in y direction, negative part */
    Matrix _AyAbs;   /*! @brief:  Flux Jacbioan in y direction, absolute part */
    Matrix _AzPlus;  /*! @brief:  Flux Jacbioan in z direction, positive part */
    Matrix _AzMinus; /*! @brief:  Flux Jacbioan in z direction, negative part */
    Matrix _AzAbs;   /*! @brief:  Flux Jacbioan in z direction, absolute part */

38
39
    Vector _scatterMatDiag; /*! @brief: diagonal of the scattering matrix (its a diagonal matrix by construction). Contains eigenvalues of the
                               scattering kernel.  */
40

Tianbai Xiao's avatar
Tianbai Xiao committed
41
42
    VectorVector _solDx; /*! @brief:  temporary storage of x-derivatives of solution */
    VectorVector _solDy; /*! @brief:  temporary storage of y-derivatives of solution */
Tianbai Xiao's avatar
Tianbai Xiao committed
43

44
45
46
47
    // ---- Member functions ----

    // IO
    /*! @brief Initializes the output groups and fields of this solver and names the fields */
48
    void PrepareVolumeOutput() override;
49
    /*! @brief Function that prepares VTK export and csv export of the current solver iteration
50
        @returns: Mass of current iteration     */
51
    void WriteVolumeOutput( unsigned idx_pseudoTime ) override;
52
53

    // Solver
54
55
    void FVMUpdate( unsigned idx_energy ) override;
    void FluxUpdate() override;
56
    void IterPreprocessing() override;
57
58
59
60
    void IterPostprocessing();

    // Helper
    void ComputeRadFlux();
61

62
    // Initialization of the Solver
63
64
65
66
    /*! @brief: parameter functions for setting up system matrix
     *  @param: degree l, it must hold: 0 <= l <=_nq
     *  @param : order k, it must hold: -l <=k <= l
     */
67
68
69
70
71
72
73
74
75
76
77
78
    double AParam( int l, int k ) const;
    double BParam( int l, int k ) const;
    double CParam( int l, int k ) const;
    double DParam( int l, int k ) const;
    double EParam( int l, int k ) const;
    double FParam( int l, int k ) const;

    double CTilde( int l, int k ) const;
    double DTilde( int l, int k ) const;
    double ETilde( int l, int k ) const;
    double FTilde( int l, int k ) const;

79
80
81
    /*! @brief: mathematical + index functions. Helper functions for setting up system matrix.
     *  @param: k: arbitrary integer
     */
82
83
84
    int Sgn( int k ) const;
    int kPlus( int k ) const;
    int kMinus( int k ) const;
steffen.schotthoefer's avatar
steffen.schotthoefer committed
85
86
87
88
89
90

    /*! @brief : computes the global index of the moment corresponding to basis function (l,k)
     *  @param : degree l, it must hold: 0 <= l <=_nq
     *  @param : order k, it must hold: -l <=k <= l
     *  @returns : global index
     */
91
    int GlobalIndex( int l, int k ) const;
92

steffen.schotthoefer's avatar
steffen.schotthoefer committed
93
94
95
96
97
    /*! @brief : Checks, if index invariant for global index holds
     *  @returns : True, if invariant holds, else false
     */
    bool CheckIndex( int l, int k ) const;

98
    /*! @brief: function for computing and setting up system matrices */
99
    void ComputeSystemMatrices();
100
    /*! @brief:  function for computing and setting up flux matrices for upwinding */
101
    void ComputeFluxComponents();
102
    /*! @brief:  fucntion for computing and setting up EV matrix for scattering kernel */
103
    void ComputeScatterMatrix();
104
    /*! @brief:  Computes Legedre polinomial of oder l at point x */
105
    double LegendrePoly( double x, int l );
106
107
108
109

    /*! @brief: Sets Entries of FluxMatrices to zero, if they are below double precision,
     *          to prevent floating point inaccuracies later in the solver
     */
110
    void CleanFluxMatrices();
111
};
112
113

#endif    // PNSOLVER_H