Estimate Continuous-Time Grey-Box Model for Heat Diffusion

This example shows how to estimate the heat conductivity and the heat-transfer coefficient of a continuous-time grey-box model for a heated-rod system.

This system consists of a well-insulated metal rod of length L and a heat-diffusion coefficient $$\kappa$$. The input to the system is the heating power u(t) and the measured output y(t) is the temperature at the other end.

Under ideal conditions, this system is described by the heat-diffusion equation—which is a partial differential equation in space and time.

$$\frac{\partial x(t,\xi )}{\partial t}=\kappa \frac{{\partial }^{2}x(t,\xi )}{\partial {\xi }^2}$$

To get a continuous-time state-space model, you can represent the second-derivative using the following difference approximation:

$$\begin{array}{l}\frac{{\partial }^{2}x(t,\xi )}{\partial {\xi }^2}=\frac{x\left(t,\xi +\Delta L\right)-2x(t,\xi )+x\left(t,\xi -\Delta L\right)}{{\left(\Delta L\right)}^2}\\ \text{where }\xi =k\cdot \Delta L\end{array}$$

This transformation produces a state-space model of order $$n=\frac{L}{\Delta L}$$, where the state variables $$x(t,k\cdot \Delta L)$$ are lumped representations for $$x(t,\xi )$$ for the following range of values:

$$k\cdot \Delta L\le \xi <\left(k+1\right)\Delta L$$

The dimension of x depends on the spatial grid size $$\Delta L$$ in the approximation.

The heat-diffusion equation is mapped to the following continuous-time state-space model structure to identify the state-space matrices:

$$\begin{array}{l}\dot{x}(t)=Fx(t)+Gu(t)+\tilde{K}w(t)\\ y(t)=Hx(t)+Du(t)+w(t)\\ x(0)=x0\end{array}$$

The state-space matrices are parameterized by the heat diffusion coefficient κ and the heat transfer coefficient at the far end of the rod h_tf. The expressions also depend upon the grid size, Ngrid, and the length of the rod L. The initial conditions x0 are a function of the initial room temperature, treated as a known quantity in this example.

1. Create a MATLAB^® file.

   The following code describes the state-space equation for this model. The parameters are κ and h_tf while the auxiliary variables are Ngrid, L and initial room temperature temp. The grid size is supplied as an auxiliary variable so that the ODE function can be easily adapted for various grid sizes.

   function [A,B,C,D,K,x0] = heatd(kappa,htf,T,Ngrid,L,temp)
   % ODE file parameterizing the heat diffusion model
   
   % kappa (first parameter) - heat diffusion coefficient
   % htf (second parameter) - heat transfer coefficient 
   %                          at the far end of rod
   
   % Auxiliary variables for computing state-space matrices:
   % Ngrid: Number of points in the space-discretization
   % L: Length of the rod
   % temp: Initial room temperature (uniform)
   
   % Compute space interval
   deltaL = L/Ngrid;
   
   % A matrix
   A = zeros(Ngrid,Ngrid);
   for kk = 2:Ngrid-1
     A(kk,kk-1) = 1;
     A(kk,kk) = -2;
     A(kk,kk+1) = 1;
   end
   
   % Boundary condition on insulated end
   A(1,1) = -1; A(1,2) = 1;
   A(Ngrid,Ngrid-1) = 1;
   A(Ngrid,Ngrid) = -1;
   A = A*kappa/deltaL/deltaL;
   
   % B matrix
   B = zeros(Ngrid,1);
   B(Ngrid,1) = htf/deltaL;
   
   % C matrix
   C = zeros(1,Ngrid);
   C(1,1) = 1;
   
   % D matrix (fixed to zero)
   D = 0;
   
   % K matrix: fixed to zero
   K = zeros(Ngrid,1);
   
   % Initial states: fixed to room temperature
   x0 = temp*ones(Ngrid,1);
   

2. Use the following syntax to define an idgrey model object based on the heatd code file:

   m = idgrey('heatd',{0.27 1},'c',{10,1,22});

   This command specifies the auxiliary parameters as inputs to the function, include the model order (grid size) 10, the rod length of 1 meter, and an initial temperature of 22 degrees Celsius. The command also specifies the initial values for heat conductivity as 0.27, and for the heat transfer coefficient as 1.

3. For given data, you can use greyest to estimate the grey-box parameter values:

   me = greyest(data,m)
   

The following command shows how you can specify to estimate a new model with different auxiliary variables:

m.Structure.ExtraArguments = {20,1,22};
me = greyest(data,m);

This syntax uses the ExtraArguments model structure attribute to specify a finer grid using a larger value for Ngrid. For more information about linear grey-box model properties, see the idgrey reference page.