An open source method of characteristics neutron transport code.
CPULSSolver Class Reference

This a subclass of the CPUSolver class for using the linear source approximation. More...

#include "src/CPULSSolver.h"

Public Member Functions

 CPULSSolver (TrackGenerator *track_generator=NULL)
 Constructor initializes array pointers for fluxes and sources. More...
 
virtual ~CPULSSolver ()
 Destructor deletes array for linear fluxes, sources and constants. CPUSolver parent class destructor handles deletion of arrays for flat fluxes and sources.
 
void initializeFluxArrays ()
 Allocates memory for boundary and scalar fluxes. More...
 
void initializeSourceArrays ()
 Allocates memory for FSR source arrays. More...
 
void initializeCmfd ()
 Initializes a Cmfd object for acceleration prior to source iteration. More...
 
void initializeExpEvaluators ()
 Initializes new ExpEvaluator objects to compute exponentials. More...
 
void initializeFSRs ()
 Initializes the FSR constant linear source component, volumes and Materials array. More...
 
void initializeFixedSources ()
 Initializes the arrays for the fixed source and its moments.
 
void setFixedSourceMomentsByCell (Cell *cell, int group, double source_x, double source_y, double source_z)
 Assign fixed source moments for a Cell and energy group. More...
 
void setFixedSourceMomentByFSR (long fsr_id, int group, double source_x, double source_y, double source_z)
 Assign fixed source moments for a FSR and energy group. More...
 
void flattenFSRFluxes (FP_PRECISION value)
 Set the scalar flux constants for each FSR and energy group to some value and the scalar flux moments to zero. More...
 
double normalizeFluxes ()
 Normalizes all FSR scalar fluxes and Track boundary angular fluxes to the total fission source (times $ \nu $). More...
 
void computeFSRSources (int iteration)
 Computes the total source (fission, scattering, fixed) in each FSR. More...
 
void tallyLSScalarFlux (segment *curr_segment, int azim_index, int polar_index, FP_PRECISION *fsr_flux, FP_PRECISION *fsr_flux_x, FP_PRECISION *fsr_flux_y, FP_PRECISION *fsr_flux_z, float *track_flux, FP_PRECISION direction[3])
 Computes the contribution to the LSR scalar flux from a Track segment. More...
 
void accumulateLinearFluxContribution (long fsr_id, FP_PRECISION weight, FP_PRECISION *fsr_flux)
 Move from buffers to global arrays the contributions of one (or several for per-stack solving) segments. More...
 
void addSourceToScalarFlux ()
 Add the source term contribution in the transport equation to the FSR scalar flux.
 
void computeStabilizingFlux ()
 Computes the stabilizing flux for transport stabilization.
 
void stabilizeFlux ()
 Adjusts the scalar flux for transport stabilization.
 
void checkLimitXS (int iteration)
 Checks to see if limited XS should be reset. More...
 
void initializeLinearSourceConstants ()
 Initialize linear source constant component and matrix coefficients.
 
double * getLinearExpansionCoeffsBuffer ()
 Returns a memory buffer to the linear source expansion coefficent matrix. More...
 
FP_PRECISION * getSourceConstantsBuffer ()
 Returns a memory buffer to the constant part (constant between MOC iterations) of the linear source. More...
 
FP_PRECISION getFluxByCoords (LocalCoords *coords, int group)
 Get the flux at a specific point in the geometry. More...
 
- Public Member Functions inherited from CPUSolver
 CPUSolver (TrackGenerator *track_generator=NULL)
 Constructor initializes array pointers for Tracks and Materials. More...
 
virtual ~CPUSolver ()
 Destructor deletes array for OpenMP mutual exclusion locks for FSR scalar flux updates, and calls Solver parent class destructor to deletes arrays for fluxes and sources.
 
int getNumThreads ()
 Returns the number of shared memory OpenMP threads in use. More...
 
void setNumThreads (int num_threads)
 Sets the number of shared memory OpenMP threads to use (>0). More...
 
void setFluxes (FP_PRECISION *in_fluxes, int num_fluxes)
 Set the flux array for use in transport sweep source calculations. More...
 
void setFixedSourceByFSR (long fsr_id, int group, FP_PRECISION source)
 Assign a fixed source for a flat source region and energy group. More...
 
void computeFSRFissionRates (double *fission_rates, long num_FSRs, bool nu=false)
 Computes the volume-averaged, energy-integrated fission or nu-fission rate in each FSR and stores them in an array indexed by FSR ID. More...
 
