Numerical Analyses and a Nonlinear Composite Controller for a Real-Time Ground Aerodynamic Heating Simulation of a Hypersonic Flying Object

Abstract

This paper contains two parts: numerical analyses and a control method. The numerical analyses of a hypersonic flying object’s aerodynamic heating environment are based on three different two-dimensional outflow fields via finite element calculations. Then, the reference temperature trajectories of a hypersonic flying object are obtained. The other one is an intelligent proportional-derivative (IPD) with a nonlinear global sliding mode control (NGSMC) based on a nonlinear extended state observer (NESO) for a real-time ground aerodynamic heating simulation of a hypersonic flying object, named a thermal-structural test with quartz lamp heaters. The composite controller is made of three sub-components: a model free frame that is independent of the system dynamic model along with an ultra-local model; a NESO for the lumped disturbances observation; and an integral sliding mode control with a nonlinear function for the observation errors compensation. The flight environment of the hypersonic flying object is from Mach number 0.6 to Mach number 5.0, with between flight altitude of 31,272 m and flight altitude of 13,577 m. The comparative results demonstrate some superiorities of the proposed composite controller in terms of tracking errors and robustness.

1. Introduction

Nowadays, the hypersonic vehicle [1] is one of new development directions in the aeronautic and astronautic fields because of its hypersonic speed. There are a variety of hypersonic vehicles, such as unmanned aerial vehicles [2], hypersonic missiles [3], space shuttles [4], and other winged/wingless aircrafts [5]. When the speed of a hypersonic vehicle exceeds that of sound, an aerodynamic heating phenomenon [6], also known as “thermal barrier”, occurs, especially if a hypersonic vehicle’s speed is more than Mach number 5.0. This serious aerodynamic heating phenomenon has two characteristics: a fast heating rate and high temperature, which not only burn through the surface of hypersonic vehicles but also affect the lifetime of the electronic devices inside the hypersonic vehicles, even resulting in flight accidents. Therefore, it is necessary to evaluate the rationality of the thermal protection systems (TPS) [7,8] of hypersonic vehicles, called thermal-structural tests.

The thermal-structural test [9,10] is developed to reproduce a real aerodynamic heating environment of hypersonic vehicles in a time sequence on the ground. So far, the heat transfer of thermal-structural test has been divided into three main ways: heat convection, heat conduction, and heat radiation. Among them, the wind tunnel experiment [11] is a kind of typical heat convection, which is realized by the high-speed relative motion between the high-speed airstream and test sample. In [12], hypersonic wind tunnels of the China Aerodynamics Research and Development Center were developed for the hypersonic aerodynamic characteristics acquisition of air-breathing hypersonic vehicles. In [13], four typical driving modes of the high enthalpy wind tunnel are introduced: air-direct-heated hypersonic wind tunnels, light-gas-heated shock tunnels, free-piston-driven shock tunnels, and detonation-driven shock tunnels. However, the flexibility of a wind tunnel’s design size is single and it cannot adapt to multi-size test samples. At the same time, running time is short, leading to limited information acquisition in thermal-structural characteristics. Different from the wind tunnel, non-convective thermal-structural tests, which usually focus on heat radiation to reproduce a real aerodynamic heating environment, are applied widely and flexibly. Common heating elements of heat radiation have many varieties: nichrome wires, silicon carbide rods, graphite rods, and quartz lamps. By contrast, the quartz lamp [14], as an infrared radiation heating element, is used in non-convective thermal-structural test due to its small thermal inertia, long lifetime, antioxidant capacity, high efficiency, and other advantages. In [15], the heat flux distributions of a quartz lamp and its array are studied via the Monte Carlo method in terms of loading power and array height. However, in view of thermal-structural test based on quartz lamp (TSTQL) control systems, its nonlinear and external disturbances lead to some problems: poor robustness, long response time, and low precision in control systems.

To solve the abovementioned problems, some control strategies, which are fuzzy control [16], proportional-integral-derivative (PID) control [17], and iterative learning control (ILC) [18] were developed for TSTQL. In [16], a fuzzy control is used for transient aerodynamic heating simulation, which is produced by the surface of a hypersonic flying object, in a quartz lamp heating device. In [17], a nonlinear PID, along with digital and analog values, are applied for aerodynamic heating tracking of quartz lamp heaters. In [18], a fractional order PD^α ILC is designed for quartz lamp infrared radiation aerodynamic heating simulation experiment. However, fuzzy control is mostly based on empirical formulae with many parameters for tuning, leading to an imbalance between high precision and efficient decision-making; the PID control just relies on tracking errors linear superposition, which poses a threat to rapidity and overshoot. ILC with some prepared calculations cannot track temperature trajectories in real time during its open loop control.

Model-free control (MFC) [19] is a model independent strategy in an ultra-local model frame, consisting of three aspects: a closed loop controller, a lumped disturbances observation, and an auxiliary controller. In general, the model-based control relies on parameters of the complex system dynamic model, but the system dynamic model itself is sometimes changeable and indeterminate. Moreover, in the absence of accurate system dynamic model information, the model-based control may have difficulty tuning parameters with uncertainties. MFC uses an ultra-local model to replace the complex system dynamic model for reducing order and simplifying the system model, and its loop is usually closed by an intelligent PID (IPID). In [19], a theoretical summary of MFC is on the basis of an ultra-local model, in which IPID and estimation techniques in model-free frames have some concrete applications from intelligent transportation systems to energy management as well. Then, an observer is employed for the lumped disturbances observation, such as time-delay estimation (TDE) and extended-state observer (ESO). In [20], an intelligent proportional-derivative neural network model-free control based on TDE is presented for a five-DOFs lower extremity exoskeleton, in which TDE is employed for system uncertainties estimation. In [21], for a permanent magnet synchronous motor (PMSM), a control strategy about intelligent proportional-integral super-twisting nonlinear fractional-order sliding mode model-free control with a linear ESO (LESO), is applied to estimate the unknown terms from the ultra-local model. As such, an IPID and an observer are integrated into an ultra-local model to achieve some goals of reducing system complexity, improving dynamics robustness. In addition, due to some observation errors and measurement noise existing, an auxiliary controller is required for real-time compensation. In [22], an unknown nonlinear system, such as an inverted pendulum plant, is replaced by a one-order ultra-local model with an intelligent proportional-derivative neural network model-free control, in which the neural network control functions as an auxiliary controller and has a function for estimation error compensation. However, neural network control may have more computational load of approximation properties.

Sliding mode control (SMC), as an auxiliary controller of MFC, is a type of the variable structure control strategy, which has two discrete structures: a reaching phase and a sliding phase. It is insensitive to parameter variation and fast dynamic response compared with neural network control functions. SMC is applied in many fields, such as the attitude stabilization of a rigid spacecraft [23], salient-pole wound rotor synchronous generators [24], and unmanned surface vehicles [25]. In [26], according to the cascaded control structure about the speed loop and the current loop from the synchronous reluctance motor drive system, a Hermite neural network-based second-order SMC is able to obtain better tracking performance and stronger robustness. In [27], SMC with a bat algorithm is designed to save energy by minimizing the forces in a two-DOF robot. In [28], a super twisting SMC with a fuzzy logic for the hybrid energy storage system of fuel cell hybrid electric vehicles. At the beginning, the system state is driven towards the sliding surface along with high-frequency switching motions, and then system state is constrained to lie in the prescribed sliding mode surface with convergence to equilibrium state. In [29], an intelligent proportional-integral sliding mode model-free control is proposed to track the maximum power point in doubly fed induction generator (DFIG) wind turbine systems, in which SMC circumvents estimation errors, whereas during the reaching phase, the system dynamics performance is sensitive to internal parametric uncertainties and external disturbances. So, a global SMC (GSMC) [30,31,32,33] without the reaching phase is proposed to suppress chattering phenomena from the high-frequency switching motions. The original system state lies in the prescribed sliding mode surface, avoiding the problem of sensitivity to internal parametric uncertainties and external disturbances. In [30], a GSMC controller is proposed for both variable order fractional systems and constant order fractional systems. In [31], due to the GSMC controller, precise locating of load position and eliminating fluctuations are obtained during anti-sway control of offshore container cranes. In [32], a continuous GSMC with a nonlinear disturbance observer is designed for the translational oscillator with rotational actuator system. In [33], an adaptive backstepping GSMC is proposed to track attitude and position of the quad-rotor unmanned aerial vehicles in finite time, and, at the same time, model uncertainties, external disturbances, and input saturation are under consideration. In [34], an adaptive global terminal sliding mode control without the reaching phase is used for the overload control system of anti-warship flying objects. These references remove the reaching phase to obtain higher performance; by contrast, the conventional first-order SMC has a reaching phase and a sliding phase leading to chattering phenomena from the high-frequency switching motions. Moreover, instead of sign function, a nonlinear function [35] is introduced into a global SMC, which alleviates steady state errors and saturation errors, achieving a goal: smaller errors corresponding to larger gains and larger errors corresponding to smaller gains.

Hence, this paper presents a control strategy: an intelligent proportional-derivative (IPD) with a nonlinear global sliding mode control (NGSMC) based on a nonlinear extended state observer (NESO) for a real-time ground aerodynamic heating simulation of a hypersonic flying object, named a TSTQL heater (TSTQLs). The composite controller has three parts: based on an ultra-local model, the system dynamic model is replaced by a model free frame; the loop of an ultra-local model is closed via IPD and NESO, which is for the lumped disturbances observation; and an NGSMC, as an auxiliary controller, consists of an integral sliding mode and a nonlinear function for the observation errors compensation. As a consequence, the control strategy, a composite controller, is proposed to track temperature trajectories from numerical analyses of a hypersonic flying object.

In this paper, we make some contributions:

• Numerical analyses of a hypersonic flying object’s aerodynamic heating environment are based on three different two-dimensional outflow fields via finite element calculation in ANSYS Workbench 2020 R2.
• The composite controller control strategy is proposed for the TSTQLs, which can provide a model free frame being independent of the system dynamic model along with an ultra-local model.
• The NESO is designed for the lumped disturbances observation and the NGSMC, an auxiliary controller of MFC, combines an integral sliding mode with a nonlinear function to achieve a stage of great tracking errors, fast response time, and strong robustness. Moreover, the NGSMC eliminates the reaching phase, suppressing chattering phenomena from the high-frequency switching motions.
• Instead of sign functions, the nonlinear function is integrated into NGSMC, alleviating steady state errors and saturation errors and achieves a goal: smaller errors corresponding to larger gains and larger errors corresponding to smaller gains.
• The comparative results demonstrate some superiorities of the proposed composite controller in terms of tracking errors and robustness.

