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

Authors: Jinghua Cao; Xianku Zhang; Guangping Yang; Xiang Zou
Journal: Complexity
 (2018)

---

## Abstract

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.


---

## Body

## 1. 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–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 roll-damping 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 network-based 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–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 mirror-mapping 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.


## 2. 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 high-duty (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 return-vapor 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.Figure 1
LNG cargo handling system plant.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.


## 3. 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.
(1)H=pdr−psrρg+hs+hd+Δh+∑fs+∑fd=pdr−psrρg+h+∑f,where H is the head of the pump. The unit is m. pdr is the pressure of the pump discharges the liquid surface, and the unit is kPa. psr is the suction surface pressure, and the unit is kPa; ρ is the density of cargo pumped cargo, and the unit is kg/m3. g is gravity acceleration, and the value is 9.8 m/s2; h=hs+hd+Δh is the total height difference between the suction level and the discharge level, and the unit is m; ∑f=∑fs+∑fd 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 second-order polynomial as shown in below.(2)H=a1Q2+a2Q+a3.Here,H is the total head of the pump, and Q is the total flow of the pump. a1, a2, and a3 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.
(3)Hpipe=Hst+∑f=h+pdr−psrρg+KQ2,where Hst is the static head which is independent of flow. K is the proportional constant depending on the resistance characteristics of the pipeline.


### 3.2. 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.
(4)hfliquid=ξldv22g,Δpfgas=ξldρv22,where hfliquid represents the head loss of the liquid along the way (m). Δpfgas represents resistance loss of the gas along the way (m); l represents pipeline length (m); d represents pipe diameter (m); v represents the average fluid velocity (m/s); g represents gravity acceleration(m/s2); ξ is the coefficient 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.


### 3.3. 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(5)VLNG=VLNG0+∫0tQdt,h=fnonVLNG,Vliquid=V−VLNG,where VLNG is the volume of LNG gas (m3). VLNG0 is the volume of LNG gas at initial state (m3). h is the height of liquid surface (m). Q is flow rate of the liquid into the cargo tank per unit time (m3/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.
(6)q=tg−tdrk=ktg−td.In the formula,q represents the heat flux (which is passed from thermal fluid outside to cold fluid inside), W/m2; rk represents the total thermal resistance of the series heat transfer process, (m2·°C)/W; k indicates the heat transfer coefficient, W/(m2·°C). tg indicates the temperature of the thermal fluid (water, air) outside the tank, °C; td indicates the temperature of the cold fluid (gas cargo, liquid cargo) in the tank, °C.


### 3.4. 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.(7)dqdt+h1dm1dt−h2dm2dt=dUdt+dWdt,where dq/dt indicates the heat transfer rate from outside into the tank (cal/s); h1dm1/dt is the heat brought into control volume in the unit time (cal/s); h2dm2/dt is the heat taken away from control volume in the unit time (cal/s); du/dt is the change of internal energy in the control volume per unit time (cal/s); dW/dt is heat consumption for external work in control volume per unit time (cal/s). Here, it is a rigid volume, dW/dt=0; let the gas temperature in the tank be T, then the relationship dq/dt and T can be expressed as dq/dt=k1−k2T.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.
(8)dqdt=k1−k2T,(9)h1dm1dt=cp1T1dm1dt−zdm1dt,(10)h2dm2dt=cpTdm2dt,where the values of k1 and k2 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). cp the specific heat of constant pressure at temperature T (cal/(mol·K)). cp1 is the specific heat of constant pressure at the temperature T1, 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.
(11)dU=dmu=mdu+udm=mcvdT+cvTdm1−cvTdm2,(12)dUdt=mcvdT+cvTdm1−cvTdm2dt,where m is the mass of natural gas in the tank (kg); cv is the constant volume specific heat of natural gas (cal/(mol·K)). Substituting (7) with (8), (9), (10), and (12), (13) is achieved.
(13)mcvdTdt=cvdm2dt−cpdm2dt−cvdm1dt−k2T+k1−zdm1dt+cp1T1dm1dt.In (13), m is the mass of the gas in the tank (kg); cp and cv meet the following relationship [24, 25] as follows.
(14)cp−cv=Rg,cp=A+BT+CT2+DT3,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
(15)mcvdTdt=k1−zdm11dt+cp11T11dm11dt+cp10T10dm10dt−k2+cpdm2dt+cvdm1dt−cvdm2dtT,where dm1/dt=dm10/dt+dm11/dt, dm10/dt is the mass of gas flowing into the control volume per unit time, and dm11/dt is the mass of vaporization gas in tank per unit time; T10 is the temperature of the gas entering the control volume; T11 is the critical temperature at the current pressure; cp10 and cp11 are the specific heat of constant pressure at temperature T10andT11.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.
(16)qs=PRTP′/8ξρπd2,ρ=PT0P0TM22.4.Equation (17) can be used to solve the current gas quality in the tank.
(17)m=m0+min−mout,where m0 denotes the gas quality in the tank at the previous moment; min denotes the gas quality entering the tank; mout is the quality of discharges gas. It should be noted that min and mout 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.


