Robust Control of Pressure for LNG Carrier Cargo Handling System via Mirror-Mapping Approach

The pressure control for a liquid nature gas (LNG) ship is vital for the cargo handling system which is a nonlinear, unstable, and controllable complex system accompanied by the dynamics of time delay. To improve the control effect, this article proposed a robust controller based on closed-loop gain shaping algorithm by mirror-mapping approach. In addition, the complete mechanism models for the system are established to predict the changes of temperature, pressure, and liquid inventory in the cargo tank. The heat exchanges and evaporation of liquid are also considered. By collecting the data from the LNG ship “Dapeng,” the system models are validated. At the same time, the comparative experiment is introduced to verify the reliability and effectiveness. The scheme has been compared with the 2-DOF structure control law of modified Smith predictor. The comparative experimental results show that the scheme proposed in this note has a strong disturbance rejection ability and steady-state performance. The controller designed in this note has advantages of simplified construction method, satisfactory control effect and robustness.


Introduction
The use of natural gas has been increasing rapidly all over the world because the gas contains low amounts of pollutants.LNG shipping is therefore an economic way of transporting large quantities of natural gas over long distances.LNG carriers are purpose-built tank vessels for transporting LNG at sea [1,2].During handling, LNG is transferred between the onshore storage tanks and ship tanks at high flow rates through single or parallel pipelines.It is always influenced by external disturbances outside.In general, these operations are energy intensive and involve stringent safety considerations [3][4][5].It is very important and critical to control the temperature and pressure in order to reduce the risk of accidents.Thus, dynamic integration models of the handling operation are needed to be established first.The model can also predict the tank level, temperature, and pressure within design limits and can help to prevent overfilling or excessive unloading [6].Temperature and pressure are connected with the gas mass of the control volume inflow or outflow in unit time.If the control performance is unsatisfactory, the pressure and flow rate will greatly vary which heavily influences the operation.In the worst case scenario, equipments such as the LNG pumps will fall below the required minimum flow rate and trip, resulting in suspended handling or LNG disaster [7,8].
Therefore, it is very important to design a robust controller to maintain pressure stability in the system.External disturbances must be taken into account in the controller design so that the closed-loop system has a good response even in such dynamics [9].Robustness is critical to the pressure control system.In [10], closed-loop gain shaping algorithm (CGSA) was put forward based on robust control theory.The advantages were obvious physical meaning and simple design procedure.With the algorithm, the controller was designed and applied to the nonlinear model of rudder rolldamping SIMO system [11].In [12], the robust adaptive control design problem was studied for a class of nontriangular nonlinear systems with unmodeled dynamics and stochastic disturbances and a robust adaptive output feedback control scheme was then proposed.An adaptive neural networkbased fin control design method was proposed in [13], which combined the adaptive backstepping technique with the neural network.Guan et al. [14] proposed a robust control scheme based on the closed-loop gain shaping algorithm that offered a simplifying controller for processes with time delay.Padhan and Majhi [15] proposed a modified Smith predictor control scheme for controlling stable and integrating and unstable processes with time delay.
Recently, there has been a plentiful interest focusing on the problem of robust control design for unstable processes or nonlinear system.Mirror-injection technology was first proposed in [16] which could map the unstable process or unstable system to a stable one.To control these mapped stable one, the same effectiveness could be achieved.In [17], robust control method was addressed for a class of pure unstable system with time delay by mirror-mapping technique.A general scheme was proposed in [18] for the high order unstable delay process with one or more positive poles, using the mirror mapping technique.Robustness analysis and design for the integrating unstable delay systems were discussed in [19].From the above studies, it is obvious to note that the controller of pressure maintenance system for the LNG carrier is not discussed.
The work in this paper is inspired by [16][17][18][19].This applies the mirror-mapping technique and CGSA to LNG carrier pressure control system that is an unstable nonlinear system with time delay.In this paper, a set of dynamic models are established by analyzing the equipment conditions in LNG cargo handling system.Then, verification of the model with real ship data is provided.In order to explain the controller design, we take the insulation space as a control plant.Linearization of the model is required.The system is converted into a regular unstable process.Consequently, stability of the system is analyzed.Then the loop-shaped unstable process is transformed into a minimum-phase system through mirror-mapping technique.Robust controller with hybrid control strategy is designed based on the mirrormapping with second-order CGSA.The proposed control scheme is compared with 2-DOF control law.The transient and steady-state performance are analyzed when different disturbance is introduced into the system.The simulation results show that the scheme proposed in this paper has better robustness and disturbance rejection.