2. Thermal-Structural Test with Quartz Lamp Heaters

In Figure 1, there is a general framework of the TSTQLs for the rationality evaluation of the hypersonic vehicle’s TPS, including aerodynamic heating data, control system, and feedback. These steps are explained in details as follows:

Step 1 (aerodynamic heating data): based on the outflow fields of a hypersonic flying object, the real flight environment needs to be added into the finite element calculation, such as flight altitude, ambient temperature, ambient sound velocity, ambient pressure, ambient density, and flight Mach number. Then, transient aerodynamic heating data in time sequence are required, which are regarded as reference temperature trajectories for the next control system’s tracking aim.

Step 2 (control system): the control system is made of five parts: temperature trajectory, control center, power regulator, quartz lamp heaters, and temperature sensor. The temperature trajectories are from the finite element calculation results, called transient aerodynamic heating. In the control center, our proposed controllers are loaded to calculate. The transistor, as a power regular, is a type of silicon-controlled rectifier (SCR), and it can obtain different power values via changing the SCR conduction angle in the AC voltage regulating circuit. Therefore, the quartz lamp heaters have a variety of heat radiation abilities, providing an adaptive heating environment for test samples. Next, in the temperature sensor module, the K_Thermocouple collects real temperature signals of quartz lamp heaters and transmits these temperature signals to the temperature transmitter isolation safety barrier. Finally, reference temperature trajectories and real temperature signals are fed back in the control center.

Step 3 (feedback): once the control system can track reference temperature trajectories successfully, the whole circle can be evaluated the TPS of hypersonic vehicles because it can reproduce the real aerodynamic heating environment on the ground.

3. Numerical Analyses

The first step of Figure 1 is aerodynamic heating data acquisition via the finite element calculation in ANSYS Workbench 2020 R2. In Figure 2, there is a three-dimensional drawing and a two-dimensional drawing of the hypersonic flying object. A hypersonic flying object is selected as the calculated object, and its specific parameters are: total length 295 mm, flying object warhead length 87 mm, flying object warhead radius 6 mm, flying object body diameter 40 mm, angle of the flying object warhead 6.31°, angle of the flying object body 1.38°, as shown in Figure 2b.

Based on the two-dimensional drawing of the hypersonic flying object, the range of the outflow field of the hypersonic flying object is 590 mm by 200 mm, in which 590 mm is 2 times length in horizontal direction and 200 mm is 5 times in longitudinal direction. In Figure 3, three different attack angles are set: 0°, 5°, and 10° and from a partially enlarged drawing of the hypersonic flying object’s warhead (Figure 3d), three different parts of the hypersonic flying object’s warhead are chosen as calculated objects: Wall 0, Wall 1, and Wall 2 because these three parts represent the whole of the hypersonic flying object’s warhead and are the maximum temperature values.

In

Table 1. The flight trajectory of the hypersonic flying object and corresponding atmospheric environment.

	H	T	a	P	R	M
a	31,272	227.922	302.6483	968	0.0148	0.61136
b	29,934	226.584	301.7587	1184	0.0182	0.768021
c	28,596	225.246	300.8664	1449	0.0224	0.924682
d	27,258	223.908	299.9715	1776	0.0276	1.081343
e	25,920	222.57	299.0738	2180	0.0341	1.238004
f	25,380	222.03	298.7108	2368	0.0372	1.394665
g	24,840	221.49	298.3473	2574	0.0405	1.551326
h	24,300	220.95	297.9834	2798	0.0441	1.707987
i	23,760	220.41	297.6191	3041	0.0481	1.864648
j	23,220	219.87	297.2543	3307	0.0524	2.021309
k	22,680	219.33	296.889	3597	0.0571	2.17797
l	22,140	218.79	296.5233	3913	0.0623	2.277316
m	21,600	218.25	296.1572	4258	0.068	2.376662
n	21,060	217.71	295.7906	4634	0.0742	2.476008
o	20,520	217.17	295.4235	5044	0.0809	2.575354
p	19,980	216.65	295.0696	5492	0.0883	2.6747
q	19,440	216.65	295.0696	5980	0.0962	2.774046
r	18,900	216.65	295.0696	6512	0.1047	2.873392
s	18,360	216.65	295.0696	7091	0.114	2.972738
t	17,820	216.65	295.0696	7721	0.1242	3.072084
u	17,280	216.65	295.0696	8407	0.1352	3.17143
v	16,861	216.65	295.0696	8981	0.1444	3.270776
w	16,450	216.65	295.0696	9583	0.1541	3.370122
x	16,040	216.65	295.0696	10,223	0.1644	3.469468
y	15,629	216.65	295.0696	10,907	0.1754	3.568814
z	15,219	216.65	295.0696	11,636	0.1871	3.66816
A	14,809	216.65	295.0696	12,413	0.1996	3.821
B	14,398	216.65	295.0696	13,244	0.213	4.212653
C	13,988	216.65	295.0696	14,129	0.2272	4.614305
D	13,577	216.65	295.0696	15,075	0.2424	5.015958

is a flight trajectory of the hypersonic flying object from flight altitude (H m) 31,272 m, flight Mach number 0.61136 to flight altitude 13,577 m, and flight Mach number 5.015958. The flight trajectory of the hypersonic flying object is a continuous flight course and that is divided into 30 Groups a–D of Table 1 with the same interval, which are different flight positions and represent the whole flight trajectory for the finite element calculation. The flight trajectories of the hypersonic flying object from Groups a–i are set with attack angle 0°; those from Groups j–x are set with attack angle 5°; those from Groups y–D are set with attack angle 10°. Based on the ATMOSCOESA (H) formula [36] and flight altitude data, the corresponding atmospheric environment values are obtained in MATLAB, including ambient temperature (T K), ambient sound velocity (a m/s), ambient pressure (P Pa), ambient density (R kg/m³). Some details about the ATMOSCOESA (H) formula are explained: SCOESA represents the U.S. Committee on Extension to the Standard Atmosphere, and the ATMOSCOESA (H) formula is a built-in MATLAB formula used for checking 1976 U.S. Standard Atmosphere. For example, if the flight altitude value of a hypersonic flying object is input into MATLAB, the corresponding atmospheric environment values (T, a, P, R) are obtained via the ATMOSCOESA (H) formula.

The two-dimensional outflow fields with different attack angles of the hypersonic flying object are imported into ANSYS Workbench 2020 R2, and its flight trajectory values with the corresponding atmospheric environment are also input into ANSYS Workbench 2020 R2. Then, all groups (a–D) adopt the triangles method to generate Mesh with element size 10.0 mm. So, the element quality of Mesh Metric is in

Table 2. The element quality of mesh metric.

	Min	Max	Average	Standard Deviation
0°	0.49732	0.99999	0.96312	0.036722
5°	0.39539	0.99997	0.96445	0.037045
10°	0.47993	0.99999	0.96397	0.037951

.

Moreover, the top and the bottom of the two-dimensional outflow fields are set as pressure-far-field; the left and the right of those are set as inlet and outlet, respectively. All boundary conditions correspond to Table 1. Some details of solution methods and solution controls are set in

Table 3. The specific parameters of solution methods.

	Pressure-Velocity Coupling				Spatial Discretization
	Scheme	Gradient	Pressure	Density	Momentum	Turbulent Kinetic Energy	Specific Dissipation Rate	Energy
0°	Coupled	Least Squares Cell Based	Second Order	Second Order Upwind	Second Order Upwind	First Order Upwind	First Order Upwind	Second Order Upwind
5°	Coupled	Green-Gauss Cell Based	Second Order	Second Order Upwind	Second Order Upwind	First Order Upwind	First Order Upwind	Second Order Upwind
10°	Coupled	Least Squares Cell Based	Second Order	Second Order Upwind	Second Order Upwind	First Order Upwind	First Order Upwind	Second Order Upwind

and

Table 4. The specific parameters of solution controls.

	Pseudo-Transient Explicit Relaxation Factors
	Pressure	Momentum	Density	Body Forces	Turbulent inetic Energy	Specific issipation Rate	Turbulent Viscosity	Energy
0°	0.5	0.5	1	1	0.75	0.75	1	0.75
5°	0.5	0.5	0.8	1	0.75	0.75	1	0.75
10°	0.5	0.5	1	0.9	0.75	0.75	1	0.75

, respectively. In addition, the solution initialization uses hybrid initialization methods to initialize the outflow fields and the number of iterations is 20. In calculation activities, the number of iterations of run calculation is 3000. In the end, all contours of groups a–D are presented in Figure 4, Figure 5 and Figure 6 corresponding to residuals in Figure 7.

According to the contours of all groups (a–D) in Figure 4, Figure 5 and Figure 6, the finite element calculation results of three chosen calculated objects (Wall 0, Wall 1, Wall 2) from the hypersonic flying object’s warhead are plotted to scatter graphs in Figure 8, Figure 9 and Figure 10. Figure 8a–c is connected to the relationship of Wall 0 between position and temperature; Figure 9a–c is connected to that of Wall 1; Figure 10a–c is connected to that of Wall 2. As a result, temperature values are positively associated with the flight Mach number.

Based on Table 1, the flight trajectory of the hypersonic flying object is divided into 30 groups: a–D, and each group has the same independent time interval, which is regarded as one second. Then, three positions of each Wall are selected as the reference points, which can represent the whole Wall 0 or Wall 1 or Wall 2. For example, Wall 0 chooses three positions: 0 mm, 3 mm, and 6 mm in Figure 8c, and the others are shown in Figure 9c and Figure 10c. Finally, the reference temperature trajectories are plotted to scatter graphs and the fitting curves of them are as follows:

Wall 0_0 mm (aim 1):

$$T^{*}\left(t\right)=4.632\ast {10}^{-4}{t}^5-0.02581t^4+0.4831t^3-3.131t^2+26.46t+214.9$$

(1)

Wall 0_3 mm (aim 2):

$$T^{*}\left(t\right)=7.184\ast {10}^{-4}{t}^5-0.04601t^4+1.008t^3-8.385t^2+40.19t+174.1$$

