HUNGARIAN JOURNAL OF INDUSTRY AND CHEMISTRY Vol. 49(1) pp. 47–58 (2021) hjic.mk.uni-pannon.hu DOI: 10.33927/hjic-2021-07 DEVELOPMENT OF COMPARTMENT MODELS FOR DIAGNOSTIC PUR- POSES BÁLINT LEVENTE TARCSAY*1, ÁGNES BÁRKÁNYI1, TIBOR CHOVÁN1, AND SÁNDOR NÉMETH1 1Research Centre for Biochemical, Environmental and Chemical Engineering, University of Pannonia, Egyetem u. 10, Veszprém, 8200, HUNGARY The importance of recognizing the presence of process faults and resolving these faults is continuously increasing par- allel to the development of industrial processes. Fault detection methods which are both robust and sensitive help to recognize the presence of faults in time to avoid malfunctions, financial loss, environmental damage or loss of human life. In the literature, the use of various model-based fault detection methods has gained a considerable degree of pop- ularity. Methods usually based on black-box models, data-based techniques or models using symbolic logic, e.g. expert systems, have become widespread. White-box models, on the other hand, have been applied less despite their consider- able robustness because of multiple reasons. Firstly, their complexity and the relatively vast amount of technological and modelling knowledge needed to construct them for industrial systems. Secondly, their large computational demand which makes them less suitable for online fault detection. In this study, the aim was to resolve these problems by developing a method to simplify the complex Computational Fluid Dynamics models employed to describe the equipment used in the chemical industry into less complex model structures. These simpler structures are Compartment Models, a type of white-box model which breaks down a complex system into smaller units with idealized behaviour. In the case of a small number of compartments, the computational load of such models is not significant, therefore, they can be employed for the purposes of online fault detection while providing an accurate representation of the system. For the purpose of identi- fying the compartmental structure, fuzzy logic was employed to create a model which approximates the real behaviour of the system as accurately as possible. Our future objective is to explore the possibility of combining this model with various diagnostic methods (expert systems, Bayesian networks, parity relations, etc.) and derive robust tools for the purpose of fault detection. Keywords: Compartment Model, Computational Fluid Dynamics, expert system, fuzzy logic 1. Introduction Fault detection and isolation has become a cardinal prob- lem in industrial practices. As the intricacy of techno- logical processes increases, the probability of faults aris- ing also grows. However, due to the stricter demands for environmentally friendly and profitable technologies, the presence of faults which could potentially cause finan- cial loss, environmental hazards or loss of life cannot be tolerated. Faults such as those defined by Venkatasubra- manian et al. [1] describe the consequences when an ob- served process variable deviates from its expected value. The source of such process faults is known as a root cause. Naturally, if left unmanaged, faults can propagate and lead to serious disasters, e.g., the Bhopal disaster or the Texas City Refinery explosion [2, 3]. To circumvent such occurrences, modern knowledge-based fault detec- tion techniques have become commonplace. These methods use the same common principle for fault detection. Firstly, some knowledge is used to cre- *Correspondence: tarcsayb@fmt.uni-pannon.hu ate a model of the system under normal operating condi- tions. Then the actual operation of the system is observed and, if the behaviour of the system does not match the assumptions drawn from the model, then an abnormality is present [4]. After identifying the presence of an ab- normality, various approaches using different techniques trace the fault signature back to the root cause depend- ing on the type of knowledge used to create the system model. Therefore, it is possible to differentiate between knowledge that originates from quantitative as well as qualitative models and that from relations based on pro- cess history [1, 5, 6]. The aforementioned categorization of knowledge used for fault detection methods stems from the early works of Venkatasubramanian [1] and is arbitrary. The expectations of fault detection and isolation methods are that they should be robust, sensitive and accurate. Robust, in this sense, means that the method should work reliably, even in the presence of noise, disturbances and changes in operating conditions. Sensitive and accurate means that the method should even be able to recognize small pro- https://doi.org/10.33927/hjic-2021-07 mailto:tarcsayb@fmt.uni-pannon.hu 48 TARCSAY, BÁRKÁNYI, CHOVÁN, AND NÉMETH cess faults and reliably identify their root causes. Usually, the bigger issue with most fault detection methods is their robustness. Commonly used models for the purposes of fault detection, e.g., a posteriori models, neural networks and statistical methods, all employ arbitrary relationships between system variables which only hold true in a cer- tain operating regime [4]. If the system itself or the operating point changes, these models often become unreliable. Therefore, they lack both robustness and flexibility. Among quantitative models used for fault detection, another type of model, the a priori or white-box model, is listed. This modelling technique uses the laws of physics, chemistry or any ap- plicable scientific discipline to describe the behaviour of various systems. Since the laws of nature cannot be vio- lated, these models are incredibly robust and hold true for a system regardless of its operating point. The problem is that these models are complex and dif- ficult to create for industrial systems since they require a vast amount of technological and modelling knowledge. Additionally, they often require a significant computa- tional load to be solved that renders them impractical for the purposes of online fault detection [1]. This is why they have not garnered widespread use in the field of di- agnostic procedures. Over the course of this work, our goal was to propose a method for simplifying the complex Com- putational Fluid Dynamics (CFD) of white-box mod- elling techniques to describe the behaviour of units used in the chemical industry. CFD methods have become widespread over the last few decades for the purpose of solving the material and energy balance equations of complex systems to describe transport phenomena as well as evaluate flow patterns, velocity fields and mixing in fluid systems. These numerical methods involve the dis- cretization of the differential equations that describe the behaviour of the unit as well as solving them through di- rect or iterative computation methods [7]. An alternative to these methods is the use of Com- partment Models (CM). CMs emulate the behaviour of complex systems by identifying individual elements of the system that can be approximated using idealized mod- els [8]. In the case of chemical technologies, the flow and transport phenomena in a system can be approximated by breaking down the system into compartments that exhibit idealized flow behaviour. CMs have become widespread in pharmacokinetic experiments to simulate drug trans- port in biological systems as well as during the evalua- tion and modelling of bioreactors and chemical systems. Savic et al. [9] used a transit CM in pharmacokinetic stud- ies involving the absorption of four different substances. They compared the results to those of a lag time model and concluded that the proposed CM produced a signifi- cantly better match to the observed experimental data and can be used as an alternative to the lag time model. Pigou and Morchain [10] employed a CM combined with a pop- ulation balance model and metabolic model to estimate the behaviour of a microbial population within industrial bioreactors. Their approach was validated by applying it to a bioreactor containing the Escherichia coli. species using experimental data. Spann et al. [11] created a risk assessment methodology for the cultivation of lactic acid bacteria based on a CM. They combined a CM with a biochemical model to monitor the species Streptococcus thermophilus online. Using Monte Carlo methods, they evaluated the risk of the biomass size not achieving its desired extent under the operating conditions of the sys- tem with regard to uncertainties in the production pro- cess. Kaur et al. [12] developed a model for the sim- ulation of a top spray fluidized bed granulator. In their work, they defined two compartments within the system, a “wet” compartment in which aggregation of the parti- cles is dominant and a “dry” compartment in which the breaking up of the particles is mainly observed. By com- bining the flow model with a population balance model, they managed to observe the effects of operational pa- rameters on the particle size distribution of the granulated product. The examples show that the CM is a contemporary and popular approach to system modelling. The main question, however, when developing a CM is how to se- lect the appropriate number and characteristics of as well as connections between compartments to obtain a model that accurately represents the flow patterns of the real system. The CM usually employs two types of idealized compartments to approximate a velocity field in a system. These are compartments which are assumed to be per- fectly stirred and homogeneous (Continuous Stirred Tank Reactor (CSTR) model) or compartments exhibiting tur- bulent plug flow (Plug Flow Reactor (PFR) model). These compartments are linked through mixer (M) and distribu- tor (D) models [13]. These latter units are usually only of mathematical significance and their volume is assumed to be zero. The component balance model for the units is showcased in the following equations. For unit CSTR: dcCSTR (t) dt = Fin VCSTR (cin (t)− cCSTR (t)) (1) FCSTR = Fin (2) For unit PFR: ∂cPFR (t, x) ∂t = −vPFR ∂cPFR (t, x) ∂x (3) FPFR = Fin (4) For unit M: cM = ∑nin i=1 cin,iFin,i∑nin i=1 Fin,i (5) FM = nin∑ i=1 Fin,i (6) Hungarian Journal of Industry and Chemistry Development of Compartment Models for Diagnostic Purposes 49 Figure 1: Step responses of ideal flow models [14] For unit D: ~cD = ~αcin (7) ~FD = ~αFin (8) The volumetric flow rate (F ) entering and exiting the units is also provided by assuming a constant temperature and dilute solutions. The models establish relationships between the volume of the unit (V ), the inlet velocity (v) and the inlet as well as actual concentrations in each ide- alized unit (c). The concentration is introduced as a func- tion of time (t) and/or Cartesian coordinates (x). In the case of the distributor unit, the outlet concentration and flow rate depend on the division rate (α). The idealized CSTR and PFR models exhibit characteristic responses to known input signals such as the Heaviside step func- tion (H) or the Dirac delta function (δ). The step response (C) of the ideal flow models as well as the serial combi- nation of CSTR units in the Tanks in Series (TIS) model are shown in Fig.1. The response functions can be used to characterize the residence time distribution (RTD) of a chemical species within the system [14]. Therefore, a popular identifica- tion method is to obtain the step response of the system and use the RTD function as a basis for identifying the structure of a compartment and the connections between the various compartments. Claudel et al. [15] employed the discipline of possibility theory to provide a method for estimating the structure of a compartment in a unit based on the RTD of the unit and the physical description of the system. They introduced “possibility” and “neces- sity” rules as well as weighing factors for these rules to assess the RTD function. They summarized the results of the rules to propose the structure of a compartment based on the characteristics of the RTD function. Egedy et al. [13] used a qualitative approach to iden- tify the structures of compartments in various systems. Their developed algorithm utilized qualitative indicators of the RTD function in an identification algorithm to pro- pose and filter out various CM structures based on their fitness. Approaches which use the CFD model as a basis for the estimation of the CM structure have also become popular in recent studies. Delafosse et al. [16] utilized CFD methods to create a CM structure for approximating the flow characteristics of a bioreactor. Fogarasi et al. [17] developed the CM of a copper leaching reactor based on experimental and CFD results. Using the CM, they opti- mized the operational parameters of the leaching reactor. Nauha et al. [18] used a hybrid CFD/CM approach to in- vestigate industrial bioreactors and predict the effects of equipment scale-up on mixing properties within the unit. Weber et al. [19] investigated the use of a hybrid CFD- based CM model for analysing a multiphase loop-reactor. They used the CM model to estimate the mass transfer and drop size distribution within the liquid-liquid extrac- tion part of the loop reactor. While all of these publications contained significant contributions towards the development of the CFD/CM hybrid model approach, they mostly suffer from the fact that the compartmentalization was conducted manually. Algorithms to automatize the identification of CM struc- tures from CFD have been proposed before. Bezzo et al. [8, 20] used spatial partitioning in their works to de- velop an algorithm that is capable of aggregating vol- umes in a unit into compartments and estimating the cross flow between the compartments based on CFD results. Tajsoleiman et al. [21] utilized a similar methodology for the classification of zones within a bioreactor that, based on CFD results, yielded distinctive flow patterns and clustered the results into a CM. Nørregaard et al. [22] developed an algorithm to identify CM structures based on CFD results by using hypothesis-driven logic. They evaluated circumferential, axial and radial bulk flows at different locations within a continuous stirred tank reac- tor. Using the obtained flow profile, they developed the CM model which resulted in a structure with 56 compart- ments. Over the course of this work, the goal was to develop an algorithm capable of identifying the CM structure of a system based on results obtained from CFD methods. The algorithm utilizes two steps to first identify the CM structure and then the flow rates between the individual compartments. Expectations for the algorithm were • Production of compartments that accurately corre- spond to the investigated volumes of the physical system • A small number of compartments to make the CM practical for online computations • Possibility to integrate empirical knowledge about the system into the compartmentalization process In future works, the possibility of using the acquired CM as a basis for a fault diagnostic system will be explored. 2. Modelling concept The proposed algorithm uses CFD results concerning the velocity field within a unit as a basis for identifying the structure of CMs. After partitioning the unit into numer- ous smaller cells, the algorithm utilizes fuzzy logic to evaluate the local velocity field within the cells and cate- gorize them into one of the idealized models. Adjacent 49(1) pp. 47–58 (2021) 50 TARCSAY, BÁRKÁNYI, CHOVÁN, AND NÉMETH Figure 2: The proposed algorithm cells that exhibit similar flow behaviour are agglomer- ated to form compartments showing idealized behaviour. The flow rate between the individual compartments is es- timated through optimization. The chart of the proposed algorithm is displayed in Fig. 2. To summarize the procedure, the compartmentaliza- tion process is conducted over six steps. 1. Calculation of a reliable steady state estimate for the velocity field of the system in question using CFD methods. 2. Definition of elementary cells (EC) by dividing the investigated volume. 3. Evaluating the flow characteristics in the individ- ual ECs through fuzzy logic and identifying them as CSTR, PFR or dead volumes. 4. Agglomerating adjacent ECs that exhibit similar flow characteristics into compartments. 5. Defining possible connections between compart- ments based on the physical geometry of the system. 6. Estimating flow rates between linked compartments through an optimization based on the step response function of the system. The steps of the proposed method and the results de- rived from it will be presented in a model system derived from the applications library of a commercial CFD simu- lator (COMSOL Multiphysics Version 5.2a). The top and side views of the system are shown in Fig. 3 with its char- acteristic dimensions. To complete the first step of the algorithm, the CFD model of the unit was developed. The system was the model of a wastewater treatment reactor which originally consisted of four baffles, two of which were removed. To model the system, the fluid within the unit was assumed to exhibit the same properties as water. To investigate the behaviour of the system, the steady-state velocity field of the unit as well as the step response function of the sys- tem to an inert tracer were necessary. To obtain these characteristics, the mass, impulse and component mass balance equations of the system have been solved numerically. The impulse balance of the sys- tem was solved assuming incompressible turbulent flow Figure 3: Top and side views of the investigated system (Remin ≈ 2 · 105 for the flow at the maximum diame- ter) within the unit using the k− ε turbulence model. The form of the Reynolds-averaged Navier-Stokes equations under these flow conditions is summarized in ρ (~v∇)~v = ∇ [ −p~I + (µ+ µT ) ( ∇ (~v) + (∇~v)T )] +F (9) ρ∇~v = 0 (10) ρ (∇~v) k = ∇ [( µ+ µT σk ) ∇k ] + Pk − ρε (11) ρ (∇~v) ε = ∇ [( µ+ µT σε ) ∇ε ] + Cε,1 ε k Pk − Cε,2 ε2 k (12) µT = ρCµ k2 ε (13) Pk = µT [ ∇~v ( ∇~v +∇~vT )−1 ] (14) The component mass balance and the step response function were calculated assuming that the velocity field within the unit had reached its steady state. No reactions were considered within the system and the change in the tracer concentration was attributed to convection and dif- fusion. The tracer was assumed to have the same phys- ical properties as the water within the tank. Under these assumptions, the component mass balance within the unit can be expressed using the partial differential equation ∂c ∂t = −~v∇c−∇ (−D∇c) (15) To solve the balance equations within the unit, a com- putational mesh was refined by taking into account both the computation time and accuracy of calculations, which was checked by calculating the relative mass balance er- ror within the unit using Emass = ( ∑nsim i=0 Ain~vin,iρin − ∑nsim i=0 Aout~vout,iρout) 2 max ( ∑nsim i=0 Ain~vin,iρin; ∑nsim i=0 Aout~vout,iρout) (16) Hungarian Journal of Industry and Chemistry Development of Compartment Models for Diagnostic Purposes 51 Table 1: Parameters of the mesh Cell type Tetrahedral Minimal size (mm) 23.1 Maximal size (mm) 77.4 Emass 5 · 10−4 Table 2: Operating conditions of the unit T (◦C) 20 Component Water vin ( ms−1 ) 0.1 pout(bar) 1 cin ( molm−3 ) 2 The mesh was defined with a construction that allows for both an effectively short computation time and ade- quate degree of accuracy. The parameters of the mesh are shown in Table 1. The steady-state velocity field and step response of the unit were calculated under the operating conditions listed in Table 2. After testing different numerical meth- ods, it was found that in the case of the steady-state impulse balance calculations, the Generalized Minimal Residuals iterative algorithm was optimal for computa- tions both with regard to computation time and computa- tional accuracy. In the case of the component mass bal- ance, by taking the same aspects into consideration, the conjugate gradient iterative method was utilized. After solving the impulse balance of the unit un- der the previously displayed conditions, the streamlined steady-state velocity profile could be observed in Fig. 4. The streamlines were defined using the internal option of COMSOL, while their densities were set proportional to the local velocity field within an observed volume. The figure shows that after the fluid enters the unit (on the left-hand side), various zones with different flow char- acteristics emerge due to the baffles. The flow entering the unit is at first unidirectional with PFR tendencies. The first baffle breaks this flow profile and creates a secluded area within which a well-mixed, circular, vortex-like flow could be observed. On the opposite side of the first baf- fle, a volume where small-scale mixing occurred could be observed. By following the main path of the flow, the area between the two baffles was a mixed regime with the presence of a significant dead volume. The area ob- structed by the second baffle also showed dead volume tendencies. The flow above the second baffle assumed a tendency that was on the verge of being mixed and PFR- like with a clear PFR tendency dominating near the out- let. During the proposed method, the velocity field was converted into a vector field which is shown in Fig. 5. The vectors on the chart represent the tangents of the veloc- ity streamlines within the unit. Their density, length and direction, as with the streamlines, was set to be propor- tional to the local magnitude and direction of the velocity Figure 4: Steady-state streamlines within the unit ( ms−1 ) Figure 5: Steady-state velocity vector field within the unit field within the reactor. The step function used to test the system and the re- sponse function of the unit are shown in Fig. 6. During the solution of the component mass balance equation dif- fusion coefficient D = 10−4 m2 s−1 was used. Based on the step response, the unit exhibited be- haviour similar to a mixed CSTR unit combined with a PFR unit. The time delay between the step and response functions (0.05 h) is indicative of PFR behaviour within the system, while the curvature of the response shows CSTR tendencies. Irregularities and jumps within the re- sponse, as can be seen at approximately 0.08 h, might be due to internal circulation flows and backflow within the unit. Figure 6: Step and response functions of the unit 49(1) pp. 47–58 (2021) 52 TARCSAY, BÁRKÁNYI, CHOVÁN, AND NÉMETH 3. Modelling results Once the CFD results had been acquired, the results were processed using MATLAB R2020b. The geometry of the system was partitioned into EC units. Over the course of this investigation, the total volume of the tank reactor (4.96 m3) was broken down into 25 rectangular ECs of equal volumes (0.2 m3) and sizes in accordance with Step 2 of the proposed algorithm (Fig. 2). Investigations using different numbers of ECs were conducted and it was con- cluded that the most accurate results can be achieved by using a number of ECs where the size of the individual cells is approximately the same as that of the cells used in the CFD simulation. After the partitioning, the vectors characterizing the velocity field displayed in Fig.5 were observed within the ECs. Since the magnitude, direction and density of the vectors are all based on the attributes of the local velocity field within the tank, key observations can be derived from them about the local flow character- istics within the unit. The following general rules were established: 1. In areas where the flow exhibits PFR characteristics, the velocity vectors are mostly unidirectional with little variance in their magnitude. The magnitude of the velocity vectors is relatively large compared to the magnitude of the inlet velocity. 2. In areas where circular flow is present that exhibits CSTR characteristics, the variance with regard to the directional components of the velocity vectors is higher. The magnitude of the velocity vectors in these areas vary but are comparable to the magni- tude of the inlet velocity. 3. In areas where little to no flow is present, also re- ferred to as dead volumes, the magnitude of the ve- locity vectors is negligible compared to that of the inlet velocity. Based on these rules, the flow characteristics within a cer- tain volume of the system can be approximately catego- rized. Therefore, the velocity vectors within the individ- ual ECs were observed and classified based on their flow characteristics. To achieve this, three measures were in- troduced. The average velocity within a cell (v1), the vari- ance in the directional components of the velocity vectors within a cell (v2) and the skewness (v3) with regard to the distribution of velocity vectors within a cell. To calculate the mean velocity within the cells, the total of the velocity vectors was calculated and averaged. Then the magnitude of the averaged vector was computed from v1 = ∣∣∣∣∑n~v i=1 ~vi (x, y, z) n~v ∣∣∣∣ (17) The variance and skewness were calculated using similar approaches. Since the coordinates determine the directions, ECs that exhibit PFR behaviour where vec- tors are unidirectional will have a high mean velocity, while those that contain circular flow with varying vec- tor directions will generally have lower mean velocities. However, the variance will increase within ECs that ex- hibit strong circular flow with velocity vectors pointing in various directions. Dead volumes have exceptionally low mean velocities. In the case of all three variables, their nu- merical values were evaluated in the individual ECs. To properly compare these parameters, the numerical value of each variable within an EC was compared to its maxi- mum value within the investigated system from ui = vi max (vi) (18) Therefore, for each variable, a quantity was produced that determines its fractional value compared to its maximum value within the unit. The values of the three parameters within the system can be seen in Fig. 7. Fig. 7 shows that the indicated variables can be used to differentiate between the flow behaviour of the unit based on the preliminary rules. The figure also explains the need for the inclusion of the skewness of the veloc- ity distribution. Given that in the area adjacent to the first baffle, the flow velocity was rather high, the variance in the magnitude and direction of the velocity vectors in that area was exceptionally high. Compared to this value, most values of variance within the system were minus- cule, except near the outlet where the circular flow met the plug-like flow exiting the reactor. The skewness pa- rameter describes the asymmetry of a given distribution function. In areas with strong circular flow patterns, the directions of the velocity vectors are more evenly dis- tributed on both sides of the mean velocity value, result- ing in a more symmetric distribution of the velocity vec- tor direction and lower skewness. However, in areas that exhibit plug flow-like behaviour, the distribution is gen- erally more skewed due to the mostly unidirectional ve- locity vectors. After obtaining the three quantities used to character- ize the flow, their values within the individual ECs were evaluated in accordance with the previously defined rules to obtain the idealized EC structure mentioned in Step 3 of the algorithm (Fig. 2). During this process, the PFR, CSTR and dead volume characteristics of all the indi- vidual ECs were examined by converting the established rules using traditional fuzzy logic [23]. Fuzzy expert sys- tems utilize empirical knowledge from subject-matter ex- perts to characterize and analyze processes [24]. They contain a working memory, inputs, a collection of rules and outputs. The working memory is the cul- mination of objects Oi, where O denotes the set of all objects. Objects consist of a fuzzy vector (~Vi) where all members of the vector are fuzzy sets or fuzzy numbers( ~Vi,1, ~Vi,2, ..., ~Vi,n ) . These fuzzy sets describe qualita- tive information about the objects, which in this case were each of the individual ECs. The inputs of the fuzzy expert system were the variables ui, which were used to describe the flow characteristics within an individual EC. Hungarian Journal of Industry and Chemistry Development of Compartment Models for Diagnostic Purposes 53 Figure 7: Characteristic velocity variables within the sys- tem A set of rules (R) was established, which in this case pertains to the flow characteristics within the unit. A fuzzy subset of R denoted as R with the member- ship function finput provides finput (R), which symbol- izes the level of belief in R for each subset [25]. The degree of belief is a number between 0 and 1 with 0 and 1 denoting disbelief and wholehearted belief in the proposed statement, respectively. The rules are obtained in the form of “IF-THEN” rules with logical operators like “OR”, “AND”, etc. In the case of our system, the real number inputs ui were converted into fuzzy numbers using monotonic sigmoidal membership functions, their general form is presented finput,j (ui) = 1 1 + e − ( ui−b a ) (19) finput,j (ui) = 1− 1 1 + e − ( ui−b a ) (20) for fuzzy numbers with index j generated from input i. The first equation is used for membership functions which are increasing in tendency, while the latter is for membership functions which are decreasing. Parameters a and b define the spread of the fuzzy numbers either side of their central value and the central values of the fuzzy numbers, respectively. The fuzzy numbers obtained through the membership functions represent the linguistic variables P belonging to the set of basic linguistic vari- ables P. In our case, these linguistic variables consisted of the set P1 = {high, low} which was used to describe whether the mean value, variance and skewness of veloc- ity vectors within an EC were high or low compared to the maximum value of the respective attribute within the system. The output actions define a fuzzy set of conclusions. The set G represents the set of possible conclusions of the expert system. A subset of this set represents the out- Table 3: Parameters of the fuzzy membership functions finput a b foutput a b Input1,low 0.01 0.05 Outputlow 0.02 0.4 Input2,low 0.01 0.05 Input3,low 0.02 0.1 Input1,high 0.03 0.2 Outputhigh 0.02 0.6 Input2,high 0.02 0.2 Input3,high 0.03 0.3 puts of our object, which in our case meant the char- acteristics of the PFR, CSTR and dead volume of in- dividual ECs [25]. Given that the linguistic variables P2 = {high, low} are represented by membership func- tion foutput,j for each characteristic, the membership functions are monotonic and sigmoidal like the input membership functions as shown in foutput,j (finput,j (R)) = 1 1 + e − ( finput,j(R)−b a ) (21) foutput,j (finput,j (R)) = 1− 1 1 + e − ( finput,j(R)−b a ) (22) In this equation, finput,j (R) represents the degree of belief in R for each subset. After defuzzification, the expert system returned crisp values that determined the PFR, CSTR and dead-volume characteristics within each EC. The input and output membership functions can be observed in Fig. 8 and Fig. 9, respectively. Since the out- put membership functions in the case of all three outputs were identical, they have only been displayed once. The first index within the legend represents the input which the membership function belongs to, while the second one corresponds to the linguistic subsets (1→ Low, 2→ High) [25]. Figure 8: Input membership functions 49(1) pp. 47–58 (2021) 54 TARCSAY, BÁRKÁNYI, CHOVÁN, AND NÉMETH Figure 9: Output membership functions Figure 10: Flow characteristics of the different EC units The parameters a and b of the membership functions are displayed in Table 3. Using the expert system, the flow characteristics of the ECs were characterized. The PFR (y1), CSTR (y2) and dead volume (y3) characteris- tics within the ECs are shown in Fig. 10. To achieve a flow map that describes the behaviour of the velocity field within each unit, the three characteris- tics were unified. The weighted average of the three out- put characteristics was computed and a function defined by Im (L) ∈ [0, 2] where L (y1, y2) = 0 corresponds to the pure dead-volume behaviour, L (y1, y2) = 1 denotes the pure CSTR behaviour and L (y1, y2) = 2 stands for the pure PFR behaviour. L (y1, y2) in each individual EC was calculated from L (y1, y2) = 2 · y1 + y2 3 . (23) Since the weighting factor of the dead-volume behaviour is 0, it was not included in the equation. The results of the process are shown in Fig. 11. The correlation between Fig. 4 and Fig. 11 is appar- ent. The model predicts clear PFR behaviour within the flow near the inlet and outlet which can also be observed from Fig. 4. The area beyond the inlet and bounded by the first baffle showed CSTR tendencies with some PFR be- haviour near the walls. A dead volume could be observed at the centre of the stirred volume which was similar to in the middle of the vortex seen in Fig. 4. The area oppo- site the volume bounded by the first baffle exhibited dead- volume tendencies with a low velocity and occasional cir- cular flow patterns, indicative of slight CSTR behaviour. The area between the two baffles showed mixed, CSTR Figure 11: Characteristic flow behaviour within the unit tendencies with a noticeable dead volume present in the middle of the vortex, as can be seen in Fig. 4. The area bounded by the second baffle mainly exhibited dead- volume tendencies with a medium degree of circulation. The volume opposite it consisted of a mixed flow. The flow near the outlet, as noted before, exhibited PFR ten- dencies. The figure presents an approximately accurate de- scription of the flow characteristics when compared to the streamlines obtained from the CFD methods. How- ever, due to fuzzification and the weighted average cal- culations to ascertain the flow tendencies within the unit, individual flow characteristics of the ECs are unclear. To clarify these values, limits were introduced to clearly cat- egorize the flow behaviour of all ECs, which were calcu- lated according to Eq. 24. The limits for this clustering process were determined by observing flow patterns ob- tained using various limits. It was found that these limits facilitated the most efficient clustering and resulted in a cluster which strongly resembled the flow characteristics observed in the CFD results. Ldisc(y1, y2) =  2, if L(y1, y2) > 1.25 1, if 1.25 ≥ L(y1, y2) ≥ 0.75 0, otherwise (24) This equation shows that the values of the corrected func- tion are discrete, that is, the ECs are categorized into three idealized classes of flow depending on which pure ideal- ized class their behaviour most resembles. The plot displaying the flowchart with ECs contain- ing the corrected values is shown in Fig. 12. The figure also indicates the locations of the inlet, outlet and baffles. Since a rectangular mesh was used to divide the system into ECs, the locations that were not part of the geome- try and the baffles were also included in the calculations. Since no flow was present at these locations, they have been categorized as clear dead volumes, which can be seen by comparing Fig. 11 and Fig. 12. For the sake of clarity, these areas were removed from the corrected fig- Hungarian Journal of Industry and Chemistry Development of Compartment Models for Diagnostic Purposes 55 Figure 12: Recognized flow behaviour within the ECs Figure 13: Defined compartments within the system ure. Fig. 12 shows a clear map indicating the idealized flow behaviour of individual ECs within the unit. Once agglomerated, adjacent ECs that exhibited similar charac- teristics could be grouped to form compartments of ideal- ized behaviour in accordance with Step 4 of the proposed algorithm (Fig. 2). The agglomeration of the ECs within this paper was conducted manually. The resulting set of compartments is shown in Fig. 13. The unit was partitioned into seven compartments, namely two PFR regions, two CSTR regions and three dead volumes. The volumes of the compartments were evaluated by totalling the volumes of the individual ECs which were assigned as part of the compartment in ques- tion. It can be seen that the resulting structure of the com- partment strongly resembles the results obtained through CFD methods seen in Fig. 4. Based on the defined struc- ture, a CM consisting of the compartments that exhibited idealized behaviour was constructed to complete Step 5 of the algorithm (Fig. 2). The CM structure is displayed Figure 14: Estimated CM structure in Fig. 14. The constructed CM structure contains the previ- ously identified compartments seen in Fig. 13. Connec- tions between the different compartments were estab- lished through the use of theoretical mixer and divider units. During the evaluation of the possible connections between compartments, the geometry of the system was taken into account as well, e.g., the placement of baffles, positions of adjacent compartments, etc. The models of the individual compartments, mixers and dividers were formulated according to the equations provided in the In- troduction. After obtaining the CM structure, Step 6 of the iden- tification algorithm (Fig. 2) involved estimating the flow rates between compartments. In the case of the structure shown in Fig.14, the flow rates between compartments were defined through the division rates (αi) of the D units. These values describe the ratio of mass flow enter- ing the ith adjacent compartment compared to the total mass flow leaving the jth reference compartment. For the identification of the αi values optimiza- tion was conducted in MATLAB R2020b through the Simulink interface. The constructed CM structure was defined in Simulink along with its connections and indi- vidual compartments. The system of differential and al- gebraic equations describing the behaviour of the system were solved numerically using Rosenbrock methods for solving stiff differential equations. The step response of the CM was observed as a function of the αi values and an objective function was formulated to estimate the α values displayed in EC(~α) = ∑nsim i=1 (CCFD,i − CCM,i) 2 max ( ∑nsim i=1 CCFD,i; ∑nsim i=1 CCM,i) (25) The objective function estimates the squared difference between the step response function acquired through CFD methods and the CM over the observed discrete- time horizon relative to the maxima of the sum of the two functions. To minimize the objective function with regard to the αi values, the interior-point method was uti- lized. The step response function of the system obtained through CFD methods and the optimized CM can be seen in Fig. 15. The resulting figure shows that the step response of the CM after optimization fits the step response of the system obtained through CFD methods well with an aver- 49(1) pp. 47–58 (2021) 56 TARCSAY, BÁRKÁNYI, CHOVÁN, AND NÉMETH Figure 15: Step response of the system with CFD methods and CM Table 4: Optimized α parameters α1 0.7 α6 0 α2 0.25 α7 1 α3 0.05 α8 0 α4 0.8 α9 0.95 α5 0.2 α10 0.05 age fit of 92%. The optimized model parameters rounded up to two decimal places are presented in Table 4. The estimated parameters in Table 4 are also in line with the general CM structure. It should be noted that the optimized α parameters of flows entering dead vol- umes (α3, α6, α8, α10) were exceptionally low, more- over, sometimes not even detectable when rounded up to two decimal places. This indicates that the step response of the system is most similar to the response obtained through CFD methods if almost no flow enters the listed volumes, indicating that they are truly dead volumes from a flow perspective. 4. Conclusions Over the course of this study, the process of creating CMs from results obtained through CFD simulations has been presented. The proposed algorithm consisted of two phases. In the first phase, the CM structure of the object was defined. For these investigations, the steady-state ve- locity field of the system was utilized. The system was partitioned into ECs in which the local velocity field was investigated. By analyzing the direction and magnitude of local velocity vectors within the ECs using fuzzy logic, a set of rules was established to categorize the behaviour of individual ECs and correlate these ECs to units exhibiting idealized behaviour. To characterize the flow behaviour, three measures were defined based on the distribution of the magnitude and direction of velocity vectors within individual ECs. Three types of ECs were considered, ones in which the flow is well-mixed and circular, similar to a CSTR, ones in which the flow is unidirectional and resembles that found in a PFR, and ones in which little to no flow is present, that is, dead-volume behaviour is exhibited. After the qualitative analysis using fuzzy logic, each EC was assigned an idealized trait. Adjacent ECs that ex- hibited similar traits were clustered into compartments showing idealized behaviour. After determining the ide- alized compartments, the flow rates between these com- partments were estimated through optimization by calcu- lating the flow rates for which the sum of the squared difference between the step response of the CFD model and CM was minimal. The estimated parameters and the observed compartments showed a strong correlation with the flow patterns which could be observed within the CFD results. The method is capable of effectively ana- lyzing CFD results and converting them into CMs in a structured framework while also allowing the addition of empirical knowledge using fuzzy logic. The proposed method could serve as a basis for the estimation of ro- bust and computationally inexpensive a priori models to describe equipment utilized within the chemical indus- try. The physical correlation between the model compart- ments and the actual flow regimes within the unit could provide opportunities to develop a priori model-based di- agnostic strategies for online fault detection among many other applications. 5. Acknowledgements This work was supported by the TKP2020-IKA-07 project financed under the 2020-4.1.1-TKP2020 The- matic Excellence Programme by the National Research, Development and Innovation Fund of Hungary. Nomenclature Fuzzy sets ~V Set of fuzzy vectors/numbers G Set of possible fuzzy conclusions O Set of fuzzy objects P Set of linguistic fuzzy variables R Set of fuzzy rules Greek letters α Distribution ratio [−] ε Turbulent dissipation rate [m2 s−3] µ Kinematic viscosity [m2 s−1] µT Turbulent viscosity [kgm−1 s−1] ρ Density [kgm−3] σε Constant for change in ε [−] σk Constant for change in k [−] τ Residence time [h] Latin letters ~F Force [N] Hungarian Journal of Industry and Chemistry Development of Compartment Models for Diagnostic Purposes 57 ~I Impulse [kgm s−1] ~v Velocity [m s−1] A Surface [m2] a Fuzzy membership function central value [−] b Fuzzy membership function variance [−] C Step response function [molm−3] c Species concentration [molm−3] Cε,1 Dissipation rate calculation constant 1 [−] Cε,2 Dissipation rate calculation constant 2 [−] Cµ Turbulent viscosity calculation constant [−] D Diffusion rate [m2 s−1] E Error function [−] F Flow rate [m3 s−1] f Fuzzy membership function H Heaviside (Step) function [molm−3] i, j Numeric indices [−] k Turbulent kinetic energy [m2 s−2] L Continuous flow characteristic function [−] Ldisc Discretized flow characteristic function [−] n Absolute number [−] p Pressure [Pa] Pk Turbulent kinetic energy due to mean velocity variations [kgm−1 s−1] T Temperature [◦C] t Time [h] u Fuzzy expert system input variable [−] V Volume [m3] x Longitudinal coordinate [m] y Fuzzy expert system output variable [−] z Lateral coordinate [m] REFERENCES [1] Venkatasubramanian, V.; Rengaswamy, R.; Yin, K.; Kavuri, S.N.: A review of process fault detection and diagnosis: Part I: Quantitative model-based methods., Comput. Chem. Eng., 2003, 27(3), 293– 311 DOI: 10.1016/S0098-1354(02)00160-6 [2] Broughton, E.: The Bhopal disaster and its after- math: a review., Environ. Health, 2005, 4(1), 1–6 DOI: 10.1186/1476-069X-4-6 [3] Isimite, J.; Rubini, P.: A dynamic HAZOP case study using the Texas City refinery explosion., J. Loss. Prev. Process. Ind., 2016, 40(1), 496–501 DOI: 10.1016/j.jlp.2016.01.025 [4] J., G.: Fault detection and diagnosis in engineering systems. (CRC press), 1998. ISBN: 9780367400439 [5] Venkatasubramanian, V.; Rengaswamy, R.; Kavuri, S.N.: A review of process fault detection and diag- nosis: Part II: Qualitative models and search strate- gies., Comput. Chem. Eng., 2003, 27(3), 313–326 DOI: 10.1016/S0098-1354(02)00161-8 [6] Venkatasubramanian, V.; Rengaswamy, R.; Kavuri, S.N.; Yin, K.: A review of process fault detection and diagnosis: Part III: Process history based meth- ods., Comput. Chem. Eng., 2003, 27(3), 327–346 DOI: 10.1016/S0098-1354(02)00162-X [7] Solin, P.; Segeth, K.; Dolezel, I.: Higher-order fi- nite element methods. (CRC press), 2003. ISBN: 9780429205279 [8] Bezzo, F.; Macchietto, S.; Pantelides, C.C.: A general methodology for hybrid multizonal/CFD models: Part I. Theoretical framework., Com- put. Chem. Eng., 2004, 28(4), 501–511 DOI: 10.1016/j.compchemeng.2003.08.004 [9] Savic, R.M.; Jonker, D.M.; Kerbusch, T.; Karls- son, M.O.: Implementation of a transit compart- ment model for describing drug absorption in phar- macokinetic studies., J. Pharmacokinet. Pharmaco- dyn., 2007, 34(5), 711–726 DOI: 10.1007/s10928-007- 9066-0 [10] Pigou, M.; Morchain, J.: Investigating the in- teractions between physical and biological heterogeneities in bioreactors using compart- ment, population balance and metabolic models., Chem. Eng. Sci., 2015, 126(1), 267–282 DOI: 10.1016/j.ces.2014.11.035 [11] Spann, R.; Gernaey, K.V.; Sin, G.: A compartment model for risk-based monitoring of lactic acid bac- teria cultivations., Biochem. Eng. J., 2019, 151(1), 2–11 DOI: 10.1016/j.bej.2019.107293 [12] Kaur, G.; Singh, M.; Matsoukas, T.; Kumar, J.and De Beer, T.; Nopens, I.: Two-compartment model- ing and dynamics of top-sprayed fluidized bed gran- ulator., Appl. Math. Model., 2019, 68(1), 267–280 DOI: 10.1016/j.apm.2018.11.028 [13] Egedy, A.; Varga, T.; Chován, T.: Compart- ment model structure identification with qualita- tive methods for a stirred vessel., Math. Com- put. Model. Dyn. Syst., 2013, 19(2), 115–132 DOI: 10.1080/13873954.2012.700939 [14] Levenspiel, O.: Chemical reaction engineering. (John Wiley & Sons), 1999. ISBN: 978-0-471-25424-9 [15] Claudel, S.; Fonteix, C.; Leclerc, J.P.; Lintz, H.G.: Application of the possibility theory to the compart- ment modelling of flow pattern in industrial pro- cesses., Chem. Eng. Sci., 2003, 58(17), 4005–4016 DOI: 10.1016/S0009-2509(03)00269-0 [16] Delafosse, A.; Collignon, M.L.; Calvo, S.; Delvi- gne, F.; Crine, M.; Thonart, P.; Toye, D.: CFD-based compartment model for description of mixing in bioreactors., Chem. Eng. Sci., 2014, 106(1), 76–85 DOI: 10.1016/j.ces.2013.11.033 [17] Fogarasi, S.; Egedy, A.; Imre-Lucaci, F.; Varga, T.; Chován, T.: Hybrid CFD-compartment approach for modelling and optimisation of a leaching reactor., Comput. Aided Chem. Eng., 2014, 33(1), 1255– 1260 DOI: 10.1016/B978-0-444-63455-9.50044-1 [18] Nauha, E.K.; Kálal, Z.; Ali, J.M.; Alopaeus, V.: Compartmental modeling of large stirred 49(1) pp. 47–58 (2021) https://doi.org/10.1016/S0098-1354(02)00160-6 https://doi.org/10.1186/1476-069X-4-6 https://doi.org/10.1016/j.jlp.2016.01.025 https://doi.org/10.1016/j.jlp.2016.01.025 https://doi.org/10.1016/S0098-1354(02)00161-8 https://doi.org/10.1016/S0098-1354(02)00162-X https://doi.org/10.1016/j.compchemeng.2003.08.004 https://doi.org/10.1016/j.compchemeng.2003.08.004 https://doi.org/10.1007/s10928-007-9066-0 https://doi.org/10.1007/s10928-007-9066-0 https://doi.org/10.1016/j.ces.2014.11.035 https://doi.org/10.1016/j.ces.2014.11.035 https://doi.org/10.1016/j.bej.2019.107293 https://doi.org/10.1016/j.apm.2018.11.028 https://doi.org/10.1080/13873954.2012.700939 https://doi.org/10.1080/13873954.2012.700939 https://doi.org/10.1016/S0009-2509(03)00269-0 https://doi.org/10.1016/j.ces.2013.11.033 https://doi.org/10.1016/B978-0-444-63455-9.50044-1 58 TARCSAY, BÁRKÁNYI, CHOVÁN, AND NÉMETH tank bioreactors with high gas volume fractions., Chem. Eng. J., 2018, 334(1), 2319–2334 DOI: 10.1016/j.cej.2017.11.182 [19] Weber, B.; von Campenhausen, M.; Maßmann, T.; Bednarz, A.; Jupke, A.: CFD based compartment- model for a multiphase loop-reactor., Chem. Eng. Sci.: X, 2019, 2, 100010 DOI: 10.1016/j.cesx.2019.100010 [20] Bezzo, F.; Macchietto, S.: A general methodology for hybrid multizonal/CFD models: Part II. Auto- matic zoning., Comput. Chem. Eng., 2004, 28(4), 513–525 DOI: 10.1016/j.compchemeng.2003.08.010 [21] Tajsoleiman, T.; Spann, R.; Bach, C.; Gernaey, K.V.; Huusom, J.K.; Krühne, U.: A CFD based au- tomatic method for compartment model develop- ment., Comput. Chem. Eng., 2019, 123(1), 236–245 DOI: 10.1016/j.compchemeng.2018.12.015 [22] Nørregaard, A.; Bach, C.; Krühne, U.; Borgb- jerg, U.; Gernaey, K.V.: Hypothesis-driven compart- ment model for stirred bioreactors utilizing com- putational fluid dynamics and multiple pH sen- sors., Chem. Eng. J., 2019, 356(1), 161–169 DOI: 10.1016/j.cej.2018.08.191 [23] Zadeh, L.: Fuzzy logic., Computer, 1988, 21(4), 83–93 DOI: 10.1109/2.53 [24] Hadjimichael, M.: A fuzzy expert system for avi- ation risk assessment., Expert Syst. Appl., 2009, 36(3), 6512–6519 DOI: 10.1016/j.eswa.2008.07.081 [25] Buckley, J.J.; Siler, W.; Tucker, D.: A fuzzy expert system., Fuzzy Sets Syst., 1986, 20(1), 1–16 DOI: 10.1016/S0165-0114(86)80027-6 Hungarian Journal of Industry and Chemistry https://doi.org/10.1016/j.cej.2017.11.182 https://doi.org/10.1016/j.cej.2017.11.182 https://doi.org/10.1016/j.cesx.2019.100010 https://doi.org/10.1016/j.compchemeng.2003.08.010 https://doi.org/10.1016/j.compchemeng.2018.12.015 https://doi.org/10.1016/j.cej.2018.08.191 https://doi.org/10.1016/j.cej.2018.08.191 https://doi.org/10.1109/2.53 https://doi.org/10.1016/j.eswa.2008.07.081 https://doi.org/10.1016/S0165-0114(86)80027-6 https://doi.org/10.1016/S0165-0114(86)80027-6 Introduction Modelling concept Modelling results Conclusions Acknowledgements