.. _section-3.9:

Steady-State and Single-Phase Transient Hydraulics
--------------------------------------------------

.. _section-3.9.1:

Introduction
~~~~~~~~~~~~

The core assembly hydraulics treatment in |SAS| includes the
calculation of coolant flow rates and pressure distributions within each
core channel. Coolant flow rates and pressures are calculated in a
number of different places in the code. They are used in the
steady-state thermal hydraulics initialization; the pre-voiding thermal
hydraulics module calculates coolant flows and pressures; and after the
onset of voiding, they are calculated in the boiling module. This
section describes only the steady-state and pre-voiding calculations.
:numref:`Chapter %s<section-12>` describes the hydraulics calculations after the onset of
voiding. :numref:`Chapter %s<section-14>` and :numref:`Chapter %s<section-16>` describe the hydraulics calculations after
pin disruption.

The coolant flow provides the heat removal from the fuel pins, and so
the main reason for calculating coolant flow rates before the onset of
voiding is to provide information for the fuel-pin temperature
calculations. Core channel flow rates are also needed for the PRIMAR-4
primary loop thermal hydraulics calculations. The pressure distribution
within a channel is calculated mainly to obtain the pressure-dependent
coolant saturation temperature used to determine the onset of boiling.

For the steady-state initialization the user specifies the initial
coolant flow rate for each channel and the outlet plenum pressure. The
code then calculates the pressure distribution in each channel, starting
from the outlet and working down to the inlet. The channel with the
largest steady-state pressure drop is used to determine the steady-state
inlet plenum pressure; and the inlet orifice coefficients in all other
channels are adjusted so that all channels have the same total
steady-state pressure drop.

In the transient hydraulic calculations, both the coolant flow rates and
the pressure distributions are calculated for each channel. Pressure
boundary conditions are used for the transient coolant flow rate
calculations. The driving pressures for these calculations are the inlet
and outlet coolant plenum pressures supplied by the PRIMAR-4 module.
Since all channels use the same inlet and outlet pressures, flow
redistribution between channels is automatically accounted for as
temperatures and flows change. Extrapolated coolant temperatures are
used to evaluate coolant properties in the transient hydraulics
calculations, since coolant flow rates for a time step are calculated
before temperatures are calculated.

:numref:`figure-3.1-2` shows the logic flow in subroutine TSCL0, the driver for
pre-voiding core channel thermal hydraulics. As indicated previously, a
multi-level time step approach is used for the transient calculations in
|SAS|. In general, the routines in a transient module start
with all quantities knows at the beginning of a time step, and a pass
through a module results in calculating the values of relevant
parameters of the module for one channel at the end of the time step.
Different time-step sizes can be used for different phenomena. In
particular, the pre-voiding coolant hydraulics time step can be smaller
than, but not larger than, the heat-transfer time step or the PRIMAR
time step. One pass through TSCL0 calculates one coolant time step for
one channel. If the coolant time step finishes a heat-transfer step,
then the temperature calculations for the heat-transfer step are also
carried out.

:numref:`figure-3.9-1` shows the logic flow in subroutine TSCNV1. This subroutine
and the routines called by it carry out the pre-viding coolant flow and
pressure calculations. Subroutine TSCNV1 also checks on whether to
switch to the boiling model.

As indicated in :numref:`figure-3.9-2`, the pre-voiding transient hydraulics
module interacts with PRIMAR and TSHTRN. PRIMAR supplies the inlet and
outlet pressures that drive the coolant flows in the core channels. In
return, TSCL0 supplies PRIMAR with the current channel flow rates, as
well as hydraulic parameters to use in estimating future flow rates for
the next PRIMAR time step. TSCL0 supplies the coolant flow rates to
TSHTRN, and TSHTRN computes the coolant temperatures used by TSCL0.

.. _figure-3.9-1:

..  figure:: media/Flowchart_3.9-1.*
	:align: center
	:figclass: align-center

	Subroutine TSCNV1, Pre-Boiling Coolant Flow Rates and Pressure Distribution

.. _figure-3.9-2:

..  figure:: media/Flowchart_3.9-2.*
	:align: center
	:figclass: align-center

	Interactions Between Pre-Voiding Transient Hydraulics and Other Modules

.. _section-3.9.2:

Basic Equations
~~~~~~~~~~~~~~~

Before the onset of voiding, the coolant is treated as incompressible,
and the basic equation used for liquid coolant flow in non-voided
channels is

.. math::
   :label: eq-3.9-1

     \frac{1}{A_{\text{c}}} \frac{ \partial \text{w}}{\partial \text{t}} + \frac{ \partial \text{p}}{\partial \text{z}} + \frac{1}{A_{\text{c}}} \frac{\partial \left( \text{wv} \right) }{\partial \text{z}} =  - \left( \frac{\partial \text{p}}{\partial \text{z}} \right)_{\text{fr}} - \left( \frac{\partial \text{p}}{\partial \text{z}} \right)_{\text{K}} - \rho_{\text{c}} g

with

:math:`w` = :math:`w(t)` - independent of :math:`z`

:math:`w` = coolant flow rate (kg/s) :math:`= \rho_{\text{c}} v A_{\text{c}}`

:math:`A_{\text{c}}` = coolant flow area

:math:`p` = pressure

:math:`\rho_{\text{c}}` = density

:math:`v` = coolant velocity (m/s)

:math:`z` = axial position

:math:`t` = time

.. math::
   :label: eq-3.9-2a

	\left( \frac{\partial p}{\partial z} \right)_{\text{fr}} = \text{friction pressure drop}

.. math::
   :label: eq-3.9-2

    \left( \frac{\partial \text{p}}{\partial \text{z}} \right)_{\text{fr}} = f\frac{w \left| w \right|}{2\rho_{\text{c}}A_{\text{c}}^{2}D_{\text{h}}}


:math:`f` = friction factor, which is approximated differently depending on the chosen model option.
When :sasinp:`IFRFAC` = 0:

.. math::
   :label: eq-3.9-3a

	f = \begin{cases}
	A_{\text{fr}}\left( \mathrm{\text{Re}} \right)^{b_{\text{fr}}} & \text{if Re} \geq \text{Re}_{\text{L}} \\
	\frac{A_{\text{fL}}}{\mathrm{\text{Re}}} & \text{if Re} < \text{Re}_{\text{L}} \\
	\end{cases}

When :sasinp:`IFRFAC` = 1:

.. math::
   :label: eq-3.9-3b

	f = A_{\text{fr}}\left( \text{Re} \right)^{b_{\text{fr}}}\  + \frac{A_{\text{fL}}}{\text{Re}}

When :sasinp:`IFRFAC` = 2:

.. math::
   :label: eq-3.9-3c

   	f = \begin{cases}
	A_{\text{fr}}\left( \mathrm{\text{Re}} \right)^{b_{\text{fr}}} & \text{if Re} > \text{Re}_{\text{TR}} \\
	\frac{A_{\text{fL}}}{\mathrm{\text{Re}}}(1-\psi)^{C_{\text{1}}}(1-\psi^{C_{\text{2}}}) + A_{\text{fr}}\left( \mathrm{\text{Re}} \right)^{b_{\text{fr}}}\psi^{C_{\text{3}}} & \text{if Re}_{\text{L}} \leq \text{Re} \leq {\text{Re}_\text{TR}}\\
	\frac{A_{\text{fL}}}{\mathrm{\text{Re}}} & \text{if Re} < \text{Re}_{\text{L}} \\
	\end{cases}

.. math::
   :label: eq-3.9-3d

	\psi = \frac{\text{log}_\text{10}(\frac{\text{Re}}{\text{Re}_\text{L}})}{\text{log}_\text{10}(\frac{\text{Re}_\text{TR}}{\text{Re}_\text{L}})}

where :math:`A_{\text{fr}}` (:sasinp:`AFR`), :math:`b_{\text{fr}}` (:sasinp:`BFR`), :math:`C_{\text{1}}` (:sasinp:`C1TRAN`),
:math:`C_{\text{2}}` (:sasinp:`C2TRAN`), :math:`C_{\text{3}}` (:sasinp:`C3TRAN`), and :math:`A_{\text{fL}}` (:sasinp:`AFLAM`) are
user-supplied correlation coefficients and the correlations are in the form given by [3-16].

:math:`\text{Re} = \frac{D_{\text{h}} \left| w \right|}{\mu A_{\text{c}}}` = Reynolds number

