.. _section-13.3:

Method of Solution
------------------

.. _section-13.3.1:

CLAP Initialization and Cladding Meltthrough
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

The temperature and melt fraction of initial cladding are computed after
the onset of boiling by the SAS fuel-pin model located in the subroutine
TSHTRV and the calculation continues after CLAP initiation. This
subroutine also contains the logical statements for CLAP initiation and
cladding melt-through.

The first cladding to move is assumed to be the first axial segment that
is completely molten; CLAP is initiated with this first melt-through.
Any other axial segments join the molten cladding region if any one of
its three radial segment nodes is molten; i.e., is above the liquidus
temperature.

When entering CLAP (which uses coolant time steps) from a heat-transfer
time step, a check is made to see whether any previously solid cladding
nodes are now to be treated as molten cladding. If so, the total mass
and energy of the initial, refrozen, and moving cladding are evaluated
and these become the mass and energy of the moving cladding.

.. _section-13.3.2:

Solution of the Cladding Energy Equations
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

.. _section-13.3.2.1:

Interface with the Fuel-pin Model
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

There are slight differences in the thermodynamic representations of
molten steel in the fuel-pin and CLAP models. In particular, the
fuel-pin model in subroutine TSHTRV has a melting band (discrete solidus
and liquidus temperatures) and a temperature-dependent specific heat,
whereas CLAP has a single melting temperature and constant specific
heats. If cladding temperatures from the fuel-pin model were used
directly in CLAP calculations (for example, in adding segments to the
molten region), then, possibly, energy would not be conserved due to
differences in the thermodynamic models. To insure energy conservation,
CLAP instead uses nodal energies rather than temperatures. A mean energy
is computed from the SAS nodal energies:

(13.3‑1)

.. _eq-13.3-1:

.. math::

	{\overline{e}}_{\text{j}} = \frac{1}{4} e_{1,\text{j}} + \frac{1}{2} e_{2,\text{j}} + \frac{1}{4} e_{3,\text{j}}

where the :math:`e_{\text{i,j}}`'s are the energies of the three radial nodes
in axial segment :math:`j`. CLAP then evaluates the intact-cladding
temperature :math:`T_{\text{i}}` and melt fraction using
:math:`\overline{e}` and Eqs. :ref:`13.2-41<eq-13.2-41>` to :ref:`13.2-43<eq-13.2-43>`.

The fuel-pin model computes fuel temperatures in the molten cladding
zone and both fuel and intact cladding temperatures for the remainder of
the fueled section. For the molten zone, the CLAP routine integrates the
fuel surface heat loss for each segment for a heat-transfer time
interval:

(13.3‑2)

.. _eq-13.3-2:

.. math::

	Q_{\text{j}} = P_{\text{r}} \Delta z_{\text{j}} \int_{\text{t}^{*}}^{t^{*} + \Delta t^{*}}{\phi_{\text{c,j}}\text{dt}}

where

:math:`\Delta z_{\text{j}} = z_{\text{j}+1} - z_{\text{j}}`

:math:`t^{*}` is time at start of current heat-transfer time step

:math:`\Delta t^{*}` is current heat-transfer time interval

:math:`\phi_{\text{c}}` is the heat flux given by Eq. :ref:`13.2-27<eq-13.2-27>`

:math:`z_{\text{j}}` is the elevation at the bottom of the ":math:`j`"th segment.

This heat loss is utilized at the outer fuel boundary condition by the
fuel-pin model. In a similar fashion, the heat loss from the outer
surface of an intact cladding segment is computed by

(13.3‑3)

.. _eq-13.3-3:

.. math::

	Q_{\text{j}} = P_{\text{e}} \Delta z_{\text{j}} \int_{\text{t}^{*}}^{t^{*} + \Delta t^{*}}{\phi_{\text{r,j}} \text{dt}}

where :math:`\phi_{\text{r}}` is the heat flux from Eq. :ref:`13.2-30<eq-13.2-30>` or :ref:`13.2-37<eq-13.2-37a>`,
depending upon the configuration. This heat loss becomes part of the
outer cladding boundary condition in the fuel-pin model, being added to
the coolant heat loss, if any.

.. _section-13.3.2.2:

Intact and Refrozen Cladding
^^^^^^^^^^^^^^^^^^^^^^^^^^^^

For this configuration, we set φ = 0 (adiabatic outer boundary) in the
refrozen cladding energy Eq. :ref:`13.2-23<eq-13.2-23>`. The energy equation is then
combined with Eqs. :ref:`13.2-30<eq-13.2-30>` and :ref:`13.2-41<eq-13.2-41>` to give

(13.3‑4)

.. _eq-13.3-4:

.. math::

	\frac{\text{dT}_{\text{s,j}}}{\text{dt}} = \xi_{1} \left( T_{\text{i}} - T_{\text{s}} \right) \bigg\rvert_{\text{j}}

where

.. math:: \xi_{1,\text{j}} = \frac{P_{\text{e}} k_{\text{c}}}{\left\lbrack A_{\text{s}}\rho_{\text{s}} c_{\text{ps}}\left( \Delta r_{\text{s}} + \Delta r_{\text{i}} \right) \right\rbrack \bigg\rvert_{\text{j}}}

Since Eq. :ref:`13.2-23<eq-13.2-23>` has been reduced from a general field equation to a
nodal equation, it becomes an ordinary differential equation (in time)
in the process. The subscript :math:`j`, of course, denotes the axial segment
number for which the parameter is evaluated.

Equation :ref:`13.3-4<eq-13.3-4>` is converted to an implicit finite difference equation
by substituting for the left side, the following:

(13.3‑5)

.. _eq-13.3-5:

.. math::

	\frac{\text{dT}_{\text{s,j}}}{\text{dt}} \approx \frac{T_{\text{s,j}}^{n}}{d}

and by evaluating :math:`T_{\text{s,j}}` on the right side at the advanced time

(13.3‑6)

.. _eq-13.3-6:

.. math::

	T_{\text{s,j}}^{n + 1} = T_{\text{s,j}}^{n} + \Delta T_{\text{s,j}}^{n}

giving (after some rearrangement):

(13.3‑7)

.. _eq-13.3-7:

.. math::

	T_{\text{s,j}}^{n} = \frac{\xi_{1}\left( T_{\text{i}} - T_{\text{s}} \right)t}{1 + \xi_{1} t} \bigg\rvert_{\text{j}}^{n}

The implicit form is used rather than an explicit solution because it is
stable for large :math:`\xi_{1}` (or small :math:`A_{\text{s}}`) and a reasonably
sized :math:`\Delta t`.

The heat flux :math:`\phi_{\text{r}}`, required for the fuel-pin boundary
condition (Eq. :ref:`13.3-3<eq-13.3-3>`) is evaluated from the energy equation:

(13.3‑8)

.. _eq-13.3-8:

.. math::

	\phi_{\text{r,j}}^{n} = \frac{A_{\text{s}}\rho_{\text{s}}C_{\text{ps}}\Delta T_{\text{s}}}{P_{\text{e}} \Delta t} \bigg\rvert_{\text{j}}^{n}

.. _section-13.3.2.3:

Evaluation of Convective Term in the Moving Cladding Energy Equation
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

The CLAP model uses a modified donor cell technique developed by Bohl
[13-2] for the evaluation of the convective term

(13.3‑9)

.. _eq-13.3-9:

.. math::

	C_{\text{v}} = A_{\text{c}} \rho_{\text{c}}v_{\text{c}}\frac{\partial \text{e}_{\text{c}}}{\partial \text{z}}

in Eq. :ref:`13.2-20<eq-13.2-20>`. This scheme uses a segment boundary mass flow
:math:`w^{*}` which is estimated by

(13.3‑10)

.. _eq-13.3-10:

.. math::

	w_{\text{j}}^{*} = \frac{1}{2}\left\lbrack A_{\text{c,j} - 1}\rho_{\text{c,j} - 1} v_{\text{c,j} - 1} + A_{\text{c,j}}\rho_{\text{c,j}}v_{\text{c,j}} \right\rbrack

The flow direction is determined by a mean flow :math:`w_{\text{m}}`:

(13.3‑11)

.. _eq-13.3-11:

.. math::

	w_{\text{m,j}} = \frac{1}{2}\left( \text{w}_{\text{j}}^{*} + w_{\text{j} + 1}^{*} \right)

The recipe for the convective term for a unblocked segment is as follows

(13.3‑12a)

.. _eq-13.3-12a:

.. math::

	\begin{align}
	C_{\text{v,j}} = w_{\text{j}}^{*} \frac{e_{\text{c,j}} - e_{\text{c,j} - 1}}{z_{\text{m,j}} - z_{\text{m,j} - 1}} && \text{for } w_{\text{m,j}} > 0, w_{\text{j}}^{*} > 0
	\end{align}

(13.3‑12b)

.. _eq-13.3-12b:

.. math::

	\begin{align}
	C_{\text{v,j}} = w_{\text{j} + 1}^{*}\frac{e_{\text{c,j} + 1} - e_{\text{c,j}}}{z_{\text{m,j} + 1} - z_{\text{m,j}}} && \text{for } w_{\text{m,j}} < 0, w_{\text{j} + 1}^{*} < 0
	\end{align}

(13.3‑12c)

.. _eq-13.3-12c:

.. math::

	\begin{align}
	C_{\text{v,j}} = 0, && \text{otherwise}
	\end{align}

The parameter :math:`z_{\text{m,j}}` is the nodal elevation, which is related
to the segment interface elevations :math:`z_{\text{j}}` by

(13.3‑13)

.. _eq-13.3-13:

.. math::

	z_{\text{m,j}} = \frac{1}{2} \left( z_{\text{j}} + z_{\text{j} + 1} \right)

If the ":math:`j`"th segment is nearly blocked, then the energy convection
computed by Eq. :ref:`13.3-12<eq-13.3-12a>` may be too large in magnitude. In this case, we
define a segment flow :math:`w` as

(13.3‑14)

.. _eq-13.3-14:

.. math::

	w_{\text{j}} = A_{\text{c,j}} \rho_{\text{c,j}} v_{\text{c,j}}

and compute the convective terms by (blocked segment only)

(13.3‑15a)

.. _eq-13.3-15a:

.. math::

	\begin{align}
	C_{\text{v,j}} = w_{\text{j}} \frac{e_{\text{c,j}} - e_{\text{c,j} - 1}}{z_{\text{m,j}} - z_{\text{m,j} - 1}} && \text{for } w_{\text{m,j}} > 0, w_{\text{j}} > 0
	\end{align}

(13.3‑15b)

.. _eq-13.3-15b:

.. math::

	\begin{align}
	C_{\text{v,j}} = w_{\text{j}}\frac{e_{\text{c,j} + 1} - e_{\text{c,j}}}{z_{\text{m,j} + 1} - z_{\text{m,j}}} && \text{for } w_{\text{mj}} < 0, w_{\text{j}} > 0
	\end{align}

(13.3‑15c)

.. _eq-13.3-15c:

.. math::

	\begin{align}
	C_{\text{v,j}} = 0, && \text{otherwise}
	\end{align}

The scheme also suppresses convection (sets :math:`C_{\text{v}} = 0`) if the
adjoining donor cell (either :math:`j-1` or :math:`j+1` depending on the sign of
:math:`w_{\text{m}}`) is blocked.

.. _section-13.3.2.4:

Moving Cladding on Bare Fuel
^^^^^^^^^^^^^^^^^^^^^^^^^^^^

The outermost radial fuel node in the fuel-pin model is duplicated in
CLAP for axial segments within the molten zone. The outer fuel
temperature in CLAP, :math:`T_{\text{f}}`, is reset to the value computed by
the pin model at the start of each heat-transfer time step. However,
between heat transfer time steps, the fuel temperature is computed in
CLAP; this is done primarily to achieve numerical stability, but it may
also improve the accuracy of results. The configuration of the SAS outer
fuel node and moving cladding is shown in :numref:`figure-13.3-1`. The transient
heat-balance equation for this outer node is

(13.3‑16)

.. _eq-13.3-16:

.. math::

	\pi\left( r_{\text{NR}}^{2} - r_{\text{NT}}^{2} \right)\rho_{\text{f}}c_{\text{pf}}\frac{\partial \text{T}_{\text{f}}}{\partial \text{t}} \
	+ 2\pi_{\text{NT}} k_{\text{f}} \frac{\partial \text{T}}{\partial \text{r}} \bigg\rvert_{\text{r}_{\text{NT}}} - \pi\left( r_{\text{NR}}^{2} \
	- r_{\text{NT}}^{2} \right) Q_{\text{NT}} = - P_{\text{r}} h \left( T_{\text{f}} - T_{\text{c}} \right)

where reference to axial core segment :math:`j` is tacitly implied.

It is assumed that the heat flow from the interior fuel nodes and the
heat generated in the outer fuel node (given by the second and third
terms in Eq. :ref:`13.3-16<eq-13.3-16>`) are relatively constant over one heat-transfer
time step. We can therefore approximate these terms by evaluating them
at the beginning of a heat-transfer time step. These terms are
represented symbolically by the constant :math:`C_{\text{f}}`, defined by

(13.3‑17)

.. _eq-13.3-17:

.. math::

	\begin{align}
	C_{\text{f}} = 2\pi r_{\text{NT}} k_{\text{f}} \frac{\partial \text{T}}{\partial \text{r}} \bigg\rvert_{\text{r}_{\text{NT}}} - \pi\left( r_{\text{NR}}^{2} - r_{\text{NT}}^{2} \right) Q_{\text{NT}}, && \text{at } t = t^{*}
	\end{align}

where :math:`t^{*}` is the beginning of the heat-transfer time step.

.. _figure-13.3-1:

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

	Schematic Showing the Interface of CLAP with the SAS Fuel Model

Because of the high heat-transfer coefficient to the molten cladding
(Eq. :ref:`13.2-28<eq-13.2-28>`), the fuel-surface temperature and molten cladding
temperature are closely coupled. For that reason, the term representing
the instantaneous heat loss to the molten steel (on the right side of
Eq. :ref:`13.3-16<eq-13.3-16>`) is retained as is, rather than incorporating it into the
constant :math:`C_{\text{f}}`. Combining Eqs. :ref:`13.3-16<eq-13.3-16>` and :ref:`13.3-17<eq-13.3-17>` and
rearranging, we have

(13.3‑18)

.. _eq-13.3-18:

.. math::

	\frac{\text{dT}_{\text{f,j}}}{\text{dt}} + \gamma_{\text{f}} = - \xi_{\text{f}} \left( T_{\text{f}} - T_{\text{c}} \right) \bigg\rvert_{\text{j}}

where

(13.3‑19)

.. _eq-13.3-19:

.. math::

	\gamma_{\text{f,j}} = \frac{C_{\text{f}}}{\left\lbrack \pi\left( r_{\text{NR}}^{2} - r_{\text{NT}}^{2} \right) \right\rbrack}\rho_{\text{f}}c_{\text{pf}} \bigg\rvert_{\text{j}}

(13.3‑20)

.. _eq-13.3-20:

.. math::

	\xi_{\text{f,j}} = \frac{P_{\text{r}} h}{\left\lbrack \pi\left( r_{\text{NR}}^{2} - r_{\text{NT}}^{2} \right)\rho_{\text{f}}c_{\text{pf}} \right\rbrack \bigg\rvert_{\text{j}}}

Combining Eqs. :ref:`13.2-20<eq-13.2-20>`, :ref:`13.2-27<eq-13.2-27>`, :ref:`13.2-42<eq-13.2-42>`, and :ref:`13.3-9<eq-13.3-9>`, one obtains the
nodal energy equation for moving cladding:

(13.3‑21)

.. _eq-13.3-21:

.. math::

	\frac{\text{dT}_{\text{c,j}}}{\text{dt}} + \gamma_{\text{c,j}} = \xi_{2} \left( T_{\text{f}} - T_{\text{c}} \right) \bigg\rvert_{\text{j}}

where

(13.3‑22)

.. _eq-13.3-22:

.. math::

	\gamma_{\text{c,j}} = \frac{C_{\text{v}}}{\left( A_{\text{c}}\rho_{\text{c}}c_{\text{pc}} \right) \bigg\rvert_{\text{j}}}

(13.3‑23)

.. _eq-13.3-23:

.. math::

	\xi_{2,\text{j}} = \frac{P_{\text{r}} h}{\left( A_{\text{c}}\rho_{\text{c}}c_{\text{pc}} \right) \bigg\rvert_{\text{j}}}

Equations :ref:`13.3-18<eq-13.3-18>` and :ref:`13.3-21<eq-13.3-21>` are converted to time-implicit difference
equations by approximating the time derivatives according to

(13.3‑24)

.. _eq-13.3-24:

.. math::

	\frac{\text{dT}_{\text{f,j}}}{\text{dt}} \simeq \frac{\Delta T_{\text{f,j}}^{n}}{\Delta t^{n}}

(13.3‑25)

.. _eq-13.3-25:

.. math::

	\frac{\text{dT}_{\text{c,j}}}{\text{dt}} \simeq \frac{\Delta T_{\text{c,j}}^{n}}{\Delta t^{n}}

and evaluating the temperatures on the right-hand side at the advanced
time :math:`t^{n+1}`. After substituting

(13.3‑26)

.. _eq-13.3-26:

.. math::

	T_{\text{c,j}}^{n + 1} = T_{\text{f,j}}^{n} + \Delta T_{\text{f,j}}^{n}

(13.3‑27)

.. _eq-13.3-27:

.. math::

	T_{\text{c,j}}^{n + 1} = T_{\text{c,j}}^{n} + \Delta T_{\text{c,j}}^{n}~ ,

the two difference equations are solved simultaneously for
:math:`\Delta T_{\text{f}}` and :math:`\Delta T_{\text{c}}` and second-order terms in :math:`\Delta t`
are discarded, with the following results

(13.3‑28)

.. _eq-13.3-28:

.. math::

	\Delta T_{\text{c,j}}^{n} = - \frac{\left\lbrack \xi_{2}\left( T_{\text{f}} - T_{\text{c}} \right) - \gamma_{\text{c}} \right\rbrack \Delta t}{1 + \left( \xi_{2} + \xi_{\text{f}} \right) \Delta t} \bigg\rvert_{\text{j}}^{n}

(13.3‑29)

.. _eq-13.3-29:

.. math::

	\Delta T_{\text{f,j}}^{n} = - \frac{\left\lbrack \xi_{\text{f}}\left( T_{\text{f}} - T_{\text{c}} \right) - \gamma_{\text{f}} \right\rbrack \Delta t}{1 + \left( \xi_{2} + \xi_{\text{f}} \right) \Delta t} \bigg\rvert_{\text{j}}^{n}

The heat flux :math:`\phi_{\text{c}}`, which is needed to evaluate the integral in
Eq. :ref:`13.3-2<eq-13.3-2>`, is computed using a finite difference form of the moving
clad energy Eq. :ref:`13.2-20<eq-13.2-20>`:

(13.3‑30)

.. _eq-13.3-30:

.. math::

	\phi_{\text{c}} = \frac{C_{\text{v}}}{P_{\text{r}}} + \frac{A_{\text{c}}\rho_{\text{c}}c_{\text{pc}}\Delta T_{\text{c}}}{P_{\text{r}} \Delta t}

.. _section-13.3.2.5:

Moving Cladding on Intact and Refrozen Cladding
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

The nodal energy equation for the moving cladding for this configuration
is obtained by combining Eqs. :ref:`13.2-20<eq-13.2-20>`, :ref:`13.2-32<eq-13.2-32>`, :ref:`13.2-42<eq-13.2-42>`, and :ref:`13.3-9<eq-13.3-9>`,
giving

(13.3‑31)

.. _eq-13.3-31:

.. math::

	\frac{\text{dT}_{\text{c,j}}}{\text{dt}} + \gamma_{\text{c}} = \xi_{2} \left( T_{\text{m}} - T_{\text{c}} \right) \bigg\rvert_{\text{j}}

where :math:`\gamma_{\text{c}}` and :math:`\xi_{2}` are defined in Eqs. :ref:`13.3-22<eq-13.3-22>` and :ref:`13.3-23<eq-13.3-23>`.

Solving this equation implicitly, i.e., replacing :math:`T_{\text{c}}` by the
updated value, approximating the time derivative by Eq. :ref:`13.3-25<eq-13.3-25>`, and
solving for :math:`\Delta T_{\text{c}}`, we get

(13.3‑32)

.. _eq-13.3-32:

.. math::

	\Delta T_{\text{c,j}}^{n} = - \frac{\left\lbrack \xi_{2}\left( T_{\text{m}} - T_{\text{c}} \right) - \gamma_{\text{c}} \right\rbrack \Delta t}{1 + \xi_{2} \Delta t} \bigg\rvert_{\text{j}}^{n}

The equation for the refrozen cladding temperature is found by
substituting Eqs. :ref:`13.2-30<eq-13.2-30>`, :ref:`13.2-31<eq-13.2-31>`, and :ref:`13.2-41<eq-13.2-41>` into the energy Eq.
:ref:`13.2-23<eq-13.2-23>`, resulting in

(13.3‑33)

.. _eq-13.3-33:

.. math::

	\Delta T_{\text{c,j}}^{n} = - \frac{\left\lbrack \xi_{2}\left( T_{\text{m}} - T_{\text{c}} \right) - \gamma_{\text{c}} \right\rbrack \Delta t}{1 + \xi_{2} \Delta t} \bigg\rvert_{\text{j}}^{n}

where

(13.3‑34)

.. _eq-13.3-34:

.. math::

	\xi_{3} = \frac{P_{\text{r}} k_{\text{c}}}{A_{\text{s}}\rho_{\text{s}}c_{\text{ps}}\Delta r_{\text{s}}}

Converting Eq. :ref:`13.3-33<eq-13.3-33>` to a finite-difference, implicit equation and
solving for :math:`\Delta T_{\text{s}}`, one obtains

(13.3‑35)

.. _eq-13.3-35:

.. math::

	\Delta T_{\text{s,j}}^{n} = - \frac{\left\lbrack \xi_{1}\left( T_{\text{i}} - T_{\text{s}} \right) - \xi_{3}\left( T_{\text{s}} - T_{\text{m}} \right) \right\rbrack \Delta t}{1 + \left( \xi_{1} + \xi_{3} \right) \Delta t} \bigg\rvert_{\text{j}}^{n}

The heat fluxes :math:`\phi_{\text{r}}` and :math:`\phi_{\text{hf}}` are required for the
fuel-pin boundary conditions (Eq. :ref:`13.3-3<eq-13.3-3>`) and to evaluate the rate of
melting (Eq. :ref:`13.2-25<eq-13.2-25>`). These are evaluated using discretized versions of
the energy Eqs. :ref:`13.2-20<eq-13.2-20>` and :ref:`13.2-23<eq-13.2-23>` after computing :math:`\phi_{\text{r}}` from
Eq. :ref:`13.2-30<eq-13.2-30>`. Then :math:`\phi_{\text{hf}}` is computed using Eq. :ref:`13.2-24<eq-13.2-24>`.

Several constraints apply to :math:`\phi_{\text{hf}}`. For positive :math:`\phi_{\text{hf}}`
(melting), more cladding cannot melt than actually exists in a frozen
state. This restriction can be expressed as

(13.3‑36)

.. _eq-13.3-36:

.. math::

	\phi_{\text{hf}} \leq \frac{A_{\text{s}}\rho_{\text{s}}\lambda_{\text{c}}}{\Delta t P_{\text{r}}}

In case :math:`\phi_{\text{hf}}` is adjusted to satisfy Eq. :ref:`13.3-36<eq-13.3-36>`, new values of
:math:`\phi`, :math:`\phi_{\text{c}}` and :math:`\phi_{\text{r}}` are computed from discretized
versions of energy Eqs. :ref:`13.2-20<eq-13.2-20>` and :ref:`13.2-23<eq-13.2-23>`, along with Eq. :ref:`13.2-24<eq-13.2-24>`, the
definition of :math:`\phi_{\text{hf}}`.

The constraint for negative :math:`\phi_{\text{hf}}` (freezing) arises from the
condition that more cladding cannot freeze then exists in the molten
state. Expressed mathematically, the limit is

(13.3‑37)

.. _eq-13.3-37:

.. math::

	\phi_{\text{hf}} \geq - \frac{A_{\text{c}}\rho_{\text{c}}\lambda_{\text{c}}}{\Delta t P_{\text{r}}}

.. _section-13.3.2.6:

Intact and Moving Cladding
^^^^^^^^^^^^^^^^^^^^^^^^^^

The molten cladding energy equation for this case is identical to that
for the case of moving cladding on intact and refrozen cladding. The
molten cladding temperature change is therefore given by Eq. :ref:`13.3-32<eq-13.3-32>`,
which is also displayed below

(13.3‑38)

.. _eq-13.3-38:

.. math::

	\Delta T_{\text{c,j}}^{n} = - \frac{\left\lbrack \xi_{2}\left( T_{\text{m}} - T_{\text{c}} \right) - \gamma_{\text{c}} \right\rbrack \Delta t}{1 + \xi_{2} \Delta t} \bigg\rvert_{\text{j}}^{n}

The quantities :math:`\phi_{1}, \phi_{2}`, and :math:`\phi_{\text{trial}}` are
computed from Eqs. :ref:`13.2-34<eq-13.2-34>` to :ref:`13.2-36<eq-13.2-36>`. Then :math:`\phi_{\text{hf}}` and
:math:`\phi_{\text{r}}` are computed according to Eq. :ref:`13.2-37<eq-13.2-37a>`, which we repeat
below

(13.3‑39a)

.. _eq-13.3-39a:

.. math::

	\begin{matrix}
	\phi_{\text{hf}} = 0, & \phi_{\text{r}} = \phi_{2} & \text{ for } \phi_{\text{trial}} \geq 0,
	\end{matrix}

(13.3‑39b)

.. _eq-13.3-39b:

.. math::

	\begin{matrix}
	\phi_{\text{hf}} = \phi_{\text{trial}}, & \phi_{\text{r}} = \phi_{1} & \text{ for } \phi_{\text{trial}} < 0.
	\end{matrix}

The case of :math:`\phi_{\text{trial}} \geq 0` corresponds to melting of the initial
cladding; the heat of fusion is included in :math:`\phi_{\text{r}}` which is the
negative of the heat flux to the pin surface.

For :math:`\phi_{\text{trial}} \leq 0`, we obtain a negative :math:`\phi_{\text{hf}}`
(freezing) which will result in the appearance of refrozen cladding
during the current time step. The constraint given by Eq. :ref:`13.3-37<eq-13.3-37>` is
checked and, if not satisfied, then (i) :math:`\phi_{\text{hf}}` is set equal to
the limit given by Eq. :ref:`13.3-37<eq-13.3-37>` and (ii) a new :math:`\phi_{\text{r}}` is computed
by

(13.3‑40)

.. _eq-13.3-40:

.. math::

	\phi_{\text{r}} = \phi_{2} + \phi_{\text{hf}}

.. _section-13.3.2.7:

Heat Loss to Structure
^^^^^^^^^^^^^^^^^^^^^^

At the onset of cladding motion, the structure is still relatively cool
and consequently constitutes a potentially significant heat sink for the
refreezing of cladding, particularly in small experimental
subassemblies. The local heat absorbed by the structure (per unit length
and unit time) is given by

(13.3‑41)

.. _eq-13.3-41:

.. math::

	Q_{\text{s,j}}^{n} = \frac{\theta P_{\text{w}}\left( T_{\text{m}} - T_{\text{w}} \right)}{\frac{\Delta r_{\text{w}}}{k} + \frac{1}{h}}

where

:math:`\theta` = multiplier on heat loss to structure (usually = 1)

:math:`P_{\text{w}}` = heated perimeter of the structure

:math:`T_{\text{w}}` = structure temperature

:math:`\Delta r_{\text{w}}` = half-thickness of structure.

The present CLAP model allows for the indirect transfer of heat form the
moving cladding to the structure by way of the frozen cladding. An
adjustment to the refrozen clad temperature :math:`\Delta T_{\text{s}}` is defined
by

(13.3‑42)

.. _eq-13.3-42:

.. math::

	T_{\text{s,j}}^{n + 1} = {\hat{T}}_{\text{s,j}}^{n + 1} + \Delta{\hat{T}}_{\text{s,j}}^{n}

where :math:`{\hat{T}}_{\text{s,j}}^{n + 1}` is the unadjusted refrozen
cladding temperature.

Normally, the temperature adjustment would be given by

(13.3‑43)

.. _eq-13.3-43:

.. math::

	\Delta {\hat{T}}_{\text{s,j}}^{n} = Q_{\text{s}}\frac{\Delta t}{A_{\text{s}}\rho_{\text{s}}c_{\text{ps}}} \bigg\rvert_{\text{j}}^{n}

However, this formulation would be unstable in the case of vanishing
:math:`A_{\text{s}}` or large :math:`\Delta t`. To achieve stability, we slightly alter
the physics; in particular, :math:`T_{\text{m}}` is replaced in Eq. :ref:`13.3-41<eq-13.3-41>` by
:math:`T_{\text{s,j}}^{n + 1}`. This adjustment in :math:`T_{\text{s}}` is made only
for segments where moving cladding is present; and, for these segments,
we would expect the frozen cladding temperature to be close to the
melting temperature, and hence the error in this approximation would be
small. With this modification the adjustment is given by

(13.3‑44)

.. _eq-13.3-44:

.. math::

	\Delta {\hat{T}}_{\text{s,j}}^{n} = \frac{\xi_{\text{w}}\Delta t^{n}\left( \Delta {\hat{T}}_{\text{s,j}}^{n + 1} - T_{\text{w}} \right)}{1 + \epsilon_{\text{w}}\Delta t^{n}}

where

(13.3‑45)

.. _eq-13.3-45:

.. math::

	\xi_{\text{w}} = \frac{\theta P_{\text{w}}}{A_{\text{s}}\rho_{\text{s}}C_{\text{ps}}\left( \frac{\Delta r_{\text{w}}}{k} + \frac{1}{h} \right)} \bigg\rvert_{\text{j}}^{n}

.. _section-13.3.3:

Refrozen Steel Area Calculation
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

The refrozen cladding mass conservation Eq. :ref:`13.2-21<eq-13.2-21>` with the source term
given by Eq. :ref:`13.2-25<eq-13.2-25>` is rewritten in the form

(13.3‑46)

.. _eq-13.3-46:

.. math::

	\rho_{\text{s}}\frac{\partial \text{A}_{\text{s}}}{\partial \text{t}} = - \frac{P_{\text{r}}\phi_{\text{hf}}}{\lambda} - A_{\text{s}}\frac{\partial\rho_{\text{s}}}{\partial \text{t}}.

The finite difference form of this equation is

(13.3‑47)

.. _eq-13.3-47:

.. math::

	A_{\text{s,j}}^{n + 1} = A_{\text{s,j}}^{n} - \frac{\left\lbrack \frac{\phi_{\text{hf}}P_{\text{r}}}{\lambda}\Delta t \bigg\rvert_{\text{j}}^{n} + A_{\text{s,j}}^{n} \left( \rho_{\text{s,j}}^{n + 1} - \rho_{\text{s,j}}^{n} \right) \right\rbrack}{\rho_{\text{s,j}}^{n + 1}}

After the density :math:`\rho_{\text{s,j}}^{n + 1}` is evaluated from
:math:`T_{\text{s,j}}^{n + 1}` using Eq. :ref:`13.2-38<eq-13.2-38>`, then the area
:math:`A_{\text{s,j}}^{n + 1}` is computed using Eq. :ref:`13.3-47<eq-13.3-47>`. The thickness of
the refrozen steel layer is computed from its area by assuming that the
refrozen steel forms an annular ring around the intact cladding.

.. _section-13.3.4:

Moving Cladding Area Calculation
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

The continuity Eq. :ref:`13.2-12<eq-13.2-12>` and source term (Eq. :ref:`13.2-25<eq-13.2-25>`) are rewritten
as

(13.3‑48)

.. _eq-13.3-48:

.. math::

	\rho_{\text{c}}\frac{\partial \text{A}_{\text{c}}}{\partial \text{t}} = - A_{\text{c}} \frac{\partial\rho_{\text{c}}}{\partial \text{t}} - C_{\text{a}} + \frac{P_{\text{r}}\varphi_{\text{hf}}}{\lambda}

where

(13.3‑49)

.. _eq-13.3-49:

.. math::

	C_{\text{a}} = \frac{\partial}{\partial \text{z}}\left( \rho_{\text{c}}A_{\text{c}}v_{\text{c}} \right)

The finite difference form of this equation is given by

(13.3‑50)

.. _eq-13.3-50:

.. math::

	A_{\text{c,j}}^{n + 1} = A_{\text{c,j}}^{n} + \frac{\left\lbrack {- A}_{\text{c,j}}^{n}\left( \rho_{\text{c,j}}^{n + 1} - \rho_{\text{c,j}}^{n} \right) - C_{\text{a}}\Delta t_{\text{n}} \
	+ \frac{{P_{\text{r}}\phi}_{\text{hf}}\Delta t}{\lambda} \bigg\rvert_{\text{j}}^{n} \right\rbrack}{\rho_{\text{c,j}}^{n + 1}} + A_{\text{mod,j}}^{n}

where :math:`A_{\text{mod}}` is a correction term, defined later.

The solution method uses donor-cell differencing to evaluate the
convective term :math:`C_{\text{a}}` in terms of the nodal fluxes :math:`w_{\text{j}}`
(see Eq. :ref:`13.3-14<eq-13.3-14>`) and nodal elevation (see Eq. :ref:`13.3-13<eq-13.3-13>`). Normally
:math:`A_{\text{mod}} = 0` and :math:`C_{\text{a}}` is computed as follows

(13.3‑51a)

.. _eq-13.3-51a:

.. math::

	\begin{align}
	C_{\text{a}} = \frac{w_{\text{j}}^{n} - w_{\text{j} - 1}^{n}}{z_{\text{m,j}} - z_{\text{m,j} - 1}} && \text{for } w_{\text{m,j}} \geq 0
	\end{align}

(13.3‑51b)

.. _eq-13.3-51b:

.. math::

	\begin{align}
	C_{\text{a}} = \frac{w_{\text{j} + 1}^{n} - w_{\text{j}}^{n}}{z_{\text{m,j} + 1} - z_{\text{m,j}}} && \text{for } w_{\text{m,j}} < 0
	\end{align}

Modifications to this formulation are made (i) to prevent cladding
transfer outside of the fuel and blanket region and (ii) to inhibit
cladding transfer to flooded or blocked nodes; the modifications consist
simply of setting the appropriate values of :math:`w_{\text{j}}` equal to zero.

Additionally, adjustments must be made for the case where the molten
cladding area exceeds the available area :math:`A_{\text{max}}`, defined by

(13.3‑52)

.. _eq-13.3-52:

.. math::

	A_{\text{max,j}}^{n + 1} = A_{\text{f,j}} - A_{\text{i,j}}^{n + 1} - A_{\text{s,j}}^{n + 1}

where

:math:`A_{\text{f}}` = total area for cladding allowed by the fuel

:math:`A_{\text{i}}` = area of the intact cladding (if any).

For the case where the computed :math:`A_{\text{c,j}}^{n + 1}` exceeds
:math:`A_{\text{max}}`, the following steps are taken: (i) the condition is
flagged by setting :math:`\text{NFULL}_{\text{j}} = 1` (initial value is 0), (ii) the
molten cladding area :math:`A_{\text{c,j}}^{n + 1}` is set equal to
:math:`A_{\text{max}}`, and (iii) the computed area for the donor segment is
later adjusted consistent with step (ii).

To accomplish step (iii), we first estimate the volume of cladding
convected into segment :math:`j`, denoted by :math:`\Gamma_{\text{j}}`, from continuity
considerations:

(13.3‑53)

.. _eq-13.3-53:

.. math::

	\Gamma_{\text{j}} \cong \left\lbrack A_{\text{max,j}}^{n + 1} - A_{\text{c,j}}^{n} - \frac{P_{\text{r}} \frac{\phi_{\text{hf}} \Delta t}{\lambda} \bigg\rvert_{\text{j}}^{n}}{\rho_{\text{c,j}}^{n + 1}} \right\rbrack \left( z_{\text{j} + 1 }z_{\text{j}} \right)

Rather than adjusting :math:`C_{\text{a}}` for the donor cell, we instead set
:math:`w_{\text{j}}` (or :math:`w_{\text{j}+1}`) equal to zero and compute a
correction :math:`A_{\text{mod}}`, which is the area change (negative) due to
the volume of fluid :math:`\Gamma_{\text{j}}` transferred to the flooded cell. For
example, for the case of the donor cell below the flooded cell we have

(13.3‑54)

.. _eq-13.3-54:

.. math::

	\begin{align}
	w_{\text{j}} = 0, A_{\text{mod,j} - 1}^{n} = - \frac{\Gamma_{\text{j}}}{\left( z_{\text{j}} - z_{\text{j} - 1} \right)} && \text{for } \text{NFULL}_{\text{j}} = 1, w_{\text{m,j} - 1} \geq 0.
	\end{align}

Similarly, for the case of the donor cell above the flooded cell we have

(13.3‑55)

.. _eq-13.3-55:

.. math::

	\begin{align}
	w_{\text{j} + 1} = 0, A_{\text{mod,j} + 1}^{n} = - \frac{\Gamma_{\text{j}}}{\left( z_{\text{j} + 2} - z_{\text{j} + 1} \right)}, && \text{for } \text{NFULL}_{\text{j}} = 1, w_{\text{m,j} + 1} < 0.
	\end{align}

.. _section-13.3.5:

Reactivity Calculation
~~~~~~~~~~~~~~~~~~~~~~

The reactivity change due to cladding relocation is computed by

(13.3‑56)

.. _eq-13.3-56:

.. math::

	\Delta K = \sum_{\text{j}}\left( m_{\text{j}}^{n} - m_{\text{j}}^{o} \right) \cdot W_{\text{j}}

where

:math:`m_{\text{j}}^{n}` = mass of cladding in the :math:`j`-th segment

:math:`m_{\text{j}}^{o}` = original mass of cladding in the :math:`j`-th
segment

:math:`W_{\text{j}}` = cladding reactivity worth distribution.

In the molten region, the mass is computed by

(13.3‑57)

.. _eq-13.3-57:

.. math::

	m_{\text{j}}^{n} = A_{\text{c}}\rho_{\text{c}}\Delta z \bigg\rvert_{\text{j}}^{n}

(13.3‑58)

.. _eq-13.3-58:

.. math::

	\text{where }\Delta z_{\text{j}} = z_{\text{j} + 1} - z_{\text{j}}

Outside the molten region the mass is evaluated using

(13.3‑59)

.. _eq-13.3-59:

.. math::

	m_{\text{j}}^{n} = m_{\text{j}}^{0} + A_{\text{s}}\rho_{\text{s}}\Delta z \bigg\rvert_{\text{j}}^{n} + A_{\text{c}}\rho_{\text{c}} \Delta z \bigg\rvert_{\text{j}}^{n} .

Due to roundoff and other possible small errors, the total mass of
cladding may not be absolutely conserved which, if left uncorrected,
would create anomalous reactivity effects. To insure mass conservation,
the computed mass distribution is renormalized by the following steps:

(13.3‑60)

.. _eq-13.3-60:

.. math::

	M^{n} = \sum_{\text{j}}{m_{\text{j}}^{n}} ,

(13.3‑61)

.. _eq-13.3-61:

.. math::

	M^{o} = \sum_{\text{j}}{m_{\text{j}}^{o}} ,

(13.3‑62)

.. _eq-13.3-62:

.. math::

	m_{\text{j}}^{n} \cdot \frac{M^{o}}{M^{n}} \rightarrow m_{\text{j}}^{n}

where the operation A → B indicates that A is substituted for B.

.. _section-13.3.6:

Moving Cladding Velocity Calculation
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

The cladding momentum equation is solved in a separate subroutine
(TSCLD2) that is called after the coolant dynamics equations have been
solved. The terms :math:`\frac{\partial \text{p}}{\partial \text{z}}` and :math:`F_{\text{v}}` in the cladding momentum
equation are evaluated at time :math:`t_{\text{n}} + 1/2 \Delta t^{n}` (by
averaging values at :math:`t^{n}` and :math:`t^{n+1}`) to
achieve maximum accuracy with respect to sodium-vapor-induced cladding
accelerations.

The pin friction term, :math:`\left( \frac{\partial \text{p}}{\partial \text{z}} \right)_{\text{fr}}`, is evaluated at the end
of the time step to obtain maximum stability under conditions of
strongly accelerating thin cladding films. Substituting

(13.3‑63)

.. _eq-13.3-63:

.. math::

	v_{\text{c,j}}^{n + 1} = v_{\text{c,j}}^{n} + \Delta v_{\text{c,j}}^{n}

into Eqs. :ref:`13.2-15<eq-13.2-15>` and :ref:`13.2-16<eq-13.2-16a>` and retaining only linear terms, we obtain
an equation of the form

(13.3‑64)

.. _eq-13.3-64:

.. math::

	F_{\text{p,j}}^{n + 1} = y_{1,\text{j}}^{n} + y_{2,\text{j}}^{n} \Delta v_{\text{c,j}}^{n}

where

(13.3‑65)

.. _eq-13.3-65:

.. math::

	\begin{align}
	y_{1,\text{j}}^{n + 1} = \frac{32 \mu_{\text{c}} v_{\text{c}}}{D_{\text{c}}^{2}} \bigg\rvert_{\text{j}}^{n} && \text{for } \text{Re}_{\text{j}}^{n} < \left( \text{Re} \right)_{\text{break}}
	\end{align}

(13.3‑66)

.. _eq-13.3-66:

.. math::

	\begin{align}
	y_{2,\text{j}}^{n} = \frac{32 \mu_{\text{c}}}{D_{\text{c}}^{2}} \bigg\rvert_{\text{j}}^{n} && \text{for } \text{Re}_{\text{j}}^{n} < \left( \text{Re} \right)_{\text{break}}
	\end{align}

(13.3‑67)

.. _eq-13.3-67:

.. math::

	\begin{align}
	y_{1,\text{j}}^{n} = \frac{b_{\text{f}}\rho_{\text{c}}v_{\text{c}} \left| v_{\text{c}} \right|}{2D_{\text{c}}} \bigg\rvert_{\text{j}}^{n} && \text{for } \text{Re}_{\text{j}}^{n} \geq \left( \text{Re} \right)_{\text{break}}
	\end{align}

(13.3‑68)

.. _eq-13.3-68:

.. math::

	\begin{align}
	y_{2,\text{j}}^{n} = \frac{b_{\text{f}}\rho_{\text{c}} \left| v_{\text{c}} \right|}{D_{\text{c}}} \bigg\rvert_{\text{j}}^{n} && \text{for } \text{Re}_{\text{j}}^{n} \geq \left( \text{Re} \right)_{\text{break}}
	\end{align}

The parameter :math:`C_{\text{m}}` is used to represent the momentum convective
term

(13.3‑69)

.. _eq-13.3-69:

.. math::

	C_{\text{m,j}} = v_{\text{c}}\frac{\partial \text{v}_{\text{c}}}{\partial \text{z}} \bigg\rvert_{\text{j}}

This parameter is evaluated by alternative formulae depending upon
whether the segment is filled with cladding. For :math:`\text{NFULL} \left( j \right)=1`,
indicating a filled segment, the term is estimated by the formula

(13.3‑70)

.. _eq-13.3-70:

.. math::

	C_{\text{m,j}}^{n} = v_{\text{c,j}}^{n}\frac{\partial \text{v}_{\text{c}}}{\partial \text{z}} \bigg\rvert_{\text{j}}^{n + 1/2}

where the derivative :math:`\frac{\partial \text{v}_{\text{c}}}{\partial \text{z}}` is
obtained by rearranging the continuity equation:

(13.3‑71)

.. _eq-13.3-71:

.. math::

	\frac{\partial \text{v}_{\text{c}}}{\partial \text{z}} = \frac{1}{\rho_{\text{c}}A_{\text{c}}}\left\lbrack \frac{\phi_{\text{hf}}P_{\text{r}}}{\lambda} - \frac{\partial\left( \rho_{\text{c}} A_{\text{c}} \right)}{\partial \text{t}} - v_{\text{c}}\frac{\left( \rho_{\text{c}}A_{\text{c}} \right)}{\partial \text{z}} \right\rbrack

and evaluating at time :math:`t^{n} + 1/2 \Delta t^{n}`. Two
alternative equations are obtained, one for each flow direction, as
given by the following:

(13.3‑72)

.. _eq-13.3-72:

.. math::

    \frac{\partial \text{v}_{\text{c}}}{\partial \text{z}_{\text{j}}} \bigg\rvert^{n + 1/2} = \frac{1}{\left( \rho_{\text{c,j}} A_{\text{c,j}} \right)^{n + 1/2}} \left\lbrack \frac{\phi_{\text{h}_{\text{f}}}P_{\text{r}}}{\lambda_{\text{c}}} \bigg\rvert_{\text{j}}^{n} - \rho_{\text{c,j}}^{n} \frac{A_{\text{c,j}}^{n + 1} - A_{\text{c,j}}^{n}}{\Delta t^{n}} \\
    - A_{\text{c,j}}^{n}\frac{\rho_{\text{c,j}}^{n + 1} - \rho_{\text{c,j}}^{n}}{\Delta t^{n}} - v_{\text{c,j}}^{n}\frac{\left( \rho_{\text{c,j} + 1}A_{\text{c,j} + 1} \right)^{n + 1/2} - \left( \rho_{\text{c,j}}A_{\text{c,j}} \right)^{n + 1/2}}{z_{\text{j} + 1} - z_{\text{j}}} \right\rbrack

for :math:`w_{\text{m,j}}^{n} < 0`

(13.3‑73)

.. _eq-13.3-73:

.. math::

    \frac{\partial \text{v}_{\text{c}}}{\partial \text{z}_{\text{j}}} \bigg\rvert^{n + 1/2} = \frac{1}{\left( \rho_{\text{c,j}}A_{\text{c,j}} \right)^{n + 1/2}} \left\lbrack \frac{\phi_{\text{h}_{\text{f}}}P_{\text{r}}}{\lambda_{\text{c}}} \bigg\rvert_{\text{j}}^{n} - \rho_{\text{c,j}}^{n}\frac{A_{\text{c,j}}^{n + 1} - A_{\text{c,j}}^{n}}{\Delta t^{n}} \\
    - A_{\text{c,j}}^{n}\frac{\rho_{\text{c,j}}^{n + 1} - \rho_{\text{c,j}}^{n}}{\Delta t^{n}} - v_{\text{c,j}}^{n}\frac{\left( \rho_{\text{c,j} + 1}A_{\text{c,j} + 1} \right)^{n + 1/2} - \left( \rho_{\text{c,j}}A_{\text{c,j}} \right)^{n + 1/2}}{z_{\text{j} + 1} - z_{\text{j}}} \right\rbrack

for :math:`w_{\text{m,j}}^{n} < 0`

where quantities on the right-hand side of Eqs. :ref:`13.3-72<eq-13.3-72>` and :ref:`13.3-73<eq-13.3-73>` with
the superscript :math:`n + 1/2` are evaluated at time :math:`t^{n} + \Delta t^{n/2}` as follows

(13.3‑74)

.. _eq-13.3-74:

.. math::

    \left( \ \ \right)^{n + 1/2} = 1/2 \left\lbrack \left( \ \ \right)^{n} + \left( \ \ \right)^{n + 1} \right\rbrack

For :math:`\text{NFULL}_{\text{j}} = 0`, indicating an unblocked segment, we use the
identity

(13.3‑75)

.. _eq-13.3-75:

.. math::

    v_{\text{c}} \frac{\partial \text{v}_{\text{c}}}{\partial \text{z}} = 1/2 \frac{\partial \text{v}_{\text{c}}^{2}}{\partial \text{z}}

and evaluate the momentum convection term as follows

(13.3‑76)

.. _eq-13.3-76:

.. math::

	\begin{align}
    C_{\text{m,j}} = \frac{\frac{1}{2} \left\lbrack \left( v_{\text{c,j} + 1}^{n} \right)^{2} - \left( v_{\text{c,j}}^{n} \right)^{2} \right\rbrack}{\left( z_{\text{m,j} + 1} - z_{\text{m,j}} \right)} && \text{for } w_{\text{m,j}}^{n} < 0
    \end{align}

(13.3‑77)

.. _eq-13.3-77:

.. math::

	\begin{align}
    C_{\text{m,j}} = \frac{\frac{1}{2} \left\lbrack \left( v_{\text{c,j}}^{n} \right)^{2} - \left( v_{\text{c,j} - 1}^{n} \right)^{2} \right\rbrack}{\left( z_{\text{m,j}} - z_{\text{m,j} - 1} \right)} && \text{for } w_{\text{m,j}}^{n} \geq 0
    \end{align}

The discretized version of he moving-cladding momentum (Eq. :ref:`13.2-19<eq-13.2-19>`),
using the linear approximation for :math:`F_{\text{p}}` (Eq. :ref:`13.3-64<eq-13.3-64>`), and the
symbolic representation for the momentum convective term (Eq. :ref:`13.3-69<eq-13.3-69>`),
is given by

(13.3‑78)

.. _eq-13.3-78:

.. math::

	\left( \rho_{\text{c,j}} \right)^{n + 1/2} \left\lbrack \frac{\Delta v_{\text{c,j}}^{n}}{\Delta t^{n}} + C_{\text{m,j}} + g \right\rbrack = - \frac{\partial \text{P}}{\partial \text{z}} \bigg\rvert_{\text{j}}^{n + 1/2} - \frac{A_{\text{v}}}{A_{\text{c}}}F_{\text{v}} \bigg\rvert_{\text{j}}^{n + 1/2} - y_{1,\text{j}}^{n} - y_{2,\text{j}}^{n} \Delta v_{\text{c,j}}^{n} -
	\begin{cases}
    0 & \text{if freezing} & \text{(13.3‑78a)} \\
    \frac{P_{\text{r,j}}^{n}{\phi}_{\text{hf,j}}^{n} \left( v_{\text{c,j}}^{n} + \Delta v_{\text{c,j}}^{n} \right)}{A_{\text{c,j}}^{n + 1/2}{\lambda}_{\text{c,j}}^{n}} & \text{if melting} & \text{(13.3‑78b)} \\
    \end{cases}

Solving this equation for :math:`\Delta v_{\text{c}}`, one obtains the form used in
the CLAP program.