## 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.
(1)H=pdr−psrρg+hs+hd+Δh+∑fs+∑fd=pdr−psrρg+h+∑f,where H is the head of the pump. The unit is m. pdr is the pressure of the pump discharges the liquid surface, and the unit is kPa. psr is the suction surface pressure, and the unit is kPa; ρ is the density of cargo pumped cargo, and the unit is kg/m3. g is gravity acceleration, and the value is 9.8 m/s2; h=hs+hd+Δh is the total height difference between the suction level and the discharge level, and the unit is m; ∑f=∑fs+∑fd 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 second-order polynomial as shown in below.(2)H=a1Q2+a2Q+a3.Here,H is the total head of the pump, and Q is the total flow of the pump. a1, a2, and a3 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.
(3)Hpipe=Hst+∑f=h+pdr−psrρg+KQ2,where Hst is the static head which is independent of flow. K is the proportional constant depending on the resistance characteristics of the pipeline.


## 3.2. 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.
(4)hfliquid=ξldv22g,Δpfgas=ξldρv22,where hfliquid represents the head loss of the liquid along the way (m). Δpfgas represents resistance loss of the gas along the way (m); l represents pipeline length (m); d represents pipe diameter (m); v represents the average fluid velocity (m/s); g represents gravity acceleration(m/s2); ξ is the coefficient 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.


## 3.3. 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(5)VLNG=VLNG0+∫0tQdt,h=fnonVLNG,Vliquid=V−VLNG,where VLNG is the volume of LNG gas (m3). VLNG0 is the volume of LNG gas at initial state (m3). h is the height of liquid surface (m). Q is flow rate of the liquid into the cargo tank per unit time (m3/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.
(6)q=tg−tdrk=ktg−td.In the formula,q represents the heat flux (which is passed from thermal fluid outside to cold fluid inside), W/m2; rk represents the total thermal resistance of the series heat transfer process, (m2·°C)/W; k indicates the heat transfer coefficient, W/(m2·°C). tg indicates the temperature of the thermal fluid (water, air) outside the tank, °C; td indicates the temperature of the cold fluid (gas cargo, liquid cargo) in the tank, °C.