:math:`Re_{\text{L}}` = Maximum Reynolds number for flow in the laminar regime (:sasinp:`RELAM`)

:math:`Re_{\text{TR}}` = Minimum Reynolds number for flow in the turbulent regime (:sasinp:`RETRAN`)

:math:`\mu` = viscosity

:math:`D_{\text{h}}` = hydraulic diameter

:math:`\left( \frac{\partial \text{p}}{\partial \text{z}} \right)_{\text{K}}` = orifice
pressure drop

:math:`g` = acceleration of gravity

:math:`\rho` = density

The axial node structure shown in :numref:`figure-3.2-3` is used for the coolant.
:numref:`figure-3.9-3` shows which variables are defined at node boundaries and
which are averages or integrals over the length of a node. With this
mesh structure, integrating :eq:`eq-3.9-1` over the length of the channel
gives

.. math::
   :label: eq-3.9-5

	I_{1}\frac{\partial \text{w}}{\partial \text{t}} + p_{\text{t}} - p_{\text{b}} + w^{2}I_{2} + A'_{\text{fL}}wI_{\text{3L}} + A'_{\text{fr}}w \left| w \right|^{1 + b_{\text{fr}}}I_{\text{3T}} + w\left| w \right| I_{4}(w) + g I_{5} = 0

where

.. math::
   :label: eq-3.9-6

    I_{1} = \int{\frac{\text{dz}}{A_{\text{c}}} = \left( \frac{\Delta z_{\text{i}}}{A} \right)_{\text{b}}{+ \left( \frac{\Delta z_{\text{i}}}{A} \right)}_{\text{t}} + \sum_{\text{JC} = 1}^{\text{MZC} - 1}{X_{\text{I}1}\left( \text{JC} \right)}}

and

.. math::
   :label: eq-3.9-7

    X_{\text{I}1}\left( \text{JC} \right) = \frac{ \Delta z \left( \text{JC} \right)}{A_{\text{c}}\left( \text{JC} \right)}

:math:`\left( \frac{\Delta z_{\text{i}}}{A} \right)_{\text{b}}` and :math:`\left( \frac{\Delta z_{\text{i}}}{A} \right)_{\text{t}}` are
effective inertial terms at the bottom and top of the subassembly.

.. math::
   :label: eq-3.9-8

    I_{2} = \sum_{\text{JC} = 1}^{\text{MZC} - 1}{X_{\text{I}2}\left( \text{JC} \right)}

.. math::
   :label: eq-3.9-9

    X_{\text{I}2}\left( \text{JC} \right) = \frac{1}{A_{\text{c}}\left( \text{JC} \right)^{2}}\left\lbrack \frac{1}{\rho_{\text{c}}\left( \text{JC} + 1 \right)} - \frac{1}{\rho_{\text{c}}\left( \text{JC} \right)} \right\rbrack

.. math::
   :label: eq-3.9-10a

    I_{\text{3L}} = {\int{\frac{\mu}{2\rho_{\text{c}}A_{\text{c}}D_{\text{h}}^{2}}}}\text{dz} = \sum_{\text{JC} = 1}^{\text{MZC} - 1}{X_{\text{I3L}}\left( \text{JC} \right)}

.. math::
   :label: eq-3.9-11a

    X_{\text{I3L}}\left( \text{JC} \right) = \frac{ \Delta z  (\text{JC} ) \overline{\mu}\left( \text{JC} \right)}{\left\lbrack \rho_{\text{c}}\left( \text{JC} \right) + \rho_{\text{c}}\left( \text{JC} + 1 \right) \right\rbrack A_{\text{c}}\left( \text{JC} \right)D_{\text{h}}\left( \text{JC} \right)^{2}}

.. math::
   :label: eq-3.9-10

    I_{\text{3T}} = {\int{\frac{1}{2\rho_{\text{c}}A_{\text{c}}^{2}D_{\text{h}}}\left\lbrack \frac{D_{\text{h}}}{\mu A_{\text{c}}} \right\rbrack}}^{b_{\text{fr}}}\text{dz} = \sum_{\text{JC} = 1}^{\text{MZC} - 1}{X_{\text{I3T}}\left( \text{JC} \right)}