(2)

Wall 0_6 mm (aim 3):

$$T^{*}\left(t\right)=5.026\ast {10}^{-6}{t}^7-5.037\ast {10}^{-4}{t}^6+1.981\ast {10}^{-2}{t}^5-0.3811t^4+3.615t^3-14.33t^2+23.84t+197.4$$

(3)

Wall 1_5 mm (aim 4):

$$T^{*}\left(t\right)=5.72\ast {10}^{-4}{t}^5-0.03267t^4+0.5924t^3-3.008t^2+10.95t+196.8$$

(4)

Wall 1_45 mm (aim 5):

$$T^{*}\left(t\right)=3.364\ast {10}^{-6}{t}^7-3.708\ast {10}^{-4}{t}^6+1.632\ast {10}^{-2}{t}^5-0.3587t^4+4.049t^3-21.44t^2+51.3t+193.7$$

(5)

Wall 1_85 mm (aim 6):

$$T^{*}\left(t\right)=3.685\ast {10}^{-5}{t}^5+1.715\ast {10}^{-3}{t}^4-0.1414t^3+2.618t^2-5.852t+226.8$$

(6)

Wall 2_5 mm (aim 7):

$$T^{*}\left(t\right)=1.407\ast {10}^{-4}{t}^5-4.506\ast {10}^{-3}{t}^4-2.05\ast {10}^{-2}{t}^3+1.799t^2-0.6687t+205.5$$

(7)

Wall 2_45 mm (aim 8):

$$T^{*}\left(t\right)=3.479\ast {10}^{-4}{t}^5-0.0224t^4+0.5084t^3-4.591t^2+22.09t+206.1$$

(8)

Wall 2_85 mm (aim 9):

$$T^{*}\left(t\right)=4.126\ast {10}^{-4}{t}^5-0.02885t^4+0.7306t^3-7.885t^2+40.4t+177.8$$

(9)

where $T^{*}\left(t\right)$ is the reference temperature trajectories of the TSTQLs system and the goodness of fit results are illustrated in

Table 5. The goodness of fit results.

	SSE	R-Square	Adjusted R-Square	RMSE
Wall 0_0 mm	12,600	0.9962	0.9954	22.91
Wall 0_3 mm	9883	0.9949	0.9939	20.29
Wall 0_6 mm	11,360	0.9911	0.9882	22.72
Wall 1_5 mm	19,000	0.9904	0.9884	28.13
Wall 1_45 mm	3780	0.9963	0.9951	13.11
Wall 1_85 mm	4263	0.9948	0.9937	13.33
Wall 2_5 mm	15,250	0.9878	0.9853	25.21
Wall 2_45 mm	2293	0.9958	0.9949	9.773
Wall 2_85 mm	7589	0.9795	0.9753	17.78

.

4. Control System

The second step of Figure 1 is the control system of TSTQLs, including the system dynamic model [9,35] and some control methods, in which the system dynamic model is derived from the energy conservation law, and IPD, NESO, MFC, NGSMC are introduced.

4.1. The System Dynamic Model of TSTQLs

In the AC voltage regulating circuit, the input voltage of quartz lamp heaters (QLH) is expressed as:

$$U\left(t\right)=\sqrt{\frac{1}{\pi }{{\displaystyle \int }}_{\alpha \left(t\right)}^{\pi }{\left[\sqrt{2}{U}_{\mathrm{I}}\sin \left(\omega t\right)\right]}^2\mathrm{d}\left(\omega t\right)}=U_{\mathrm{I}}\sqrt{\frac{\sin \left[2\alpha \left(t\right)\right]}{2\pi }+\frac{\pi -\alpha \left(t\right)}{\pi }}$$

(10)

where $U\left(t\right)\in {ℝ}^{+}$ is the input voltage of QLH (V); $\alpha \left(t\right)\in {ℝ}^{+}$ is the SCR conduction angle (rad); $U_{\mathrm{I}}\in {ℝ}^{+}$ is the supply voltage (V); $\omega t\in \left[\begin{array}{cc}0,& \pi \end{array}\right]$ is the phase angle (rad).

Based on the electric power equation, the input electric power of QLH is calculated:

$$P\left(t\right)=\frac{U^2\left(t\right)}{R}=\frac{U_{\mathrm{I}}^2}{R}\left\{\frac{\sin \left[2\alpha \left(t\right)\right]}{2\pi }+\frac{\pi -\alpha \left(t\right)}{\pi }\right\}$$

(11)

where $P\left(t\right)\in {ℝ}^{+}$ is the input electric power of QLH (W); $R\in {ℝ}^{+}$ is the total resistance of QLH (Ω).

According to the theories of thermodynamics and heat transfer, the output electrothermal energy of QLH is given:

$$Q\left(t\right)=cm\left[T\left(t\right)-T\left(t-\Delta t\right)\right]+A\left\{\beta \left[T\left(t\right)-T\left(t-\Delta t\right)\right]+\lambda \left[T\left(t\right)-T\left(t-\Delta t\right)\right]+\epsilon \sigma F{T}^4\left(t\right)\Delta t\right\}$$

(12)

where the left of (12) is the output electrothermal energy $Q\left(t\right)\in {ℝ}^{+}$ of QLH (J); the right of (12) is the internal energy $cm\left[T\left(t\right)-T\left(t-\Delta t\right)\right]$, heat convection $A\beta \left[T\left(t\right)-T\left(t-\Delta t\right)\right]$, heat conduction $A\lambda \left[T\left(t\right)-T\left(t-\Delta t\right)\right]$, and heat radiation $A\epsilon \sigma F{T}^4\left(t\right)\Delta t$, respectively. $c\in {ℝ}^{+}$ is the specific heat capacity of quartz lamp filament ($\mathrm{J}/\mathrm{kg}·\mathrm{K}$); $m\in {ℝ}^{+}$ is the mass of quartz lamp filament ($\mathrm{kg}$); $T\left(t\right)\in {ℝ}^{+}$ is the current temperature of QLH (K); $T\left(t-\Delta t\right)\in {ℝ}^{+}$ is the QLH’s previous temperature (K) of the time interval $\Delta t\in {ℝ}^{+}$ (s); $A\in {ℝ}^{+}$ is the surface area of quartz lamp tube (${\mathrm{m}}^2$); $\beta \in {ℝ}^{+}$, $\lambda \in {ℝ}^{+}$, and $\epsilon \in {ℝ}^{+}$ are the coefficients of heat convection ($\mathrm{W}/{\mathrm{m}}^2·\mathrm{K}$), heat conduction ($\mathrm{W}/\mathrm{m}·\mathrm{K}$), and heat radiation blackness, respectively; $\sigma \in {ℝ}^{+}$ is the Stephen Boltzmann’s constant ($\mathrm{W}/{\mathrm{m}}^2·{\mathrm{K}}^4$); and $F\in {ℝ}^{+}$ is the angle coefficient.

On the basis of the energy conservation law, combining (11) with (12), the system dynamic model of TSTQLs is obtained:

$$\frac{U_{\mathrm{I}}^2}{R}\left\{\frac{\sin \left[2\alpha \left(t\right)\right]}{2\pi }+\frac{\pi -\alpha \left(t\right)}{\pi }\right\}\Delta t=cm\left[T\left(t\right)-T\left(t-\Delta t\right)\right]+A\left\{\beta \left[T\left(t\right)-T\left(t-\Delta t\right)\right]+\lambda \left[T\left(t\right)-T\left(t-\Delta t\right)\right]+\epsilon \sigma F{T}^4\left(t\right)\Delta t\right\}$$

(13)

where the system dynamic model of TSTQLs is established between input variable $\alpha \left(t\right)$ and output variable $T\left(t\right)$.

All parameters are explained in detail in

Table 6. All parameters of QLH.

Symbol	Unit	Description
$U\left(t\right)$	V	the input voltage of QLH
$\alpha \left(t\right)$	rad	the SCR conduction angle
$U_{\mathrm{I}}$	V	the supply voltage of QLH
$\omega t$	rad	the phase angle
$P\left(t\right)$	W	the input electric power of QLH
$R$	Ω	the total resistance of QLH
$c$	$\mathrm{J}/\mathrm{kg}·\mathrm{K}$	the specific heat capacity of quartz lamp filament
$m$	$\mathrm{kg}$	the mass of quartz lamp filament
$T\left(t\right)$	K	the current temperature of QLH
$T\left(t-\Delta t\right)$	K	the previous QLH’s temperature of the time interval
$\Delta t$	s	the time interval
$A$	${\mathrm{m}}^2$	the surface area of quartz lamp tube
$\beta$	$\mathrm{W}/{\mathrm{m}}^2·\mathrm{K}$	the heat convection coefficient of QLH
$\lambda$	$\mathrm{W}/\mathrm{m}·\mathrm{K}$	the heat conduction coefficient of QLH
$\epsilon$		the heat radiation blackness coefficient of QLH
$\sigma$	$\mathrm{W}/{\mathrm{m}}^2·{\mathrm{K}}^4$	the Stephen Boltzmann’s constant
$F$		the angle coefficient

.

4.2. Control Methods

Both sides of (13) are divided by $\Delta t$, and (13) can be further calculated:

$$\frac{U_{\mathrm{I}}^2}{R}\left\{\frac{\sin \left[2\alpha \left(t\right)\right]}{2\pi }+\frac{\pi -\alpha \left(t\right)}{\pi }\right\}=cm\frac{T\left(t\right)-T\left(t-\Delta t\right)}{\Delta t}+A\left[\beta \frac{T\left(t\right)-T\left(t-\Delta t\right)}{\Delta t}+\lambda \frac{T\left(t\right)-T\left(t-\Delta t\right)}{\Delta t}+\epsilon \sigma F{T}^4\left(t\right)\right]$$

(14)

Considering the system dynamic model of TSTQLs with uncertainties, (14) can be written:

$$\begin{gathered} \frac{U_{\mathrm{I}}^2}{R}\left\{\frac{\sin \left[2\alpha \left(t\right)\right]}{2\pi }+\frac{\pi -\alpha \left(t\right)}{\pi }\right\} \\ =cm\frac{T\left(t\right)-T\left(t-\Delta t\right)}{\Delta t}+A\left[\beta \frac{T\left(t\right)-T\left(t-\Delta t\right)}{\Delta t}+\lambda \frac{T\left(t\right)-T\left(t-\Delta t\right)}{\Delta t}+\epsilon \sigma F{T}^4\left(t\right)\right]+\Delta G\left(t\right) \end{gathered}$$