## 3.4. 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.(7)dqdt+h1dm1dt−h2dm2dt=dUdt+dWdt,where dq/dt indicates the heat transfer rate from outside into the tank (cal/s); h1dm1/dt is the heat brought into control volume in the unit time (cal/s); h2dm2/dt is the heat taken away from control volume in the unit time (cal/s); du/dt is the change of internal energy in the control volume per unit time (cal/s); dW/dt is heat consumption for external work in control volume per unit time (cal/s). Here, it is a rigid volume, dW/dt=0; let the gas temperature in the tank be T, then the relationship dq/dt and T can be expressed as dq/dt=k1−k2T.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.
(8)dqdt=k1−k2T,(9)h1dm1dt=cp1T1dm1dt−zdm1dt,(10)h2dm2dt=cpTdm2dt,where the values of k1 and k2 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). cp the specific heat of constant pressure at temperature T (cal/(mol·K)). cp1 is the specific heat of constant pressure at the temperature T1, 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.
(11)dU=dmu=mdu+udm=mcvdT+cvTdm1−cvTdm2,(12)dUdt=mcvdT+cvTdm1−cvTdm2dt,where m is the mass of natural gas in the tank (kg); cv is the constant volume specific heat of natural gas (cal/(mol·K)). Substituting (7) with (8), (9), (10), and (12), (13) is achieved.
(13)mcvdTdt=cvdm2dt−cpdm2dt−cvdm1dt−k2T+k1−zdm1dt+cp1T1dm1dt.In (13), m is the mass of the gas in the tank (kg); cp and cv meet the following relationship [24, 25] as follows.
(14)cp−cv=Rg,cp=A+BT+CT2+DT3,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
(15)mcvdTdt=k1−zdm11dt+cp11T11dm11dt+cp10T10dm10dt−k2+cpdm2dt+cvdm1dt−cvdm2dtT,where dm1/dt=dm10/dt+dm11/dt, dm10/dt is the mass of gas flowing into the control volume per unit time, and dm11/dt is the mass of vaporization gas in tank per unit time; T10 is the temperature of the gas entering the control volume; T11 is the critical temperature at the current pressure; cp10 and cp11 are the specific heat of constant pressure at temperature T10andT11.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.
(16)qs=PRTP′/8ξρπd2,ρ=PT0P0TM22.4.Equation (17) can be used to solve the current gas quality in the tank.
(17)m=m0+min−mout,where m0 denotes the gas quality in the tank at the previous moment; min denotes the gas quality entering the tank; mout is the quality of discharges gas. It should be noted that min and mout 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.


## 4. Model Validation

