# HydDown - User guide and technical reference (v0.30.0)

Anders Andreasen

## ► To cite this version:

Anders Andreasen. HydDown - User guide and technical reference (v0.30.0). ORS Consulting. 2025. [https://hal.science/hal-04858235v3](https://hal.science/hal-04858235v3)

HAL Id: hal-04858235

[https://hal.science/hal-04858235v3](https://hal.science/hal-04858235v3)

Submitted on 10 Dec 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# HydDown

User guide and technical reference

Anders Andreasen
# Contents

<table><thead><tr><th></th><th></th><th></th></tr></thead><tbody><tr><td>1</td><td>Introduction</td><td>1</td></tr><tr><td>1.1</td><td>Citing HydDown</td><td>2</td></tr><tr><td>1.2</td><td>Background</td><td>2</td></tr><tr><td>1.3</td><td>Getting the software</td><td>3</td></tr><tr><td>1.4</td><td>Requirements</td><td>4</td></tr><tr><td>1.5</td><td>Testing</td><td>4</td></tr><tr><td>1.6</td><td>Units of measure</td><td>5</td></tr><tr><td>1.7</td><td>Credit</td><td>6</td></tr><tr><td>1.8</td><td>License</td><td>6</td></tr><tr><td>2</td><td>Usage</td><td>7</td></tr><tr><td>2.1</td><td>Basic usage</td><td>7</td></tr><tr><td>2.2</td><td>Demos</td><td>7</td></tr><tr><td>2.3</td><td>Calculation methods</td><td>8</td></tr><tr><td>2.4</td><td>Script</td><td>9</td></tr><tr><td>2.5</td><td>Module import</td><td>10</td></tr><tr><td>2.6</td><td>Input file examples</td><td>10</td></tr><tr><td>2.7</td><td>Input fields and hierarchy</td><td>11</td></tr><tr><td>2.7.1</td><td>Calculation</td><td>12</td></tr><tr><td>2.7.2</td><td>Vessel</td><td>12</td></tr><tr><td>2.7.3</td><td>Initial</td><td>12</td></tr><tr><td>2.7.4</td><td>Valve</td><td>13</td></tr><tr><td>2.7.5</td><td>Heat transfer</td><td>13</td></tr><tr><td>2.7.6</td><td>Validation</td><td>14</td></tr><tr><td>3</td><td>Theory</td><td>17</td></tr><tr><td>3.1</td><td>Thermodynamics</td><td>17</td></tr><tr><td>3.1.1</td><td>Equation of state</td><td>17</td></tr><tr><td>3.1.2</td><td>First law for flow process</td><td>19</td></tr><tr><td>3.2</td><td>Flow devices</td><td>22</td></tr><tr><td>3.2.1</td><td>Restriction Orifice</td><td>22</td></tr></tbody></table>
<table><tbody><tr><td>3.2.2</td><td>Pressure safety valve / Relief valve</td><td>23</td></tr><tr><td>3.2.3</td><td>Control Valve</td><td>26</td></tr><tr><td>3.3</td><td>Heat transfer</td><td>28</td></tr><tr><td>3.3.1</td><td>Natural convection</td><td>28</td></tr><tr><td>3.3.2</td><td>Mixed convection</td><td>29</td></tr><tr><td>3.3.3</td><td>Nucleate boiling heat transfer</td><td>30</td></tr><tr><td>3.3.4</td><td>Conduction</td><td>30</td></tr><tr><td>3.3.5</td><td>Fire heat loads</td><td>32</td></tr><tr><td>3.4</td><td>Vessel geometry</td><td>33</td></tr><tr><td>3.5</td><td>Model implementation</td><td>34</td></tr><tr><td>3.5.1</td><td>Isothermal process</td><td>35</td></tr><tr><td>3.5.2</td><td>Isenthalpic process</td><td>35</td></tr><tr><td>3.5.3</td><td>Isentropic process</td><td>35</td></tr><tr><td>3.5.4</td><td>Isenergetic process</td><td>36</td></tr><tr><td>3.5.5</td><td>Energy balance</td><td>36</td></tr><tr><td>3.6</td><td>Multicomponent mixtures</td><td>37</td></tr><tr><td><b>4</b></td><td><b>Validation</b></td><td><b>41</b></td></tr><tr><td>4.1</td><td>Hydrogen discharge (Type I)</td><td>41</td></tr><tr><td>4.2</td><td>Hydrogen filling (Type I)</td><td>44</td></tr><tr><td>4.3</td><td>Nitrogen discharge (Type I)</td><td>47</td></tr><tr><td>4.4</td><td>Multi-component discharge (Type I)</td><td>48</td></tr><tr><td>4.5</td><td>1-D transient heat transfer</td><td>50</td></tr><tr><td>4.5.1</td><td>Validation against commercial simulation tool</td><td>50</td></tr><tr><td>4.5.2</td><td>Validation against KIT experiment (Type IV)</td><td>51</td></tr><tr><td>4.6</td><td>Validation against GasTeF experiments (Type IV)</td><td>53</td></tr><tr><td>4.7</td><td>Validation against Type III cylinder discharge experiments</td><td>55</td></tr><tr><td>4.8</td><td>Validation against Type III cylinder filling experiments</td><td>60</td></tr><tr><td><b>5</b></td><td><b>Similar software</b></td><td><b>63</b></td></tr><tr><td><b>6</b></td><td><b>Future development of HydDown</b></td><td><b>65</b></td></tr><tr><td><b>7</b></td><td><b>Example use cases</b></td><td><b>67</b></td></tr><tr><td></td><td><b>References</b></td><td><b>69</b></td></tr></tbody></table>
# 1 Introduction

HydDown is an open source python3 tool for calculation of hydrogen (or other pure component) vessel/container depressurization and filling. The HydDown logo shown in Fig. 1.1 visualizes the key parameters and transport phenomena during gas vessel filling or discharging. The thermodynamic state inside the vessel changes over time as seen from immediately observable variables temperature (T) and pressure (P). This is caused by change in fluid inventory (density) due to flow of gas either in or out of the vessel. Further, heat is transferred from or to the surroundings via convective heat transfer on the in- and outside of the vessel with heat being conducted through the vessel wall.

<figure><img src="image_2.png" /><figcaption>**Figure 1.1:** HydDown logo</figcaption></figure>

Running the code is as simple as:

1 python main.py input.yml

where `main.py` is the main script and `input.yml` is the input file in Yaml syntax.
## 1.1 Citing HydDown

If you use HydDown please cite the following reference (Andreasen 2021):

Andreasen, A., (2021). HydDown: A Python package for calculation of hydrogen (or other gas) pressure vessel filling and discharge. Journal of Open Source Software, 6(66), 3695, https://doi.org/10.21105/joss.03695

1 @article{Andreasen2021,
2 doi = {10.21105/joss.03695},
3 url = {https://doi.org/10.21105/joss.03695},
4 year = {2021},
5 publisher = {The Open Journal},
6 volume = {6},
7 number = {66},
8 pages = {3695},
9 author = {Anders Andreasen},
10 title = {HydDown: A Python package for calculation of hydrogen (or other gas) pressure vessel filling and discharge},
11 journal = {Journal of Open Source Software}
12 }

This manual can be cited as (Andreasen 2024):

Anders Andreasen. HydDown - User guide and technical reference. 2024. (hal-04858235)

1 @report{Andreasen2024,
2 url = {https://hal.science/hal-04858235},
3 year = {2024},
4 publisher = {HAL open science},
5 author = {Anders Andreasen},
6 title = {HydDown -- user guide and technical reference},
7 }

## 1.2 Background

HydDown started as a small spare-time project for calculation of vessel filling and depressurization behaviour. This was mainly to demonstrate that although perceived as a very tedious and difficult task to write your own code for such an apparently complex problem, actually only a fairly limited amount of code is necessary if you have a good thermodynamic backend. The code has evolved to a point where additional features will increase the complexity to the point where it is no longer a simple tool.

A few choices has been made to keep things simple:

* Coolprop is used as thermodynamic backend
* Only pure substances are considered (limited multi-component capabilities are included)

* No temperature stratification in the gas phase

* A default of no temperature gradient through vessel wall.

* Heat transfer is modelled as constant or simplified using empirical correlations

These choices make the problem much simpler to solve: First of all, the pure substance Helmholtz energy based equation of state (HEOS) in coolprop offers a lot of convenience in terms of the property pairs/state variables that can be set independently. Using only a single gas phase species also means that component balances is redundant and 2 or 3-phase flash calculations are not required. That being said, the principle used for a single component is more or less the same, even for multicomponent mixtures with potentially more than one phase.

In the latest revision 1-D transient heat conduction through the vessel wall is now an option if required for low thermal conductivity materials and e.g. type III/IV vessels (1-D heat transfer not yet implemented for fire heat load). Further, rigorous two-phase calculations are now possible by specifying a liquid level in the vessel (1-D heat transfer not yet implemented, only 0-D). Rigorous two-phase single component calculations is also supported in the latest release implementing an equilibrium assumption i.e. gas and liquid phase is in thermal equilibrium, although the vessel wall in contact with gas and liquid, respectively, is allowed to differ.

In case multi-component two-phase behaviour is required, the HydDown sibling ORS *openthermo* is recommended. This software package also implements a partial/non-equilibrium assumption for the gas and liquid phase (Andreasen and Stegelmann 2025).

1. @article{https://doi.org/10.1002/prs.70035,
2.   author = {Andreasen, Anders and Stegelmann, Carsten},
3.   title = {Open source pressure vessel blowdown modeling under partial phase equilibrium},
4.   journal = {Process Safety Progress},
5.   doi = {https://doi.org/10.1002/prs.70035},
6. }

A preprint of the paper is also available.

## 1.3 Getting the software

The source code can be obtained either from GitHub (via git or via the latest tar-ball release) or via pip. No packaged releases have currently been planned for conda.

The main branch is located here:

https://github.com/andr1976/HydDown
Clone the repo by:

1 git clone https://github.com/andr1976/HydDown.git

Running from source via git the dependencies must be installed manually from the repo root dir:

1 pip install -r requirements.txt

Installation of latest release via **pip** also install dependencies automatically:

1 pip install hyddown

In case `pip` links to a v2.7 of python you will get an error. If so try the following:

1 python3 -m pip install hyddown

where `python3` is the symlink or full path to the `python3` executable installed on your system.

## 1.4 Requirements

* Python (Check CoolProp for latest supported python version)

* Numpy

* matplotlib

* Coolprop (6.4.1)

* cerberus

* PyYaml

* pandas

* Scipy

* tqdm

The script is running on Windows 10 x64, with stock python installation from python.org and packages installed using `pip`. It should also run on Linux (it does on an Ubuntu image on GitHub) or in any conda environment as well, but this hasn't been checked.

## 1.5 Testing

Although testing is mainly intended for automated testing (CI) during development using github actions, testing of the installed package can be done for source install by:

1 python -m pytest
run from the root folder. In writing 33+ test should pass.

For the package installed with pip navigate to the install directory `../site-packages/hyddown` and run:

```
1 python -m pytest
```

## 1.6 Units of measure

The SI unit system is adapted for this project. The following common units are used in the present project and this also applies to the units used in the input files:

<table><caption>**Table 1.1:** Unit system</caption><thead><tr><th>Property</th><th>Unit</th><th>Comment</th></tr></thead><tbody><tr><td>Temperature</td><td>K</td><td>°C is used in plots</td></tr><tr><td>Pressure</td><td>Pa</td><td>bar is used in plots</td></tr><tr><td>Mass</td><td>kg</td><td></td></tr><tr><td>Volume</td><td>m<sup>3</sup></td><td></td></tr><tr><td>Time</td><td>s</td><td></td></tr><tr><td>Energy</td><td>J</td><td></td></tr><tr><td>Duty/power</td><td>W</td><td></td></tr><tr><td>Length</td><td>m</td><td></td></tr><tr><td>Area</td><td>m<sup>2</sup></td><td></td></tr><tr><td>Heat flux</td><td>W/m<sup>2</sup></td><td></td></tr><tr><td>Heat transfer coefficient</td><td>W/(m<sup>2</sup> K)</td><td></td></tr><tr><td>Thermal conductivity</td><td>W/(m K)</td><td></td></tr><tr><td>Density</td><td>kg/m<sup>3</sup></td><td></td></tr><tr><td>Heat capacity</td><td>J/(kg K)</td><td></td></tr></tbody></table>

As will be noted when presenting the equations implemented in the code, some of the equations utilise different units than the ones listed in tbl. 1.1. However, it is important to note that unit conversions are built in to the methods implemented, so the user shall not worry about unit conversion.
## 1.7 Credit

In the making of this document a great deal of material has been sourced (and modified) from a good colleague's M.Sc. thesis (Eriksen and Bjerre 2015), from co-published papers [Bjerre et al. (2017)](Andreasen et al. 2018) and from on-line material published under permissive licenses (with proper citation). Further, the making of this project would not have possible without the awesome CoolProp library (I. H. Bell et al. 2014). I am also thankful for enlightening discussions with colleague Jacob Gram Iskov Eriksen (Ramboll Energy, Denmark) and former Ramboll Energy colleague Carsten Stegelmann (ORS Consulting) in relation to vessel depressurization, nozzle flow and heat transfer considerations.

The present document is typeset using Markdown + pandoc with the Eisvogel template.

## 1.8 License

### MIT License

Copyright (c) 2021-2025 Anders Andreasen

Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
# 2 Usage

## 2.1 Basic usage

From the source repository, running the code is as simple as:

1. `python src/main.py input.yml`

where `main.py` is the main script and `input.yml` is the input file in Yaml syntax.

The Yaml input file is edited to reflect the system of interest.

Further, a usable copy of a main script is installed in the python installation Scripts/ folder:

1. `/python3x/Scripts/hyddown_main.py`

or from GitHub/in the source folder:

1. `HydDown/Scripts/hyddown_main.py`

Various example files are included for inspiration under site-packages:

1. `/site-packages/hyddown/examples/`

or in the `git/source` directory tree

1. `HydDown/src/hyddown/examples/`

## 2.2 Demos

A few demonstrations of the codes capability are available as `streamlit` apps:

* Vessel gas pressurisation/depressurization for various gases at https://share.streamlit.io/andr1976/hyddown/main/scripts/streamlit_app.py

* Response to external fire of vessel equipped with a PSV using the Stefan-Boltzmann fire equation at https://share.streamlit.io/andr1976/hyddown/main/scripts/streamlitSBapp.py
* Vessel depressurization including slow opening of valve at https://share.streamlit.io/andr1976/hyddown/main/scripts/Streamlit_cvapp.py

The various streamlit apps are also included in the scripts folder for running the scripts on a local machine.

## 2.3 Calculation methods

The following methods are implemented:

* Isothermal: constant temperature of the fluid during depressurization (for a very slow process with a large heat reservoir)

* Isenthalpic: constant enthalpy of the fluid, no heat transfer with surroundings, no work performed by the expanding fluid for depressurization, no work added from inlet stream

* Isentropic: constant entropy of the fluid, no heat transfer with surroundings, PV work performed by the expanding fluid

* Isenergetic: constant internal energy of the fluid

* Energy balance: this is the most general case and is based on the first law of thermodynamics applied to a flow process.

* Relief: Relief valve dimensioning for gas filled vessels subject to fire. This method provides a dynamic approach to relief valve dimensioning providing realistic orifice size, compared to the very conservative API521 approach.

For `isothermal/isenthalpic/isentropic/isenergetic` calculations the minimal input required are:

* Initial conditions (pressure, temperature)

* vessel dimensions (ID, length)

* valve parameters (Cd, diameter, back-pressure)

* Calculation setup (time step, end time)

* Type of gas

If heat transfer is to be considered the calculation type `energybalance` is required. A few options are possible:

* Fixed U: U-value, i.e. overall heat transfer coefficient, is required and ambient temperature required. Based on the ambient (external) temperature, the fluid temperature, and the U-value, the heat flux Q is calculated for each time step. In this calculation method the temperature of the vessel wall is not calculated:
$$Q = UA(T_{ambient} - T_{fluid})$$

* Fixed Q: The external heat flux, Q, applied to the fluid is specified and constant. The ambient temperature is not required. Further, the vessel wall temperature is redundant and not calculated.

* Specified h: The external heat transfer coefficient must be specified and either the internal heat transfer coefficient is provided or calculated from the assumption of natural convection from a vertical cylinder at high Gr number. Ambient temperature is required. Using this method the wall temperature is calculated from an energy balance over the vessel wall taking in and out flux to/from the external ambient plenum as well as heat flux to/from the fluid inventory into account.

* Fire: The Stefan-Boltzmann equation is applied for estimating the external heat duty. The fire heat flux depends on the vessel wall temperature, and the wall temperature is continuously updated as in the method with specified h.

More elaborate description of the required input for the different calculation types are provided in Sec. 2.7.

## 2.4 Script

HydDown comes with a script which can be used as a command-line tool to start calculations. If an input filename (path) is given as the first argument, this input file will be used. If no arguments are passed, the script will look for an input file with the name `input.yml`. The content of the main script is shown below:

```python
1 import yaml
2 import sys
3 from hyddown import HydDown
4
5 if __name__ == "__main__":
6     if len(sys.argv) > 1:
7         input_filename = sys.argv[1]
8     else:
9         input_filename = "input.yml"
10
11     with open(input_filename) as infile:
12         input = yaml.load(infile, Loader=yaml.FullLoader)
13
14 hdown=HydDown(input)
15 hdown.run()
16 hdown.verbose=1
17 hdown.plot()
18
```
## 2.5 Module import

To use HydDown simple import the main calculation class HydDown.

```python
1 from hyddown import HydDown
```

## 2.6 Input file examples

When using HydDown a dictionary holding all relevant input in order for HydDown to do vessel calculations shall be provided when the class is initialized. One way is to read an input file. For HydDown a Yaml format is chosen, but if JSON is a preference this should also work with the Yaml parser being substituted with a JSON parser.

An example of a minimal file for an isentropic vessel depressurization (no heat transfer) is shown below

```yaml
1 vessel:
2   length: 1.524
3   diameter: 0.273
4 initial:
5   temperature: 388.0
6   pressure: 15000000.
7   fluid: "N2"
8 calculation:
9   type: "isentropic"
10  time_step: 0.05
11  end_time: 100.
12 valve:
13   flow: "discharge"
14   type: "orifice"
15   diameter: 0.00635
16   discharge Coef: 0.8
17   back pressure: 101300.
```

A more elaborate example which includes heat transfer and with validation data (some data points dropped for simplicity) included:

```yaml
1 vessel:
2   length: 1.524
3   diameter: 0.273
4   thickness: 0.025
5   heat_capacity: 500
6   density: 7800.
7   orientation: "vertical"
8 initial:
9   temperature: 288.0
```
10 pressure: 15000000.
11 fluid: "N2"
12 calculation:
13 type: "energybalance"
14 time_step: 0.05
15 end_time: 100.
16 valve:
17 flow: "discharge"
18 type: "orifice"
19 diameter: 0.00635
20 discharge Coef: 0.8
21 back pressure: 101300.
22 heat transfer:
23 type: "specified_h"
24 temp ambient: 288.
25 h_outer: 5
26 h_inner: 'calc'
27 validation:
28 temperature:
29 gas high:
30 time: [0.050285, ..., 99.994]
31 temp: [288.93, ..., 241.29]
32 gas low:
33 time: [0.32393, ..., 100.11]
34 temp: [288.67, ..., 215.28]
35 wall low:
36 time: [0.32276, ..., 100.08]
37 temp: [288.93, ..., 281.72]
38 wall high:
39 time: [0.049115, ..., 100.06]
40 temp: [289.18, ..., 286.09]
41 pressure:
42 time: [0.28869, ..., 98.367]
43 pres: [150.02, ..., 1.7204]

## 2.7 Input fields and hierarchy

In the following the full hierarchy of input for the different calculation types is summarised.

At the top level the following fields are accepted, with the last being optional and the second last
dependant on calculation type:

1 initial: mandatory

2 vessel: mandatory

3 calculation: mandatory

4 valve: mandatory

5 heat transfer: depends on calculation type

6 validation: optional
### 2.7.1 Calculation

The subfields under `calculation`, with value formats and options are:

1 calculation:

2 type: "isothermal", "isentropic", "isenthalpic", "constantU", "energybalance"

3 time_step: number

4 end_time: number

The simulation end time is specified as well as the fixed time step used in the integration of the differential equations to be solved. The four main calculation types are shown as well.

### 2.7.2 Vessel

1 vessel:

2 length: number, mandatory

3 diameter: number, mandatory

4 thickness: number, required when heat transfer is calculated

5 heat_capacity: number, required when heat transfer is calculated

6 density: number, required when heat transfer is calculated

7 thermal_conductivity: number, required for 1-D transient heat transfer

8 liner_thickness: number, required only for bi-material 1-D transient heat transfer

9 liner_heat_capacity: number, required only for bi-material 1-D transient heat transfer

10 liner_density: number, required only for bi-material 1-D transient heat transfer

11 liner_thermal_conductivity: number, required only for bi-material 1-D transient heat transfer

12 orientation: string, required when heat transfer is calculated

13 liquid_level: number, optional for two-phase

14 type: string, optional for heads other than flat-end, "DIN", "ASME F&D", or "Hemispherical"

### 2.7.3 Initial

1 initial:

2 temperature: number, mandatory

3 pressure: number, mandatory

4 fluid: string, mandatory, e.g. "N2", "H2"
### 2.7.4 Valve

The `valve` field determines the mode of the process if it for depressurization/discharge from the vessel using the value `discharge` or if mass flow is entering the vessel, filling it up/increasing pressure with the value `filling`.

Different types of mass flow devices can be specified:

* Restriction `orifice`

* Relief valve/pressure safety valve (only for `discharge` not `filling`)

* Control valve (`controlvalve`)

* A specified mass flow (`mdot`)

For the physical devices a `backellaner` is required for the flow calculations. The value of the `backellaner` is also used to specify the reservoir pressure when the vessel is filled. See also Sec. 3.2 for details about the calculation of flow rates.

1. **valve:**
2. **flow:** `string`, mandatory `"discharge"` or `"filling"`
3. **type:** `string`, mandatory `"orifice"`, `"controlvalve"`, `"psv"`, `"mdot"`
4. **back风pressure:** `number`, required **for** type `"orifice"`, `"controlvalve"` and `"psv"`
5. **diameter:** `number`, required **for** `"orifice"` and `"psv"`
6. **discharge_coef:** `number`, required **for** `"orifice"` and `"psv"`
7. **Cv:** `number`, required **for** `"control_valve"`
8. **characteristic:** `string`, optional **for** `"control_valve"`
9. **timeConstant:** `number`, optional **for** `"control_valve"`

For the calculations using a control valve as mass transfer device a linear rate for the actuator can be specified using the `timeConstant` field. An associated valve characteristic is associated in order to translate the linear actuator rate to an opening `Cv` (of max. `Cv`). Three characteristics can be specified:

* Linear

* Equal percentage (exponential)

* Quick opening/fast (square root)

### 2.7.5 Heat transfer

For more information about the actual estimation of heat transfer see also Sec. 3.3. For the case of `fire` heat input predetermined parameters for the Stefan-Boltzmann fire equation are used to calculate different background heat loads cf. Tbl. 2.1
<table><caption>**Table 2.1:** Fire heat loads</caption><thead><tr><td>Source</td><td>Fire type</td><td>Heat load (kW/m²)</td></tr></thead><tbody><tr><td>API521</td><td>Pool</td><td>60</td></tr><tr><td>API521</td><td>Jet</td><td>100</td></tr><tr><td>Scandpower</td><td>Pool</td><td>100</td></tr><tr><td>Scandpower</td><td>Jet</td><td>100</td></tr></tbody></table>

1. `heat_transfer`:

2. `type: string, mandatory, "specified_h", "specified_Q", "specified_U", "s-b"`

3. `temp_ambient: number, required for type "specified_h", "specified_U"`

4. `h_outer: number, required for type "specified_h"`

5. `h_inner: number or 'calc', required for type "specified_h"`

6. `U_fix: number, required for type "specified_U"`

7. `Q_fix: number, required for type "specified_Q"`

8. `fire: string, required for type "s-b", 'api_pool', 'api.jet', 'scandpower_pool', 'scandpower.jet'`

9. `D_throat: number, required for flow type "filling", set to vessel ID as a starting point`

### 2.7.6 Validation

In order to plot measured data against simulated data the field `validation` is included.

The following arrays (one or more) are supported in the sub-field `temperature`:

* `gas_high`: highest measured values of the bulk gas

* `gas_low`: lowest measured values of the bulk gas

* `gas_mean`: average measured values of the bulk gas

* `wall_mean`: average measured values of the vessel (inner) wall

* `wall_high`: i.e. highest measured values of the vessel (inner) wall

* `wall_low`: i.e. lowest measured values of the vessel (inner) wall

* `wall_outer`: i.e. external measured values of the vessel (outer) wall

* `wall_inner`: i.e. internal measured values of the vessel (inner) wall

For each of the above fields arrays for `time` and `temp` shall be supplied with matching length.

A field for measured vessel pressure is also possible, where `time` and `pres` shall be identical length arrays. See also example below.
1 validation:
2 temperature:
3 gas_high:
4 time: [0.050285, ..., 99.994]
5 temp: [288.93, ..., 241.29]
6 gas_low:
7 time: [0.32393, ..., 100.11]
8 temp: [288.67, ..., 215.28]
9 wall_low:
10 time: [0.32276, ..., 100.08]
11 temp: [288.93, ..., 281.72]
12 wall_high:
13 time: [0.049115, ..., 100.06]
14 temp: [289.18, ..., 286.09]
15 pressure:
16 time: [0.28869, ..., 98.367]
17 pres: [150.02, ..., 1.7204]
<figure><img src="image_3.png" /></figure>
# 3 Theory

In this chapter the basic theory and governing equations for the model implementation in HydDown is presented. The following main topics are covered:

* thermodynamics

* mass transfer

* heat transfer

## 3.1 Thermodynamics

### 3.1.1 Equation of state

The equation of state used by HydDown is the Helmholtz energy formulation as implemented in CoolProp (I. H. Bell et al. 2014). Most of the text in the present section has been sourced from the CoolProp documentation to be as accurate and true to the source as possible. The Helmholtz energy formulation is a convenient construction of the equation of state because all the thermodynamic properties of interest can be obtained directly from partial derivatives of the Helmholtz energy.

It should be noted that the EOS are typically valid over the entire range of the fluid, from subcooled liquid to superheated vapor, to supercritical fluid. In general, the EOS are based on non-dimensional terms $\delta$ and $\tau$, where these terms are defined by

$$ \delta = \rho/\rho_c $$

$$ \tau = T_c/T $$

where $\rho_c$ and $T_c$ are the critical density of the fluid if it is a pure fluid. For pseudo-pure mixtures, the critical point is typically not used as the reducing state point, and often the maximum condensing temperature on the saturation curve is used instead.

The non-dimensional Helmholtz energy of the fluid is given by

$$ \alpha = \alpha^0 + \alpha^r $$
where $\alpha^0$ is the ideal-gas contribution to the Helmholtz energy, and $\alpha^r$ is the residual Helmholtz energy contribution which accounts for non-ideal behavior. For a given set of $\delta$ and $\tau$, each of the terms $\alpha^0$ and $\alpha^r$ are known. The exact form of the Helmholtz energy terms is fluid dependent, but a relatively simple example is that of Nitrogen, which has the ideal-gas Helmholtz energy of

$$ \alpha^0 = \ln \delta + a_1 \ln \tau + a_2 + a_3 \tau + a_4 \tau^{-1} + a_5 \tau^{-2} + a_6 \tau^{-3} + a_7 \ln[1 - \exp(-a_8 \tau)] $$

and the non-dimensional residual Helmholtz energy of

$$ \alpha^r = \sum_{k=1}^{6} N_k \delta^{i_k} \tau^{j_k} + \sum_{k=7}^{32} N_k \delta^{i_k} \tau^{j_k} \exp(-\delta^{l_k}) + \sum_{k=33}^{36} N_k \delta^{i_k} \tau^{j_k} \exp(-\phi_k(\delta - 1)^2 - \beta_k(\tau - \gamma_k)^2) $$

and all the terms other than $\delta$ and $\tau$ are fluid-dependent correlation parameters.

The other thermodynamic parameters can then be obtained through analytic derivatives of the Helmholtz energy terms. For instance, the pressure is given by

$$ p = \rho RT \left[ 1 + \delta \left( \frac{\partial \alpha^r}{\partial \delta} \right)_\tau \right] $$

and the specific internal energy by

$$ \frac{u}{RT} = \tau \left[ \left( \frac{\partial \alpha^0}{\partial \tau} \right)_\delta + \left( \frac{\partial \alpha^r}{\partial \tau} \right)_\delta \right] $$

and the specific enthalpy by

$$ \frac{h}{RT} = \tau \left[ \left( \frac{\partial \alpha^0}{\partial \tau} \right)_\delta + \left( \frac{\partial \alpha^r}{\partial \tau} \right)_\delta \right] + \delta \left( \frac{\partial \alpha^r}{\partial \delta} \right)_\tau + 1 $$

which can also be written as

$$ \frac{h}{RT} = \frac{u}{RT} + \frac{p}{\rho RT} $$

The specific entropy is given by

$$ \frac{s}{R} = \tau \left[ \left( \frac{\partial \alpha^0}{\partial \tau} \right)_\delta + \left( \frac{\partial \alpha^r}{\partial \tau} \right)_\delta \right] - \alpha^0 - \alpha^r $$
and the specific heats at constant volume and constant pressure respectively are given by

$$ \frac{c_v}{R} = -\tau^2 \left[ \left( \frac{\partial^2 \alpha^0}{\partial \tau^2} \right)_\delta + \left( \frac{\partial^2 \alpha^r}{\partial \tau^2} \right)_\delta \right] $$

$$ \frac{c_p}{R} = \frac{c_v}{R} + \frac{\left[ 1 + \delta \left( \frac{\partial \alpha^r}{\partial \delta} \right)_\tau - \delta \tau \left( \frac{\partial^2 \alpha^r}{\partial \delta \partial \tau} \right) \right]^2}{\left[ 1 + 2\delta \left( \frac{\partial \alpha^r}{\partial \delta} \right)_\tau + \delta^2 \left( \frac{\partial^2 \alpha^r}{\partial \delta^2} \right)_\tau \right]} $$

The EOS is set up with temperature and density as the two independent properties, but often other inputs are known, most often temperature and pressure because they can be directly measured. As a result, if the density is desired for a known temperature and pressure, it can be obtained iteratively.

### 3.1.2 First law for flow process

The control volume sketched in Fig. 3.1, separated from the surrounding by a control surface, is used as a basis for the analysis of an open thermodynamic system with flowing streams (fs) in and out, according to (Smith, Van Ness, and Abbott 1996)

A general mass balance or continuity equation can be written:
<figure><img src="image_4.png" /><figcaption>**Figure 3.1:** Control volume with one entrance and one exit. The image has been sourced from (WikiMedia 2008,).</figcaption></figure>

$$ \frac{m_{cv}}{dt} + \Delta (\dot{m})_{fs} = 0 \qquad (3.1) $$

The first term is the accumulation term i.e. the rate of change of the mass inside the control volume, $m_{cv}$, and the $\Delta$ in the second term represents the difference between the outflow and the inflow

$$ \Delta (\dot{m})_{fs} = \dot{m}_2 - \dot{m}_1 $$

An energy balance for the control volume, with the first law of thermodynamics applied, needs to account for all the energy modes with can cross the control surface. Each stream has a total energy

$$ U + \frac{1}{2}u^2 + zg $$

where the first term is the specific internal energy, the second term is the kinetic energy and the last
term is the potential energy. The rate at which each stream transports energy in or out of the control volume is given by

$$ \dot{m}(U + \frac{1}{2}u^2 + zg) $$

and in total

$$ \Delta \left[ \dot{m} \left( U + \frac{1}{2} u^2 + zg \right) \right]_{fs} $$

Furthermore, work (not to be confused with shaft work) is also associated with each stream in order to move the stream into or out from the control volume (one can think of a hypothetical piston pushing the fluid at constant pressure), and the work is equal to *PV* on the basis of the specific fluid volume. The work rate for each stream is

$$ \dot{m}(PV) $$

and in total

$$ \Delta [\dot{m}(PV)]_{fs} $$

Further, heat may be transferred to (or from) the control volume at a rate $\dot{Q}$ and shaft work may be applied, $\dot{W}_{shaft}$. Combining all this with the accumulation term given by the change in total internal energy the following general energy balance can be written:

$$ \frac{d(mU)_{cv}}{dt} + \Delta \left[ \dot{m} \left( U + \frac{1}{2} u^2 + zg \right) \right]_{fs} + \Delta [\dot{m}(PV)]_{fs} = \dot{Q} + \dot{W}_{shaft} $$

Applying the relation $H = U + PV$, setting $\dot{W}_{shaft} = 0$ since no shaft work is applied to or extracted from the vessel, and assuming that kinetic and potential energy changes can be omitted the energy balance simplifies to:

$$ \frac{d(mU)_{cv}}{dt} + \Delta [\dot{m}H]_{fs} = \dot{Q} $$

The equation can be further simplified if only a single port acting as either inlet or outlet is present:

$$ \frac{d(mU)_{cv}}{dt} + \dot{m}H = \dot{Q} \qquad (3.2) $$
where the sign of $\dot{m}$ determines if the control volume is either emptied or filled. The continuity equation Eq. 3.1 and the energy balance Eq. 3.2 combined with the equation of state are the key equations that shall be solved/integrated in order to calculate the change in temperature and pressure as a function of time.

## 3.2 Flow devices

### 3.2.1 Restriction Orifice

When a fluid flows through a restriction or opening such as an orifice, the velocity will be affected by conditions upstream and downstream. If the upstream pressure is high enough, relative to the downstream pressure, the velocity will reach the speed of sound (Ma = 1) and the flow rate obtained will be the critical flow rate. This condition is referred to as choked flow. The maximum downstream pressure for the flow to still be sonic (Ma = 1), is when $P_d = P_c$. The ratio of the critical and upstream pressure is defined by equation Eq. 3.3.

$$ \frac{P_c}{P_u} = \left( \frac{2}{k+1} \right)^{\frac{k}{k-1}} \qquad (3.3) $$

* $P_c$ is the critical pressure. [kPa]

* $P_d$ is the downstream pressure. [kPa]

* $P_u$ is the upstream pressure. [kPa]

* k is the isentropic coefficient, approximated by the ideal gas heat capacity ratio $C_p/C_v$. [-]

In order to calculate the mass flow rate through an orifice equations are used based on literature from the Committee for the Prevention of Disasters (Bosch and Weterings 2005).

To account for the difference in choked and non-choked flow a set limit pressure is introduced as in equation Eq. 3.4. If the downstream pressure, $P_{down}$, is below the pressure limit, $P_{\{limit}\}$, then the flow is choked, and the pressure used, $P_{used}$, in equation Eq. 3.5 should be the pressure limit, $P_{limit}$. Otherwise if the downstream pressure, $P_{down}$, is greater than or equal to the pressure limit, $P_{limit}$, the flow is no longer choked and the pressure used should be the downstream pressure, $P_{down}$ (Bosch and Weterings 2005).

$$ P_{limit} = P_{up} \cdot \left( \frac{2}{k+1} \right)^{\frac{k}{k-1}} \qquad (3.4) $$

$$ \dot{m}_{flow} = C_d \cdot A \cdot \sqrt{\left(\frac{2k}{k-1}\right) \cdot P_{up} \cdot \rho \cdot \left(\frac{P_{used}}{P_{up}}\right)^{\frac{2}{k}} \left(1 - \left(\frac{P_{used}}{P_{up}}\right)^{\frac{k-1}{k}}\right)} \qquad (3.5) $$
* $\rho$ is the density of the gas upstream. [kg/m³]

* $P_{limit}$ is the pressure limit of the upstream absolute pressure. [bara]

* $P_{up}$ is the absolute pressure upstream of the orifice. [bara]

* $k$ is the ratio of the heat capacities at constant pressure, $C_p$, and at constant volume, $C_v$.

* $P_{down}$ is the absolute pressure downstream of the orifice. [bara]

* $P_{used}$ is the pressure used in the mass flow equation based on choked or non-choked conditions. [bara]

* $\dot{m}_{flow}$ is the mass flow through the orifice. [kg/s]

* $C_d$ is the discharge coefficient of the orifice opening. [–]

* $A$ is the cross sectional area of the orifice. [m²]

### 3.2.2 Pressure safety valve / Relief valve

A PSV / relief valve is a mechanical device actuated by the static pressure in the vessel and a conventional PSV is often used for gas/vapor systems. A conventional PSV is a spring-loaded device which will activate at a predetermined opening pressure, and relieve the vessel pressure until a given reset pressure has been reached. Both the opening pressure and the reset pressure is above the vessel operating pressure and the PSV will remain closed until the pressure inside the vessel increases to the opening pressure. The operation of a conventional spring-loaded PSV is based on a force balance. A conventional PSV can be seen in Fig. 3.2. A spring exerts a force on a disc blocking the inlet of the PSV. When the pressure inside the vessels reaches the opening pressure, the force exerted on the disc by the gas will be larger than the force exerted by the spring and the PSV will open and allow the gas to flow out of the vessel. The flow of gas out of the vessels will lower the pressure and thereby also the force exerted on the disc. When the pressure in the vessels is reduced to the reset pressure, the PSV will close and the disc will again hinder the gas flow.
<figure><img src="image_5.png" /><figcaption>**Figure 3.2:** Conventional/pop action PSV adapted from (Eriksen and Bjerre 2015) and (API 2014a)</figcaption></figure>

The relief valve model implemented in HydDown is the API 520 equations (API 2014a) for gas relief for both sonic/critical as well as subcritical flow. No corrections factors are implemented in HydDown.

For sonic flow (critical flow), as indicated in equation, the mass flow through the PSV can be determined by equation eq. 3.6.

$$W = \frac{AC \cdot K_d \cdot K_b \cdot K_c \cdot P_1}{\sqrt{\frac{T \cdot Z}{M}}} \qquad (3.6)$$

* A is the effective discharge area. [mm²]

* W is the mass flow through the device. [kg/h]

* C is a coefficient, as a function of k, as defined in equation eq. 3.7.

* K<sub>d</sub>, K<sub>b</sub>, and K<sub>c</sub> are correction factors.

* P<sub>1</sub> is the allowable upstream absolute pressure. [kPa]

* T is the temperature of the inlet gas at relieving conditions. [K]

* M is the molecular mass of the gas at relieving conditions. [kg/kmol]

* Z is the compressibility factor for the gas.

K<sub>d</sub> is the effective coefficient of discharge, with a typical value of 0.975, for an installed PSV. K<sub>b</sub> is a
back-pressure correction factor between 0 and 1, assumed to be 1. $K_c$ is a correction factor used when a rupture disk is installed upstream, otherwise it is 1. In the present implementation a value of 1 is assumed.

$$C = 0.03948 \sqrt{k \left(\frac{2}{k+1}\right)^{\left(\frac{k+1}{k-1}\right)} \quad (3.7)$$

For subsonic flow (subcritical flow), the effective discharge area of the PSV is determined by equation eq. 3.8.

$$W = \frac{A \cdot F_2 \cdot K_d \cdot K_c}{17.9 \sqrt{\frac{T \cdot Z}{M \cdot P_1 \cdot (P_1 - P_2)}} \quad (3.8)$$

$F_2$ is the coefficient of subcritical flow which can be determined from eq. 3.9.

$$F_2 = \sqrt{\left(\frac{k}{k-1}\right) r^{\left(\frac{2}{k}\right)} \left(\frac{1 - r^{\left(\frac{k-1}{k}\right)}}{1 - r}\right) \quad (3.9)$$

where $r$ is the ratio of backpressure to upstream relieving pressure, $P_2/P_1$.

When modelling a pop action PSV/relief valve under dynamic conditions, the valve will go from closed to fully open in a short period of time when the set pressure, $P_{set}$, is reached. The pop action is illustrated in Fig. 3.3 which shows the opening and closing hysteresis of the PSV as a function of pressure. In order to close the pressure shall be reduced below the reseat pressure.

<figure><img src="image_6.png" /><figcaption>**Figure 3.3:** Relief valve hysteresis adapted from (Eriksen and Bjerre 2015)</figcaption></figure>
When specifying PSVs it common to use standard API sizes as shown in tbl. 3.1

<table><caption>**Table 3.1:** Standard PSV orifice sizes according to API</caption><thead><tr><th>Size</th><th>Area [in²]</th><th>Area [m²]</th></tr></thead><tbody><tr><td>D</td><td>0.110</td><td>7.09676 · 10<sup>-5</sup></td></tr><tr><td>E</td><td>0.196</td><td>1.26451 · 10<sup>-4</sup></td></tr><tr><td>F</td><td>0.307</td><td>1.98064 · 10<sup>-4</sup></td></tr><tr><td>G</td><td>0.503</td><td>3.24515 · 10<sup>-4</sup></td></tr><tr><td>H</td><td>0.785</td><td>5.06450 · 10<sup>-4</sup></td></tr><tr><td>J</td><td>1.287</td><td>8.30320 · 10<sup>-4</sup></td></tr><tr><td>K</td><td>1.838</td><td>1.18580 · 10<sup>-3</sup></td></tr><tr><td>L</td><td>2.853</td><td>1.84064 · 10<sup>-3</sup></td></tr><tr><td>M</td><td>3.600</td><td>2.32257 · 10<sup>-3</sup></td></tr><tr><td>N</td><td>4.340</td><td>2.79999 · 10<sup>-3</sup></td></tr><tr><td>P</td><td>6.380</td><td>4.11612 · 10<sup>-3</sup></td></tr><tr><td>Q</td><td>11.050</td><td>7.12901 · 10<sup>-3</sup></td></tr><tr><td>R</td><td>16.000</td><td>1.03225 · 10<sup>-2</sup></td></tr><tr><td>T</td><td>26.000</td><td>1.67741 · 10<sup>-2</sup></td></tr></tbody></table>

### 3.2.3 Control Valve

For calculating the mass flow through a control valve, the ANSI/ISA [Borden (1998)](ISA 1995) methodology also described in IEC 60534 (IEC 2011) is applied.

The flow model for a compressible fluid in the turbulent regime is:

$$W = CN_6F_pY\sqrt{x_{sizing}p_1\rho_1}$$

or equivalently:

$$W = CN_8F_p p_1 Y \sqrt{\frac{x_{sizing} M}{T_1 Z_1}}$$
* C is the flow coefficient ($C_v$ or $K_v$)

* $N_8$ is a unit specific constant, 94.8 for $C_v$ and bar as pressure unit

* $F_P$ is a piping geometry factor [-]

* Y is the expansion factor [-]

* $x_{sizing}$ is the pressure drop used for sizing[-]

* $p_1$ is the upstream pressure [bar]

* $\rho_1$ is the upstream density [kg/m³]

* M is the molecular weight [kg/kmol]

* $T_1$ is the upstream temperature [K]

* $Z_1$ is the upstream compressibility [-]

In HydDown, the piping geometry factor is not yet implemented and assumed to be 1. The pressure drop ratio $x_{sizing}$ used for sizing is determined as the lesser of the actual pressure drop ratio, $x$, and the choked pressure drop ratio $x_{choked}$. The actual pressure drop ratio is given by:

$$x \frac{\Delta p}{p_1}$$

The pressure drop ratio at which flow no longer increases with increased value in pressure drop ratio is the choked pressure drop ratio, given by the following equation:

$$x_{choked} = F_{\gamma}x_{TP}$$

The factor $x_T$ is based on air near atmospheric pressure as the flowing fluid with a specific heat ratio of 1.40. If the specific heat ratio for the flowing fluid is not 1.40, the factor $F_\gamma$ is used to adjust $x_T$. Use the following equation to calculate the specific heat ratio factor:

$$F_{\gamma} = \frac{\gamma}{1.4}$$

where $\gamma$ is the ideal gas $C_p/C_v$. It should be noted that the above equation has been derived from perfect gas behaviour and extension of an orifice model with $\gamma$ in the range of 1.08 to 1.65. If used outside the assumptions flow calculations may become inaccurate.

The expansion factor Y accounts for the change in density as the fluid passes from the valve inlet to the vena contracta. It also accounts for the change in the vena contracta area as the pressure differential is varied.

$$Y = 1 - \frac{x_{sizing}}{3x_{choked}}$$
## 3.3 Heat transfer

### 3.3.1 Natural convection

Experiments have indicated that the internal heat transfer mechanism for a vessel subject to depres-surization can be well approximated by that of natural convection as found from measured Nusselt numbers being well correlated with Rayleigh number, with no apparent improvement in model performance by inclusion of the Reynolds number in the model (Woodfield, Monde, and Mitsutake 2007).

To determine the heat transfer for the gas-wall interface, the following is applied cf. equation Eq. 3.10:

$$ \frac{dQ}{dt} = \dot{Q} = hA(T_s - T_{gas}) \qquad (3.10) $$

* $dQ$ is the change in thermal energy due to convective heat transfer. [J]

* $dt$ is the change in time during the heat transfer. [s]

* $h$ is the convective heat transfer. [W/m²·K]

* $A$ is the area normal to the direction of the heat transfer. [m²]

* $T_s$ is the inner surface temperature of the geometry. [K]

* $T_{gas}$ is the temperature of the bulk gas inside the vessel. [K]

The convective heat transfer will need to be estimated for the gas-wall interface, by the use of empirical relations for the Nusselt number. The Nusselt number describes the ratio of convective heat transfer to conductive heat transfer, normal to a surface area, as given in equation eq. 3.11.

$$ Nu = \frac{hL}{k} \qquad (3.11) $$

* $Nu$ is the Nusselt number. [-]

* $h$ is the convective heat transfer. [W/m²·K]

* $L$ is a characteristic length of the geometry. [m]

* $k$ is the thermal conductivity of the gas. [W/m·K]

The characteristic length $L$ used is the height of the gas volume i.e. the length of a vertical vessel or the diameter of a horizontal vessel.

The empirical correlations used to calculate the Nusselt number of the gas-wall interface is a function of the Rayleigh number, which can be defined by the Grashof number and Prandtl number, as in equation Eq. 3.12:

$$ Ra = Gr \cdot Pr \qquad (3.12) $$

* $Ra$ is the Rayleigh number. [-]
* $Gr$ is the Grashof number. [-]

* $Pr$ is the Prandtl number. [-]

The Grashof number is a dimensionless number which approximates the ratio of the buoyancy forces to viscous forces, as given in equation [Eq. 3.13]:

$$Gr = \frac{\beta g \rho^2 L^3 \Delta T}{\mu^2} \qquad (3.13)$$

The Prandtl number is a dimensionless number defined as the ratio of the momentum diffusivity to thermal diffusivity, as given in equation Eq. 3.14:

$$Pr = \frac{c_p \mu}{k} \qquad (3.14)$$

* $\beta$ is the coefficient of volume expansion. [1/K]

* $g$ is the standard acceleration of gravity. [m/s²]

* $\rho$ is the gas density. [kg/m³]

* $L$ is the characteristic length. [m]

* $\Delta T$ is the temperature difference of the surface and gas. [K]

* $\mu$ is the dynamic viscosity. [kg/m·s]

* $c_p$ is the heat capacity of gas. [J/kg·K]

* $k$ is the thermal conductivity of gas. [J/m·K]

It is important to note that the properties in the above equations shall be evaluated at the fluid film temperature which can be approximated by the average of the the fluid bulk temperature and the vessel wall temperature (Geankoplis 1993).

The mechanism of heat transfer on the outside of the vessel at ambient conditions is similar to the above. In HydDown the external heat transfer is not modelled currently. A heat transfer coefficient shall be provided.

### 3.3.2 Mixed convection

Experiments have indicated that the internal heat transfer mechanism for a vessel subject to filling can be well approximated by that of combined natural convection and forced convection as found from measured Nusselt numbers being well correlated with Rayleigh and Reynolds number (Woodfield, Monde, and Mitsutake 2007).

For mixed convection the effective Nusselt number, $Nu$, can be approximated by

$$Nu = (Nu_{forced} + Nu_{natural})^{\frac{1}{n}}$$
During charging with different gases (H₂, N₂ and Argon), Woodfield et al. demonstrated that in order to provide a good fit to the experimentally determined Nusselt number a correlation based on both Reynolds and Rayleigh number was necessary. The following formula was found to provide a good fit (n=1):

$$Nu = Nu_{forced} + Nu_{natural} = 0.56Re_{d}^{0.67} + 0.104Ra_{H}^{0.352}$$

### 3.3.3 Nucleate boiling heat transfer

For heat transfer between vessel wall and liquid, this is treated as nucleate boiling and estimated via the Rohsenow correlation (Rohsenow 1951) as implemented in the `ht` python library (C. Bell et al. 2024).

$$h = \mu_L \Delta H_{vap} \left[ \frac{g(\rho_L - \rho_v)}{\sigma} \right]^{0.5} \left[ \frac{C_{p,L} \Delta T_e^{2/3}}{C_{sf} \Delta H_{vap} Pr_L^n} \right]^3$$

* $\rho_L$ is the density of the liquid [kg/m³]

* $\rho_v$ is the density of the produced gas [kg/m³]

* $\mu_L$ is the viscosity of liquid [Pa s]

* $k_L$ is the thermal conductivity of liquid [W/m K]

* $C_{p,L}$ is the heat capacity of liquid [J/kg K]

* $H_{vap}$ is the heat of vaporization of the fluid [J/kg]

* $\sigma$ is the surface tension of liquid [N/m]

* $T_e$ is the excess wall temperature, [K]

* $C_{sf}$ is Rohsenow coefficient specific to fluid and metal [-]

* n constant, 1 for water, 1.7 (default) for other fluids usually [-]

In the present study a value of the Rohsenow coefficient of 0.013 is applied and the exponent n is set to 1.7. See also (Pioro 1999; Jabardo et al. 2004).

### 3.3.4 Conduction

For accurate prediction of the outer and especially the inner wall temperature for correct estimation of internal convective heat transfer and the average material temperature, the general equation of 1-D unsteady heat transfer shall be solved:

$$\frac{\delta T}{\delta t} = \frac{k}{\rho c_p} \frac{\delta^2 T}{\delta x^2}$$
* T is temperature

* x is the spatial (1-D) coordinate

* k is the thermal conductivity

* ρ is the material density

* Cₚ is the specific heat capacity

To be solved, the initial values and boundary values must be specified. In its default state (if thermal cobductivity is not applied for the vessel), HydDown does not include the unsteady heat transfer model, i.e. the assumption is that the temperature from outer to inner surface is uniform and equal to the average temperature. This is obviously a crude approximation, but might be justified depending in the Biot number:

$$Bi = \frac{hL}{k}$$

The Biot number is a simple measure of the ratio of the heat transfer resistances at the surface of a body to the inside of a body. The ratio gives an indication to which extent the temperature will vary in space (gradient) when the body is subject to a displacement in temperature at the surface boundary layer. Striednig et al. (Striednig et al. 2014) concluded that for a type I (steel) cylinder the Biot number was approx. 0.03 and hence the error in assuming a uniform temperature in the vessel wall was low.

With a typical thermal conductivity of 45 W/mK for steel and a heat transfer coefficient up to 600 W/m²K (Woodfield, Monde, and Mitsutake 2007) the Biot number for a vessel with a wall thickness of 2 cm is 0.27. This is significantly higher than that approximated by (Striednig et al. 2014). However, the Biot number is significantly lower than 1, and the assumption of a uniform temperature is reasonable. However, for increased wall thickness, and/or for different materials with lower thermal conductivity, the error may grow to an unacceptable level.

Especially for vessels with low conductivity materials (or very thick walls) accurate estimation of the vessel wall temperatures requires the 1-D transient heat transfer problem to be solved. HydDown incorporates the *thermesh* code provided under an MIT license by Wouter Grouve (Grouve, n.d.). The implemented model applies Cartesian coodinates as also applied in the h2fills program by NREL (Kuroki et al. 2021), i.e. the curved vessel wall is assumed a flat plate. This may to some extent be justified by the relatively large radius of a cylindrical storage container compared to the vessel wall thickness (see also justification references in (Kuroki et al. 2021)).

In HydDown, the 1-D transient heat conductivity problem is solved with *thermesh* assuming temperature independent vessel material density, heat capacity and thermal conductivity applying the Crank-Nicholson scheme using Neumann boundary conditions:
$$ -k_z \left. \frac{\partial T}{\partial z} \right|_{z=0} = q_{\text{left}}(t), \quad -k_z \left. \frac{\partial T}{\partial z} \right|_{z=L} = q_{\text{right}}(t) $$

The heat flux at the outer (left) and inner (right) is calculated/updated for each major time step for the explicit solver for mass and energy balances. For each major time step thermesh is used to update the outer and inner wall temperature by solving the 1-D transient heat conductivity using a minor time step equal to 1/10th of the major time step and 11 nodes. For each time thermesh is called during major time steps the mesh is initialised with the temperature profile at the end of the former major time step.

### 3.3.5 Fire heat loads

The heat transfer from the flame to the shell is modelled using the recommended approach from Scandpower (Hekkelstrand and Skulstad 2004) and API (API 2014b). The heat transfer from the flame to the vessel shell is divided into radiation, convection, and reradiation as seen in equation Eq. 3.15:

$$ q_f = \underbrace{\alpha_s \cdot \varepsilon_f \cdot \sigma \cdot T_f^4}_{\text{Radiation}} + \underbrace{h_f \cdot (T_f - T_s(t))}_{\text{Convection}} - \underbrace{\varepsilon_s \cdot \sigma \cdot T_s(t)^4}_{\text{Reradiation}} \quad (3.15) $$

* $q_f$ is the flame heat flux. [W/m²]

* $\alpha_s$ is the vessel surface absorptivity. [-]

* $\varepsilon_f$ is the flame emissivity. [-]

* $\sigma$ is the Stefan-Boltzmann constant, $\sigma = 5.67 \cdot 10^{-8} [W/m^2 \cdot K^4]$

* $T_f$ is the flame temperature. [K]

* $h_f$ is the convection heat transfer coefficient between the flame and the surface. [W/m²·K]

* $T_s(t)$ is the time dependent surface temperature. [K]

* $\varepsilon_s$ is the surface emissivity. [-]

This model assumes that the pressure vessel is fully engulfed by the flame. This means that the view factor for the radiation is unity and is therefore not taken into consideration. The convective heat transfer coefficients for a jet fire and a pool fire, and recommended values for the emissivity and absorptivity, are given by Scandpower as (Hekkelstrand and Skulstad 2004):

* $h_{jet \ fire} = 100 [W/m^2 \cdot K]$

* $h_{pool \ fire} = 30 [W/m^2 \cdot K]$

* $\alpha_s = 0.85$

* $\varepsilon_s = 0.85$

* $\varepsilon_f = 1.0$ (optical thick flames, thickness > 1 m)
The flame temperature is found by solving equation Eq. 3.16 for the incident heat flux in relation to the ambient conditions. The flame temperature is kept constant throughout the simulation:

$$q_{total} = \sigma \cdot T_f^4 + h_f \cdot (T_f - T_{amb}) \quad (3.16)$$

* $q_{total}$ is the incident flame heat flux as given in table Tbl. 3.2. [W/m²]

* $T_{amb}$ is the ambient temperature ≈ 293 K (20° C)

The heat flux used to calculate the flame temperature is given in table tbl. 3.2.

<table><caption>**Table 3.2:** Incident heat fluxes for various fire scenarios given by Scandpower (Hekkelstrand and Skulstad 2004)</caption><thead><tr><th></th><th>Small jet fire [kW/m²]</th><th>Large jet fire [kW/m²]</th><th>Pool fire [kW/m²]</th></tr></thead><tbody><tr><td>Peak heat load</td><td>250</td><td>350</td><td>150</td></tr><tr><td>Background heat load</td><td>0</td><td>100</td><td>100</td></tr></tbody></table>

## 3.4 Vessel geometry

All vessels modelled are assumed of cylindrical shape. The following shapes are available in openthermo all provided by the Python fluids library (C. Bell 2025).

* Flat-end vessel

* ASME F&D

* DIN (28011)

* 2:1 Semi-elliptical

* Hemispherical

Both ASME F&D, DIN and 2:1 semielliptical are variants of a torispherical vessel. See also the document Calculating Tank Volume. The hemispherical ends are half-spheres extending one radius out.
<table><caption>**Table 3.3:** Vessel geometry details. For torispherical tank heads, the following *f* and *k* parameters are used in standards (C. Bell 2025). *f* is the dish-radius parameter for tanks with torispherical heads or bottoms, *k* is the knuckle-radius parameter for tanks with torispherical heads or bottoms</caption><thead><tr><th>Vessel geometry</th><th>f</th><th>k</th></tr></thead><tbody><tr><td>2:1 semi-elliptical</td><td>0.9</td><td>0.17</td></tr><tr><td>ASME F&amp;D</td><td>1</td><td>0.06</td></tr><tr><td>DIN 28011</td><td>1</td><td>0.1</td></tr></tbody></table>

Using the *fluids* library partial volumes, surface area (full and partial) and liquid level (from partial volume) can be calculated and used internally in *openthermo*.

## 3.5 Model implementation

A simple (naive) explicit Euler scheme is implemented to integrate the mass balance over time, with the mass rate being calculated from an orifice/valve equation as described in Sec. 3.2:

$$m_{cv}(i + 1) = m_{cv}(i) + \dot{m}(i)\Delta t \quad (3.17)$$

$$\dot{m}(i) = f(P, T, i)$$

For each step, the mass relief/ left in the vessel is known. Since the volume is fixed the mass density is directly given.

For the calculation methods (isentropic, isenthalpic, isenergetic, etc.), CoolProp allows specifying density and either H, S, or U directly - this is very handy and normally only TP, PH, or TS property pairs are implemented and you would need to code a second loop to make it into am UV, VH, or SV calculation.

In the following, the high-level steps in the solution procedure are outlined for the different calculation types. To illustrate calls to the equation of state in CoolProp, the notation *EOS(D, T)* corresponds to calculation of the state at known density and temperature, for example.
### 3.5.1 Isothermal process

For an isothermal process, the solution procedure for each calculation step is the following with the mass in the control volume being calculated from Eq. 3.17:

$$ 
\begin{aligned}
T(i+1) &= T(i) \\
D(i+1) &= \frac{m_{cv}(i+1)}{V} \\
P(i+1) &= EOS(D(i+1), T(i+1))
\end{aligned}
$$

The notation EOS(T, P) refers to a call to CoolProp either via PropsSI or the AbstractState solving the equation of state for specified temperature and pressure. Different state pairs are applied including combinations of density (D), specific enthalpy (H), specific entropy (S), specific internal energy (U).

### 3.5.2 Isenthalpic process

For an isenthalpic process:

$$ 
\begin{aligned}
H(i+1) &= H(i) \\
D(i+1) &= \frac{m_{cv}(i+1)}{V} \\
P(i+1) &= EOS(D(i+1), H(i+1)) \\
T(i+1) &= EOS(D(i+1), H(i+1))
\end{aligned}
$$

### 3.5.3 Isentropic process

For an isentropic process:

$$ 
\begin{aligned}
S(i+1) &= S(i) \\
D(i+1) &= \frac{m_{cv}(i+1)}{V} \\
P(i+1) &= EOS(D(i+1), S(i+1)) \\
T(i+1) &= EOS(D(i+1), S(i+1))
\end{aligned}
$$
### 3.5.4 Isenergetic process

For an isenergetic process:

$$ 
\begin{aligned}
U(i+1) &= U(i) \\
D(i+1) &= \frac{m_{cv}(i+1)}{V} \\
P(i+1) &= EOS(D(i+1), U(i+1)) \\
T(i+1) &= EOS(D(i+1), U(i+1))
\end{aligned}
$$

### 3.5.5 Energy balance

The general first law applied to a flow process as outlined in Sec. 3.1.2 subject to an explicit Euler scheme is:

$$ 
\begin{aligned}
D(i+1) &= \frac{m_{cv}(i+1)}{V} \\
U_{cv}(i+1) &= \frac{m_{cv}(i)U_{cv}(i) - (\dot{m}(i)H(i) + \dot{Q}(i))\Delta t}{m_{cv}(i+1)}
\end{aligned}
\quad (3.18) $$

The above assumes that mass flow is positive when leaving the control volume and heat rate is positive when transferred to the control volume. $H(i)$ is the specific enthalpy of the fluid in the control volume for a discharging process and it is equal to the energy of the entering stream for a filling process. The heat rate is calculated as outlined in [#Sec:heat].

For the vessel wall a simple heat balance is also made:

$$ T_{wall}(i+1) = T_{wall}(i) + \frac{\dot{Q}_{outer} - \dot{Q}_{inner}}{c_p m_{vessel}} \Delta t $$

where $\dot{Q}_{outer}$ is the convective heat transfer to or from the ambient surroundings from the outer surface of the vessel, with positive values indicating that heat is transferred from the surroundings to the vessel wall. This is either a fixed heat transfer coefficient with a specified ambient temperature (outer surface) or a calculated fire heat load. $\dot{Q}_{inner}$ is the internal heat transfer, either a fixed number or calculated as natural convection (for discharge) or mixed natural and forced convection (filling).

When the temperature profile of the vessel wall is calculated using the 1-D transient conductivity equation the vessel temperature as calculated above represent the mean vessel temperature and addition to this the outer and inner wall temperatures are calculated.
If a fixed $\dot{Q}$ is provided or a fixed overall heat transfer coefficient is provided the vessel wall temperature heat balance is not solved.

The remaining steps are update of temperature and pressure:

$$P(i + 1) = EOS(D(i + 1), U(i + 1))$$

$$T(i + 1) = EOS(D(i + 1), U(i + 1))$$

## 3.6 Multicomponent mixtures

Although not an initial requirement, the code can handle multi-component gas mixtures. When calculations are done for multicomponent mixtures, the code runs significantly slower. Please note, that in case that liquid condensate is formed during discharge calculations and even if the calculations does not stop, the results cannot be considered reliable. This is because component balances are not made and it is always assumed that the discharge composition is the same as the global composition inside the vessel. When liquid condensate is formed the composition of the vapour phase will differ from the global composition. Further, the only gas properties are estimated during multicomponent calculations, and no liquid heat transfer is considered, although a liquid phase may be present. This will make the estimation of wall temperature non-conservative if low temperature excursion is investigated.

There are a few examples of multicomponent mixtures included with HydDown. In order to specify multicomponent mixtures the below example can be used as guidance:

1 "Methane[9.01290129e-01]&Ethane[6.35063506e-02]&N2[7.80078008e-03]&C02[2.34023402e-02]&Propane[3.50035004e-03]&Butane[5.00050005e-04]"

For component names please refer to the Coolprop manual.

An example calculation for the above mixture is shown in Fig. 3.4.

The pressure and temperature trajectory is visualised along with the fluid mixture phase envelope in Fig. 3.5. As seen from the figure, this case is borderline and the pressure/temperature trajectory just coincides with the dew line on the phase envelope. This plot is included in the example main script that comes with HydDown and serves as an important quality control.

For multi-component mixtures CoolProp does not provide solver for PH and UV-problems and these solvers have been included in HydDown by wrapping the CoolProp PT solver with a Nelder-Mead algorithm for solving for internal energy and density or pressure and enthalpy.
<figure><img src="image_7.png" /><figcaption>**Figure 3.4:** Gas discharge of multicomponent mixture</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>**Figure 3.5:** Multicomponent mixture phase envelope and pressure and temperature trajectory of the vessel inventory during discharge</figcaption></figure>
<figure><img src="image_9.png" /></figure>
# 4 Validation

The code is provided as-is. However, comparisons have been made to some experiments from the literature.

The following gases and modes are considered:

* High pressure hydrogen discharge

* High pressure hydrogen filling

* High pressure nitrogen discharge

* High pressure multi-component gas discharge

These cases are all for Type I (steel) cylinders and the simple heat transfer / heat balance model has been applied, ignoring the thermal gradient in the cylinder wall. It has been checked by applying the 1D transient heat conduction model that the results are basically unaffected by this assumption. In addition to these, a few cases for Type IV cylinders have also been included in order to validate the 1D transient heat conduction model applied. This includes

* High hydrogen discharge - comparison with commercial software.

* High pressure helium discharge experiment (Karlshuhe Institute of Technology)

* High pressure hydrogen discharge (GasTeF experiment)

Details on Type III experiments are sparse but some general observations are made for the difference between type III and type IV cylinders and a single type III filling experiment is modelled.

Furthermore, HydDown has also been benchmarked against external code by Ruiz and Moscardelli (Ruiz Maraggi and Moscardelli 2023), who compared their GeoH2 app against HydDown and found excellent agreement for both pressure, mass flow rate and gas temperature for three different cases of discharge and filling.

## 4.1 Hydrogen discharge (Type I)

Calculations with HydDown for high pressure hydrogen vessel discharge have been compared to three experiments from (Byrnes, Reid, and Ruccia 1964). Byrnes et al. (Byrnes, Reid, and Ruccia 1964) have reported a number of rapid depressurization experiments using both nitrogen and hydrogen, with
blowdown times being varied between 14 seconds to 33 minutes. The vessel used was a type "K" gas cylinder with a volume of 0.0499 m³, an ID of 8.56 inches and shell length of 55 inches. The vessel is placed in a vertical position. A control valve is mounted on the neck of the cylinder with a fixed opening percentage during each experiment, set to provide the required depressurization rate. It is mentioned that the gas cylinder is placed inside a weather proof cabinet with perlite insulation in order to emulate adiabatic conditions. However the simulations run here are with a finite heat transfer coefficient on the outside of the cylinder. In their paper, detailed pressure and temperature traces are provided for three experimental runs: 7, 8, and 9, with depressurization from 2000 psia in 30 seconds, 480 seconds and 14 seconds, respectively. These experiments are simulated in Fig. 4.1, Fig. 4.2 and Fig. 4.3.

<figure><img src="image_10.png" /><figcaption>**Figure 4.1:** Calculations of hydrogen discharge experiment run 7 from (Byrnes, Reid, and Ruccia 1964). The figure shows calculated gas and wall temperature (full lines) compared to experiments (upper left), calculated and experimental pressure (upper right), specific thermodynamic state variables (lower left), and the calculated vent rate (lower right).</figcaption></figure>

The figure layout is general for all the validation studies conducted and will be described in brief for eased readability. Each experiment figure is divided in 4 subplots: the upper left shows the simulated bulk gas temperature and vessel wall (average) temperature with measured values included with stipulated lines (if available), the upper right shows the simulated and measured pressure (if available, shown with points), the lower left shows the specific enthalpy, specific internal energy, and the specific entropy of the bulk gas content as a function of time, the lower right shows the simulated mass rate
leaving or entering the vessel.

As seen from Fig. 4.1, Fig. 4.2, and Fig. 4.3, the calculations generally compare well with the experimental results. It appears as though the faster depressurization times are predicted the best, as seen from run 7 and 9. The much slower run 8 does not have as good a fit for the bulk gas temperature. However, on an absolute scale the results are still within 5°C of the experimentally determined gas temperature. The trend in *goodness of fit* could indicate that the heat transfer from the outside of the vessel to the gas lacks some details, which become detectable at the very slow depressurization. For the fast depressurization the heat transfer from the outside is of less importance and the inside heat transfer coefficient may be predicted better.

<figure><img src="image_11.png" /><figcaption>**Figure 4.2:** Simulation of run 8 from (Byrnes, Reid, and Ruccia 1964).</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>**Figure 4.3:** Simulation of run 9 from (Byrnes, Reid, and Ruccia 1964).</figcaption></figure>

## 4.2 Hydrogen filling (Type I)

In order to compare HydDown to experimental values of hydrogen filling operation the experiments reported by Striednig *et al.* (Striednig et al. 2014) are used. Experimental results were obtained using a type I tank (steel) with a volume of 0.0235 m³. The filling of hydrogen of gas into the vessel was carried out with an upstream reservoir kept at 350 bar and the flow was controlled with an electronically controlled dispenser. Vessel details are provided below:

<table><caption>**Table 4.1:** Key hydrogen vessel data (Striednig et al. 2014)</caption><thead><tr><td>Parameter</td><td>Value</td></tr></thead><tbody><tr><td>Mass</td><td>69.5 kg</td></tr><tr><td>Length</td><td>0.61 m</td></tr><tr><td>OD</td><td>0.280 m</td></tr><tr><td>Wall thickness</td><td>0.0129 m</td></tr></tbody></table>
<table><thead><tr><td>Parameter</td><td>Value</td></tr></thead><tbody><tr><td>Density</td><td>7740 kg/m³</td></tr><tr><td>Heat capacity</td><td>470 J/(kg K)</td></tr></tbody></table>

(Striednig et al. 2014) reports many different experiments. In this comparison we will use experiments applying different pressurisation rates / durations using identical initial conditions. The experiments are reported in Fig. 6 and Fig. 7 in (Striednig et al. 2014).

In its current stage HydDown does not directly support specifying a fixed pressurisation rate expressed as a constant rate of change in pressure / pressure ramp rate (PRR). In order to emulate the experiments, an orifice model is used with an artificially high reservoir pressure in order to ensure that the flow is critical during the entire pressurisation. This gives a constant mass rate and translates to a fairly constant pressure ramp rate, although with real gas effects being clearly observable.

The results of the simulations with HydDown and the experiments from (Striednig et al. 2014) are shown in Fig. 4.4, Fig. 4.5, Fig. 4.6. As seen from the figures the experimental gas temperature appears to be simulated very accurately. The simulation with the fastest pressurisation under-predicts the experimental temperature slightly. The main features of the pressurisation temperature profile is a rapid increase in temperature, at the beginning of the pressurisation, where the effects of filling dominate. At approximately 20-30 s, the temperature increase levels off to a plateau, where the length of the plateau depends on the pressurisation rate. The slower the depressurization, the longer the plateau. During this temperature plateau, the heat transfer from the gas to the vessel and surroundings balance the temperature increase that would occur if performed at adiabatic (isolated) conditions. When the gas flow has ceased upon reaching the final pressure, the temperature starts to drop due to cooling of the gas by the colder vessel wall driven by convective heat transfer. The temperature stabilises at a temperature which is higher than the initial ambient temperature due to the heating of the steel vessel. The cooling of the vessel by ambient air is slower and would require a longer run time to be clearly visible.
<figure><img src="image_13.png" /><figcaption>**Figure 4.4:** Simulation of *H*₂ pressurization using 5 MPa/min (Striednig et al. 2014).</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>**Figure 4.5:** Simulation of *H*₂ pressurization using 10 MPa/min (Striednig et al. 2014).</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>**Figure 4.6:** Simulation of H₂ pressurization using 30 MPa/min (Striednig et al. 2014).</figcaption></figure>

## 4.3 Nitrogen discharge (Type I)

Calculations with HydDown are compared to experiment I1 from ref. (Haque et al. 1992). The experiment is a blowdown of a vertically oriented cylindrical vessel with flat ends. The vessel length is 1.524 m, the inside diameter is 0.273 m, and the wall thickness is 25 mm. The vessel is filled with N₂ at 150 bar and 15°C. Ambient temperature is 15°C. The blowdown orifice diameter is 6.35 mm. The results are shown in Fig. 4.7. The discharge coefficient of the orifice has been set to 0.8 in order to match the vessel pressure profile. The back pressure is set to atmospheric conditions.

As seen from Fig. 4.7, the calculations compare well with the experimental results. The calculated temperature of the bulk vapor is within the experimental range of measured temperature at all times during the simulation. It is also noted that the minimum temperature is reached at approximately the same time as in the experiments. The calculated vessel inner wall temperature does not decline as rapidly as the experiments, but from around a calculation time of 60 s, the temperature is within the experimentally observed inner wall temperature.
<figure><img src="image_16.png" /><figcaption>**Figure 4.7:** Calculations of nitrogen discharge emulating experiment I1 from (Haque et al. 1992). The figure shows calculated gas and wall temperature (full lines) compared to experiments (upper left), calculated and experimental pressure (upper right), specific thermodynamic state variables (lower left), and the calculated vent rate (lower right).</figcaption></figure>

**4.4 Multi-component discharge (Type I)**

The experiments were conducted by Haque *et al.* (Haque et al. 1992) and the experimental data has been extracted from the Ph.d. thesis of Wong (Wong 1998).

The experimental vessel is a full-size suction scrubber for a gas compressor which has a total volume of 2.78 m3, with length and inside diameter of 3.240 m (2.75 m tan-to-tan) and 1.130 m respectively. The
wall thickness is 59 mm. For modelling the length of a flat-ended cylinder has been adjusted to give a total volume of 2.78 m3. An orifice diameter of 6.35 mm was used for discharge and for HydDown modelling the discharge coefficient was adjusted to 0.97 to match the experimental pressure. The results of HydDown and comparison with the experimental results are shown in fig. 4.8.

<figure><img src="image_17.png" /><figcaption>**Figure 4.8:** Calculations of discharge of a non-condensable mixture (0.91/0.09 molefraction of methane/ethane) emulating experiment I1 from (Haque et al. 1992). The figure shows calculated gas and wall temperature (full lines) compared to experiments (upper left), calculated and experimental pressure (upper right), specific thermodynamic state variables (lower left), and the calculated vent rate (lower right).</figcaption></figure>

The minimum gas temperature simulated by HydDown is -6°C, compared to the minimum measured average gas temperature of -10.3°C. The difference in measured vs calculated wall temperature is 2.3°C, with a lower measured wall temperature.
## 4.5 1-D transient heat transfer

A few notes about vessels with poor thermal conductivity and composite materials. When modelling systems with high Biot number and in particular composite materials the complexity increases significantly. Not just because of the more difficult numerical problem, but even more so because of uncertainty in key parameters such as thermal conductivity, density and heat capacity. Composite materials such as carbon fibre or glass fibre reinforced epoxy systems can be manufactured in many different ways (fibre orientation etc.) which effects the previously mentioned properties. These properties, in particular the thermal conductivity, will influence the results significantly. If these properties are not accurately informed for the system to be analysed, sourcing data from literature shall be done with caution.

### 4.5.1 Validation against commercial simulation tool

The first validation case is made against the commercial tool Honeywell Unisim Design in dynamics mode. The initial conditions and vessel details are summarised below cf. tbl. 4.2.

<table><caption>**Table 4.2:** Key vessel data and initial conditions</caption><thead><tr><td>Parameter</td><td>Value</td></tr></thead><tbody><tr><td>Length</td><td>6.4 m</td></tr><tr><td>ID</td><td>0.619 m</td></tr><tr><td>Wall thickness</td><td>0.024 m</td></tr><tr><td>Density</td><td>2000 kg/m³</td></tr><tr><td>Heat capacity</td><td>962 J/(kg K)</td></tr><tr><td>Thermal conductivity</td><td>0.5 W/(m K)</td></tr><tr><td>Outside h</td><td>2 W/(m² K)</td></tr><tr><td>Initial pressure</td><td>182 bar</td></tr><tr><td>Initial temperature</td><td>6 °C</td></tr><tr><td>Discharge mass flow</td><td>0.02 kg/s</td></tr><tr><td>Gas</td><td>Hydrogen</td></tr></tbody></table>

The comparison between HydDown results and the corresponding simulations using Unisim Design Dynamics are shown in Figure 4.9: As seen from the results the results obtained using HydDown
closely resembles the Unisim Design results. The outer wall temperature is matched perfectly and the final gas temperature and inner vessel temperature deviates marginally by 1.5 and 1.9 °C, respectively. Unisim predicts a slightly lower gas temperature and HydDown predicts a slightly lower inner wall temperature.

<figure><img src="image_18.png" /><figcaption>**Figure 4.9:** Calculations of vessel wall temperature (inner/outer) with 1D transient heat conduction during hydrogen discharge for comparison with simulation performed with Unisim Design Dynamics. The figure shows calculated gas and wall temperature (lines) compared to Unisim simulations (left) with Unisim results shown with points, calculated and Unisim simulation pressure (denoted "Experimental") (right).</figcaption></figure>

### 4.5.2 Validation against KIT experiment (Type IV)

Using data from Molkov *et al.* [Dadashzadeh and Molkov (2017)](Molkov et al. 2021) experiments performed at HYKA-HyJet research facility at Karlsruhe Institute of Technology (KIT) for a 19 liter type IV pressure vessel are simulated. The storage vessel was initially charged to 700 bar with helium gas and then cooled down to a normal room temperature (293 K) before start of the depressurisation. Tank characteristics were not available and the required parameters were extracted from a similar tank by Molkov *et al.*, see references within refs. [Dadashzadeh and Molkov (2017)](Molkov et al. 2021). Discharge was through 1 mm nozzle and a discharge coefficient of 0.9 was applied in HydDown as in the study by Molkov *et al.*. The initial conditions and vessel details are summarised below cf. tbl. 4.3.
<table><caption>**Table 4.3:** Key vessel data and initial conditions</caption><thead><tr><th colspan="2">Type IV tank</th></tr></thead><tbody><tr><td>External length</td><td>904 mm</td></tr><tr><td>ID</td><td>180 mm</td></tr><tr><td>OD</td><td>228 mm</td></tr><tr><th colspan="2">HDPE liner</th></tr><tr><td>Thickness</td><td>7 mm</td></tr><tr><td>Density</td><td>945 kg/m³</td></tr><tr><td>Heat capacity</td><td>1584 J/(kg K)</td></tr><tr><td>Thermal conductivity</td><td>0.385 W/(m K)</td></tr><tr><th colspan="2">CFRP shell</th></tr><tr><td>Thickness</td><td>17 mm</td></tr><tr><td>Density</td><td>1360 kg/m³</td></tr><tr><td>Heat capacity</td><td>1020 J/(kg K)</td></tr><tr><td>Thermal conductivity</td><td>0.5 W/(m K)</td></tr><tr><th colspan="2">Initial conditions</th></tr><tr><td>Outside h</td><td>8 W/(m² K)</td></tr><tr><td>Initial pressure</td><td>700 bar</td></tr><tr><td>Initial temperature</td><td>20 °C</td></tr><tr><td>Gas</td><td>Helium</td></tr></tbody></table>

In HydDown the type IV pressure vessel geometry as assumed that of a flat ended cylinder, with the length adjusted to give a total inventory volume of 19 liter. In order to simulate the system a lumped heat capacity and density for the type IV cylinder is made. The thermal conductivity is set to that of the liner material, and the outer heat transfer coefficient is set to 8 W/(m²K) as applied in the H2fills software (Kuroki et al. 2021). For a type IV vessel, the application of an average value for density, heat capacity and thermal conductivity can be an acceptable approximation since the liner and composite shell material are not too dissimilar.
<figure><img src="image_19.png" /><figcaption>**Figure 4.10:** Calculations of vessel wall temperature (inner/outer) with 1D transient heat conduction during helium discharge for comparison with KIT experiments. The figure shows calculated gas and wall temperature (lines) compared to experimental gas temperature (left) and calculated and measured pressure (right).</figcaption></figure>

HydDown simulations are compared with the KIT experiment in Figure 4.10. As seen from the results the depressurisation pressure is matched very well. The experimental gas temperature is also matched fairly well. The lowest simulated gas temperature is 182.3 K and the lowest experimental temperature is 177.5 K a difference of 4.8 K. The minimum gas temperature occurs at an earlier time compared to the experimental results. The former at approx. 75 sec. and the latter at approx. 100 sec. This is also in agreement with the simulation model presented by Molkov *et al.* (Dadashzadeh and Molkov 2017). The final gas temperature of the simulations is 237 K compared to the experimental value of 216 K.

Molkov *et al.* (Dadashzadeh and Molkov 2017) also made a thermal analysis of the thermocouple arrangement used in KIT experiment in order to estimate thermal lag in the temperature measurement. Incorporating the thermal model of the thermocouple arrangement displayed an improved prediction of the time of the minimum measured gas temperature as well as the final measured temperature.

## 4.6 Validation against GasTeF experiments (Type IV)

The experimental setup is described in more detail in refs. [Acosta et al. (2014)](de Miguel et al. 2015). Different tank types where tested both for filling and discharge including a 29 liter type IV tank and a 40 liter type III tank. In the following the type IV tank is modelled. Key details such as the initial conditions and vessel details are summarised below cf. tbl. 4.4. The geometry and thermal properties of the liner and composite shell material was not informed in the cited references and it is assumed to be similar
to the KIT vessel described in the previous section, with the liner and composite shell thickness scaled to match the total thickness of the GasTeF type IV vessel.

<table><caption>**Table 4.4:** Key vessel data and initial conditions. The mass flow during discharge is constant at 1.8 g/s until the pressure drops below 5 MPa after which it drops</caption><tbody><tr><th colspan="2">Type IV tank</th></tr><tr><td>External length</td><td>827 mm</td></tr><tr><td>ID</td><td>230 mm</td></tr><tr><td>OD</td><td>279 mm</td></tr><tr><td>Liner</td><td>HDPE</td></tr><tr><td>Composite shell</td><td>G&CRPE</td></tr><tr><th colspan="2">Initial conditions</th></tr><tr><td>Initial pressure</td><td>700 bar</td></tr><tr><td>Initial temperature</td><td>25±2°C</td></tr><tr><td>Gas</td><td>Hydrogen</td></tr><tr><td>Discharge rate</td><td>1.8 g/s</td></tr></tbody></table>

The simulations with HydDown and comparison against the GasTeF experiment is shown in Figure 4.11. The experimental setup included several thermocouples mounted in different positions inside the test vessel, both at the center line and nearer the top and bottom. The experimental setup did not include a direct measurement of the internal liner temperature interface towards the gas nor the composite shell. The experimental points for the gas temperature as shown in the figure includes the lowest measured temperatures (near the bottom), the gas temperature measured in the middle and towards the top (highest temperatures). The experiments display significant temperature stratification during discharge experiments.
<figure><img src="image_20.png" /><figcaption>**Figure 4.11:** Calculations of vessel wall temperature (inner/outer) with 1D transient heat conduction during hydrogen discharge for comparison with GasTeF experiments. The figure shows calculated gas and wall temperature (lines) compared to experimental gas temperature.</figcaption></figure>

As seen the prediction of the external surface temperature of the composite shell is matched very well. The average gas temperature as calculated with HydDown matches the temperatures recorded in the lower half of the test vessel. Further, the minimum in gas temperature between 450 and 500 seconds are also matched well. The measured pressure was not reported in ref. (de Miguel et al. 2015)

**4.7 Validation against Type III cylinder discharge experiments**

As shown previously type IV cylinder thermal behaviour is matched fairly well when the liner/composite thermal material properties are lumped. For a Type III vessel which has very dissimilar thermal properties of the liner and the composite shell, this approach is not advisable. Luckily, the implemented 1D transient heat conduction code *thermesh* allows modelling bi-materials and this is also implemented in HydDown.

Unfortunately, it has not been possible to find reported discharge experiments which include detailed information about the exact liner and shell composite thickness. In lack of a good and complete
validation cases, an example is made using the KIT cylinder in the previous section and replacing the HPDE liner with an aluminum liner.

<table><caption>**Table 4.5:** Key vessel data and initial conditions for Type III cylinder simulations. The mass flow during discharge is constant at 1.8 g/s until the pressure drops below 5 MPa after which it drops.</caption><thead><tr><th colspan="2">Type III tank</th></tr></thead><tbody><tr><td>External length</td><td>904 mm</td></tr><tr><td>ID</td><td>180 mm</td></tr><tr><td>OD</td><td>228 mm</td></tr><tr><th colspan="2">Aluminum liner</th></tr><tr><td>Thickness</td><td>7 mm</td></tr><tr><td>Density</td><td>2700 kg/m³</td></tr><tr><td>Heat capacity</td><td>900 J/(kg K)</td></tr><tr><td>Thermal conductivity</td><td>237 W/(m K)</td></tr><tr><th colspan="2">CFRP shell</th></tr><tr><td>Thickness</td><td>17 mm</td></tr><tr><td>Density</td><td>1360 kg/m³</td></tr><tr><td>Heat capacity</td><td>1020 J/(kg K)</td></tr><tr><td>Thermal conductivity</td><td>0.5 W/(m K)</td></tr><tr><th colspan="2">Initial conditions</th></tr><tr><td>Outside h</td><td>8 W/(m² K)</td></tr><tr><td>Initial pressure</td><td>700 bar</td></tr><tr><td>Initial temperature</td><td>20 °C</td></tr><tr><td>Gas</td><td>Helium</td></tr></tbody></table>

The results of HydDown simulation with the artificial type III cylinder is shown in Figure fig. 4.12. As seen the immediate effect of replacing the HDPE liner with aluminum (compared to Figure fig. 4.10) is a lower temperature drop in the gas phase and a liner temperature which approached the gas temperature very closely. This is in general agreement with the observations in ref. (de Miguel et al. 2015).
<figure><img src="image_21.png" /><figcaption>**Figure 4.12:** Calculations of vessel wall temperature (inner/outer) with 1D transient heat conduction during helium discharge for a hypothetical type III cylinder with dimensions similar to the KIT type IV cylinder.</figcaption></figure>

In Figures fig. 4.13 and fig. 4.14 the simulated temperature across the vessel liner/composite wall is shown for the real type IV cylinder from the KIT experiment and for the *artificial* type III cylinder, respectively. As seen there is a significant difference between the two bi-materials. The type III cylinder has almost zero temperature gradient across the aluminum liner, which can be rationalised by the very high thermal conductivity. Thus, the entire thermal gradient is over the composite shell material. These observations are in-line with simulation results from ref. (Melideo et al. 2017).
<figure><img src="image_22.png" /><figcaption>**Figure 4.13:** Calculations of type IV vessel wall temperature profile with 1D transient heat conduction. z=0 is the bonding interface between liner and composite, z<0 is the liner and z>0 is the composite shell</figcaption></figure>
<figure><img src="image_23.png" /><figcaption>**Figure 4.14:** Calculations of type III vessel wall temperature profile with 1D transient heat conduction. z=0 is the bonding interface between liner and composite, z<0 is the liner and z>0 is the composite shell</figcaption></figure>
## 4.8 Validation against Type III cylinder filling experiments

For HydDown type III cylinder filling experimental validation, the work of Dicken and Mérida (Dicken and Mérida 2007) is used. They conducted experiments of filling a 74 liter type III cylinder and also compared measurements with CFD calculations. Details of type III cylinder tested is shown in Table tbl. 4.6. For HydDown simulation the length of a flat-ended cylinder has been adjusted to give a total cylinder volume of 74 liter. The time-dependent filling mass flow has been sourced from ref. (Dicken and Mérida 2007) and scaled in order to match the final experimental pressure.

<table><caption>**Table 4.6:** Key vessel data and initial conditions for Type III cylinder simulations.</caption><thead><tr><th colspan="2">Type III tank</th></tr></thead><tbody><tr><td>External length</td><td>893 mm</td></tr><tr><td>ID</td><td>358 mm</td></tr><tr><td>OD</td><td>396 mm</td></tr><tr><th colspan="2">Aluminum liner</th></tr><tr><td>Thickness</td><td>4 mm</td></tr><tr><td>Density</td><td>2700 kg/m<sup>3</sup></td></tr><tr><td>Heat capacity</td><td>900 J/(kg K)</td></tr><tr><td>Thermal conductivity</td><td>167 W/(m K)</td></tr><tr><th colspan="2">CFRP shell</th></tr><tr><td>Thickness</td><td>15 mm</td></tr><tr><td>Density</td><td>938 kg/m<sup>3</sup></td></tr><tr><td>Heat capacity</td><td>1494 J/(kg K)</td></tr><tr><td>Thermal conductivity</td><td>1.0 W/(m K)</td></tr><tr><th colspan="2">Initial conditions</th></tr><tr><td>Outside h</td><td>8 W/(m<sup>2</sup> K)</td></tr><tr><td>Initial pressure</td><td>93 bar</td></tr><tr><td>Initial temperature</td><td>20.25 °C</td></tr><tr><td>Gas</td><td>H2</td></tr></tbody></table>

Simulation results from HydDown is compared to the experimental data of Dicken and Mérida in Figure
fig. 4.15. As seen from the results the gas temperature calculated with HydDown is generally higher than the experimental results. The final calculated temperature is 74.2°C, compared to the final measured temperature of 68.9°C. The CFD calculations performed by Dicken and Mérida (Dicken and Mérida 2007) was also higher than the experimental value (71.3°C). Another CFD analysis of the investigated system by Hall and Ramasamy (Hall and Ramasamy 2024) revealed a final mean gas temperature of 71.4°C, and a zero-dimensional model by the same authors gave an end temperature of 72.4°C.

<figure><img src="image_24.png" /><figcaption>**Figure 4.15:** Calculations of vessel pressure, gas temperature and wall temperature (inner/outer) with 1D transient heat conduction model during hydrogen filling. Comparison is made against experimental results from Dicken and Mérida (Dicken and Mérida 2007) for pressure and gas temperature. Inner wall (liner) temperature is compared to CFD simulations (Dicken and Mérida 2007)</figcaption></figure>

In lack of measured vessel inner wall temperature (liner) the values calculated at the end of filling by Dicken and Mérida using their CFD model is used for comparison. The CFD simulations revealed very large variations in the heat transfer coefficient and resulting wall temperatures. For comparison with HydDown the average CFD calculated temperature is used. This is taken as the arithmetic average of
all seven positions from front to back of the cylinder. As seen from Figure fig. 4.15 the average value
from HydDown (52.5°C)compares very well with the average CFD results at the end of filling (52.1°C).
However, estimated inner wall temperature spanned temperatures from 26 to 67 °C, with the highest
temperatures recorded opposite to the entry of hydrogen.
# 5 Similar software

To be written
<figure><img src="image_25.png" /></figure>
# 6 Future development of HydDown

To be written
<figure><img src="image_26.png" /></figure>
# 7 Example use cases

To be written
<figure><img src="image_27.png" /></figure>
References

Acosta, B., P. Moretto, N. de Miguel, R. Ortiz, F. Harskamp, and C. Bonato. 2014. "JRC Reference Data from Experiments of on-Board Hydrogen Tanks Fast Filling." *International Journal of Hydrogen Energy* 39 (35): 20531-37. https://doi.org/10.1016/j.ijhydene.2014.03.227.

Andreasen, Anders. 2021. "HydDown: A Python Package for Calculation of Hydrogen (or Other Gas) Pressure Vessel Filling and Discharge." *Journal of Open Source Software* 6 (66): 3695. https://doi.org/10.21105/joss.03695.

———. 2024. "HydDown – User Guide and Technical Reference." HAL open science. https://hal.science/ hal-04858235v1.

Andreasen, Anders, Filippo Borroni, Marcos Zan Nieto, Carsten Stegelmann, and Rudi P. Nielsen. 2018. "On the Adequacy of API 521 Relief-Valve Sizing Method for Gas-Filled Pressure Vessels Exposed to Fire." *Safety* 4 (1). https://doi.org/10.3390/safety4010011.

Andreasen, Anders, and Carsten Stegelmann. 2025. "Open Source Pressure Vessel Blowdown Modeling Under Partial Phase Equilibrium." *Process Safety Progress*. https://doi.org/10.1002/prs.70035.

API. 2014a. *API 520 Sizing, Selection, and Installation of Pressure-Relieving Devices*. American Petroleum Institute.

API. 2014b. "Pressure-relieving and Depressuring Systems, API Standard 521." API Standard 521, Sixth Edition, January. American Petroleum Institute.

Bell, Caleb. 2025. "Fluids: Fluid Dynamics Component of Chemical Engineering Design Library (ChEDL)." https://github.com/CalebBell/fluids.

Bell, Caleb, Pierre Lesouhaitier, VikramGovindarajan, and jeremy rutman. 2024. "CalebBell/Ht: 1.0.8." Zenodo. https://doi.org/10.5281/zenodo.14321704.

Bell, Ian H., Jorrit Wronski, Sylvain Quoilin, and Vincent Lemort. 2014. "Pure and Pseudo-Pure Fluid Thermophysical Property Evaluation and the Open-Source Thermophysical Property Library Cool-Prop." *Industrial & Engineering Chemistry Research* 53 (6): 2498-508. https://doi.org/10.1021/ie4033999.

Bjerre, Michael, Jacob G. I. Eriksen, Anders Andreasen, Carsten Stegelmann, and Matthias Mandø. 2017. "Analysis of pressure safety valves for fire protection on offshore oil and gas installations." *Process Safety and Environmental Protection* 105: 60-68. https://doi.org/10.1016/j.psep.2016.10.008.

Borden, Guy, ed. 1998. *Control Valves – Practical Guides for Measurement and Control*. Instrument Society of America.

Bosch, C. J. H. van den, and R. A. P. M. Weterings, eds. 2005. *Methods for the Calculation of Physical*
*Effects (Yellow Book) – CPR 14E.* 3rd ed. Committee for the Prevention of Disasters.

Byrnes, W. R., R. C. Reid, and F. E. Ruccia. 1964. "Rapid Depressurization of Gas Storage Cylinder." *Industrial & Engineering Chemistry Process Design and Development* 3 (3): 206-9. https://doi.org/10.1021/i260011a004.

Dadashzadeh, Makarov, M., and V. Molkov. 2017. "Non-Adiabatic Blowdown Model: A Complimentary Tool for the Safety Design of Tank-TPRD System." In *Proceedings of the Internationsl Conference on Hydrogen Safety*, 1-18. Hamburg, Germany: International Association for Hydrogen Safety. https://hysafe.info/uploads/papers/2017/186.pdf.

de Miguel, N., B. Acosta, P. Moretto, and R. Ortiz Cebolla. 2015. "The Effect of Defueling Rate on the Temperature Evolution of on-Board Hydrogen Tanks." *International Journal of Hydrogen Energy* 40 (42): 14768-74. https://doi.org/10.1016/j.ijhydene.2015.06.038.

Dicken, C. J. B., and W. Mérida. 2007. "Modeling the Transient Temperature Distribution Within a Hydrogen Cylinder During Refueling." *Numerical Heat Transfer, Part A: Applications* 53 (7): 685-708. https://doi.org/10.1080/10407780701634383.

Eriksen, Jacob Gram Iskov, and Michael Skov Bjerre. 2015. "Analysis of the Application and Sizing of Pressure Safety Valves for Fire Protection on Offshore Oil and Gas Installations." Aalborg University. https://projekter.aau.dk/projekter/files/213880237/Analysis_of_the_application_and_sizing_of>_pressure_safety_valves_for_fire_protection_on_offshore_oil_and_gas_installations._By_PEC_T10_1_F15.pdf.

Geankoplis, C. J. 1993. *Transport Processes and Unit Operations*. Chemical Engineering. PTR Prentice Hall.

Grouve, Walter. n.d. "Thermesh – Finite Element Code for Transient One-Dimensional Heat Conduction in Python." https://github.com/wjbg/thermesh.

Hall, C., and V. Ramasamy. 2024. "Modelling the Conjugate Heat Transfer During the Fast-Filling of High-Pressure Hydrogen Vessels for Vehicular Transport." *International Journal of Thermofluids* 21: 100527. https://doi.org/10.1016/j.ijft.2023.100527.

Haque, M. A., M. Richardson, G. Saville, G. Chamberlain, and L. Shirvill. 1992. "Blowdown of pressure vessels. Pt. II: Experimental validation of computer model and case studie." *Trans. IChemE B* 70: 10-17.

Hekkelstrand, B., and P. Skulstad. 2004. *Guidelines for the Protection of Pressurised Systems Exposed to Fire*. Scandpower Risk Management AS.

IEC. 2011. "IEC 60534-2-1 Industrial-Process Control Valves – Part 2-1: Flow Capacity – Sizing Equations for Fluid Flow Under Installed Conditions." 2nd ed. European Committee for Electrotechnical Standardization.

ISA. 1995. *Flow Equations for Sizing Control Valves, ANSI/ISA S75.01 Standard, ISA-S75.01-1985 (r 1995)*. Instrument Society of America.

Jabardo, J. M. Saiz, E. Fockink da Silva, G. Ribatski, and S. F. de Barros. 2004. "Evaluation of the Rohsenow Correlation Through Experimental Pool Boiling of Halocarbon Refrigerants on Cylindrical
Surfaces.” *Journal of the Brazilian Society of Mechanical Sciences and Engineering* 26: 218-30.
https://doi.org/10.1590/S1678-58782004000200015.

Kuroki, Taichi, Kazunori Nagasawa, Michael Peters, Daniel Leighton, Jennifer Kurtz, Naoya Sakoda, Masanori Monde, and Yasuyuki Takata. 2021. "Thermodynamic Modeling of Hydrogen Fueling Process from High-Pressure Storage Tank to Vehicle Tank." *International Journal of Hydrogen Energy* 46 (42): 22004-17. https://doi.org/10.1016/j.ijhydene.2021.04.037.

Melideo, D., D. Baraldi, B. Acosta-Iborra, R. Ortiz Cebolla, and P. Moretto. 2017. "CFD Simulations of Filling and Emptying of Hydrogen Tanks." *International Journal of Hydrogen Energy* 42 (11): 7304-13. https://doi.org/10.1016/j.ijhydene.2016.05.262.

Molkov, V., M. Dadashzadeh, S. Kashkarov, and D. Makarov. 2021. "Performance of Hydrogen Storage Tank with TPRD in an Engulfing Fire." *International Journal of Hydrogen Energy* 46 (73): 36581-97. https://doi.org/10.1016/j.ijhydene.2021.08.128.

Pioro, I. l. 1999. "Experimental Evaluation of Constants for the Rohsenow Pool Boiling Correlation." *International Journal of Heat and Mass Transfer* 42 (11): 2003-13. https://doi.org/10.1016/S0017-9310(98)00294-4.

Rohsenow, Warren M. 1951. "A Method of Correlating Heat Transfer Data for Surface Boiling of Liquids." Technical {Report}. Cambridge, Mass. : M.I.T. Division of Industrial Cooporation, [1951]. https://dspace.mit.edu/handle/1721.1/61431.

Ruiz Maraggi, Leopoldo M., and Lorena G. Moscardelli. 2023. "Modeling Hydrogen Storage Capacities, Injection and Withdrawal Cycles in Salt Caverns: Introducing the GeoH2 Salt Storage and Cycling App." *International Journal of Hydrogen Energy* 48 (69): 26921-36. https://doi.org/10.1016/j.ijhydene.2023.03.293.

Smith, J. M., H. C. Van Ness, and M. M. Abbott. 1996. *Introduction to Chemical Engineering Thermodynamics*. Fifth. McGraw-Hill.

Striednig, Michael, Stefan Brandstätter, Markus Sartory, and Manfred Klell. 2014. "Thermodynamic Real Gas Analysis of a Tank Filling Process." *International Journal of Hydrogen Energy* 39 (16): 8495-8509. https://doi.org/10.1016/j.ijhydene.2014.03.028.

WikiMedia. 2008,. https://commons.wikimedia.org/wiki/File:First_law_open_system.svg.

Wong, Shan Meng Angela. 1998. "DEVELOPMENT OF a MATHEMATICAL MODEL FOR BLOWDOWN OF VESSELS CONTAINING MULTICOMPONENT HYDROCARBON MIXTURES." Department of Chemical Engineering, University College London. https://discovery.ucl.ac.uk/id/eprint/1317912/1/299471.pdf.

Woodfield, Peter Lloyd, Masanori Monde, and Yuichi Mitsutake. 2007. "Measurement of Averaged Heat Transfer Coefficients in High-Pressure Vessel During Charging with Hydrogen, Nitrogen or Argon Gas." *Journal of Thermal Science and Technology* 2 (2): 180-91. https://doi.org/10.1299/jtst.2.180.
<figure><img src="image_28.png" /></figure>