Mathematical Modelling of a Propellent Gauging System: A Case Study on PRISMA

Abstract

Propellant gauging is crucial for a spacecraft approaching the end of its lifespan. Current gauging systems for satellites typically have an accuracy rate of a few months to a year at the end of their operational life. Therefore, it is essential to determine the appropriate gauging system for mission operations. This research focuses on modeling the propellant gauging system for PRISMA, an Earth Observation (EO) satellite of the Italian Space Agency. The analysis centers on implementing algorithms that calibrate the remaining propellant mass in the satellite tank using traditional methods such as bookkeeping (BKP) and pressure-volume-temperature (PVT). To enhance accuracy in quantification, an unconventional approach called thermal propellant gauging (TPG) has been considered. Preliminary computations were conducted using data obtained from the PRISMA thermal model to understand the calibration accuracy of the three methods. At the end of its operational life, the BKP and PVT methods exhibited error rates of 4.6% and 4.8%, respectively, in calculating the mass. In contrast, the TPG method demonstrated a significantly higher precision with an error rate of 1.86%. However, at the beginning of the satellite’s operational life, the PVT and TPG methods showed error rates of 1.0% and 1.3%, respectively, while the BKP technique reported an error rate of 0.1%. Based on these findings, it has been concluded that combining the BKP and TPG approaches yields superior results throughout the satellite’s lifespan. Furthermore, the researchers have determined the specific time duration for which each of these distinct approaches can be effectively utilized.

1. Introduction

One of the important parameters to consider when planning the operational life of a satellite is the residual propellant on board. The amount of leftover propellant, especially for commercial satellites, indicates how long the satellite will remain operational. Towards the end of their service life, the available onboard propellant, which requires precise calibration, becomes crucial in determining satellite de-orbiting, relocation, and replacement [1]. Inaccurate estimations of the available propellant had an impact on the famous Apollo moon landing and presented challenges in extending mission life before its end. Precise and cost-effective quantification of propellant always poses a challenging investment aspect, as indicated in NASA’s technological roadmaps [2]. It is obvious that precise propellant gauging techniques onboard are necessary for earth observation, interplanetary, and other commercial satellites. Besides, a better gauging strategy in an application delivers uncertainties of 2–5% that combines the sensor’s inaccuracy, cost, and design constraints. Satellite engineers usually consider relevant margins for the propellant onboard to ensure a continuous operation of around 10–15 years. In addition, the inclusion of accurate gauging strategies will improve the spacecraft’s life. Furthermore, the presence of precise propellant monitoring will assist the mission operators in improving the orbital life by lowering station-keeping corrections and other regular maintenance in emergency situations. Eventually, the attention rises to propellant monitoring from the past decade, and much research is being carried out on this problem. Usually, the satellite industry follows classical techniques such as bookkeeping (BKP), thruster usage accounting, and the pressure−volume−temperature (PVT) thermodynamic method. Lately, an unconventional method, i.e., the thermal propellant gauging (TPG) method has come under the spotlight due to its intriguing approach to quantifying the available propellant onboard [3].

The traditional BKP technique measures the consumed propellant by predicting the flow rate of the thruster. The remaining propellant onboard will be derived by comparing the predicted mass to the loaded initial mass of the propellant. Since it requires data on consumed propellant, the level of complexity is less in gauging the available propellant. Accurate estimation of the thruster’s flow rate and the exact information on initial propellant mass navigates the precise measurement of the leftover propellant. Though the calculation process appears simple, the possibility of uncertainty in the results is significant because the inaccurate flow rate estimation and inaccurate information on the propellant mass at the beginning of life (BoL) lead to greater uncertainties in the overall outcome. To predict the correct flow rate, the data on pulse widths and thruster on-times are necessary. At the same time, thruster temperature, feed pressure, and duty cycle will affect the thruster’s flow rate. Commonly, the consumed propellant, which will be estimated using the BKP strategy, consists of a 2.5–3.5% error [4,5,6]. However, the variation in thruster orifice diameter with respect to time, as well as the efficiency differences between continuous and pulse firing, implicate the error rate and pose criticalities to consider a constant error for the calibrated fuel. Eventually, the efficiency of the BKP technique decreases and becomes the least accurate at the end of life (EoL). In addition, this technique is incapable of quantifying the individual tank’s propellant mass in a pair configuration. All these factors in the classical BKP technique have forced researchers to focus on other methods.

The gas law technique/PVT strategy is another traditional method for calculating propellant volume using temperature and pressure readings from telemetry data, and it follows the ideal gas law in quantification [4,7,8,9,10,11]. The EoL of a spacecraft will be calculated depending on the helium’s conservation (typically) which is a pressurant gas; thus, gas law equations will be applied to it. Like the BKP approach, it does not require great computation efforts, but the quantification strategy is different since it does not require previous firing history. Since the PVT method requires pressure and temperature readings, the precision of the calibrated propellant volume depends on mounted temperature and pressure sensors, the tank’s volume, the tank’s stretch, the pressurizing gas mass, etc. Temperature differences will be ignored during propellant gauging because the tank is assumed to be isothermal. From the literature, it has been revealed that the pressure sensor’s accuracy directly affects the estimated propellant volume [8]. Furthermore, the inaccurate information on pressurant mass during the calibration time also alters the precision of the measured fuel since the mass of the pressurant may change as the propellant mass decreases. This is because the pressurizing gas may undergo re-pressurization, be absorbed by the propellant, or be due to leakage. Eventually, inaccurate information on pressurant mass has a direct impact on the accuracy of the quantified propellant. This results in the accuracy decrement in PVT as the propellant decreases [8,12,13]. The PVT technique’s accuracy decreases with time because the pressure is inversely proportional to volume, which means that the square root of the pressure and sensitivity are inversely proportional to each other.

Despite the classical techniques such as BKP and PVT, an unorthodox strategy known as the thermal gauging method has grabbed the attention of the space industry based on its different approach to measuring the on-board propellant. Like the other methods, it does not require previous firing history or pressurant temperature and pressure readings to estimate the life of the mission but is rather a characteristic of the material (a standalone/independent technique). The propellant tank’s thermal inertia will be determined to derive the quantity of available propellant [14,15,16,17,18]. The temperature rise will be observed to determine the heat capacity of the tank by applying a known heat load; then, the recorded temperature rise values will be compared against the tank’s thermal model given simulation results to estimate the precise mass of the propellant. Anyway, temperature sensors will be placed on the selected locations of the tank to generate some temperature points that will aid in predicting the tank’s temperature. Since it is possible to overcome the drawbacks of the other two classical techniques, the TPG Method necessitates greater computational efforts. In addition, this technique allows the satellite manufacturer to achieve optimum design to place the TPG system with feasible allocation of the parts like a heater, temperature sensors, the tank’s thermal connection, and a suitable temperature limit to apply on the tanks. As time goes on, the precision of the TPG increases because the sensitivity of the temperature rises with respect to tank mass increases as tank mass decreases, making it more promising to gauge tank load at the EoL. The accumulation of errors in the BK and PVT methods makes TPG more sophisticated as the mission’s end of life approaches. The thermal gauging method can measure the individual tank loads in a multitank system and bipropellant systems. The TPG has been implemented for the ABS1A satellite and GEOStartm1 satellites to compare the estimations carried out by different schemes, and the thermal gauging technique was able to estimate accurately with an uncertainty of less than 1 kg [16,19]. However, the nonuniform heat distribution when heat is applied to the tanks and nonuniform distribution of propellant incorporates the temperature gradients that need to be considered as TPG assumes no temperature gradients during the quantification.

Nevertheless, each strategy is unique in propellant-load quantification. BK and PVT are more precise initially but become less accurate as time goes on. Furthermore, TPG is less accurate in the beginning and delivers precise quantification as the spacecraft approaches its EoL. This trend of the three discussed techniques leads to the concept of hybridization [20], in which any of the techniques, in addition to TPG, will be included to achieve a sophisticated system of measuring propellant load throughout the mission’s life. ASTRIUM-SAS, an on-orbit telecommunication satellite, employs the thermal propellant gauging technique and BKP technique to quantify the correct propellant mass onboard [21]. ABS-1A, Lockheed Martin’s telecommunication satellite, used the BKP and PVT approach in conjunction with the thermal gauging method and attained an accurate estimate of available propellant. In addition, the incorporation of this hybrid technique resulted in a 6-month extension to the mission’s life before de-orbiting [16]. The following are the main contributions of the paper.

• Mathematical modelling of propellent gauging systems
• Errors analysis of the propellant gauging systems
• Performing case study on PRISMA satellite.

These things considered, the current research employs the three above-mentioned methods to assess the accuracy of the propellant load estimate on PRISMA, an Earth Observation (EO) hyperspectral satellite of the Italian Space Agency. The open literature reveals that the inclusion of TGM onboard the spacecraft yields an accurate estimate of the propellant load as the mission’s lifetime is decreasing. Initially, the plan was to employ the BK and PVT methods in the mission. However, based on open literature and previous missions’ results, TPG has also been incorporated. The point of focus for the current analysis is the implementation of algorithms that calibrate the residual propellant mass in the satellite tank using BK, PVT, and TPG. The effectiveness of each technique has been parametrized in terms of error modelling, mass estimates, and their comparisons which could be useful for future satellite manufacturers.

2. PRISMA Mission and Mathematical Modelling

On 22 March 2019, the PRISMA research and demonstration mission was launched aboard the VEGA rocket. The HyperSpectral Earth Observer (HypSEO) project, which was a collaboration between ASI and the Canadian Space Agency, served as the cornerstone for the initial concept investigations. PRISMA plays a vital role in the current and future global EO setting for both the scientific community and end users due to its ability to record data internationally with very high spectral resolution across a wide range of wavelengths [22,23]. PRISMA has the capacity to collect, transmit, and store imagery from all panchromatic/hyperspectral bands, covering a total area of 200,000 km² daily across the entire globe. This results in 30 km × 30 km square Earth tiles. The PRISMA mission operates primarily in two modes: primary and secondary. The primary mode involves collecting panchromatic data, while the secondary mode focuses on gathering hyperspectral data from specific targets as requested by end users. The mission is designed to operate in the “background” mode to fully utilize the system’s resources. The PRISMA space segment consists of one modest class spacecraft. Its payload includes a hyperspectral/panchromatic camera equipped with visible near-infrared (VNIR) and short-wave infrared (SWIR) detectors. It comprises a medium-resolution panchromatic camera (PAN) with a resolution of 5 m (ranging from 400 nm to 700 nm) and an image spectrometer with a spatial resolution of 30 m. The image spectrometer can acquire data across a range of spectral bands from 400 nm to 2505 nm (i.e., from 400 nm to 700 nm in VNIR and from 920 nm to 2505 nm in SWIR). The PRISMA hyperspectral sensor uses a prism to quantify the dispersion of incoming radiation on 2-D matrix detectors, enabling the capture of numerous spectral bands from a single ground strip. Numerous works based on PRISMA EO data can be found in the literature [24,25,26,27,28].

The following section particularly describes three propellant gauging techniques, namely BKP, PVT, and TPG, that have been considered in this study for estimating the residual mass of the PRISMA satellite. Additionally, it provides a description of the inputs that will be acquired from PRISMA telemetry data. A typical satellite propellant gauging system comprises several components, including [13,17,29].

• Fuel tank: It is where the satellite stores its propellant.
• Pressure transducers: They measure the pressure of the propellant inside the tank.
• Temperature sensors: They measure the temperature of the propellant.
• Level sensors: They measure the level of the propellant inside the tank.
• Control electronics: They receive data from the pressure transducers, temperature sensors, and level sensors, and process this data to calculate the amount of propellant remaining in the tank.
• Communication interface: It enables the control electronics to transmit the calculated propellant level data to the satellite’s telemetry system.

Together, these components form a complete propellant gauging system that provides accurate and reliable measurements of the propellant remaining in the satellite’s fuel tank. The PRISMA propulsion system is shown in Figure 1.

2.1. Input Data from PRISMA

PRISMA is an Earth observation hyperspectral satellite developed by the Italian Space Agency (ASI) and built at OHB-Italia.

Table 1. Input parameters of PRISMA.

Parameters	Input Values	Values Type	Associated Error	Units
$M_{hyd}^{load}$	53.70	Measured at fueling	$M_{hyd}^a$	Kg
$V_{pipe}$	0.109	Calculated [CAD & supplier data]	$V_{pipe}^a$	L
$T_0$	293.15	Measured at fueling	$T_0^a$	K
$\dot{m}$	1.102 (BoL)	Calculated at initial helium pressure and burnt mass	${\dot{m}}^a$	kg/s
	0.567 (EoL)
$P_{He}$	21.6 (BoL)	Measured at fueling (BoL)	$P_{He}^a$	bar
	11.0 (EoL)	Calculated (EoL)
t	2	Data sampling time	$t^a$	s
$P_0$	21.59	Measured at fueling	$P_0^a$	Bar
$V_{tank}$	103.2	From supplier data	$V_{tank}^a$	L
T	283.15–323.15	Allowed range in orbit	$T_{hyd}^a$	K
$P_e$	10	From supplier data		W

presents the input values and associated errors used for estimating residual mass and conducting error analysis. It is worth noting that the avionics driver can maintain a constant power consumption of 10 W throughout the mission.

2.2. BKP Approach

A mathematical model used for the gauging analysis based on the BKP method is shown Equation (1). As already described in the introduction, the BKP method is based on the accurate logging of the thrust duration and parameters. To evaluate the residual propellant in the tank, the following model (Equation (1)) was considered for the analysis.

$$M_{hyd}^{\prime }=M^{\prime }\left(t_i\right)=M_{hyd}^{load}-V_{pipe}{\rho }_{hyd\_in}-{\sum }_{j=0}^i{\int }_{t_j}^{t_{j+1}}\dot{m}\left(P_{He},M_{hyd}\right)dt$$

(1)

The propellant mass ejected from the satellite is the function of mass flow rate ($m$), which again is dependent on the pressure of helium $(P_{He})$ and its initial mass ($M_{hyd}$). The temperature is not taken into consideration because the pressure is directly measured using a pressure transducer. This model also considers the degradation of thrusters, i.e., an increase in orifice diameter, which affects the burnt mass through a linear dependency in the mass flow rate model.

The residual propellant mass is calculated by subtracting the sum of all the contributions of the ejected mass. In a similar way, errors are also accumulated over time. The mass flow rate through the thruster was calculated using the expression given by the manufacturer, which is as in Equation (2).

$$\dot{m}=\left(-0.00024P_{He}^2+0.02986P_{He}+0.01804\right).\left(1-0.0006564M_b\right)$$

(2)

Moreover, the maximum deviation mentioned by the manufacturer is 9%. The above model gives a mass flow rate through one thruster, but PRISMA satellite’s configuration has two thrusters, and so the mass flow rate needs to be multiplied by two.

Error Analysis for BKP Approach

The model provided in Equation (2) is for calculating the residual propellant mass as the function of parameters namely, , $M_{hyd}^{load}$, $V_{pipe}$, ${\rho }_{hyd}\left(T_0\right)$ and $\dot{m}$ where the mass flow rate is again a function of pressure and burnt mass. To find the errors occurring in the residual propellant mass estimation, it is required to differentiate and derive the error on residual mass for all the parameters.

Since the mass flow rate depends on pressure $(P_{He})$ and burnt mass ($M_b)$, two models were obtained by differentiating mass flow rate for the two parameters ($P_{He}, M_b$) individually, and the expressions are as in Equations (3) and (4).

$$\frac{d\dot{m}}{d{P}_{He}}=\left(-0.00048P_{He}+0.02986\right).\left(1-0.0044551M_b\right)$$

(3)

$$\frac{d\dot{m}}{d{M}_b}=-0.0044551\left(-0.00024P_{He}^2+0.02986P_{He}+0.01804\right)$$

(4)

Burnt mass $\left(M_b\right)$ used in the above expression can be simply calculated by multiplying the mass flow rate $\left(g/s\right)$ with the time burnt. Errors occurring due to different parameters are tabulated as in

Table 2. Error models of BKP method.

	Parameters	Relative Error Model	Absolute Error Models
Constant errors	Propellant loaded mass ($M_{hyd}^{load})$	$M_{load}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{M}_{hyd}^{load}}=1$	$M_{load}^{abs}=M_{load}^{rel}.M_{hyd}^a$
	Volume of the pipe ($V_{pipe})$	$M_{Vpipe}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{V}_{pipe}}=-{\rho }_{hyd\_in}$	$M_{Vpipe}^{abs}=M_{Vpipe}^{rel}.V_{pipe}^a$
	Initial Temperature ($T_0)$	$M_{T_0}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{T}_0}=-V_{pipe}\frac{d{\rho }_{hyd\_in}}{d{T}_0}$	$M_{T_0}^{abs}=M_{T_0}^{rel}.T_0^a$
Variable errors	Mass flow rate ($\dot{m})$	$M_{\dot{m}}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d\dot{m}}=\frac{\delta }{2}\frac{\dot{m}}{\dot{m}}$	$M_{\dot{m}}^{abs}=M_{\dot{m}}^{rel}.{\dot{m}}^a$
	Helium Pressure ($P_{He})$	$M_{P_{He}}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{P}_{He}}=-\frac{d\dot{m}}{d{P}_{He}}\delta$	$M_{P_{He}}^{abs}=M_{P_{He}}^{rel}.P_{He}^a$
	Burnt mass ($M_b$)	$M_{M_b}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{M}_b}=-\frac{d\dot{m}}{d{M}_b}\delta$	$M_{M_b}^{abs}=M_{M_b}^{rel}.M_b^a$
	Sampling time (t)	$M_t^{rel}=\frac{d{M}_{hyd}^{\prime }}{dt}=-\dot{m}$	$M_t^{abs}=M_t^{rel}.t^a$

. Relative error models tabulated in Table 2 are obtained by differentiating for its related parameters, whereas absolute error models were obtained by multiplying with variable associated errors to relative error models.

To find out the combined error due to all parameters, the root sum squared method, which uses the variations to find out the total error, is applied. In this case, two major types of errors are involved viz., constant errors and variable errors, so it is required to find root sum squares separately. Because variable error increases over time, it is cumulated over time and added to the constant error using the RSS method. in addition, it has been quantified as:

$${(M_{hyd}^{cum})}_{i+1}=\sqrt{{(M_{\dot{m}}^{abs})}^2+{(M_{P_{He}}^{abs})}^2+{(M_{M_b}^{abs})}^2+{(M_t^{abs})}^2+{({(M_{hyd}^{cum})}_i)}^2}$$

(5)

For the initial sequence of errors, the cumulative error will be zero. The total error will be calculated by implementing the RSS method to constant errors summing up with the cumulative error for the present sampling time (t), and the expression is as given in Equation (6).

$$M_{hyd}^{Err}=\sqrt{{(M_{load}^{abs})}^2+{(M_{V_{pipe}}^{abs})}^2+{(M_{T_0}^{abs})}^2+{(M_{hyd}^{cum})}^2}$$

(6)

2.3. Pressure−Volume−Temperature Approach

The PVT method is based on the pressure and temperature measurements of the systems, considering its evolution with time due to the blow-down configuration of the system. By measuring the pressure and temperature of the pressurant and hydrazine, respectively, available in the tank and pipe in stationary conditions, volume (change) is calculated, whereas the residual propellant volume results from the difference between the initial volume, expansion volume of the gas, and eventually the thermal conditions of the propellant in the tank and pipes. The model considered for the analysis is as shown in Equations (7) and (8):

$$M_{hyd}^{\prime }={\rho }_{hyd\_T}{V}_{hyd}^{\prime }={\rho }_{hyd\_T}\left(V_{hyd}^0-\Delta V_{ull}-\Delta V_{hyd\_pipe}\right)$$

(7)

This expression can be further expanded as

$$M_{hyd}^{\prime }={\rho }_{hy{d}_T}\left(V_{hyd}^0-\left(V_{ull}^{\prime }-V_{ull}^0\right)-\left(V_{hy{d}_{pipe}}^{\prime }-V_{hy{d}_{pipe}}^0\right)\right)$$

(8)

The volumes given in the above expression can be re-written as in Equations (9)–(11)

$$V_{hyd}^0=\frac{M_{hyd}^0}{{\rho }_{hy{d}_{in}}}=\frac{M_{hyd}^{load}-M_{hy{d}_{pipe}}}{{\rho }_{hy{d}_{in}}}=\frac{M_{hyd}^{load}-V_{hy{d}_{pipe}}^0 {\rho }_{hy{d}_{in}}}{{\rho }_{hy{d}_{in}}}$$

(9)

$$V_{ull}^{\prime }-V_{ull}^0=\left(\frac{V_{ull}^{\prime }}{V_{ull}^0}-1\right)V_{ull}^0=\left(\frac{T_{He}{P}_0}{P_{He}{T}_0}-1\right)V_{ull}^0=\left(\frac{T_{He}{P}_0}{P_{He}{T}_0}-1\right)\left(V_{tank}-V_{hyd}^0\right)$$

(10)

$$V_{hyd\_pipe}^{\prime }-V_{hy{d}_{pipe}}^0=\left(\frac{V_{hy{d}_{pipe}}^{\prime }}{V_{pipe}}-\frac{V_{hy{d}_{pipe}}^0}{V_{pipe}}\right)V_{pipe}=\left(\frac{ {\rho }_{hy{d}_{in}}}{{\rho }_{hyd\_Tpipe}}-1\right)V_{pipe}$$

(11)

The temperature dependent density (${\rho }_{hyd\_T}$) of hydrazine can be written as shown in Equation (12).

$${\rho }_{hyd\_T}={\rho }_{hyd}^{avg}+\frac{d{\rho }_{hyd}^{avg}}{d{T}_{avg}}{T}_{hyd}$$

(12)

Average density (${\rho }_{hyd}^{avg}$) of the hydrazine can be found considering the maximum and minimum operating temperature of the satellite and its respective densities and similarly change in average density ($\frac{d{\rho }_{hyd}^{avg}}{d{T}_{avg}}$) with respect to average temperature can be found. Similarly, ${\rho }_{hyd\_Tpipe}$ can also be calculated considering the temperature of the pipe.

Error Analysis on PVT Approach

The expression of residual propellant mass is a function of different parameters viz., loaded propellant mass ($M_{hyd}^{load})$, the volume of the pipe ($V_{pipe})$, initial temperature ($T_0)$, initial pressure $(P_0)$, the pressure of helium ($P_{He})$, the volume of the tank ($V_{tank})$, the temperature of the helium ($T_{He})$, the temperature of the hydrazine ($T_{hyd})$ and the temperature of the pipe ($T_{pipe})$. So, to find out the errors occurring with residual propellant mass estimation, it is required to partially differentiate and derive for measuring the errors, which are shown in

Table 3. Error models of PVT method.

Parameters	Relative Error Model
	Relative Error	Absolute Error
$M_{hyd}^{load}$	$M_{load}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{M}_{hyd}^{load}}=\frac{{\rho }_{hyd\_T}}{{\rho }_{hyd\_in}}\frac{T_{He}{P}_0}{P_{He}{T}_0}$	$M_{load}^{abs}=M_{load}^{rel}.M_{hyd}^a$
$V_{pipe}$	$M_{Vpipe}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{V}_{pipe}}=\left(\frac{T_{He}{P}_0}{P_{He}{T}_0}-1\right)\left({\rho }_{hyd\_Tpipe}\right)$	$M_{Vpipe}^{abs}=M_{Vpipe}^{rel}.V_{pipe}^a$
$T_0$	$\begin{array}{cc}\hfill M_{T_0}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{T}_0}& ={\rho }_{hy{d}_T}(\frac{M_{hyd}^{load}}{{\rho }_{hy{d}_{in}}{}^2} \frac{d{\rho }_{hyd}^{avg}}{d{T}_{avg}}\hfill \\ \hfill & -\frac{T_{He}{P}_0}{P_{He}{T}_0}\left(V_{tank}-\frac{\left(M_{hyd}^{load}-V_{pipe} {\rho }_{hy{d}_{in}}\right)}{{\rho }_{hy{d}_{in}}}\right)\hfill \\ \hfill & +\left(\frac{T_{He}{P}_0}{P_{He}{T}_0}-1\right)\left(\frac{M_{hyd}^{load}}{{\rho }_{hy{d}_{in}}{}^2} \frac{d{\rho }_{hyd}^{avg}}{d{T}_{avg}}\right)-\frac{1}{{\rho }_{hy{d}_{Tpipe}}} \frac{d{\rho }_{hyd}^{avg}}{d{T}_{avg}})\hfill \end{array}$	$M_{T_0}^{abs}=M_{T_0}^{rel}.T_0^a$
$P_0$	$M_{P_0}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{P}_0}={\rho }_{hy{d}_T}\left(\frac{T_{He}}{\raisebox{1ex}{$P_{He}$}\!\left/ \!\raisebox{-1ex}{$T_0$}\right.}\left(V_{tank}-\frac{\left(M_{hyd}^{load}-V_{pipe} {\rho }_{hy{d}_{in}}\right)}{{\rho }_{hy{d}_{in}}}\right)+\left(\frac{T_{He}{P}_0}{P_{He}{T}_0}-1\right)a_1\right)$	$M_{P_0}^{abs}=M_{P_0}^{rel}.P_0^a$
$P_{He}$	$M_{P_{He}}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{P}_{He}}={\rho }_{hy{d}_T}\left(-\frac{T_{He}{P}_0}{P_{He}^2{T}_0}\right)\left(V_{tank}-\frac{\left(M_{hyd}^{load}-V_{pipe} {\rho }_{hy{d}_{in}}\right)}{{\rho }_{hy{d}_{in}}}\right)$	$M_{P_{He}}^{abs}=M_{P_{He}}^{rel}.P_{He}^a$
$V_{tank}$	$M_{Vtank}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{V}_{tank}}={\rho }_{hy{d}_T}\left(\frac{T_{He}{P}_0}{P_{He}{T}_0}-1\right)$	$M_{Vtank}^{abs}=M_{Vtank}^{rel}.V_{tank}^a$
$T_{He}$	$M_{T_{He}}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{T}_{He}}={\rho }_{hy{d}_T}\left(\frac{P_0}{P_{He}{T}_0}\right)\left(V_{tank}-\frac{\left(M_{hyd}^{load}-V_{pipe} {\rho }_{hy{d}_{in}}\right)}{{\rho }_{hy{d}_{in}}}\right)$	$M_{T_{He}}^{abs}=M_{T_{He}}^{rel}.T_{He}^a$
$T_{hyd}$	$M_{T_{hyd}}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{T}_{hyd}}=\frac{M_{hyd}^{\prime } \frac{d{\rho }_{hyd}^{avg}}{d{T}_{avg}}}{{\rho }_{hyd\_T}}$	$M_{T_{hyd}}^{abs}=M_{T_{hyd}}^{rel}.T_{hyd}^a$
$T_{pipe}$	$M_{Tpipe}^{rel}=\frac{d{M}_{hyd}^{\prime }}{d{T}_{pipe}}=\frac{{\rho }_{hyd\_T} V_{pipe} {\rho }_{hyd\_in } \frac{d{\rho }_{hyd}^{avg}}{d{T}_{avg}}}{{\rho }_{hyd\_Tpipe}{}^2}$	$M_{T_{pipe}}^{abs}=M_{T_{pipe}}^{rel}.T_{pipe}^a$

.

Similar to the BKP approach, the root sum squared method can be implemented to find out the total error causing the deviation in residual propellant mass and is expressed in Equation (13).

$$M_{hyd}^{Err}=\sqrt{{(M_{load}^{abs})}^2+{(M_{V_{pipe}}^{abs})}^2+{(M_{T_0}^{abs})}^2+{(M_{P_0}^{abs})}^2+{(M_{P_{He}}^{abs})}^2+{(M_{Vtank}^{abs})}^2+{(M_{T_{He}}^{abs})}^2+{(M_{T_{hyd}}^{abs})}^2+{(M_{T_{pipe}}^{abs})}^2}$$

(13)

The error obtained by the formula above (Equation (13)) gives the combined error due to all parameters for pressure and volume at one moment.

2.4. Thermal Propellant Gauging Approach

The TPG method is based on estimating the propellant tank’s thermal mass using electrical heaters attached to it and a couple of temperature sensors. The tank design should be such that it is thermally insulated from the spacecraft by using blankets wrapped around the tank to limit IR radiation and poor conductive supports. The heater power must be sufficient for the overall heat leakages to be able to raise the tank temperature in a short time (say a couple of minutes) and to allow for uniform spreading of the heat within the tank. It must be investigated in detail if, under zero-gravity conditions, the heat is being distributed efficiently enough into the liquid propellant; for instance, by using capillary material at the walls. From the slope of the temperature response after heating up, the heat capacity of the system is estimated as follows, which is the sum of the tank and propellant heat capacity. There are three possible ways available to calculate the residual propellant using thermal inertia. They are:

• Newton’s law of heating;
• By the thermal response of the system (in our case, systems can be considered as a combination of tank, propellant, and pressurant);
• Using mass fuel ratios;

In this case, the last one (using mass fuel ratios) was used to estimate the residual propellant mass. The temperature response as function $T\left(t\right)$ when applying a power $(P_e$) to a system with heat capacity $(H_s$) and effective conduction to the spacecraft $C$ can be calculated with the differential Equation (14).

$$P_e-H_s\frac{\mathsf{\Delta }T}{\mathsf{\Delta }t}-C \mathsf{\Delta }T=0$$

(14)

where the heat capacity of the system can be elaborated as in Equation (15)

$$Q_s=H_s\frac{\mathsf{\Delta }T}{\mathsf{\Delta }t}=H_{t\_up}\frac{\mathsf{\Delta }{T}_{t\_up}}{\mathsf{\Delta }t}+H_{t\_dn}\frac{\mathsf{\Delta }{T}_{t\_dn}}{\mathsf{\Delta }t}+H_{prop}\frac{\mathsf{\Delta }{T}_{prop}}{\mathsf{\Delta }t}+H_{He}\frac{\mathsf{\Delta }{T}_{He}}{\mathsf{\Delta }t}$$

(15)

$H_{t\_up}=M_{t\_up}C{V}_{t\_up}$ for upper tank

$H_{t\_dn}=M_{t\_dn}C{V}_{t\_dn}$ for lower tank

$H_{prop}=M_{prop}C{V}_{prop}$ for propellant

$H_{He}=M_{He}C{V}_{He}$ for pressurant

Whereas heat loss can also be elaborated as

$$\begin{gathered} Q_{loss}=C\mathsf{\Delta }T=k_{i{f}_{up}}\left(T_{t_{up}}-T_{i{f}_{up}}\right)+k_{i{f}_{dn}}\left(T_{t_{dn}}-T_{i{f}_{dn}}\right)+C_{t_{up}}\left(T_{t_{up}}^4-T_{mli\_up}^4\right)+C_{t_{dn}}\left(T_{t_{dn}}^4-T_{mli\_dn}^4\right) \\ +C_{mli\_up}\left(T_{mli\_up}^4-T_{ext}^4\right)+C_{mli\_dn}\left(T_{mli\_dn}^4-T_{ext}^4\right) \end{gathered}$$

(16)

$C_{t_{up}}=\frac{\sigma A_{t_{up}}}{\left(\frac{1}{{\in }_t}+\frac{1}{{\in }_{mli}}-1\right)}$ for upper tank

$C_{t_{dn}}=\frac{\sigma A_{t_{dn}}}{\left(\frac{1}{{\in }_t}+\frac{1}{{\in }_{mli}}-1\right)}$ for lower tank

$C_{mli\_up}=\frac{\sigma A_{mli\_up}}{\left(\frac{1}{{\in }_{mli}}+\frac{1}{{\in }_{ext}}-1\right)}$ for upper thermal blanket

$C_{mli\_dn}=\frac{\sigma A_{mli\_dn}}{\left(\frac{1}{{\in }_{mli}}+\frac{1}{{\in }_{ext}}-1\right)}$ for lower thermal blanket

$C\mathsf{\Delta }T$ considers the heat losses to spacecraft due to all heat transfer processes, namely, conduction, convection, and radiation. To find this, a detailed thermal analysis of the propulsion system is required.

$$Q_s^{\prime }=Q_s-Q_{loss}$$

(17)

Using the balance equation, it is possible to find the residual propellant mass by rearranging the terms. Heat capacity to be considered in the balance is shown in the following equation.

$$H_s^{\prime }\frac{\mathsf{\Delta }T}{\mathsf{\Delta }t}=H_{t\_up}\frac{\mathsf{\Delta }{T}_{t\_up}}{\mathsf{\Delta }t}+H_{t\_dn}\frac{\mathsf{\Delta }{T}_{t\_dn}}{\mathsf{\Delta }t}+H_{He}\frac{\mathsf{\Delta }{T}_{He}}{\mathsf{\Delta }t}$$

(18)

Subtracting actual heat capacity (current heat capacity) from true heat capacity (initial heat capacity) gives the actual propellant heat capacity and dividing it with the specific heat constant of propellant gives the remaining propellant mass.

$$M_{prop}=\frac{Q_s^{\prime }-H_s^{\prime }}{PWR}$$

(19)

where $PWR=C{V}_{prop}\ast \frac{\mathsf{\Delta }{T}_{prop}}{\mathsf{\Delta }t}$

$$\beta =\frac{M_{\mathrm{prop}}}{M_{\mathrm{prop\_ini}}}$$

(20)

Error Analysis of the TPG Approach

An error analysis for the TPG approach has been carried out using first-order equations without considering gas compressibility, propellant vapor pressure, pressurant solvability, tank stretch, thermal expansion of the tank, and the propellant. The heat capacity of the system considered for this error analysis is in

Table 4. Error models of the TPG method.

Parameter	Consideration	Error
$Q_s^{\prime }$	$Q_s^{\prime } \frac{\mathsf{\Delta }t}{\mathsf{\Delta }T}$	$d{H}_s=\left(\frac{d{Q}_s^{\prime } }{Q_s^{\prime } }+\frac{d\mathsf{\Delta }t }{\mathsf{\Delta }t }+\frac{d\mathsf{\Delta }Ts }{\mathsf{\Delta }Ts }\right)ß$
$H_s^{\prime }$	$M_{tank} C{V}_{tank}$	$d{H}_s^{\prime }= =\left(\frac{d{Q}_s^{\prime } }{Q_s^{\prime } }+\frac{d\mathsf{\Delta }t }{\mathsf{\Delta }t }+\frac{d\mathsf{\Delta }Ts }{\mathsf{\Delta }Ts }\right)ß$
PWR	$C{V}_{prop}\ast \frac{\mathsf{\Delta }{T}_{prop}}{\mathsf{\Delta }t}$	$dPWR=\left(\frac{d{Q}_s^{\prime } }{Q_s^{\prime } }+\frac{d\mathsf{\Delta }t }{\mathsf{\Delta }t }+\frac{d\mathsf{\Delta }Ts }{\mathsf{\Delta }Ts }\right)ß$

.

The equations given above (Table 4) are the simple model, like the above methods. In this also, the total error was found by the root sum square of errors that is occurred due to different parameters.

2.5. Algorithm Definition

Using three different approaches described above and considering the errors obtained from analyses, an algorithm to measure remaining residual mass was developed. This algorithm works in the following procedural steps. To begin with, it reads the CSV file, which contains the telemetry data obtained from the satellite when it is onboard. As the first step, it checks whether the software status indicates the orbit maneuvering condition during the firing.

• Once the status has been confirmed, it checks whether the two thrusters were working or only one was working during the firing.
• After this check, it reads the date and time of that operation to verify whether this analysis has already been carried out.
• According to the active number of thrusters, it calculates the total mass flow rate from the tank., i.e., in case two thrusters were working, the mass flow rate will be doubled.
• Then it considers the pressure and temperature values, which are in alignment with the above two conditions, and then calculates residual propellant mass and relative errors using the formulae.
• The two propellant masses and errors will be calculated applying two different approaches viz BKP and PVT method.
• Finally, it prints the results into two txt files; one is for the total history of firing with firing time and the other one will have a summary file, which also provides input for the next analysis.
• A graph will be plotted using the data from the summary file after every firing.

The software status, the indication of active thrusters, and the pressure and temperature of the tank can be obtained from telemetry data when the satellite is onboard.

3. Software Implementation

3.1. Analysis of Propellant Gauging

The mathematical model is implemented in high-level numerical programming software ‘OCTAVE’ [30,31] as a script. The code is divided into three phases, as shown in Figure 2.

3.1.1. Input Phase

In the input phase, the data required to do analysis was inserted by reading the text file, which contains the data about the propulsion system.

3.1.2. Main Phase

The software architecture consists of a primary script that utilizes sampling time to perform a step-by-step analysis until the propellant is fully processed. In the BKP approach, the condition for completing the propellant is used in a while loop, which performs calculations for each time step of 2 seconds according to the given condition. On the other hand, the PVT method conducts computations for uniform pressure and temperature intervals, taking into account the maximum and minimum values. In the main phase, several significant functions are employed. One of them is the “Polyfit” function, which is utilized to determine the average temperature ($T_{avg}$)) and the rate of change of density concerning the average temperature $\left(d{\rho }_{hyd}^{avg}\right)/\left(d{T}_{avg}\right)$, considering the maximum and minimum temperatures. Specifically, the function “$( \mathrm{p}=polyfit\left(\mathrm{x},\mathrm{y},\mathrm{n}\right)$” provides the coefficients for a polynomial, p(x), of degree n, which serves as the best fit (according to the least-squares method) for the data in y. The coefficients within p are arranged in decreasing powers, and the length of p is n + 1. Here, p(x) is represented as

$$p\left(x\right)=p_1{x}_n+p_2{x}_{n-1}+\dots +p_{n}x+p_{n+1}$$

(21)

3.1.3. Output Phase

In the output phase, plots have been generated to analyze various aspects of residual propellant and errors in quantifying it. These plots aim to determine the rate of propellant usage and the duration for which it lasts. The next section presents the plotted results.

3.2. Telemetry Data Analysis for Propellant Gauging

The OCTAVE software was utilized to analyze telemetry data for measuring residual propellant gauging. It generates error graphs of BKP, PVT, and TPG in a single figure as output.

4. Results and Discussion

Since overall research has been carried out on the 3 propellant gauging strategies, a relevant discussion is incorporated on each one in this section accompanied with obtained results. And this explanation includes the rationalities which are used to choose the optimum techniques among the 3 for the PRISMA spacecraft.

4.1. BKP Approach

The first three error sources reported in the Table 2 are constant errors during the mission whereas the other quantities are evolving over time and are additive (after each burn the error must be computed and summed to the previous error in a cumulative way). The software OCTAVE 4.4.1 was used to estimate the residual mass and the total error accumulated depending on parameters considering the sampling time as 2 s. In the following graph, the result of the analysis was plotted, assuming a temperature of 20° (293.15 K) and for a single firing with a sampling data time of 2 s. It was observed from the trend plot that the mass of 53.70 kg will last for 64,456 s as shown in Figure 3.

In Figure 4a,b, the error for the same conditions was computed and plotted with respect to firing time. It can be observed that errors are increasing over time and accumulated to maximum at EoL. At the BoL, the error is limited and strictly related to the uncertainty in the loaded mass (0.1 kg) plus the trapped mass due to the uncertainty in the piping volume (0.008 kg) and the initial temperature through the hydrazine density model (0.00004 kg). As the mission evolves, the mass uncertainty increases after every burn. At EoL, the uncertainty reaches a value of 2.50 kg, about 4.6% of the initially loaded mass.

4.2. PVT Approach

The uncertainty in the input parameters has been expanded in partial derivatives, whereas the numerical values depend on instantaneous parameters. The partial derivative values are provided for BoL and EoL conditions (21.59 and 11.0 bar), both at a standard temperature of 293.15 K (system assumed to be isothermal).

Figure 5 depicts the behaviour of the residual usable hydrazine mass provided (in the case of an isothermal system). The family of pressure-mass curves is parameterized according to the varying temperatures. Temperature ranges considered here are 278.15–323.15 K.

The absolute error in the estimation is provided in the following plot, as shown in Figure 6.

As a reference, the error in the estimation for a given case (20 °C) is provided in the Figure 7. The driving error source is constituted by the quantified pressure of the system, which with a precision of 1.2% in the acquisition chain, creates an error between 0.69 and 2.6 kg (4.6%).

4.3. Thermal Propellant Gauging Approach

For the propellant’s density to remain constant, its temperature must be maintained at 293.15 K, and the heater must be controlled accordingly. When the firing is OFF, and the heater is ON, thermal propellant gauging can be conducted. The thermal model already developed was helpful in modelling the thermal propellant gauging considering all losses. The residual mass for increasing time (Temperature slopes) is shown in Figure 8. The graph below was obtained by using the equations shown in Section 2.4.

Figure 8 shows the propellant mass remaining after every second up to the remaining 5 kg and the temperature slopes were calculated by considering sampling time of 2 s. Figure 9a demonstrates the presence of the error in the above calculation for the whole duration. It is evident that the error in TPG decreases at the EoL. The model considered for plotting the results below is a linear error model, and the model could be further improved by performing the sensitivity analysis. From Figure 9b, it is observed that the error in the estimated mass using TPG method is 0.70 kg at the BoL whereas it decreases to 1.02 kg at the EoL when compared to the other two methods BKP (2.5 kg) and PVT (2.6 kg). So, it is suggestible to use these readings after a certain point to increase the life of the satellite, which is explained in the following section.

4.4. Comparative Study of Errors on Different Methods

The estimated and error masses of various methodologies at the BoL and EoL are displayed in

Table 5. Error values of three different methods.

Method	Estimated Error (kg)	Error Percentage (%)
	BoL EoL	BoL EoL
BKP	0.1 2.5	0.1 4.6
PVT	0.6 2.6	1.0 4.8
TPG	0.7 1.0	1.3 1.8

. At the BoL, the error is observed to be 0.1 kg for the BKP method and 0.7 kg for the TPG method. However, at the end, the error is 1 kg for thermal propellant estimation compared to 2.5 kg for the BKP method. In addition, it is evident that the PVT method is subpar compared to the other approaches.

In addition, we have evaluated only the TPG and BKP due to their advantages over the PVT method. Figure 10 depicts the utilization of propellant mass measured with BKP and thermal methodologies from which it is possible to observe the difference in propellant mass between the two distinct methods. During the operational phase, as reported in Table 5, the error of the TPG method decreases over time, whereas the error bars of the BKP method increase.

Figure 11 depicts the duration after which the BKP method error becomes greater than the thermal propellant measurement. After approximately 29,616 s, it is recommended to use thermal propellant gauging, which provides a more accurate estimate of the remaining propellant. By incorporating this hybrid strategy, the mission duration will be extended, and residual propellent will be utilized more efficiently.

The three different models have shown three different results, and those are as follows:

• Firstly, the BKP method is one of the best methods for gauging propellant from BoL up to a certain limit. After this, the error in the BKP values increases since it calculates by a cumulative method. In this case, the error was 0.1 kg (0.1%) in the beginning, and it increased to 2.5 kg (4.6%) till the end of life (EoL).
• Secondly, the PVT method has shown a moderate result (error) which is 0.69 kg (1%) and increased up to 2.6 kg (4.8%) till the EoL. Since it uses two parameters, i.e., pressure and temperature values for calculation, both parameters contribute to the error. The BKP and PVT method are passive which means it is not required to send any commands to satellite especially for this, but the data acquired by default is enough to calculate the residual mass.
• Finally, the TPG method, which was the focus in this research work calculated the residual propellant mass from input and output heat transfer values from heaters and temperature sensors, respectively. This is an active/passive method of gauging which means heaters could be switched ON/OFF manually by giving a command from the ground station if necessary, and the change in temperature with respect to time for the tank should be acquired from telemetry data. Using temperature slopes, the residual mass can be calculated, and the models for this were developed in the project. The error contribution at the beginning was 0.7 kg (1.3%), whereas at the end, it was 1.02 kg (1.86%).
• Usually, the PVT method will be used by most of the satellites, but TPG exclusively developed here is giving better results than the other two methods at end of life. So, it was concluded to use two different methods for gauging the propellant in the satellite tank throughout its life.
• The duration up to which the BKP method can be used was also found. It was observed as 29,616 s up to which the BKP method provides better accuracy, and after this is recommended to shift to the TPG method for the tank fill level.
• Hence, developing the TPG model was more effective and will allow us to save 1.5 kg of propellant at the end of life which can add more ~5 months to the operational life of the satellite.

5. Conclusions

Satellite propellant gauging systems are one of the significant components of the satellite technology since these provide accurate propellant measuring, prevent overfilling, premature shutdown, and enable efficient fuel utilization. These techniques serve a key role in assuring satellite safety, stability, and performance, as well as ensures that satellites can fulfil their intended missions while conserving fuel resources. With the growing demands for space exploration and satellite operations, the requirement for reliable and accurate propellant gauging systems is more vital than ever, and continued improvements in this technology are anticipated to continue to drive innovation and progress in the field. In this research work, different types of propellant gauging methods, namely Bookkeeping (BKP), Pressure−Volume−Temperature (PVT), and Thermal Propellant Gauging (TPG) were analyzed for implementation on the PRISMA satellite and the respective error models were also developed. The challenging task in this study was to develop the TPG method, as it was rarely implemented on satellites because of its complexity in modelling. Using these models, it is possible to measure the residual propellant mass using the obtained telemetry data from the satellite. Most satellites will typically employ the PVT approach; however, the TPG method developed exclusively here, outperforms the other two techniques at the end of life. As a result, it has been decided to employ two distinct ways of measuring propellant in the satellite tank during its lifetime. It was also discovered how long the BKP approach can be used. The calculations have been predicted that the BKP method provides superior accuracy up to 29,616 s, after which it is recommended to switch to the TPG method for monitoring the tank fill level. As a result, building the TPG model was more efficient, allowing 1.5 kg of propellant to be saved at the end of life, potentially adding 5 months to the operational life of the satellite.