Predictive Control Strategy for Continuous Production Systems: A Comparative Study with Classical Control Approaches Using Simulation-Based Analysis

Abstract

Due to today’s technological development and information progress, an increasing number of physical systems have become interconnected and linked together through communication networks, thus resulting in Cyber-Physical Systems (CPSs). Continuous manufacturing, which involves the manufacture of products without interruption, has become increasingly important in many industries, including the pharmaceutical and chemical industries. CPSs can be used to control and monitor the production process, which is essential in enabling continuous manufacturing. This paper is focused on the modeling and control of physical systems required in tablet production using dry granulation. Tablets are a type of oral dosage form that is commonly used in the pharmaceutical industry. They are solid, compressed forms of medication that are formulated to release the active ingredients in a manner that allows for optimal absorption and efficacy. Thus, a model predictive control (MPC) strategy is applied to a plant model to test the designed controller and to analyze the obtained performances. The simulation results are compared with those obtained using other control algorithms, linear quadratic regulator (LQR) and proportional-integral-derivative (PID), applied to the same plant model. The results showed that the predictive control strategy performed significantly better than the other two control strategies.

1. Introduction

Advances in data acquisition systems, information technology (IT), and networking technologies have led to the development of Cyber-Physical Systems (CPSs). These systems are complex infrastructures that integrate digital technology into physical reality in order to automate and digitalize physical processes [1]. Research and advances in CPS technologies have increasingly become part of emerging trends in IT fields such as the Internet of Things [2], Big Data [3], Cloud Computing, and Artificial Intelligence [4]. Cyber-physical systems offer numerous benefits to technology, industries, and organizations. They improve efficiency, productivity, and safety in manufacturing, healthcare, transportation, and other sectors. By integrating physical devices and software systems, cyber-physical systems enable real-time monitoring, analysis, and control of complex processes, resulting in enhanced performance, reduced costs, and better decision-making capabilities. As a result, cyber-physical systems have the potential to revolutionize multiple industries and transform the way organizations operate. In the context of Industry 4.0, next-generation technologies have propelled the development of smart manufacturing due to the many benefits they offer [5]. The economic benefits are among the most important for manufacturers as they can reduce their production costs by reducing manpower due to process automation, increasing production, and reducing raw material losses. Increasing production volume and quality of the yielded products are important factors in the industry’s technological development due to the increasing demand on the market.

Because CPSs use both hardware and software resources to integrate computing, communication, and control processes, they are used in many fields such as automotive, aerospace, transportation, energy, and manufacturing. Among the many applications of CPS is the pharmaceutical industry, which needs to be upgraded due to the increasing demand for products, highlighted during the COVID pandemic.

The first step towards the modernization of industry is the transition to continuous production, necessary for the emergence of smart factories [6]. Continuous production is a production flow method used to manufacture, produce, or process materials without interruption. During continuous processing, materials are transported from one processing unit to another in a continuous flow to produce a finished product. Regardless of the industry in which it is applied, continuous production offers the major advantage of a higher production rate as well as an increased standard of product quality. It also decreases production time, the possibility of human error, and per unit costs, while production flexibility increases and the ability to scale up is improved. This type of process is used by companies involved in oil refining, metal smelting, power plants, sanitary wastewater treatment, and many more. Other industries, such as pharmaceuticals, are beginning to understand the benefits of a continuous production process, whereby products can be manufactured in a shorter time, at the lowest possible cost, while minimizing waste [7].

At the moment, most pharmaceutical companies are focused on batch production because it is cheaper compared to continuous manufacturing (CM) [8]. This approach requires each batch to go through the same steps at the same time. Unlike continuous production, where both raw material and finished product are charged/discharged continuously throughout the process, in batch production, all raw material is charged at the start of the process and the products are discharged at once [9]. This is a disadvantage because quality checks and inspection of products are made on samples from an already manufactured batch, compared to continuous production where products can be inspected in advance at each stage.

Pharmaceutical products can be found in many forms, depending on the purpose of use and how they are administered. Thus, they can be classified into oral (pills/tablets, capsules, syrups, solutions), cutaneous (unguents, gels, creams, solutions), and inhalation (solution, spray) pharmaceutical forms. Depending on the classification of the products, the method of manufacture also differs. This paper is focused on solid dosage forms, in which case the research carried out so far on the continuous production of these types of drugs still requires improvement [8]. In [10], the risks and economics of adopting continuous production for pharmaceutical tablets compared to conventional batch production are presented. The results of the studies demonstrate the potential of CM to make internal production of pharmaceuticals more economically attractive than external production of these products. A case study on the sustainability of using CM in pharmaceutical production is also illustrated in [11]. The examples analyzed in this article show that in most situations there are positive sustainability factors for the use of this method in the manufacturing process of chemical compounds needed for pharmaceutical products.

The transition to continuous manufacturing of pharmaceuticals is supported by a report of the Subcommittee on Advanced Manufacturing, Committee on Technology of the National Science and Technology Council, US, in October 2022 [12], which recommends “advanced continuous manufacturing, in-line process monitoring and control, integrated AI-assisted systems, and novel cell culture techniques".

Solid dosage forms can be classified according to the manufacturing process as follows: direct compression (DC), dry granulation (DG), or wet granulation (WG) [7]. Continuous direct compression includes multiple loss-in-weight feeders, a blender required for mixing active pharmaceutical ingredients (APIs) with powdered excipients, and a rotary tablet press. Optionally, a mill (for grounding ingredients from feeders) and a feeder lubricant can be added before the blender. Additionally, in the case of the WG route, a twin screw granulator is used to grind the primary powder. After that, the dryer reaches the desired moisture level and a mill ensures that the granules are the desired size. In addition to the processing units used in DC, a roller compactor and a mill are required in DG after blending the ingredients. The manufacturing process should be chosen according to the physical properties of the tablets (weight variation, thickness, hardness, diameter) and tablet performance (disintegration and dissolution).

Propelled by market requirements and technological advancement, the pharmaceutical industry is encouraged to adopt continuous manufacturing to increase the efficiency and flexibility of drug production by optimizing costs and production times and enhancing product quality. In this work, the focus is on developing and testing a novel model predictive control strategy for such processes. The proposed algorithm is designed for a simulator (developed in Simulink) of the continuous dry granulation manufacturing process of pharmaceutical tablets. To better illustrate the performance of the proposed strategy, the results obtained from the simulations are compared with those provided by two classical control approaches: PID and LQR. Thus, it can be noticed that MPC provides significantly improved results in terms of stability of the production process and increased quality of the final products.

Considering all the aspects mentioned above, the main contributions of this paper, in comparison with our previous work and the current literature, are the following:

• A new predictive control strategy, adapted for a continuous pharmaceutical tablet manufacturing plant using dry granulation, was developed.
• The control performance of the algorithm was tested and analyzed by simulations performed on a benchmark simulator designed based on the models that are available in the literature.
• A comparison between the results provided by the proposed predictive control strategy and those obtained using PID and LQR algorithms was performed.