To validate the above model, data was collected from a LNG ship named “Dapeng.” Substitute the values of parameter in Table1 and Table 2 into the equations above, the liquid level, temperature, pressure, and liquid volume of tanks can be obtained as functions of time.Table 1
Nomenclature of the parameters appearing in the formula.
Parameter
Value
Parameter
Value
Rg
0.514 (kj/(kg·K))
z
121.7 (cal/g)
P0
101.325kPa
R
8.314 (J/(mol·K))
T0
273.5 K
P′
121.40kPa
M
16
g
9.8m/s2Table 2
Parameters of liquid cargo tanks for thermal load calculation.
Name
Value
Name
Value
Liquid capacityV0
41743.5m3
Surface area of cargo tankS1
7028m2
Tanker waterline external surface areaS2
3968m2
The surface area of the liquid cargo shipS3
3060m2
Surface area of fuel channelS4
93.6m2
Surface area of ballast tankS5
3536m2
Isolated space side areaS6
230m2
Main shielding thicknessδ1
0.7 mm
Main shielding thermal conductivityλ1
45 W/(m⋅K)
Main insulation thicknessδ2
230 mm
Thermal conductivity of the main insulation layerλ2
0.04 W/(m⋅K)
Secondary shielding thicknessδ3
0.7 mm
Secondary shielding thermal conductivityλ3
45 W/(m⋅K)
Secondary insulation thicknessδ4
300 mm
Secondary insulation thermal conductivityλ4
0.04 W/(m⋅K)
Resin thicknessδ5
10 mm
Thermal conductivity of resin layerλ5
0.302 W/(m⋅K)
Thickness of tank inner shellδ6
20 mm
Thermal conductivity of tank inner shellλ6
45 W/(m⋅K)
Fuel temperatureTf
40°C
Sea water heat exchange coefficientα1
2908 W/(m⋅K)
Fuel heat transfer coefficientα2
2368 W/(m2⋅K)All the constant parameters that appear in the formula are listed in Table1. 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 Figures2(a)–2(c). The curves of gas temperature, tank pressure, and liquid volume in the tank are achieved below.Figure 2
The comparison with real ship.
(a)
Gas temperature verification with real data.
(b)
Tank pressure verification with real data.
(c)
Liquid volume verification with real data.It is shown from Figures2(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 Figure2(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.


## 5. 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.(18)ρN2gP¨e=Fr−κu2Ṗe2,WherePe is the pressure difference of the primary and secondary insulating layers. u is nitrogen supply rate. ρN2 is the nitrogen density supplied, kg/m3; g is the gravitational acceleration; Fr 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 Pe0 = 2.0 mbar, u0 = 42 m3/h. A system transfer function of linearization is given below.
(19)Gs=ΔPesΔus=3.25ss−0.78e−0.25s.The time delaye−0.25s of the insulating space pressure difference system Gs is less than 0.5 s, so the 1st order Padé approximation of the time delay can be obtained as follows.
(20)G1s=e−0.25s=10.25s+1.Take the nominal model below.(21)G0s=3.25ss−0.78.SoGs is expressed as
(22)Gs=ΔPesΔus=G0sG1s=3.25ss−0.78.10.25s+1.Letss−0.780.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.Figure 3
Root locus of zero-pole for the unstable system.


### 5.2. 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–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. 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.Definition 2.
The process that zero-poles in the open right half plane of unstable process are substituted for their mirror-images and then a minimum-phase system is constructed is called mirror-injection process.In order to obtain better results, normalizing (23) is as follows.
(23)Gs=G0sG1s=3.25/0.78s1/0.78s−1.10.25s+1,After mirror-mapping, (24) is acquired as follows.
(24)Gs=G'0sG1s=3.25/0.78s1/0.78s+1.10.25s+1.From the bode plots of the two system models in Figure4, it can be seen that their spectrums are similar. The root locus of Gs 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.Figure 4
Bode plots of the two system models.Figure 5
Root locus of zero-pole of stable system.Figure 6
Step response of the stable system.


### 5.3. Controller Design

The controller is designed with closed-loop gaining shaping algorithm (CGSA). CGSA is a simplifiedH∞ 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 closed-loop 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 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/T1 (here T1 = 0.2).The controller is designed through shaping the frequency spectrum curve ofT which can be approximately regarded as the frequency spectrum curve of a second-order 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.
(25)1T1s+12=GK1+GK,K=1GT1sT1s+2.TheS&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.Figure 7
S&T singular value curves.To eliminate the steady-state error, a minor constantɛ (0.001) is added into the denominator of the G0s Then (26) is obtained.
(26)G'0s=3.25/0.78s1/0.78s+1+ε.Thus, substitute (26) into (24) and (25), according to closed-loop gain shaping algorithm, a linear proportional-integral-differential (PID) plus the first order filter controller is obtained.
(27)K1=1T1s+21/0.783.25/0.78T1s+13.25/0.78T1+ε3.25/0.78T1s.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 G1s, (28) is obtained.
(28)K=0.25s+1T1s+21/0.783.25/0.78T1s+13.25/0.78T1+ε3.25/0.78T1s.The block diagram of pressure difference control for insulating spaces is shown in Figure8. Gs in (19) is substituted into Figure 8.Figure 8
The block diagram of pressure control for the insulation spaces.


## 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.(18)ρN2gP¨e=Fr−κu2Ṗe2,WherePe is the pressure difference of the primary and secondary insulating layers. u is nitrogen supply rate. ρN2 is the nitrogen density supplied, kg/m3; g is the gravitational acceleration; Fr 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 Pe0 = 2.0 mbar, u0 = 42 m3/h. A system transfer function of linearization is given below.
(19)Gs=ΔPesΔus=3.25ss−0.78e−0.25s.The time delaye−0.25s of the insulating space pressure difference system Gs is less than 0.5 s, so the 1st order Padé approximation of the time delay can be obtained as follows.
(20)G1s=e−0.25s=10.25s+1.Take the nominal model below.(21)G0s=3.25ss−0.78.SoGs is expressed as
(22)Gs=ΔPesΔus=G0sG1s=3.25ss−0.78.10.25s+1.Letss−0.780.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.Figure 3
Root locus of zero-pole for the unstable system.


## 5.2. 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–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. 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.Definition 2.
The process that zero-poles in the open right half plane of unstable process are substituted for their mirror-images and then a minimum-phase system is constructed is called mirror-injection process.In order to obtain better results, normalizing (23) is as follows.
(23)Gs=G0sG1s=3.25/0.78s1/0.78s−1.10.25s+1,After mirror-mapping, (24) is acquired as follows.
(24)Gs=G'0sG1s=3.25/0.78s1/0.78s+1.10.25s+1.From the bode plots of the two system models in Figure4, it can be seen that their spectrums are similar. The root locus of Gs 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.Figure 4
Bode plots of the two system models.Figure 5
Root locus of zero-pole of stable system.Figure 6
Step response of the stable system.


## 5.3. Controller Design

The controller is designed with closed-loop gaining shaping algorithm (CGSA). CGSA is a simplifiedH∞ 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 closed-loop 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 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/T1 (here T1 = 0.2).The controller is designed through shaping the frequency spectrum curve ofT which can be approximately regarded as the frequency spectrum curve of a second-order 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.
(25)1T1s+12=GK1+GK,K=1GT1sT1s+2.TheS&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.Figure 7
S&T singular value curves.To eliminate the steady-state error, a minor constantɛ (0.001) is added into the denominator of the G0s Then (26) is obtained.
(26)G'0s=3.25/0.78s1/0.78s+1+ε.Thus, substitute (26) into (24) and (25), according to closed-loop gain shaping algorithm, a linear proportional-integral-differential (PID) plus the first order filter controller is obtained.
(27)K1=1T1s+21/0.783.25/0.78T1s+13.25/0.78T1+ε3.25/0.78T1s.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 G1s, (28) is obtained.
(28)K=0.25s+1T1s+21/0.783.25/0.78T1s+13.25/0.78T1+ε3.25/0.78T1s.The block diagram of pressure difference control for insulating spaces is shown in Figure8. Gs in (19) is substituted into Figure 8.Figure 8
The block diagram of pressure control for the insulation spaces.


## 6. 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: Kffs and Ks.
(29)Kffs=1.9820.20s2+0.77s+15.3420.25s3+0.99s2+1.64s,Ks=5.02s2+1.28s+1.980.12s2+1.282s.The block of 2-DOF structure control laws based on modified Smith predictor is shown in Figure9. Fs is a filter.Figure 9
The block diagram of pressure control for the insulation space.The comparative results are given in Figures10–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 (Tr) 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 (Ts) to reach its output final value is obviously shorter than the later one. The compared results are shown in Table 3.Figure 10
System response comparison without disturbance.Figure 11
System response comparison with step disturbance.Figure 12
System response comparison with ramp signal.Table 3
Parameters of liquid cargo tanks for thermal load calculation.
Control law
Mpt
Css
Overshoot (σp) (%)
Tr
Ts
Tp
The algorithm in this note
2.133
2.002
6.54
1.913
7.292
2.634
2-DOF
2.181
2
9.05
2.911
19.824
3.743
The rise time of response is defined asTr. The peak value of the step response is denoted by Mpt. Css is the value of stability sate. The time to reach peak value is Tp. Overshoot is defined by the equation: σp=Mpt−Css/Css×100. These values are calculated and given 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 Figure11. 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.


## 7. 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 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.



---

*Source: 7465391-2018-07-15.xml*