The temperature model of the thermal re-radiation model in multilayer insulation systems

Abstract In this manuscript, the temperature model of multilayer insulation (MLI) on a spacecraft solar panel in thermal radiation is considered. The thermal insulation performance of MLI and the stress of the spacecraft are also considered. The thermal equilibrium equation of the surface temperature of MLI spacecraft in the presence of current transfer, and the iterative temperature acquisition algorithm of the front and back surface of solar panel MLI are derived. The Gaussian elimination method is used to invert the high order sparse matrix composed of error coefficients, thus realizing the real-time calculation of thermal re-radiation(TRR). The thermodynamic model was applied to the motion simulation of GPS Block IIR satellites, and the theoretical model was verified based on the actual temperature measurement data of several GPS Block IIR satellites. The simulation results have certain reference value for studying the thermodynamic modeling of MLI spacecraft and TRR perturbation calculation.


Introduction
In recent years, great progress has been made in the fields of gravity, the sea surface measurement, global land surface, and glacier measurement. A large number of scientific measurement satellites support significant progress in these areas. In order to guarantee the accuracy of the above works, the accuracy of the satellite orbit determination would be ensured primarily (Ziebart et al. 2005). One of the problems is to improve the precision of the perturbation models of spacecraft by analyzing and modeling some perturbation terms which have not been given enough attention before.
These effects include solar radiation pressure (SRP), thermal re-radiation (TRR), earth radiation pressure (ERP) and effects due to satellite transmission. These effects are difficult to model because they rely on knowledge of instantaneous radiation flux, satellite attitude, material, as well as optical and thermal properties. Among them, TRR is one of the least examined models because it was of least concern in the past. In fact, it has an important impact on the medium and long term orbital motion of satellites.
According to special relativity, mass and energy are equivalent. Therefore the radiation of the photon has momentum. According to the law of conservation of momentum, the pressure exerted on a spacecraft is equal to the momentum leaving the surface every second in the opposite direction, which is thermal re-radiation force. In order to calculate the magnitude of this force, the direct distribution of the emitted radiation must be known. This is the theoretical basis of TRR perturbation. Rubincam (1987) showed that the reflector embedded on the surface of LAGEOS would explain the observed drag behavior of the satellite, whose thermal power was indeed large enough to disrupt the motion of the satellite. Beutler et al (1994) of the Centre for Orbit Determination in Europe (CODE), Vigue et al. (1994) and Gregorius (1996) of The Jet Propulsion Laboratory (JPL) independently used a set of parameters to describe the effects of solar pressure and thermal effects on the orbit determination of a GPS satellite. Fliegel and Gallini (1996) explained the importance of the spacecraft thermodynamic model and pointed out that it accounts for approximately 5% of the total solar perturbation. It is called the T20 model in the satellite model of GPS Block II. Adhya et al. (2004) further pointed out that the thermal component of the nonconservative force field accounts for about 5% of the solar radiation pressure generally, but in some special satellites, such as the GPS Block IIR, the thermal effect accounts for 10% of the total solar perturbation. Marquis and Krier (2000) outlined the mechanism of partial nonconservative force, including SRP and TRR perturbations of solar panels, radiators, antennas and other spacecraft surfaces covered by multilayer insulation (MLI). Ferrier et al. (2002) of The Centre Nationale des Etudes Spatiale (CNES) described a simplified surface mechanical model based on the perturbation of SRP and TRR in the shape of the Jason1 satellite. Ziebart et al. (2003) studied and further demonstrated some systematic biases in this model. In the modified model, the perturbation calculation is based on the spacecraft temperature difference in the MLI structure. Grey et al. (2014) analyzed the effects of four kinds of disturbing force: SRP, TRR, ERP and antenna thrust on the reduction of GPS satellites orbit prediction, and specifically compared the effects of four kinds of forces. Darugna et al. (2018) established a ray-tracing model for Japan QZS-1 satellite and discussed the deviation from orbit determination caused by ignoring the thermal effect in this model.

