Example: Bilinear Estimator for a Plate Heat Exchanger

In this example, we will develop a bilinear estimator for a plate heat exchanger (PHEX) which is inspired by the work of Grunert et al. from 2015 (DOI). You can see the scheme of the whole heating appliance in the following picture.

The state we want to estimate is the outlet temperature of the plate heat exchanger, which is the Secondary Heat Exchanger in the picture. Designing an appropriate estimator would make the omission of the outlet temperature sensor possible. In order to achieve that, it is necessary to design an observer which can estimate all states accurately.

The estimator is based on a system model, which we need to set up first. In the interests of simplification, we can use a free body model which consists only of the PHEX, the sensors and the tubes which are located between the PHEX and the sensors.

Creation of the Components

The PHEX is already modeled as an ODESCA_Component. You can find it in the folder Examples->Components with the name OCLib_PlateHex. If you need help creating a custom component, please refer to our Getting Started Section.

Having modeled the PHEX, we are able to create an instance of it.

PHEX_comp = OCLib_PlateHex('myPHEX')

As you should see in the Command Window, the properties of the created component are still quite empty. The reason is a so called construction parameter which needs to be predisposed here. This parameter is the number of nodes in which the PHEX is divided. In every node, the internal states such as temperature and massflow are considered locally constant. The number of nodes affects the number of states, therefore, this parameter needs to be set before MATLAB can create any equations. We will divide the PHEX in 3 nodes.

PHEX_comp.setConstructionParam('Nodes',3)
PHEX_comp % looking at the properties again

Now, looking at the properties again, you should see that the equations have been calculated for 6 states (3 temperatures on the primary and 3 on the secondary side). For more information about the properties, please refer to our Class Documentation Section.

Take a look at all the properties of the component. The only missing fields are inside the property 'param'. In the next step we will complete the component by setting all the parameters to the values in the mentioned paper.

PHEX_comp.setParam('cHex', 500);
PHEX_comp.setParam('mHex', 1.154);
PHEX_comp.setParam('Volume1', 0.216*10^(-3));
PHEX_comp.setParam('Volume2', 0.24*10^(-3));
PHEX_comp.setParam('HexArea', 0.216);
PHEX_comp.setParam('RhoFluid', 998);
PHEX_comp.setParam('cFluid', 4182);

As you can read in the original paper, the thermal transmittance between primary and secondary side is given by a characteristic curve with three optimized parameters. This behavior is also already modeled as an ODESCA_Component in the Examples folder with the name OCLib_HeatTransferFit. We will create this component now and set the parameters. Later, we will take care of the connection between the components.

HTF_phex = OCLib_HeatTransferFit('myHeatTransferFit');
HTF_phex.setParam('c1', 6529);
HTF_phex.setParam('c2', 4029);
HTF_phex.setParam('c3', 3.29);

The model of the system also includes two temperature sensors (for outlet and return temperature) and a pipe between the PHEX and the return temperature sensor. You can find those components in the Examples folder, as well: OCLib_TSensor and OCLib_Pipe. The creation as well as the parameterization is done in the following.

OTSensor = OCLib_TSensor('myOTSensor');
OTSensor.setParam('TimeConst',1.8);
OTSensor.setParam('Gain',1);

RTSensor = OCLib_TSensor('myRTSensor');

PHEXtoRTSensor_Pipe = OCLib_Pipe('myPipe');
PHEXtoRTSensor_Pipe.setConstructionParam('Nodes',1);
PHEXtoRTSensor_Pipe.setParam('VPipe',0.55*1e-3);
% Note: There is no information given about the parameters cPipe and mPipe but only about their product C_t. 
%       Since these two parameters only appear as C_t in the equations, we can simply set one of these parameters to 1 and the other to the value of C_t.
PHEXtoRTSensor_Pipe.setParam('cPipe',1);
PHEXtoRTSensor_Pipe.setParam('mPipe',150);

As you might have noticed, we didn't set all the parameters of the new components yet. The reason is that those parameters already appeared in the some of the other components. We will use an ODESCA function to equalize those parameters later to avoid mistakes and reduce computation time.

Creation of the System

After we created and parameterized all the components, we can now create an ODESCA_System for our further work. At first, we initialize the system with the PHEX and add all the components to it.

PHEX_sys = ODESCA_System('myPHEXSystem',PHEX_comp);
PHEX_sys.addComponent(HTF_phex);
PHEX_sys.addComponent(OTSensor);
PHEX_sys.addComponent(RTSensor);
PHEX_sys.addComponent(PHEXtoRTSensor_Pipe);

As already mentioned, we will equalize the parameters now.

PHEX_sys.param % take a look at the parameter list before equalizing
PHEX_sys.equalizeParam('myPHEX_cFluid',{'myPipe_cFluid'})
PHEX_sys.equalizeParam('myPHEX_RhoFluid',{'myPipe_RhoFluid'})
PHEX_sys.equalizeParam('myOTSensor_Gain',{'myRTSensor_Gain'})
PHEX_sys.equalizeParam('myOTSensor_TimeConst',{'myRTSensor_TimeConst'})
PHEX_sys.param % take a look at the parameter list after equalizing