System Description
LNG carrier cargo operations includes the following steps: (1) Insulation space inerting; (2) cargo tank drying; (3) cargo tank inerting; (4) cargo gassing-up; (5) cargo tank cooling-down; (6) losading; (7) loaded voyage; (8) ballast voyage; (9) unloading; (10) cargo tank warm-up; (11) and cargo tank aeration.A complete description of the process can be found at [6].Here, it only focuses on the loading and unloading process.During the loading process, the boil-off gas is forced to return to the shore tank by highduty (HD) compressors.HD heaters heat the LNG vapor delivered by the HD compressors to the specified temperature.Low-duty (LD) compressors are provided in the cargo machinery room for maintaining constant cargo tank pressure and delivering boil-off gas according to the boiler demand.It is necessary to keep the cargo tank pressure between 90~120 mbar.While during the discharging process, LNG vapor is supplied from the shore to maintain pressure in the cargo tanks.If the shore does not supply returnvapor to the cargo tank, the LNG vaporizer produces vapor by bleeding LNG from the main line and supplies it to the cargo tanks.It need not to start HD compressors in this process.So the plants in the LNG carrier cargo handling system are mainly composed of storage tanks, pumps, compressors, gas pipe line, and valves.The outline of the whole plants is shown in the flowing part of Figure 1.Some duplicate plants of the process are removed for clarity.
For different type of the LNG ship, the cabin pressure range will be different.Usually, it is important to maintain the tank's pressure at least 100 mbar in order to avoid cavitations and ensure having good suction at the pumps.If the tank's pressure falls to 60 mbar, the returned gas is requested to increase from the shore.If two main cargo pumps are in use in a tank, the discharge valve on one pump should adjust to 40% and pump should be stopped when the level reaches 1.1 m.This is in order to reduce turbulence around the pump suction.

Dynamic Model
In this section, a dynamic model for the cargo handling operation is developed.The dynamic model for the cargo handling operation involves modeling the heat transfer and the dynamics of the liquid and vapor in the storage tank.Several main models and connections are shown in [20,21].The development of the dynamic model for each of these components is described below.
3.1.Pump.The liquid pump head depends on the characteristics of the pump itself and the characteristics of the piping.The head of liquid pump can be solved by the Bernoulli equation between the center of the pump inlet and the suction liquid surface.Equation (1) gives the formula for direct solution.
where H is the head of the pump.The unit is m.p dr is the pressure of the pump discharges the liquid surface, and the unit is kPa.p sr is the suction surface pressure, and the unit is kPa; ρ is the density of cargo pumped cargo, and the unit is kg/m 3 .g is gravity acceleration, and the value is 9.8 m/s 2 ; h = h s + h d + Δh is the total height difference between the suction level and the discharge level, and the unit is m; ∑f = ∑f s + ∑f d is the total resistance of the pipeline between suction and discharge, and the unit is m.Using polynomial curve fitting technique, the HQ characteristic curve of the liquid cargo pump is a secondorder polynomial as shown in below.
Here, H is the total head of the pump, and Q is the total flow of the pump.a 1 , a 2 , and a 3 are the constant determined by the characteristics of the cargo pump.
Equation (3) gives the characteristics of the supply pipeline of a single liquid cargo pump.
where H st is the static head which is independent of flow.K is the proportional constant depending on the resistance characteristics of the pipeline.

Pipeline Resistance Calculation Model.
According to the length and diameter distribution of the pipeline system, ( 4) is used to distribute the gas pressure difference between the gas supply and tank pressure of the liquid pump head reasonably.3 Complexity of resistance along the way which is related to the viscosity of the fluid, the flow rate, the inner diameter of the pipe, and the roughness of the pipe wall.