The thermal theory of MLI spacecraft
For a spacecraft, the sun is the main source of heat radiation, secondary sources include the earth and the moon. Different spacecraft components absorb different amounts of radiation, depending on whether they face the source and how much they are shielded from other components. In general, some parts of a spacecraft absorb long periods of thermal radiation while others do not due to the attitude constraints. For navigation satellites, such as GPS and BEI-DOU, the +Z antenna must always point to the earth, while the solar panel must always point to the sun. In this case, the total energy absorbed and emitted by each component depends on its material composition and optical properties, so the surface temperature may vary greatly. The outer surface of a spacecraft radiates heat directly proportional to temperature and emissivity (Adhya 2005). This situation provides a way for us to simplify the problem when analyzing TRR perturbation. Since the irradiation area of the solar panel of spacecraft is much larger and longer than those of other components, TRR perturbation of the solar panel of spacecraft can be mainly considered in TRR calculations, and other components with smaller irradiation area can be treated as secondary, or they may even be ignored under some circumstance. Based on the observations and calculations of the Solar and Heliospheric Observatory (SOHO), and considering all the wavelengths of solar radiation, the sun provides an average energy flux of 1367 W m −2 (Frolich and Lean 1998) or 1350 W m −2 (Pisacane 2005). Within an AU, the rate of change in an astronomical unit of distance is about 3.3% per annum due to the eccentricity of the earth orbit. Variation in solar illumination reflected solar radiation from the earth and the amount of heat emitted by the earth within a day can cause large temperature fluctuations in orbiting spacecraft. The uneven heating of different parts of the spacecraft will lead to a large temperature difference in the spacecraft. Accessories, such as the boom and sun sensor, must be protected from high-temperature propulsion thrusters, plumes, and external loads. For these reasons, a multilayer insulator (MLI) has been applied to spacecraft as a thermal isolation system.
The principle of MLI is that if there are two parallel planes, their spacing size is much smaller than the lengthwidth size in a vacuum, and they have the properties of a gray body, the radiant heat transfer heat flux q between them is Where σ is the Stefan-Boltzmann constant (5.6699 × 10 −9 W m −2 K −4 ), A is the surface area of plane 1 or plane 2, T 1 and T 2 are the thermodynamic temperatures of plane 1 and plane 2 respectively, ε 1 and ε 2 are the emissivity of plane 1 and plane 2 respectively (Keller 1971;Min and Guo 1998). When ε 1 = ε 2 = ε, the above formula can be simplified as If a parallel plane with the same size and the same surface emissivity is placed between the two parallel planes, then the radiant heat transfer heat flux q between the two parallel planes is If N reflection screens are placed between the two parallel planes, then the radiant heat transfer heat flux q between the two parallel planes is By comparing equations (2) and (4), it can be seen that the heat radiation rate of the two planes can be reduced to 1 N+1 of original after adding N layer reflection screen of the two parallel planes, which plays an obvious role in heat insulation.

The calculation of temperature in MLI spacecraft
According to Stefan-Boltzmann law: Φ = σT 4 , consider the radiation of objects, the formula can be revised to: Φ = εσT 4 , where T is the absolute temperature of the surface and σ is the Stefan-Boltzmann constant, and ε is the surface emissivity of the object.
In the analysis, suppose there are n(n ≥ 3) surfaces in this simulation model, there are n − 1 inner interfaces between these layers and 2 outer surfaces. The outer surfaces may be exposed to radiation from the sun or the earth. Energy can be channeled through a surface into the satellite energy system. The inner surface can transmit heat to either of the two panels, but some of it will be converted into electricity to power the spacecraft.
The general idea of the temperature model based on the thermal radiation effect of the spacecraft MLI is as follows: 1) List the thermal balance equation of each layer of the spacecraft MLI. MLI surface temperature t 1 , t 2 , . . . , tn, t n+1 to be estimated, where the t 1 is the absolute temperature of the outer surface of layer 1, t 2 is the absolute temperature between the layer 1 and layer 2. . . t n+1 is the absolute temperature of the outer surface of layer n; 2) List the equations between measurement errors e 1 , e 2 , . . . , e n+1 and temperature errors ∆t 1 , . . . , ∆t n+1 , and get the variance matrix D; 3) Calculate the measurement error e 1 , e 2 , . . . , e n+1 and variance matrix D, get the temperature errors ∆t 1 , ∆t 2 , . . . ; 4) When the temperature error in this iteration is larger than the error requirement, the temperature in this iteration can be calculated by adding the estimated value of temperature and the error, and the temperature t 1 , t 2 , . . . , tn, t n+1 can be used as the estimated value of the temperature in the next iteration, return to 2) and proceed the next iteration of the algorithm; 5) When the temperature error less than the error requirement, according to the estimated value of the temperature in the previous iteration t 1 , t 2 , . . . , tn, t n+1 and the temperature error of this iteration process ∆t 1 ,∆t 2 , . . . , ∆t n+1 , get the final temperature of MLI layers by adding them.
According to Fourier's law: energy conducted through the panel per unit time is k(t1−tn+1) ∆x , where t 1 is the absolute temperature of the front surface and t n+1 is the absolute temperature of the back surface, k is the absorption and ∆x is the thickness of the material.
In unit time, for each layer, there is thermal equilibrium as shown in the following equation, for front surface: For surface 1: For surface 2: For surface 3: . . .

For surface N:
Where, αn is the surface absorptivity of layer n, E 1 and E n+1 are the solar irradiance in unit time of outer surface in layer 1 and n respectively, ε 1 and ε n+1 are the surface emissivity of the outer surface in layer 1 and n respectively, σ is Stefan-Boltzmann constant, P is the power draw in unit time, A is the area of the panel. In the above equation, the current is derived from surface 2. Since this term is constant, it does not matter which side of the system this term is in when differentiates both sides of this equation. Suppose the error of each of the above equations are e 1 , e 2 , . . . , e n+1 , we can calculate the equations corrections as follows: Where And The temperature errors can be get as follows: Finally, The final perturbation force is (Adhya 2005): c is the speed of light in the vacuum.