.. math::
   :label: eq-3.9-11

    X_{\text{I3T}}\left( \text{JC} \right) = \frac{ \Delta z (\text{JC})}{\left\lbrack \rho_{\text{c}}\left( \text{JC} \right) + \rho_{\text{c}}\left( \text{JC} + 1 \right) \right\rbrack A_{\text{c}}\left( \text{JC} \right)^{2}D_{\text{h}}\left( \text{JC} \right)}\ \left\lbrack \frac{D_{\text{h}}\left( \text{JC} \right)}{\overline{\mu}\left( \text{JC} \right)A_{\text{c}}\left( \text{JC} \right)} \right\rbrack^{b_{\text{fr}}}

.. math::
   :label: eq-3.9-12

    I_{4}(w) = \sum_{\text{JC} = 1}^{\text{MZC} - 1}{X_{\text{I}4}\left( \text{JC} \right)} + \sum_{\text{i} = 1}^{\text{NZONE} + 1}{\frac{1}{2\rho_{\text{c}}\left( JC_i \right) A_i^2} K_{i}\left( Re_i \right) }

.. math::
   :label: eq-3.9-12b

    X_{\text{I}4}\left( \text{JC} \right) = \frac{1}{\left\lbrack \rho_{\text{c}}\left( \text{JC} \right) + \rho_{\text{c}}\left( \text{JC} + 1 \right) \right\rbrack A_{\text{c}}\left( \text{JC} \right)^{2}}K_{\text{or}}\left( \text{JC} \right)

.. math::
   :label: eq-3.9-13

    I_{5} = \int{\rho_{\text{c}}} \text{dz} = \sum_{\text{JC} = 1}^{\text{MZC} - 1}{X_{\text{I}5}\left( \text{JC} \right)}

.. math::
   :label: eq-3.9-14

    X_{\text{I}5}\left( \text{JC} \right) = .5\ \left\lbrack \rho_{\text{c}}\left( \text{JC} \right) + \rho_{\text{c}}\left( \text{JC} + 1 \right) \right\rbrack \Delta z \left( \text{JC} \right)

:math:`p_{\text{b}}` = pressure at bottom of channel, and :math:`p_{\text{t}}` = pressure
at top of channel. :math:`A'_{fL}` and :math:`A'_{fr}` are defined as:

.. math::
   :label: AflAfr

    A'_{\text{fL}} = 0 \text{ and } A'_{\text{fr}} = A_{\text{fr}} & \quad \text{if Re} > \text{Re}_{\text{TR}} \\
    A'_{\text{fL}} = A_{\text{fL}}(1-\psi)^{C_{\text{1}}}(1-\psi^{C_{\text{2}}}) \text{ and } A'_{\text{fr}} = A_{\text{fr}}\psi^{C_{\text{3}}} & \quad \text{if Re}_{\text{L}} \leq \text{Re} \leq {\text{Re}_\text{TR}} \\
    A'_{\text{fL}} = A_{\text{fL}} \text{ and } A'_{\text{fr}} = 0 & \quad \text{if Re} < \text{Re}_{\text{L}} \\

.. _figure-3.9-3:

..  figure:: media/image15.png
	:align: center
	:figclass: align-center

	Definition of Variables on the Axial Coolant Mesh

.. _section-3.9.3:

Flow Orifices
~~~~~~~~~~~~~

As mentioned in the section on the |SAS| channel treatment in
:numref:`section-3.2`, the channel is divided axially into a number of zones. One
zone represents the fuel-pin region. Other zones represent regions above
and below the fuel pins. At the bottom of each zone, flow direction-dependent
orifice coefficients can be specified by the user. Additionally,
orifice coefficients can be used at the top of the upper zone.
These axial locations are hereafter called zone interfaces.
The orifice coefficients at zone interfaces can represent orifice blocks,
subassembly inlet and exit losses, and the pin support grid at the bottom of the pins.
Users may define two different types of orifice coefficients at zone interfaces,
a constant average coefficient: :math:`K_{\text{or}} \left( \text{JC} \right)`
and a Re-dependent interface coefficient, :math:`K_{i}\left( Re_i \right)`.
For each type of orifice, the user supplies one orifice coefficient for positive flow and
another for negative flow (anisotropic values).
The anisotropic constant average coefficients are defined using :sasinp:`XKORI`.
The anisotropic Re-dependent coefficients are defined using :sasinp:`ZoneLossCoefTableID`
and a Zone Interface Loss Coefficient Tables, described in :numref:`section-LossCoeff-C`.