In this paper, a particular application, i.e., oral dosage form production in a continuous fashion is used as a study case to test and evaluate the fundamental algorithms that have been developed. The purpose of this paper is to illustrate the advantage of modeling and control techniques in the continuous production of pharmaceuticals. The established methods used in the manufacture of solid dosage forms are compared in Section 2, while Section 3 details the dry granulation process model. In this paper, a predictive control strategy is applied to a novel pilot plant of a tablet manufacturing process, which considers system constraints at the design stage. Section 4 focuses on the control architecture of the system and the predictive control strategy applied to the process. The method proposed in the current paper has significant improvements, which will be illustrated and compared with other control strategies in Section 5. The conclusions of these results and future directions of work are presented in Section 6.

2. Direct Compression vs. Dry/Wet Granulation

Either dry granulation (DG) or wet granulation (WG), the use of a granulation method in tablet manufacture is preferable due to the issues of poor fluidity, segregation, and content uniformity. The main advantages of DG compared to WG are the simplified operating process and the absence of the additional energy-consuming drying step. DG is also suitable for processing temperature-sensitive and easily hydrolyzed compounds [13].

In batch production, direct compression (DC) is considered a risky method mainly due to the risk of segregation of the mixture. Therefore, formulations that tend to segregate are often granulated to ensure consistent uniformity of active ingredients contained in the final product. However, it has been observed that when production is performed in a continuous mode, even using DC, increased homogeneity can be achieved for products with high segregation tendencies due to the minimization of internal gaps and elimination of semi-continuous steps [14]. The strong tendency of APIs to agglomerate is a strong argument for the use of granulation in batch manufacturing. This problem is completely solved by the use of a properly designed continuous system in which the material lacks the time and the environment to re-agglomerate [15]. Compared to DC and DG, WG is one of the most difficult pharmaceutical processes to model due to the increased number of steps and variables used in system design [13].

In terms of scale-up, the feasibility of direct compression is negatively influenced due to speed-sensitive compactibility or flow limitations compared to WG. Dry granulation based on roller compaction increases density and flow [16]. The simplicity of the DC process compared with (dry/wet) granulation reduces the complexity, risks, and costs of processing materials. Technically, any additional operations introduced by the granulation process come with validation, yield, cleaning, and documentation issues, in addition to time and effort issues in high-cost good manufacturing practice (GMP) containment facilities [15]. In addition, the heat and humidity problems that can occur in WG may be unacceptable for labile active substances [16]. In direct compression, product quality and consistency are essential, an advantage offered by the possibility of direct expression of input material properties [16]. However, the simplicity of the DC process requires additional constraints on the formulation components compared to production processes involving granulation [17].

Another aspect of the differences between the three routes of manufacture of tablets is the coloring of the tablets. DC and DG are unable to provide the same color intensity compared to WG [18]. In dry granulation, the powder mixture is compressed into granules without the use of a liquid binder. Since no liquid is added to the mixture, the risk of color variations due to moisture or reaction with the binder is minimized, but in its absence, the resulting mixture may have an inconsistent distribution of colorants in the final tablets. Wet granulation involves the addition of a liquid binder to the powder mixture, which can help to distribute the colorant evenly throughout the granules. In the case of direct compression, the homogeneity of the mixture can be more challenging to achieve, which can result in uneven distribution of the colorant in the final tablets.

In [16], the case of sodium starch glycolate is exemplified whose insolubility property does not exclude incompatibilities associated with sodium salts, especially in the presence of moisture. The reaction rate of the finished product is directly proportional to the water activity in it. In this situation, it is not recommended to use DC and DG methods because water activity is probably inhomogeneous in the tablet matrix and local microambitates with high water activity may exist. Another concrete example also illustrated in [16] is for vitamin tablets for which both WG and DG are used. However, in this case, it is recommended to use the DG method which offers a higher yield and more stable products due to the removal of water.

A summary of the advantages/disadvantages of the three production routes for pharmaceutical tablets is presented in

Table 1. Comparison of Wet/Dry Granulation and Direct Compression [13,14,15,16,17,18,19,20].

Property	Wet Granulation	Dry Granulation	Direct Compression
Compactability	Harder tablets in case of hard compactable substances.	Influenced by powder particle size and shape.	Potential problem for high loading of poorly compactable API substances.
Flow	The granules formed are slightly more spherical than powders and have better flowability.	There may be some issues with powder flow.	Raw materials must have proper flowability and mixed with APIs, sometimes they may need lubricants before compression.
Particle size	Greater with a longer range.		Narrower with a narrower range.
Content uniformity	It ensures better uniformity of content.	The resulting granulation increased confidence in uniformity.	It is at risk because it is difficult to accurately mix a small amount of API into a large volume of excipients.
Mixing	Prevents segregation of components.	Segregation of components may occur after mixing.	A high shear can reduce particle size.
Lubrication	Not so sensitive.	The compression step becomes easier and not sticky.	Reduces mixing time.
Disintegration	Increased intragranular levels are required because of the negative impact of wet granulation on disintegrants.	They have an improved disintegration time because the dry binder used has a lower adhesive effect and therefore a quicker disintegration.	It allows them to disintegrate into API particles rather than granules.
Dissolution	Providing hydrophilic properties to the surface of the granules can improve the dissolution rate.	The slowness of dissolution from granules during storage, particularly if an intragranular disintegrant is not used.	Difference in dissolution speeds up the process and allows better absorption for API tablets that are poorly soluble.
Cost	Higher investments costs because of the time, labour, energy, and equipment.	Lower equipment costs than wet granulation.	DC has an economic advantage over granulation as it requires fewer resources.
Sensitivity to raw material variability	Raw material wetting is influenced more by changes in raw material properties.	The properties of the raw material matter, the characteristics of API powders, and excipients are important.	Precise selection of excipients is needed as raw materials must have adequate flowability and compressibility for a successful operation.
Stability	Not suitable for use on heat or moisture-sensitive materials.	Suitable for using on heat or moisture-sensitive materials.
Tableting speed	Higher		Decreased speed if the flow is low.

.

3. Process Description

In the context of Cyber-Physical Systems, the pharmaceutical industry needs to take a new approach to the manufacturing process. Currently, the manufacturing process is performed in batches, where the materials are stored before being sent to the next processing stage and their quality is tested off the production line [9]. In the case of continuous manufacturing, the products obtained in each stage are sent directly to the next processing unit. Thus, the design of a continuous manufacturing system involves a complex control strategy to ensure the required quality standards, as well as proper process control. All this is then quantified in increased productivity and reduced production and manufacturing times. Continuous production, unlike batch production, is a flow production strategy. This means that the process never stops, and there are no breaks between the different steps of product creation.

3.1. Process Structure

Before being released for distribution and consumption, tablets are subjected to quality tests, such as the content of active ingredient test, weight uniformity test, dissolution test, or hardness test. The most common defects in compressing table processes are weight variation, low hardness/low mechanical strength of tablets, poor mixing, or layered tablet splitting [21].