(15)

$$\Delta G\left(t\right)=\mathsf{\Delta }cm\frac{T\left(t\right)-T\left(t-\Delta t\right)}{\Delta t}+A\left[\Delta \beta \frac{T\left(t\right)-T\left(t-\Delta t\right)}{\Delta t}+\Delta \lambda \frac{T\left(t\right)-T\left(t-\Delta t\right)}{\Delta t}+\epsilon \sigma \Delta F{T}^4\left(t\right)\right]$$

(16)

where $\mathsf{\Delta }c$, $\Delta \beta$, $\Delta \lambda ,$ and $\Delta F$ are denoted as internal parametric uncertainties, and $\Delta G\left(t\right)$ is the lumped parametric uncertainties.

According to the definition of a derivative, the derivative of $T\left(t\right)$ is:

$$\frac{dT\left(t\right)}{dt}=\underset{\Delta t\to 0}{\lim }\frac{T\left(t\right)-T\left(t-\Delta t\right)}{\Delta t}$$

(17)

Substituting (17) to (15), the system dynamic model of TSTQLs can be written:

$$\frac{U_{\mathrm{I}}^2}{R}\left\{\frac{\sin \left[2\alpha \left(t\right)\right]}{2\pi }+\frac{\pi -\alpha \left(t\right)}{\pi }\right\}=cm\frac{dT\left(t\right)}{dt}+A\left[\beta \frac{dT\left(t\right)}{dt}+\lambda \frac{dT\left(t\right)}{dt}+\epsilon \sigma F{T}^4\left(t\right)\right]+\Delta G\left(t\right)$$

(18)

Then, (18) can be further calculated:

$$\frac{dT\left(t\right)}{dt}=G\left(t\right)-\frac{U_I^2}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}\alpha \left(t\right)$$

(19)

$$G\left(t\right)=\frac{1}{cm+A\left(\beta +\lambda \right)}\left\{\frac{U_{\mathrm{I}}^2\sin \left[2\alpha \left(t\right)\right]}{2R\pi }+\frac{U_{\mathrm{I}}^2}{R}-A\epsilon \sigma F{T}^4\left(t\right)-\Delta G\left(t\right)\right\}$$

(20)

where $G\left(t\right)$ is the lumped disturbances of TSTQLs and it contains periodic oscillations of the trigonometric function and strong nonlinearity of the high-order term, viewed as input variable disturbances and output variable disturbances.

4.2.1. NESO

A NESO [29,37] is designed for the lumped disturbances observation in the absence of the accurate system dynamic model, which is defined:

$$\left\{\begin{array}{c}{e}_1\left(t\right)=z_1\left(t\right)-T\left(t\right)\\ e_2\left(t\right)=z_2\left(t\right)-G\left(t\right)\\ \frac{d{z}_1\left(t\right)}{dt}=z_2\left(t\right)-{\beta }_1{e}_1\left(t\right)+b\alpha \left(t\right)\\ \frac{d{z}_2\left(t\right)}{dt}=-{\beta }_2{\left|e_1\left(t\right)\right|}^{1/2}\mathrm{sign}\left[e_1\left(t\right)\right]\\ z_2\left(t\right)=\widehat{G}\left(t\right)\end{array}\right.$$

(21)

where $z_1\left(t\right)$ is the observation of the current temperature $T\left(t\right)$; $z_2\left(t\right)$ is the observation of the lumped uncertainties $G\left(t\right)$ and $\widehat{G}\left(t\right)$ is an another expression form $z_2\left(t\right)$; $e_1\left(t\right)$ is the observation error of the observer state $z_1\left(t\right)$; $e_2\left(t\right)$ is the observation error of the observer state $z_2\left(t\right)$. ${\beta }_1$, ${\beta }_2$, $b$ are gain parameters and they satisfy:

$$\begin{gathered} {\beta }_1>0 \\ {\beta }_2>0 \\ b=-\frac{U_I^2}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]} \\ \mathrm{sign}\left[e_1\left(t\right)\right]=\left\{\begin{array}{cc}1& e_1\left(t\right)>0\\ 0& e_1\left(t\right)=0\\ -1& e_1\left(t\right)<0\end{array}\right. \end{gathered}$$

(22)

Define:

$$\tilde{G}\left(t\right)=G\left(t\right)-\widehat{G}\left(t\right)$$

(23)

where $\tilde{G}\left(t\right)$ is the observation error of $G\left(t\right)$.

Proof. 

The stability of NESO is proven [38].

Substituting (22) to (19), the system dynamic model of TSTQLs is further given:

$$\frac{dT\left(t\right)}{dt}=G\left(t\right)+b\alpha \left(t\right)$$

(24)

Substituting (24) to (21), then (21) is expressed:

$$\frac{d{e}_1\left(t\right)}{dt}=e_2\left(t\right)-{\beta }_1{e}_1\left(t\right)$$

(25)

where $\mathrm{d}{e}_1\left(t\right)/\mathrm{d}t$ is the derivative of $e_1\left(t\right)$ and satisfies:

$$\frac{d{e}_1\left(t\right)}{dt}=\frac{d{z}_1\left(t\right)}{dt}-\frac{dT\left(t\right)}{dt}$$

(26)

Then, the derivative of $e_2\left(t\right)$ is:

$$\frac{d{e}_2\left(t\right)}{dt}=-{\beta }_2{\left|e_1\left(t\right)\right|}^{1/2}\mathrm{sign}\left[e_1\left(t\right)\right]-\frac{dG\left(t\right)}{dt}$$

(27)

Define a Lyapunov function:

$$V={{\displaystyle \int }}_0^{e_1\left(t\right)}2{\beta }_2{\left|x\right|}^{1/2}\mathrm{sign}\left(x\right)\mathrm{d}x+e_{11}^2\left(t\right)$$

(28)

where $e_{11}\left(t\right)$ is:

$$e_{11}\left(t\right)=e_2\left(t\right)-{\beta }_1{e}_1\left(t\right)$$

(29)

In addition, the derivative of (29) is:

$$\frac{d{e}_{11}\left(t\right)}{dt}=-{\beta }_2{\left|e_1\left(t\right)\right|}^{1/2}\mathrm{sign}\left[e_1\left(t\right)\right]-\frac{dG\left(t\right)}{dt}-{\beta }_1\left[e_2\left(t\right)-{\beta }_1{e}_1\left(t\right)\right]$$

(30)

Based on the mean value theorem for integrals, the (28) is calculated:

$$V=2{\beta }_2{\left|ϵ\right|}^{1/2}\mathrm{sign}\left(ϵ\right)e_1\left(t\right)+e_{11}^2\left(t\right)$$

(31)

where $\exists ϵ\in \left[\begin{array}{cc}0,& e_1\left(t\right)\end{array}\right]$, $e_1\left(t\right)>0$ is for assumption and $V>0$.

Substituting (29) and (30), then the derivative of (28) is:

$$\begin{gathered} \frac{dV}{dt}=2{\beta }_2{\left|e_1\left(t\right)\right|}^{1/2}\mathrm{sign}\left[e_1\left(t\right)\right]e_{11}\left(t\right)+2e_{11}\left(t\right)\left\{-{\beta }_2{\left|e_1\left(t\right)\right|}^{1/2}\mathrm{sign}\left[e_1\left(t\right)\right]-\frac{dG\left(t\right)}{dt}-{\beta }_1{e}_{11}\left(t\right)\right\} \\ =2e_{11}\left(t\right)\left[-\frac{dG\left(t\right)}{dt}-{\beta }_1{e}_{11}\left(t\right)\right] \end{gathered}$$

(32)

where $\mathrm{d}V/\mathrm{d}t<0$ satisfies:

If

$$\begin{gathered} \frac{dG\left(t\right)}{dt}>0 \\ \left|e_{11}\left(t\right)\right|>\frac{dG\left(t\right)}{dt}\frac{1}{{\beta }_1} \end{gathered}$$

(33)

Or if

$$\begin{gathered} \frac{\mathrm{d}G\left(t\right)}{\mathrm{d}t}<0 \\ \left|e_{11}\left(t\right)\right|>-\frac{dG\left(t\right)}{dt}\frac{1}{{\beta }_1} \end{gathered}$$

(34)

So,

$$\left|e_{11}\left(t\right)\right|>\left|\frac{dG\left(t\right)}{dt}\right|\frac{1}{{\beta }_1}$$

(35)

Hence, when $\exists ϵ\in \left[\begin{array}{cc}0,& e_1\left(t\right)\end{array}\right]$, $e_1\left(t\right)>0$ is for assumption and $\left|e_{11}\left(t\right)\right|>\left|\frac{dG\left(t\right)}{dt}\right|\frac{1}{{\beta }_1}$, the stability of NESO has been proven. □

4.2.2. Controller Design

Define the system error:

$$e\left(t\right)=T^{*}\left(t\right)-T\left(t\right)$$

(36)

The derivative of (36) is:

$$\frac{de\left(t\right)}{dt}=\frac{d{T}^{*}\left(t\right)}{dt}-\frac{dT\left(t\right)}{dt}$$

(37)

Definition 1.

The controller 1 with NESO is designed by an integral sliding mode control (ISMC) [39] called the ISMCNESO controller, which is illustrated in Figure 11.

Define the ISMC:

$$s\left(t\right)=e\left(t\right)+{{\displaystyle \int }}_0^t\left(a+bk\right)e\left(\tau \right)d\tau$$

(38)

where (38) is an integral sliding mode surface.$k$ is a gain parameter, and $a$, $b$ are from the system dynamic model parameters, which satisfy:

$$\begin{gathered} a=-\frac{A\epsilon \sigma F}{cm+A\left(\beta +\lambda \right)} \\ b=-\frac{U_I^2}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]} \end{gathered}$$

(39)

The derivative of (38) is:

$$\frac{ds\left(t\right)}{dt}=\frac{de\left(t\right)}{dt}+\left(a+bk\right)e\left(t\right)$$

(40)

Define the ISMCNESO controller:

