Sage Journals: Discover world-class research

Abstract

This article addresses a long-term thermal dynamic analysis problem for aircraft fuel tank and the related solution. The new method combines the thermal node network analysis and the finite element numerical calculation based on the OpenFOAM software. Compared to the thermal node network analysis algorithm, the computational precision is significantly improved by numerical calculation on more delicate grid in the solid region. On the other hand, compared to conventional finite element numerical calculation method, the computational efficiency is significantly improved since the liquid level at every moment is determined by the volume integral of every grid along the gravity direction instead of volume of fluid method on the Fluent software and two thermal nodes are used to represent the air and the fuel, respectively, and solve the two-phase heat transfer with the fuel tank by the aid of grid computation of the heat exchange area. Taken a typical fuel tank heated by aerodynamic heat as an example, the results comparing to the volume of fluid calculation by Fluent software show that the computation time is reduced to 10.3%, while the max deviation of the tank wall temperature is kept within 6.5%, which meets the requirements of engineering calculation precision.

Keywords

Thermal dynamic analysis fuel tank aircraft thermal node network OpenFOAM software

Introduction

The demands for thermal management in wide range of air vehicle systems are growing rapidly along with the increased mission power, vehicle survivability, flight speeds, and so on. The heat generated by different equipments must be collected, transported, and rejected to ram air or fuel, two principal heat sinks. In addition, increased emphasis on high stealth performance has led to abolish ram air opening and created greater reliance on fuel as the primary heat sink. The heat returned to fuel tank will result in thermally induced operating limits by its temperature increasing greatly during mission profiles. So, system-level thermal management design is requested and the dynamic thermal analysis of fuel tank temperature along the whole trajectory turns to be more significant.

It is difficult to compute conjugate heat transfer (CHT) between the fuel–air two-phase mixture and the solid region of the fuel tank by finite element numerical calculation algorithm used in the traditional computational fluid dynamics (CFD) field. Han et al.^1,2 proposed a sharp interface model to establish the governing equations in different phase regions. Since it is difficult to catch the exact phase interface in complex case, Ettrich et al.,³ Zheng et al.,⁴ and others^5,6 proposed a new phase-field algorithm for the simulation of multi-phase flows, which assumes that there is a transition layer of finite thickness between two phases and therefore establishes the governing equations in the whole multi-phase region. Most CFD software uses this algorithm called as the volume of fluid (VOF) model, but it costs too much time for the too small time step about 0.01s.

In order to cut down the computation time, Lan et al.⁷ introduced an engineering algorithm based on thermal node network model. In this algorithm, fuel–air mixture and solid region are simplified as thermal nodes network with a certain temperature and quality. The finite difference energy equations are established for these thermal nodes. Obviously, it can only got a very sparse temperature distribution since the thermal nodes network is often much looser than the element grid in CFD. Xu and Zhuang⁸ and Cao et al.⁹ used this method to predict the temperature at each node of the fuel system of a fighter, but they both took the fuel–air mixture as a single phase flow and got very rough results. It showed that the model is too simplified and the computation efficiency is very difficult to meet engineering requirement since too many parameters of the heat balance equations need to be determined by experience.

Concerning the specific demand on predicting the transient temperature progress of the aircraft fuel tank with long-term mission profile, Peng et al.¹⁰ of our research group tried the dynamic thermal node network analysis method on the Sinda/Fluint software to solve the two-phase heat transfer problem of the fuel tank. It showed that compared to the VOF method on the Fluent software, the method can obtain a fairly accurate results while the calculation efficiency is improved greatly at the same time.

In this article, a new algorithm further combining the thermal node network model with the finite element numerical calculation is designed. It uses two thermal nodes with certain quality and temperature to respectively represent the two phases and solve the two-phase heat transfer problem in the fuel tank by aid of grid computation. The solid region is partitioned meshes and the whole dynamic thermal temperature distribution is calculated by the finite element numerical calculation.

Afterwards, the above algorithm is implemented on the OpenFOAM software, which is the free CFD software with open-source code and widely applied in recent years.^11–13 Taken a typical fuel tank heated by aerodynamic heat as an example, the calculation results comparing to the method of VOF by Fluent software show that the computation time is reduced to 10.3% which means the whole computation time for a typical trajectory of the magnitude of a thousand seconds is controlled within a week, while the max deviation of the tank wall temperature is kept within 6.5%, which meets the requirements of engineering calculation precision.

Description to overall algorithm

Because the typical fuel tank with the gravity-driven two-phase flow does not have obvious phase changing, the new algorithm combines thermal network model and numerical algorithm. According to the algorithm, the fuel and air in the fuel tank are abstracted as thermal nodes and solid fields are solved by finite element numerical algorithm.

Since fuel and air in the fuel tank are separately fully mixed by natural convection, convective heat transfer between fluid and solid region can be computed with the temperature of fluid nodes and the reasonable convection heat transfer coefficient. It means that the new algorithm does not calculate fluid region by finite element numerical algorithm. The radiation model of fluid–solid surface cannot rely on the fluid field. So, surface-to-surface (S2S) radiation model is chosen in this algorithm. The flow diagram of the overall algorithm is shown in Figure 1. The time interval for twice update of fluid–solid interface is identified by $Δ t_{update}$ .

Figure 1.

The flow chart of the new algorithm.

There are two sections of computation need to be executed in every step. In the first section, update the heat flux of fluid–solid surface:

Compute the VOF nodes and get the liquid height (phase interface) by gravity direction volume integral.

Choose the convection heat transfer coefficient according to the fluid touch state of fluid–solid surface and compute the convection heat flux on it.

Calculate the radiation heat flux on fluid–solid surface by the S2S model.

Add the radiation and the convection heat fluxes on fluid–solid surface to the total heat flux.

In the second section, calculate the temperature of fluid nodes and solid regions:

Calculate the total heat get by fluid nodes in this $Δ t_{update}$ and update those temperatures.

Set the boundary condition of fluid–solid surface according to the total heat flux and compute the transient temperature field of solid regions.

Compute phase interface

In order to efficiently calculate phase interface for the fuel tank with complex geometry and gravity changing, the new algorithm sorts fluid grid units from low to high and calculates the phase interface according to the volume of the residual fuel. The specific method is as follows:

Generate tetrahedron grid for the air–fuel region.

Calculate the volume and the center point coordinate of every fluid grid unit and sort them from low to high according to their center point coordinate and gravity.

Obtain the volume of the residual fuel.

Count the volume of grid units from low to high until the total volume of the first i grid units more than the volume of the residual fuel. The height of the ith unit is the liquid level and the location of phase interface.

If the size of grid units is under millimeter level, the result of this method can meet the engineering design requirement.

Thermal network model of fluid nodes

Temperature of the air or the fuel nodes can be derived through the heat balance equation. Taking the heat convection and the heat of inflow and outflow into account, the heat balance equation can be written as follows

$\frac{d (mCpT)}{dt} = {(mCpT)}_{in} - {(mCpT)}_{out} + Q$ (1)

where the subscripts in and out imply the inflow and the outflow, respectively.

After the heat transfer during the update time step of $Δ t_{update}$ , the new temperature of the fuel node can be calculated as follows

${T_{f}}^{'} = T_{f} + \frac{Q_{f}}{ρ_{f} \cdot C p_{f} \cdot V_{f}}$ (2)

where T, V, $ρ$ , C_p, and Q are the parameters for temperature, volume, density, specific heat, and the total heat through the fuel–solid coupling interface. The subscript f implies the fuel node.

With the fuel flows out the fuel tank, the air of same volume flows into the fuel tank and mixes with the original air. After the heat transfer during the update time step of $Δ t_{update}$ , the new temperature of the air node can be calculated as follows

${T_{a}}^{'} = \frac{(T_{a} + \frac{Q_{a}}{ρ_{a} \cdot C p_{a} \cdot V_{a}}) \cdot V_{a} + T_{ain} \cdot v_{ain} \cdot Δ t_{update}}{V_{a} + v_{ain} \cdot Δ t_{update}}$ (3)

where the subscript a implies the air node and the v_ain represents the volume flow rate of air into the fuel tank.

The total heat flux from the fuel tank to the air or the fuel node during the update time step of $Δ t_{update}$ can be obtained by integrating the convective heat flux of each face unit of the fluid–solid coupling interface, which can be calculated as follows

$φ_{f / a} = h_{f / a} (T_{f / a} - T_{wi}) \cdot Δ t_{update} \cdot Δ S_{i}$ (4)

where T_wi and ΔS_i are respectively the temperature and the area of the face unit, and h is the convective heat transfer coefficient, which can be assigned by engineering experience formula or numerical analysis of the flow field.

The total heat flux of each face unit of the fluid–solid coupling interface contains two parts of the convection heat flux and the radiation heat flux, which can be written as

$ϕ = ϕ_{f / a} + ϕ_{rad}$ (5)

where $ϕ_{rad}$ is the radiation heat flux of each face unit.

Then, the temperature gradient on the face unit can be set as $ϕ / K$ , where K is the thermal conductivity of the fuel tank.

Extension of thermal network model

A typical aircraft fuel system has one supply tank and several storage tanks for a fighter with single engine.^14,15 The whole flow diagram of the fuel system can be illustrated as one tree structure as shown in Figure 2.

Figure 2.

Typical aircraft fuel system flow diagram.

In order to solve CHT problem for aircraft fuel system with structure of tree, the air–fuel thermal network model is extended. Every fuel tank with one air node and one fuel node is regarded as a thermal node of this thermal network. The supply tank is the root node of the network. Storage tanks are branch node or leaf node according to whether it has an upstream node or not. The air is allowed to flow into the thermal network only at the leaf nodes and the fuel to flow out at the root node. When a fuel tank has run out of fuel, it supplies air to its downstream nodes.

All the fuel tanks of a fuel system have the same gas material and fuel material, which is defined on the root node. The volume flow rate of outflow is defined on the leaf nodes. The flow rate of the branch nodes and the root node is decided by their upstream nodes.

In order to calculate temperature for this thermal network model, calculation of the heat and mass transfer between fuel tank nodes needs to be added to the original calculation process:

Downstream search from leaf nodes to compute the volume flow rate for every tank node.

Calculate the volume and the average temperature of the inflow’s air and fuel for every tank node.

Mix the inflow fluid and the original fluid. Update the temperature and the volume of air node and fuel node for every tank node.

These steps can be carried out after transient thermal analysis of the fluids nodes and the solid region. With the quality of heat and mass transfer in $t ~ t + Δ t_{update}$ , these steps update the fluid nodes’ temperature and volume to $t + Δ t_{update}$ .

Implementation on OpenFOAM software

The new algorithm has three main steps: divide the air–fuel region according to the drawdown; calculate the radiation heat flux on fluid–solid surface; obtain the temperature of the fluid nodes and the solid region. In order to implement these steps, geometrical information of the grid must be calculated beforehand. The overall process is shown in Figure 3.

Figure 3.

Flow chart of the implementation process.

Divide the air–fuel region according to the drawdown

Get the grid information of air–fuel mixture and sort the grid units from low to high according to the gravity. Get the serial number of the grid units on either side of the face unit from the file named as “owner” and “neighbor”; get the serial number of the grid units’ vertexes from the file named “face”; get the coordinate of the vertexes from the file named “point”; calculate the volume and center point coordinate of the grid units and sort them.

Get the information of the fluid–solid surface’s face unit. Get the serial number of the fluid–solid surface’s face unit from the file named as “boundary”; get the serial number of the face units’ vertexes from the file named as “face”; get the coordinate of the vertexes from the file named as “point”; calculate the area and center point coordinate of the face units.

Read radiation information

Set radiation calculating parameter. Set the quantity ceiling of surface bunches in the file named as “viewFactorsDict” and initialize the radiation parameter file named as “Qr.”

Calculate the surface bunches with the application named as “viewFactiosDict” and get the radiation angle factors with the application named as “viewFactiosGen.”

Get the temperature of the fluid–solid surface’s face unit from the file named as “T”; get information for surface bunches and compute their temperature; calculate the effective radiation heat flux with the effective radiation function; get the radiation angle factors from the files named as “F” and “globalFaceFaces”; get the radiation heat flux.

Compute transient temperature field

Set the end time of transient as $t + Δ t_{update}$ in the file named as “controlDict.”

Calculate the volume of the residual fuel and compute the fluid level according to the grid information of air–fuel mixture.

Calculate the heat convection between fluid nodes and solid region according to the grid information of the fluid–solid surface.

Set temperature gradient on the fluid–solid surface in the file named as “T.”

Calculate the volume and the temperature of the fluid nodes at $t + Δ t_{update}$ .

Calculate the temperature field of the solid region with the application named as “chtmultiRegionFoam.”

Dynamic thermal analysis

Geometric model

The 1/2 axis symmetry geometric model of a typical storage fuel tank over the engine inlet is shown in Figure 4. The region respectively labeled as top, horizontal, vertical and bottom from top to bottom represents the thermal insulation layer, in which the top arc is heated by aerodynamic heating and the other lines are heated by the engine inlet air.

Figure 4.

Geometric model of the typical fuel tank.

There are two equipments labeled as X-EQ1 and X-EQ2, which geometric models are simplified as two cylinders with 1.5-mm-thick shell inside the tank and contacted with the tank shell through several fins.

Simulation conditions

The fuel tank is heated from the outside of the thermal insulation layer by the aerodynamic heat of the air flow with the recovery temperature of 1730 K. The natural convection heat transfer coefficients on each inner wall are calculated using the following formula separately. The results are shown in Table 1

$Nu = C {(Gr \cdot \Pr)}^{n}$ (6)

where C and n are the coefficients; $Nu, Gr$ , and Pr are Nusselt number, Grashof number, and Prandtl number, respectively.

Table 1.

Convection heat transfer coefficient (W/m² K).

Top	Bottom	Horizontal	Vertical	Tank/fuel	Tank/air
55.2	137.9	344.8	55	100	20

At the initial time, the tank is full of fuel. The initial temperature of the fuel and the solid region is 303 K. As time goes on, the fuel outflows go into the tank with constant speed and the empty volume is full of air.

Analysis results

This typical fuel tank heat transfer problem is solved by the new algorithm. The results of the temperature rise of the tank shell and the equipment are shown in Figure 5.

Figure 5.

The temperature rise curves of different locations calculated by the new algorithm.

The dynamic thermal characteristic of the fuel tank can be found from these results as follows:

In the first 360 s, as little heat flux through the thermal insulation layer, the fluid and the solid region of the fuel tank do not have obvious temperature rise.

After 360 s, the temperature of the tank shell and the two equipments starts to significantly rise. The top shell of the tank has the quickest temperature rise since it is always exposed to the air. Oppositely, the bottom shell of the tank has the lowest temperature rise since it is always submerged in the fuel. Similarly, the temperature rise of the horizontal shell and the vertical shell has a significantly increase when they suddenly expose to the air with the drop of the liquid level. So, it shows that the fuel in the fuel tank is a very good heat sink since its specific heat is much more than the air.

The results by the VOF method on the Fluent software are shown in Figure 6 and the temperature contrast results of these two methods are shown in Table 2. These results show that the following:

The temperature rise curves of the VOF method and the new algorithm have a good similarity. At 1500 s, the relative temperature difference is less than 6.5%. So, the accuracy of the new algorithm can meet the engineering design requirement.

According to the results of Fluent software, the temperature rise rate of the horizontal shell and the vertical shell also has an sudden increase when they are exposed to the air, but the extend of the rise is less than the result by the new algorithm and the temperature rise curves do not have an obvious turning point. The reason is that VOF method is a phase-field algorithm and there is a transition layer between the two phase fluids. But the transition layer is ignored in order to simplify the computation of the thermal network analysis method in this article.

In terms of efficiency, it takes 168 h by Fluent software and only 17.3 h by the new algorithm, which means the computation time is reduced nearly 90% by the new algorithm.

Figure 6.

The temperature rise curves of different locations by Fluent software.

Table 2.

The temperature contrasts at different moments by these two methods.

Time (s)	Temperature (°C)	Top	Bottom	Horizontal	Vertical	X_EQ1	X_EQ2
600	Fluent	35.7	33.1	33.0	32.6	33.4	33.9
	New method	37.6	32.7	32.1	31.7	33.3	33.8
	Deviation (%)	5.2	−1.2	−2.7	−2.8	−0.4	−0.2
900	Fluent	44.1	35.4	37.7	35.5	40.1	40.9
	New method	46.8	35.0	35.7	33.9	39.6	40.2
	Deviation (%)	6.1	−1.1	−5.4	−4.6	−1.2	−1.7
1200	Fluent	53.3	37.5	48.5	40.1	48.1	49.4
	New method	56.2	37.1	51.0	38.2	47.1	48.3
	Deviation (%)	5.4	−1.0	5.1	−4.6	−2.2	−2.3
1500	Fluent	63.0	40.7	61.3	49.6	56.3	59.0
	New method	66.3	40.4	65.4	53.0	56.0	58.4
	Deviation (%)	5.3	−0.8	6.7	6.9	−0.6	−1.0

Conclusion

Aiming at solving the long-term thermal dynamic analysis problem for aircraft fuel tank, a new method combining the thermal node network analysis and the finite element numerical calculation based on the OpenFOAM software is designed and implemented. A case study results showed that the computational precision is significantly improved compared to the thermal node network analysis method and the computational efficiency is significantly improved comparing to the VOF calculation by Fluent software.

The fuel has a good heat sink effect which can protect the temperature rise of the equipment and the tank shell much lower. When the solid region is exposed to the air in the fuel tank with the drop of the liquid level, its temperature will start to rise significantly.

In this article, only one pair of thermal nodes are used to represent the air phase and the fuel phase respectively to solve the heat transfer problem in the fuel tank by aid of grid computation. If there is more than one cycle of the natural convection among the air or the fuel with the solid region, more thermal nodes method with accurate determination of the heat transfer coefficient need to be studied in future research.

Footnotes

Handling Editor: Ito Kazuhisa

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research,authorship,and/or publication of this article.

Funding

The author(s) received no financial support for the research,authorship,and/or publication of this article.

ORCID iDs

Sujun Dong

Fanchao Meng

Dechun Guo

References

Han

Adams

NA.

Adaptive multi-resolution method for compressible multi-phase flows with sharp interface model and pyramid data structure. J Comput Phys 2014; 262: 131–152.

Han

Adams

NA.

Scale separation for multi-scale modeling of free-surface and two-phase flows with the conservative sharp interface method. J Comput Phys 2015; 280: 387–403.

Ettrich

Choudhury

Tschukin

et al . Modelling of transient heat conduction with diffuse interface methods. Model Simul Mater Sci Eng 2014; 22: 085006.

Zheng

Babaee

Dong

et al . A phase-field method for 3D simulation of two-phase heat transfer. Int J Heat Mass Tran 2015; 82: 282–298.

Mirzabeygi

Zhang

Three-dimensional numerical model for the two-phase flow and heat transfer in condensers. Int J Heat Mass Tran 2015; 81: 618–637.

Ren

Shi

Cai

et al . Influence of secretion on airflow dynamics of mechanical ventilated respiratory system. IEEE/ACM Trans Comput Biol Bioinform 2018; 15: 1660–1668.

Lan

Zhu

Zhao

JQ.

Modeling and simulation of general fuel system model. J Aerospace Power 2014; 29: 1623–1630.

Zhuang

DM.

Research of heat management for aircraft fuel system. J Aerospace Power 2007; 11: 1833–1837.

Cao

Zhuang

Ning

et al . Numerical calculation of fluid network system of fighter. J Aircraft Des 2002; 2002: 37–41.

10.

Peng

Dong

Guo

. A dynamic node thermal analysis method of fuel tank based on Sinda/Fluint software. In: Proceedings of the 2016 IEEE/CSAA international conference on aircraft utility systems, Beijing, China, 10–12 October 2016, pp.979–983. New York: IEEE.

11.

Gandhi

Sathe

Joshi

et al . Two phase natural convection: CFD simulations and PIV measurement. Chem Eng Sci 2011; 66: 3152–3171.

12.

Fontoura

DVR

Matos

Nunhez

JR.

A three-dimensional two-phase flow model with phase change inside a tube of petrochemical pre-heaters. Fuel 2013; 110: 196–203.

13.

Habla

Marschall

Hinrichsen

et al . Numerical simulation of viscoelastic two-phase flows using Open-Foam. Chem Eng Sci 2011; 66: 5487–5496.

14.

Feng

Xiugan

Fuel thermal management system for high-performance fighter aircraft. J Beijing Univ Aeronaut Astronaut 2009; 35: 1353–1356.

15.