In addition, in the fuel-pin region, constant average orifice coefficients representing equally-spaced grid spacers can be specified.
These orifice coefficients are used as the values of :math:`K_{\text{or}} \left( \text{JC} \right)` in :eq:`eq-3.9-12b` for the appropriate axial nodes.
Grid spacer loss coefficients are defined using :sasinp:`NGRDSP` and :sasinp:`XKORGD`.
Spacer grid loss coefficients are assumed to be the same in the forward and reverse directions.

The effect of each constant average orifice is distributed evenly over an axial node,
unlike the Re-dependent interface coefficients which are concentrated at a point.

.. _section-LossCoeff-C:

Anisotropic Re-dependent Loss Coefficients
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

Anisotropic Re-dependent loss coefficients can be defined at zone interfaces using Interface Loss Coefficient Tables.
Zone interfaces are defined as the inlet of a channel axial zone and the outlet of a channel.
By default, all zone interfaces are assumed to have a loss coefficient of zero.
Only zone interfaces identified in an Interface Loss Coefficient Table will overwrite the default interface loss coefficient.
The pressure drop at zone interface :math:`i` caused by the loss coefficient is defined as

.. math::
   :label: int-Pdrop

	P_{i+} - P_{i-} = K_i(Re_i) \frac{w|w|}{2\rho_i A_i^2}

where

.. math::
   :label: int-Re

	Re_i = \frac{Dh_i |w|}{A_i \mu_i}

:math:`{variable}_i` is the variable at interface :math:`i`, and :math:`P_{i+}`/:math:`P_{i-}` are the pressure downstream and upstream the interface, respectively.
In the core channels, the pressure reported at a channel zone interface is the upstream pressure, :math:`P(JC_i) = P_{i-}`.
At the channel outlet, the reported pressure is the downstream pressure, :math:`P(MZC) = P_{(NZONE+1)+}`.

Two functional forms are currently available to capture the dependence of the loss coefficient on the Reynolds number:

.. math::
   :label: Re-loss-form

	K_i(Re_i) = \begin{cases}
	C_1 + C_2 Re_i^{C_3} & \text{Kind} = 1 \\
	C_1 + C_2 exp\left(Re_i{C_3}\right) & \text{Kind} = 2 \\
	\end{cases}

.. math::
   :label: Re-loss-input

    \begin{split}
        C_{1} & \geq 0 \\
        C_{2} & \geq 0 \\
        C_{3} & \ge -2
    \end{split} & \text{Kind} = 1 \\

    \begin{split}
           C_{1} & \geq 0 \\
           C_{1} + C_{2} & \geq 0
    \end{split} & \text{Kind} = 2

Two sets of coefficients may be defined for a given zone interface; one set for forward flow and one set for reverse flow.
Additionally, an interface area, :math:`A_i`, and hydraulic diameter, :math:`Dh_i`, may be defined for the evaluation of the interface pressure drop and interface Re.
If an interface area and interface hydraulic diameter are not provided or are defined as zero, the area and hydraulic diameter of the node upstream or downstream the interface are used.
If the upstream node has a smaller or equal area than the downstream zone, the upstream node's area and hydraulic diameter are used for the interface.
If the upstream node has a bigger area than the downstream zone, the downstream node's area and hydraulic diameter are used for the interface.

At very low flow rates, the derivative of the loss coefficients may approach infinity.
This introduces numerical difficulty when SAS attempts to reverse flow through a channel.
In order to improve the robustness of the semi-implicit solution scheme, an optional lower limit has been introduced for Re-dependent loss coefficients, :math:`ReL`.
When the Reynold number drops below :math:`ReL`, :eq:`Re-loss-form` is evaluated at :math:`ReL`.
Additionally, the derivative of :math:`K_i(Re_i)` is assumed to be zero.

.. list-table::  Input in a Core Channel Interface Loss Coefficient Table
    :header-rows: 1
    :align: center
    :widths: 1,4,2

    * - Column Label
      - Description
      - Range
    * - `iZone`
      - Index of the channel zone where the coefficient will be defined at the inlet.
        `iZone` = 1 is the inlet of the channel. `iZone` = `NZONE` + 1 is the outlet of the channel.
      - 1 :math:`\leq` `iZone` :math:`\leq` `NZONE` + 1
    * - `Kind`
      - Kind of functional form for calculating the anisotropic loss coefficient. See :eq:`Re-loss-form`
      - 1 or 2
    * - `C1f`, `C2f`, `C3f`
      - Function coefficients for forward flow
      - See :eq:`Re-loss-input`
    * - `C1b`, `C2b`, `C3b`
      - Function coefficients for backward flow
      - See :eq:`Re-loss-input`
    * - `A` (optional)
      - Area used to evaluate the Reynolds number and minor pressure loss at the
        interface. Optional area and hydraulic diameter must be provided as a pair.
      - :math:`A \geq 0`
    * - `Dh` (optional)
      - Hydraulic Diameter used to evaluate the Reynolds number at the
        interface. Optional area and hydraulic diameter must be provided as a pair.
      - :math:`Dh \geq 0`
    * - `ReL` (optional)
      - Lower Reynolds number limit for evaluating the functional form.
        For Reynolds numbers below the lower limit, the functional form is evaluated at `ReL` and the
        derivative is assumed to be 0.0.
      - :math:`ReL \geq 0`


.. _ZoneLossInput:

Example Channel Anisotropic Re-dependent Loss coefficient table input::

	TABLE <id> CHAN1
	iZone Kind  C1f    C2f       C3f   C1b    C2b     C3b   A            Dh   ReL
	1     1     0.5 100000        -1   1.5      2       0  2.2          6.5   1.0
	3     2     1       10 -0.000001     2 0.0001 0.00001 0.021        0.08   0.0
	4     1     1   100000        -1   2.2  20000      -1 0.002 0.003183099  10.0
	END


.. _section-3.9.4:

Finite Difference Equations - Coolant Flow Rates
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

If
:math:`w_{1} = w\left( t \right),\ w_{2} = w\left( t +  \Delta t  \right),\  \Delta w  = w_{2} - w_{1}`,
then we approximate :math:`\Delta w` as

.. math::
   :label: eq-3.9-15

    \Delta w = \left\lbrack \theta_{1}\left. \frac{\partial \text{w}}{\partial \text{t}} \right|_{t} + \theta_{2}\left. \frac{\partial \text{w}}{\partial \text{t}} \right|_{t +  \Delta t }  \right\rbrack \Delta t

and

.. math::
   :label: eq-3.9-16

    \theta_{1} + \theta_{2} = 1.0

in which :math:`\theta_2` is the degree of implicitness. For small time
steps a semi-implicit calculation is used, with :math:`\theta_1` =
:math:`\theta_2 = .5~`. For a fully implicit calculation for large time
steps, :math:`\theta_1 = 0` and :math:`\theta_2 = 1.0~`. The value used for
:math:`\theta_2` is

.. math::
   :label: eq-3.9-17

    \theta_{2} = \frac{6.12992 +  2.66054x +  x^{2}}{2\ x\ 6.12992 +  3.56284x +  x^{2}}

where

.. math::
   :label: eq-3.9-18

    x = \Delta t/\tau

and :math:`\tau` is the time constant for flow rate changes given by

.. math::
   :label: eq-3.9-19

    \tau = I_{1}/d_{1}

The definition of :math:`d_1` is given by :eq:`eq-3.9-25` below. This
correlation for :math:`\theta_2` is discussed in :numref:`section-A3.1`. The
additional approximations are made that

.. math::
   :label: eq-3.9-20

    \left( w_{1} + \Delta w \right)^{2} \simeq w_{1}^{2} + 2\Delta w\ w_{1}

.. math::
   :label: eq-3.9-21

    \left( w_{1} + \Delta w \right) \left| w_{1} + \Delta w \right| \simeq w_{1} \left| w_{1} \right| + 2\left| w_{1} \right| \Delta w

.. math::
   :label: eq-3.9-22

    \left( w_{1} + \Delta w \right)\left| w_{1} + \Delta w \right|^{1 + b_{\text{fr}}} \simeq w_{1}\left| w_{1} \right|^{1 + b_{\text{fr}}} + \left( 2 + b_{\text{fr}} \right)\left| w_{1} \right|^{1 + b_{\text{fr}}}\ \Delta w

.. math::
   :label: I4-approx

	I_{4}(w_1 + \Delta w) \simeq I_{4}(w_1) + \frac{\partial I_{4}}{\partial Re}\frac{\partial Re}{\partial w}\Delta w

These equations are combined to give

.. math::
   :label: eq-3.9-23
   :nowrap:

    \begin{multline}
    I_{1}\frac{ \Delta w }{ \Delta t } + \theta_{1}\left( p_{\text{t}1} - p_{\text{b}1} \right) + \theta_{2}\left( p_{\text{t}2} - p_{\text{b}2} \right) + I_{2}\left( w_{1}^{2} + 2\theta_{2}w_{1} \Delta w  \right) \\
     + A'_{\text{fL}}I_{\text{3L}}( w_{1} + \theta_{2} \Delta w ) + A'_{\text{fr}}I_{\text{3T}}\left\lbrack w_{1}\left| w_{1} \right|^{1 + b_{\text{fr}}} + \theta_{2}\left( 2 + b_{\text{fr}} \right) \left| w_{1} \right|^{1 + b_{\text{fr}}} \Delta w  \right\rbrack \\
    + \left( w_1 |w_1| + \theta_{2} 2 |w_1| \Delta w \right) I_{4}(w_1) + \theta_{2} w_1 |w_1| \frac{\partial I_{4} }{\partial Re} \frac{\partial Re}{\partial w} \Delta w + gI_{5} = 0
    \end{multline}

where :math:`p_{\text{t}1}` and :math:`p_{\text{t}2}` are the pressures at the top of
the subassembly at the beginning and end of the time step. Similarly,
:math:`p_{\text{b}1}` and :math:`p_{\text{b}2}` are the pressures at the bottom of the
subassembly at the beginning and end of the time step, respectively.

Solving for the coolant flow rate change, :math:`\Delta w`, gives

.. math::
   :label: eq-3.9-24

    \Delta w = \Delta t \frac{\left\lbrack \theta_{1}\left( p_{\text{b}1} - p_{\text{t}1} \right) + \theta_{2}\left( p_{\text{b}2} - p_{\text{t}2} \right) - I_{2}w_{1}^{2} - A'_{\text{fL}}I_{\text{3L}}w_{1} - A'_{\text{fr}}I_{\text{3T}}w_{1}\left| w_{1} \right|^{1 + b_{\text{fr}}}  - I_{4}(w_1) w_{1} \left| w_{1} \right| - g I_{5} \right\rbrack}{\left( I_{1} + \Delta t\ \theta_{2}d_{1} \right)}

where

.. math::
   :label: eq-3.9-25

    d_{1} = 2 w_{1} I_{2} + A'_{\text{fL}}I_{\text{3L}} + \left( 2 + b_{\text{fr}} \right){ \left| w_{1} \right|}^{1 + b_{\text{fr}}} A'_{\text{fr} }I_{\text{3T}} + 2 |w_1| I_{4}(w_1) + w_1 |w_1| \frac{\partial I_{4}}{\partial Re}\frac{\partial Re}{\partial w}

The temperatures used to evaluate :math:`\rho_{\text{c}}` and
:math:`\overline{\mu}` in  :eq:`eq-3.9-9`, :eq:`eq-3.9-11a`, :eq:`eq-3.9-11`, :eq:`eq-3.9-12b` and :eq:`eq-3.9-14` are
extrapolated coolant temperatures, extrapolated to the end of the
coolant step.

.. _section-3.9.5:

Coolant Pressures
~~~~~~~~~~~~~~~~~

After the flow rate has been calculated, the coolant pressures in a
channel at the end of a time step are calculated. First, the pressure at
node :math:`\text{MZC}`, the last node at the top of the subassembly, is calculated as

.. math::
   :label: eq-3.9-26

    p\left( \text{MZC} \right) = p_{\text{t}2} + \left\lbrack \frac{\Delta z_{\text{i}}}{A} \right\rbrack_{\text{t}}\frac{\partial \text{w}}{\partial \text{t}}

where

.. math::
   :label: eq-3.9-27

    p_{\text{t}2} = p_{\text{x}2} + \rho_{\text{cout}} g \left\lbrack z_{\text{plu}} - z_{\text{c}}\left( \text{MZC} \right) \right\rbrack

:math:`z_{\text{plu}}` = outlet plenum elevation

:math:`\rho_{\text{cout}}` = density at outlet

:math:`p_{\text{x}2}` = coolant outlet plenum pressure at the end of the time
step

Note that :math:`p_{\text{x}2}`, the outlet plenum pressure supplied by the
primary loop module, is defined at an elevation :math:`z_{\text{plu}}`, the
elevation of the upper plenum, whereas
:math:`p\left( \text{MZC} \right)` is calculated at
:math:`z\left( \text{MZC} \right)`, so the gravity head
term occurs in :eq:`eq-3.9-27`.

After :math:`p\left( \text{MZC} \right)` has been calculated,
the other pressures are calculated, starting at node :math:`\text{MZC}`-1 working down,
using

.. math::
    :label: eq-3.9-28
    :nowrap:

    \begin{multline}
    p \left( \text{JC} \right) = p \left( \text{JC} + 1 \right) + X_{\text{I}1}\left( \text{JC} \right)\frac{\partial \text{w}}{\partial \text{t}} + w^{2} X_{\text{I}2}\left( \text{JC} \right) + A'_{\text{fL}}wX_{\text{I3L}}\left( \text{JC} \right) \\
    + A'_{\text{fr}} w \left| w \right|^{1 + b_{\text{fr}}}X_{\text{I3T}}\left( \text{JC} \right) + \left| w \right| w X_{I4} \left( \text{JC} \right) + g\ X_{I5} \left( \text{JC} \right)
    \end{multline}

and

.. math::
   :label: eq-3.9-28-b

    p \left( \text{JC} \right) = \begin{cases}
    p \left( \text{JC} \right) + \left| w \right| w K_{i}/\left(2\rho_{\text{c}} A_i^2 \right)  & \text{JC} = \text{JC}_i \\
    p \left( \text{JC} \right)  & \text{JC} \ne \text{JC}_i \\
    \end{cases}

Note :math:`\text{JC}_i = \text{MZC}-1` for the channel outlet, :math:`i=\text{NZONE}+1`.

In :eq:`eq-3.9-26` and :eq:`eq-3.9-28`, the value used for :math:`\frac{\partial w}{\partial t}` is calculated at
the end of the time step for the transient calculation as

.. math::
   :label: eq-3.9-29

    \frac{\partial \text{w}}{\partial \text{t}} = \frac{\left\lbrack p_{\text{b}2} - p_{\text{t}2} - w_{2}^{2} I_{2} - A'_{\text{fL}}w_{2}I_{\text{3L}} - A'_{\text{fr}}w_{2}\left| w_{2} \right|^{1 + b_{\text{fr}}} I_{\text{3T}} - w_{2}\left| w_{2} \right|I_{4}(w_2) - g I_{5} \right\rbrack}{I_{1}}

The pressure, :math:`p_{\text{b}2}` at the bottom of the subassembly
at the end of the time step is calculated as

.. math::
   :label: eq-3.9-30

    p_{\text{b}2} = p_{\text{in}2} - \rho_{\text{cin}} g \left\lbrack z_{\text{c}}\left( 1 \right) - z_{\text{pll}} \right\rbrack

:math:`z_{\text{pll}}` = plenum location. Note that
:math:`p_{\text{in}2}`, the inlet plenum pressure supplied by the
primary loop module, is defined at an elevation
:math:`z_{\text{pll}}`, the elevation of the lower plenum.

The pressure distribution in the steady state is calculated the same as
in the transient, except that :math:`\frac{\partial w}{\partial t}` is zero in the steady-state
calculation.

At the end of the steady state calculation, SAS will adjust the values of :sasinp:`XKORI` such that the pressure drop
across each channel is consistent, given the user provided flow rates.
These adjustments are reported in the standard output file.
The magnitude of these adjustments can be limited using :sasinp:`ADJKTOL`.
If :sasinp:`ADJKTOL` is set to a positive value and SAS determines an adjustment to :sasinp:`XKORI` that is larger than the
defined tolerance is required, SAS will report an error and the simulation will be terminated.