Tank Mechanism Model.
A nonlinear function curve of the liquid surface height and the compartment volume can be obtained by polynomial curve fitting or piecewise linear interpolation.The mathematical model of tank volume and level can be written as where V LNG is the volume of LNG gas (m 3 ).V LNG0 is the volume of LNG gas at initial state (m 3 ).h is the height of liquid surface (m).Q is flow rate of the liquid into the cargo tank per unit time (m 3 /h).The leakage of external heat into the cargo tank is accompanied by three heat exchanging stages: conducted heat transfer, convection heat transfer, and radiation heat transfer.The heat transfer equation in the cargo tank is shown below.
In the formula, q represents the heat flux (which is passed from thermal fluid outside to cold fluid inside), W/m 2 ; r k represents the total thermal resistance of the series heat transfer process, (m 2 • °C)/W; k indicates the heat transfer coefficient, W/(m 2 • °C).t g indicates the temperature of the thermal fluid (water, air) outside the tank, °C; t d indicates the temperature of the cold fluid (gas cargo, liquid cargo) in the tank, °C.

Mathematical Model of Tank Temperature and Pressure.
During loading and discharging, the change of pressure and temperature is unstable with the heat exchanging.In order to model the pressure and temperature in the cargo tank, the upper vapor (ideal gas) space of the cargo tank is taken as the control volume.According to the hypothesis [22], the heat transferred to the cargo tanks is used for the evaporation of the surface fluids in the tanks.According to the first law of thermodynamics, energy is conserved.Thus, the energy equation of the control volume is expressed in (7).
Therefore, in the loading and unloading operation, the heat from the external into the cabin can also be regarded as the heat exchange with the cabin air. where In addition, the vaporization of liquid cargo in cargo tanks will absorb a part of the heat.According to the literature [23], the amount of heat transferred into the cargo tanks in the unit time can be expressed in ( 8), (9), and (10) as follows.
where the values of k 1 and k 2 are determined by the design parameters of the cargo tank and the external environmental conditions.z is the latent heat of vaporization of LNG (cal/g).c p the specific heat of constant pressure at temperature T (cal/(mol•K)).c p1 is the specific heat of constant pressure at the temperature T 1 , cal/(mol•K).Equation (11) indicates the energy change in the control volume per unit time.According to the assumption that the gas in the control volume is ideal gas, ( 12) is established by the basic theory of engineering thermodynamics.
In (13), m is the mass of the gas in the tank (kg); c p and c v meet the following relationship [24,25] as follows.

14
where A, B, C, and D are constants in the literature [25].Considering that the average daily vaporization rate of the tank is 0.15%, the critical temperature under the tank pressure still needs to be considered.The tank temperature mathematical model can be written as The mathematical model of exhaust flow is shown in (16).That means that the pressure relief valve will work when the pressure in the tank exceeds the set value.
Equation ( 17) can be used to solve the current gas quality in the tank.
where m 0 denotes the gas quality in the tank at the previous moment; m in denotes the gas quality entering the tank; m out is the quality of discharges gas.It should be noted that m in and m out have different sources in different cargo handling operations.With (8), ( 13), (15), and ( 17), the temperature and gas quality can be obtained as function of time.The mass of gas flowing in and out of the tank can be expressed by the equation of state PV = mrT and r = R/M.According to the ideal gas state equation at any moment, the pressure can be obtained as function of time.

Model Validation
To validate the above model, data was collected from a LNG ship named "Dapeng."Substitute the values of parameter in Table 1 and Table 2 into the equations above, the liquid level, temperature, pressure, and liquid volume of tanks can be obtained as functions of time.
All the constant parameters that appear in the formula are listed in Table 1.Other parameters depend on the material properties of different ships.These parameters for specific conditions from "Dapeng" are listed in Table 2.
The actual data of pressure, temperature, and liquid volume during loading and unloading are also collected to verify model outputs in this note.The results of verification are shown in Figures 2(a)-2(c).The curves of gas temperature, tank pressure, and liquid volume in the tank are achieved below.
It is shown from Figures 2(a) and 2(b) that the gas temperature and volume are consistent with the actual situation.At the beginning, the temperature rises, because the ingress heat is introduced into the cargo tank.Subsequently, evaporation is caused by the large return-gas; thus, the temperature tends to drop slowly.This is suitable for large LNG ship.The results is consistent with the actual situation.
For the pressure shown in Figure 2(c), it tends to be consistent at the beginning.After that, the pressure will gradually increase.This is because the dynamic models are the uncontrolled model.If the controller is added, the pressure will be within a certain range.This will be studied in the following part of this note.There are several reasons for this inaccuracy of model predictions such as compressor performance curve used in compressor model is not exact.Furthermore, a single heat transfer coefficient is used in a particular heat exchanger.However, this difference in prediction is acceptable since this is not huge.