void printInputParamsSummary ()
 A function that prints a summary of the input parameters.
 
void tallyScalarFlux (segment *curr_segment, int azim_index, int polar_index, FP_PRECISION *fsr_flux, float *track_flux)
 Computes the contribution to the FSR scalar flux from a Track segment. More...
 
void accumulateScalarFluxContribution (long fsr_id, FP_PRECISION weight, FP_PRECISION *fsr_flux)
 Move the segment(s)' contributions to the scalar flux from the buffer to the global scalar flux array. More...
 
void tallyCurrent (segment *curr_segment, int azim_index, int polar_index, float *track_flux, bool fwd)
 Tallies the current contribution from this segment across the the appropriate CMFD mesh cell surface. More...
 
void transferBoundaryFlux (Track *track, int azim_index, int polar_index, bool direction, float *track_flux)
 Updates the boundary flux for a Track given boundary conditions. More...
 
void getFluxes (FP_PRECISION *out_fluxes, int num_fluxes)
 Fills an array with the scalar fluxes. More...
 
void initializeFixedSources ()
 Populates array of fixed sources assigned by FSR.
 
void printFSRFluxes (std::vector< double > dim1, std::vector< double > dim2, double offset, const char *plane)
 A function that prints the source region fluxes on a 2D mesh grid. More...
 
void printFluxesTemp ()
 A function that prints fsr fluxes in xy plane at z=middle. More...
 
void printNegativeSources (int iteration, int num_x, int num_y, int num_z)
 A function that prints the number of FSRs with negative sources in the whole geometry subdivided by a 3D lattice. The number of negative sources per energy group is also printed out. More...
 
- Public Member Functions inherited from Solver
 Solver (TrackGenerator *track_generator=NULL)
 Constructor initializes an empty Solver class with array pointers set to NULL. More...
 
virtual ~Solver ()
 Destructor deletes arrays of boundary angular fluxes, scalar fluxes and sources for each FSR and energy group. More...
 
void setGeometry (Geometry *geometry)
 Sets the Geometry for the Solver. More...
 
GeometrygetGeometry ()
 Returns a pointer to the Geometry. More...
 
TrackGeneratorgetTrackGenerator ()
 Returns a pointer to the TrackGenerator. More...
 
FP_PRECISION getFSRVolume (long fsr_id)
 Returns the calculated volume for a flat source region. More...
 
int getNumPolarAngles ()
 Returns the number of angles used for the polar quadrature. More...
 
int getNumIterations ()
 Returns the number of source iterations to converge the source. More...
 
double getTotalTime ()
 Returns the total time to converge the source (seconds). More...
 
double getKeff ()
 Returns the converged eigenvalue $ k_{eff} $. More...
 
double getConvergenceThreshold ()
 Returns the threshold for source/flux convergence. More...
 
FP_PRECISION getMaxOpticalLength ()
 Get the maximum allowable optical length for a track segment. More...
 
bool isUsingDoublePrecision ()
 Returns whether the solver is using double floating point precision. More...
 
bool isUsingExponentialInterpolation ()
 Returns whether the Solver uses linear interpolation to compute exponentials. More...
 
bool is3D ()
 Returns whether the Solver is tackling a 3D problem. More...
 
void initializeSolver (solverMode solver_mode)
 Initializes most components of Solver. Mostly needed from the Python side.
 
void printFissionRates (std::string fname, int nx, int ny, int nz)
 Prints fission rates to a binary file. More...
 
void setKeffFromNeutronBalance ()
 Sets computation method of k-eff from fission, absorption, and leakage rates rather than from fission rates. keff = fission/(absorption + leakage)
 
void setResidualByReference (std::string fname)
 Sets residuals to be computed a error relative to a reference. More...
 
void dumpFSRFluxes (std::string fname)
 Prints scalar fluxes to a binary file. More...
 
void loadInitialFSRFluxes (std::string fname)
 Load the initial scalar flux distribution from a binary file. More...
 
void loadFSRFluxes (std::string fname, bool assign_k_eff=false, double tolerance=0.01)
 Load scalar fluxes from a binary file. More...
 
double getFlux (long fsr_id, int group)
 Returns the scalar flux for some FSR and energy group. More...
 
double getFSRSource (long fsr_id, int group)
 Returns the source for some energy group for a flat source region. More...
 