For tablet manufacturing, the most common processes are wet granulation, direct compression, or dry granulation. The first two continuous manufacturing methods are shown in Figure 1, where, with red, the route of wet granulation is exemplified and, with green, the route via direct compaction.

This work is focused on the dry granulation process of tablets, modeled in Matlab-Simulink [22]. The model designed for the process considered in this paper is based on numerical data found in the literature. The main execution elements included in the designed model are the feeders for active pharmaceutical ingredients and excipients, the blender and the hoppers associated with the roller compactor (RC) and tablet press (TP). Figure 2 shows the order of the production steps and the dependencies between the processing units. As shown in the above-mentioned figure, the feeders of dry ingredients are connected to a continuous blender that mixes all these ingredients. This results in a powder mixture that feeds a roller compactor that converts the powder blend into the ribbon. A mill is integrated with the roller compactor that breaks down the ribbon and forms granules. These granules are fed to the tablet press via a rotary feed frame. The powder blend fills a die and is subsequently compressed to create a tablet; after compaction, the tablets are coated [23].

3.2. Process Modelling

The main purpose of the unit operation of feeding is to be able to supply APIs, excipients, and lubricants to the next unit, considering a certain constant flow rate, to match the desired mixture composition and total material flow rate. For the considered continuous process, the screw feeder is mostly used. It consists of a feed hopper, one (or more) bridge-breaking systems, and a single (or twin) configuration of conveying screws through which material is distributed out of the unit [24]. The great advantage of this type of power supply is the ability to manage very low flow rates, providing good accuracy. In this operation, an important parameter is the feed factor, which is defined as the amount of powder present in a screw pitch at each rotation of the feed screw. It can be calculated as follows:

$$ff=\frac{\dot{m}}{\omega },$$

(1)

where $\dot{m}$ is the mass flow rate and $\omega$ is the screw speed.

Note that, the weight in the feeder hopper is assumed to be constant. Using Heckel’s equation [25], the changes in density inside the feeder unit are characterized as:

$${\rho }_{feeding}\left(W\right)={\rho }_{max}-({\rho }_{max}-{\rho }_{min})e^{-K\frac{Wg}{A_{feeder}}},$$

(2)

where W is the weight of the material inside the feeder, ${\rho }_{feeding}$ is the perceived density that is present at the inlet of the screw, ${\rho }_{max}$ is the density at which the perceived density remains constant while more materials are being added to the feeder, $q_{min}$ is the minimum density and occurs when $W\to 0$, meaning that the feeder is nearly empty, K is the modified pressure decay constant, g is the gravity acceleration, and $A_{feeder}$ is the cross-sectional area of the feeder.

Since the density is linked with the feed factor, it can be used to characterize this as well. In [26,27], a feed factor model based on Heckel’s equation was developed and analyzed, in the following form:

$$f{f}_{aparent}\left(W\right)=f{f}_{max}-(f{f}_{max}-f{f}_{min})e^{-\beta W},$$

(3)

where $f{f}_{max}$ corresponds to the maximum relative bulk density within the screws, and $f{f}_{min}$ is the minimum relative bulk density or powder compressibility within the screws. The constant parameter $\beta$ inside the exponential represents the rate at which the amount of powder present in the screw pitch decreases with respect to a decrease in feeder weight [28] and is defined as:

$$\beta =\frac{Kg}{A_{feeder}}.$$

(4)

In order to capture the feed factor variability that is a consequence of the noise, a feed factor noise term is added to Equation (3), resulting in:

$$ff\left(W\right)=f{f}_{aparent}+f{f}_{var}.$$

(5)

A truncated cosine Fourier series is used to realize this variability and the following equation results:

$$f{f}_{var}\left(W\right)=\sum _{j=1}^3{A}_{j}cos(W-P_j),$$

(6)

where $A_j$ and $P_j$ are constant material dependent parameters.

The blender unit guarantees uniformity of concentration and homogenization. The construction of the mathematical model for a continuous blending process is based on the population equation [29]:

$$\frac{\partial F(x,t)}{\partial t}+\frac{\partial }{\partial x}\left(F(x,t)\frac{dx}{dt}\right)=R_{form}(x,t)-R_{dep}(x,t)+{\dot{F}}_{in}(x,t)-{\dot{F}}_{out}(x,t),$$

(7)

where F is the population distribution function, x is the vector of internal coordinates used to express the particle size and t is the time. $R_{form}$ and $R_{dep}$ indicate the particles that are formed and depleted, respectively, while ${\dot{F}}_{in}$ and ${\dot{F}}_{out}$ are the rates at which the components are fed/neglected to/by the system.

Due to the complexity of the system and the lack of numerical parameters, the blender has been divided into several compartments. These local homogeneous compartments are created by discretizing the mixing vessel space in both axial and radial dimensions. Applying Equation (7) to the process under consideration results in a simple mass balance of a compartment containing the powder flows:

$$\frac{\partial m_{i,j}}{\partial t}=F_f[m_{i-1,j}-m_{i,j}]+F_b[m_{i+1,j}-m_{i,j}]+F_r[m_{i,j+1}+m_{i,j-1}-m_{i,j}],$$

(8)

where $m_{i,j}$ is the mass holdup in the compartment and the forward $\left(F_f\right)$, backward $\left(F_b\right)$, and radial $\left(F_r\right)$ fluxes are obtained using:

$$\begin{array}{cc}\hfill & F_f=a{\omega }_{blender}+b,\hfill \\ \hfill & F_b=c{\omega }_{blender}+d,\hfill \\ \hfill & F_r=e,\hfill \end{array}$$

(9)

in which ${\omega }_{blender}$ is the speed of the blender, while constant parameters $a,b$ (for forward flux), $c,d$ (for backward flux), and e (for radial flux) are estimated from experimental data, depending on the type of material (API or excipient).

Considering the following definition for the flow at the outlet of the blender:

$$F_{blender}=\sum _{j=1}^{N_r}{F}_{f,API}{m}_{API,i=N_a,j}+\sum _{j=1}^{N_r}{F}_{f,Excipient}{m}_{Excipient,i=N_a,j},$$

(10)

where $m_{API}$ and $m_{Excipient}$ are the powder holdups of the last compartments in the axial direction, $N_r$ represents the number of radial discretization, and $N_a$ is the number of axial compartments, the outlet concentration $\left(C_{API}\right)$ is resolved:

$$C_{API}=\frac{{\sum }_{j=1}^{N_r}{F}_{f,API}{m}_{API,i=N_a,j}}{F_{blender}}.$$

(11)

The RC hopper and the TP hopper used in this paper are working on the same principle. The hoppers modeled in this work consist of a conical section at the bottom and a cylindrical section at the top. The hopper model described below will finally provide the resulting height and mass flow rate. In steady-state, the mass flow through all units is the same. Thus, both hoppers will operate under the same operating conditions. This fact causes a delay in the density coming from the previous unit (blender or roller compactor).

