Main solver routines. More...
Functions/Subroutines | |
subroutine, public | get_precision () |
find number of significant places in double precision word compute the precision, iplace, of double precision word on this machine More... | |
subroutine, public | matrex_solve (trx, itot, non, n0, n, delta_t, diag, a, b, zero_flux_step, is_adjoint, nterm, nshrt, cutoff) |
subroutine | decay (kd, non0, loc, itot, non, n_final, t, diag, aa, long, xp, n_initial, xend, source_term, zero_flux_step, is_adjoint, nshrt, cut) |
DECAY determines which nuclides have half-lives that are less than about 10% of the time interval and solves for their concentrations using the Bateman solution, Eq. (F7.2.15). The end-of-interval concentrations are stored in the array N_FINAL. XDECAY also augments the concentrations of long-lived daughter products and updates the augmented beginning-of-interval concentrations in the array N_INITIAL. More... | |
subroutine | term (kd, non0, loc, itot, non, concentrations, t, diag, a, long, xp, xtemp, xend, source_term, zero_flux_step, is_adjoint, nterm, nshrt, cut) |
XTERM forms the reduced coefficient matrix such that it satisfies the criteria inherent in Eq. (F7.2.11). This involves the elimination of the transitions involving those nuclides with short half-lives, the computation of matrix elements required in this elimination, and the use of the augmented beginning-of interval concentrations calculated in XDECAY. XTERM goes on to solve the reduced system of equations by the matrix exponential method. This method uses the recursion relation, Eqs. (F7.2.9) or (F7.2.24), to obtain the solution, Eq. (F7.2.23). The solution is contained in the array XNEW. MATREX Performs the matrix exponential method recursion computation, after arrays initialized in XTERM. More... | |
subroutine | equil (kd, non0, loc, itot, non, n_final, diag, a, long, xp, n_initial, source_term, zero_flux_step) |
Get number of terms needed for matrix expansion. More... | |
subroutine | matrex (nlarge, non, itot, t, err, cim0, cimn, long, nonp, locp, ap, d, csum) |
Variables | |
real(c_double), parameter | akdjq_treshold = 1.0d-6 |
real(c_double), parameter | decay_long_treshold = -50.0d0 |
real(c_double), parameter | equal_d_treshold = 1.0d-7 |
real(c_double) | chek4 |
precision-related limit More... | |
real(c_double) | ratio4 |
precision-related limit More... | |
real(c_double) | qxn |
precision-related limit More... | |
real(c_double) | axn |
precision-related limit More... | |
real(c_double), parameter | xlow = tiny(1.0)*2.0 |
value close to lowest, xlow 2.00d-38 More... | |
character(256), dimension(256), private | rbuffer |
buffer for writing stuff More... | |
Detailed Description
Main solver routines.
Function/Subroutine Documentation
subroutine, public kernel_matrex::get_precision | ( | ) |
find number of significant places in double precision word compute the precision, iplace, of double precision word on this machine
References axn, chek4, qxn, ratio4, rbuffer, and origen_iofunctions_m::rstop().
Referenced by matrex_solve().
subroutine, public kernel_matrex::matrex_solve | ( | type(origen_transitionmatrixp), intent(in) | trx, |
integer, intent(in) | itot, | ||
integer, intent(in) | non, | ||
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(in) | diag, | ||
real(c_double), dimension(non), intent(in) | a, | ||
real(c_double), dimension(itot), intent(in) | b, | ||
logical, intent(in) | zero_flux_step, | ||
logical, intent(in) | is_adjoint, | ||
integer, intent(in) | nterm, | ||
integer, intent(in) | nshrt, | ||
real(c_double), intent(in) | cutoff | ||
) |
- Parameters
-
[in] trx Library data [in] itot Total number of nuclides [in] non Number of non-zero transition matrix elements [in] diag Diagonal elements of the transition matrix [in] a Nondiagonal elements of the transition matrix [in] b Right side of the equation (source term) [in] n0 initial concentrations [out] n final concentrations [in] delta_t Time step length in seconds [in] zero_flux_step is this a zero-flux (=decay) step? [in] is_adjoint is this an adjoint calculation? [in] nterm solver options [in] nshrt solver options [in] cutoff Results for nuclides with atomic fraction below this are zeroed.
which decay elements are long-lived
EXP(-diag*t) for long-lived, 0 otherwise
References decay(), equil(), get_precision(), and term().
Referenced by depletion::depletion_step().
|
private |
DECAY determines which nuclides have half-lives that are less than about 10% of the time interval and solves for their concentrations using the Bateman solution, Eq. (F7.2.15). The end-of-interval concentrations are stored in the array N_FINAL. XDECAY also augments the concentrations of long-lived daughter products and updates the augmented beginning-of-interval concentrations in the array N_INITIAL.
- Parameters
-
[in] kd Library data [in] non0 Library data [in] loc Library data [in] itot Total number of nuclides [in] non Number of non-zero transition matrix elements [out] n_final Current n_final [in] t Step length [in] diag Diagonal elements [in] aa Nondiagonal elements [out] long Long-live decay nuclide map [out] xp EXP(-diag*t) for long-lived, 0 otherwise [in,out] n_initial Initial concetrations (will be augmented). [out] xend MRR 00-012 p. 7 eq. 10 - the additional term to be added to matrix exponential solution [in] source_term Source term [in] zero_flux_step is this a zero-flux step? [in] is_adjoint is this an adjoint calculation?
dummy; diag(nchain)*t
the second sum in F7.2.15 accumulator
large sum accumulators
Diag(c)
Diag(b)
Diag(a)
A(b->c) / Diag(a) * que or que
just for user output
nuclide (main loop)
parent of c
parent of b
index of reaction b -> c
index of reaction a -> b
boundaries for t_bc-cycle
boundaries for t_ab-cycle
chain size
chain: diag
chain: xp
Bateman – chain loop variable (outer)
Bateman – chain loop variable (inner)
dDiag(u)
XP(c)
dXP(u)
dc / du relative difference
short-lived chain
pointer to the current position in chain
accumulator for the first product in F7.2.15
References akdjq_treshold, axn, decay_long_treshold, equal_d_treshold, origen_definitions_m::err, and origen_iofunctions_m::rstop().
Referenced by TransitionMatrixP::is_mts_properly_ordered(), matrex_solve(), FakeFactory::newspLibrary_10nuclide(), Origen::tabulateMaterialLoss(), and Origen::tabulateNuclideInfo().
|
private |
XTERM forms the reduced coefficient matrix such that it satisfies the criteria inherent in Eq. (F7.2.11). This involves the elimination of the transitions involving those nuclides with short half-lives, the computation of matrix elements required in this elimination, and the use of the augmented beginning-of interval concentrations calculated in XDECAY. XTERM goes on to solve the reduced system of equations by the matrix exponential method. This method uses the recursion relation, Eqs. (F7.2.9) or (F7.2.24), to obtain the solution, Eq. (F7.2.23). The solution is contained in the array XNEW. MATREX Performs the matrix exponential method recursion computation, after arrays initialized in XTERM.
- Parameters
-
[in] kd Library data [in] non0 Library data [in] loc Library data [in] itot Total number of nuclides [in] non Number of non-zero transition matrix elements [out] concentrations Initial concentrations [in] diag diagonal terms (removal) [in] a Nondiagonal elements [in] long which isotopes are long-lived enough? [in] is_adjoint is this an adjoint calculation?
References akdjq_treshold, equal_d_treshold, origen_definitions_m::err, and matrex().
Referenced by matrex_solve().
|
private |
Get number of terms needed for matrix expansion.
XEQUIL sets short-lived daughters of long-lived precursors into secular equilibrium with their parents. This results in the set of equations, Eq. (F7.2.18), which is solved through the Gauss-Seidel iterative technique. The iterative sequence for the concentrations is given in Eq. (F7.2.19).
- Parameters
-
[in] kd Library data [in] non0 Library data [in] loc Library data [in] itot Total number of nuclides [in] non Number of non-zero transition matrix elements [in] a Nondiagonal elements [in,out] n_final Final concentrations [in] diag Diagonal elements of the transition matrix [in] long Which nuclides are long-lived [in] xp EXP(-diag*t) for long-lived, 0 otherwise [in] source_term Fixed source [in] n_initial Initial concentrations [in] zero_flux_step is this a zero-flux step?
convergence buffers
i - daughter nuclide
j->i transition pointer
j - parent nuclide
bounds for parent iteration
iteration counter
summation term
-diag(i)
-diag(j)
summing over the transition matrix row
vector of guesses of the short-lived concentration productions
N0(i) * XP(i) for long-lived nuclides
References qxn, origen_iofunctions_m::rstop(), and const::xlow.
Referenced by matrex_solve().
|
private |
- Parameters
-
[in] non Number of compressed transition matrix elements [in] itot Total number of nuclides, = library.itot [in] t Length of time step [in] err "small number", = data_C.err [in] nonp Pointer to the last parent of a nuclide for the reduced matrix (eq. to NON0) [in] long Long-lived decay nuclides map
References origen_iofunctions_m::rstop().
Referenced by term().
Variable Documentation
real(c_double), parameter decay_long_treshold = -50.0d0 |
Referenced by decay().
real(c_double) chek4 |
precision-related limit
Referenced by get_precision().
real(c_double) ratio4 |
precision-related limit
Referenced by get_precision().
real(c_double) qxn |
precision-related limit
Referenced by equil(), and get_precision().
real(c_double) axn |
precision-related limit
Referenced by decay(), and get_precision().
real(c_double), parameter xlow = tiny(1.0)*2.0 |
value close to lowest, xlow 2.00d-38
character(256), dimension(256), private rbuffer |
buffer for writing stuff
Referenced by get_precision().