The predictor-corrector mechanism and solver encapsulation. More...
Functions/Subroutines | |
subroutine | deplete (solver_opts, lib, itot, n0, n, delta_t, source_term, flux, power, input_by_flux, is_adjoint, fixed_fission_energy, reprocessing_rate, feed_rate, adjoint_source_term_reaction, adjoint_source_term_decay) |
Perform a single depletion step (including the predictor/corrector mechanism). More... | |
subroutine | depletion_step (solver_opts, lib, itot, n0, n, delta_t, corrector, powrst, fluxst, source_term, flux, power, input_by_flux, is_adjoint, fixed_fission_energy, reprocessing_rate, feed_rate, adjoint_source_term_reaction, adjoint_source_term_decay) |
Run the solver once. More... | |
subroutine | powerflux (lib, itot, power, flux, concentrations, input_by_flux, fixed_fission_energy, is_adjoint) |
Convert power to flux or vice versa. More... | |
subroutine | fluxo (corrector, powrst, fluxst, flux, power, zero_flux_step, input_by_flux) |
Computes the average flux over a time interval using predictor corrector method, build the diagonal of the transition matrix (DIAG array), convert flux to power or vice versa, multiply A matrix off-diagoal elements by flux. More... | |
subroutine | powerf (lib, itot, factor, fission_rate, concentrations, zero_flux_step) |
Calculate the power-flux conversion factor, based on given fission rate and given concentrations. More... | |
Detailed Description
The predictor-corrector mechanism and solver encapsulation.
This module encapsulates the predictor-corrector mechanism of depletion and provides the parts of the depletion mechanism which are common to all solvers.
The deplete
subroutine executes two solver runs for predictor and corrector steps.
The depletion_step
subroutine calculates the flux (or power) integrated over the time step, builds the transition matrix A for the solver (by folding the reaction and decay data together using the evaluated flux), builds up the source vector B (consisting of the user-specified constant feed term and from the terms arising from the adjoint solution, if applicable) and runs the solver itself.
Function/Subroutine Documentation
subroutine depletion::deplete | ( | type(solver_options), intent(in) | solver_opts, |
type(origen_librarywrapper), intent(in) | lib, | ||
integer, intent(in) | itot, | ||
real(c_double), dimension(itot), intent(in) | n0, | ||
real(c_double), dimension(itot), intent(out) | n, | ||
real(c_double), intent(in) | delta_t, | ||
real(c_double), dimension(itot), intent(inout) | source_term, | ||
real(c_double), intent(inout) | flux, | ||
real(c_double), intent(inout) | power, | ||
logical, intent(in) | input_by_flux, | ||
logical, intent(in) | is_adjoint, | ||
logical, intent(in) | fixed_fission_energy, | ||
real(c_double), dimension(itot), intent(in) | reprocessing_rate, | ||
real(c_double), dimension(itot), intent(in) | feed_rate, | ||
real(c_double), dimension(itot), intent(in) | adjoint_source_term_reaction, | ||
real(c_double), dimension(itot), intent(in) | adjoint_source_term_decay | ||
) |
Perform a single depletion step (including the predictor/corrector mechanism).
- Parameters
-
[in] solver_opts solver settings [in] lib Library data [in] itot Number of nuclides in the library [in] n0 initial concentrations [out] n final concentrations [in,out] source_term Continuous processing feed rate / source term from the adjoint calculation [in] delta_t Time step length in seconds [in] input_by_flux flux is used for input? or power? [in] is_adjoint is this an adjoint calculation? [in] fixed_fission_energy is the fission energy fixed or problem dependent? [in] reprocessing_rate continuous reprocessing rates [in] feed_rate continuous feed rates [in] adjoint_source_term_reaction reaction component of the adjoint source term [in] adjoint_source_term_decay decay component of the adjoint source term
buffers for predictor-corrector power and flux
References depletion_step(), and origen_iofunctions_m::rstop().
Referenced by origen_casewrapper_m::flip_time().
subroutine depletion::depletion_step | ( | type(solver_options), intent(in) | solver_opts, |
type(origen_librarywrapper), intent(in) | lib, | ||
integer, intent(in) | itot, | ||
real(c_double), dimension(itot), intent(in) | n0, | ||
real(c_double), dimension(itot), intent(inout) | n, | ||
real(c_double), intent(in) | delta_t, | ||
logical, intent(in) | corrector, | ||
real(c_double), intent(inout) | powrst, | ||
real(c_double), intent(inout) | fluxst, | ||
real(c_double), dimension(itot), intent(inout) | source_term, | ||
real(c_double), intent(inout) | flux, | ||
real(c_double), intent(inout) | power, | ||
logical, intent(in) | input_by_flux, | ||
logical, intent(in) | is_adjoint, | ||
logical, intent(in) | fixed_fission_energy, | ||
real(c_double), dimension(itot), intent(in) | reprocessing_rate, | ||
real(c_double), dimension(itot), intent(in) | feed_rate, | ||
real(c_double), dimension(itot), intent(in) | adjoint_source_term_reaction, | ||
real(c_double), dimension(itot), intent(in) | adjoint_source_term_decay | ||
) |
Run the solver once.
- Parameters
-
[in] solver_opts solver type identifier [in] lib Library data [in] itot Total number of nuclides [in] n0 initial concentrations [in,out] n final concentrations [in] delta_t Time step length in seconds [in] corrector is this a corrector step? [in,out] powrst saved power for predictor-corrector [in,out] fluxst saved flux for predictor-corrector [in,out] source_term Continuous processing feed rate / source term from the adjoint calculation [in] input_by_flux flux is used for input? or power? [in] is_adjoint is this an adjoint calculation? [in] fixed_fission_energy is the fission energy fixed or problem dependent? [in] reprocessing_rate continuous reprocessing rates [in] feed_rate continuous feed rates [in] adjoint_source_term_reaction reaction component of the adjoint source term [in] adjoint_source_term_decay decay component of the adjoint source term
Diagonal elements of the transition matrix
Nondiagonal elements of the transition matrix
is this a zero-flux (=decay) step?
References origen_definitions_m::err, fluxo(), kernel_matrex::matrex_solve(), powerflux(), kernel_radau5::radau5_solve(), and origen_iofunctions_m::rstop().
Referenced by deplete().
subroutine depletion::powerflux | ( | type(origen_librarywrapper), intent(in) | lib, |
integer, intent(in) | itot, | ||
real(c_double), intent(inout) | power, | ||
real(c_double), intent(inout) | flux, | ||
real(c_double), dimension(itot), intent(in) | concentrations, | ||
logical, intent(in) | input_by_flux, | ||
logical, intent(in) | fixed_fission_energy, | ||
logical, intent(in) | is_adjoint | ||
) |
Convert power to flux or vice versa.
- Parameters
-
[in] lib Library data [in] itot Total number of nuclides [in,out] flux user-given flux (if index is 0; otherwise not used) [in,out] power user-given power (if index is 1; otherwise not used) [in] concentrations Concentrations to be used [in] input_by_flux depletion by flux or by power? [in] fixed_fission_energy is the fission energy fixed or problem-dependent? [in] is_adjoint is this an adjoint calculation?
Local holder for fiss array
References origen_definitions_m::err, powerf(), and origen_iofunctions_m::rstop().
Referenced by depletion_step().
subroutine depletion::fluxo | ( | logical, intent(in) | corrector, |
real(c_double), intent(inout) | powrst, | ||
real(c_double), intent(inout) | fluxst, | ||
real(c_double), intent(inout) | flux, | ||
real(c_double), intent(inout) | power, | ||
logical, intent(in) | zero_flux_step, | ||
logical, intent(in) | input_by_flux | ||
) |
Computes the average flux over a time interval using predictor corrector method, build the diagonal of the transition matrix (DIAG array), convert flux to power or vice versa, multiply A matrix off-diagoal elements by flux.
Given an input value for the specific power, XFLUXO solves for the average flux over the time interval by using Eq. (F7.2.30) and energy per fission data. Given an input value for the neutron flux, XFLUXO solves for the specific power by Eq. (F7.2.32). Having obtained a value for the neutron flux, XFLUXO generates the first-order rate constants for neutron-induced reactions. Then the diagonal matrix element is formed for each nuclide by combining the disintegration constant, the rate constant for neutron-induced reactions, and the rate constant for removal processes that are proportional to the nuclide concentration. The diagonal matrix elements are stored in the D array.
- Parameters
-
[in] corrector false => save flux (fluxst/powrst), true => average saved and actual flux [in,out] powrst saved power from predictor to corrector [in,out] fluxst saved flux from predictor to corrector [in,out] flux user-given flux (if index is 0; otherwise not used) [in,out] power user-given power (if index is 1; otherwise not used) [in] zero_flux_step is this a zero-flux step? [in] input_by_flux depletion by flux or by power?
Referenced by depletion_step().
subroutine depletion::powerf | ( | type(origen_librarywrapper), intent(in) | lib, |
integer, intent(in) | itot, | ||
real(c_double), intent(out) | factor, | ||
real(c_double), intent(in) | fission_rate, | ||
real(c_double), dimension(itot), intent(in) | concentrations, | ||
logical, intent(in) | zero_flux_step | ||
) |
Calculate the power-flux conversion factor, based on given fission rate and given concentrations.
- Parameters
-
[in] lib Library data [in] itot Total number of nuclides [out] factor power-flux factor to be calculated [in] fission_rate normalization factor [in] concentrations used concentrations (usually XTEMP) [in] zero_flux_step is this a zero-flux step?
Cycle variables
Temporary local copies of library data objects
Energy library size
Nuclide IDs in energy library
Fission energy in energy library
Prompt gamma capture energy (mev) fission energy in energy library
Populate local array copies
Garbage collection
References origen_definitions_m::err.
Referenced by powerflux().