void setTrackGenerator (TrackGenerator *track_generator)
 Sets the Solver's TrackGenerator with characteristic Tracks. More...
 
void setConvergenceThreshold (double threshold)
 Sets the threshold for source/flux convergence. More...
 
virtual void setFixedSourceByFSR (long fsr_id, int group, double source)
 Assign a fixed source for a flat source region and energy group. More...
 
void setFixedSourceByCell (Cell *cell, int group, double source)
 Assign a fixed source for a Cell and energy group. More...
 
void setFixedSourceByMaterial (Material *material, int group, double source)
 Assign a fixed source for a Material and energy group. More...
 
void setMaxOpticalLength (FP_PRECISION max_optical_length)
 Set the maximum allowable optical length for a track segment. More...
 
void setExpPrecision (double precision)
 Set the precision, or maximum allowable approximation error, of the the exponential interpolation table. More...
 
void useExponentialInterpolation ()
 Informs the Solver to use linear interpolation to compute the exponential in the transport equation.
 
void useExponentialIntrinsic ()
 Informs the Solver to use the exponential intrinsic exp(...) function to compute the exponential in the transport equation.
 
void setSolverMode (solverMode solver_mode)
 Choose between direct and adjoint mode. More...
 
void setRestartStatus (bool is_restart)
 Informs the Solver that this is a 'restart' calculation and therefore k_eff, track angular and region scalar fluxes should not be reset. More...
 
void allowNegativeFluxes (bool negative_fluxes_on)
 Informs the Solver that this calculation may involve negative fluxes for computing higher eigenmodes for example. More...
 
void correctXS ()
 Directs OpenMOC to correct unphysical cross-sections. More...
 
void stabilizeTransport (double stabilization_factor, stabilizationType stabilization_type=DIAGONAL)
 Directs OpenMOC to use the diagonal stabilizing correction to the source iteration transport sweep. More...
 
void setInitialSpectrumCalculation (double threshold)
 Instructs OpenMOC to perform an initial spectrum calculation. More...
 
void setCheckXSLogLevel (logLevel log_level)
 Determines which log level to set cross-section warnings. More...
 
void setChiSpectrumMaterial (Material *material)
 Sets the chi spectrum for use as an inital flux guess. More...
 
void resetMaterials (solverMode mode)
 Returns the Material data to its original state. More...
 
void fissionTransportSweep ()
 This method performs one transport sweep using the fission source. More...
 
void scatterTransportSweep ()
 This method performs one transport sweep using the scatter source. More...
 
void computeFlux (int max_iters=1000, bool only_fixed_source=true)
 Computes the scalar flux distribution by performing a series of transport sweeps. More...
 
void computeSource (int max_iters=1000, double k_eff=1.0, residualType res_type=TOTAL_SOURCE)
 Computes the total source distribution by performing a series of transport sweep and source updates. More...
 
void computeEigenvalue (int max_iters=1000, residualType res_type=FISSION_SOURCE)
 Computes keff by performing a series of transport sweep and source updates. More...
 
void printBGQMemory ()
 
float * getBoundaryFlux (long track_id, bool fwd)
 Returns the boundary flux array at the requested indexes. More...
 
void setVerboseIterationReport ()
 Sets the solver to print extra information for each iteration.
 
void printTimerReport ()
 Prints a report of the timing statistics to the console.
 
FP_PRECISION * getFluxesArray ()
 A function that returns the array of scalar fluxes. More...
 
void limitXS ()
 Limits cross-sections so that there are no negative cross-sections. More...
 
void setLimitingXSMaterials (std::vector< int > material_ids, int reset_iteration)
 Instructs MOC to limit negative cross-sections for early iterations. More...
 
void checkLimitXS (int iteration)
 Checks to see if limited XS should be reset. More...
 
void setOTFTransport ()
 Activate On-The-Fly transport, to OTF ray-trace and propagate the track angular fluxes at the same time.
 

Protected Attributes

double * _FSR_lin_exp_matrix
 
FP_PRECISION * _FSR_source_constants
 
FP_PRECISION * _scalar_flux_xyz
 
FP_PRECISION * _reduced_sources_xyz
 
FP_PRECISION * _stabilizing_flux_xyz
 
std::vector< std::vector< double > > _fixed_sources_xyz
 
std::map< std::pair< int, int >, std::vector< double > > _fix_src_xyz_FSR_map
 