$${\alpha }_{\mathrm{ISMCNESO}}\left(t\right)={\alpha }_{\mathrm{eq}\_\mathrm{ISMCNESO}}\left(t\right)+{\alpha }_{\mathrm{re}\_\mathrm{ISMCNESO}}\left(t\right)$$

(41)

where ${\alpha }_{\mathrm{ISMCNESO}}\left(t\right)$ is the input of the ISMCNESO controller, consisting of the equivalent control ${\alpha }_{\mathrm{eq}\_\mathrm{ISMCNESO}}\left(t\right)$ and the reaching law ${\alpha }_{\mathrm{re}\_\mathrm{ISMCNESO}}\left(t\right)$.

Substituting (19) and (37) to (40) and let $\mathrm{d}s\left(t\right)/\mathrm{d}t=0$, the equivalent control ${\alpha }_{\mathrm{eq}}\left(t\right)$ is calculated:

$${\alpha }_{\mathrm{eq}\_\mathrm{ISMCNESO}}\left(t\right)=\left[\widehat{G}\left(t\right)-\frac{d{T}^{*}\left(t\right)}{dt}-\left(a+bk\right)e\left(t\right)\right]\ast \frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}$$

(42)

Define the reaching law ${\alpha }_{\mathrm{re}\_\mathrm{ISMCNESO}}\left(t\right)$:

$$\begin{gathered} {\alpha }_{\mathrm{re}\_\mathrm{ISMCNESO}}\left(t\right)=-\epsilon \mathrm{sign}\left[s\left(t\right)\right] \\ \mathrm{sign}\left[s\left(t\right)\right]=\left\{\begin{array}{cc}1& s\left(t\right)>0\\ 0& s\left(t\right)=0\\ -1& s\left(t\right)<0\end{array}\right. \end{gathered}$$

(43)

where the ${\alpha }_{\mathrm{re}\_\mathrm{ISMCNESO}}\left(t\right)$ is traditional constant reaching law and $\epsilon >0$ is a gain parameter.

Then, (41) is further calculated:

$${\alpha }_{\mathrm{ISMCNESO}}\left(t\right)=\left[\widehat{G}\left(t\right)-\frac{d{T}^{\ast }\left(t\right)}{dt}-\left(a+bk\right)e\left(t\right)\right]\ast \frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}-\epsilon \mathrm{sign}\left[s\left(t\right)\right]$$

(44)

Proof. 

The stability of the ISMCNESO controller is proven.

Define a Lyapunov function:

$$V=\frac{1}{2}{s}^2\left(t\right)$$

(45)

The derivative of (45) is:

$$\frac{dV}{dt}=s\left(t\right)\frac{ds\left(t\right)}{dt}$$

(46)

Substituting (19), (37), and (40) to (46), (46) is further calculated:

$$\begin{gathered} \frac{dV}{dt}=s\left(t\right)\left[\frac{de\left(t\right)}{dt}+\left(a+bk\right)e\left(t\right)\right] \\ \frac{dV}{dt}=s\left(t\right)\left[\frac{d{T}^{*}\left(t\right)}{dt}-G\left(t\right)+\frac{U_I^2}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}\alpha \left(t\right)+\left(a+bk\right)e\left(t\right)\right] \end{gathered}$$

(47)

Substituting (44) to (47), (47) is further calculated:

$$\begin{gathered} \frac{dV}{dt}=s\left(t\right)\left\{\frac{d{T}^{*}\left(t\right)}{dt}\right.-G\left(t\right) \\ \left.\left.+\frac{U_I^2}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}\left\{\left[\widehat{G}\left(t\right)-\frac{d{T}^{*}\left(t\right)}{dt}-\left(a+bk\right)e\left(t\right)\right]\right.\frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}-\epsilon \mathrm{sign}\left[s\left(t\right)\right]\right\}+\left(a+bk\right)e\left(t\right)\right\} \\ \frac{dV}{dt}=-\tilde{G}\left(t\right)s\left(t\right)-\frac{\epsilon U_I^2}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}\left|s\left(t\right)\right| \end{gathered}$$

(48)

Assumption 1.

In practical applications, the lumped disturbances are bounded.

Define a boundary condition:

$$\left|\tilde{G}\left(t\right)\right|<\frac{\epsilon U_I^2}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}$$

(49)

Hence, according to Assumption 1, $\mathrm{d}V/\mathrm{d}t<0$ and the ISMCNESO controller is stable. □

Definition 2.

The controller 2 with NESO is designed by theNGSMC, named the NGSMCNESO controller, which is depicted in Figure 12.

Define the NGSMC:

$$s\left(t\right)=fal\left[e\left(t\right),\gamma ,\eta \right]+{{\displaystyle \int }}_0^t\left(a+bk\right)fal\left[e\left(\tau \right),\gamma ,\eta \right]\mathrm{d}\tau -\frac{e\left(0\right)}{{\eta }^{1-\gamma }}$$

(50)

where the (50) is a nonlinear global sliding mode surface. $e\left(0\right)$ is the initial state of system error, and $\gamma$, $\eta$,$a$, $b$ are gain parameters, which satisfy:

$$\begin{gathered} fal\left[e\left(t\right),\gamma ,\eta \right]=\left\{\begin{array}{cc}{\left|e\left(t\right)\right|}^{\gamma }\mathrm{sign}\left[e\left(t\right)\right]& \left|e\left(t\right)\right|>\eta \\ \frac{e\left(t\right)}{{\eta }^{1-\gamma }}& \left|e\left(t\right)\right|<\eta \end{array}\right. \\ 0<\gamma <1 \\ \eta >0 \\ \mathrm{sign}\left[e\left(t\right)\right]=\left\{\begin{array}{cc}1& e\left(t\right)>0\\ 0& e\left(t\right)=0\\ -1& e\left(t\right)<0\end{array}\right. \\ a=-\frac{A\epsilon \sigma F}{cm+A\left(\beta +\lambda \right)} \\ b=-\frac{U_I^2}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]} \end{gathered}$$

(51)

The derivative of (50) is:

$$\frac{ds\left(t\right)}{dt}=\left\{\begin{array}{cc}\gamma {\left|e\left(t\right)\right|}^{\gamma -1}\frac{de\left(t\right)}{dt}+\left(a+bk\right){\left|e\left(t\right)\right|}^{\gamma }\mathrm{sign}\left[e\left(t\right)\right]& \left|e\left(t\right)\right|>\eta \\ \frac{1}{{\eta }^{1-\gamma }}\frac{de\left(t\right)}{dt}+\left(a+bk\right)\frac{e\left(t\right)}{{\eta }^{1-\gamma }}& \left|e\left(t\right)\right|<\eta \end{array}\right.$$

(52)

Define the NGSMCNESO controller:

$${\alpha }_{\mathrm{NGSMCNESO}}\left(t\right)={\alpha }_{\mathrm{eq}\_\mathrm{NGSMCNESO}}\left(t\right)+{\alpha }_{\mathrm{re}\_\mathrm{NGSMCNESO}}\left(t\right)$$

(53)

where ${\alpha }_{\mathrm{NGSMCNESO}}\left(t\right)$ is the input of the NGSMCNESO controller, consisting of the equivalent control ${\alpha }_{\mathrm{eq}\_\mathrm{NGSMCNESO}}\left(t\right)$ and the reaching law ${\alpha }_{\mathrm{re}\_\mathrm{NGSMCNESO}}\left(t\right)$.

Substituting (19) and (37) to (52) and let $ds\left(t\right)/dt=0$, the equivalent control ${\alpha }_{\mathrm{eq}\_\mathrm{NGSMCNESO}}\left(t\right)$ is calculated:

$$\begin{gathered} \left|e\left(t\right)\right|>\eta \\ {\alpha }_{\mathrm{eq}\_\mathrm{NGSMCNESO}}\left(t\right)=-\frac{1}{\gamma {\left|e\left(t\right)\right|}^{\gamma -1}}\frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}\left(a+bk\right){\left|e\left(t\right)\right|}^{\gamma }\mathrm{sign}\left[e\left(t\right)\right] \\ +\frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}\left[\widehat{G}\left(t\right)-\frac{d{T}^{*}\left(t\right)}{dt}\right] \\ \left|e\left(t\right)\right|<\eta \\ {\alpha }_{\mathrm{eq}\_\mathrm{NGSMCNESO}}\left(t\right)=-{\eta }^{1-\gamma }\frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}\left(a+bk\right)\frac{e\left(t\right)}{{\eta }^{1-\gamma }}+\frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}\left[\widehat{G}\left(t\right)-\frac{d{T}^{*}\left(t\right)}{dt}\right] \end{gathered}$$

(54)

Then, the reaching law ${\alpha }_{\mathrm{re}\_\mathrm{NGSMCNESO}}\left(t\right)$ of the NGSMCNESO controller is the same as that of the ISMCNESO controller. So, (53) is further calculated:

$${\alpha }_{\mathrm{NGSMCNESO}}\left(t\right)=\left\{\begin{array}{c}-\frac{1}{\gamma {\left|e\left(t\right)\right|}^{\gamma -1}}\ast \frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}\ast \left(a+bk\right){\left|e\left(t\right)\right|}^{\gamma }\mathrm{sign}\left[e\left(t\right)\right]\\ +\frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}\left[\widehat{G}\left(t\right)-\frac{d{T}^{*}\left(t\right)}{dt}\right]-\epsilon sign\left[s\left(t\right)\right] \left|e\left(t\right)\right|>\eta \\ -{\eta }^{1-\gamma }\frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}\left(a+bk\right)\frac{e\left(t\right)}{{\eta }^{1-\gamma }}\\ +\frac{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}{U_I^2}\left[\widehat{G}\left(t\right)-\frac{d{T}^{*}\left(t\right)}{dt}\right]-\epsilon sign\left[s\left(t\right)\right] \left|e\left(t\right)\right|<\eta \end{array}\right.$$

(55)

where, for the NGSMCNESO controller, on the basis of the traditional integral sliding mode surface (38), a nonlinear global sliding mode surface (50) with a nonlinear Function (51) has the integral component with the exponential stability for the convergence rate in finite time $t_s$.

The calculated process of $t_s$:

If $\left|e\left(t\right)\right|>\eta$, (52) can be rewritten:

$$\gamma {\left|e\left(t\right)\right|}^{\gamma -1}\frac{de\left(t\right)}{dt}+\left(a+bk\right){\left|e\left(t\right)\right|}^{\gamma }\mathrm{sign}\left[e\left(t\right)\right]=0$$