If you take a look at the parameter list in the PHEX system, you should see that the parameters in the second arguments of the equalizing function have vanished. Those parameters have been replaced inside of the equations with the parameters in the first arguments of the equalizing function, respectively.

At this point, the system is completely parameterized, however, the components are not connected in any way but exist parallel and independently of one another inside the system. We will connect them by substituting in- and outputs of the respective components.

PHEX_sys.connectInput('myPipe_TempIn','myPHEX_Temperature1Out')
PHEX_sys.connectInput('myPipe_mDotIn','myPHEX_Massflow1Out')
PHEX_sys.connectInput('myOTSensor_TempIn','myPHEX_Temperature2Out')
PHEX_sys.connectInput('myRTSensor_TempIn','myPipe_TempOut')
PHEX_sys.connectInput('myPHEX_Massflow1In','myHeatTransferFit_Massflow1Out')
PHEX_sys.connectInput('myPHEX_Massflow2In','myHeatTransferFit_Massflow2Out')
PHEX_sys.connectInput('myPHEX_kHex','myHeatTransferFit_k_Hex')

The four remaining inputs of the system are the two massflows and the two temperatures of the primary and secondary side of the PHEX. However, we can still see all of the outputs of the components, though, in the real system, we are only able to measure the return and the outlet temperature using the two sensors. To get a model of the real system, we can remove all unwanted outputs as follows.

PHEX_sys.removeOutput('myPHEX_Temperature1Out')
PHEX_sys.removeOutput('myPHEX_Temperature2Out')
PHEX_sys.removeOutput('myPHEX_Massflow1Out')
PHEX_sys.removeOutput('myPHEX_Massflow2Out')
PHEX_sys.removeOutput('myHeatTransferFit_k_Hex')
PHEX_sys.removeOutput('myHeatTransferFit_Massflow1Out')
PHEX_sys.removeOutput('myHeatTransferFit_Massflow2Out')
PHEX_sys.removeOutput('myPipe_TempOut')
PHEX_sys.removeOutput('myPipe_mDotOut')
PHEX_sys % take a look at the system now

If you take a look at the system now, you should see 9 states, 4 inputs, 2 outputs and 15 parameters.

Bilinearization of the System

For the creation of a bilinear observer, we need a bilinear approximation of the nonlinear system. The approximation is only accurate close to the center of the series. We want to bilinearize the system at the most common steady state of the system, which we assume as the following input u0 along with the resulting state x0.

u0 = [60; 10; 0.25; 0.15];

The resulting state x0 will be calculated in the following using the system equations. We also want to calculate y0 for setting up the bilinear observer later.

[PHEX_sys_f,PHEX_sys_g] = PHEX_sys.calculateNumericEquations();
steadystate = vpasolve(subs(PHEX_sys_f,PHEX_sys.u,u0),PHEX_sys.x);
x0 = double([steadystate.x1; steadystate.x2; steadystate.x3; steadystate.x4; steadystate.x5; steadystate.x6; steadystate.x7; steadystate.x8; steadystate.x9]);
y0 = double(subs(PHEX_sys_g,PHEX_sys.x,x0));

Now, we can create an ODESCA_SteadyState which is the base for the creation of a system approximation such as the bilinearization. If (x0,u0) is no valid steady state, the creation of an ODESCA_SteadyState will display a warning. The function bilinearize finally creates an instance of the class ODESCA_Bilinear.

ss = PHEX_sys.createSteadyState(x0,u0,'mySteadyState');
PHEX_sys_bilin = ss.bilinearize();

All the different classes (ODESCA_Component, ODESCA_System, ODESCA_SteadyState, ODESCA_Bilinear) are still connected to each other so that it is always distinct, which system belongs to which approximation et cetera. For example

PHEX_sys_bilin.steadyState.system.components

will give you the names of all components linked to the system.

Creation of the Estimator

For the common LQE (linear quadratic estimation) approach, only the linear parts of the system are considered to calculate the feedback matrix L. The weights Q and R are chosen as unity matrices with decreased values of R to achieve higher observer gains. The lqr function of matlab can be applied to the dual system for calculating the observer gain.

Q = eye(9,9);
R = eye(2,2)*0.01;
L = (lqr(PHEX_sys_bilin.A',PHEX_sys_bilin.C',Q,R))';

To finally see the results, you can create a Simulink model with the nonlinear system and the bilinear observer. The creation of the nonlinear system is done by the following command.

PHEX_sys.createNonlinearSimulinkModel

However, we already created the whole test system for you in Examples->Systems->Estimation_PHEX. Note, that you need to load the Example_PHEX.mat or to execute the previous commands in order to run the Simulink model.

As you can see, for different initial states the simple bilinear observer estimates all states accurately when the temperatures perform steps. Even when the sensors signals are disturbed, all temperatures inside the PHEX are nearly stationary accurate. The stationary accuracy of the outputs could be achieved by using a reduced observer. The only problem of this observer can be detected when the massflows perform steps. This issue is tackled in the mentioned paper.