std::map< std::pair< Cell *, int >, std::vector< double > > _fix_src_xyz_cell_map
 
bool _stabilize_moments
 
bool _fixed_source_moments_on
 
- Protected Attributes inherited from CPUSolver
int _num_threads
 
omp_lock_t * _FSR_locks
 
- Protected Attributes inherited from Solver
int _num_azim
 
int _num_groups
 
long _num_FSRs
 
long _num_fissionable_FSRs
 
FP_PRECISION * _FSR_volumes
 
Material ** _FSR_materials
 
Material_chi_spectrum_material
 
TrackGenerator_track_generator
 
Geometry_geometry
 
int _num_materials
 
Quadrature_quad
 
int _num_polar
 
int _fluxes_per_track
 
Track ** _tracks
 
std::vector< int > _limit_xs_materials
 
int _reset_iteration
 
bool _limit_xs
 
std::map< int, Material * > _original_materials
 
std::map< int, Material * > _limit_materials
 
int *** _tracks_per_stack
 
bool _SOLVE_3D
 
solverMode _solver_mode
 
bool _is_restart
 
bool _user_fluxes
 
bool _fixed_sources_on
 
bool _correct_xs
 
bool _stabilize_transport
 
bool _verbose
 
bool _calculate_initial_spectrum
 
double _initial_spectrum_thresh
 
bool _load_initial_FSR_fluxes
 
bool _calculate_residuals_by_reference
 
bool _negative_fluxes_allowed
 
std::string _initial_FSR_fluxes_file
 
std::string _reference_file
 
logLevel _xs_log_level
 
segmentationType _segment_formation
 
long _tot_num_tracks
 
float * _boundary_flux
 
float * _start_flux
 
float * _boundary_leakage
 
FP_PRECISION * _scalar_flux
 
FP_PRECISION * _old_scalar_flux
 
FP_PRECISION * _reference_flux
 
FP_PRECISION * _stabilizing_flux
 
FP_PRECISION * _fixed_sources
 
std::vector< FP_PRECISION * > _groupwise_scratch
 
double * _regionwise_scratch
 
std::map< std::pair< int, int >, FP_PRECISION > _fix_src_FSR_map
 
std::map< std::pair< Cell *, int >, FP_PRECISION > _fix_src_cell_map
 
std::map< std::pair< Material *, int >, FP_PRECISION > _fix_src_material_map
 
FP_PRECISION * _reduced_sources
 
double _k_eff
 
bool _keff_from_fission_rates
 
int _num_iterations
 
double _converge_thresh
 
double _stabilization_factor
 
stabilizationType _stabilization_type
 
ExpEvaluator *** _exp_evaluators
 
int _num_exp_evaluators_azim
 
int _num_exp_evaluators_polar
 
Timer_timer
 
Cmfd_cmfd
 
std::string _source_type
 
bool _OTF_transport
 

Additional Inherited Members

- Protected Member Functions inherited from CPUSolver
void zeroTrackFluxes ()
 Zero each Track's boundary fluxes for each energy group (and polar angle in 2D) in the "forward" and "reverse" directions.
 
void copyBoundaryFluxes ()
 Copies values from the start flux into the boundary flux array for both the "forward" and "reverse" directions.
 
void tallyStartingCurrents ()
 Computes the total current impingent on boundary CMFD cells from starting angular fluxes.
 
void flattenFSRFluxesChiSpectrum ()
 Set the scalar flux for each FSR to a chi spectrum.
 
void storeFSRFluxes ()
 Stores the FSR scalar fluxes in the old scalar flux array.
 
void computeFSRFissionSources ()
 Computes the total fission source in each FSR. More...
 
void computeFSRScatterSources ()
 Computes the total scattering source in each FSR. More...
 
void transportSweep ()
 This method performs one transport sweep of all azimuthal angles, Tracks, Track segments, polar angles and energy groups. More...
 
void computeKeff ()
 Compute $ k_{eff} $ from successive fission sources.
 
double computeResidual (residualType res_type)
 Computes the residual between source/flux iterations. More...
 
- Protected Member Functions inherited from Solver
void initializeMaterials (solverMode mode)
 Initializes the Material's production matrices. More...
 
void countFissionableFSRs ()
 Counts the number of fissionable flat source regions. More...
 
void checkXS ()
 All material cross-sections in the geometry are checked for consistency. More...
 
void calculateInitialSpectrum (double threshold)
 Performs a spectrum calculation to update the scalar fluxes. More...
 