(56)

When $e\left(t\right)>0$, the (56) is given:

$$\begin{gathered} \gamma e^{\gamma -1}\left(t\right)\frac{de\left(t\right)}{dt}+\left(a+bk\right)e^{\gamma }\left(t\right)=0 \\ \frac{de\left(t\right)}{e\left(t\right)}=-\frac{\left(a+bk\right)}{\gamma }\mathrm{d}t \end{gathered}$$

(57)

Then, $t_s$ is obtained:

$$t_s=-\frac{\gamma \ln e\left(t\right)}{\left(a+bk\right)}+\zeta$$

(58)

where $\zeta \in ℝ$ is an arbitrary constant.

When $e\left(t\right)<0$, (56) is given:

$$\begin{gathered} \gamma {\left[-e\left(t\right)\right]}^{\gamma -1}\frac{de\left(t\right)}{dt}-\left(a+bk\right){\left[-e\left(t\right)\right]}^{\gamma }=0 \\ \frac{de\left(t\right)}{e\left(t\right)}=-\frac{\left(a+bk\right)}{\gamma }\mathrm{d}t \end{gathered}$$

(59)

Then, $t_s$ is obtained:

$$t_s=-\frac{\gamma \ln \left[-e\left(t\right)\right]}{\left(a+bk\right)}+\zeta$$

(60)

If $\left|e\left(t\right)\right|<\eta$, (52) can be rewritten:

$$\begin{gathered} \frac{1}{{\eta }^{1-\gamma }}\frac{de\left(t\right)}{dt}+\left(a+bk\right)\frac{e\left(t\right)}{{\eta }^{1-\gamma }}=0 \\ \frac{de\left(t\right)}{e\left(t\right)}=-\left(a+bk\right)\mathrm{d}t \end{gathered}$$

(61)

When $e\left(t\right)>0$, (61) is given:

$$\ln e\left(t\right)=-\left(a+bk\right)t$$

(62)

Then, $t_s$ is obtained:

$$t=-\frac{\ln e\left(t\right)}{\left(a+bk\right)}+\zeta$$

(63)

When $e\left(t\right)<0$, (61) is given:

$$\ln \left[-e\left(t\right)\right]=-\left(a+bk\right)t$$

(64)

Then, $t_s$ is obtained:

$$t_s=-\frac{\ln \left[-e\left(t\right)\right]}{\left(a+bk\right)}+\zeta$$

(65)

As a whole, $t_s$ is obtained:

$$t_s=\left\{\begin{array}{cc}-\frac{\gamma \ln \left|e\left(t\right)\right|}{\left(a+bk\right)}+\zeta & \left|e\left(t\right)\right|>\eta \\ -\frac{\ln \left|e\left(t\right)\right|}{\left(a+bk\right)}+\zeta & \left|e\left(t\right)\right|<\eta \end{array}\right.$$

(66)

Moreover, the NGSMC (50), a global sliding mode surface, lets the initial system state be trapped on the prescribed sliding mode surface (50), eliminating chattering phenomena from the high-frequency switching motions and avoiding the problem of sensitivity to internal parametric uncertainties and external disturbances; a nonlinear Function (51) can provide: smaller errors corresponding to larger gains and larger errors corresponding to smaller gains.

Proof. 

The stability of the NGSMCNESO controller is proven.

Define a Lyapunov function:

$$V=\frac{1}{2}{s}^2\left(t\right)$$

(67)

The derivative of (67) is:

$$\frac{dV}{dt}=s\left(t\right)\frac{ds\left(t\right)}{dt}$$

(68)

Substituting (19), (37), and (52) to (68), (68) is further calculated:

$$\frac{dV}{dt}=s\left(t\right)\left\{\begin{array}{cc}\gamma {\left|e\left(t\right)\right|}^{\gamma -1}\left\{\frac{d{T}^{*}\left(t\right)}{dt}-G\left(t\right)+\frac{U_I^2}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}\alpha \left(t\right)\right\}+\left(a+bk\right)fal\left[e\left(t\right),\gamma ,\eta \right]& \left|e\left(t\right)\right|>\eta \\ \frac{1}{{\eta }^{1-\gamma }}\left\{\frac{d{T}^{*}\left(t\right)}{dt}-G\left(t\right)+\frac{U_I^2}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}\alpha \left(t\right)\right\}+\left(a+bk\right)fal\left[e\left(t\right),\gamma ,\eta \right]& \left|e\left(t\right)\right|<\eta \end{array}\right.$$

(69)

Substituting (55) to (69), (69) is further calculated:

$$\frac{dV}{dt}=s\left(t\right)\left\{\begin{array}{cc}\gamma {\left|e\left(t\right)\right|}^{\gamma -1}\left\{-\tilde{G}\left(t\right)+\frac{-\epsilon U_I^2\mathrm{sign}\left[s\left(t\right)\right]}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}\right\}& \left|e\left(t\right)\right|>\eta \\ \frac{1}{{\eta }^{1-\gamma }}\left\{-\tilde{G}\left(t\right)+\frac{-\epsilon U_I^2\mathrm{sign}\left[s\left(t\right)\right]}{R\pi \left[cm+A\left(\beta +\lambda \right)\right]}\right\}& \left|e\left(t\right)\right|<\eta \end{array}\right.$$

(70)

Hence, according to Assumption 1, $dV/dt<0$ and the NGSMCNESO controller is stable. □

Definition 3.

The controller 3 with NESO is designed by an IPD and a NGSMC in an ultra-local model frame of MFC, called the composite controller, as is shown in Figure 13.

Based on the system dynamic model (19), define:

$$\frac{dT\left(t\right)}{dt}=G\left(t\right)+\chi \alpha \left(t\right)$$

(71)

where (71) is a one-order ultra-local model [19]. $G\left(t\right)$ are the lumped disturbances and $\chi$ is a gain parameter without any real physical meaning.

Define an IPD as a closed-loop controller:

$${\alpha }_{\mathrm{IPD}}\left(t\right)=\frac{1}{\chi }\left[-\widehat{G}\left(t\right)+\frac{d{T}^{*}\left(t\right)}{dt}+K_{\mathrm{P}}e\left(t\right)+K_{\mathrm{D}}\frac{de\left(t\right)}{dt}\right]$$

(72)

Substituting (71) and (72) to (37), the system error equation is given:

$$\begin{gathered} \frac{de\left(t\right)}{dt}=\frac{d{T}^{*}\left(t\right)}{dt}-G\left(t\right)-\chi \frac{1}{\chi }\left[-\widehat{G}\left(t\right)+\frac{d{T}^{*}\left(t\right)}{dt}+K_{P}e\left(t\right)+K_D\frac{de\left(t\right)}{dt}\right] \\ \frac{de\left(t\right)}{dt}+\tilde{G}\left(t\right)+\left[K_{P}e\left(t\right)+K_D\frac{de\left(t\right)}{dt}\right]=0 \end{gathered}$$

(73)

According to Laplace transform $F\left(s\right)$, (73) is further calculated:

$$\begin{gathered} sE\left(s\right)+K_{\mathrm{P}}E\left(s\right)+K_{\mathrm{D}}sE\left(s\right)+\tilde{G}\left(s\right)-\tilde{G}\left(0^{+}\right)=0 \\ \underset{t\to \infty }{\lim }e\left(t\right)=\underset{s\to 0}{\lim }sE\left(s\right)=\frac{s}{s+K_{\mathrm{P}}+K_{\mathrm{D}}s}\tilde{G}\left(0^{+}\right)-\tilde{G}\left(s\right) \end{gathered}$$

(74)

where (74) reveals a fact that the system error $e\left(t\right)$ relies on some gain parameters $K_{\mathrm{P}}$ and $K_{\mathrm{D}}$ and the lumped disturbances observer NESO, but the system error $e\left(t\right)$ has difficulty trending to zero in the presence of bounded $\tilde{G}\left(t\right)$ and measurement noise. Hence an auxiliary controller is added into the IPDNGSMC controller for the observation errors compensation.

Use an NGSMC as an auxiliary controller and define the IPDNGSMC controller:

$${\alpha }_{\mathrm{IPDNGSMC}}\left(t\right)=\frac{1}{\chi }\left[-\widehat{G}\left(t\right)+\frac{d{T}^{*}\left(t\right)}{dt}+K_{\mathrm{P}}e\left(t\right)+K_{\mathrm{D}}\frac{de\left(t\right)}{dt}\right]+{\alpha }_{\mathrm{NGSMC}}\left(t\right)$$

(75)

where ${\alpha }_{\mathrm{NGSMC}}\left(t\right)$ is an auxiliary controller.

Substituting (37) and (75) to (71), the system error is:

$$\tilde{G}\left(t\right)+\frac{de\left(t\right)}{dt}+K_{\mathrm{P}}e\left(t\right)+K_{\mathrm{D}}\frac{de\left(t\right)}{dt}+\chi {\alpha }_{\mathrm{NGSMC}}\left(t\right)=0$$

(76)

where ${\alpha }_{\mathrm{NGSMC}}\left(t\right)$ is made of the equivalent control ${\alpha }_{\mathrm{eq}\_\mathrm{NGSMC}}\left(t\right)$ and the reaching law ${\alpha }_{\mathrm{re}\_\mathrm{NGSMC}}\left(t\right)$.

Substituting (76) to (52) and letting $\mathrm{d}s\left(t\right)/\mathrm{d}t=0$, the equivalent control of an auxiliary controller is calculated:

$${\alpha }_{\mathrm{eq}\_\mathrm{NGSMC}}\left(t\right)=\left\{\begin{array}{cc}\frac{1}{\chi \gamma {\left|e\left(t\right)\right|}^{\gamma -1}}\left(a+bk\right){\left|e\left(t\right)\right|}^{\gamma }\mathrm{sign}\left[e\left(t\right)\right]-\frac{1}{\chi }\left[\tilde{G}\left(t\right)+K_{\mathrm{P}}e\left(t\right)+K_{\mathrm{D}}\frac{de\left(t\right)}{dt}\right]& \left|e\left(t\right)\right|>\eta \\ \frac{{\eta }^{1-\gamma }}{\chi }\left(a+bk\right)\frac{e\left(t\right)}{{\eta }^{1-\gamma }}-\frac{1}{\chi }\left[\tilde{G}\left(t\right)+K_{\mathrm{P}}e\left(t\right)+K_{\mathrm{D}}\frac{de\left(t\right)}{dt}\right]& \left|e\left(t\right)\right|<\eta \end{array}\right.$$

(77)

where ${\alpha }_{\mathrm{eq}\_\mathrm{NGSMC}}\left(t\right)$ is the equivalent control of an auxiliary controller.

Define the reaching law ${\alpha }_{\mathrm{re}}\left(t\right)$ [40]:

$${\alpha }_{\mathrm{re}\_\mathrm{NGSMC}}\left(t\right)=-k_1{\left|s\left(t\right)\right|}^p\mathrm{sign}\left[s\left(t\right)\right]-k_{2}f{\left[s\left(t\right)\right]}^q$$

(78)

where $k_1$, $k_2$, $p$, $q$ are gain parameters. $q$ is a positive odd integer and $r$ is a positive integer. All parameters satisfy:

$$\begin{gathered} k_1>0 \\ {k}_2>0 \\ 0\le p<1 \\ f\left(s\right)=\left\{\begin{array}{cc}s\left(t\right)& \left|s\left(t\right)\right|<r\\ \mathrm{sign}\left[s\left(t\right)\right]& \left|s\left(t\right)\right|\ge r\end{array}\right. \end{gathered}$$

(79)

So, the whole auxiliary controller is expressed:

$$\begin{gathered} \left|e\left(t\right)\right|>\eta \\ {\alpha }_{\mathrm{NGSMC}}\left(t\right)=\frac{1}{\chi \gamma {\left|e\left(t\right)\right|}^{\gamma -1}}\left(a+bk\right){\left|e\left(t\right)\right|}^{\gamma }\mathrm{sign}\left[e\left(t\right)\right]-\frac{1}{\chi }\left[\tilde{G}\left(t\right)+K_{\mathrm{P}}e\left(t\right)+K_{\mathrm{D}}\frac{de\left(t\right)}{dt}\right] \\ -k_1{\left|s\left(t\right)\right|}^p\mathrm{sign}\left[s\left(t\right)\right]-k_{2}f{\left[s\left(t\right)\right]}^q \\ \left|e\left(t\right)\right|<\eta \\ {\alpha }_{\mathrm{NGSMC}}\left(t\right)=\frac{{\eta }^{1-\gamma }}{\chi }\left(a+bk\right)\frac{e\left(t\right)}{{\eta }^{1-\gamma }}-\frac{1}{\chi }\left[\tilde{G}\left(t\right)+K_{\mathrm{P}}e\left(t\right)+K_{\mathrm{D}}\frac{de\left(t\right)}{dt}\right]-k_1{\left|s\left(t\right)\right|}^p\mathrm{sign}\left[s\left(t\right)\right]-k_{2}f{\left[s\left(t\right)\right]}^q \end{gathered}$$

(80)

where ${\alpha }_{\mathrm{NGSMC}}\left(t\right)$ is the whole of an auxiliary controller.

Then, the input of the composite controller is written:

$$\begin{gathered} \left|e\left(t\right)\right|>\eta \\ {\alpha }_{\mathrm{composite} \mathrm{controller}}\left(t\right)=\frac{1}{\chi }\left[-\widehat{G}\left(t\right)+\frac{d{T}^{*}\left(t\right)}{dt}\right]+\frac{1}{\chi \gamma {\left|e\left(t\right)\right|}^{\gamma -1}}\left(a+bk\right){\left|e\left(t\right)\right|}^{\gamma }\mathrm{sign}\left[e\left(t\right)\right] \\ -k_1{\left|s\left(t\right)\right|}^p\mathrm{sign}\left[s\left(t\right)\right]-k_{2}f{\left[s\left(t\right)\right]}^q \\ \left|e\left(t\right)\right|<\eta \\ {\alpha }_{\mathrm{composite} \mathrm{controller}}\left(t\right)=\frac{1}{\chi }\left[-\widehat{G}\left(t\right)+\frac{d{T}^{*}\left(t\right)}{dt}\right]+\frac{{\eta }^{1-\gamma }}{\chi }\left(a+bk\right)\frac{e\left(t\right)}{{\eta }^{1-\gamma }}-k_1{\left|s\left(t\right)\right|}^p\mathrm{sign}\left[s\left(t\right)\right]-k_{2}f{\left[s\left(t\right)\right]}^q \end{gathered}$$

(81)

where ${\alpha }_{\mathrm{composite} \mathrm{controller}}\left(t\right)$ is the composite controller including IPD, NGSMC and NESO. For the composite controller, IPD (72), as a closed-loop controller, NESO, as a lumped disturbances observation, NGSMC (50), as an auxiliary controller are integrated into a one-order ultra-local model frame of MFC (71). The composite controller maintains the advantages of the NGSMCNESO controller (55) and introduces a reaching law (78) to further suppress chattering phenomena from the high-frequency switching motions, which is the combination of an exponential reaching law and a power reaching law. From the reaching law (78), when the system error state is far away from the sliding mode surface (50), $1\le \left|s\left(t\right)\right|$, then ${\left|s\left(t\right)\right|}^p<\left|s\left(t\right)\right|$ means that the reaching rate (78) alleviates chattering phenomena; when the system error state is near the sliding mode surface (50), $\left|s\left(t\right)\right|<1$, then $\left|s\left(t\right)\right|<{\left|s\left(t\right)\right|}^p$ means that the reaching rate (78) is faster than the traditional exponential reaching law, and $r$ limits the amplitude from the high-frequency switching motions, so the reaching rate (78) has variable reaching speed.

Proof. 

The stability of the composite controller is proven.

Define a Lyapunov function:

$$V=\frac{1}{2}{s}^2\left(t\right)$$

(82)

The derivative of (82) is:

$$\frac{dV}{dt}=s\left(t\right)\frac{ds\left(t\right)}{dt}$$

(83)

Substituting (37), (52), and (71) to (83), (83) is further calculated:

$$\frac{dV}{dt}=s\left(t\right)\left\{\begin{array}{cc}\gamma {\left|e\left(t\right)\right|}^{\gamma -1}\left\{\frac{d{T}^{*}\left(t\right)}{dt}-G\left(t\right)-\chi \alpha \left(t\right)\right\}+\left(a+bk\right){\left|e\left(t\right)\right|}^{\gamma }\mathrm{sign}\left[e\left(t\right)\right]& \left|e\left(t\right)\right|>\eta \\ \frac{1}{{\eta }^{1-\gamma }}\left\{\frac{d{T}^{*}\left(t\right)}{dt}-G\left(t\right)-\chi \alpha \left(t\right)\right\}+\left(a+bk\right)\frac{e\left(t\right)}{{\eta }^{1-\gamma }}& \left|e\left(t\right)\right|<\eta \end{array}\right.$$

(84)

Substituting (81) to (84), (84) is further calculated:

$$\frac{dV}{dt}=s\left(t\right)\left\{\begin{array}{cc}\gamma {\left|e\left(t\right)\right|}^{\gamma -1}\left\{-\tilde{G}\left(t\right)+\chi \left\{k_1{\left|s\left(t\right)\right|}^p\mathrm{sign}\left[s\left(t\right)\right]+k_{2}f{\left[s\left(t\right)\right]}^q\right\}\right\}& \left|e\left(t\right)\right|>\eta \\ \frac{1}{{\eta }^{1-\gamma }}\left\{-\tilde{G}\left(t\right)+\chi \left\{k_1{\left|s\left(t\right)\right|}^p\mathrm{sign}\left[s\left(t\right)\right]+k_{2}f{\left[s\left(t\right)\right]}^q\right\}\right\}& \left|e\left(t\right)\right|<\eta \end{array}\right.$$

(85)

If $\left|s\left(t\right)\right|<r$, (85) is further calculated:

$$\frac{dV}{dt}=\left\{\begin{array}{cc}\gamma {\left|e\left(t\right)\right|}^{\gamma -1}\left\{-\tilde{G}\left(t\right)s\left(t\right)+\chi \left[k_1{\left|s\left(t\right)\right|}^{p+1}+k_2{s}^{q+1}\left(t\right)\right]\right\}& \left|e\left(t\right)\right|>\eta \\ \frac{1}{{\eta }^{1-\gamma }}\left\{-\tilde{G}\left(t\right)s\left(t\right)+\chi \left[k_1{\left|s\left(t\right)\right|}^{p+1}+k_2{s}^{q+1}\left(t\right)\right]\right\}& \left|e\left(t\right)\right|<\eta \end{array}\right.$$

(86)

Then, $k_1{\left|s\left(t\right)\right|}^{p+1}+k_2{s}^{q+1}\left(t\right)>0.$

If $\left|s\left(t\right)\right|\ge r$,(85) is further calculated:

$$\frac{dV}{dt}=\left\{\begin{array}{cc}\gamma {\left|e\left(t\right)\right|}^{\gamma -1}\left\{-\tilde{G}\left(t\right)s\left(t\right)+\chi \left[k_1{\left|s\left(t\right)\right|}^{p+1}+k_2\left|s\left(t\right)\right|\right]\right\}& \left|e\left(t\right)\right|>\eta \\ \frac{1}{{\eta }^{1-\gamma }}\left\{-\tilde{G}\left(t\right)s\left(t\right)+\chi \left[k_1{\left|s\left(t\right)\right|}^{p+1}+k_2\left|s\left(t\right)\right|\right]\right\}& \left|e\left(t\right)\right|<\eta \end{array}\right.$$

(87)

Then, $k_1{\left|s\left(t\right)\right|}^{p+1}+k_2\left|s\left(t\right)\right|>0.$

$q$ is a positive odd integer and $r$ is a positive integer, $r\ge 1$, $0\le p<1$ and Let $\chi <0.$

when $\left|s\left(t\right)\right|<1$, then ${\left|s\left(t\right)\right|}^{p+1}\le \left|s\left(t\right)\right|$, $s^{q+1}\left(t\right)<\left|s\left(t\right)\right|$, so, $\left|\chi \right|$, $\left|k_1\right|$, $\left|k_2\right|$ are large enough for $\mathrm{d}V/\mathrm{d}t<0$;

when $1\le \left|s\left(t\right)\right|<r$, then $\left|s\left(t\right)\right|\le {\left|s\left(t\right)\right|}^{p+1}$, $\left|s\left(t\right)\right|\le s^{q+1}\left(t\right)$, so, $\left|\tilde{G}\left(t\right)\right|<\left|\chi k_1\right|$ or $\left|\tilde{G}\left(t\right)\right|<\left|\chi k_2\right|$ or $\left|\tilde{G}\left(t\right)\right|<\left|\chi k_1\right|+\left|\chi k_2\right|$ for $\mathrm{d}V/\mathrm{d}t<0$;

when $r\le \left|s\left(t\right)\right|$, then $\left|s\left(t\right)\right|\le {\left|s\left(t\right)\right|}^{p+1}$, $\left|s\left(t\right)\right|\le s^{q+1}\left(t\right)$, so, $\left|\tilde{G}\left(t\right)\right|<\left|\chi k_1\right|$ or $\left|\tilde{G}\left(t\right)\right|<\left|\chi k_2\right|$ or $\left|\tilde{G}\left(t\right)\right|<\left|\chi k_1\right|+\left|\chi k_2\right|$ for $\mathrm{d}V/\mathrm{d}t<0$;

Hence, the boundary condition: $\left|\chi \right|$, $\left|k_1\right|$, $\left|k_2\right|$ are large enough or $\left|\tilde{G}\left(t\right)\right|<\left|\chi k_1\right|$ or $\left|\tilde{G}\left(t\right)\right|<\left|\chi k_2\right|$ or $\left|\tilde{G}\left(t\right)\right|<\left|\chi k_1\right|+\left|\chi k_2\right|$ for $\mathrm{d}V/\mathrm{d}t<0$ and the composite controller is stable. □

5. Simulation Results

To further validate some superiorities of the proposed composite controller, some specific parameters are set, including the TSTQLs system’s parameters and some controllers’ parameters in

Table 7. All parameters of TSTQLs system and controllers.

${\mathit{U}}_{\mathit{I}}$	$\mathit{R}$	$\mathit{c}$	$\mathit{m}$		$\mathit{A}$
$220 \mathrm{V}$	$3.08 \mathsf{\Omega }$	$130 \mathrm{J}/\mathrm{kg}·\mathrm{K}$	$1.46\ast {10}^{-2} \mathrm{kg}$		$2.9\ast {10}^{-3} {\mathrm{m}}^2$
$\beta$	$\lambda$	$\epsilon$	$\sigma$		$F$
$11.6 \mathrm{W}/{\mathrm{m}}^2·\mathrm{K}$	$174 \mathrm{W}/\mathrm{m}·\mathrm{K}$	$0.97$	$5.67\ast {10}^{-8} \mathrm{W}/{\mathrm{m}}^2·{\mathrm{K}}^4$		$1$
${\beta }_1$	${\beta }_2$	$k$	$\epsilon$	$\gamma$	$\eta$
$70$	$20$	$20$	$0.5$	$0.8$	$4$
$k_1$	$k_2$	$p$	$q$	$r$	$\chi$
$8$	$200$	$0.9$	$5$	$0.3$	$-3000$

.

To obtain proper parameters of the composite controller, some steps are explained. The sampling period used in the simulation tests is 30 s. Step 1, tune $\chi$ via increasing its value from negative to positive until the control performance is degrade; Step 2, maintain $\chi$ and tune $k_1$, $k_2$ via increasing their values and satisfy $k_1>0$ and $k_2>0$; Step 3, tune $p$, $q$ by satisfying $0\le p<1$ and $q$ is a positive odd integer; Step 4, tune $k$, $\eta$, $\gamma$ by satisfying $0<\gamma <1$ and $\eta >0$ until the tracking errors decrease and the chattering trend declines; Step 5, if these above steps are successful, maintain the same parameters of the ISMCNESO controller and NGSMCNESO controller, and tune the other parameters to achieve a great performance as comparison.

As is shown in Figure 14, Figure 15, Figure 16, Figure 17, Figure 18, Figure 19, Figure 20, Figure 21 and Figure 22, there are three comparative control methods, including the ISMCNESO controller, the NGSMCNESO controller, and the composite controller. The ISMCNESO controller is based on a traditional integral sliding mode surface (38) with a NESO (21). The NGSMCNESO controller and the composite controller have the same nonlinear sliding mode surface (50). In Figure 14a, the fitting curve of the reference temperature trajectory from Wall 0_0 mm, named aim 1, is from 247.596 K to 1595.98 K between 0 s and 30 s and the sequence of other aims is from aim 2 to aim 9 corresponding to Figure 15a, Figure 16a, Figure 17a, Figure 18a, Figure 19a, Figure 20a, Figure 21a and Figure 22a. From Figure 14b, Figure 15b, Figure 16b, Figure 17b, Figure 18b, Figure 19b, Figure 20b, Figure 21b and Figure 22b, the reference temperature tracking performances among the ISMCNESO controller, NGSMCNESO controller, and composite controller have a similar trend. However, in terms of tracking errors and their partial enlarged figures, the composite controller has the least chattering phenomena compared with the ISMCNESO controller and NGSMCNESO controller because of IPD (72) and the reaching law (78). The MFC provides the model free nature for independence of the system dynamic model and the reaching law (78) suppresses chattering phenomena from the high-frequency switching motions. In addition, the less chattering phenomena occurs in the NGSMCNESO controller than that in ISMCNESO owing to the nonlinear Function (51).

To analyze tracking performance of the composite controller quantitatively, the calculation results of RMSE and MAX are in

Table 8. The calculation results of RMSE and MAX.

	ISMCNESO Controller		NGSMCNESO Controller		Composite Controller
	RMSE	MAX	RMSE	MAX	RMSE	MAX
Wall 0_0 mm	0.486	1.000700	0.161	0.410	0.00682	0.174
Wall 0_3 mm	0.491	1.00300	0.132	0.431	0.00634	0.192
Wall 0_6 mm	0.257	0.924	0.109	0.400	0.00235	0.0628
Wall 1_5 mm	0.534	1.0100	0.117	0.402	0.00532	0.135
Wall 1_45 mm	0.136	1.0531	0.103	0.400	0.00204	0.0993
Wall 1_85 mm	0.116	1.0519	0.0968	0.383	0.00202	0.0580
Wall 2_5 mm	0.162	1.0678	0.111	0.406	0.00291	0.0799
Wall 2_45 mm	0.118	1.00392	0.0887	0.346	0.00277	0.0574
Wall 2_85 mm	0.0832	0.922	0.0759	0.393	0.00283	0.0520

, which are root mean square error $\mathrm{RMSE}=\sqrt{{{\displaystyle \sum }}_{i=1}^N{e}_i^2/N}$ and maximum error $\mathrm{Max}=Max\left|e_i\right|\left|i=1~N\right.$, respectively. In Table 7, The RMSE values of these three controllers: composite controller < NGSMCNESO controller < ISMCNESO controller and the MAX values are in the same sequence. The quantitative calculation results can effectively demonstrate the superiorities of the nonlinear sliding mode surface (50) with the nonlinear Function (51), MFC (71), and the reaching law (78).

To study the robustness of the composite controller, the time-varying resistance is chosen as a kind of external disturbance: R1, R2, R3. As is shown in Figure 23, the value of R1 is set as a step signal which is 100 Ω from 14 s to 16 s and 3.08 Ω in the other time. The values of R1 and R2 are also set as a step signal with 300 Ω and 500 Ω, respectively, and that of the other time is 3.08 Ω. From Figure 24, Figure 25, Figure 26, Figure 27, Figure 28, Figure 29, Figure 30, Figure 31 and Figure 32, the larger the value of the step signal is, the lager the amplitude of chattering will be. For example, the overshoot of R3 for the composite controller is 0.25 K in Figure 24b and the system state reaches the steady state after 1 s. At 16 s, the value of R3 is from 500 Ω to 3.08 Ω, and the overshoot of R3 for the composite controller is about 2.5 K in Figure 24a. Although the external disturbances have two sudden changes which have an adverse impact on tracking performance, the composite controller shows a strong robustness. Considering the effect of some measurement noises, the external signal is set as a step signal which is 100 Ω from 4 s to 5 s and from 25 s to 26 s on the basis of external disturbances before (R1, R2, R3). The tracking errors of Wall 2_45 mm and Wall 2_85 mm for the composite controller are shown in Figure 33. At the beginning of 4 s and 25 s, the tracking error trend has a slight chattering and converges to zero quickly. Hence, the composite controller is insensitive to measurement noises.

6. Conclusions

In this paper, a general framework for the TSTQLs for the rationality evaluation of a hypersonic vehicle’s TPS is introduced, including two key roles: numerical analyses and control methods. Numerical analyses of a hypersonic flying object’s aerodynamic heating environment are based on three different two-dimensional outflow fields via finite element calculations in ANSYS Workbench 2020 R2. The reference temperature trajectories are plotted to scatter graphs and the fitting curves of scatter graphs are obtained from three calculated objects: Wall 0, Wall 1, and Wall 2. In order to track these reference temperature trajectories, a composite controller is proposed for a real-time ground aerodynamic heating simulation of a hypersonic flying object, named the TSTQLs. The proposed composite controller is made of IPD, NGSMC, and NESO. In an ultra-local model of MFC frame, the system dynamic model of the TSTQLs is replaced by IPD controller, and the NESO is as the observation for the lumped disturbances from internal parametric uncertainties, periodic oscillations, strong nonlinearity, and external disturbances. The NGSMC with a nonlinear function is an auxiliary controller to achieve a stage of great tracking errors and strong robustness. Moreover, the NGSMC eliminates the reaching phase suppressing chattering phenomena from the high-frequency switching motions. In addition, the reaching law is the combination of an exponential reaching law and a power reaching law and further improves the dynamic performance of chattering and robustness. Finally, the comparative simulation results demonstrate some superiorities of the proposed composite controller for the TSTQLs.

For the proposed composite controller strategy, this control method brings so many parameters for tuning. In the future, we will focus on the adaptive control. During a real flight process of the hypersonic vehicle, the hypersonic vehicle will also face serious high-frequency vibration and high-intensity sound fields, and hence aerodynamic heating, vibration, and sound will be considered synthetically, and we will also focus on the interaction mechanism among sound, vibration, and heating. Moreover, the heating environment is in time sequence, and then the group control needs to be considered for the TSTQLs with the circular quartz lamp group.