An open source method of characteristics neutron transport code.
|
This a subclass of the Solver class for multi-core CPUs using OpenMP multi-threading. More...
#include "src/CPUSolver.h"
Public Member Functions | |
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... | |
Geometry * | getGeometry () |
Returns a pointer to the Geometry. More... | |
TrackGenerator * | getTrackGenerator () |
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 . 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 Member Functions | |
virtual void | initializeFluxArrays () |
Allocates memory for Track boundary angular flux and leakage and FSR scalar flux arrays. More... | |
virtual void | initializeSourceArrays () |
Allocates memory for FSR source arrays. More... | |
virtual void | initializeFSRs () |
Initializes the FSR volumes and Materials array. More... | |
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. | |
virtual void | flattenFSRFluxes (FP_PRECISION value) |
Set the scalar flux for each FSR and energy group to some value. More... | |
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. | |
virtual double | normalizeFluxes () |
Normalizes all FSR scalar fluxes and Track boundary angular fluxes to the total fission source (times ). | |
void | computeFSRFissionSources () |
Computes the total fission source in each FSR. More... | |
void | computeFSRScatterSources () |
Computes the total scattering source in each FSR. More... | |
virtual void | computeFSRSources (int iteration) |
Computes the total source (fission, scattering, fixed) 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... | |
virtual void | computeStabilizingFlux () |
Computes the stabilizing flux for transport stabilization. | |
virtual void | stabilizeFlux () |
Adjusts the scalar flux for transport stabilization. | |
virtual void | addSourceToScalarFlux () |
Add the source term contribution in the transport equation to the FSR scalar flux. | |
void | computeKeff () |
Compute from successive fission sources. | |
double | computeResidual (residualType res_type) |
Computes the residual between source/flux iterations. More... | |
Protected Member Functions inherited from Solver | |
virtual void | initializeExpEvaluators () |
Initializes new ExpEvaluator object to compute exponentials. | |
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... | |
virtual void | initializeCmfd () |
Initializes a Cmfd object for acceleration prior to source iteration. 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... | |
This a subclass of the Solver class for multi-core CPUs using OpenMP multi-threading.
CPUSolver::CPUSolver | ( | TrackGenerator * | track_generator = NULL | ) |
Constructor initializes array pointers for Tracks and Materials.
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.
track_generator | an optional pointer to the TrackGenerator |
void CPUSolver::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.
fsr_id | the id of the fsr |
weight | the quadrature weight (only for 3D ray tracing) |
fsr_flux | the buffer containing the segment(s)' contribution |
|
virtual |
Computes the volume-averaged, energy-integrated fission or nu-fission rate in each FSR and stores them in an array indexed by FSR ID.
This is a helper method for SWIG to allow users to retrieve FSR fission or nu-fission rates as a NumPy array. An example of how this method can be called from Python is as follows:
fission_rates | an array to store the fission rates (implicitly passed in as a NumPy array from Python) |
num_FSRs | the number of FSRs passed in from Python |
nu | whether return nu-fission (true) or fission (default, false) rates |
Implements Solver.
|
protectedvirtual |
Computes the total fission source in each FSR.
This method is a helper routine for the openmoc.krylov submodule.
Implements Solver.
|
protectedvirtual |
Computes the total scattering source in each FSR.
This method is a helper routine for the openmoc.krylov submodule.
Implements Solver.
|
protectedvirtual |
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.
Implements Solver.
Reimplemented in CPULSSolver, and VectorizedSolver.
|
protectedvirtual |
Computes the residual between source/flux iterations.
res_type | the type of residuals to compute (SCALAR_FLUX, FISSION_SOURCE, TOTAL_SOURCE) |
Implements Solver.
|
protectedvirtual |
Set the scalar flux for each FSR and energy group to some value.
value | the value to assign to each FSR scalar flux |
Implements Solver.
Reimplemented in CPULSSolver.
|
virtual |
Fills an array with the scalar fluxes.
This class method is a helper routine called by the OpenMOC Python "openmoc.krylov" module for Krylov subspace methods. Although this method appears to require two arguments, in reality it only requires one due to SWIG and would be called from within Python as follows:
out_fluxes | an array of FSR scalar fluxes in each energy group |
num_fluxes | the total number of FSR flux values |
Implements Solver.
int CPUSolver::getNumThreads | ( | ) |
Returns the number of shared memory OpenMP threads in use.
|
protectedvirtual |
Allocates memory for Track boundary angular flux and leakage and FSR scalar flux arrays.
Deletes memory for old flux arrays if they were allocated for a previous simulation.
Implements Solver.
Reimplemented in CPULSSolver, and VectorizedSolver.
|
protectedvirtual |
Initializes the FSR volumes and Materials array.
This method allocates and initializes an array of OpenMP mutual exclusion locks for each FSR for use in the transport sweep algorithm.
Reimplemented from Solver.
Reimplemented in CPULSSolver.
|
protectedvirtual |
Allocates memory for FSR source arrays.
Deletes memory for old source arrays if they were allocated for a previous simulation.
Implements Solver.
Reimplemented in CPULSSolver, and VectorizedSolver.
void CPUSolver::printFluxesTemp | ( | ) |
A function that prints fsr fluxes in xy plane at z=middle.
Recommend deletion, since redundant printFluxes
void CPUSolver::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.
dim1 | coordinates of the mesh grid in the first direction |
dim2 | coordinates of the mesh grid in the second direction |
offset | The location of the mesh grid center on the perpendicular axis |
plane | 'xy', 'xz' or 'yz' the plane in which the mesh grid lies |
void CPUSolver::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.
iteration | the current iteration |
num_x | number of divisions in X direction |
num_y | number of divisions in Y direction |
num_z | number of divisions in Z direction |
void CPUSolver::setFixedSourceByFSR | ( | long | fsr_id, |
int | group, | ||
FP_PRECISION | source | ||
) |
Assign a fixed source for a flat source region and energy group.
Fixed sources should be scaled to reflect the fact that OpenMOC normalizes the scalar flux such that the total energy- and volume-integrated production rate sums to 1.0.
fsr_id | the flat source region ID |
group | the energy group |
source | the volume-averaged source in this group |
|
virtual |
Set the flux array for use in transport sweep source calculations.
This is a helper method for the checkpoint restart capabilities, as well as the IRAMSolver in the openmoc.krylov submodule. This routine may be used as follows from within Python:
NOTE: This routine stores a pointer to the fluxes for the Solver to use during transport sweeps and other calculations. Hence, the flux array pointer is shared between NumPy and the Solver.
in_fluxes | an array with the fluxes to use |
num_fluxes | the number of flux values (# groups x # FSRs) |
Implements Solver.
void CPUSolver::setNumThreads | ( | int | num_threads | ) |
Sets the number of shared memory OpenMP threads to use (>0).
num_threads | the number of threads |
void CPUSolver::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.
void CPUSolver::tallyScalarFlux | ( | segment * | curr_segment, |
int | azim_index, | ||
int | polar_index, | ||
FP_PRECISION * | fsr_flux, | ||
float * | track_flux | ||
) |
void CPUSolver::transferBoundaryFlux | ( | Track * | track, |
int | azim_index, | ||
int | polar_index, | ||
bool | direction, | ||
float * | track_flux | ||
) |
|
protectedvirtual |
This method performs one transport sweep of all azimuthal angles, Tracks, Track segments, polar angles and energy groups.
The method integrates the flux along each Track and updates the boundary fluxes for the corresponding output Track, while updating the scalar flux in each flat source region.
Implements Solver.
|
protected |
OpenMP mutual exclusion locks for atomic FSR scalar flux updates
|
protected |
The number of shared memory OpenMP threads