Robust Control Studies
5.1.Control Plant Description.Take the model of insulation space for example to discuss the pressure control.A controller is designed to control the pressure difference in the insulating layer.The pressure in insulation space is required between 0.2 kPa and 0.4 kPa.The pressure is maintained by inputting nitrogen or withdrawing nitrogen from the space.The relationship between pressure difference and nitrogen supply rate is expressed as follows.
Where P e is the pressure difference of the primary and secondary insulating layers.u is nitrogen supply rate.ρ N 2 is the nitrogen density supplied, kg/m 3 ; g is the gravitational acceleration; F r is the damping force under nonideal conditions; κ is a constant, determined by the material properties of the insulating layer.According to [26], (18) can be linearized near the equilibrium point P e0 = 2.0 mbar, u 0 = 42 m 3 /h.A system transfer function of linearization is given below.
The time delay e −0 25s of the insulating space pressure difference system G s is less than 0.5 s, so the 1st order Padé approximation of the time delay can be obtained as follows.
G 0 s = 3 25 s s − 0 78 21 Let s s − 0 78 0 25s + 1 = 0.The root of characteristic equation 0.78 > 0. There is the right half plane (RHP) pole at s = 0.78 in (22); hence, the system is unstable.It is a nonminimum-phase plant (NMP) system with an unstable controlled plant, which is shown in Figure 3.

System Stabilization Process.
It is necessary to investigate stability and limitations on performance imposed by RHP zeros and poles.Mirror-mapping technique provides a novel design tool in [16][17][18][19].The basic principle is that the largest singular value curve of the mirror-mapping process is similar to one of the original unstable process.Hence, the results from mirror-mapping can also stabilize the original system.Definition 2. The process that zero-poles in the open right half plane of unstable process are substituted for their mirrorimages and then a minimum-phase system is constructed is called mirror-injection process.
In order to obtain better results, normalizing ( 23) is as follows.
From the bode plots of the two system models in Figure 4, it can be seen that their spectrums are similar.The root locus of G s is shown in Figure 5.There are no poles in the RHP.The closed-loop system step response curve of the controlled plant without perturbation is shown in Figure 6.The system steady-state error is 0, so the system is stable.

Controller Design.
The controller is designed with closed-loop gaining shaping algorithm (CGSA).CGSA is a simplified H ∞ mixed sensitivity algorithm [10] by shaping directly the singular value curves of S (the sensitivity function) and T (the complementary sensitivity function).According to the H ∞ robust control theory, the closedloop frequency spectrum (i.e., equal to T) of a typical control system has a low pass characteristics to guarantee robust performance.The frequency asymptote of the frequency spectrum determines how much the system is sensitive to the frequency outside the valid frequency bandwidth, that is, the disturbance frequency.In order to    7 Complexity obtain better robustness and appropriate order form, the slope −40 dB/dec is adopted which can be derive from Figure 4, so bandwidth frequency is 1/T 1 (here T 1 = 0.2).
The controller is designed through shaping the frequency spectrum curve of T which can be approximately regarded as the frequency spectrum curve of a secondorder inertial system with the largest singular value of 1.This is equal to the condition that the damping ratio is 1 compared to a standard oscillating second-order system, so it guarantees that there is no peak value in the frequency spectrum of T. According to the closed-loop gain shaping algorithm in [10], then (25) is obtained.
The S&T singular value curves of the controlled plant are shown in Figure 7.It is shown that the designed controller has better robustness from the aspect of H ∞ control mix sensitivity algorithm.
To eliminate the steady-state error, a minor constant ɛ (0.001) is added into the denominator of the G 0 s Then ( 26) is obtained.
Thus, substitute ( 26) into ( 24) and ( 25), according to closed-loop gain shaping algorithm, a linear proportionalintegral-differential (PID) plus the first order filter controller is obtained.
Actually, ɛ < 0.001 is the smallest constant; ɛ can derive the integral effect for the designed controller and it also reproduces the effect of uncertain constant disturbance upon the closed-loop system.In order to eliminate the effect of G 1 s , (28) is obtained.

