Abstract
In this paper, we present a system identification (SI) procedure that enables building linear time-dependent fractional-order differential equation (FDE) models able to accurately describe time-dependent behavior of complex systems. The parameters in the models are the order of the equation, the coefficients in it, and, when necessary, the initial conditions. The Caputo definition of the fractional derivative, and the Mittag-Leffler function, is used to obtain the corresponding solutions. Since the set of parameters for the model and its initial conditions are nonunique, and there are small but significant differences in the predictions from the possible models thus obtained, the SI operation is carried out via global regression of an error-cost function by a simulated annealing optimization algorithm. The SI approach is assessed by considering previously published experimental data from a shell-and-tube heat exchanger and a recently constructed multiroom building test bed. The results show that the proposed model is reliable within the interpolation domain but cannot be used with confidence for predictions outside this region. However, the proposed system identification methodology is robust and can be used to derive accurate and compact models from experimental data. In addition, given a functional form of a fractional-order differential equation model, as new data become available, the SI technique can be used to expand the region of reliability of the resulting model.
1. Introduction
Complex systems are common in mechanical engineering applications, for example, automobiles, washing machines, and thermal power plants. As illustrated schematically in Figure 1, each of these systems is composed of many interconnected subsystems. In principle each subsystem can be studied in isolation and its input and output connected to others. However, for purposes of control we may think of the overall system as being of the single-input single-output (SISO) type, i.e., the interconnections between subsystems are all internal, and there is only a single input and a single output of the overall system, where is time. The control of such systems may either be model based or not [1]. Control procedures that are not based on a mathematical model include techniques such as PID [2], fuzzy logic [3], or other procedures that work on the error signal. For model based control, on the other hand, one should have an approximate mathematical model of the time-dependent behavior of the system to be able to predict what is going to happen and then control it (approximate because otherwise there would be no need for control). For simplicity we will only consider linear systems. There are two types of inputs that are of interest in complex SISO systems: one is relaxation and the other is periodic. In the former, which is what we will be concerned with here, a step input is applied to the system and the output goes in some dynamic fashion from one constant value to another.
One approach to modeling is from the ground up, i.e., to use first principles to model each component and their interactions to create a so-called white box model. If we think of each subsystem as being governed by a single ordinary differential equation (ODE), then the overall system is governed by a large set of coupled ODEs. If, however, there is a governing partial differential equation (PDE), then that may be considered to be equivalent to an infinite set of ODEs. Modeling is possible only when there is physical and mathematical understanding of the behavior of each subsystem. In any case, for control of the overall system one ends up working with a large set of ODEs that must be solved in real time. This may be computationally undesirable for control purposes.
In many cases it is advantageous to use existing or continually acquired experimental data to make predictions [4]. This is a black-box approach and its main advantage is that the input-output response of the overall system can be directly obtained, but a disadvantage is that there is no physical understanding of the dynamic behavior compared to the ground-up approach. For example, one can have a procedure based on artificial neural networks which can either be trained beforehand or continuously trained [5]. Alternatively one can fit an analytical model to the data, and this is system identification (SI). For control purposes the data-based model should be good enough to approximate reality but also be solvable in real time without intensive computation. The SI process leads to a mathematical model that can make predictions that differ little from reality in the region of interest and is compact enough to be easily solved and used for control purposes. The most common SI procedure is to propose a mathematical model of the relaxation dynamics of the system, usually in the form of differential equations, with certain free constants that are adjusted to best fit measurements. In this manifestation SI is really determination of parameters in the proposed model.
Two new aspects will be explored here in relation to SI.
(a) Fractional-Order Model. Usually the proposed mathematical models involve integer-order derivatives, but this lacks the generality provided by derivatives of fractional order [6–8]. In addition, the required mathematical description of the system is directly linked to its complexity, and when the latter increases so does the former. Though any system representation can be achieved using integer-order conservation equations, their solution, in the context of system control, requires extremely large CPU times, and it becomes necessary to develop accurate and compact models that can be used to determine the corresponding behavior in a timely and reliable manner. Since it has been recently demonstrated that fractional-based models are able to describe systems that are complex [9], and a physical understanding of the subsystem processes is not of concern for control of the overall system, the use of suitably defined fractional derivatives is appropriate [10, 11]. There are many different definitions of fractional derivatives available, though there are some common features between them [12]. Initial conditions for time-dependent systems is an important issue for fractional-order systems: unlike integer-order derivatives where the number of initial conditions corresponds to the order of the differential equation, many of the definitions of fractional derivatives make the equation infinite dimensional so that an infinite number of initial conditions are needed. Here we use the Caputo definition of the fractional derivative since it enables the prescription and application of initial conditions of integer order typically encountered in physical systems, and analytical solutions are often possible.
(b) Error Minimization by Global Regression. The mathematical model with the best approximation to experimental data is obtained through a process of error minimization. On the one hand this can be done through an algorithm that searches locally in the neighborhood of an initial guess, a procedure that is effective if there is a single minimum. On the other hand there are procedures that search globally within the permitted range of parameters for one or more minima. One of these techniques is that of genetic algorithms [4] and another is simulated annealing, which is what we will use here [13, 14].
In the present investigation, we propose a methodology of analysis to derive accurate models – based on fractional-order differential equations – that approximate relaxation processes of complex systems. To this end, the paper presents first a brief description of the background information about fractional calculus. Next, the system identification technique, based on global regression of an error-cost function, is introduced in detail. Application to previously published experimental data derived from (1) a heat exchanger and (2) a recently constructed building test facility is later carried out with special emphasis on the accuracy of the model within the region of interest. A vital conclusion is the nonuniqueness in the resulting parameter set that defines both the fractional-order model and the corresponding initial conditions. Although the proposed model is accurate within the interpolation region, it is unreliable outside it.
2. Background on Fractional Calculus
Fractional (noninteger) calculus, which can be thought of as a generalization of the well-known integer calculus, has a long history dating back to Leibniz [9]. Nevertheless, it is not as celebrated by the scientific community as its integer counterpart, mainly due to insufficient enthusiasm about its potential usefulness in applications. Although significant progress in its theoretical basis has been achieved, primarily due to contributions from 18th- and 19th-century mathematicians (e.g., Laplace, Fourier, Abel, and Liouville), it is only in the last two decades that the subject has significantly broadened its applications in physics and engineering. This is reflected by the increased number of publications (i.e., books [9, 15–17], review monograms [10], and scientific articles) devoted to the topic in a variety of fields, like anomalous diffusion [18], semi-infinite tree networks of mechanical, electrical, and hydrodynamic equipment [19, 20], and viscoelastic systems [21], among others. In the context of modeling and control of thermal devices and phenomena, studies, such as those of Aoki et al. [22], Pineda et al. [23], Gabano and Poinot [24], and Caponetto et al. [25], demonstrate that fractional derivatives provide good approximations for describing the dynamic behavior of heat transfer processes.
In principle, the field of fractional calculus entails the generalization of the concept of a derivative beyond that of the common integer order. It is concerned with the meaning of the “in-between” derivatives. For instance, the 5/4-derivative of a function ,can be thought of as a derivative in-between the first and second derivatives. Though it seems straightforward to think of a derivative whose order is between two consecutive integers, the definition of a fractional derivative is not unique and various definitions are possible as long as it satisfies a required set of mathematical rules [12]. Two of the most common definitions are those of Riemann-Liouville [16] and Caputo [26], the latter being the type of derivative used in the present study.
Let us start with the Cauchy formula for repeated integrals of integer order , for the function , which is defined aswhere is the integral operator. By defining the Gamma function asthe term in (2) can be replaced with , to obtainwhere is an arbitrary positive real number and is a dummy variable of integration. Note that ; thus, for a value of such that which establishes the Riemann-Liouville fractional derivative of order , asThe Caputo definition of the fractional derivative of order , on the other hand, is defined as [26]so that is explicitly given asIn (6) and (8), the subscripts RL and C, respectively, denote Riemann-Liouville and Caputo fractional derivatives.
Regardless of the definition used, the Laplace transform is a powerful technique to solve equations involving fractional derivatives. For instance, the Riemann-Liouville definition in the Laplace domain iswhere it can be seen that the transformation requires the following initial conditions: for all . It is important to note that, since does not need to be an integer number, the terms are all of fractional order. This poses a problem for applications in physical systems, where the initial conditions are typically given in terms of integer orders. On the other hand, the Laplace transform for the Caputo definition iswhere now the initial conditions necessary to evaluate (10) are prescribed as integer order, , rather than fractional-order derivatives. Therefore, in the subsequent sections of the paper, we use the Caputo definition of a fractional-order derivative.
3. System Identification
The first step in the SI process is to propose the model and then to determine the parameters therein that will approximate – as closely as possible – its output to the target data. For our purposes, the proposed model is the linear fractional differential equation, given aswhere and are constants, is the order (noninteger, , with ) of the equation, is its input, is the corresponding output, and is a time-like independent variable. Since the interest is on the conditions of operation that would take the system from one state to another, i.e., relaxation processes, taking to be the Heaviside function, the Laplace transform of (11) iswhich can then be solved if the initial conditions for , and its corresponding derivatives, are provided. The inverse Laplace transform of the above equation leads to the following analytical form for where is the typical Mittag-Leffler function [17]. Note that the output from (13) is completely defined if the constants are known. This implies that one must provide derivatives at (); i.e., for , , and for , , and are necessary to determine the solution.
To find the unknowns , , , , for , in (11) and (13), from known data (e.g., values ), a least-squares method is used to minimize the difference between the sample- and predicted-values of . This is equivalent to minimizing an Euclidean norm, i.e., the variance of the error, given bywhere is the number of data points, usually from experimental measurements; , for , are values predicted by the model and , for , refer to the experiments. It is to be noted that , with being the vector of unknown parameters, is a smooth manifold in a -dimensional space. The goal is now to search for the values of the parameter set , such that is a minimum. This process can be carried out either using local – such as gradient-based – methods or global optimization techniques like evolutionary or deterministic algorithms, each with advantages and drawbacks [27–31]. However, as we will show in a later section, under certain cases (14) has multiple local minima, stemming from its nonlinear nature with respect to the arguments, thus leading us to seek the global minimum which will provide the optimal values of the parameter set required in the model (11), and the corresponding initial conditions. This optimization procedure is carried out here by the simulated annealing (SA) technique.
Description. The SA technique is inspired by the molecular calculation of the cooling of a physical system in which random agitation is used to avoid entrapment in local extrema [29]. A starting point in the space of unknowns, , is randomly selected, and a cycle of random moves along each coordinate direction is then performed. The new point is accepted if it gives a better value of . If it is worse, it is accepted only with a certain probability, , where is the change in value of and is a dynamic parameter that is analogous to the temperature of a system being cooled [32]. The process is repeated with decreasing , and step size until convergence within a certain tolerance is reached. The procedure described above [33] has been successfully implemented and used in the context of heat transfer correlations by Pacheco-Vega et al. [34], and it is the one followed here.
It is to be noted that other global optimization algorithms can also be used in the search for the set of parameters . These include the popular genetic algorithm (GA) [28] and the interval method (IM) [35], among several alternatives, which have been successfully used in a variety of applications. Though each technique has advantages over the others (particularly in terms of CPU time necessary to find a solution), the GA will only find the region – with a high degree of certainty – where the global optimum is located. From this perspective, only IMs mathematically guarantee that the global optimum has been found, while the SA is only probabilistically – but not deterministically – guaranteed to find such an optimum solution.
SI using the fractional model in (11), herein referred to as FOSI, is applied to experimental data from two complex thermal systems: (1) shell-and-tube heat exchanger data analyzed by Mayes [36] and (2) a multiroom building test facility. Preliminary analysis on the applicability of the FOSI methodology to a set of analytical problems and to heat exchanger data has been recently reported in Li et al. [37]. In what follows, the search domain for , in (11), is restricted to with .
4. Shell-and-Tube Heat Exchanger Data
4.1. Description
Details about the experimental setup in the thermal systems laboratory at the University of Notre Dame, and the corresponding data obtained with it, have been reported in Mayes [36]; here we provide only a brief description of the problem at hand. The heat exchanger, a schematic of which is shown in Figure 2, corresponds to the common shell-and-tube configuration. The figure illustrates flow directions of the cold and hot fluids, and , respectively, along with their corresponding terminal inlet, and , and outlet, and , bulk temperatures. The reported tests [36] were carried out under time-dependent conditions until thermal equilibrium in the device was achieved. Subsequently, while maintaining constant values of , , and , a step change in , i.e., , was applied. A data acquisition system (DAQ) and a personal computer (PC) were used to record and store the time-dependent values of the inlet and outlet temperatures of the two fluids.
For purposes of analysis, we take the approach of [36] and work with normalized variables. The dimensionless outlet temperature of the hot fluid and time variable are defined aswhere and are, respectively, the values of hot-fluid outlet temperatures at the initial and steady states; is the rise time, defined as the time required for to reach of . Thus, , with . The experiments showing the temporal evolution of are illustrated graphically in Figure 3(a), where it is clear that some type of relaxation process from one state (with having a small positive value) to another (where is close to 0.8) took place. Since the goal here is to find the parameters in the fractional model of the form given by (11) which fit the experimental data the best, for clarity we define and , so that the fractional-order differential equation becomeswhere and are constants, is the order in the search domain , , and is the input resulting from the step change in the mass flow rate of the cold fluid . Once the model has been defined (i.e., Eq. (16)), the FOSI technique is applied to the normalized data.
(a) Experimental data and model approximation
(b) Section of surface . A = local minimum; B = global minimum
In this case, the variance of the error is given aswhere and correspond to predictions and experiments, respectively, and is the number of data sets. The set of unknown parameters that will be considered in the forthcoming sections will be either or . Notice that, when seeking the optimum values in the parameter set, either one or two initial conditions and are included; i.e., we have relaxed the restriction regarding the number of unknown parameters and now included either one or two initial conditions in the search. The inclusion of the initial conditions in the parameter search brings the advantage of finding the most appropriate values for them directly from the operating conditions of the system and, therefore, provides a better approximation from the model thus obtained. However, as it will be shown in the next section, this generates multiplicity of solutions of the parameter set, some of which may not be physically feasible.
4.2. Initial Conditions and Nonuniqueness
Initial conditions provide the set of values of the system at the point of departure towards a new set of states. Mathematically, these conditions can be easily set up so that the model can be solved. Experimentally, however, this may not be the case since it is difficult to measure the value of the function (e.g., temperature) at a specific time (e.g., ), and even more difficult to establish its rate of change, i.e., the derivative. Sensitivity of the sensors, and even where they are placed, plays a crucial role and increases the uncertainty of the corresponding measurement.
As an example of this situation, let us focus on the experimental measurements of the shell-and-tube heat exchanger described before. The system was analyzed in [36] for dynamic conditions of the outlet temperatures after some thermal equilibrium in the device was achieved. Figure 3(a) shows the experimental data (digitized from [36]) corresponding to the evolution of the dimensionless temperature . From the figure it can be observed that the general trend in the values of the experimental data near () is to decrease in magnitude before showing a sustained relaxation-type increase (with some variation about the mean values) to about at . This seems to indicate that, unless sudden changes in temperature occurred in the laboratory, the heat-exchange device might not have actually reached the expected thermal equilibrium, and it raises the question about the actual level of uncertainty that may be present while attempting to accurately establish the initial conditions for complex systems since the values and , in this case, are not feasible.
An example is the case of a sensor measuring temperature that either is not appropriately calibrated or that has a slow response time-constant, then the initial states for the temperature or its rate of change will not be correctly established. In other cases, the system may be so complex that establishing actual initial conditions is virtually impossible. A car comprising a large number of subsystems may have initial states for each of those components that are different than those of the overall system, and deciding which initial state is the most appropriate (and thus, the one that should be used) may be extremely difficult. A similar situation occurs in a fluid, where initial conditions must be specified at all locations. From this perspective, the inclusion of the initial conditions in the unknown parameter set may be advantageous since the optimal values of these states can be directly computed from the data via SI. On the other hand, mathematically, by increasing the number of adjustable parameters it may be possible to obtain a better approximation from the model. This is what we have done within the context of the two thermal systems considered here.
In what follows, we apply the global-regression-based FOSI to the heat exchanger experimental data to obtain the parameters in (16) and – when necessary – the corresponding initial conditions that minimize the variance of the error [Eq. (17)]. The search for the parameters in is constrained to the following ranges: , , and . The initial condition is kept as , as established in [36]. For purposes of comparison, a gradient-based local optimization algorithm (LOA) is also applied to seek the best-possible parameters in . In applying the LOA technique, the search is constrained only for the order , while the domain for all other parameters remains unconstrained.
The results from the two fractional-order models are shown in Figure 3(a) and Table 1. From the figure it can be seen that the model obtained by global regression (model B) provides a much better approximation to the experimental data than that from the local-regression procedure (model A). Although both models are somewhat closer in the middle regions, i.e., , their predictions deviate significantly at the two ends of the relaxation process, with the SA-based model (also referred to as a four-parameter model in Figure 3(a) since the vector comprised four parameters) being much closer to the original data. On the other hand, from the table it can be seen that the two fractional models have very different values of the parameters , , and , and particularly the fact that only appears in the fractional-order model obtained via global regression (model B) since its value of lies between and . This outcome is indicative of the nonuniqueness of the resulting fractional differential equation model, along with the initial condition considered in the parameter set .
Figure 3(b) shows a section of the hypersurface that passes through the two minima and . For clarity, a location coordinate is defined such that , where and , corresponding to the two minima. The value of the cost function from the LOA-based model (model A) is 3 times higher than that of the SA model (model B in Table 1). The above result leads us to ask which model is better, and a possible answer may be the model equation that provides the minimum error, i.e., the global minimum. This is accomplished in the present study through the implementation of SA in the SI methodology.
By using the heat exchanger data described before, and relaxing the restriction on the number of parameters in to now include the additional initial condition , , i.e., a five-parameter fractional-order model, then it is clear as illustrated in Figure 4 that not only the overall solution but also the value of is much closer to the target data than that of the four-parameter model. It is to be noted that, although not shown explicitly, there is also a multiplicity of solutions in , the global optimum being found by the SA algorithm. The multiplicity in the parameter set seems to be due to the fact that the initial condition(s) is(are) included in the search rather than the model itself. Finally, Table 2 shows, quantitatively, the improved accuracy of the five-parameter fractional model with a value of , which is an improvement of 4.5% with respect to that of the model based on four parameters.
4.3. Reliability of the Model
The use of fractional-order approximation models stems from the fact that as the complexity of a system increases, so does its required mathematical representation. Though this can be achieved starting from conservation equations, for purposes of control it is important that the plant model is compact, efficient, and sufficiently accurate. For instance, a model based on the Navier-Stokes equations may be more accurate than an approximate model, but would be computationally intensive for control (for design and prediction, on the other hand, the mathematical model needs to be accurate). An approximate model, however, may also serve the purpose of bringing the output close to the desired state from where some other form of nonmodel control (e.g., PID) can take over. From this perspective, the fractional-order model given in (11) and (16) can be regarded simply as an approximation to the governing equations (PDEs), and certainly it is not a model based on first principles (i.e., derived from the laws of physics for a specific problem), but is developed based on a curve-fitting process from a set of experimental data in a specific domain. Thus, it follows that the resulting mathematical model can only be expected to perform well in the domain of the data from which it was derived (i.e., interpolation). In other words, the model may not work – and it does not have to work – beyond the domain in which the data apply (i.e., within the cloud of training data). Making predictions outside the domain of the data from which the model was derived would correspond to extrapolation, and therefore unreliable.
This data-driven model interpolation/extrapolation issue is particularly problematic for the case of complex systems, an example of which is heat exchangers, and the phenomena associated with them. In these thermal systems, complexity arises from geometrical configurations, the large number of parameters involved in their operation, and the nonlinear nature of the system. In the latter case, for example, there is nonlinearity due to variation of properties (e.g., density, viscosity, and thermal conductivity), with temperature. Nonlinearity in the phenomena then generates the possibility of bifurcations, including instability and transition to turbulence. Finally, the conservation equations themselves, which provide the most accurate description of the system, are nonlinear (e.g., in the advective terms). It is apparent that nonlinearities associated with the system may be a main reason for constraining predictions from data-driven models to within the interpolation region only.
The issue of interpolation vs. extrapolation is prevalent in models derived from experimental data. For instance, techniques such as artificial neural networks and correlation equations, among others, are unreliable for extrapolation; yet, they are extremely useful for predictions within the region of training data of the system behavior at operating points different from those used to derive the model. This is illustrated in Figure 5 in the context of the current investigation, which shows the results from the fractional-order model given in (16) and an ordinary second-order linear oscillator given bywith , zero initial conditions, and fitted constants , , and obtained from the same data. The figure shows the evolution of for the ranges and . For completeness, the results from a first-order model are also included. From Figure 5(a), the results show that, in the range , both models predict the behavior of the system very well, with the fractional-order model being more accurate ( and , respectively), and both show (see Figure 5(b)) some discrepancy with respect to the expected behavior in the range , with the second-order model being qualitatively more accurate than the fractional-order. However, in this range the actual values of the expected behavior from both models are incorrect, indicating a lack of reliability.
(a) Predictions in interpolated region
(b) Predictions in extrapolated region
The issue of reliability of mathematical models is an active research topic [38], and several studies have attempted to establish the region of reliability by defining it either as a hyperbox or by using the convex hull of the data. Others have concluded that the number and distribution of the experimental data are important factors in the reliability of empirical models [39, 40]. In fact, in the context of modeling with artificial neural networks, Pacheco-Vega et al. [41] discovered that even inside the cloud of data – and its corresponding convex hull – extrapolation could exist if the number and distribution of data points in a high-dimensional parameter space were limited in number and location, and proposed a methodology to establish the upper limits in the prediction error from such models.
4.4. Range of Applicability
In practice, experimental uncertainty may be present due to nonideal sensor measurement, thus it may be necessary to establish the applicability of the resulting model for the system under analysis. In the present case, we use the experimental data, along with a definition based on the root-mean-square (RMS) error, to establish the range of applicability for the five-parameter fractional-order model, i.e., (16), with fitted and .
The first step is to define the RMS error asThis error provides the baseline value for the range of applicability from the original predictions. The next step is to shift the data up and down by a value , asThe global-regression-based SI is then applied to the shifted data and the optimal parameters obtained. The results are given in Table 2 and Figure 6. Table 2 presents a quantitative account of the results showing the range of applicability of the fractional differential equation model, where the order of the equation remains essentially unaltered, while the other - and -based shifted parameters provide the values corresponding to the upper and lower bounds. Figure 6 supplies a graphical account of the same results, where it is clear that the model is able to accurately predict the behavior of the system.
5. Multiroom Building Data
The second test deals with a more complex multiroom building system. A subscaled two-floor building facility, designed specifically to be a scaled model of an actual building for energy-related studies, is located in the thermofluids laboratory at the California State University, Los Angeles, and is used to conduct experiments for thermal analysis and control [42]. A schematic of the experimental facility, with overall dimensions of 1.2 m 0.92 m 1.1 m, and its photographic depiction are shown in Figure 7. The test bed has a two-floor configuration with four rooms in each floor. Its structure is made of wood with interior walls covered with drywall and insulation (additional details are in [42]). Average temperatures of air in the rooms are measured by type-K thermocouples, while incandescent light bulbs serve as internal heat sources. Supply and return vents are installed in each room, and an external cooling unit is used to provide cold air through a set of ducts, each connected to the corresponding supply vent. Flow rates of air delivered into each room can be modified by a set of dampers (valves) through Arduino microcontrollers. Time-dependent information of the room temperatures, , is collected in two DAQ boards and stored in a PC for further analysis. LabVIEW serves to interface the controller and the experimental system. The facility has been recently studied by Baltazar et al. [43] in the context of thermal control with fuzzy logic.
(a) Schematic view
(b) Photograph
The average air temperature measurements for the eight rooms of the building are illustrated in Figure 8 which shows that each room is initially heated to a temperature above 27°C and then cooled down for 30 min [43]. From the figure it can be seen that all temperature curves depict three features: (1) the cooling rate in each room is different due to uneven delivery of airflow by the dampers; (2) there is a nearly exponential decay in the average temperatures, which is expected due to the convection process that takes place in each room; and (3) the decay in the temperatures is altered by two peaks, at = 17 min and = 25 min, resulting from the influence that the laboratory temperature has on the operation of the cooling unit.
For the analysis that follows, room 1 is selected because of its smooth transition from a maximum temperature, , to a temperature considered as minimum, , just before the peak occurs at . As before, it is useful to normalize its mean temperature and the time variable , as where is the dimensionless average temperature for air, the dimensionless time, C, C, and min. With the aforementioned definitions, the linear fractional-order model is again defined aswhere results from the sudden change in the delivery of cooling air flow from zero to a maximum.
The FOSI methodology is now applied to the experimental temperature measurements of room 1 of the test facility. Again, the variance of the error is , the parameter set is , and the parameter search is constrained to the following ranges: , , , and . Using the same conditions, LOA is also used for purposes of comparison. The results from the two methods are shown in Figure 9 and Table 3, where it is clear that, once again, multiplicity of solutions for the parameter set , and , occurs. From the figure it can be observed that the fractional-order model obtained by the FOSI, SA technique, accurately predicts the dynamic behavior of the temperature dynamics in room 1 for the entire time range. It is particularly remarkable that the conditions of and at are very close to those of the data. On the other hand, the results from the LOA-based model are close to the measurements in the mid-section, but deviate significantly at the extremes, particularly near the initial times.
These results are confirmed quantitatively in Table 3, where it is noticed that the fractional-order model obtained via global regression outperforms the corresponding LOA model by close to 40% improvement in the accuracy of the predictions. Moreover, the values of the parameters obtained from global regression seem to be consistent with the experimental data, particularly those of the rate of change in temperature , at times close to zero, whereas those computed by the LOA technique can be associated with physically unfeasible conditions, as mentioned in the previous section.
Finally, it is important to note that, though commonly done with integer-order calculus, a physical interpretation of the fractional derivative in the modeling of the two thermal systems analyzed here is not necessary since the corresponding fractional-order model – in each case – is merely a stand-in for the true mathematical model. Some efforts have been made in the physical interpretation of fractional-order derivatives (see [44, 45], for example), but at the moment they do not have the acceptance that terms such as velocity and acceleration have in mechanics, for instance. The issue of the physical interpretation of fractional-order models would become more important if initial conditions were being assigned based on the initial state of the plant.
6. Concluding Remarks
Complex systems, of which heat exchangers and multiroom buildings are examples, can be found in widespread applications. For their performance control it is necessary to obtain an accurate description of their dynamic behavior. In this regard, fractional-order differential equations, derived via SI, may provide compact and efficient models of complex systems. In this work, a fractional-order-based system SI procedure grounded on global regression has been proposed to build accurate models from data. The SI methodology seeks the optimum values of a parameter set that includes the fractional order of the differential equation, its parametric constants, and, when necessary, the initial conditions. This last case arises when accurate knowledge of the initial conditions, whether the value of the function or its rate of change, is not possible, and such conditions have to be included as part of the search. It is to be noted, however, that there may be cases when some values of these conditions may not be physically feasible.
The application of the SI methodology to experimental measurements from a shell-and-tube heat exchanger and a multiroom building has confirmed that the approach is accurate and robust. The fractional-based model obtained via global regression provides better approximations to the data than those obtained by SI with local optimization algorithms. The results show that not only the apparent multiplicity in the parameter set that includes the order of the proposed model equation along with its constants and the initial conditions, as provided by the results from the global regression analysis, is due to the fact that the initial condition(s) is(are) included in the search, rather than in the model itself, but also the proposed fractional-order differential equation model is reliable within the interpolation region, but cannot be used with confidence for predictions outside this region. The SI methodology proposed here is robust and can be used to derive accurate and compact models from experimental data. In addition, given a functional form of a fractional-order differential equation model, as new data become available the SI technique can be used to expand the region of reliability (interpolation) of the resulting model. The methodology described here can be extended to model other complex physical systems.
Nomenclature
| : | Fractional integral operator of order |
| : | Coefficient for the linear fractional differential equation |
| : | Riemann-Liouville fractional derivative of order |
| : | Caputo fractional derivative of order |
| : | Two-parameter Mittag-Leffler function |
| : | -th integer order |
| : | Mass flow rate |
| : | -th initial condition |
| : | Fractional order |
| : | Variance of the error |
| : | Number of experimental data points |
| : | Vector of search parameters |
| : | Fluid bulk temperature |
| : | Time |
| : | Input variable |
| : | Output variable |
| : | rms error |
| : | Gamma function |
| : | Dimensionless temperature |
| : | Independent variable |
| : | Dimensionless time |
| : | Cold fluid |
| : | Experimental value |
| : | Hot fluid |
| in: | Inlet |
| out: | Outlet |
| : | Predicted value |
| DAQ: | Data acquisition system |
| FDE: | Fractional-order differential equation |
| FOSI | Fractional-order system identification |
| LOA: | Local optimization algorithm |
| ODE: | Ordinary differential equation |
| SA: | Simulated annealing algorithm |
| SI: | System identification |
| PDE: | Partial differential equation. |
Data Availability
The data used in this study can be made available, after publication, via a hyperlink to the corresponding author’s personal website.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
Kin M. Li was the recipient of a CREST-CEaS fellowship, for which we are grateful. The comments from an anonymous referee have greatly influenced the final version of this paper and are very much appreciated. This work has been supported by an NSF HRD-1547723 grant.