Calculation of matrix D and its inverse matrix
In equation (11), In this matrix, Based on the values of the above matrix elements, D is a matrix in the following form: As it can be seen from the above matrix, it is a sparse matrix whose form is similar to the identity matrix. In the case that MLI has multiple layers, for example n ≥ 4, it is a high-order sparse matrix. In order to realize the realtime calculation of temperature, it is necessary to invert the above matrix quickly. Based on the analysis of the inverse matrix method, the fast inverse of the high order sparse matrix whose form is similar to the identity matrix can be divided into two steps. The first step is to confirm that it is a non-singular matrix, which is equivalent to confirm that it does have an inverse matrix. The second step is to use the Gaussian Elimination method to quickly reverse the matrix.
According to the Gaussian elimination method, d i,i−1 (2 ≤ i ≤ n) and need to be transfer to zero in the matrix while .

Test in GPS Block IIR satellites
The satellite data used for model simulation was acquired from GPS Block IIR. In this type of satellite, current is derived from the interior or outer surface of solar panels. When the current is derived from the solar panel, its rated power is 90. According to the attitude requirements of GPS Block IIR, its solar panel always faces the sun, so the direction and area of vertical irradiation are basically constant, which will lead to stable TRR perturbation force of the MLI due to the temperature difference.
The Table 1 (Fliegel and Gallini 1996) is the area and optical parameters of the main components of GPS Block IIR. It can be seen from the table that the total area of solar panels among the components with high reflectivity are much larger than the total area of panel beams, and the reflectivity of other components are extremely low.
In addition to solar panels, almost all surfaces of the GPS Block IIR bus is also covered by MLI, including the plume shield and the radome. In general, the GPS Block IIR is a typical MLI structure satellite. The large structure of GPS Block IIR surface is W-sensor high and low band antenna and NAP UHF antenna. They are composed of thin cylindrical aluminum components covered in black tape (Fliegel and Gallini 1996). Since aluminum has high thermal conductance, it can be considered that the temperature distribution of the axis of each cylinder is uniform. In this case, since there is no temperature difference, the TRR perturbation force generated by the main part of the satellite is equal to 0. Therefore, only TRR perturbation of solar panels is considered, while the force on the spacecraft body is ignored.
We used the thermal characteristics of the GPS Block IIR solar panel in this experiment (Fliegel and Gallini 1996): The component materials and the initial temperatures are (Fliegel and Gallini 1996, the initial temperature in the table is obtained by random estimation) (Table 3).

Test results
The following table shows the final temperature calculation results. If the specification maximum error is less than 0.0005 K, the iteration is processed for six times. The final temperature of the front layer is 318.8658 K and the final temperature of the back layer is 311.1424 K. The perturbation force of TRR is 9.393705227 × 10 −7 N.
The following table shows the temperature errors calculated in the six iterations. As it can be seen from the table, the temperature changes greatly in the initial stage of iteration, and then gradually decrease until it shrinks to within the error range. Figure 1 is the matrix diagram in iterative computation. Due to the large span of matrix data, the illustration in this paper adopts a normalized way to show its data relations more clearly. It can be seen from the figure that the temperature error of MLI monolayer has a strong correlation with the temperature of this layer and the temperature of the adjacent two layers, while the correlation with other layers of MLI significantly decreases.
In Figure 2, it can be seen that with the increase of iteration times, the temperature of MLI rapidly converges from the initial estimate of about 500K to around 300K, and the convergence is fast in the early iteration stage and gradually slows down in the later stage.
In Figure 3, it can be seen that with the increase of iteration times, the temperature error rapidly converges from −120 degrees to nearly 0, and the convergence is fast in the  early stage of iteration and gradually slows down in the later stage. It can be seen that even when the initial temperature estimation differs greatly from the actual temperature, the algorithm has the characteristics of rapid convergence of the actual temperature.
If the MLI of GPS Block IIR solar panel is taken as a whole, according to equation (1), it can be seen that the heat flux q is 57.74113937.
The temperature measured by the sensor on the front of the GPS Block IIR solar panel can be obtained in the satellite telemetry data. The measured temperature results are as follows (Adhya 2005): In SVN41, SVN43, SVN44, SVN46, SVN51 and SVN54 satellite, the average temperature is 322.708K, compared with the results of simulation model calculation, the maximum error rate is −2.7%, the minimum error is 0.3%, the average temperature error rate is 1.17%. Therefore, the model accuracy in this paper is enough to support the orbit calculation of general satellites.

Conclusions
A method is proposed to calculate the temperature of MLI through initial values estimation and iterative algorithm. Based on the thermodynamic equilibrium equation of MLI surface temperature, when current is derived from the solar panel to the spacecraft bus, this method deduces the iterative acquisition algorithm of front surface and back surface temperature of MLI of the solar panel and then obtains TRR perturbation force. Compared with the actual temperature data measured by GPS Block IIR, in the case where the initial temperature error of the initial estimation is 200 K, the final temperature error can be controlled within the range of 0.0005 K through six iterations, and the average error is about 1.17%.