28
The block diagram of pressure difference control for insulating spaces is shown in Figure 8. G s in ( 19) is substituted into Figure 8.

Simulation and Analysis
In order to verify the effectiveness of the controller in this paper, a comparative experiment is conducted.The compared control scheme comes from the literature [27].The compared law is derived from (29) which is designed as 2-DOF structure based on modified Smith predictor and includes two parts: K f f s and K s .K f f s = 1 982 0 20s 2 + 0 77s + 1 5 342 0 25s 3 + 0 99s 2 + 1 64s , The block of 2-DOF structure control laws based on modified Smith predictor is shown in Figure 9. F s is a filter.
The comparative results are given in Figures 10-12.The real line represents control performance by using the proposed method.The dashed line describes the control performance in [27].Figure 10 shows that the rise time of response 8 Complexity (T r ) is shorter than the later one.Also, the percent overshoot in the proposed scheme is 2.51% lower than that in 2-DOF scheme.The setting time (T s ) to reach its output final value is obviously shorter than the later one.The compared results are shown in Table 3.
In order to verify the robustness of the proposed scheme, disturbances are added at 30 s. Disturbance of the step signal −1 is added at 30 s as shown in Figure 11.The disturbance of ramp signal is added in Figure 12.The comparison of output responses is presented in Figures 11 and 12.
From the results of the simulation above, it is obvious that the proposed scheme in this paper can tackle the load disturbances along with the control input more efficiently.In addition, the time delay is treated by one order Padé approximation.But the simulation uses the original time delay term.That means the controlled plant has a certain perturbation.The results show that the proposed control scheme has better control performance and robustness.The designed controller order use lower order and easy to be realized.The control law proposed in this paper has a quick adjust in time and strong anti-interference ability and can quickly recover to the initial control state after the disturbance is introduced.The control performance is satisfactory and the proposed control method is robust.

Conclusion
This paper establishes a mechanism model based on thermodynamics and dynamics for LNG carrier cargo handling system.Models are validated with actual data from LNG ship "Dapeng."Results show that the model outputs are basically consistent with actual situation, and the system is controllable.To control the gas pressure, a robust controller is designed with the closed-loop gain shaping algorithm by  9 Complexity mirror-mapping method.Insulation space is taken as controlled plant.Comparative simulations are conducted by comparing with the 2-DOF structure based on modified Smith predictor.The results show that the scheme proposed in this note has stronger robustness and steady-state performance.The mathematical models established in this paper can be applied to the LNG carrier handling simulator.And it also helpful to shipbuilding and optimizing the process of handling.The control scheme proposed in this paper can be used for other system models and provides a good reference for LNG control system design.The rise time of response is defined as T r .The peak value of the step response is denoted by M pt .C ss is the value of stability sate.The time to reach peak value is T p .Overshoot is defined by the equation: σ p = M pt − C ss /C ss × 100.These values are calculated and given in Table 3.
Liquid volume verification with real data.

Figure 2 :
Figure 2: The comparison with real ship.
Definitions 1 and 2 give the brief description on the mirror mapping technique.Definition 1.The symmetric values of the zero-poles in the open right half plane of unstable process about the imaginary axis are referred to as their mirror-images.

Figure 3 :
Figure 3: Root locus of zero-pole for the unstable system.

Figure 5 :
Figure 5: Root locus of zero-pole of stable system.

Figure 4 :
Figure 4: Bode plots of the two system models.

Figure 9 :Figure 10 :
Figure 9: The block diagram of pressure control for the insulation space.
Complexitywhere dm 1 /dt = dm 10 /dt + dm 11 /dt , dm 10 /dt is the mass of gas flowing into the control volume per unit time, and dm 11 /dt is the mass of vaporization gas in tank per unit time; T 10 is the temperature of the gas entering the control volume; T 11 is the critical temperature at the current pressure; c p10 and c p11 are the specific heat of constant pressure at temperature T 10 and T 11.

Table 1 :
Nomenclature of the parameters appearing in the formula.

Table 2 :
Parameters of liquid cargo tanks for thermal load calculation.Tanker waterline external surface area S 2 3968 m 2The surface area of the liquid cargo ship S 3 3060 m 2

Table 3 :
Parameters of liquid cargo tanks for thermal load calculation.