The powder that leaves the blender comes into a hopper as well as the material that exits the roller compactor [22]. This unit is modeled as a conical, spherical hopper-type tank. The dynamics of the system are guided by the mass balance equation, defined as the difference between the inlet mass flow rate that is coming from the blender $\left({\dot{m}}_{in}\right)$ and the outlet mass flow rate $\left({\dot{m}}_{out}\right)$:

$${\rho }_b\frac{dV}{dt}={\dot{m}}_{in}-{\dot{m}}_{out},$$

(12)

where ${\rho }_b$ is the density at the outlet of the hopper and V is the volume of material in the hopper.

The dynamics of the hopper are modeled by the change in volume in the hopper, depending on the powder’s height:

$$\begin{array}{cc}\hfill h<H_1& :\frac{dV}{dt}=\frac{1}{3}\pi \left[{\left(\frac{R_2}{H_1}h\right)}^2+R_1\left(\frac{R_2}{H_1}h\right)+R_1^2\right]\frac{dh}{dt},\hfill \\ \hfill h\ge H_1& :\frac{dV}{dt}=\pi R_2^2\frac{dh}{dt},\hfill \end{array}$$

where $R_1$ is the radius of the conic section in the bottom and $R_2$ is the radius of the cylindrical section at the top of the hopper, h is the height of the powder in the hopper and $H_1$ is the height of the considered section.

The roller compaction process involves turning a powder mixture into a ribbon using a set of counter-rotating rolls that will compress the material. When the particles are fed to the rolls, they are initially considered to be in the slip region of the process [30]. The nip region is the area where the rollers are closest together. The last area, namely the release region, is where the ribbon is already released from the rolls and transported forward. This process is a complex one, depending on a high number of parameters and operating conditions.

An important parameter in this model is the nip angle, defined as the angular position $(\theta =\alpha )$, at which the material transforms from slip to non-slip conditions. It depends only on material-related parameters and the following equation is used to calculate it:

$$\frac{4\left(\frac{\pi }{2}-\alpha -\nu \right)tan{\delta }_E}{cot(B-\mu )-cot(B+\mu )}=\frac{K\left(2cos\alpha -1-\frac{h_0}{R}\right)tan\alpha }{cos\alpha },$$

(13)

where $\nu$ is a variable calculated in relation to the effective angle of friction $\left({\delta }_E\right)$ and the angle of wall friction, K is the compressibility factor, $h_0$ is half of the roll gap width, and R is the radius of the roll. The variables B and $\mu$ are calculated as:

$$\begin{array}{cc}\hfill & B=\frac{\alpha +\mu +\frac{\pi }{2}}{2},\hfill \\ \hfill & \mu =\frac{\pi }{4}-\frac{{\delta }_E}{2}.\hfill \end{array}$$

(14)

In the nip region, the compression behavior can be calculated using an empirical model, $\sigma =C_1{\rho }^K$, where $\sigma$ is the material stress and $\rho$ is the compact density, while $C_1$ is the pre-exponential coefficient of material compression.

Based on Johanson’s model [31] and mass balance equation, the roller compactor model is computed using the compact density at the exit point $\left({\rho }_{exit}\right)$ and the inlet density $\left({\rho }_{in}\right)$:

$$\begin{array}{cc}\hfill & {\int }_0^{{\theta }_{in}}\rho \left(\theta \right)cos\theta d\theta =\underset{\mathrm{nip\; region}}{\underbrace{{\int }_0^a\rho \left(\theta \right)cos\theta d\theta }}+\underset{\mathrm{slip\; region}}{\underbrace{{\int }_a^{{\theta }_{in}}\rho \left(\theta \right)cos\theta d\theta }}=\hfill \\ \hfill & ={\rho }_{exit}\left(\frac{h_0}{R}\right)\left\{\frac{2\left(1+\frac{h_0}{R}\right)}{\sqrt{\frac{h_0}{R}\left(2+\frac{h_0}{R}\right)}}arctan\left[\sqrt{1+\frac{2R}{h_0}}tan\frac{\alpha }{2}\right]-\alpha \right\}+{\rho }_{in}(cos\nu -sin\alpha ).\hfill \end{array}$$

(15)

In Equation (15), ${\rho }_{in}$ is a variable equal to $\pi /2-\nu$ and $\rho \left(\theta \right)$ is density profile, relative to angular position $\theta$. Finally, the outlet mass flow rate can be calculated, ${\dot{m}}_{out}=RWhW{\rho }_{exit}$, where W is the roll width and $h=2h_0$.

The tableting process involves compacting powder mixtures to form a solid dosage form. The performance of tableting processes can be assessed based on tablet active ingredient content, content uniformity, weight variability, and physical properties such as friability, hardness, and dissolution performance [30]. In this paper, the final goal of the design of the tableting process model is to obtain the tablet’s hardness and the mass flow out of the tablet. For the first requirement, the pre-compression $\left({C\_P}_{pre}\right)$ and main compression $\left({C\_P}_{main}\right)$ pressures are formulated based on the Kawakita equation [32]:

$$\begin{array}{cc}\hfill & C\_P_{pre}=\frac{V_0-V_{pre}}{b_{pre}\left(V_0\left({ϵ}_0+1\right)+V_{pre}\right)},\hfill \\ \hfill & C\_P_{main}=\frac{V_{pre}-V_{tablet}}{b_{main}\left(V_{pre}\left({ϵ}_{main}+1\right)+V_{tablet}\right)},\hfill \end{array}$$

(16)

where $b_{pre}$ and $b_{main}$ are material dependent parameters. In the following equations, the required volumes ($V_{tablet}$ for the volume of the tablet, pre-compression volume $V_{pre}$ and the feed volume $V_0$) are defined:

$$\begin{array}{cc}\hfill & V_{tablet}=L_{tablet}{A}_{tablet},\hfill \\ \hfill & V_{pre}=L_{pre}{A}_{tablet},\hfill \\ \hfill & V_0=L_{depth}{A}_{tablet},\hfill \end{array}$$

(17)

where $L_{tablet}$, $L_{pre}$ and $L_{depth}$ are the heights of the tablet, of the powder in the die, and of the powder after compression in the die, respectively, while $A_{tablet}$ represents the area of the tablet. It is notable that, in the considered process model, the pre-compression and main compression forces are calculated as multiplication by ${10}^6$ of the associated pressures.

Finally, the hardness of the tablet is obtained as follows:

$$H_{tablet}=H_{max}\left(1-e^{\frac{(1-{ϵ}_0)V_0}{V_{tablet}}-{\rho }_{rc}+log\frac{1-(1-{ϵ}_0)V_0}{(14-{\rho }_{rc})V_{tablet}}}\right),$$

(18)

where maximum tablet hardness $\left(H_{max}\right)$ and critical relative density $\left({\rho }_{rc}\right)$ are known variables.

In the case of the mass flow out of the tablet press, this is determined as the product of the production rate $\left(R_{tablet}\right)$ and the mass of a single tablet $\left(m_{tablet}\right)$: $\dot{m}=R_{tablet}{m}_{tablet}$, where $R_{tablet}=T_{turret}{n}_{die}60$, considering the speed of the turret $T_{turret}$ and the number of dies on the turret $n_{die}$.

4. Closed-Loop Control System Framework

As the manufacturing processes of tablets is a complex system, the need for advanced control systems grows increasingly important. To achieve the desired product quality and consistency, closed-loop control systems have been widely adopted in dry granulation processes. The use of closed-loop control systems allows for the real-time monitoring and adjustment of process parameters to ensure the final product meets the required specifications.

4.1. Control System Architecture

The process described in the previous section is modeled in Matlab as a plant, which is the $8\times 8$ input-output system presented in Figure 3. The nominal values (Value) of these variables and the associated constraints (lower bounds—LB and upper bounds—UB) are given in

Table 2. Nomina values and constraints of inputs variables.

Input	Value	LB	UB	Unit
Screw speed excipient feeder $\left({\omega }_{Exc}\right)$	207.6	0	240	rpm
Screw speed API feeder $\left({\omega }_{API}\right)$	37.4	0	240	rpm
Hydraulic pressure RC $\left(P_{h}ydraulic\right)$	1	1	10	MPa
Feed speed $\left(u_d\right)$	2.017	1	5	cm/s
Angular velocity rolls RC $\left({\omega }_{RC,rolls}\right)$	5	1	10	rpm
Turret speed TP $\left({\omega }_{TP,turret}\right)$	45	40	50	rpm
Height tablet $\left(L_{tablet}\right)$	0.004	0.0038	0.005	m
Feed volume $\left(V_{feed}\right)$	9.6 × 10${}^{-7}$	9 × 10${}^{-7}$	11 × 10${}^{-7}$	m${}^3$

and

Table 3. Nominal values and constraints of outputs variables.

Output	Value	LB	UB	Unit
Mass flow rate outlet blender $\left({\dot{m}}_{blender}\right)$	20	17	23	kg/h
Concentration API $\left(C_{API}\right)$	0.15	0.10	0.20	-
Density outlet RC $\left({\rho }_{RC}\right)$	1.057	0.8	1.2	g/cm${}^3$
Roller gap RC $\left(h_{0,RC}\right)$	1.6	1	5	mm
Mass flow rate outlet RC $\left({\dot{m}}_{out,RC}\right)$	20	17	23	kg/h
Mass flow rate outlet TP $\left({\dot{m}}_{TP}\right)$	20	17	23	kg/h
Hardness tablet $\left(H_{tablet}\right)$	5.433	4	6	MPa
Mass tablet $\left(M_{tablet}\right)$	0.4566	0.44	0.47	g

.

The two tables are divided into three sections, indicated by the white and gray regions. The output of each section is influenced by the inputs it receives. These three sections are internally connected through density, and it is essential that the mass flow rate is equal across all three parts: ${\dot{m}}_{blender}={\dot{m}}_{out,RC}={\dot{m}}_{TP}$.

The architecture for the implementation of the predictive control strategy is resumed in Figure 4, where the Plant block is illustrated in Figure 3 with 8 inputs (manipulated by the MPC controller) and 8 outputs. A centralized control solution involving a single MPC controller was used to manage all processing units involved in the continuous dry granulation process for the manufacture of solid dosage forms.

4.2. Model Predictive Control Strategy

Starting from a linear discrete-time system with m states, n inputs, and p outputs described by the state-space equation in Equation (19), the state $x_d\in {\Re }^m$, control input $u_d\in {\Re }^n$ and output $y_d\in {\Re }^p$ vectors of the system at discrete time t are defined:

$$\begin{array}{cc}\hfill & x_d[t+1]=A_d{x}_d\left[t\right]+B_d{u}_d\left[t\right],\hfill \\ \hfill & y_d\left[t\right]=C_d{x}_d\left[t\right]+D_d{u}_d\left[t\right],\hfill \end{array}$$

(19)

where $A_d\in {\Re }^{m\times m}$, $B_d\in {\Re }^{m\times n}$, $C_d\in {\Re }^{p\times m}$, $D_d\in {\Re }^{p\times n}$ are the discretized system matrices obtained using a specified sampling time $T_s$.

Supposing that the matrix pairs $(A_d,B_d)$ and $(C_d,A_d)$ are stabilizable and detectable, respectively, the goal is to find a state-feedback control law, expressed based on the model predictive control law gain $K_{MPC}$, of the form:

$$u_d\left[t\right]=K_{MPC}x\left[t\right],$$

(20)

so that the origin of the closed-loop system $x_d[t+1]=(A_d+B_d{K}_{MPC})x_d\left[t\right]$ is globally asymptotically stable, which implies that the matrix $(A_d+B_d{K}_{MPC})$ is stable.

Thus, the predicted state can be written in matrix form as $\widehat{X}\left[t\right]=Mx\left[t\right]+CU\left[t\right]$ [33,34,35], considering the following:

$$U\left[t\right]=\left(\begin{array}{c}u\left[t\right|t]\\ u[t+1|t]\\ ⋮\\ u[t+hp-1|t]\end{array}\right),$$

(21)

with $u[t+i|t]$, $i=0,\cdots ,hp-1$, represents the future control sequence with $hp$ the prediction horizon,

$$\widehat{X}\left[t\right]=\left(\begin{array}{c}\widehat{x}[t+1|t]\\ \widehat{x}[t+2|t]\\ ⋮\\ \widehat{x}[t+hp|t]\end{array}\right),$$

(22)

where $\widehat{x}[t+i|t]$ is the predicted value of the state vector $x[t+i+1|t]=A_d\widehat{x}[t+i|t]+B_{d}u[t+i|t],i=0,1,2,\cdots$ with the initial condition described by $\widehat{x}\left[t\right|t]=x\left[t\right]$,

$$M=\left[\begin{array}{c}{A}_d\\ A_d^2\\ ⋮\\ A_d^{hp}\end{array}\right],C=\left[\begin{array}{cccc}{B}_d& 0& \cdots & 0\\ A_d{B}_d& B_d& \cdots & 0\\ ⋮& ⋮& \ddots & ⋮\\ A_d^{hp-1}{B}_d& A_d^{hp-2}{B}_d& \cdots & B_d\end{array}\right].$$

(23)

To design an optimal predictive control law expressed as:

$$U^{*}\left[t\right]=arg\underset{U}{min}{J}_{MPC}\left[t\right],$$

(24)

it is necessary to minimize a quadratic cost function of the following form:

$$J_{MPC}\left[t\right]=U_0^T\left[t\right]H{U}_0\left[t\right]+2U_0^T\left[t\right]F{x}_0\left[t\right]+x_0^T\left[t\right]G{x}_0\left[t\right],$$

(25)

with

$$\begin{array}{cc}\hfill & U_0\left[t\right]=U\left[t\right]-U_{ss}\left[t\right],\hfill \\ \hfill & x_0\left[t\right]=x\left[t\right]-x_{ss}\left[t\right],\hfill \\ \hfill & H=C^T{\tilde{Q}}_{MPC}C+{\tilde{R}}_{MPC},\hfill \\ \hfill & F=C^T{\tilde{Q}}_{MPC}M,\hfill \\ \hfill & G=M^T{\tilde{Q}}_{MPC}M+Q_{MPC},\hfill \end{array}$$

(26)

where at discrete-time instant t, $U_{ss}\left[t\right]$ and $x_{ss}\left[t\right]$ correspond to the steady-state values determined from a given reference. The matrices ${\tilde{Q}}_{MPC}$ and ${\tilde{R}}_{MPC}$ are defined as:

$${\tilde{Q}}_{MPC}=\left[\begin{array}{cccc}{Q}_{MPC}& 0& \cdots & 0\\ 0& \ddots & & ⋮\\ ⋮& & Q_{MPC}& 0\\ 0& \cdots & 0& {\overline{Q}}_{MPC}\end{array}\right],{\tilde{R}}_{MPC}=\left[\begin{array}{cccc}{R}_{MPC}& 0& \cdots & 0\\ 0& \ddots & & ⋮\\ ⋮& & R_{MPC}& 0\\ 0& \cdots & 0& R_{MPC}\end{array}\right].$$

(27)

The matrices $Q_{MPC}$, $R_{MPC}$, and ${\overline{Q}}_{MPC}$ are positive definite matrices, with $Q_{MPC}$ possibly being positive semi-definite. These matrices are utilized to weigh the states, inputs, and terminal state appropriately.

Assuming there are no constraints, the optimization problem can be solved offline to obtain $U^{*}\left[t\right]=-H^{-1}Fz\left[t\right]$ and in line with the receding horizon principle, only the initial control command in $U^{*}$ is implemented on the system, as described by $u\left[t\right]=u^{*}\left[t\right|t]=K_{MPC}x\left[t\right]$, with $K_{MPC}=-\left[I_{hp}0\cdots 0\right]H^{-1}F$.

In practical applications, system variables, such as inputs, outputs, and states, are limited by physical constraints, which can include input constraints (actuator limits) or state/output constraints (that can be active during transients or in steady-state). Predictive controllers offer a significant benefit in that they can determine an optimal nonlinear control law for systems with state, input, and/or output constraints in real-time at each sampling period.

The constraints on system variables appear usually due to physical limitations and can be expressed as linear inequalities, which are defined as:

$$\begin{array}{ccc}\hfill u_{min}\le & u\left[t\right]& \le u_{max}\hfill \\ \hfill x_{min}\le & x\left[t\right]& \le x_{max}\hfill \\ \hfill y_{min}\le & y\left[t\right]& \le y_{max}.\hfill \end{array}$$

(28)

Given the quadratic cost function $J_{MPC}\left[t\right]$ defined in Equation (25) and the constraints specified earlier, the minimization of the cost function involves solving a quadratic problem, which is expressed as [33,34,35]: ${min}_U{U}^{T}HU+2x^T\left[t\right]F^{T}U$, subject to the constraints that were previously stated.

5. Simulation Results

In this section, the MPC methodology proposed for continuous manufacturing of solid dosage forms using the dry granulation route is tested and compared by simulation with other classical control approaches (described in Section 5.2) applied to the same simulation setup presented in Section 5.1. The obtained illustrative results are presented in Section 5.3 and are completed by numerical analysis in Section 5.4.

5.1. Simulation Setup

The Matlab-Simulink environment was used to model and control the continuous dry granulation tablet production process described in Section 3. A benchmark simulator is created to evaluate the feasibility of developing control strategies that are intended to control the real physical plant. The simulator considered in this work consists of a feeder, a blender, a roller compactor with an RC hopper, a co-mill, and a tablet press with a TP hopper. The Simulink implementation of the fixed part involves a complex system where each processing unit mentioned previously is modeled individually.

This study simulates a powder mixture consisting of API and excipient material. The excipients and APIs feeders get the input speed of the associated screws and provide the respective flow rates at the output. It should be noted that all parameters of the feed factor model used in this work are material dependent. The material properties that have an impact on the feeder performance are: conditioned bulk density, compressibility, permeability, cohesion, flow function coefficient, and the angle of internal friction. The blender is a crucial unit in the pharmaceutical industry since it guarantees concentration uniformity and homogeneity. Several modeling techniques are being used to simulate the dynamics of the powder in the blender. The outputs of the feeders model and the speed of the blender are required for the blender model to provide the outlet concentration, the relative standard deviation as well as the mass flow and the density at the outlet of the blender. In the case of the roller compactor, the key operating variables are the feed speed, roll speed, roll pressure, and the gap between the rolls. The gap and the roller force—the actively controlled process parameters—are the inputs for the model, whereas the screw speed—the changing process parameter—is neglected. The throughput can be predicted once the ribbon density is known. The auger filling may fluctuate from one revolution to another. This fluctuation is modeled by introducing a disturbance in the inlet bulk density. Hereafter, the co-mill is designed based on the outlet density and the outlet mass flow rate from the RC with the co-mill speed. This model will provide the necessary inputs to the hopper associated with the tablet press which is modeled similarly to that of the roller compactor. Finally, the tablet press involves a complex model where it is necessary to have the powder solid density determined previously, turret speed, feed volume, and thickness of the tablet to establish the tablet’s rate, hardness, and mass of the tablet.

5.2. Comparison Control Strategies

Two control strategies, Proportional-Integral-Derivative (PID) [22] and Linear Quadratic Regulator (LQR) [36], are presented and applied to the same simulator to analyze and compare the performance of the proposed strategy for end-to-end continuous pharmaceutical production. PID controllers have been designed for the individual control of the process units involved in the manufacturing process, while for the LQR-based strategy has been proposed a unitary control strategy for the entire process, using a single controller, which is adapted to the system requirements.

Decentralized control strategies offer an alternative to simplify the control of multiple-input multiple-output (MIMO) processes by decomposing the complex system into a series of single-input single-output (SISO) subsystems [37]. However, decentralized control systems have become increasingly popular due to their ability to handle large-scale systems with multiple loops. In such systems, each loop has its own controller, and decisions are made independently of other loops. In [38], a novel approach to designing decentralized controllers based on nonlinear optimization, with the aim of minimizing the impact of disturbances in coupled loops. The proposed methodology achieves improved control performance by taking into account the interdependence between loops and optimizing the control gains for each individual loop.

Figure 5 illustrates the decentralized control strategy applied to the system in which one or more PID controllers are designed for each processing unit involved in the continuous tablet production process, as detailed in Section 5.2.1.

As in the case of predictive control, a centralized control strategy was also applied in the case of LQR, as can be observed in Figure 6, where the Plant block is illustrated in Figure 3 with 8 inputs (manipulated by the LQR controller) and 8 outputs. In addition, an observer has been employed to estimate the system’s unknown states, thereby expanding the control architecture, as detailed in Section 5.2.2.

5.2.1. Proportional-Integral-Derivative Controller

A PID controller is a control loop feedback algorithm used to automatically adjust an output variable based on an input variable and a set point. Thus, it calculates the error as the difference between the desired output and the system’s response. This type of controller is widely used in industrial control systems due to its robust performance and functional simplicity. Several elements can degrade PID controller performance, such as process interactions, process dead time, non-linearity, and process constraints [39]. The control function is defined as the sum of the proportional $\left(P\right)$, integral $\left(I\right)$, andderivative $\left(D\right)$ terms [40]:

$$u\left[t\right]=P\left[t\right]+I\left[t\right]+D\left[t\right],$$

(29)

with

$$\begin{array}{cc}\hfill & P\left[t\right]=K_{p}e\left[t\right],\hfill \\ \hfill & I\left[t\right]=K_i{\int }_0^{t}e\left[\tau \right]d\tau ,\hfill \\ \hfill & D\left[t\right]=K_d\frac{d}{dt}e\left[t\right],\hfill \end{array}$$

(30)

where $K_{p,i,d}$ is the proportional/integral/derivative gain, e represents the tracking error at time t and $\tau$ is the variable of integration.

Tuning the parameters means adjusting the values of the gains to achieve optimal control system performance. In this context, two tuning methods were applied to the simulator, specifically the Frequency Response tool (FRtool), which is a controller design tool, and an optimization-based method (Integral of Time Absolute Error).

FRtool [41] is a program that allows a graphical design of the controller. It is able to tackle systems with time-delays, which is not the case when using the standard Matlab controller design tool, since it operates with frequency diagrams (Nichols charts). The dead time present in the system is not approximated and is treated in the right way.

The Integral of Time Absolute Error (ITAE) [42] is one of the possible criteria of an optimization tuning method. The method involves analytically calculating the gain of the controller $\left(K_p\right)$, integral time $\left(K_i\right)$, and the derivative time $\left(K_d\right)$ for PID-controlled systems. A Matlab program is used to find the optimum values for the PID controller to achieve the best performance to satisfy the system requirements: reducing the overshoot, achieving a good load disturbance rejection, and maintaining a high system response. The corresponding objective function is defined as follows:

$$ITAE={\int }_0^{\infty }t\left|e\left[t\right]\right|dt,$$

(31)

where t is the time and $e\left[t\right]$ is the difference between set point and controlled variable.

Since the dry granulation route consists of three parts (feeder + blender, roller compactor, and tablet press), the evaluation will also be done three ways. So, after comparing the different control strategies and tuning methods in the first part (feeder + blender), the best one will be taken to evaluate the second one. The same will be done for the third part.

In the discussion of the first part, the best structure for the feedforward will be chosen first based on the evaluation, which is used in the other parts as well. In the first part of the granulation route, the cascade control strategy combined with a static feedforward structure, where the PID controllers were tuned using the FRtool, performed better. It was compared to open-loop, only cascade, and other feedforward structures. This output will be fed into the next unit, i.e., the roller compactor (part II of the dry granulation route). Next up is the evaluation of the control strategy for the roller compactor. In this case, the ITAE tuning method is preferred with the feedback combined with feedforward as control strategy. This strategy performs very well and keeps the output variables within the control bounds. This result (part I tuned with FRtool consisting of cascade feedback and feedforward control as a control strategy, part II tuned using the ITAE optimization method consisting of feedback/feedforward control as control strategy) is now used to be fed into the third and last part of the dry granulation route (i.e., the tablet press). The same conclusions as in part II can be drawn, specifically that the combined feedback and feedforward (static structure) control strategy is the most effective.

5.2.2. Linear Quadratic Regulator

The LQR control algorithm [43] is based on the state-space model of the process and can handle multiple feedback variables. This leads to a more robust control algorithm as there is more information about the system available. The approach is based on minimizing a quadratic cost function [44] and can be extended to Linear Quadratic Integral control (LQI) or Linear Quadratic Gaussian regulator (LQG).

A control structure with state feedback is implemented using a mathematical model in matrix form (19). The purpose is to find an optimal state-variable feedback control law, which guarantees the desired closed-loop performances and is defined as:

$$u_d\left[t\right]=-K_{LQR}{x}_d\left[t\right],$$

(32)

where $K_{LQR}$ is the state feedback gain. Also, a quadratic cost function $J_{LQR}$ is defined to quantify the controller performance and to design the optimal state-variable feedback gain [44]:

$$J_{LQR}=\sum _{t=0}^{\infty }\{x_d^T\left[t\right]Q_{LQR}{x}_d\left[t\right]+u_d^T\left[t\right]R_{LQR}{u}_d\left[t\right]+2x_d^T\left[t\right]N_{LQR}{u}_d\left[t\right]\},$$

(33)

where the weighting parameters, $Q_{LQR}$ and $N_{LQR}$ represent semi-positive definite symmetric matrices of dimension $(m\times m)$ and $R_{LQR}$ is a positive definite symmetric matrix of dimension $(n\times n)$. To ensure efficient control of each state of the system with minimum effort, as per the aforementioned performance index, it is necessary to select appropriate weight matrices $Q_{LQR}$ and $R_{LQR}$, which determine the closed-loop performance.

Optimizing the performance index $J_{LQR}$ and solving for $K_{LQR}$ implies the solution $P_{LQR}$ of the Discrete Algebraic Riccati Equation [45]:

$$\begin{array}{cc}\hfill P_{LQR}=& Q_{LQR}+A_d^T{P}_{LQR}{A}_d-\hfill \\ \hfill & -(A_d^T{P}_{LQR}{B}_d+N_{LQR}){(B_d^T{P}_{LQR}{B}_d+R_{LQR})}^{-1}(B_d^T{P}_{LQR}{A}_d+N_{LQR}^T).\hfill \end{array}$$

(34)

By iteratively solving the Riccati equation of the related matrix, the optimal feedback gain can be achieved:

$$K_{LQR}={(B_d^T{P}_{LQR}{B}_d+R_{LQR})}^{-1}(B_d^T{P}_{LQR}{A}_d+N_{LQR}^T).$$

(35)

In order to achieve the best performance of feedback control, an integral feedback control strategy is used to account for disturbances and noises from the external environment. This is achieved by designing a Linear Quadratic Regulator with Integral action (LQI), which provides zero steady-state error using an integrator. Thus, the Linear Quadratic Regulator with Integral action [46] is designed to provide zero steady-state error using an integrator, added as an additional state $\left(x_i\right)$ to the state vector $\left(x_d\right)$. This represents the integral of the error between the desired output r, and the actual output $y_d$. Thus, the extended state vector is given as:

$$x_e=\left[\begin{array}{c}{x}_d\\ x_i\end{array}\right],$$

(36)

where $x_i[t+1]=x_i\left[t\right]+T_s(r\left[t\right]-y_d\left[t\right])$, with sampling time $T_s$.

The optimal gain matrix $K_{LQI}$ required to determine the control law given by $u_i\left[t\right]=-K_{LQI}{x}_e$ is obtained similarly to the method presented above.

When the full state variable $x_d\left[t\right]$ cannot be measured in the implementation of a controller, an observer is typically used. The Luenberger observer [47] employs a correction term based on the measured output $y_d$ and its estimate ${\widehat{y}}_d$ to model the current system. Assuming observability of the pair $(A_d,C_d)$ in the discretized system from Euqation (19), a full-order observer can be expressed using the following equation [35]:

$$\begin{array}{cc}\hfill {\widehat{x}}_d[t+1]=A_d{\widehat{x}}_d\left[t\right]& +B_d{u}_d\left[t\right]+K_{obs}(y_d\left[t\right]-{\widehat{y}}_d\left[t\right])\hfill \\ \hfill & {\widehat{y}}_d=C_d{\widehat{x}}_d,\hfill \end{array}$$

(37)

where $K_{obs}$ is the constant observer gain matrix. A conventional method used to model the observer is pole assignment, whereby the selection of gain matrix $K_{obs}$ can result in various states for the system $(A_d,C_d)$. If $K_{obs}$ is chosen to place the eigenvalues of $A_d-K_{obs}{C}_d$ at any desired locations in the plane, then the system can be fully observed. The observer is designed with the goal of estimating the actual state $x_d\left[t\right]$. Based on the definition of the estimation error:

$$e_d\left[t\right]=x_d\left[t\right]-{\widehat{x}}_d\left[t\right],$$

(38)

the observer dynamics lead to the following equation for the estimation error:

$$e_d[t+1]=x_d[t+1]-{\widehat{x}}_d[t+1]=(A_d-K_{obs}{C}_d)e_d\left[t\right].$$

(39)

5.3. Illustrative Results

This section will review the performance of the controlled tablet manufacturing process using a predictive control strategy, which was implemented using the Model Predictive Control Toolbox for Matlab/Simulink, in comparison to the results obtained in [22] using PID, as well as those based on LQR in [36]. Both the performance of disturbance rejection and set-point tracking are analyzed and discussed. The comparative results obtained from the simulations are shown in Figure 7, Figure 8, Figure 9, Figure 10, Figure 11, Figure 12, Figure 13, Figure 14, Figure 15 and Figure 16. The time interval $[0,500]$ s is considered the time until the system sets up and reaches a normal operating mode. Thus, the variations of the signal occurring in this interval do not impact the performance of the proposed strategy.

Defining specific boundaries is essential to assess the efficacy of the control strategy. To evaluate the performance of the controller, upper and lower limits are established for the controlled variables, according to the numerical values presented in Table 3. These values are plotted on the graphs in Figure 7, Figure 8 and Figure 9, Figure 11 and Figure 13, Figure 14, Figure 15 and Figure 16 with red dashed lines. If the controlled variable remains within these predetermined boundaries, it can be concluded that the controller is performing well, the yielded control law, i.e., the manipulated variable, being feasible.

The control inputs of the associated outputs were represented in order to evaluate the performance of the proposed control strategy. The speed of the excipient feeder (Figure 7) is paired with the mass flow rate at the outlet of the blender (Figure 8). The speed of the API feeder (Figure 9) required to provide the concentration at the outlet of the blender (Figure 11Figure 12) can be observed in detail in Figure 10. The hardness of the tablet (Figure 14) is determined by the height of the tablet (Figure 13) and the mass of the tablet (Figure 16) is paired with the feed volume (Figure 15) at the inlet side.

Two kinds of disturbances are introduced into the system, namely on the one hand structural disturbances that represent the inherent feed rate variability and on the other hand measurement errors (or disturbances) which represent measurement noise coming for example from the equipment.

Some minor random disturbances have been added to certain variables to account for unknown disturbances in the plant that can relate to the measurements and feeds. Their influence is more visible in Figure 11 and Figure 16, where the controllers are not able to completely reject the introduced disturbances. The mass of the tablet, plotted in Figure 16, is controlled considerably well around its reference with the disturbance not as intense as before for the PID controller. The step disturbance added at time $T=1400$ s affects the output, regardless of the used control strategy.

5.4. Numerical Analysis

A numerical analysis was also made considering the performance of the three control methods applied to a continuous pharmaceutical tablet production system. Normalized integral-square-error (NISE) is a performance metric that is commonly used to evaluate the accuracy of a control system. It is a variant of the integral-square-error (ISE) [42] metric that takes into account the magnitude of the input signal, as well as the magnitude of the error signal.

The ISE is the integral of the square of the error signal over the time interval of interest. This measures the cumulative deviation of the system output from the set point over time. The NISE metric then normalizes the ISE by dividing it by the integral of the square of the input signal over the same time interval. This takes into account the fact that a larger input signal will result in a larger error signal, and therefore a larger ISE. By normalizing the ISE, the NISE metric provides a more meaningful measure of the system’s accuracy. A lower NISE value indicates better performance of the control system.

The results are illustrated in

Table 4. Normalized integral-square-error (NISE).

	Out1 (${10}^{-3}$)	Out2 (${10}^{-7}$)	Out7 (${10}^{-3}$)	Out8 (${10}^{-6}$)
MPC	3.519	6.263	2.655	2.361
LQR	4.279	7.223	2.717	3.158
PID	168	14,612	1021	175

. It can be noticed that the predictive solution proposed in the current work gives the best results. As it could be observed in the results illustrated in the previous section, the PID control strategy presents significantly degraded performance, leading to the manufacturing of out-of-specification products.

6. Conclusions

In the context of cyber physical systems, continuous production can be considered a necessity because this technology can offer a lot of benefits for companies as well as for the whole economy. First of all, continuous manufacturing allows companies to produce higher quality products more efficiently and at lower cost. This is because CM is a continuous, automated process that eliminates the necessity to stop and restart the production process and reduces the time and costs associated with it. In addition, continuous production allows companies to produce more products in the same period of time, which can lead to increased productivity and profitability. Also, by eliminating the risk of human error, continuous production can increase product quality and reduce the costs associated with eliminating damaged products. The cyber-physical system can also be configured and controlled remotely, which offers increased flexibility and agility in production. Through continuous monitoring and data analysis, companies can identify and eliminate problems before they have a negative impact on the production process.

In this study, an advanced control strategy, i.e., MPC, was implemented to control a multiple-input multiple-output (MIMO) system that represents a continuous pharmaceutical plant. The results obtained from the experiments performed on a simulator modeled in Simulink were compared with those resulting from the application of two classical control methods, PID and LQR. The performance of each strategy in terms of stability of the production process and the achievement of high-quality products was analysed. MPC demonstrated better performances in terms of the continuous production process of solid dosage forms by dry granulation because it is a control strategy based on complex mathematical models, which can be adapted to different system operating conditions. This predictive capability makes MPC more effective in situations where quick decisions need to be made and the production process adjusted in real-time.

However, in the case of continuous production in the pharmaceutical industry, there is still potential for improvement and future research directions in this area. One of these is the development of hybrid control strategies that combine several types of control techniques to achieve even better performance. Further research is also needed to assess the impact of variations in material quality and environmental conditions on the performance of cyber-physical continuous manufacturing systems.