void clearTimerSplits ()
 Deletes the Timer's timing entries for each timed code section code in the source convergence loop. More...
 

Detailed Description

This a subclass of the CPUSolver class for using the linear source approximation.

Constructor & Destructor Documentation

◆ CPULSSolver()

CPULSSolver::CPULSSolver ( TrackGenerator track_generator = NULL)

Constructor initializes array pointers for fluxes and sources.

The constructor retrieves the number of energy groups and FSRs and azimuthal angles from the Geometry and TrackGenerator if passed in as parameters by the user. The constructor initalizes the number of OpenMP threads to a default of 1.

Parameters
track_generatoran optional pointer to the TrackGenerator

Member Function Documentation

◆ accumulateLinearFluxContribution()

void CPULSSolver::accumulateLinearFluxContribution ( long  fsr_id,
FP_PRECISION  weight,
FP_PRECISION *  fsr_flux 
)

Move from buffers to global arrays the contributions of one (or several for per-stack solving) segments.

Parameters
fsr_idregion index
weightthe quadrature weight (only for 3D ray tracing)
fsr_fluxbuffer storing contribution to the region's scalar flux

◆ checkLimitXS()

void CPULSSolver::checkLimitXS ( int  iteration)

Checks to see if limited XS should be reset.

For the linear source, the linear expansion coefficients should also be reset, to use the non-limited cross sections.

Parameters
iterationThe MOC iteration number

◆ computeFSRSources()

void CPULSSolver::computeFSRSources ( int  iteration)
virtual

Computes the total source (fission, scattering, fixed) in each FSR.

This method computes the total source in each FSR based on this iteration's current approximation to the scalar flux and its moments. Fixed source moments are currently not supported.

Reimplemented from CPUSolver.

◆ flattenFSRFluxes()

void CPULSSolver::flattenFSRFluxes ( FP_PRECISION  value)
virtual

Set the scalar flux constants for each FSR and energy group to some value and the scalar flux moments to zero.

Parameters
valuethe value to assign to each FSR scalar flux

Reimplemented from CPUSolver.

◆ getFluxByCoords()

FP_PRECISION CPULSSolver::getFluxByCoords ( LocalCoords coords,
int  group 
)

Get the flux at a specific point in the geometry.

Parameters
coordsThe coords of the point to get the flux at
groupthe energy group

◆ getLinearExpansionCoeffsBuffer()

double * CPULSSolver::getLinearExpansionCoeffsBuffer ( )

Returns a memory buffer to the linear source expansion coefficent matrix.

Returns
_FSR_lin_exp_matrix a pointer to the linear source coefficient matrix

◆ getSourceConstantsBuffer()

FP_PRECISION * CPULSSolver::getSourceConstantsBuffer ( )

Returns a memory buffer to the constant part (constant between MOC iterations) of the linear source.

Returns
_FSR_source_constants a pointer to the linear source constant part

◆ initializeCmfd()

void CPULSSolver::initializeCmfd ( )
virtual

Initializes a Cmfd object for acceleration prior to source iteration.

For the linear source solver, a pointer to the flux moments is passed to the Cmfd object so that they can be updated as well in the prolongation phase.

Reimplemented from Solver.

◆ initializeExpEvaluators()

void CPULSSolver::initializeExpEvaluators ( )
virtual

Initializes new ExpEvaluator objects to compute exponentials.

Using the linear source incurs calculating extra exponential terms.

Reimplemented from Solver.

◆ initializeFluxArrays()

void CPULSSolver::initializeFluxArrays ( )
virtual

Allocates memory for boundary and scalar fluxes.

Deletes memory for old flux arrays if they were allocated for a previous simulation. Calls the CPUSolver parent class initializeFluxArrays for track anguar fluxes and flat scalar fluxes.

Reimplemented from CPUSolver.

◆ initializeFSRs()

void CPULSSolver::initializeFSRs ( )
virtual

Initializes the FSR constant linear source component, volumes and Materials array.

This method calls parent class CPUSolver's initalizeFSRs to allocate and initialize an array of OpenMP mutual exclusion locks for each FSR for use in the transport sweep algorithm.

Reimplemented from CPUSolver.

◆ initializeSourceArrays()

void CPULSSolver::initializeSourceArrays ( )
virtual

Allocates memory for FSR source arrays.

Deletes memory for old source arrays if they were allocated for a previous simulation.

