Commit 2e81926b authored by Steffen Schotthöfer's avatar Steffen Schotthöfer
Browse files

clean up solvers... replaced temp variable psiNew with _solNEw


Former-commit-id: fe81a279
parent ce972d95
...@@ -52,8 +52,8 @@ class MNSolver : public Solver ...@@ -52,8 +52,8 @@ class MNSolver : public Solver
void WriteOutputFields( unsigned idx_pseudoTime ) override; void WriteOutputFields( unsigned idx_pseudoTime ) override;
// Solver // Solver
void FVMUpdate( VectorVector& psiNew, unsigned idx_energy ) override; void FVMUpdate( unsigned idx_energy ) override;
void FluxUpdate( VectorVector& psiNew ) override; void FluxUpdate() override;
void IterPreprocessing() override; void IterPreprocessing() override;
void IterPostprocessing(); void IterPostprocessing();
/*! @brief : Construct flux by computing the Moment of the sum of FVM discretization at the interface of cell /*! @brief : Construct flux by computing the Moment of the sum of FVM discretization at the interface of cell
......
...@@ -48,8 +48,8 @@ class PNSolver : public Solver ...@@ -48,8 +48,8 @@ class PNSolver : public Solver
void WriteOutputFields( unsigned idx_pseudoTime ) override; void WriteOutputFields( unsigned idx_pseudoTime ) override;
// Solver // Solver
void FVMUpdate( VectorVector& psiNew, unsigned idx_energy ) override; void FVMUpdate( unsigned idx_energy ) override;
void FluxUpdate( VectorVector& psiNew ) override; void FluxUpdate() override;
void IterPreprocessing() override; void IterPreprocessing() override;
void IterPostprocessing(); void IterPostprocessing();
......
...@@ -25,8 +25,8 @@ class SNSolver : public Solver ...@@ -25,8 +25,8 @@ class SNSolver : public Solver
void WriteOutputFields( unsigned idx_pseudoTime ) override; void WriteOutputFields( unsigned idx_pseudoTime ) override;
// Solver // Solver
void FVMUpdate( VectorVector& psiNew, unsigned idx_energy ) override; void FVMUpdate( unsigned idx_energy ) override;
void FluxUpdate( VectorVector& psiNew ) override; void FluxUpdate() override;
void IterPreprocessing() override; void IterPreprocessing() override;
void IterPostprocessing() override; void IterPostprocessing() override;
......
...@@ -89,9 +89,9 @@ class Solver ...@@ -89,9 +89,9 @@ class Solver
/*! @brief Performs postprocessing for the current solver iteration */ /*! @brief Performs postprocessing for the current solver iteration */
virtual void IterPostprocessing() = 0; virtual void IterPostprocessing() = 0;
/*! @brief Constructs the flux update for the current iteration and stores it in psiNew*/ /*! @brief Constructs the flux update for the current iteration and stores it in psiNew*/
virtual void FluxUpdate( VectorVector& psiNew ) = 0; virtual void FluxUpdate() = 0;
/*! @brief Computes the finite Volume update step for the current iteration */ /*! @brief Computes the finite Volume update step for the current iteration */
virtual void FVMUpdate( VectorVector& psiNew, unsigned idx_energy ) = 0; virtual void FVMUpdate( unsigned idx_energy ) = 0;
// Helper // Helper
/** /**
......
...@@ -152,16 +152,16 @@ void MNSolver::ComputeRadFlux() { ...@@ -152,16 +152,16 @@ void MNSolver::ComputeRadFlux() {
} }
} }
void MNSolver::FluxUpdate( VectorVector& psiNew ) { void MNSolver::FluxUpdate() {
// Loop over the grid cells // Loop over the grid cells
for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) { for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) {
// Dirichlet Boundaries stay // Dirichlet Boundaries stay
if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) continue; if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) continue;
psiNew[idx_cell] = ConstructFlux( idx_cell ); _solNew[idx_cell] = ConstructFlux( idx_cell );
} }
} }
void MNSolver::FVMUpdate( VectorVector& psiNew, unsigned idx_energy ) { void MNSolver::FVMUpdate( unsigned idx_energy ) {
// Loop over the grid cells // Loop over the grid cells
for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) { for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) {
// Dirichlet Boundaries stay // Dirichlet Boundaries stay
...@@ -169,14 +169,14 @@ void MNSolver::FVMUpdate( VectorVector& psiNew, unsigned idx_energy ) { ...@@ -169,14 +169,14 @@ void MNSolver::FVMUpdate( VectorVector& psiNew, unsigned idx_energy ) {
for( unsigned idx_system = 0; idx_system < _nTotalEntries; idx_system++ ) { for( unsigned idx_system = 0; idx_system < _nTotalEntries; idx_system++ ) {
psiNew[idx_cell][idx_system] = _sol[idx_cell][idx_system] - _solNew[idx_cell][idx_system] = _sol[idx_cell][idx_system] -
( _dE / _areas[idx_cell] ) * psiNew[idx_cell][idx_system] /* cell averaged flux */ ( _dE / _areas[idx_cell] ) * _solNew[idx_cell][idx_system] /* cell averaged flux */
- _dE * _sol[idx_cell][idx_system] * - _dE * _sol[idx_cell][idx_system] *
( _sigmaT[idx_energy][idx_cell] /* absorbtion influence */ ( _sigmaT[idx_energy][idx_cell] /* absorbtion influence */
+ _sigmaS[idx_energy][idx_cell] * _scatterMatDiag[idx_system] ); /* scattering influence */ + _sigmaS[idx_energy][idx_cell] * _scatterMatDiag[idx_system] ); /* scattering influence */
} }
psiNew[idx_cell][0] += _dE * _Q[0][idx_cell][0]; _solNew[idx_cell][0] += _dE * _Q[0][idx_cell][0];
} }
} }
......
...@@ -70,7 +70,7 @@ void PNSolver::ComputeRadFlux() { ...@@ -70,7 +70,7 @@ void PNSolver::ComputeRadFlux() {
} }
} }
void PNSolver::FluxUpdate( VectorVector& psiNew ) { void PNSolver::FluxUpdate() {
// Loop over all spatial cells // Loop over all spatial cells
for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) { for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) {
...@@ -79,7 +79,7 @@ void PNSolver::FluxUpdate( VectorVector& psiNew ) { ...@@ -79,7 +79,7 @@ void PNSolver::FluxUpdate( VectorVector& psiNew ) {
// Reset temporary variable psiNew // Reset temporary variable psiNew
for( unsigned idx_sys = 0; idx_sys < _nTotalEntries; idx_sys++ ) { for( unsigned idx_sys = 0; idx_sys < _nTotalEntries; idx_sys++ ) {
psiNew[idx_cell][idx_sys] = 0.0; _solNew[idx_cell][idx_sys] = 0.0;
} }
// Loop over all neighbor cells (edges) of cell j and compute numerical fluxes // Loop over all neighbor cells (edges) of cell j and compute numerical fluxes
...@@ -87,35 +87,36 @@ void PNSolver::FluxUpdate( VectorVector& psiNew ) { ...@@ -87,35 +87,36 @@ void PNSolver::FluxUpdate( VectorVector& psiNew ) {
// Compute flux contribution and store in psiNew to save memory // Compute flux contribution and store in psiNew to save memory
if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::NEUMANN && _neighbors[idx_cell][idx_neighbor] == _nCells ) if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::NEUMANN && _neighbors[idx_cell][idx_neighbor] == _nCells )
psiNew[idx_cell] += _g->Flux( _solNew[idx_cell] += _g->Flux(
_AxPlus, _AxMinus, _AyPlus, _AyMinus, _AzPlus, _AzMinus, _sol[idx_cell], _sol[idx_cell], _normals[idx_cell][idx_neighbor] ); _AxPlus, _AxMinus, _AyPlus, _AyMinus, _AzPlus, _AzMinus, _sol[idx_cell], _sol[idx_cell], _normals[idx_cell][idx_neighbor] );
else else
psiNew[idx_cell] += _g->Flux( _AxPlus, _solNew[idx_cell] += _g->Flux( _AxPlus,
_AxMinus, _AxMinus,
_AyPlus, _AyPlus,
_AyMinus, _AyMinus,
_AzPlus, _AzPlus,
_AzMinus, _AzMinus,
_sol[idx_cell], _sol[idx_cell],
_sol[_neighbors[idx_cell][idx_neighbor]], _sol[_neighbors[idx_cell][idx_neighbor]],
_normals[idx_cell][idx_neighbor] ); _normals[idx_cell][idx_neighbor] );
} }
} }
} }
void PNSolver::FVMUpdate( VectorVector& psiNew, unsigned idx_energy ) { void PNSolver::FVMUpdate( unsigned idx_energy ) {
// Loop over all spatial cells // Loop over all spatial cells
for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) { for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) {
// Dirichlet cells stay at IC, farfield assumption // Dirichlet cells stay at IC, farfield assumption
if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) continue; if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) continue;
// Flux update // Flux update
for( unsigned idx_sys = 0; idx_sys < _nTotalEntries; idx_sys++ ) { for( unsigned idx_sys = 0; idx_sys < _nTotalEntries; idx_sys++ ) {
psiNew[idx_cell][idx_sys] = _sol[idx_cell][idx_sys] - ( _dE / _areas[idx_cell] ) * psiNew[idx_cell][idx_sys] /* cell averaged flux */ _solNew[idx_cell][idx_sys] = _sol[idx_cell][idx_sys] - ( _dE / _areas[idx_cell] ) * _solNew[idx_cell][idx_sys] /* cell averaged flux */
- _dE * _sol[idx_cell][idx_sys] * - _dE * _sol[idx_cell][idx_sys] *
( _sigmaT[idx_energy][idx_cell] /* absorbtion influence */ ( _sigmaT[idx_energy][idx_cell] /* absorbtion influence */
+ _sigmaS[idx_energy][idx_cell] * _scatterMatDiag[idx_sys] ); /* scattering influence */ + _sigmaS[idx_energy][idx_cell] * _scatterMatDiag[idx_sys] ); /* scattering influence */
} }
psiNew[idx_cell][0] += _dE * _Q[0][idx_cell][0]; // Source Term
_solNew[idx_cell][0] += _dE * _Q[0][idx_cell][0];
} }
} }
......
...@@ -42,7 +42,7 @@ void SNSolver::ComputeRadFlux() { ...@@ -42,7 +42,7 @@ void SNSolver::ComputeRadFlux() {
} }
} }
void SNSolver::FluxUpdate( VectorVector& psiNew ) { void SNSolver::FluxUpdate() {
// Legacy Code form second order reconstruction: // Legacy Code form second order reconstruction:
...@@ -104,12 +104,12 @@ void SNSolver::FluxUpdate( VectorVector& psiNew ) { ...@@ -104,12 +104,12 @@ void SNSolver::FluxUpdate( VectorVector& psiNew ) {
// Loop over all ordinates // Loop over all ordinates
for( unsigned idx_quad = 0; idx_quad < _nq; ++idx_quad ) { for( unsigned idx_quad = 0; idx_quad < _nq; ++idx_quad ) {
// Reset temporary variable // Reset temporary variable
psiNew[idx_cell][idx_quad] = 0.0; _solNew[idx_cell][idx_quad] = 0.0;
// Loop over all neighbor cells (edges) of cell j and compute numerical fluxes // Loop over all neighbor cells (edges) of cell j and compute numerical fluxes
for( unsigned idx_neighbor = 0; idx_neighbor < _neighbors[idx_cell].size(); ++idx_neighbor ) { for( unsigned idx_neighbor = 0; idx_neighbor < _neighbors[idx_cell].size(); ++idx_neighbor ) {
// store flux contribution on psiNew_sigmaS to save memory // store flux contribution on psiNew_sigmaS to save memory
if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::NEUMANN && _neighbors[idx_cell][idx_neighbor] == _nCells ) if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::NEUMANN && _neighbors[idx_cell][idx_neighbor] == _nCells )
psiNew[idx_cell][idx_quad] += _solNew[idx_cell][idx_quad] +=
_g->Flux( _quadPoints[idx_quad], _sol[idx_cell][idx_quad], _sol[idx_cell][idx_quad], _normals[idx_cell][idx_neighbor] ); _g->Flux( _quadPoints[idx_quad], _sol[idx_cell][idx_quad], _sol[idx_cell][idx_quad], _normals[idx_cell][idx_neighbor] );
else { else {
/* switch( reconsOrder ) { /* switch( reconsOrder ) {
...@@ -146,10 +146,10 @@ void SNSolver::FluxUpdate( VectorVector& psiNew ) { ...@@ -146,10 +146,10 @@ void SNSolver::FluxUpdate( VectorVector& psiNew ) {
// default: first order solver // default: first order solver
default: default:
*/ */
psiNew[idx_cell][idx_quad] += _g->Flux( _quadPoints[idx_quad], _solNew[idx_cell][idx_quad] += _g->Flux( _quadPoints[idx_quad],
_sol[idx_cell][idx_quad], _sol[idx_cell][idx_quad],
_sol[_neighbors[idx_cell][idx_neighbor]][idx_quad], _sol[_neighbors[idx_cell][idx_neighbor]][idx_quad],
_normals[idx_cell][idx_neighbor] ); _normals[idx_cell][idx_neighbor] );
// } // }
} }
} }
...@@ -157,31 +157,31 @@ void SNSolver::FluxUpdate( VectorVector& psiNew ) { ...@@ -157,31 +157,31 @@ void SNSolver::FluxUpdate( VectorVector& psiNew ) {
} }
} }
void SNSolver::FVMUpdate( VectorVector& psiNew, unsigned idx_energy ) { void SNSolver::FVMUpdate( unsigned idx_energy ) {
for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) { for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) {
// Dirichlet cells stay at IC, farfield assumption // Dirichlet cells stay at IC, farfield assumption
if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) continue; if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) continue;
// loop over all ordinates // loop over all ordinates
for( unsigned idx_quad = 0; idx_quad < _nq; ++idx_quad ) { for( unsigned idx_quad = 0; idx_quad < _nq; ++idx_quad ) {
// time update angular flux with numerical flux and total scattering cross section // time update angular flux with numerical flux and total scattering cross section
psiNew[idx_cell][idx_quad] = _sol[idx_cell][idx_quad] - ( _dE / _areas[idx_cell] ) * psiNew[idx_cell][idx_quad] - _solNew[idx_cell][idx_quad] = _sol[idx_cell][idx_quad] - ( _dE / _areas[idx_cell] ) * _solNew[idx_cell][idx_quad] -
_dE * _sigmaT[idx_energy][idx_cell] * _sol[idx_cell][idx_quad]; _dE * _sigmaT[idx_energy][idx_cell] * _sol[idx_cell][idx_quad];
} }
// compute scattering effects // compute scattering effects
psiNew[idx_cell] += _dE * _sigmaS[idx_energy][idx_cell] * _scatteringKernel * _sol[idx_cell]; // multiply scattering matrix with psi _solNew[idx_cell] += _dE * _sigmaS[idx_energy][idx_cell] * _scatteringKernel * _sol[idx_cell]; // multiply scattering matrix with psi
// Source Term // Source Term
if( _Q.size() == 1u ) { // constant source for all energies if( _Q.size() == 1u ) { // constant source for all energies
if( _Q[0][idx_cell].size() == 1u ) // isotropic source if( _Q[0][idx_cell].size() == 1u ) // isotropic source
psiNew[idx_cell] += _dE * _Q[0][idx_cell][0]; _solNew[idx_cell] += _dE * _Q[0][idx_cell][0];
else else
psiNew[idx_cell] += _dE * _Q[0][idx_cell]; _solNew[idx_cell] += _dE * _Q[0][idx_cell];
} }
else { else {
if( _Q[0][idx_cell].size() == 1u ) // isotropic source if( _Q[0][idx_cell].size() == 1u ) // isotropic source
psiNew[idx_cell] += _dE * _Q[idx_energy][idx_cell][0]; _solNew[idx_cell] += _dE * _Q[idx_energy][idx_cell][0];
else else
psiNew[idx_cell] += _dE * _Q[idx_energy][idx_cell]; _solNew[idx_cell] += _dE * _Q[idx_energy][idx_cell];
} }
} }
} }
......
...@@ -164,10 +164,10 @@ void Solver::Solve() { ...@@ -164,10 +164,10 @@ void Solver::Solve() {
IterPreprocessing(); IterPreprocessing();
// --- Compute Fluxes --- // --- Compute Fluxes ---
FluxUpdate( _solNew ); FluxUpdate();
// --- Finite Volume Update --- // --- Finite Volume Update ---
FVMUpdate( _solNew, idx_energy ); FVMUpdate( idx_energy );
// --- Postprocessing --- // --- Postprocessing ---
IterPostprocessing(); IterPostprocessing();
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment