Harhara and Hasan BMC Chemical Engineering

https://doi.org/s42480-020-0029-1
(2020) 2:5
BMC Chemical Engineering

RESEARCH ARTICLE Open Access

Dynamic modeling of heat exchanger

tube rupture
Ahmed Harhara and M.M. Faruque Hasan*

Abstract
One fault that occurs with heat exchangers is a tube rupture, an overpressure scenario in which high pressure fluid
flows into the low pressure region. It is a serious safety concern that may lead to significant damage. Accurate
prediction of the pressure build-up after a rupture is critical to determine the appropriate size of a relief device and
avoid exceeding allowable pressure limits. This paper describes a model-based step-by-step methodology to predict
dynamic pressure profiles during tube rupture for liquid-liquid, vapor-liquid, and flashing liquid-liquid systems. The
transient effects of the relief valve are considered. The effects of choked flow must also be considered for accurate
maximum pressure predictions. Using a dimensionless analysis, the pressure ratio and density ratio are shown to
significantly impact the severity of this incident. Results show that vapor-liquid systems result in the highest pressure
surges.
Keywords: Process safety, Modeling shell and tube heat exchanger, Tube rupture

Introduction fluid enters the shell-side, increasing the system pressure.

In process industries, pressure relief valves are commonly At the set point, a pressure relief valve opens, venting
referred to as the last line of defense. Code requires that shell-side fluid to the flare. With heat exchangers in the
when a process fails, regardless of its built-in redundan- offshore industry having tube-side pressures into the hun-
cies, there exists an independent protection system (e.g., dreds of bar, dynamic modeling is critical to selecting an
pressure relief valve) powered by the medium it pro- appropriate relief device. An additional benefit of accu-
tects. [1] Therefore, correct pressure relief valve sizing rately sizing overpressure scenarios via dynamic analysis
is critical to ensure the protection of equipment, prop- is the potential to reduce the weight of equipment, low-
erty, and life. Fortunately, for most overpressure scenarios, ering capital costs. [3] This is especially important in the
pressure relief standards set by API 521 help guide engi- offshore industry, where topside weight reduction is an
neers on how to tackle various overpressure scenarios. important design criterion.
For heat exchanger tube ruptures, an overpressure event Before modeling the tube rupture phenomena, it is help-
that has been documented in industry, API 521 recom- ful to understand its relevance to industry. Tube ruptures
mends a dynamic analysis be performed when there exist can occur in the chemical, oil, gas, offshore, and even the
large differences in pressure. [2] However, a framework for nuclear industry. Knowing where tube ruptures are more
dynamically modeling tube rupture scenarios is not given. common allows one to scrutinize one set of systems over
A sketch of the tube rupture scenario is given in Fig. 1. others. For reasons including the need to protect propri-
This scenario involves a heat exchanger with high pres- etary process information, reliability data on the failure
sure on the tube-side and low pressure on the shell-side. of equipment in industry is hard to come by. However,
Exchanger tubes over time can weaken (shown as a zigzag there are efforts to counter this such as the Mary Kay
stress pattern in Fig. 1). If unnoticed or not replaced in O’Connor Process Safety Center’s Instrument Reliability
time, this weak point may fail. High pressure tube-side Network [4].
Because of this, it is often difficult to find statistics on
*Correspondence: [email protected] the frequency of tube ruptures. Moreover, some failure
1
Artie McFerrin Department of Chemical Engineering, Texas A&M University, rate databases may record a tube rupture failure under a
TX 77843-3122 College Station, USA

© The Author(s). 2020 Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0
International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and
reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the
Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver
(http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 2 of 21

Fig. 1 Sketch of heat exchanger tube rupture

different name, making it even harder to estimate the like- capture the pressure-time relationship. This was followed
lihood of this event. Fortunately, there are some databases by a stress analysis to determine if shell failure was likely to
that give us a clue as to where tube ruptures occur. The occur. Simpson developed a different shell-average pres-
United Kingdom’s Hydrocarbon Release Database col- sure model based on the assumption that a tube failure
lected data on offshore incidents. [5] By filtering out will result in an expanding translating spherical explo-
exchanger related releases, Table 1 shows an estimate of sion [7]. The models by Fowler et al. and Simpson both
which areas in a plant experienced a tube rupture. From accounted for the inertia of the shell-side fluid and relief
the table, it is evident that most offshore tube ruptures line. They did not however simulate the sudden pres-
occur in gas compression systems. This is expected since sure increase following a tube break (commonly referred
gas compression systems contain a series of intercool- to as water-hammer). Sumaria et al. modeled the tran-
ers, with high pressure gases on the tube-side, and low sient pressure surge by using a lumped parameter model
pressure cooling water on the shell-side. to incorporate the compressibility of the shell-side fluid
The model by Fowler et al. was one of the earliest and relief line. [8] This model was applied to a high
attempts to describe the tube rupture transient process. pressure gas entering a low pressure liquid. Ennis et al.
[6] A shell-average pressure model was developed to extended this work by considering the case of a tube-side
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 3 of 21

Table 1 Recorded Offshore Systems Experiencing Tube Rupture adequate protection of the heat exchanger. There are some
Incidents issues that exist with this method however. The first is the
Severity concern that the shell-side of the exchanger can fail before
System Major Significant Minor Total the relief valve can open. The reasoning behind this is
Flare/HP 1 1 that tube rupture scenarios result in fast-acting transient
Gas compression 4 21 13 38 shockwaves which can reach peak pressure in millisec-
onds. [16] This is especially the case for high pressure
Processing/gas/dehydration 3 5 8
differences between the shell and tube-side.
Processing/gas/LPG/condensate 1 2 3
Another concern with this technique is that the orifice
Processing/gas/sour gas 1 1 calculation method ignores the shell-side throughout the
Processing/oil/oil treatment 1 1 calculation. This is a problem because the tube rupture
Separation/gas production 1 1 failure scenario is entirely focused on the shell-side failing
Utilities/gas/fuel gas 4 4 since it is the low pressure side of the exchanger. Orifice-
style calculations overlook how the shell-side pressure
Utilities/oil/heat transfer oil 1 1
profile develops throughout the tube failure — the effects
Total 5 28 25 58
of mixing between the shell and tube-side fluids are
ignored. Lastly, orifice calculations are more appropriate
for tube rupture calculations that do not have a sub-
flashing liquid and accounted for the spatial and tempo- cooled liquid on the tube-side. Subcooled liquids will flash
ral aspects of the attached piping system. [9] Furthermore, once entering the low pressure shell-side, further com-
Botros addressed the importance of accounting for plicating the calculation of the required relief rate. For
piping systems and specifying appropriate boundary heat exchanger tube ruptures, the goal is to prevent the
conditions [10]. shell-side of the exchanger from failing. Thus, a dynamic
Cassata et al used commercial transient analysis soft- simulation gives a more complete picture of the overpres-
ware to model a tube rupture by performing a series of sure scenario, increasing a plant’s confidence in its ability
simulations with the burst occurring at the U bend. [11] to withstand a tube rupture.
Pressure profiles at various points along the exchanger The novel contributions of this paper include a step-by-
were generated. Ewan et al calculated pressure profiles step methodology on dynamically sizing a relief system
when varying the location of the tube burst. [12] For for a tube rupture. The model will be useful in answering
relief devices close to the tube burst, pulse widths were critical safety related questions such as what is the max-
higher than bursts far from the relief device. Nagpal eval- imum pressure the shell-side reaches, how long after the
uated several tube rupture cases using dynamic simula- tube rupture is maximum pressure reached, and how long
tions and provides a chart for when dynamic simulations is the shell design pressure exceeded. Examples on how
are recommended. [13] The simulation results show that to apply this methodology to a liquid-liquid, vapor-liquid,
pressure relief devices can be used to protect liquid- and flashing liquid-liquid system are given. This paper
filled exchangers. Aside from simulating the severity of explains how to use process simulation software to gen-
this event through consequence modeling, Acosta and erate the necessary thermophysical properties required
Siu examined the risk of this scenario through the use of for accurate pressure profile predictions. Finally, a dimen-
dynamic event trees. [14] Lastly, a Joint Industry Project sionless analysis of properties that affect tube ruptures is
conducted a full-scale heat exchanger tube rupture. [15] performed.
These experimental results were compared against blind The article is organized as follows. “Background”
simulations performed by consultants. The results sup- section provides a background on this overpressure event,
port using one dimensional dynamic models to model the in addition to examining the common causes of heat
pressure profiles during tube ruptures. exchanger tubes failing. “Background” section also pro-
One technique on estimating the relief rate during a vides a brief overview of ASME VIII pressure vessel code
tube rupture is to apply a steady state calculation. This and its relation to the tube rupture scenario. This helps
is done by first assuming a complete tube rupture has to specify an acceptable pressure limit during a tube rup-
occurred (referred to in the literature as a guillotine rup- ture event. Lastly, “Background” section also describes the
ture). A standard hydraulic orifice calculation is then differences in relief devices used to protect exchangers
performed to determine what is the rate of fluid entering and their associated response times. “Methods” section
the shell-side. This number is then doubled due to the fact models the tube rupture. “Background” section first lists
that fluid is entering from both ends of the ruptured tube. the assumptions in this model, followed by the modeling of
Finally, this rate is then inputted into a relief valve calcu- each component that is fed into a tube rupture dynamic
lation in order to determine what size valve is required for pressure profile algorithm. “Results & Discussion” section
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 4 of 21

illustrates how the equations in “Methods” section are scenarios, pressure vessels must vent excess pressure to
used in a practical setting, with examples on liquid- a flare, recycled to the process, or to the environment
liquid, vapor-liquid, and flashing liquid-liquid systems. (depending on the fluid being vented). For the most part,
A discussion on the differences between these sys- pressure relief valves are set at or below the MAWP. There
tems is included. “Dimensionless Analysis of Tube Rupt- are exceptions to this when dealing with multiple relief
ure” section describes a dimensionless analysis performed valves. Note that the MAWP and the design pressure often
for a liquid-liquid tube rupture. “Mitigating Tube Rup- refer to the same pressure and are used interchangeably.
ture Scenarios” section includes strategies on mitigating a This paper makes no distinction between the two terms.
tube rupture. This is followed by this paper’s conclusion, For a single pressure relief valve protecting a system,
in “Conclusion” section. a 10% accumulation in pressure beyond the MAWP is
allowed [22]. For multiple relief valves providing protec-
Background tion, the accumulation is allowed to increase to 16% above
A tube in an exchanger can fail from fatigue, corrosion, the MAWP [22]. The last common overpressure scenario
increased temperature, and conflict with baffling. [17] A is for dealing with external fire, in which case, accumu-
case study from SABIC Europe shows baffle hammering lation can increase up to 21% beyond the MAWP [22].
as having the potential to damage tubes. [18] Depending For pressure vessels, the hydrotest pressure is tested any-
on the tolerance between the baffle and tubes, excessive where between 1.3 to 1.5 times the MAWP [23]. This work
movement by the tubes will result in continued strik- assumes a hydrotest pressure of 1.5 times the design pres-
ing against the baffle. Over time, this has the possibility sure. This work can easily be reformulated to apply to a
to form a cavity. This creates a weak point in the tube, hydrotest pressure of 1.3 times the design pressure, how-
effectively reducing the design pressure of the tube. Simi- ever. Lastly, because of the low probability of tube ruptures
larly, Shahrani et al studied the failure mechanism of heat occurring in comparison with other more common relief
exchanger tubes in a plant. [19] It was concluded that fret- events, and because a tube rupture occurs over a short
ting between the baffles and tubes, in combination with period of time, the shell-side hydrotest pressure is set as
the corrosion due to the presence of sulfur, led to the tube the upper bound. Some plants may not even consider
failing. Vibration from equipment such as compressors an overpressure scenario up to the hydrotest pressure to
can also cause tube failures via baffle hammering. One be a layer of protection analysis (LOPA) scenario. [24]
way to counter this is to reduce the velocity to minimize Figure 2 shows an overview of where the different pres-
vibrations. For new exchangers, materials of construction sures are located in relation to one another. The range
and/or adding design features such as baffle spacing may of pressures for a tube rupture scenario can theoretically
also allow for increased vibration. vary from the shell-side operating pressure to the tube-
Metal erosion from high velocity particulates over time side hydrotest pressure. A normal operating window for
can also thin tubes and result in a rupture. [17] The ero- the shell-side is between its operating pressure and design
sion rate primarily depends on the velocity, diameter, pipe pressure. Somewhere between these pressures, a pressure
thickness, and mass of particulates. The impact angle also relief device is set to open. Similarly, the tube-side also has
affects the severity of erosion. Because of the need for the an operating window between its operating pressure and
fluid to change direction, in a piping system, pipe bends design pressure, as well as a relief device set between these
frequently exhibit higher rates of erosion. [20] Another pressures. The next threshold after the shell-side design
case study by Khilnaney examines the failures of tube rup- pressure is its hydrotest pressure. When modeling the
tures at a nuclear plant cooling water exchanger. [21] In dynamic pressure profile, any pressure between the oper-
this study, hydrogen sulfide was detected in the cooling ating pressure and hydrotest pressure will be considered
tower, thus suspecting a tube leak had occurred. It was acceptable. Beyond this pressure, the severity of a tube
found that in 50% of these tubes, fluid elastic instability rupture and the likelihood of the shell-side failing increase
across the tube may have played a role in their failure. with pressure (shown as the region in yellow-to-red region
Lastly, corrosion from harsh chemicals such as sulfur and in Fig. 2).
caustic chemicals can weaken tubes to the point of failure. For most overpressure scenarios, relief valves are able to
ASME VIII code governs the performance of pressure open in time for pressure increases. With heat exchanger
vessels. It specifically applies to pressure vessels operating tube ruptures however, there is some concern that tube
at or exceeding 15 psig. [22] Most exchangers in a plant ruptures (which can generate peak pressures in millisec-
have operating pressures that fall under this criterion onds) can overpressure the shell-side of an exchanger
and are consequently governed by ASME VIII code. The before a pressure relief valve can take effect. A 2002 study
code also defines a maximum allowable working pressure prepared by Pipeline Simulation and Integrity Ltd for the
(MAWP) based on an equipment’s material, thickness, United Kingdom’s Health and Safety Executive tested the
and other mechanical properties. During overpressure respons times of pressure relief valves and rupture disks.
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 5 of 21

Fig. 2 Range of pressures for tube rupture scenario

[16] The study found that fast-acting relief valves do exist. liquid) has a path to exit the exchanger. Phases II, III, and
Note that if a vendor includes a datasheet that lists the IV are modeled with the following assumptions:
response time, that should be used in the model. An if-
• Only one tube is assumed to rupture.
statement can be included that doesn’t allow any fluid to
• The rupture is assumed to be a full-bore rupture.
vent through the relief valve from the initial rupture until
• A tube rupture exposes both ends of the severed
the response time has been met.
tube, allowing fluid to enter the shell-side at twice the
Methods cross-sectional area of one tube.
• The tube-side fluid enters the shell-side as an
This paper examines the three most common tube
rupture scenarios: the liquid-liquid system, the isentropic nozzle flow.
• The effects of temperature are ignored (e.g., increased
vapor-liquid system, and the flashing liquid-liquid
system. For all shell and tube exchanger configura- heating resulting in the buildup of pressure).
• An infinite reservoir of tube-side fluid exists. This is a
tions, tube ruptures can be divided into four distinct
phases [25]: conservative assumption. Without this assumption,
one would need to determine what is the quantity of
• Phase I: Sudden rupture tube-side fluid available (from upstream and
• Phase II: Development of transient shell-side pressure downstream units). One can also determine what is
wave the rate of tube-side pressure loss (and decrease in
• Phase III: Shell-side liquid discharges through velocity) as the fluid is being lost to the shell-side.
pressure relief system • Pressure relief systems have an instantaneous
• Phase IV: Shell-side vapor discharges through response time (i.e., zero milliseconds).
pressure relief system • No credit is taken for outflow via inlet or outlet
piping. The phrase "taking credit for outflow" allows
The first phase generates a shock wave from the rup- one to reduce the relieving requirements for
ture of the tube. This is, however, considered to be nearly overpressure scenarios by accounting for fluid leaving
instantaneous and is not modeled. With an open path, the system. This however comes with a set of
the high pressure side fluid enters the low pressure side. requirements to ensure that an open path exists (and
Because the tube-side fluid acts against the shell-side is adequately sized) for fluid to escape. For this
fluid, the shell-side fluid is pushed into the relief system. model, all excess pressure must be vented through
Once this is complete, the tube-side fluid (either vapor or available relief systems.
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 6 of 21

• All areas of the exchanger’s shell-side are in function that yields a mass flux for all pressures between
equilibrium. The effects of localized pressure are the shell and tube-side pressure.
ignored.
• The exchanger shell-side can safely handle pressures
up to the hydrotest pressure. G = A1 × P(t)2 + A2 × P(t) + A3 (4)

What follows describes how one can model and predict Vapor Flow Rate
the liquid and vapor flow rates, flashing liquid flow rates, The methodology for the calculation of vapor flow
vapor densities, mass fluxes from the tube to shell-side through the ruptured tube is very similar to that of
and from the shell-side to outside through a pressure relief the liquid flow rate. A series of isentropic flashes from
valve, and the dynamic pressure profiles as a consequence high pressure to low pressure conditions are performed.
of tube rupture in a shell and tube heat exchanger. Using the relationship between pressure and specific vol-
ume, Equations 1 - 3 are solved to yield a mass flux
Liquid Flow Rate as a function of pressure. Up to this point, the solution
The volumetric energy balance for isentropic nozzle flow method for the liquid flow rate and the vapor flow rate
is given by Eq. 1.[26] are identical. However, for most vapor-liquid tube rup-
tures, a choked flow condition will be present. The mass
 P2 flux function, G, when calculated, will show whether or
2
2× P1 v dP not a choked condition exists. If that is the case, the
G = (1)
v2t mass flux versus pressure must be adjusted to reflect
the choked condition. This is fitted using a cubic fit in
Where G is the tube-side mass flux entering the shell the following format (Eq. 5), where α1, α2, α3, and α4
(kg/s-m2 ), v is the specific volume of the tube-side fluid are the coefficients that give the best fit (e.g., method of
(m3 /kg), P1 and P2 are the pressure (Pa) the tube-side fluid least squares). An example of how this is performed is
experiences at the inlet and outlet, respectively, and t is the presented later.
fluid condition at the throat of the nozzle where the cross-
sectional area is minimized. Assuming a constant pressure
step size, the integration portion in Eq. 1 can be expressed G = α1 × P(t)3 + α2 × P(t)2 + α3 × P(t) + α4 (5)
via Eq. 2.
Flashing Liquid Flow Rate
⎛ ⎞ The flow rate of a flashing liquid entering the shell-side
 P2  n−1 incorporates techniques used in calculating the liquid
h
v dP = ⎝v1 + vn + 2 × vj ⎠ (2) flow rate and vapor flow rate. First, a series of isentropic
P1 2
j=2 flashes from the tube-side pressure to the shell-side pres-
sure are performed. Equations 1 - 3 are solved to yield
Where v is the specific volume of the tube-side fluid
a mass flux as a function of pressure. If a choked condi-
(m3 /kg), h is the constant pressure step size chosen for
tion does result, then the mass flux should be fitted as
summation purposes (Pa), n is the index for fluid con-
a cubic equation (Eq. 5). These steps are identical to the
ditions at the assumed throat pressure (i.e., the assumed
ones performed in the “Vapor Flow Rate” section. The
endpoint pressure for the integral), and j is the incre-
main difference, however, is that the downstream phase
ment counter used for summation purposes. Combining
in a flashing liquid may take the form of a vapor, two-
Equations 1 and 2, we are able to obtain the mass flux
phase, or liquid, depending on its pressure. Thus, the
equation for liquid flow through a tube (Eq. 3).
vapor phase fraction at each pressure interval needs to be
calculated in addition to the total mass flux. The vapor
 phase fraction, yt , is expressed as a piecewise function

h v1 + vn + 2 × n−1
j=2 vj
as follows:
G2 = (3)
v2t

C1 × P(t) + C2 P(t) ≤ Ptbp
Once G is determined for different pressure intervals, yt = (6)
0 otherwise
the mass flux should be fitted against pressure by using
a quadratic fit in the following format (Eq. 4), where A1, Where C1 and C2 are the coefficients that give the best
A2, and A3 are the coefficients that give the best fit (e.g., fit, and Ptbp is the tube-side fluid’s bubble point pressure.
method of least squares). Smaller pressure intervals will The vapor phase fraction can then be multiplied with the
increase the accuracy of G. The end result is a continuous total mass flux, G, to yield the tube-side vapor mass flux,
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 7 of 21

Gv , and one minus the vapor phase fraction multiplied by Note that the equations used for flow through a pres-
G yields the tube-side liquid mass flux, Gl : sure relief valve can be modified to model spring loaded,
balanced bellows, and pilot valves. An example of this
performed is by Singh who developed a one-dimensional
Gv = G × y t (7) dynamic model for a spring loaded pressure safety valve
(PSV). [27] Moreover, any equations of flow through a
PSV provided by a vendor can also be substituted here.
Gl = G × (1 − yt ) (8) The expression can also be modified to reflect the reseat
pressure of a relief valve.
Density of Tube-Side Vapor, ρtv
With liquid-liquid systems, the density is assumed to Bulk Modulus of Elasticity
be constant. However, with vapor-liquid systems, density Interpolation should be used to determine the bulk mod-
must be changed as a function of pressure. If the vapor ulus for a fluid, B, shown in Eq. 14. While using the same
density is treated as an ideal gas, it would be inaccurate temperature as that of the system, the high and low side
at high pressures and low temperatures. If density ver- pressures (P1 and P0 ) and their respective volumes (V1 and
sus pressure is fitted using experimental data, the results V0 ) are chosen to calculate the bulk modulus. Although
will be accurate. The problem with this, however, is that the bulk modulus does vary, it is reasonable to assume it is
the number of gases available to deploy this method are constant throughout the entire range of pressures during
limited. A general way is to use an equation of state to the overpressure event.
determine the density versus pressure relationship. This
can be accomplished by using widely available process
simulation software (e.g., Aspen HYSYS). Density and P1 − P 0
B= V1 −V0
(14)
pressure are expressed via a linear relationship as follows: V0

Tube-Side Vapor and Liquid Volumes Present in Shell, Vtv

ρtv = B1 × P(t) + B2 (9)
and Vtl
Where B1 and B2 are the coefficients that give the best Vtv and Vtl represent the tube-side vapor and tube-side
fit. liquid volumes that are present in the exchanger shell.
At the beginning of a tube rupture, no tube-side fluid is
Rate of Tube-Side Fluid Entering Shell-Side
present in the shell. As fluid starts entering the shell-side,
With the tube-side mass flux, G, determined, the area Vtv and Vtl , are calculated by the following Equations:
of a single tube and the tube-side mass flow rate can be
obtained by Equations 10, 11, and 12. Note that the mass ˙ tv
dV ṁtv
flow rate is listed as two times the flux multiplied by = (15)
dt ρtv
the area. This is because the flow entering the shell-side
will be entering from upstream of one ruptured tube and
downstream of the other ruptured tube. ˙ tl
dV ṁtl
= (16)
dt ρtl
2
Atube = πrtube (10) Where ṁtv and ρtv are the tube-side vapor mass flow
rate (kg/s) and density (kg/m3 ), respectively. And ṁtl and
ρtl are the tube-side liquid mass flow rate (kg/s) and
m˙tl = 2Gl Atube (11) density (kg/m3 ), respectively.

m˙tv = 2Gv Atube (12) Shell Pressurization

A conservation of volume balance on the exchanger shell
Rate of Shell-Side Fluid Exiting Relief Device yields Eq. 17. [9] This allows one to determine the change
Assuming a pressure relief valve can be modeled as flow in pressure over time during a tube rupture. For liquid-
through an orifice, the following gives the mass flow rate liquid systems, the vapor associated terms ṁtv , ρtv , Vtv ,
exiting the relief device: and ctv0 are removed from the equation. Similarly, for
vapor-liquid systems, the terms ṁtl , ρtl , Vtl , and Btl are
  removed from the equation. In the case of a flashing tube-
Apsv C0 2ρsl gc P(t) P(t) ≥ Pset side liquid with liquid on the shell-side, all terms in Eq. 17
m˙psv = (13)
0 otherwise would be present.
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 8 of 21

  ṁtv
+ ṁtl ṁpsv
ρtl − ρsl
dP ρtv
for the liquid-liquid system solved in “Liquid-Liquid Sys-
= Vsl Vshell
(17) tem” section. The integral term in Table 2 refers to the
dt + Bsl + Bshell solution of Eq. 2. These properties are used to determine
⎧ Vtv
⎫ the tube-side mass flux which can be used to calcu-
⎪
⎪ vapor ⎪
⎪ late the mass flow rate entering the shell-side. The bulk
⎨ c2tv0 ρtv ⎬
Vtl
= liquid (18) modulus is also calculated at this stage. Once the above
⎪
⎪
Btl ⎪
⎩ Vtv
c2tv0 ρtv
+ Vtl
Btl flashing liquid ⎪
⎭ information are obtained, the equations relating the mass
flux and pressurization of the shell are used to calcu-
late the shell-side pressure. An appropriate time-step is
Where ṁtv and ρtv are the tube-side vapor mass flow selected to solve Eq. 17 numerically. If the time-step is
rate (kg/s) and density (kg/m3 ), respectively. ṁtl and ρtl too large, then the resolution will not be detailed enough.
are the tube-side liquid mass flow rate (kg/s) and den- If this occurs, the dynamic pressure profile obtained
sity (kg/m3 ), respectively. ṁpsv and ρsl are the PSV mass may not display the true peak of the tube rupture phe-
flow rate (kg/s) and shell-side liquid density (kg/m3 ), nomena. Too small a time-step, on the other hand,
respectively.  represents the tube-side volume and bulk increases the number of iterations needed to achieve
modulus terms. Vsl and Bsl are the volume of shell-side a solution.
liquid (m3 ) remaining in the shell-side and the shell-side
liquid’s bulk modulus (Pa), respectively. Vshell and Bshell Liquid-Liquid System
are the shell volume (m3 ) and the bulk modulus (Pa) Consider the following example. A carbon steel shell and
of the shell material of construction (e.g., carbon steel), tube exchanger has ethylene glycol in the tube-side and
respectively. The three cases for  are given by Eq. 18. cooling water in the shell-side. The exchanger contains
Vtv and ctv0 are the volume of tube-side vapor (m3 ) that 150 tubes each 2 meters long with an inner diameter of 15
has entered the shell-side and the tube-side vapor’s speed mm. The shell’s volume is 7.5 m3 . The high and low side
of sound (m/s), respectively. Vtl and Btl are the volume operating pressures are 10 bar and 1 bar, respectively. The
of tube-side liquid (m3 ) that has entered the shell-side high and low side operating temperatures are 100 °C and
and the tube-side liquid’s bulk modulus (Pa), respectively. 20 °C, respectively. The shell-side has a PSV set pressure of
This differential equation is the main mathematical model 1.2 bar (equal to its design pressure) and a hydrotest pres-
that will be solved and is valid for liquid-liquid, vapor- sure of 1.8 bar. The following steps are used to determine
liquid, and flashing liquid-liquid systems. The previous the required PSV size in order to adequately protect this
variables, described in earlier subsections all feed into this system.
model. Step 1: Obtain the physical properties of the liquid-
Since the initial shell-side pressure in Eq. 17 is known, liquid system.
this differential model can be viewed as an initial value Using the Peng-Robinson equation of state, at 1 bar and
problem. Furthermore, this equation is a first order 20 °C, the density of water is 1011 kg/m3 . Similarly, at
differential equation. A simple explicit method can be 10 bar and 100 °C, the density of ethylene glycol is 1055
used to solve the model. Euler’s method is used for kg/m3 . We then determine the bulk modulus for each
this paper. Using Matlab, a custom algorithm based fluid. It is determined that water has a bulk modulus of
on Euler’s method is developed [28]. The accuracy of 3.44931 ×109 Pa, and ethylene glycol has a bulk modulus
this method increases with a smaller step size. Because of 0.89769 ×109 Pa. The bulk modulus of common met-
tube ruptures can exhibit peak pressures in millisec- als can be found in engineering reference manuals [29].
onds, a step size of no larger than 0.1 millisecond is Carbon steel has a bulk modulus (Bshell ) of 159 ×109 Pa.
appropriate. Step 2: Use process simulation software to calculate
thermodynamic properties.
Results & Discussion Using Aspen HYSYS [30], a property table is created for
Figure 3 lists the major steps in generating a tube rup- the tube-side fluid. Isentropic flash calculations (in incre-
ture pressure profile. The first step in modeling a tube ments of 1 bar) are performed for ethylene glycol from the
rupture is collecting the appropriate data about the sys- tube-side pressure (10 bar) to the shell-side pressure (1
tem. These include, but are not limited to, equipment bar). The results are listed in Table 2.
diagrams, process simulation parameters (pressure, tem- Step 3: Calculate the mass flow rate through the tube,
perature, density, etc.), and relief valve datasheets. The entering the shell-side.
second step is to obtain the necessary thermodynamic Using the mass flux equation and property table, the
properties. One way to achieve this is by creating a table mass flux can be calculated (shown in Table 2). Once
of isentropic flashes from the tube-side conditions to the the mass flux for each pressure is determined, the data
shell-side conditions. Table 2 lists the isentropic flashes is then fitted to a quadratic equation in order to obtain
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 9 of 21

Fig. 3 Steps needed to generate dynamic pressure profile during tube rupture

the mass flux as a function of pressure. For the ethylene- the calculations presented shall assume a zero second
glycol water system described, Eq. 4 is generated with response time for relief valves.
coefficients A1, A2, A3 equal to -434.4, 526.4, and 41854.5, The pressure profiles for the liquid-liquid system are
respectively. Note that the pressure in Eq. 4 refers to the shown in Fig. 4. The pressure (bar) is plotted against time
shell-side pressure. At a pressure of 1 bar and 10 bar, (milliseconds) for various standard PSV orifice sizes. This
this equation yields the maximum and minimum tube- figure uses a time axis of 500 milliseconds, enough to cap-
side flow rates. This is reasonable since the largest driving ture the pseudo steady state pressures for all orifice sizes.
force exists when the shell-side pressure is 1 bar. Similarly, A dotted line representing the hydrotest pressure (1.8 bar)
when the shell-side is at 10 bar, no driving force exists. is shown. As stated earlier, we do not wish to exceed the
Step 4: Using an appropriate time-step, determine hydrotest pressure for this scenario. At time is equal to
the pressure at time t. 0 milliseconds, the tube rupture begins. This period in
At time is equal to zero, the shell-side pressure is set time has the largest difference in pressure between the
to 1 bar. Using Matlab [28], Equations 4, 9-13, and 15-16 two exchanger sides.For orifice sizes D, E, F, and G, the
are then simultaneously solved using a time-step of 1 mil- hydrotest pressure is exceeded at 25 milliseconds. An 10
lisecond. The required data are shown in Table 3. Note additional miliseconds shows orifice size H exceeding the
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 10 of 21

Table 2 Properties of ethylene-glycol isentropic flashes and (resulting from its smaller size), this PSV would most
mass flux at different pressure intervals likely be preferred in a plant.
Pressure Mass Density Specific Volume Integral Term Mass Flux
bar kg/m3 1e4 × m3 /kg m2 /s2 kg/s-m2 Vapor-Liquid System
10 1055 9.482 Consider the following example. A carbon steel shell and
9 1054 9.483 95 14522
tube exchanger has methane in the tube-side and cooling
water in the shell-side. The exchanger contains 100 tubes
8 1054 9.484 190 20536
each 3 meters long with an inner diameter of 10 mm. The
7 1054 9.485 285 25149 shell’s volume is 7.5 m3 . The high and low side operat-
6 1054 9.486 379 29037 ing pressures are 5 bar and 1 bar, respectively. The high
5 1054 9.487 474 32462 and low side operating temperatures are 100 °C and 20 °C,
4 1054 9.488 569 35558 respectively. The shell-side has a PSV set pressure of 1.2
3 1054 9.489 664 38404
bar (equal to its design pressure) and a hydrotest pressure
of 1.8 bar. We would like to determine the required PSV
2 1054 9.490 759 41052
size in order to adequately protect this system following
1 1054 9.491 854 43539 the steps:
Molar Entropy = -131 kJ/kgmole-°C Step 1: Obtain the physical properties of the vapor-
liquid system.
hydrotest pressure. These orifice sizes reach a pseudo Using the Peng-Robinson equation of state, at 1 bar and
steady state condition at approximately 300 miliseconds, 20 °C, the density of water is 1011 kg/m3 . Between 5 bar
with shell-side pressures ranging from 3 bar to above 9 and 1 bar, the density of methane in kg/m3 is given by
bar. Thus, for orifice sizes D, E, F, G, and H, the shell-side Eq. 9 with coefficients B1 and B2 equal to 0.4747 and 0.58,
failing is almost certain due to the hydrotest pressure respectively. As previously stated, water has a bulk modu-
being significantly exceeded. Orifice sizes J and K are the lus of 3.44931 ×109 Pa. Similarly, carbon steel has a bulk
only PSV sizes that adequately protect against this over- modulus (Bshell ) of 159 ×109 Pa.
pressure event. While orifice size J does exceed the design Step 2: Use process simulation software to calculate
pressure, it is comfortably below the hydrotest pressure. thermodynamic properties.
Orifice size K is even more conservative with shell-side Using Aspen HYSYS [30], a property table is created for
pressures not exceeding the design pressure during this the tube-side fluid. Isentropic flash calculations (in incre-
tube rupture. Due to the lower cost of the J orifice PSV ments of 0.4 bar) are performed for methane from the

Table 3 Parameters used in liquid-liquid, vapor-liquid, and flashing liquid-liquid systems

Ethylene-glycol & Water Methane & Water Propane & Water
Parameter Value Parameter Value Parameter Value
rtube 0.0075 m rtube 0.005 m rtube 0.005 m
ρtl 1055 kg/m3 ctv0 505.2 m/s ctv0 228 m/s
Btl 8.9769 ×108 Pa ρsl 1011 kg/m3 ρsl 1011 kg/m3
ρtl 1011 kg/m3 Bsl 3.4493 x 109 Pa Bsl 3.4493 x 109 Pa
Bsl 3.4493 ×109 Pa Vsl 7.5 m3 Vsl 7.5 m3
Vshell 7.5 m3 Vshell 7.5 m3 Vshell 7.5 m3
Vsl 7.5 m3 Bshell 159 x 109 Pa Bshell 159 x 109 Pa
Bshell 159 ×109 Pa α1 -34.219 α1 -8.131
A1 -434.4 α2 219.62 α2 323.33
A2 526.4 α3 -439.53 α3 -3295.7
A3 41854.5 α4 997.29 α4 27649
B1 0.4747 B1 2.32
B2 0.58 B2 -1.5468
C1 -0.025
C2 0.5285
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 11 of 21

Fig. 4 Shell-side pressure profiles for ethylene-glycol water systems for different PSV sizes

tube-side pressure (5 bar) to the shell-side pressure (1 bar). the calculations presented shall assume a zero second
The results are listed in Table 4. response time for relief valves.
Step 3: Calculate the mass flow rate through the tube, Figure 5 shows the mass flux profile for methane at
entering the shell-side. various downstream pressures. Two mass flux profiles
Using the mass flux equation and property table, for methane are included, one which takes into account
the mass flux can be calculated and is given in choked flow (listed as red circles) and one that excludes
Table 4. Note that the methane will exhibit choked the effects of choked flow (listed as blue triangles). The
flow, and the mass flux must reflect this (see Fig. 5). choked flow region is roughly between downstream pres-
Thus, the mass flux versus pressure data is fitted by sures of 1 bar and 2.6 bar. In other words, from an initial
a cubic polynomial (Eq. 5) with the coefficients α1, pressure and temperature of 5 bar and 100 °C, the mass
α2, α3, and α4 equal to -34.219, 219.62, -439.53, and flow rate of methane transported through an isentropic
997.29, respectively. Between 1 bar and 5 bar, this nozzle (as is the case in a tube rupture) increases with
equation yields the maximum and minimum tube-side decreasing downstream pressure until 2.6 bar. Between
flow rates. downstream pressures of 1 bar and 2.6 bar, the mass
Step 4: Using an appropriate time-step, determine flow rate of methane is constant and at its maximum
the pressure at time t. (731 kg/s-m2 ). Neglecting the effects of choked flow can
At time is equal to zero, the shell-side pressure is set result in a misleading mass flux, as can be seen by the
to 1 bar. Using Matlab [28], Equations 5, 9-14, and 16 are negative parabolic shaped trendline in Fig. 5. Accepting
then simultaneously solved using a time-step of 1 mil- this erroneous mass flux may lead to an undersized relief
lisecond. The required data are shown in Table 3. Note system.
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 12 of 21

Table 4 Properties of methane isentropic flashes and mass flux at different pressure intervals
Pressure Mass Density Specific Volume Integral Term Mass Flux Corrected Mass Flux
bar kg/m3 m3 /kg m2 /s2 kg/s-m2 kg/s-m2
5.0 2.598 0.3849
4.6 2.433 0.4110 7698 302 302
4.2 2.265 0.4414 25287 509 509
3.8 2.095 0.4774 43663 619 619
3.4 1.921 0.5207 63625 685 685
3.0 1.742 0.5740 85518 721 721
2.6 1.559 0.6414 109826 731 731
2.2 1.370 0.7299 137252 718 731
1.8 1.174 0.8518 168886 682 731
1.4 0.968 1.0328 206578 622 731
1.0 0.749 1.3348 253931 534 731
Molar Entropy = 178.5 kJ/kgmole-°C

The pressure profiles for the methane water system are pressure (1.8 bar) is included. For this tube rupture, orifice
shown in Fig. 6. By plotting the pressure (bar) versus size Q appears to provide adequate overpressure protec-
time (milliseconds), the effects of selecting different PSV tion, without being conservative. One notable difference
orifice sizes (D, J, N, P, Q, R, and T) to protect this sys- between this system and the liquid-liquid system is the
tem are shown. A dotted line representing the hydrotest speed at which pseudo steady state pressure is reached.

Fig. 5 Tube-side mass flux vs. downstream pressure for methane with an initial pressure of 5 bar and initial temperature of 100 °C
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 13 of 21

Fig. 6 Shell-side pressure profiles for methane-water systems for different PSV sizes

Figure 6 took approximately 30 milliseconds as opposed These pressure oscillations can be broken into a series of
to 300 milliseconds in Fig. 4. This highlights the severity steps. First, the system pressure builds up due to an over-
of a tube rupture in vapor-liquid systems. pressure event, causing the PSV to open. Because the PSV
is sufficiently large to handle the influx of fluid entering,
Comparison of Liquid-Liquid and Vapor-Liquid Systems the PSV vents all of the excess fluid along with additional
Assuming the same tube and shell-side pressures, for an fluid. At this point, the system pressure decreases, even-
equivalent mass entering the shell-side, a vapor-liquid sys- tually dropping below the set pressure. Since the system
tem reaches peak pressure significantly faster (and as a pressure is below the set pressure, the PSV reseats to
result is more severe) than that of a liquid-liquid sys- its original position. This cycle repeats for the duration
tem. This is true even though the fluid flow (on a mass of the overpressure event. More sophisticated modeling
basis) that enters the shell-side for vapors is significantly can include the effects of the PSV reseating. The rapid
less than that of liquids. The vapor-liquid simulations opening and closing of the relief valve is known as chatter-
indicate that the peak pressures are reached in less than ing and has been documented as affecting the reliability
50 ms as opposed to the liquid-liquid cases which took of a PSV. [31] Chattering can occur due to oversizing a
between 100 and 300 ms, depending on the size of the relief valve, which appears to be the case for orifices T
PSV orifice. Shell-side transient peak pressures are able to and R. For more common overpressure scenarios such
significantly exceed their hydrotest pressures. For systems as blocked outlet, control valve failures, and fires, plants
left with inadequate overpressure protection, these simu- may desire to select a relief valve that avoids chattering.
lations demonstrate the consequences of tube ruptures. However, because of the low frequency of tube ruptures
The T and R orifices, as shown in Fig. 6, resulted in in plants, chattering is not considered to be a primary
the pressure oscillating above and below the set pressure. concern.
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 14 of 21

Another distinction between vapor-liquid and liquid- Using the Peng-Robinson equation of state, the density
liquid systems is the presence of choked flow. For vapor- of 60 °C liquid propane between 30 bar and 6 bar is on
liquid systems, ignoring choked flow can result in a relief average 446 kg/m3 . The density of vapor propane between
system being undersized. Figure 6 shows the tube-side 21 bar and 6 bar is given by Eq. 9 with coefficients B1
mass flux when including and excluding choked flow. For and B2 equal to 2.32 and -1.5468, respectively. As pre-
downstream pressures that are close to the upstream pres- viously stated, water has a density of 1011 kg/m3 and a
sures, it is possible the fluid may not enter the choked bulk modulus of 3.44931 ×109 Pa. Liquid propane has a
region. This would yield the same mass flow rate whether bulk modulus of 0.1536 ×109 Pa. Carbon steel has a bulk
choked flow is included or excluded. This is seen in the modulus (Bshell ) of 159 ×109 Pa.
pressure interval between 5 bar and 2.6 bar in both Table 4 Step 2: Use process simulation software to calculate
and Fig. 5. This should not be left to chance however. thermodynamic properties.
Testing for choked flow should always be performed. Using Aspen HYSYS [30], a property table is created
for the tube-side fluid. Isentropic flash calculations (in
Flashing Liquid-Liquid System increments of 1.5 bar) are performed for propane from
Consider the following example. A carbon steel shell and the tube-side pressure (30 bar) to the shell-side pressure
tube exchanger has propane in the tube-side and cooling (6 bar). Because vapor is generated, the vapor fraction is
water in the shell-side. The exchanger contains 100 tubes included in the property table. The results are listed in
each 4 meters long with an inner diameter of 10 mm. The Table 5.
shell’s volume is 7.5 m3 . The high and low side operat- Step 3: Calculate the mass flow rate through the tube,
ing pressures are 30 bar and 6 bar, respectively. The high entering the shell-side.
and low side operating temperatures are 60 °C and 6 °C, Using the mass flux equation and property table, the
respectively. The shell-side has a PSV set pressure of 7.2 mass flux can be calculated and is given in Table 5.
bar (equal to its design pressure) and a hydrotest pressure Methane exhibits choked flow with a maximum mass flux
of 10.8 bar. We would like to determine the required PSV of 27170 kg/s-m2 . The corrected mass flux versus pressure
size in order to adequately protect this system following is fitted against a cubic polynomial (Eq. 5) with the coeffi-
the steps: cients α1, α2, α3, and α4 equal to -8.131, 323.33, -3295.7,
Step 1: Obtain the physical properties of the flashing and 27649, respectively. The vapor fraction is expressed
liquid-liquid system. as Eq. 6 with coefficients C1 and C2 equal to -0.025 and

Table 5 Properties of propane isentropic flashes and mass flux at different pressure intervals
Pressure Mass Density Integral Term Mass Flux Corrected Mass Flux Vapor Fraction
bar kg/m3 m2 /s2 kg/s-m2 kg/s-m2
30.0 434.9
28.5 434.5 172 8069 8069 0
27.0 434.1 589 14900 14900 0
25.5 433.6 935 18751 18751 0
24.0 433.2 1281 21928 21928 0
22.5 432.8 1627 24692 24692 0
21.0 432.4 1974 27170 27170 0
19.5 330.7 2374 22786 27170 0.0367
18.0 248.5 2903 18936 27170 0.0791
16.5 192.2 3595 16294 27170 0.1194
15.0 151.0 4482 14294 27170 0.1581
13.5 119.6 5606 12662 27170 0.1954
12.0 94.8 7024 11242 27170 0.2318
10.5 74.9 8816 9943 27170 0.2675
9.0 58.4 11102 8707 27170 0.3029
7.5 44.7 14064 7493 27170 0.3384
6.0 33.0 18014 6267 27170 0.3747
Molar Entropy = 104.5 kJ/kgmole-°C
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 15 of 21

0.5285, respectively. The vapor fraction equation is only mass flux is present between the downstream pressures 21
valid between pressures 6 bar and 21 bar (the tube-side bar and 30 bar.
bubble point pressure). The pressure profiles for the propane water system are
Step 4: Using an appropriate time-step, determine shown in Fig. 8. Compared to the previous two cases,
the pressure at time t. this system requires roughly between 100 to 200 mil-
At time is equal to zero, the shell-side pressure is set to liseconds before reaching pseudo steady state. A dot-
6 bar. Using Matlab [28], Equations 6-16, and 16 are then ted line line at 10.8 bar represents the hydrotest pres-
simultaneously solved using a time-step of 1 millisecond. sure. The smallest PSV orifice size that does not exceed
The required data are shown in Table 3. Note the calcula- this pressure is K. One interesting feature of this sys-
tions presented shall assume a zero second response time tem is that it functions as both a liquid-liquid and
for relief valves. vapor-liquid tube rupture. In other words, as the tube
Figure 7 describes the vapor and liquid mass fluxes rupture begins, the high pressure liquid flashes when
for propane with changes in pressure. The x-axis is the exposed to the shell-side. This phase of the tube rup-
downstream pressure in bar. This downstream pressure is ture is closely associated with a vapor-liquid rupture.
equivalent to the shell-side pressure. Thus, Fig. 7 can be Over time, pressure builds up and the amount of flash-
interpreted as beginning from the left of the graph (at a ing decreases. Eventually, the liquid exiting the tube-side
downstream pressure of 6 bar), and shifting to the right will fully remain a liquid in the shell-side, exhibiting char-
as the tube rupture develops. If no PSV is present, this acteristics similar to that of a liquid-liquid rupture. Ori-
would continue until the downstream pressure is in equi- fice sizes D and E are both piece-wise functions since
librium with the upstream pressure (30 bar). In addition to around 20 bar, the shell-side is entirely a liquid due to the
the x-axis, two y-axis are shown, one being the mass flux high pressure.
(kg/s-m2 ) and the other being the vapor fraction. Between
downstream pressures of 21 bar and 30 bar, the vapor frac- Dimensionless Analysis of Tube Rupture
tion is zero. Between pressures 6 bar and 21 bar, the vapor
fraction increases with decreasing downstream pressure. The severity of a tube rupture can be affected by a myr-
The vapor mass flux is generated by multiplying the total iad of properties. These include the fluid’s density, phase,
mass flux with the vapor fraction. Therefore, no vapor pressure, and temperature. Exchanger properties such as

Fig. 7 Tube-side mass flux vs. downstream pressure for propane with an initial pressure of 30 bar and initial temperature of 60 °C
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 16 of 21

Fig. 8 Shell-side pressure profiles for propane-water systems for different PSV sizes

the shell volume and tube diameter affect the flux and severity of a tube rupture, while a larger shell volume slows
the rate of pressure increase. Lastly, pressure relief prop- the accumulation of pressure. For Eq. 19, characteristic
erties (e.g. valve size and opening times) affect the rate lengths of 1e-4 m and 1 m3 for Dc and Vc were selected.
of fluid exiting the shell. To better understand the differ- The pressure ratio (Eq. 20) is investigated because it is
ent parameters that can affect a tube rupture, properties most likely the primary driver in the severity of a tube
were grouped together to form dimensionless numbers. rupture. In fact, API 521 recommends a dynamic sim-
The following numbers are proposed: ulation when pressure differences exceed 1000 psi. The
dimensionless number Pmax in Eq. 21 is the maximum
pressure the shell-side experiences during a tube rupture
Dtube /Dc
τ= (19) divided by the shell design pressure. Lastly, the effects of
Vshell /Vc the density ratio’s (Eq. 22) impact on tube rupture will be
studied. A series of liquid-liquid tube rupture simulations
Ptube
Pratio = (20) were performed. In addition to the dimensionless num-
Pshell bers proposed, these simulations also measure the time
P(t)shell max it takes to reach the maximum pressure. The simulations
Pmax = (21)
Pshell design contained ethylene-glycol and water on the tube and shell-
ρtube side, respectively. In Fig. 9, for a constant pressure ratio of
ρratio = (22) 10, the maximum pressure was plotted against the density
ρshell
ratio (which was adjusted from 0.1 to 10). The dimension-
less number, τ , was set at 13.3, 20, and 40. These values for
τ in Eq. 19 evaluates the size of the tube in relation to τ are obtained from a 7.5 m3 shell with 1 cm, 1.5 cm, and
the exchanger shell. A larger tube diameter increases the 3 cm diameter tubes, respectively.
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 17 of 21

Fig. 9 Shell-side peak pressure vs. density ratio for varying τ with a pressure ratio of 10

From Fig. 9, it is seen that the density ratio can signif- capturing the severity of a tube rupture. A scenario that
icantly affect the ability of an exchanger to pressurize. A reaches peak pressure quickly may result in the shell-side
high density ratio reduces the ability of a tube rupture to failing due to the PSV not opening in time. This is in
result in overpressure. There exists for each system a den- contrast to a longer time to reach Pmax , which allows the
sity ratio threshold, in which no increase in pressure is PSV to open, limiting our concern to the maximum pres-
present for this overpressure scenario. This relationship sure the system reaches. The same pressure ratio (held
can be explained by thinking of an infinitely large density constant at 10), density ratios, and values for τ as Fig. 9
ratio. This would mean that the tube-side fluid is infinitely were used. The results are shown in Fig. 10. All else being
small and cannot pressurize the system. In mathematical equal, it is preferred to have a tube rupture with a longer
terms, an infinite tube-side density forces the dP dt term in time to reach maximum pressure. For an equivalent den-
Eq. 17 to zero. Alternatively, low density ratios were able to sity ratio, a faster time to reach Pmax is observed for higher
reach 10 bar regardless of the tube diameter and shell vol- τ values. For a fixed τ , a higher density ratio increases the
ume. When looking at the effects of τ , for the maximum time it takes for a system to reach its maximum pressure.
pressure a system experiences, a low value for τ results in In addition, for all values of τ , the maximum pressure is
a lower peak pressure. This is due to the fact that a low achieved in less than 350 milliseconds.
τ value results in a small tube diameter in relation to the The third set of dimensionless analysis simulations
shell volume. Thus, the cross sectional area of flow will be (Fig. 11) study the impact the pressure ratio has on the
less, allowing a pressure relief device to mitigate the tube severity of a tube rupture. The pressure ratio was varied,
rupture’s impact. in increments of 2, from 2 bar to 10 bar. Two values of τ
A second set of simulations were performed with the were used, 20 and 40. Similarly, three values for the den-
aim of studying the time it takes to reach the system peak sity ratio were used, 0.49, 1.09, and 5.04. From Fig. 11, a
pressure, Pmax . This parameter is important to accurately higher pressure ratio leads to a higher maximum pressure.
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 18 of 21

Fig. 10 Time from initial rupture to peak pressure vs. density ratio for varying τ with a pressure ratio of 10

This is not surprising. What is interesting to note, how- upgraded to comply with the two-thirds rule, then a pres-
ever, is how strong the density ratio can affect this peak sure relief valve would not be required for that system.
pressure. As an example, a pressure ratio of 10 with a τ Another strategy to mitigate tube ruptures is preventing
of 20, can yield a maximum pressure of 1.2 bar, 3 bar, or exchanger tubes from failing in the first place. Nowadays,
7 bar. This wide variation in maximum pressure is due to nondestructive technology such as eddy current testing
the systems having density ratios of 5.04, 1.09, and 0.49, can help estimate the remaining life of exchanger tubes.
respectively. An increase in τ also shows an increase in a [32] This technology can be used to more frequently mon-
tube rupture’s maximum pressure, although not as severe itor exchangers that are more susceptible to erosion and
as the density ratio. A pressure ratio of 8, with a density corrosion. Proactive prevention can also focus on detect-
ratio of 0.49 yields a maximum pressure of 5 bar (τ equal ing future operational hazards, such as by incorporating
to 20) and 8 bar (τ equal to 40). alarms within a process model. [33]
One can also select an appropriately sized relief valve
Mitigating Tube Rupture Scenarios to handle a tube rupture scenario. For the ethylene-glycol
Depending on how far along in the lifecycle of an water example in “Liquid-Liquid System” section, without
exchanger a plant is, mitigating the severity of tube rup- a pressure relief valve, the shell-side would be pressur-
ture incidents can be done through one of several ways. ized from 1 bar to 9 bar in 300 ms. However, by installing
First, one should flag exchangers that will have large differ- either relief valve sizes J or K, the shell-side would exhibit
ences in pressure. These exchangers are most susceptible peak pressures of 1.43 and 1.27 bar, respectively. More-
to failure in the event of a tube rupture. For exchangers over, while relief size K exhibits chattering during this tube
that are still in the design phase, increasing the strength rupture, relief size J does not. Furthermore, the results
of the exchanger is one way of mitigating the effects of a from this paper indicate vapor-liquid and flashing liquid-
tube rupture. This design change does increase the invest- liquid ruptures to be more severe than liquid-liquid tube
ment cost of the exchanger, however. If the exchanger is ruptures. Thus, more resources (e.g., increased design
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 19 of 21

Fig. 11 Shell-side peak pressure vs. pressure ratio relationship for varying τ and density ratios

pressure and larger PSV orifice sizes) can be dedicated to paper covered the main steps in generating a pressure
these higher risk exchangers. Lastly, it was found that the profile, and a step-by-step calculation was performed for
density ratio and pressure ratio can significantly affect the an ethylene-glycol water, methane water, and propane
ability of a shell to pressurize. As a result, systems with a water system. A comparison of the PSV sizes reveals that
high density ratio and low pressure ratio are preferred for vapor-liquid and liquid-liquid systems are the most and
this overpressure scenario. least severe cases, respectively. Flashing liquid-liquid tube
rupture scenarios performed in between these two cases.
Conclusion The examples covered in this paper can serve as a basis
Shell and tube heat exchangers are commonly used in for approaching liquid-liquid, vapor-liquid, and flashing
the oil, gas, chemical, and nuclear industries. Although liquid-liquid systems.
rare, tube rupture overpressure events may compromise
the mechanical integrity of an exchanger and can lead to Abbreviations
the equipment’s failure. This has the potential to result Apsv : Relief valve orifice area; Bshell : Shell material of construction bulk modulus;
in catastrophic failures and should be modeled with rig- Bsl : Shell-side liquid bulk modulus; Btl : Tube-side liquid bulk modulus; C0 :
Coefficient of discharge; Dc : Characteristic tube diameter; Dtube : Tube
orous sizing methods. This paper points out the chal- diameter; P(t): Transient shell pressure; Pmax,norm : Normalized maximum
lenges in modeling tube ruptures. The importance of transient shell-side pressure; Pratio : Dimensionless number relating tube and
accurately modeling tube rupture scenarios increases with shell pressures; rtube : Tube radius; Vc : Characteristic shell volume; Vshell : Shell
volume; Vtl : Volume of tube-side liquid present in shell-side; Vtv : Volume of
large differences in pressure between the shell and tube. tube-side vapor present in shell-side; ṁpsv : Relief valve mass flow rate; ṁtl :
In addition, low tube-to-shell density ratios increase the Mass flow rate of tube-side liquid entering shell-side; ṁtv : Mass flow rate of
severity of this event. By switching from an orifice-style tube-side vapor entering shell-side; ρratio : Dimensionless number relating tube
and shell densities; ρsl : Shell-side liquid density; ρtl : Tube-side liquid density;
calculation into a more rigorous dynamic simulation, one ρtv : Tube-side vapor density; τ : Dimensionless number relating tube size to
can better predict the effects of a tube rupture. This shell volume; MAWP: Maximum allowable working pressure; PSV: Pressure
Harhara and Hasan BMC Chemical Engineering (2020) 2:5 Page 20 of 21

safety valve or pressure relief valve; G: Mass flux of tube-side fluid entering 16. Pipeline Simulation and Integrity Ltd. Testing and analysis of relief device
shell-side; h: Pressure step size; set: Shell-side set pressure opening times. Tech Rep, Health Safety Executive. 2002.
17. Schwartz MP. Four types of heat exchanger failures. Plant Eng. 1982;23:
Acknowledgements 45–50.
The authors gratefully acknowledge support from the Mary Kay O’Connor 18. Megens B. Tube rupture study of a 300 bar heat exchanger. In: ASME 2014
Process Safety Center, and the Artie McFerrin Department of Chemical Pressure Vessels and Piping Conference. New York: American Society of
Engineering at Texas A&M University. Mechanical Engineers Digital Collection; 2014.
19. Shahrani S, Al-Subai S. Failure analysis of heat exchanger tubes. J Fail Anal
Authors’ Contributions Prevent. 2014;14(6):790–800.
AH and MMFH wrote the manuscript. AH developed the method with support 20. Njobuenwu DO, Fairweather M. Modelling of pipe bend erosion by dilute
from MMFH. MMFH supervised the project. All authors read and approved the particle suspensions. Comput Chem Engineer. 2012;42:235–47.
final manuscript. 21. Khilnaney V. Heat exchanger vibrations - a case study (Paper No 512).
National symposium on commissioning and operating experiences in
Funding
heavy water plants and associated chemical industries. 1992;24(16):501.
This research was funded by the Mary Kay O’Connor Process Safety Center at
India: Bhabha Atomic Research Centre.
Texas A&M University.
22. Boiler A, Code PV. Section viii division 1. UG-126 Pressure Relief Valves to
Availability of data and materials UG-129 Marking, ASME International, New York. 2010.
The datasets generated and analysed during the current study are available 23. Stewart M, Lewis OT. Heat Exchanger Equipment Field Manual: Common
from the corresponding author on reasonable request. Operating Problems and Practical Solutions. Houston: Gulf Professional
Publishing; 2012.
Competing interests 24. American Institute of Chemical Engineers, Center for Chemical Process
The authors declare that they have no competing interests. Safety. Guidelines for Initiating Events and Independent Protection Layers
in Layer of Protection Analysis. Hoboken: John Wiley & Sons; 2014.
Received: 25 September 2019 Revised: 22 January 2020 Accepted: 23 25. Deddis C. Shell & tube heat exchanger overpressure protection from tube
January 2020 rupture. Tech Rep. 2011.
26. Standard 520: Sizing API Selection. Part I — Sizing and Selection. 2008.
27. Singh A. An analytical study of the dynamics and stability of a spring
loaded safety valve. Nuclear Engineer Design. 1982;72(2):197–204.
References
28. MathWorks. MATLAB 9.2. MathWorks. 2017.
1. Hellemans M. The Safety Relief Valve Handbook: Design and Use of
29. Krauss G, Hudok D. Properties and selection: Iron, steel and high
Process Safety Valves to ASME and International Codes and Standards.
performance alloys. ASM Metals Hand Book, Materials. 1991;1:.
Oxford: Elsevier; 2009.
30. AspenTech. Aspen HYSYS V10. AspenTech. 2017.
2. American Petroleum Institute. API Standard 521: Pressure-relieving and
31. Hagey E. The psv that did not fail—misconceptions about psvs. Process
Depressuring Systems, 6th ed. , Washington, DC: American Petroleum
Safety Progress. 2013;32(1):84–9.
Institute; 2016.
32. Golovin V, Pechnikov N, Shchelkov V, Tsivadze AY. Determination of the
3. Feng J, Aggarwal A, Dasgupta S, Shariat H. Using dynamic analysis to
life cycle of heat-exchange tubes of vapor condensers on the basis of
reduce weight of offshore installations. 2009.
statistical analysis of local pitting corrosion according to data of eddy
4. Ahammad M, Hatanaka L, Hui J, Lapeyrouse K, Roche E. Learnings from
current testing. Protect Metals Phys Chem Surfaces. 2018;54(6):1221–32.
mary kay o’connor process safety center (mkopsc) instrument reliability
33. Ahooyi TM, Soroush M, Arbogast JE, Seider WD, Oktem UG.
network’s project on pressure transmitter maintenance data collection.
Model-predictive safety system for proactive detection of operation
In: Annu Int Symp; 2016. p. 523–32.
hazards. AIChE J. 2016;62(6):2024–42.
5. United Kingdom’s Health and Safety Executive Offshore Hydrocarbon
Releases. 2016.
6. Fowler D, Herndon T, Wahrmund R. An analysis of potential overpressure Publisher’s Note
of a heat exchanger shell due to a ruptured tube. In: ASME Petroleum Springer Nature remains neutral with regard to jurisdictional claims in
Mechanical Engineering Conference, Dallas TX. New York: American published maps and institutional affiliations.
Society of Mechanical Engineers; 1968.
7. Simpson LL. Dynamics of underwater gas tubing failure. Master’s thesis,
West Virginia University. 1971.
8. Sumaria V, Rovnak J, Heitner I, Herbert R. Model to predict transient
consequences of a heat exchanger tube rupture. Proceedings of the
American Petroleum Institute Division of Refining. 1976;55(41):631–654.
9. Ennis C, Botros K, Patel C. Dynamic model for a heat exchanger tube
rupture discharging a high-pressure flashing liquid into a low-pressure
liquid-filled shell. J Loss Prevent Proc Industries. 2011;24(1):111–21.
10. Botros KK. Importance of accounting for the piping system and boundary
conditions in determining the maximum surge pressure following
heat-exchanger tube-rupture. J Loss Prevent Proc Industries. 2015;37:
63–73.
11. Cassata J, Feng Z, Dasgupta S, Samways R. Prevent overpressure failures
on heat exchangers. Hydrocarbon Process. 1998;77(11):123–8.
12. Ewan B, Moatamedi M. Design aspects of chemical plant exposed to
transient pressure loads. Chem Engineer Res Design. 2000;78(6):866–70.
13. Nagpal S. Evaluate heat-exchanger tube-rupture scenarios using dynamic
simulation. Chem Engineer. 2015;122(2):48.
14. Acosta C, Siu N. Dynamic event trees in accident sequence analysis:
application to steam generator tube rupture. Reliability Engineer Syst
Safety. 1993;41(2):135–54.
15. Thyer A, Wilday A, Bankes G. The experimental study and simulation of
tube rupture in shell-and-tube heat exchangers. In: Institution of
Chemical Engineers Symposium Series. Marylebone: Institution of
Chemical Engineers; 1999; 2000.