Reimplemented from CPUSolver.

◆ normalizeFluxes()

double CPULSSolver::normalizeFluxes ( )
virtual

Normalizes all FSR scalar fluxes and Track boundary angular fluxes to the total fission source (times $ \nu $).

Returns
norm_factor the normalization factor on the scalar fluxes and moments

Reimplemented from CPUSolver.

◆ setFixedSourceMomentByFSR()

void CPULSSolver::setFixedSourceMomentByFSR ( long  fsr_id,
int  group,
double  src_x,
double  src_y,
double  src_z 
)

Assign fixed source moments for a FSR and energy group.

Parameters
fsr_idthe id of the FSR
groupthe energy group (starts at 1)
src_xthe volume-averaged source x-moment in this group
src_ythe volume-averaged source y-moment in this group
src_zthe volume-averaged source z-moment in this group

◆ setFixedSourceMomentsByCell()

void CPULSSolver::setFixedSourceMomentsByCell ( Cell cell,
int  group,
double  src_x,
double  src_y,
double  src_z 
)

Assign fixed source moments for a Cell and energy group.

This routine will add fixed source moments to all instances of the Cell in the geometry (e.g., all FSRs for this Cell).

Parameters
cella pointer to the Cell of interest
groupthe energy group (starts at 1)
src_xthe volume-averaged source x-moment in this group
src_ythe volume-averaged source y-moment in this group
src_zthe volume-averaged source z-moment in this group

◆ tallyLSScalarFlux()

void CPULSSolver::tallyLSScalarFlux ( segment curr_segment,
int  azim_index,
int  polar_index,
FP_PRECISION *  fsr_flux,
FP_PRECISION *  fsr_flux_x,
FP_PRECISION *  fsr_flux_y,
FP_PRECISION *  fsr_flux_z,
float *  track_flux,
FP_PRECISION  direction[3] 
)

Computes the contribution to the LSR scalar flux from a Track segment.

This method integrates the angular flux for a Track segment across energy groups (and polar angles in 2D), and tallies it into the scalar flux buffers, and updates the Track's angular flux.

Parameters
curr_segmenta pointer to the Track segment of interest
azim_indexazimuthal angle index for this 3D Track
polar_indexpolar angle index for this 3D Track
fsr_fluxbuffer to store segment contribution to region scalar flux
fsr_flux_xbuffer to store contribution to the x scalar flux moment
fsr_flux_ybuffer to store contribution to the y scalar flux moment
fsr_flux_zbuffer to store contribution to the z scalar flux moment
track_fluxa pointer to the Track's angular flux
directionthe segment's direction

Member Data Documentation

◆ _fix_src_xyz_cell_map

std::map< std::pair<Cell*, int>, std::vector<double> > CPULSSolver::_fix_src_xyz_cell_map
protected

A map of fixed sources moments keyed by the pair (Cell*, energy group)

◆ _fix_src_xyz_FSR_map

std::map< std::pair<int, int>, std::vector<double> > CPULSSolver::_fix_src_xyz_FSR_map
protected

A map of fixed sources moments keyed by the pair (FSR ID, energy group)

◆ _fixed_source_moments_on

bool CPULSSolver::_fixed_source_moments_on
protected

Whether fixed linear source moments have been provided

◆ _fixed_sources_xyz

std::vector<std::vector<double> > CPULSSolver::_fixed_sources_xyz
protected

Optional user-specified fixed linear source in each FSR & energy group

◆ _FSR_lin_exp_matrix

double* CPULSSolver::_FSR_lin_exp_matrix
protected

The FSR linear expansion matrix values for each FSR

◆ _FSR_source_constants

FP_PRECISION* CPULSSolver::_FSR_source_constants
protected

The FSR source constants for each FSR and energy group

◆ _reduced_sources_xyz

FP_PRECISION* CPULSSolver::_reduced_sources_xyz
protected

An array of the reduced source x, y, and z terms

◆ _scalar_flux_xyz

FP_PRECISION* CPULSSolver::_scalar_flux_xyz
protected

An array of the scalar flux x, y, and z terms

◆ _stabilize_moments

bool CPULSSolver::_stabilize_moments
protected

Whether to stabilize the flux moments

◆ _stabilizing_flux_xyz

FP_PRECISION* CPULSSolver::_stabilizing_flux_xyz
protected

The stabilizing flux for each energy group in each FSR


The documentation for this class was generated from the following files: