Applied Science, Faculty of
Engineering, School of (Okanagan)
DSpace
UBCO
Hosseinionari, Hadi
Doctor of Philosophy - PhD
University of British Columbia
The importance of thermoforming in manufacturing has increased due to the growing demand for plastic products across various industries. The large-scale production of thermoformed parts has necessitated the development of precise and cost-effective process controllers. The thermoforming industry lags behind manufacturing processes, like injection molding, in terms of automation, control, and process understanding. As such, the thermoforming industry is currently facing the critical task of increasing its degree of process control. The thermoforming process consists of two primary phases: the heating phase and the vacuum and cooling phase. The heating phase stands out as the most critical step, as the outcomes of the subsequent phase hinge on its effectiveness. Despite the development of control techniques for the heating phase of the process, oven heater temperatures in the thermoforming industry are still mostly regulated by experienced operators using trial and error. The designed control algorithm for the heating phase in thermoforming must be easy to implement, cost-effective, robust to environmental disturbances, and should minimize control time to maximize production efficiency. To bridge these existing gaps, the first objective of this thesis is to develop a precise and computationally efficient model for the heating phase in thermoforming. This 2D model, by utilizing the experimental radiation pattern, achieves a level of accuracy comparable to a 3D model (decreasing the MSE by four-fold compared to traditional 2D models) while maintaining the computational efficiency characteristic of a 2D model (ten times faster than 3D models). The complexity of this model makes it difficult for traditional controller development. Thus, it is used to collect data for data-driven methods or to train various reinforcement learning (RL) control approaches. The controller is strategically designed and trained based on the system’s input and output instead of facing the complexities of solving mathematical equations involved in the heating phase model, which includes radiation heat transfer and actuators’ dynamics. Because it uses robust search algorithms, RL is especially effective at optimizing controllers in systems with highly uncertain, possibly stochastic, or nonlinear dynamics. They are compared to PI and Model Predictive Control (MPC), an advanced model-based process control technique to Benchmark the effectiveness of the developed RL control strategies.
Physics-based Simulation andReinforcement Learning Control of theHeating Phase in the ThermoformingProcessbyHadi HosseinionariB.Sc., The University of Mohaghegh Ardabili, 2011M.Sc., K. N. Toosi University of Technology, 2013A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFDOCTOR OF PHILOSOPHYinThe College of Graduate Studies(Mechanical Engineering)THE UNIVERSITY OF BRITISH COLUMBIA(Okanagan)July 2024© Hadi Hosseinionari, 2024The following individuals certify that they have read, and recommendto the College of Graduate Studies for acceptance, a dissertation entitled:Physics-based Simulation and Reinforcement Learning Con-trol of the Heating Phase in the Thermoforming Processsubmitted by Hadi Hosseinionari in partial fulfilment of the requirementsof the degree of Doctor of PhilosophyDr. Rudolf Seethaler, School of Engineering.SupervisorDr. Abbas S. Milani, School of Engineering.Supervisory Committee MemberDr. Bhushan Gopaluni, The University of British Columbia (Vancouver).Supervisory Committee MemberDr. Ian Foulds, School of Engineering.University ExaminerDr. Syed Imtiaz, Memorial University of Newfoundland.External ExamineriiAbstractThe importance of thermoforming in manufacturing has increased dueto the growing demand for plastic products across various industries. Thelarge-scale production of thermoformed parts has necessitated the develop-ment of precise and cost-effective process controllers. The thermoforming in-dustry lags behind manufacturing processes, like injection molding, in termsof automation, control, and process understanding. As such, the thermo-forming industry is currently facing the critical task of increasing its degreeof process control.The thermoforming process consists of two primary phases: the heatingphase and the vacuum and cooling phase. The heating phase stands out asthe most critical step, as the outcomes of the subsequent phase hinge on itseffectiveness. Despite the development of control techniques for the heatingphase of the process, oven heater temperatures in the thermoforming indus-try are still mostly regulated by experienced operators using trial and error.The designed control algorithm for the heating phase in thermoforming mustbe easy to implement, cost-effective, robust to environmental disturbances,and should minimize control time to maximize production efficiency.To bridge these existing gaps, the first objective of this thesis is to de-velop a precise and computationally efficient model for the heating phasein thermoforming. This 2D model, by utilizing the experimental radiationpattern, achieves a level of accuracy comparable to a 3D model (decreasingthe MSE by four-fold compared to traditional 2D models) while maintainingthe computational efficiency characteristic of a 2D model (ten times fasterthan 3D models). The complexity of this model makes it difficult for tradi-tional controller development. Thus, it is used to collect data for data-drivenmethods or to train various reinforcement learning (RL) control approaches.The controller is strategically designed and trained based on the system’sinput and output instead of facing the complexities of solving mathematicalequations involved in the heating phase model, which includes radiation heattransfer and actuators’ dynamics. Because it uses robust search algorithms,RL is especially effective at optimizing controllers in systems with highlyuncertain, possibly stochastic, or nonlinear dynamics. They are comparediiito PI and Model Predictive Control (MPC), an advanced model-based pro-cess control technique to Benchmark the effectiveness of the developed RLcontrol strategies.ivLay SummaryThe advent of the fourth industrial revolution is reshaping manufac-turing processes, rendering factories more efficient, customizable, and en-vironmentally sustainable. At the forefront of this transformation are in-telligent decision-making tools driven by artificial intelligence (AI), whichhave become integral to the evolving landscape of modern manufacturing.Manufacturing enterprises are progressively harnessing the power of AI tooptimize various facets of their operations, spanning production, marketing,and enhancing customer satisfaction.This dissertation aims to develop and implement a precise, robust, andscaleable temperature control algorithm during the thermoforming processheating stage. This work supports using reinforcement learning (RL) toovercome the difficulties of controlling the heating phase in thermoforming.The control algorithm proposed is easy to implement while functioning effec-tively with different sizes of heater banks. Most importantly, it is designedto function seamlessly without requiring much computational effort, eitheroffline or online.vPrefaceThis dissertation results from my research under the supervision of Dr.Rudolf Seethaler at the University of British Columbia. I was responsiblefor designing the solution framework, constructing the lab setup, conductingexperiments, developing an accurate simulator for the heating phase, design-ing and evaluating control strategies using the simulator, implementing thedesigned controllers on the lab-scale setup, and analyzing the obtained re-sults. Notably, this research was conducted in collaboration with KohlerCompany in Canada, providing valuable support and contributing thermo-plastic materials for experimentation. No generative AI was used for anyaspects of this work. The original research carried out in this dissertationhas led to multiple journal publications, a summary of which is providedbelow:1. A portion of Chapter 3 has been published in the Journal of Man-ufacturing and Materials Processing. I was responsible for developing themodeling framework, conducting the experiments, and analyzing the results,as well as writing the manuscript:Hosseinionari, Hadi, Milad Ramezankhani, Rudolf Seethaler, and Ab-bas S. Milani. 2023. ”Development of a Computationally Efficient Modelof the Heating Phase in Thermoforming Process Based on the Experimen-tal Radiation Pattern of Heaters,” Journal of Manufacturing and MaterialsProcessing 7, no. 1: 48. https: // doi. org/ 10. 3390/ jmmp70100482. A portion of Chapters 3, 5, and 6 has been published in the Journalof Manufacturing Processing. I was responsible for developing the modelingframework, developing the Deep reinforcement Learning and Model Predic-tive Control Algorithms for the heating phase and integrating them, andanalyzing the results, as well as writing the manuscript:Hosseinionari, Hadi, and Rudolf Seethaler. ”The integration of ModelPredictive Control and deep Reinforcement Learning for efficient thermalvicontrol in thermoforming processes.” Journal of Manufacturing Processes115 (2024): 82-93. https: // doi. org/ 10. 1016/ j. jmapro. 2024. 01.0853. A portion of Chapters 3 and 6 has been submitted to the Journalof Manufacturing Processing. My contributions to this submission includedformulating the modeling framework, designing the Multi-agent Deep Rein-forcement Learning algorithms for the heating phase, deploying the devel-oped controller on the real-world lab-scale thermoforming setup, analyzingthe outcomes, and composing the manuscript:Hosseinionari, Hadi, Rudolf Seethaler, Rohith Jayaraman Krishnamurthy,Iman Jalilvand, Abbas S. Milani. ”A Scalable Multi-Agent Deep Reinforce-ment Learning in Thermoforming: An Experimental Evaluation of Ther-mal Control by Infrared Camera-Based Feedback.” Journal of ManufacturingProcesses. (Under Review)4. A portion of Chapters 3 and 5 has been submitted to the Journalof Robotics and Computer-Integrated Manufacturing. My contributions tothis submission included formulating the modeling framework, deploying thedeveloped controller on the real-world lab-scale thermoforming setup, ana-lyzing the outcomes, and composing the manuscript:Hosseinionari, Hadi, Mohammad Bajelani, Klaske van Heusden, AbbasS. Milani, Rudolf Seethaler. ”Towards Indirect Data-Driven Predictive Con-trol in Thermoforming Process: Radiation Thermal Control Stage.” Journalof Robotics and Computer-Integrated Manufacturing. (Under Review)5. During my PhD studies, I have also been involved in some side projectsas part of the collaborations with other research labs at UBC as follows:Jalilvand, I., Jiyoung, J., Hosseinionari, H., Seethaler, R., Gopaluni,B., Milani, A.S. (2023). An Interactive Digital Twin of a Composite Man-ufacturing Process for Training Operators via Immersive Technology. In:Stephanidis, C., Antona, M., Ntoa, S., Salvendy, G. (eds) HCI Interna-tional 2023 Posters. HCII 2023. Communications in Computer and Infor-mation Science, vol 1836. Springer, Cham. https: // doi. org/ 10. 1007/978-3-031-36004-6_ 28Hosseinionari, Hadi, Reza Azarmanesh. (2023). Septic Tank SensorviiDesign. MITACS project in collaboration with Mountain Side EarthworksCompanyviiiTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiiLay Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . ixList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiiList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xivAcknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . xxvDedication . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .xxviChapter 1: Introduction . . . . . . . . . . . . . . . . . . . . . . . 11.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31.2 Major Contributions and Structure of the Proposed Research 7Chapter 2: Literature Review . . . . . . . . . . . . . . . . . . . 112.1 Thermoforming Process Modeling . . . . . . . . . . . . . . . . 112.2 Temperature Sensing . . . . . . . . . . . . . . . . . . . . . . . 142.3 Control System Approaches . . . . . . . . . . . . . . . . . . . 162.3.1 Proportional-Integral-Derivative (PID) . . . . . . . . . 172.3.2 Model Predictive Control (MPC) . . . . . . . . . . . . 182.3.3 Reinforcement Learning . . . . . . . . . . . . . . . . . 19Chapter 3: Modeling . . . . . . . . . . . . . . . . . . . . . . . . . 223.1 Thermal Dynamics of the Heating Element . . . . . . . . . . 223.2 Thermal Dynamics of the Thermoplastic Sheet . . . . . . . . 253.3 UBC Lab-scale Thermoforming Setup . . . . . . . . . . . . . 31ixTABLE OF CONTENTS3.4 Experiment to Determine Heating Element’s Radiation Pat-tern . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 343.5 Physics-based Heating Phase Model Verification . . . . . . . . 383.6 Chapter Summary . . . . . . . . . . . . . . . . . . . . . . . . 42Chapter 4: Stability Analysis and PI Control of the HeatingPhase in Thermoforming Process . . . . . . . . . . 444.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 444.2 Control Objectives for Thermoforming Systems . . . . . . . . 444.2.1 Heating Beyond Glass Transition Temperature . . . . 454.2.2 Operating Time Efficiency . . . . . . . . . . . . . . . . 454.2.3 Terminal Control Error . . . . . . . . . . . . . . . . . 454.2.4 Robustness . . . . . . . . . . . . . . . . . . . . . . . . 464.2.5 Scalability . . . . . . . . . . . . . . . . . . . . . . . . . 464.3 Assessing Feasibility of Reference Temperature Profiles onThermoplastic Sheets . . . . . . . . . . . . . . . . . . . . . . . 474.4 Open-loop Control Results . . . . . . . . . . . . . . . . . . . . 534.5 PI Control Design and Stability Analysis . . . . . . . . . . . . 554.5.1 State Feedback Stability . . . . . . . . . . . . . . . . . 564.5.2 Sheet Temperature Control Loop . . . . . . . . . . . . 614.5.3 PI Control Loop Stability . . . . . . . . . . . . . . . . 624.5.4 Heater Temperature Control Loop . . . . . . . . . . . 624.5.5 PI Results on the Simulator . . . . . . . . . . . . . . . 634.5.6 Chapter Summary . . . . . . . . . . . . . . . . . . . . 68Chapter 5: Model Predictive Control of the Heating Phasein Thermoforming Process . . . . . . . . . . . . . . 695.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 695.2 Adaptive MPC for the Heating Phase . . . . . . . . . . . . . 705.2.1 Physics Based Model Linearization . . . . . . . . . . . 715.2.2 Optimization Problem . . . . . . . . . . . . . . . . . . 735.2.3 Simulation Result . . . . . . . . . . . . . . . . . . . . 735.3 NARX-MPC . . . . . . . . . . . . . . . . . . . . . . . . . . . 765.3.1 Data-Driven Modeling: Nonlinear Autoregressive Ex-ogenous (NARX) . . . . . . . . . . . . . . . . . . . . . 765.3.2 NARX-MPC for Heating Phase in Thermoforming Sys-tems . . . . . . . . . . . . . . . . . . . . . . . . . . . . 785.3.3 Data-driven Modeling for Thermoforming Systems . . 805.3.4 NARX Model Evaluation . . . . . . . . . . . . . . . . 815.3.5 NARX-MPC Results . . . . . . . . . . . . . . . . . . . 86xTABLE OF CONTENTS5.4 Chapter Summary . . . . . . . . . . . . . . . . . . . . . . . . 97Chapter 6: Reinforcement Learning Control of the HeatingPhase in Thermoforming Process . . . . . . . . . . 996.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 996.1.1 Markov Decision Process (MDP) . . . . . . . . . . . . 996.1.2 Reward Function . . . . . . . . . . . . . . . . . . . . . 1006.1.3 Policy . . . . . . . . . . . . . . . . . . . . . . . . . . . 1016.1.4 Value Functions . . . . . . . . . . . . . . . . . . . . . 1016.1.5 Bellman Equations in Reinforcement Learning . . . . 1026.1.6 Learning Optimal Policies . . . . . . . . . . . . . . . . 1046.2 Deep Reinforcement Learning . . . . . . . . . . . . . . . . . . 1056.2.1 Taxonomy of Reinforcement Learning Algorithms . . . 1066.2.2 Proximal Policy Optimization (PPO) . . . . . . . . . 1096.3 Deep Reinforcement Learning for Thermoforming System . . 1116.4 MPC-Guided DRL Approach in Thermoforming System . . . 1136.4.1 Observation and Action Space . . . . . . . . . . . . . 1156.4.2 Reward Function . . . . . . . . . . . . . . . . . . . . . 1166.4.3 Agent Training . . . . . . . . . . . . . . . . . . . . . . 1186.4.4 Simulation Results . . . . . . . . . . . . . . . . . . . . 1226.5 Multi-Agent DRL Approach in Thermoforming System . . . . 1266.5.1 Observation and Action Space . . . . . . . . . . . . . 1286.5.2 Reward Function . . . . . . . . . . . . . . . . . . . . . 1316.5.3 Agent Training . . . . . . . . . . . . . . . . . . . . . . 1316.5.4 Simulation Results . . . . . . . . . . . . . . . . . . . . 1346.5.5 Cycle-to-cycle Control Analysis . . . . . . . . . . . . . 1396.5.6 Robustness Analysis . . . . . . . . . . . . . . . . . . . 1416.5.7 Experimental Results . . . . . . . . . . . . . . . . . . 1436.6 Chapter Summary . . . . . . . . . . . . . . . . . . . . . . . . 148Chapter 7: Conclusions and Future Works . . . . . . . . . . . 1507.1 Research Contributions and Conclusions . . . . . . . . . . . . 1507.2 Content of this thesis . . . . . . . . . . . . . . . . . . . . . . . 1527.3 Future Works . . . . . . . . . . . . . . . . . . . . . . . . . . . 154Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 156Appendix . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 177Appendix A: Evaluation of the PI gain matrices and closed-loopsystem poles . . . . . . . . . . . . . . . . . . . . . . . 178xiList of TablesTable 1.1 Applications of standard thermoformable amorphousthermoplastics [Das22]. . . . . . . . . . . . . . . . . . . 2Table 2.1 Comparison of previous heating phase modeling . . . . 15Table 3.1 Average surface temperature variation models’ coeffi-cients (Equation (3.1)) in two modes for Trough Ce-ramic 500 W. . . . . . . . . . . . . . . . . . . . . . . . 24Table 3.2 System parameters. . . . . . . . . . . . . . . . . . . . . 34Table 3.3 Heater bank condition during the test scenario. . . . . 38Table 3.4 Comparison of mean square errors along six cut-linesin two methods. . . . . . . . . . . . . . . . . . . . . . . 41Table 4.1 Adjacent zones temperatures feasibility. . . . . . . . . 51Table 4.2 PI performance on the simulator for temperature con-trol thermoforming heating phase. . . . . . . . . . . . 66Table 5.1 A comparison between PI and AMPC on thermoform-ing heating phase simulator with 15 heaters for refer-ence temperature profiles in 4.8. . . . . . . . . . . . . . 75Table 5.2 Accuracy of the proposed NARX model for the testdata: NRMSE for each Zone v.s. N-step closed-loopsimulation prediction horizon. . . . . . . . . . . . . . . 86Table 5.3 A comparison between PI, AMPC, and NARX-MPCon thermoforming heating phase simulator with 15heaters for reference temperature profiles in 4.8. . . . . 88Table 5.4 A comparison between open-loop, PI, and NARX-MPCin cycle-to-cycle control approach. . . . . . . . . . . . 90Table 5.5 NARX-MPC validation. . . . . . . . . . . . . . . . . . 97Table 6.1 Estimated DRL agent training time versus the numberof the heaters in the heater bank. . . . . . . . . . . . . 113xiiLIST OF TABLESTable 6.2 Probability function parameters. . . . . . . . . . . . . 115Table 6.3 Network Structures of the implemented artificial neu-ral networks for PPO. . . . . . . . . . . . . . . . . . . 119Table 6.4 Hyperparameters. . . . . . . . . . . . . . . . . . . . . . 119Table 6.5 A comparison between PI, AMPC, NARX-MPC, andMPC-guided DRL on thermoforming heating phasesimulator with 15 heaters for reference temperatureprofiles in 4.8. . . . . . . . . . . . . . . . . . . . . . . . 125Table 6.6 Network Structures of the implemented artificial neu-ral networks for PPO in multi-agent RL approach. . . 128Table 6.7 State vector in zones around the sheet’s periphery. . . 131Table 6.8 Hyperparameters. . . . . . . . . . . . . . . . . . . . . . 132Table 6.9 A comparison between PI, AMPC, NARX-MPC, MPC-guided DRL, and Multi-agent DRL on thermoformingheating phase simulator with 15 heaters for referencetemperature profiles in 4.8. . . . . . . . . . . . . . . . 139Table 6.10 A comparison between open-loop, PI, NARX-MPC,and Multi-agent DRL in cycle-to-cycle control approach.140Table 6.11 Robustness comparison between NARX-MPC and Multi-agent DRL for three system parameter variation. . . . 143Table 6.12 Multi-agent DRL validation. . . . . . . . . . . . . . . . 148Table 7.1 A comparison between NARX-MPC and Multi-agentDRL experimental results for uniform and non-uniformreference temperature distributions. . . . . . . . . . . . 152xiiiList of FiguresFigure 1.1 Schematic of basic vacuum thermoforming. . . . . . . 1Figure 1.2 Theoretical variation of Young modulus according totemperature for amorphous and semi-crystalline poly-mers . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3Figure 1.3 Thermoformed plastic market size estimation in bil-lion USD [pre]. . . . . . . . . . . . . . . . . . . . . . . 4Figure 1.4 Thermoplastic on the mold. . . . . . . . . . . . . . . . 4Figure 1.5 Thickness distribution [T-S, Coma]. . . . . . . . . . . 6Figure 1.6 Defects in output parts during the thermoformingprocess. T1 and T2 show two spots on the productwith different thicknesses. . . . . . . . . . . . . . . . . 7Figure 1.7 Integration of the IR camera in the control loop . . . 9Figure 2.1 Control algorithm’s taxonomy . . . . . . . . . . . . . 17Figure 3.1 Two different heater technologies are mainly beingused in thermoforming [Cer]. . . . . . . . . . . . . . . 23Figure 3.2 Actual surface temperature variation over time in dif-ferent power levels for trough-ceramic-500W vs. modeloutput [HS24]. Used with permission from ©. Else-vier.2024. . . . . . . . . . . . . . . . . . . . . . . . . . 24Figure 3.3 FVM vs. FDM vs. FEM [MA] . . . . . . . . . . . . . 25Figure 3.4 Heat transfer between thermoplastic sheet, radiantheaters, and ambient [HRSM23]. Used under open-access license ©.MDPI.2023. . . . . . . . . . . . . . . 26Figure 3.5 Heat transfer coefficient for free convection betweenthe cold surrounding and hot plate [PTvZ+13] [HRSM23].Used under open-access license ©.MDPI.2023. . . . . 29Figure 3.6 Plate arrangements (parallel plates) [ES93] [HRSM23].Used under open-access license ©.MDPI.2023. . . . . 30xivLIST OF FIGURESFigure 3.7 Heaters numbering in laboratory-scale thermoform-ing heater bank[HRSM23]. Used under open-accesslicense ©.MDPI.2023. . . . . . . . . . . . . . . . . . . 32Figure 3.8 Laboratory-scale thermoforming setup employed. . . . 33Figure 3.9 Electrical components diagram for 8 heaters in two-circuit AC source. Line A: heaters 1, 2, 3, 4. Line B:heaters 12, 13, 14, 15. . . . . . . . . . . . . . . . . . . 33Figure 3.10 Radiation pattern in different types of ceramic heat-ing elements [Cer]. . . . . . . . . . . . . . . . . . . . . 35Figure 3.11 Laboratory-scale setup: temperature distribution onthe sheet using a trough ceramic heater with 500 Wmaximum power consumption[HRSM23]. Used underopen-access license ©.MDPI.2023. . . . . . . . . . . . 35Figure 3.12 Reflected and absorbed radiation on opaque surface[HRSM23].Used under open-access license ©.MDPI.2023. . . . . 36Figure 3.13 Extracted experimental α×Fh→m,n for each elementof the sheet[HRSM23]. Used under open-access li-cense ©.MDPI.2023. . . . . . . . . . . . . . . . . . . 37Figure 3.14 Heating elements number 4 and 8 with surface tem-perature of 603 K and 803 K, respectively[HRSM23].Used under open-access license ©.MDPI.2023. . . . . 38Figure 3.15 Temperature distribution on the sheet[HRSM23]. Usedunder open-access license ©.MDPI.2023. . . . . . . 39Figure 3.16 Actual and simulated temperature variations alongthree cut lines in the y direction of the thermoplas-tic sheet[HRSM23]. Used under open-access license©.MDPI.2023. . . . . . . . . . . . . . . . . . . . . . . 40Figure 3.17 Actual and simulated temperature variations alongthree cut lines in the x direction of the thermoplas-tic sheet[HRSM23]. Used under open-access license©.MDPI.2023. . . . . . . . . . . . . . . . . . . . . . . 41Figure 3.18 Temperature distribution on the sheet after 120 sec,using COMSOL[HRSM23]. Used under open-accesslicense ©.MDPI.2023. . . . . . . . . . . . . . . . . . . 42Figure 4.1 An industrial example of large sheet sagging duringheating phase [CMS]. . . . . . . . . . . . . . . . . . . 46Figure 4.2 The distance separating the heater bank and the sheetand its impact on the interaction between system in-put and output. . . . . . . . . . . . . . . . . . . . . . 47xvLIST OF FIGURESFigure 4.3 Steady-state temperature distribution on the thermo-plastic sheet at different gaps between heater and sheet 48Figure 4.4 Center zone, longitudinal adjacent 1 zone, longitudi-nal adjacent 2 zone, and transverse adjacent zone onthe thermoplastic sheet . . . . . . . . . . . . . . . . . 49Figure 4.5 Reachable temperature range in the adjacent zones. . 50Figure 4.6 A sample reference temperature distribution on thethermoplastic sheet. . . . . . . . . . . . . . . . . . . . 51Figure 4.7 Reference temperature feasibility evaluation in thedifferent gaps between the heater bank and the ther-moplastic sheet. . . . . . . . . . . . . . . . . . . . . . 52Figure 4.8 Temperature distribution from open-loop control (d=15 cm). . . . . . . . . . . . . . . . . . . . . . . . . . . 53Figure 4.9 Cycle-to-cycle open-loop control (d= 15 cm) . . . . . 54Figure 4.10 Temperature distribution from open-loop control (d=1.5 cm). . . . . . . . . . . . . . . . . . . . . . . . . . . 54Figure 4.11 Cycle-to-cycle open-loop control (d= 1.5 cm) . . . . . 55Figure 4.12 a cascaded block diagram of the system . . . . . . . . 57Figure 4.13 PI control (d= 1.5 cm) . . . . . . . . . . . . . . . . . 63Figure 4.14 PI control (d= 15 cm) . . . . . . . . . . . . . . . . . . 64Figure 4.15 Cycle-to-cycle PI control (d= 15 cm) . . . . . . . . . 65Figure 4.16 Robustness evaluation of the designed controller onthe simulator. (a) The average error of zones undervarying convection heat transfer coefficients and gapsbetween the heater bank and the sheet. (b) Maxi-mum overshoot under varying convection heat trans-fer coefficients and gaps between the heater bank andthe sheet. (c) The average error of zones under vary-ing convection heat transfer coefficients and sheet ab-sorptivity coefficients. (d) Maximum overshoot undervarying convection heat transfer coefficients and sheetabsorptivity coefficients. . . . . . . . . . . . . . . . . . 67Figure 5.1 Model Predictive Control [HS24]. Used with permis-sion from ©.Elsevier.2024. . . . . . . . . . . . . . . . 69Figure 5.2 Temperature profile on the sheet on the simulator. . . 74Figure 5.3 Zone errors on the thermoplastic sheet. . . . . . . . . 74Figure 5.4 Control input to 15 heating elements. . . . . . . . . . 75Figure 5.5 a single-output NARX model block diagram. . . . . . 77Figure 5.6 MPC block diagram . . . . . . . . . . . . . . . . . . . 79xviLIST OF FIGURESFigure 5.7 MPC algorithm flowchart in thermoforming process. . 80Figure 5.8 Experiment Design for Data Collection: ExcitationInputs. . . . . . . . . . . . . . . . . . . . . . . . . . . 81Figure 5.9 1-step temperature prediction results refer to compar-ing the temperature prediction and the actual tem-perature in the fifteen zones (blue line: dataset, redline: NARX model’s output, yellow line: error). . . . 83Figure 5.10 Sample Auto-correlation of output yj , j ∈ {1, ..., 15}with 99% Confidence Intervals: Test dataset . . . . . 84Figure 5.11 Sample Cross-correlation of input uj and output yj , j ∈{1, ..., 15} with 99% Confidence Intervals: Test dataset 85Figure 5.12 Temperature profile on the sheet at 1000 sec in thesimulator. . . . . . . . . . . . . . . . . . . . . . . . . . 87Figure 5.13 Zone errors on the thermoplastic sheet . . . . . . . . 87Figure 5.14 Control input to 15 heating elements during controlprocess . . . . . . . . . . . . . . . . . . . . . . . . . . 87Figure 5.15 The cycle-to-cycle control performance of the Narx-MPC in the heating phase of the thermoforming pro-cess. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89Figure 5.16 Robustness evaluation of the designed controller onthe simulator (NARX model is trained with initialvalues of h = 5,W/m2K, d = 15, cm, and α = 0.8).(a) The average error of zones under varying con-vection heat transfer coefficients and gaps betweenthe heater bank and the sheet. (b) Maximum over-shoot under varying convection heat transfer coef-ficients and gaps between the heater bank and thesheet. (c) The average error of zones under vary-ing convection heat transfer coefficients and sheet ab-sorptivity coefficients. (d) Maximum overshoot undervarying convection heat transfer coefficients and sheetabsorptivity coefficients. . . . . . . . . . . . . . . . . . 91Figure 5.17 Reference temperature profiles for experiments . . . . 92Figure 5.18 Experimental temperature profile measured by the IRcamera: (a) t = 250 [sec], (b) t = 300 [sec], (c) t = 350[sec], (d) t = 400 [sec] . . . . . . . . . . . . . . . . . . 93Figure 5.19 Experimental zone errors on the thermoplastic sheetfor a uniform reference temperature distribution. . . . 94xviiLIST OF FIGURESFigure 5.20 Control input to 15 heating elements during the con-trol process for uniform reference temperature distri-bution. . . . . . . . . . . . . . . . . . . . . . . . . . . 94Figure 5.21 Experimental temperature profile measured by the IRcamera: (a) t = 125 [sec], (b) t = 250 [sec], (c) t = 375[sec], (d) t = 500 [sec] . . . . . . . . . . . . . . . . . . 95Figure 5.22 Experimental zone errors on the thermoplastic sheetfor a nonuniform reference temperature distribution. . 96Figure 5.23 Control input to 15 heating elements during the con-trol process. . . . . . . . . . . . . . . . . . . . . . . . 96Figure 6.1 Reinforcement Learning framework. . . . . . . . . . . 100Figure 6.2 Growth in search results of reinforcement learning onGoogle from 2007 to 2022 [Rl-]. . . . . . . . . . . . . 105Figure 6.3 Map of reinforcement learning algorithms (DP (Dy-namic Programming), TD (Temporal Difference), MC(Monte Carlos), I2A (Imagination-Augmented Agent)[WRR+17], DQN (Deep Q-Network) [MKS+13], TRPO(Trust Region Policy Optimization) [SLA+15], ACKTR(Actor-Critic using Kronecker-Factored Trust Region),AC (Actor-Critic), A2C (Advantage Actor-Critic), A3C(Asynchronous Advantage Actor Critic) [MBM+16],DDPG (Deep Deterministic Policy Gradient) [LHP+15],TD3 (Twin Delayed DDPG) [FHM18], SAC (SoftActor-Critic) [HZAL18], PPO (Proximal Policy Op-timization) [SWD+17].) . . . . . . . . . . . . . . . . . 107Figure 6.4 Value-based and policy-based methods [ZY20]. . . . . 107Figure 6.5 Monte Carlo methods and temporal-difference meth-ods [ZY20]. . . . . . . . . . . . . . . . . . . . . . . . . 108Figure 6.6 Actor–Critic PPO algorithm process flow chart [HS24].Used with permission from ©.Elsevier.2024. . . . . . 109Figure 6.7 Flowchart of MPC-guided approach for DRL training[HS24]. Used with permission from ©.Elsevier.2024. . 114Figure 6.8 Distribution of Action selection probability duringtraining process [HS24]. Used with permission from©.Elsevier.2024. . . . . . . . . . . . . . . . . . . . . . 115Figure 6.9 Error bands and overheating area[HS24]. Used withpermission from ©.Elsevier.2024. . . . . . . . . . . . 116Figure 6.10 Reward output[HS24]. Used with permission from©.Elsevier.2024. . . . . . . . . . . . . . . . . . . . . . 118xviiiLIST OF FIGURESFigure 6.11 DRL episode reward during training process when thesystem comprises only one heater controlling the tem-perature of a specific zone . . . . . . . . . . . . . . . 120Figure 6.12 Episode reward during traditional DRL (PPO) train-ing process for 15 heaters[HS24]. Used with permis-sion from ©.Elsevier.2024. . . . . . . . . . . . . . . . 121Figure 6.13 Episode reward during DRL (PPO) training processwith MPC guided approach for 15 heaters[HS24]. Usedwith permission from ©.Elsevier.2024. . . . . . . . . 122Figure 6.14 Episode reward during retraining process[HS24]. Usedwith permission from ©.Elsevier.2024. . . . . . . . . 123Figure 6.15 Temperature profile on the sheet at 1000 sec in thesimulator. . . . . . . . . . . . . . . . . . . . . . . . . . 124Figure 6.16 Zone errors on the thermoplastic sheet. . . . . . . . . 124Figure 6.17 Control input to 15 heating elements. . . . . . . . . . 124Figure 6.18 Single and Multi-Agent RL in heating phase of ther-moforming . . . . . . . . . . . . . . . . . . . . . . . . 126Figure 6.19 The RL architecture in radiation thermal control basedon PPO with Actor-Critic styleThe RL architecturein radiation thermal control based on PPO with Actor-Critic style . . . . . . . . . . . . . . . . . . . . . . . . 127Figure 6.20 Central and adjacent zones definition; ZCe: centralzone, Zu: up neighboring zone, Zd: down neighboringzone, Zr: right neighboring zone, Zl: left neighboringzone. . . . . . . . . . . . . . . . . . . . . . . . . . . . 129Figure 6.21 Zones around the sheet’s periphery . . . . . . . . . . 130Figure 6.22 Random input power to adjacent heaters during agenttraining . . . . . . . . . . . . . . . . . . . . . . . . . . 133Figure 6.23 DRL training during only one episode . . . . . . . . . 133Figure 6.24 DRL agent’s episode reward during the training pro-cess for only one heater. . . . . . . . . . . . . . . . . . 134Figure 6.25 RL agent’s performance on the simulator . . . . . . . 135Figure 6.26 Temperature profile on the sheet after 1000 sec in thesimulator. . . . . . . . . . . . . . . . . . . . . . . . . . 136Figure 6.27 Zone errors on the thermoplastic sheet. . . . . . . . . 136Figure 6.28 Control input to 15 heating elements during the con-trol process. . . . . . . . . . . . . . . . . . . . . . . . 137Figure 6.29 Maximum overshoot obtained from simulator in dif-ferent uniform reference temperature distribution onthe sheet. . . . . . . . . . . . . . . . . . . . . . . . . . 137xixLIST OF FIGURESFigure 6.30 The cycle-to-cycle control performance of the multi-agent DRL in the heating phase of the thermoformingprocess . . . . . . . . . . . . . . . . . . . . . . . . . . 140Figure 6.31 Robustness evaluation of the designed controller onthe simulator (NARX model is trained with initialvalues of h = 5,W/m2K, d = 15, cm, and α = 0.8).(a) The average error of zones under varying con-vection heat transfer coefficients and gaps betweenthe heater bank and the sheet. (b) Maximum over-shoot under varying convection heat transfer coef-ficients and gaps between the heater bank and thesheet. (c) The average error of zones under vary-ing convection heat transfer coefficients and sheet ab-sorptivity coefficients. (d) Maximum overshoot undervarying convection heat transfer coefficients and sheetabsorptivity coefficients. . . . . . . . . . . . . . . . . . 142Figure 6.32 Experimental temperature profile measured by the IRcamera: (a) t = 200 [sec], (b) t = 250 [sec], (c) t = 300[sec], (d) t = 750 [sec] . . . . . . . . . . . . . . . . . . 144Figure 6.33 Experimental zone errors on the thermoplastic sheetfor a uniform reference temperature distribution. . . . 145Figure 6.34 Control input to 15 heating elements during the con-trol process. . . . . . . . . . . . . . . . . . . . . . . . 145Figure 6.35 Experimental temperature profile measured by the IRcamera: (a) t = 150 [sec], (b) t = 200 [sec], (c) t = 250[sec], (d) t = 845 [sec] . . . . . . . . . . . . . . . . . . 146Figure 6.36 Experimental zone errors on the thermoplastic sheetfor a nonuniform reference temperature distribution. . 147Figure 6.37 Control input to 15 heating elements during the con-trol process. . . . . . . . . . . . . . . . . . . . . . . . 147xxGlossary of NotationS(Xt) Sum of dilated and translated scaling functions. 77W (Xt) Sum of dilated and translated wavelets. 77A System matrix. 57Ah Surface area of the heaters. 29B Input matrix. 57C Output transformation matrix. 58E Vector of zones’ average error. 59Fh→i View factor between the hth heater and the ith element on the sheet.29Gt Cumulative reward (Return). 101H Number of heaters. 26K Full state feedback matrix. 59L Linear mapping vector. 77LV (µ) Critic Objective Function. 110LCLIP(θ) Clipped Surrogate Objective Function. 110N Number of samples. 24Nc Control horizons. 73Np Prediction horizons. 73Nx Number of sheet divisions in the x direction. 26Ny Number of sheet divisions in the y direction. 26xxiGlossary of NotationP Projection matrix. 77Q∗π(s, a) Optimal Action-Value Function. 103Qπ(s, a) Action-Value Function. 102T Vector of sheet elements’ temperature. 26Tg Glass transition temperature. 3Tm Melting temperature. 3Tm,n Temperature of element in row m and column n on the sheet. 26Top Vector of sheet elements’ temperature at operating point. 57V ∗π (s) Optimal State-Value Function. 103Vπ(s) State-Value Function. 102Xt Regressor vector. 77∆x Sheet element’s length. 28∆y Sheet element’s width. 28∆z Sheet thickness. 34Ω Weight matrix. 73Π Weight matrix. 73α Absorbed emission. 36αac Learning rate for actor-network. 119αcr Learning rate for critic-network. 119αd Thermal diffusivity. 27β Entropy loss weight. 119Q˙conduction Conduction heat transfer. 28Q˙convection Convective heat transfer. 28Q˙radiation Radiation heat transfer. 28xxiiGlossary of Notationq˙ Energy generation rate per unit volume of the medium. 27ϵ Clipping Parameter. 110γ Discount factor. 99, 100Aˆt Advantage Function. 109Rˆ Random number generator function. 114|B| Batch size. 119x Entire states of the system, including sheet element temp. and heaters’temp. 71A Action space. 99H Experience horizon. 119M Maximum of zones’ absolute error at each time. 117N Number of training epochs per update. 119P State transition probability. 99R Reward function. 99S State space. 99T di Zone’s desired temperature on the thermoplastic sheet. 45Ti Zone’s average temperature on the thermoplastic sheet. 45Tover Allowable overshoot temperature in a zone determined by the materialproperties. 45U Input admissible set. 73Y Output admissible set. 45dU Sheet element’s internal energy rate. 27µ Parameters of the value function (critic) network. 109ϕ Average surface temperature of the heater. 23ϕop Average surface temperature of the heater at operating point. 57xxiiiGlossary of Notationπ Policy. 101ρ Density. 27ρr Reflected emission. 36σ Stefan-Boltzmann constant. 29θ Parameters of the policy (actor) network. 109ε Offset. 77εe Effective emissivity of the sheet-heater. 29εh Heater emissivity. 29εs Sheet emissivity. 29εt exploration noise. 113at Action at time t. 99cp Specific heat capacity. 27d Disrtance between the heater and the sheet. 30e Vector of elements’ error on the sheet. 58h History of events. 100hc Convective heat transfer coefficient. 28nx Number of zone divisions in the x direction. 57ny Number of zone divisions in the y direction. 57rt Reward at time t. 100rt(θ) Likelihood Ratio. 110st Current state. 99u Heaters input heating power in watts. 23v Volume of the element on the sheet. 28xxivAcknowledgementsI began my Ph.D. journey under unusual circumstances, trapped in thelimbo of a lengthy study permit process. My visa application to Canadatook an unprecedented 28 months, forcing me to defer my start semestertwice. I started my Ph.D. remotely from Iran during COVID-19, and thosedays were incredibly challenging. I was working to earn a living and wasdedicating the rest of my time to university research until midnight. Thisschedule continued for 16 months, compounded by the constant fear of losingmy position despite all my efforts. I am grateful to Dr. Rudolf Seethalerfor trusting me during those difficult times and retaining me in his researchgroup. Working with him throughout this time has been an honor, andknowing him has been the most rewarding aspect of my Ph.D. experience.His unwavering support, expert guidance, and constant encouragement havebeen instrumental in completing my PhD research. Without his guidance,this work would not have been possible.I am indebted to my thesis committee members, Dr. Abbas S. Milaniand Dr. Bhushan Gopaluni, for guiding me through my Ph.D. with theirinvaluable feedback. I am also deeply thankful to Dr. Milad Ramezankhani,Mohammad Bajelani, and Rohith Jayaraman Krishnamurthy, who were al-ways present to provide their input and share their expert knowledge.I thank my friends Saeed Askariani, Navid Keshmiri, Reza Azarmanesh,Mohammad Ramezankhani, Arman Rajayi, and Ali Davoodi for making thisresearch journey delightful and memorable. I also wish to express my grat-itude to my off-campus friends, especially Majid Kavianfar, who supportedand encouraged me throughout my Ph.D.Lastly, and most importantly, I thank my family for their love and sup-port. Thank you, Mom and Dad, for allowing me to pursue my academicdreams far from home and unconditionally supporting me. I am proud tohave amazing siblings like you, Ziba, Mehdi, Mina, and Sina. Finally, myspecial thanks to my lovely spouse, Negin, for ceaselessly inspiring me topursue what is meaningful in life. Your presence in my life is my greatestmotivation to strive and succeed.xxvDedicationTo my devoted motherFereshte Mohammadi Khaneghah,my lovely wifeNegin Ghiaeimoghaddam,and to the memory of my beloved uncleDr. Hossein Mohammadi Khaneghah.Your unwavering support and love have been my greatest inspiration.xxviChapter 1IntroductionThere are three main distinct phases to the thermoforming process. Theplastic sheet is heated in the first step to becoming ductile, which makesit softer and more moldable. Different types of plastic require differenttemperatures to be ductile. After the sheet is heated to the proper tem-perature, it moves to the second stage, which involves cooling and molding.The plastic is guaranteed to take on the required shape by several methods,including creating a vacuum between the sheet and the mold. After that,the molded part cools and takes on its shape inside the mold. Excess plasticis cut away during the third phase (trimming) to achieve the final productshape. Thermoforming processes are categorized based on the thickness ofthe thermoplastic sheets used: thin gauge (< 0.25 mm), medium-weightsheet thermoforming (0.25 mm to 1.5 mm), and heavy-sheet gauge (> 1.5mm) [AAM21]. Figure 1.1 shows the schematic of this process [LCRMC+18].Figure (1.1) Schematic of basic vacuum thermoforming. (a) heating; (b)sealing or pre-stretch; (c) forming and cooling; (d) trimming [LCRMC+18].In many contemporary applications, amorphous thermoplastic polymerslike polystyrene (PS), High Impact Polystyrene (HIPS), Polyvinyl Chloride(PVC), Acrylonitrile-Butadiene-Styrene (ABS), and polymethyl methacry-late (PMMA) serve as fundamental constituents. These thermoplastics areoften combined to create complex multi-layer products, catering to variousrequirements in different sectors. Table 1.1 presents an overview of the keyapplication sectors for several thermoformable amorphous thermoplastics.1Table (1.1) Applications of standard thermoformable amorphous thermoplastics [Das22].Thermoformable amorphousthermoplasticsHIPS PS PMMA PVC ABS PCTg(oC) 96 90 100 90 90-120 150Forming window (oC) 160-205 150-190 150-190 110-150 150-205 175-230DailylifeproductsFood packaging[Thr08, MM10, GLM+04][KMK13]*Disposable containers,Serving trays, Plastic cups[KST13]* [TMAS17]* [PV10]*Healthcare and medical[Thr08], [BDGS18]*Cosmetic cases, Trays,Tablet packaging[KST13]* [Cha94]*Consumer goods [SPG+16]Masks, Textured surfacesInner liner of refrigerators[KST13]* [TMAS17]* [PV10]*HightechnologyproductsElectronic components [SSA+17],[VPY+15, CHC+21, WZZ+19]2.5 D shaped circuits,Encapsulation of printedcircuit board assemblies[KST13]*Transport sectors [Thr08, Ros79]Car dashboards,Compartment doors,Aircraft cabin wall panels,[KST13]* [Thr08] [Thr08] [Arn13]*(*) Patetnts21.1. MOTIVATIONFrom a thermo-mechanical perspective, the appeal of amorphous ther-moplastics lies in their wide forming window (FW), which stems from theirrubber-like behavior. This behavior occurs above the glass transition tem-perature (Tg) and extends over various temperatures. The thermoforma-bility of these materials hinges on finding a balance between withstandingsignificant deformations and maintaining a certain level of rigidity to preventexcessive sheet flow during heating. The forming window can be character-ized by conducting temperature sweeps via Dynamic Mechanical Analysis(DMA) and monitoring the evolution of storage moduli, as depicted in Fig-ure 1.2. typically, the forming window ranges from 30 to 60 °C above Tg formost amorphous thermoplastics. In contrast, semi-crystalline thermoplas-tics transition from rigid to rubbery states near their melting temperature(Tm). This report concentrates on the thermoforming of amorphous ther-moplastics.Figure (1.2) Theoretical variation of Young modulus according to tem-perature for amorphous and semi-crystalline polymers [Ern05]1.1 MotivationPlastic consumption has increased globally in recent decades. The in-crease in developing nations such as China and India demonstrates a clear re-lationship between a country’s economic development and its per capita plas-tic consumption. Many industries use plastic extensively, including packag-ing, construction, healthcare, automotive, aerospace, and recreation. Thefavorable properties of plastics, such as their low density, impermeability toliquids and gases, strong mechanical performance in tensile and impact resis-tance, and appealing aesthetic qualities like color, gloss, and tactile appeal,31.1. MOTIVATIONare responsible for their widespread adoption. Moreover, their popularityacross industries is partly attributed to their ease of conversion into manyuseful products. The global thermoformed plastics market was valued atUSD 41 billion in 2021 and is expected to reach around USD 62.1 billion by2030. North America (United States, Canada, and Mexico) dominated theglobal market with a share of around 54% in revenue in 2023 [pre].Figure (1.3) Thermoformed plastic market size estimation in billion USD[pre].Studies have demonstrated that the temperature distribution of the ther-moplastic sheet significantly affects the final product’s thickness distributionin the thermoforming process. To achieve a product with a uniform thicknessdistribution, the temperature profile created on the sheet must be compati-ble with the geometry of the mold. For example, for a parabolic mold shape,in order to avoid having a higher thickness on the floor, the sheet shouldbe heated in the center to a higher temperature compared to other areas,according to Figure 1.4.Figure (1.4) Thermoplastic on the mold.Figure 1.5 illustrates an example of a thermoformed sheet with two differ-ent initial temperature profiles. As shown, when an optimized temperatureprofile is applied to the thermoplastic sheet, the thickness gradient in the fi-41.1. MOTIVATIONnal product decreases from 1.9 mm (when a uniform temperature is applied)to 0.5 mm. This highlights the importance of precise temperature controlin producing high-quality thermoformed products.For a given optimized temperature profile on the sheet, the feasibility ofachieving this reference temperature on the thermoplastic must be assessed,considering the hardware configuration of the setup. Once the referencetemperature is feasible, an appropriate control algorithm is needed to cre-ate this profile on the thermoplastic sheet. Real-time control architecturesfor manufacturing processes must usually satisfy robustness. Despite this,operators often manually regulate the heating bed in the thermoformingprocess, leading to operational inefficiencies and faulty products (see figure1.6). Empirical evidence suggests that this manual regulation produces 5-10% of wasted raw materials and energy.This directly conflicts with the United Nation’s Sustainable DevelopmentGoals that encourage us to save natural resources and minimize the energyused in production processes [tUN15]. Traditional, centralized architecturesbased on PID or H-infinity control mechanisms are insufficiently flexible torespond effectively to these types of highly dynamic processes [ABD+06].These controllers cannot handle process constraints such as the maximumand minimum heater temperatures and the heaters’ limited heating andcooling rates [CBH11]. Due to the limited heating and cooling rate, it isoften observed that the heaters do not track the controller’s control inputs[Aje04]. In MPC, the computation of the optimization problem at everysampling instant may require complex calculations demanding a very fastprocessor [ML99]. Due to the online optimization at every sampling instant,MPC has not been an effective technique for dealing with large multivari-able constrained systems that increase computational complexity in solvingthe optimization problem. Moreover, although several issues like stability,feasibility, and performance of linear MPC control are well developed andunderstood [MRRS00, BMDP02], much work must be done in nonlinearMPC. Thus, the large system size and the presence of nonlinearities in theheating phase of thermoforming seem to have discouraged control engineersfrom using optimal control techniques for this kind of system. This studyproposes reinforcement learning control as a model-free approach for over-coming the above-mentioned heating phase control difficulties. Instead ofdealing with complex mathematical equations of the heating phase model inactuators and radiation heat transfer, the controller is designed and trainedusing command data applied to the system and the system’s response. Re-inforcement learning control uses potent search algorithms for optimal con-trollers of systems with nonlinear, possibly stochastic, unknown, or highly51.1. MOTIVATION(a) Mold geometry.(b) Thickness distribution for a initial uniform temperature profile on the sheet.(c) Thickness distribution for an optimized initial temperature profile on the sheet.Figure (1.5) Thickness distribution [T-S, Coma].61.2. MAJOR CONTRIBUTIONS AND STRUCTURE OF THE PROPOSED RESEARCHFigure (1.6) Defects in output parts during the thermoforming process. T1and T2 show two spots on the product with different thicknesses.uncertain dynamics.1.2 Major Contributions and Structure of theProposed ResearchThe core aim of this thesis is to design and deploy a scaleable, precise,and robust control algorithm tailored for overseeing the heating phase in thethermoforming process. Temperature management is a critical determinantin this phase, where deviations from excessively high or low temperatures canpotentially compromise the integrity of thermoformed products, manifestedin issues like wrinkling and irregular thickness (see Figure 1.6). Moreover,optimizing the heating duration for each sheet is essential to producing ef-ficiency. Therefore, the ultimate goal is to develop a control algorithm thatis easy to implement and scalable for integration into industrial setups. Ad-ditionally, the algorithm should exhibit robustness against environmentaldisturbances commonly encountered in industrial environments. The accu-racy of the temperature distributions achieved on the thermoplastic sheetby the end of the control process is paramount. The quality of the formedproducts is heavily contingent upon precise thermal control, underscoringthe importance of achieving accurate temperature distributions. Two dis-tinct approaches are proposed to achieve this objective: model-based andmodel-free control methodologies.Chapter 1 provides an introduction to the presented thesis. Chapter71.2. MAJOR CONTRIBUTIONS AND STRUCTURE OF THE PROPOSED RESEARCH2 delves into the background of the research presented. It outlines thevarious types of modeling employed in thermoforming, examines existingtemperature-sensing methodologies utilized in thermoforming manufactur-ing, and explores the diverse control strategies previously investigated. Ad-ditionally, it introduces two potential approaches that will be explored inthis thesis: Model Predictive Control and Reinforcement Learning.In Chapter 3, the thermal dynamics of the thermoforming process’s heat-ing phase are modeled in detail. This chapter develops a novel 2D physics-based model to investigate the temperature distribution on the thermoplas-tic sheet during the heating phase of the thermoforming process. An exper-imental view factor matrix was introduced into the 2D model for the firsttime. With this method, the 2D model with low computing effort achievedthe accuracy of a 3D model [HRSM23]. A lab-scale thermoforming setupwas designed and built to study process optimization, and the details of thedeveloped hardware are provided in this chapter. Integrating an IR camerain the control loop to evaluate the thermoplastic sheet temperature feedbackin the thermoforming systems is introduced and implemented for the firsttime. This lab-scale setup was used to design the experiment to producethe experimental view factor matrix between the heater and thermoplasticsheets. Also, this lab-scale setup was utilized to verify the performance ofthe created heating phase simulator.Chapter 4 of this thesis assesses the system’s stability and explores its dy-namic behavior. The control objectives for thermoforming systems are thendiscussed. A novel method is introduced to assess the feasibility and reach-ability of a reference temperature distribution on the thermoplastic sheet,considering the configuration of the setup. The open-loop heating time ofthe sheet for a particular reference temperature distribution is studied onthe simulator in the cycle-to-cycle. The results of the open-loop approach,which is used in the industry, are compared to closed-loop control resultsin the next Chapters to show the importance of automation in this phase.The effectiveness of a proportional-integral (PI) controller is examined asa fundamental control strategy. This chapter explains the shortcomings ofusing only basic control techniques through careful analysis. It emphasizesthe importance of having more complex control strategies adapted to thespecifics of the system.Chapter 5 investigates the application of Model Predictive Control (MPC)as a robust model-based control technique in the heating phase of the ther-moforming process. The chapter initially focuses on designing an adaptiveMPC system based on the linearized model derived from the physics-basedsimulator. Even though the results are better than PI’s on the physics-based81.2. MAJOR CONTRIBUTIONS AND STRUCTURE OF THE PROPOSED RESEARCHFigure (1.7) Integration of the IR camera in the control loopsimulator, hardware limitations in implementation and the high computa-tional cost associated with online execution incentivized the exploration ofalternative control approaches. Thus, in the second step, a novel data-drivenmodel of the thermoforming system is developed using Nonlinear AutoRe-gressive with eXogenous inputs (NARX) techniques. This model is thenintegrated into the MPC framework and deployed on the lab-scale setup.The results showcase a substantial enhancement in MPC performance facil-itated by integrating the NARX model.Chapter 6 studies the realm of Deep Reinforcement Learning (DRL) asa transformative approach for industrial control systems, particularly focus-ing on its application in the heating phase of the thermoforming process.The chapter underscores the time-consuming and expensive nature of theoffline training process for DRL when dealing with systems with large ac-tion and observation spaces, prompting the need for innovative solutions.Two distinct approaches are proposed to address this challenge. Firstly, ahybrid method is introduced, seamlessly integrating MPC with DRL. Thisnovel approach harnesses the predictive capabilities of MPC to accelerate thelearning performance of RL agents, resulting in a remarkable 100-fold im-provement in average episode rewards over traditional DRL methods. Thisapproach is named MPC-guided DRL in this study [HS24]. Secondly, amulti-agent DRL method is presented, wherein each agent controls a specificheating element within the thermoforming setup. This approach showcasessignificant advantages over conventional single-agent RL, particularly re-garding training efficiency and control precision. The offline training processis expedited 110 times compared to traditional DRL, leading to an 8 times91.2. MAJOR CONTRIBUTIONS AND STRUCTURE OF THE PROPOSED RESEARCHreduction in the final error margin. Experimental trials on the laboratory-scale setup validate the applicability and effectiveness of the multi-agent RLapproach, demonstrating robustness in achieving and maintaining desiredtemperature profiles across the thermoplastic sheet.Chapter 7 concludes the thesis by summarizing the key findings. Itelucidates the significant contributions of this research and emphasizes itsimpact and relevance in advancing the understanding and application of con-trol strategies in thermoforming processes. A qualitative comparison tableoffers an insightful comparison of two selected and implemented techniqueson the experimental setup. The overall content of the thesis is explained,and finally, future research directions stemming from the proposed method-ologies are discussed, providing valuable insights for further investigation inthe field.10Chapter 2Literature ReviewComposite manufacturing is a process that can potentially be impactedby industry 5.0 [RFLP22] due to the large number of manufacturing stepsstill performed manually or with minimal automation. Any controller designrequires an accurate model that describes the system’s behavior. For tradi-tional model-based control, it determines the controller structure and allowsfor optimizing the control system’s performance [FPHG15a]. For model-freecontrollers, a precise model can be used to replace expensive physical exper-iments during the training phase of the control systems. Many studies havebeen conducted to reduce the infeasible requirement of large data availabil-ity for developing accurate learning models [RNSM21, RCN+21]. However,since model-free control methods such as Reinforcement Learning (RL) re-quire a large number of “trial and error” episodes or interactions within anenvironment, utilizing simulators is a cost-effective and time-efficient strat-egy for achieving results [NLH20a, MLL+21, MLB+]. Once the algorithmhas been sufficiently trained offline through interaction with a simulated en-vironment, offline learning would stop, and the trained algorithm would bedeployed on the target hardware, similar to a control law. At this stage,the need for continued training on target hardware highly depends on thesimulator’s accuracy during the offline training. Moreover, a precise modelis also necessary to simulate the developed control technique to ensure thestability and robustness of the proposed technique before its implementationin the real world.2.1 Thermoforming Process ModelingThermoforming is the process of heating and deforming a thermoplas-tic sheet into the desired shape [Thr08, Thr96]. Typically, radiant heatersare positioned on one or both sides of the plastic sheet. The length ofthe heating cycle required to soften the sheet adequately depends on thephysical, chemical, and mechanical properties of the polymer, its thick-ness, and its color. The final product quality in thermoforming is influ-enced by various factors, including material characteristics, heating con-112.1. THERMOFORMING PROCESS MODELINGsistency, pressure levels, and surrounding moisture. These variables canintroduce unexpected variations in the thermoforming process, making trou-bleshooting complex. Some defects that affect the final product in multipleways can be introduced during the process, both structurally and aesthet-ically. Among them, the undesired presence of wrinkles and the variationof thickness of the formed component are the most frequently observed inprevious works [DGN+22, EE14, WHB13, BHM16]. Developing reliablepredictive tools for thermoformed product quality assessment involves im-plementing complex nonlinear models, such as finite element simulations[WWJ+03, XHB19, BLAV22] or deep learning models [BBX+21, YH04].Data-intensive deep learning models demand extensive datasets for practicaltraining and may suffer from sub-optimal generalization performance. Someresearch has been done to develop a sim-to-real transfer-learning frameworkto address the data scarcity and heterogeneity in smart manufacturing ap-plications [RHSM23]. These models help predict the required sheet temper-ature distributions and vacuum pressures.Traditionally, the thermoforming industry has relied on operators man-ually adjusting heating elements to minimize defects. Feedback controlschemes could improve this approach but require an optimal temperatureprofile on the sheet for a particular mold to prevent defects in the finalproduct. Thermoplastic sheet deformation models during the thermoform-ing process can be used in an optimization algorithm to achieve this optimaltemperature profile. In these models, the temperature profile on the sheetimpacts the physical and thermomechanical properties in different zones,resulting in different sheet thicknesses as a final product. Hosseini et al.[HBMZ06] created a model that describes how a pre-stretched sheet de-forms during any stage of the vacuum thermoforming process. Additionally,they developed a method based on the wall-thickness distribution criterionfor predicting and improving the quality of the finished products. Wang etal. [WHB13] presented a simulation method for multi-layer thermoplasticdeformation during thermoforming. Semi-discrete shell elements were usedto model each prepreg layer. Xiong et al. [XHB19] studied the consolida-tion behavior of thermoplastic composite prepregs during the thermoformingprocess based on a generalized Maxwell approach. Wang et al. [WLPL22]presented a fully coupled thermomechanical computational model based onthe hot optimal transportation mesh-free method for resin-based frictioncomposites. Wanger et al. [WKM+22] developed a thermoforming modelthat approximatively includes the deforming plastic sheet and the pressurebox on the deforming plastic sheet. Bean et al. [BLAV22] developed a toolthat employs the forming simulation’s local fiber alignment and thickness122.1. THERMOFORMING PROCESS MODELINGpredictions to generate a more precise structural model. This is anticipatedto provide a more accurate prediction of the final part’s performance than astructural model that utilizes initial ply thicknesses and fiber orientations.Consequently, the utility of a forming simulation is increased.The first task in controlling any system is to build a proper model.This model must include the physical interactions between the heating el-ements and the plastic sheet and the behavior of actuators. It is also nec-essary to design efficient controllers for the thermoforming process. Someresearchers have developed simple models to understand heat transfer duringthe heating phase of the thermoforming process, focusing on various heattransfer mechanisms for different polymer processing techniques. Throne[Thr96] developed a model for heat transfer in semitransparent polymersduring the heating phase, considering the wavelength dependency of sheetabsorptivity and heater emissivity. Duarte and Covas [DC02, DC03] inves-tigated the inverse heating problem for thermoforming machines that usethin sheet materials. They created a model based on an analytical viewfactor to determine the optimal heating element temperatures, resulting ina uniform sheet surface temperature. In some research, a 2D finite dif-ferential analysis was combined with an energy balance method to modelthe heating process during the thermoforming of thin thermoplastic sheets[BMBD17, BvMvBD17, WED+17, SLX+10, LMX+10a]. In [Moo02], Mooreintroduces a model for the sheet reheat phase, considering conduction, con-vection, and radiation. This model discretizes the plastic sheet across itsthickness to account for heat propagation only through the thickness. In[Kum05], the dependency of absorptivity and emissivity on the wavelengthof transmitted heat is analyzed, influencing the model development for theheating phase. Researchers in [Sch03] established the spectral propertiesof infrared heat emitters. Additionally, [YBD02] analyzes the uncertaintyof model parameters related to the input-output relationship of the heat-ing phase, focusing on parameter uncertainty in the simulation of the sheetreheat phase in the thermoforming process.In [MSMD98], researchers emphasize optimizing the reheating stage inblow molding and the thermoforming process. They develop a model topredict transient temperature distributions over both thin and thick-gaugepolypropylene thermoformed sheets, using an effective radiative heat trans-fer coefficient and the effective bulk temperature. Ajersch [Aje04] improvesthis model by incorporating heat absorption within the sheet as it heats.Building on Ajersch’s model, Gauthier [Gau08] proposes enhancements tothe heating phase of a thermoforming machine by considering the transmis-sion factor of incident heat over the plastic surface. Gauthier uses the Beer-132.2. TEMPERATURE SENSINGLambert Law to evaluate the portions of heat energy transmitted throughand absorbed by the plastic sheet layer. Chy et al. [CB11] improved Gau-thier‘s model by more accurately evaluating input parameters for convectionheat transfer, investigating other heating sources, and developing a modelfor the air temperature inside the oven. Erchiqui (in [Erc18]) considered aspecific 3D volumetric enthalpy-based finite element method. They also useda technique based on the method of contour proposed in [EN07] to calculatethe view factor, which avoids the problem of singularity while reducing thenumber of numerical arithmetic operations.Although these models effectively simulate the heating phase of the ther-moforming process, discrepancies still exist between the simulation and ex-perimental results. Table 2.1 compares the previously developed heatingphase models in terms of four criteria: considering conduction heat transferin the thermoplastic sheet, considering convection heat transfer on the sheetsurface, including the heater model, and assumption in calculating the viewfactor matrix. A significant weakness in these models is their reliance onmathematical relationships for the view factor, which assumes the heatingelement is a flat surface with a constant temperature. The heater’s out-put radiation pattern varies significantly and depends on its geometry. Inthis thesis, an improved model for the heating phase of the thermoformingprocess is developed by incorporating experimental view factor matrices,resulting in superior-quality predictions.2.2 Temperature SensingTo control the heating phase and achieve the desired temperature distri-bution on the sheet’s surface, a controller requires continuous feedback fromthe measured surface temperature to adjust the emitted infrared heat flux.When performing a contactless temperature measurement, several factorsmust be considered: the sheet material’s temperature-dependent emissivity,the emission and sensor’s wavelength dependence, the sensor’s directionaldependence, and the experimental difficulties associated with filtering outinterfering reflectance [MSMD98, HS94a, HS94b].Infrared temperature measurement has been a focal point in optimiz-ing the thermoforming process. The use of IR heaters, coupled with tem-perature sensors for each heating zone, ensures even heating of the plasticsheet, which is critical for achieving optimal control and quality in the fi-nal product. IR sensors, including cameras and pyrometers, are utilizedextensively in the industry for their ability to measure temperature without142.2. TEMPERATURE SENSINGTable (2.1) Comparison of previous heating phase modelingYear Authors Cond. Conv. HeaterModelingView FactorAssumption2002 Moore [Moo02] × ✓ ✓ Flat Heater2002 A. Yousefi[YBD02]× ✓ × Flat Heater2004 Ajresh [Aje04] × ✓ ✓ Flat Heater2005 Gauthier et al.[GAB+05]× ✓ × Flat Heater2008 Gauthier [Gau08] × ✓ ✓ Flat Heater2011 Chy et al. [CB11] × ✓ ✓ Flat Heater2011 Erchiqui et al.[ESY07]× × × Flat Heater2017 Buffel et al.[BMBD17]× ✓ × Flat Heater2018 Erchiqui [Erc18] ✓ ✓ ✓ Flat Heater2023 Hosseinionari etal.[HRSM23]✓ ✓ ✓ -direct contact [AORADdJRT19]. This capability is particularly beneficialfor monitoring the temperature of moving, extremely hot, or difficult-to-reach objects. The temperature range required for forming is governed bythe type and grade of thermoplastic material used, as well as the intricacyof the product. Several researchers examined the measurement validity andaccuracy of infrared temperature sensors, similar to those used to detectsheet surface temperatures within a thermoforming oven, and found themreliable. Initially, Lappin and Martin [LM96] examined the feasibility ofimplementing these sensors on a thermoforming machine and the potentialcalibration and noise concerns that would arise. Bugbee et al. [BDMT98]classified IR sensors into three categories and distinguished them based ontheir calibration and sensor body temperature errors. Baker et al. [BNK01]devised and prototyped a successful approach to calibrating infrared ther-mometers. Benjamin Neubig [NB19] proposed a new approach by combiningthree different infrared sensors for determining the temperature distributionalong the sheet thickness at any time of the heating phase. Jacob Wilson[WED+17] formulated a temperature distribution control of a rod utilizingthe spatial data from the thermal camera as feedback. Mekid [Mek10] useda thermal infrared camera to acquire spatial temperature on machines forerror compensation via open architecture controllers. Using an IR camera152.3. CONTROL SYSTEM APPROACHEScan be advantageous in two ways. First, instead of dealing with numerousIR sensors that each detect the temperature of a small portion of the sheet,an IR camera can monitor the entire surface. In addition, using an IR cam-era allows calculating the average error for an area instead of using the errorof a single point on the sheet, which can improve controller performance[WED+17].2.3 Control System ApproachesThe taxonomy of control algorithms can be categorized into several do-mains, each focusing on different aspects of control systems. Figure 2.1presents a brief taxonomy of control algorithms that can be employed.Numerous control algorithms have been proposed to address thermal con-trol problems, each with unique advantages and limitations. Notable ap-proaches include model predictive control [CBH11, LVL19], Fuzzy-PID con-trol [YCH+08], robust control [VBHP23], adaptive control [HCLB18], andlearning-based algorithms [YCH+23, OF21, MJ16, BPM+17]. Model-basedcontrollers rely on accurate mathematical models of the system dynamicsand can be sensitive to discrepancies between the assumed model and theactual system. If the model used for controller design does not perfectly rep-resent the real system, or if there are uncertainties in the model parameters,the controller’s performance may degrade [FvHK14]. Furthermore, manymodel-based controllers are designed based on linear system models. Whilelinearization can be useful for certain applications, it may not accuratelycapture the full dynamics of nonlinear systems. Nonlinearities, especially inhighly dynamic systems, can lead to modeling errors and challenges in con-trol design [PS02]. Conversely, data-driven [SAW+23] and learning-basedcontrol methods face challenges related to sample efficiency and algorithmtraining time [HW13].The thermoforming heating process is designed to heat the thermoplasticsheet rapidly to the desired forming temperature. To optimize productionefficiency in terms of production rate, minimizing the time required to cre-ate each part is crucial. Theoretically, the quickest way to obtain a uniformthrough-thickness temperature is to heat the sheet’s outside surface to theforming temperature and maintain it constant. The material’s thermal diffu-sivity impacts the heat transmission rate through it and the time required forthe sheet’s core to reach the forming temperature in this procedure [Buf16].Many thermoforming industries traditionally rely on trial-and-error meth-ods to adjust heater power, which can be inefficient in achieving desired162.3. CONTROL SYSTEM APPROACHESFigure (2.1) Control algorithm’s taxonomysheet temperature setpoints [ABD+06]. Few studies have investigated theoptimal heater input over time [LMX+10b]. These methods are sensitive toenvironmental disturbances and can take a long to obtain the desired tem-perature across the sheet since they perform as open-loop systems. Thus,adopting closed-loop sheet temperature control becomes essential for thethermoforming process as it automates the adjustments of heater power, ef-ficiently achieving the desired sheet temperature profiles [NB19, PNMS21].2.3.1 Proportional-Integral-Derivative (PID)The control of the heating phase of the thermoforming process is chal-lenging due to its nonlinear and time-varying nature. Additionally, the highlevel of uncertainty in various material parameters and the system’s multi-input multi-output (MIMO) nature, with significant coupling between inputsand outputs, further complicate control. Moreover, several hard constraintsmust be considered, including the maximum allowable temperatures for thesheet and the heaters, critical for ensuring material integrity and processsafety. Traditional control methods cannot accommodate these process con-straints [BBM17].Currently, standard proportional-integral-derivative (PID) controllers aresometimes used to control the heating phase of the thermoforming process[CBH11, CB10, HARFCLCA18, Yan08]. However, PID controllers face sig-nificant challenges due to the long time constants inherent in the process172.3. CONTROL SYSTEM APPROACHES[AM99]. These long-time constants result from the low heating and coolingrates of the heaters and the plastic sheet, causing the sheet temperature toovershoot the set point and respond slowly to disturbances. Moreover, PIDcontrollers cannot effectively handle process constraints such as maximumand minimum heater temperatures and limited heating and cooling rates.The limited heating and cooling rates often prevent heaters from ac-curately tracking the control inputs from the PID controller. Additionally,PID controllers typically use minimal model information to calculate controlinputs. Due to the nonlinearities in the heating phase, PID controllers strug-gle to achieve the desired performance across different operating conditions.This inadequate control performance can lead to part quality variations thatfail to meet the industry’s stringent quality standards. Furthermore, PIDcontrollers are unsuited for controlling multi-input multi-output (MIMO)systems.2.3.2 Model Predictive Control (MPC)MPC has emerged as a robust control strategy in various industrialapplications due to its ability to handle multivariable systems with con-straints [DHMDK14]. MPC requires accurate dynamic models, and theprocesses of modeling and identification represent the most demanding andtime-intensive aspects of any industrial MPC project, as highlighted by[Ric93, Zhu98, ZPWZ13]. For MPC maintenance, re-identifying the modelis the primary task. Beyond the challenge of model identification, compre-hending the intricacies of MPC theory, tuning methodologies, and evaluatingcontrol performance proves to be complex. This complexity contributes tothe rarity of highly skilled MPC engineers. Integrating Nonlinear Autore-gressive with Exogenous Inputs (NARX) models within MPC frameworks,known as NARX MPC, has further enhanced the capability to manage com-plex nonlinear systems. NARX MPC offers several advantages, includingimproved prediction accuracy for nonlinear systems, enhanced disturbancerejection, and better handling of multivariable systems with constraints.These benefits lead to more efficient and reliable control systems in indus-trial applications, from process control in bioethanol fermentation to tem-perature regulation in thermoforming and optimization in manufacturingprocesses [LNK+18, ZFLL24].In manufacturing, NARX MPC has been applied to stabilize systems un-der various conditions, achieving close to ideal control performance [SGM+21].NARX MPC has been successfully applied to the bioethanol fermentationprocess, where a NARX-based MPC controller demonstrated superior per-182.3. CONTROL SYSTEM APPROACHESformance in disturbance rejection and setpoint tracking compared to linearMPC controllers. The NARX-MPC could handle variations in operatingconditions with faster response times and more minor errors, proving its ef-fectiveness in maintaining the desired product concentration and fermentertemperature setpoints [MA15].In [NB23], NARXMPC has been applied in chemical processes to controla continuous stirred-tank reactor with a heated jacket. The methodologyproposed for identifying sparse polynomial NARX models using elastic netfor real-time implementable NMPC demonstrated good closed-loop perfor-mance under setpoint changes, unmeasured disturbances, and measurementnoise. In smart buildings, a learning-based MPC approach for thermal con-trol has been proposed, utilizing NARX models to optimize the thermalenvironment efficiently in [EA19]. Hu et al. [HFZ+24] established an accu-rate combustion process prediction model based on NARX, deployed with anMPC control strategy, and realized the uniform combustion in the boiler bycontrolling the opening travel of secondary windgates. Jung and Lee [JL23]applied NARX-MPC on the post-combustion CO2 capture (PCC) processto balance model complexity with prediction accuracy. The researchers in[BXFS22] proposed a method to design a control strategy for NARX mod-els. This strategy aims to ensure offset-free tracking of constant referenceswhile maintaining the nominal stability of the closed-loop system.Physics-based Model Predictive Control (MPC) [CBH11] and IterativeLearning Model Predictive Control (ILMPC) [CB12] were investigated forthe heating phase of the thermoforming process. While these methodsshowed promise, they resulted in high computational costs and requireda full-state feedback solution, which was only evaluated through simulationstudies for nominal cases. Additionally, these methods are challenging toscale to new setups because the physics-based model must be adapted tothe specific hardware configuration.To address these challenges, this study proposes data-driven methods toreduce computational costs. These methods involve using input-output datato create a model, which can be used to design a control algorithm (e.g.,NARX-MPC) or directly develop control algorithms based on input-outputdata (e.g. Reinforcement Learning). This approach simplifies the modelingprocess and enhances scalability and adaptability to different setups.2.3.3 Reinforcement LearningReinforcement Learning (RL) has emerged as a transformative approachin industrial control systems, particularly in applications requiring thermal192.3. CONTROL SYSTEM APPROACHEScontrol and manufacturing process optimization. The evolution of RL fromits theoretical foundations to practical industrial applications marks a sig-nificant trajectory in the advancement of intelligent control systems. Theinception of RL can be traced back to the early works on optimal control anddecision-making processes within operations research. The seminal work byWatkins and Dayan on Q-learning laid the groundwork for RL algorithmsthat could learn optimal policies through interactions with the environmentwithout requiring a model of the system dynamics [PRL+21]. This model-free characteristic of RL is particularly advantageous in complex industrialsettings where accurate models are difficult to obtain or computationally in-feasible to use in real-time control. RL has been applied to various systemsin thermal control, ranging from Heating, Ventilation, and Air Conditioning(HVAC) to more specialized manufacturing processes. For instance, stud-ies have demonstrated the applicability of RL to HVAC systems, where thegoal is to minimize energy consumption while maintaining occupant comfort[KKPP23, SIV22]. These applications often involve using Deep Reinforce-ment Learning (DRL), which combines RL with deep neural networks tohandle high-dimensional state spaces, such as those encountered in buildingtemperature control [NSHJ21, AMM21]. The manufacturing sector has alsowitnessed the integration of RL in process control. RL has been employedin thermal manufacturing processes to optimize the temperature controlof Laser Powder Bed Fusion, a technique used in additive manufacturing[OF21]. The RL framework developed for this application could dynamicallyalter process parameters to minimize defects and ensure product consistency.The versatility of RL is further exemplified in its application to materialhandling within manufacturing systems. By learning optimal policies forthe allocation and movement of materials, RL contributes to the overall effi-ciency of the manufacturing process [GP21]. Moreover, RL has been appliedto optimize the scheduling and dispatching of manufacturing jobs, address-ing complex constraints and dynamic changes in production requirements[Win20]. RL has been utilized to manage the heating phase in hot stamp-ing processes in the specific area of thermal control for manufacturing. AnRL-based system was trained to produce high-quality components with dy-namic cycle time optimization, demonstrating the potential of RL to adaptto varying process conditions and improve production efficiency [NPBB+22].Despite the promising results, implementing RL in industrial settings is chal-lenging. One of the primary hurdles is the need for extensive data collectionand training, which can be time-consuming and resource-intensive. Ad-ditionally, deploying RL policies in real-world systems must consider safetyand robustness to ensure that the learned policies do not lead to undesirable202.3. CONTROL SYSTEM APPROACHESor hazardous outcomes [LDJ23, PB22]. Researchers have explored variousstrategies to address these challenges, such as using simulation environmentsfor training RL agents before deployment [SSR23]. Transfer learning tech-niques have also been investigated to leverage knowledge gained in one taskto accelerate learning in another, potentially reducing the training time anddata requirements [dRTAOR+22]. In summary, the application of RL in in-dustrial thermal control and manufacturing has shown significant promise,with numerous studies demonstrating its potential to optimize complex pro-cesses and reduce operational costs. The continuous development of RLalgorithms and techniques and the increasing availability of computationalresources will likely enhance RL’s capabilities and adoption in the industrialsector.21Chapter 3ModelingThis chapter delves into modeling the sheet reheating process integral tothermoforming, prefaced by an overview of the process and introducing keyterms used throughout the modeling phase.3.1 Thermal Dynamics of the Heating ElementIn the thermoforming process, various heaters ensure that plastic sheetsare heated to the appropriate temperatures for molding. The primary typesof heaters used in thermoforming include ceramic and quartz heaters, asshown in Figure 3.11. Ceramic Heaters: These heaters function as electrical resistors, wherethe current flowing through them generates radiant heat. The power outputof ceramic heaters is regulated through the current’s Pulse Width Modula-tion (PWM). Ceramic heaters are known for their durability and ability toemit long-wavelength infrared radiation, making them suitable for efficientlyheating thick materials.2. Quartz Heaters: Quartz heaters pass electricity through a heatingcoil encased within a quartz tube. The quartz material facilitates the rapidtransmission of infrared heat, making these heaters ideal for applicationsrequiring quick heating. Quartz heaters are valued for their fast responsetimes and are particularly effective for materials requiring precise tempera-ture control.223.1. THERMAL DYNAMICS OF THE HEATING ELEMENT(a) Ceramic HeatingElement.(b) Quartz Heating Ele-ment.Figure (3.1) Two different heater technologies are mainly being used inthermoforming [Cer].Obtaining a precise model of the surface temperature variation of heatingelements in relation to the input electrical power is an essential componentof our thermal model. Since the Trough Ceramic 500 W heating element hasbeen used in our lab-scale thermoforming setup, an experiment was designedto obtain a mathematical model of the heating element to simulate its aver-age surface temperature dynamics in relation to input power. The surfacetemperature variation of the heating element was measured over time usingan infrared (IR) camera at different power consumption levels. The heatingand cooling rates of the heater depend highly on the ambient temperatureand heat transfer via convection. The experiments were conducted in alab setting where the environmental temperature was at 21°C. This studyintroduces a predictive model designed to estimate the fluctuations in thesurface temperature of a heater, as denoted by Equation (3.1). This modeloperates as a first-order system, assuming that the steady-state gain andtime constant change depending on the system’s input power:ϕ˙(u) =1c1ec2u + c3ec4u(c5ec6u + c7ec8u − ϕ− 273) (3.1)Here, ’u’ signifies the input power in watts [W], while ’ϕ’ representsthe average surface temperature of the heater in Kelvin [K]. Specific con-ditions, covering both initial and boundary conditions, must be enforcedto determine the coefficients regulating this equation. At the beginningof the scenario (t = 0), when the input power (u) is nil (u = 0), the233.1. THERMAL DYNAMICS OF THE HEATING ELEMENTheater’s surface temperature corresponds to the ambient temperature, im-plying ϕambient = c5+c7 . Furthermore, as time extends to infinity (t =∞),the heater’s surface temperature attains its steady state point, yieldingϕ˙(u) = 0 and ϕss = c5ec6u + c7ec8u .Table (3.1) Average surface temperature variation models’ coefficients(Equation (3.1)) in two modes for Trough Ceramic 500 W.Model Coefficients HeatingmodelCoolingmodeltime constantc1 355.7 104.7c2 -308e-5 -0.02905c3 8.133 268.5c4 4162e-6 -2321e-6gainc5 311.5 311.5c6 1063e-6 1063e-6c7 -290.5 -290.5c8 -1147e-5 -1147e-50 200 400 600 800 1000 1200Time [sec]0200400600Temperature[o C] 55 W107 W155 W197 W244 W287 W357 W407 W501 WModel(a) heating up.200 400 600 800Time [sec]0200400600Temperature[o C]0 W26 W67 W129 W204 W300 W377 W413 W440 WModel(b) cooling down.Figure (3.2) Actual surface temperature variation over time in differentpower levels for trough-ceramic-500 W vs. model output [HS24]. Used withpermission from ©. Elsevier.2024.The coefficients that define this heating model under downward installa-tion for the Trough Ceramic 500 W heater are detailed in Table 3.1. A com-parison of the model and real-world heater dynamics (Figure 3.2) demon-strates the model’s accuracy. The relative average error for the model isapproximately 4 percent, according to its definition in Equation 3.2, whereN is the number of samples throughout the experiment. It should be notedthat even though the model was originally developed for the heating process,it also performs well for the cooling process.243.2. THERMAL DYNAMICS OF THE THERMOPLASTIC SHEETe¯r =[1NN∑i=1( |ϕactual − ϕmodel|ϕactual)i]× 100 (3.2)3.2 Thermal Dynamics of the ThermoplasticSheetThe heat transfer equations can be solved using several numerical meth-ods, including the Finite Difference Method (FDM), the Finite VolumeMethod (FVM), and the Finite Element Method (FEM). Figure 3.3 high-lights the key differences between these methods: With the FDM, the dif-ferential equation is approximated by a different equation, making use of anequidistant structured calculating mesh. This method is straightforward toimplement. In contrast, the FEM approximates the solution space of thePDE by interpolating the unknown function using a set of basis functions.This method offers more flexibility in decomposing the solution domain butis more complex than FDM. The FVM begins by integrating the PDE overa control volume. Applying the divergence theorem, the integral over thevolume is substituted with an integral over the control volume’s boundary.While simpler than FEM, FVM retains the flexibility of decomposing thesolution domain, albeit with specific volume decomposition requirements.Figure (3.3) FVM vs. FDM vs. FEM [MA]In this section, a model of the heating system with 15 heating elementsabove a thermoplastic sheet is formulated using a finite difference analysis253.2. THERMAL DYNAMICS OF THE THERMOPLASTIC SHEETand an energy balance method. In this model, radiation heat transfer fromheating coils to the thermoplastic sheet and from the thermoplastic sheetto the environment, conduction heat transfer inside the thermoplastic sheet,and convection heat transfer between the thermoplastic sheet and the ambi-ent are considered, as shown in Figure 3.4(a). As depicted in Figure 3.4(a),a zone corresponds to a specific sheet section under a heater, primarily ex-posed to the radiation emitted by the corresponding heating element. Eachzone can receive radiation from all the heaters in the system. The numberof zones is equal to the total number of heating elements.(a) Zone definition. (b) Elements on the sheet.Figure (3.4) Heat transfer between thermoplastic sheet, radiant heaters,and ambient [HRSM23]. Used under open-access license ©.MDPI.2023.The thermoplastic sheet is divided into Nx columns and Ny rows, cre-ating NxNy equal 2D elements with the thickness of d. The temperatureof each element on the sheet is represented by the vector T , where the ithcomponent corresponds to Tm,n, with m and n denoting the row and columnindices of the element. The relationships between i, m, and n are given byEquation (3.3):m =⌊i− 1Nx⌋+ 1, n = i−⌊i− 1Nx⌋Nx, ∀i = 1, 2, . . . , NxNy, (3.3)where ⌊.⌋ represents the floor operator. The system has H heating elements,and the vector u = [u1, u2, ..., uH ]T represents their power consumptions,263.2. THERMAL DYNAMICS OF THE THERMOPLASTIC SHEETwhile ϕ = [ϕ1, ϕ2, ..., ϕH ]T corresponds to the average surface temperatureof each heater. The heating phase of the thermoforming process can beapproximated by a nonlinear system, described by Equation (3.4):[T˙ϕ˙]=[f1(T, ϕ)f2(ϕ, u)], (3.4)where T˙ and ϕ˙ are temperature derivatives of sheet elements and heaters,respectively, and f1 and f2 are nonlinear functions that f1 : DT × Dϕ →RNxNy , DT ⊂ RNxNy , Dϕ ⊂ RH and according to Equation 3.1, f2 : Dϕ ×Du → RH , Du ⊂ RH .Figure 3.4(b) visually illustrates the elements of the thermoplastic sheet.Since thin sheets are being investigated in this study and the sheet’s thick-ness is very small compared to its length and width, just one layer of elementsis considered through the thickness of the sheet. Each element’s four faces(±x,±y) are considered for heat transfer between elements via conduction.According to Figure 3.4(a), the radiation heat transfer from the heating coilsto the sheet occurs through the top face (−z), and every exposed face ofthe element on the sheet undergoes heat transfer via convection. Equation3.5 represents the general heat conduction expression for transient condi-tions with two-dimensional effects, constant properties, and internal heatgeneration [Hol09].∂2T∂x2+∂2T∂y2+q˙k=1αd∂T∂t(3.5)In this equation, αd =kρcpis the thermal diffusivity of the material, k isthe material’s conductivity [W/(m·K)], ρ is the material’s density [kg/m3],cp is the material’s specific heat capacity [J/(kg·K)], and q˙ is the rate atwhich energy is generated per unit volume of the medium [W/m3]. Thecentral-difference form of ∂2T∂x2and ∂2T∂y2according to Equation 3.6 can beused to obtain the finite-difference form of Equation 3.5.∂2T∂x2∣∣∣∣(m,n)=T(m+1,n) + T(m−1,n) − 2T(m,n)(∆x)2,∂2T∂y2∣∣∣∣(m,n)=T(m,n+1) + T(m,n−1) − 2T(m,n)(∆y)2,(3.6)Equation (3.7) illustrates the element’s temperature derivative basedon the internal energy rate (dU). By determining the heat flow on the273.2. THERMAL DYNAMICS OF THE THERMOPLASTIC SHEETsheet from each element in all directions, these values can be aggregated tocalculate the net heat flow, representing the change in the internal energyof the element according to Equation (3.8).T˙ =dUρvcp, (3.7)dU = Q˙radiation + Q˙convection + Q˙conduction, (3.8)In this equation, ρ represents the density of the thermoplastic material, vdenotes the element’s volume, and cp signifies the heat capacity of the ther-moplastic sheet material. To ascertain the temperature change, we examinethe alteration in internal energy for each element on the sheet.We discretize in time using the integer p as t = p∆t. Hence, the timederivative is in terms of the temperature difference at new time (p+ 1) andprevious time (p), separated by the time interval ∆t according to Equation3.9.∂T∂t∣∣∣∣m,n=T p+1m,n − T pm,n∆t(3.9)In this scenario, three distinct types of heat transfer are considered: radi-ation, conduction, and convection(Q˙radiation, Q˙conduction, Q˙convection). First,we have conduction heat transfer between adjacent elements in the ther-moplastic sheet according to Equation (3.10). The element is treated as aplane with a constant temperature for conduction heat transfer in a specificdirection. Consequently, ∆L signifies the width of the element in the x or ydirections, Ae is the cross-sectional surface area, and k denotes the material’sthermal conductivity. The temperature change (∆T ) can be determined bycomparing the temperatures at the previous time step.Q˙conduction = kAe∆T∆L, (3.10)Furthermore, convective heat transfer to the surrounding air is consid-ered according to Equation (3.11). Here, ∆x∆y refers to the element’ssurface area exposed to airflow, while Tambient denotes the air temperature.hc is the convective heat transfer coefficient, which is different for the topand bottom faces of the thermoplastic sheet. For convection heat transfer,the surface temperature of the element is assumed to be the temperaturecalculated at the previous time step.Q˙convection = hc∆x∆y(T − Tambient), (3.11)283.2. THERMAL DYNAMICS OF THE THERMOPLASTIC SHEETThe heat transfer coefficient for free convection between a cold surround-ing and a horizontal hot plate differs in the top and bottom faces. Thiscoefficient depends on the plate surface’s temperature and varies accordingto Figure 3.5.Figure (3.5) Heat transfer coefficient for free convection between the coldsurrounding and hot plate [PTvZ+13] [HRSM23]. Used under open-accesslicense ©.MDPI.2023.Finally, this model incorporates heat transfer through radiation from theheaters. The heat added to each element is proportional to its position inrelation to the heaters. The Equation (3.12) represents the calculation ofradiation heat transfer, where Ah signifies the surface area of the heaters[m2], εe denotes the effective emissivity of the sheet-heater, σ representsthe Stefan-Boltzmann constant, H is the number of heaters, Fh→i is theview factor between the hth heater and the ith element on the sheet, ϕhis the absolute temperature of the hth emitter [K], and Ti is the absolutetemperature of the ith element on the sheet, [K] [HMDS20, Mod03].Q˙radiation = AhεeσNxNy∑i=1H∑h=1Fh→i(ϕ4h − T 4i ), (3.12)If all of the radiation leaving one plane reaches the other, the effectiveemissivity of the sheet-heater system is defined in Equation 3.13, where εhand εs are the emissivities for the heater and sheet, respectively [Hol09].εe =[1εh+1εs− 1]−1(3.13)293.2. THERMAL DYNAMICS OF THE THERMOPLASTIC SHEETThe relationship between the energy received by the second surface(sheet) and the energy emitted by the first surface (heater) is known as theview factor. Mathematical equations for view factors are available for therelative positions of the surfaces [EN07, IM16]. The analytical view factorfor two parallel, unequal, and non-coaxial rectangle plates is shown in Figure3.6 and is given by Equation 3.14 [ES93]. In this equation, (x1, y1), (x1, y2),(x2, y1), and (x2, y2) represent the coordinates of the vertices of plane A1,and (ξ1, η1), (ξ1, η2), (ξ2, η1), and (ξ2, η2) represent the coordinates of thevertices of plane A2. These two planes are spaced apart by d.Figure (3.6) Plate arrangements (parallel plates) [ES93] [HRSM23]. Usedunder open-access license ©.MDPI.2023.F(1→2) =1(x2 − x1)(y2 − y1)2∑l=12∑k=12∑j=12∑i=1(−1)(i+j+k+l)G(xi, yj , ηk, ξl)(3.14)G(xi, yj , ηk, ξl) =12π((yj − ηk)√(xi − ξl)2 + d2 tan−1(yj − ηk√(xi − ξl)2 + d2))+12π((xi − ξl)√(yj − ηk)2 + d2 tan−1(xi − ξl√(yj − ηk)2 + d2))− d24πln((xi − ξl)2 + (yj − ηk)2 + d2)(3.15)303.3. UBC LAB-SCALE THERMOFORMING SETUPEmploying these heat transfer mechanisms, the model provides a com-prehensive understanding of temperature changes and their effects withinthe system. Equation 3.16 can be used to obtain new nodal temperaturesat a new time step (p+1) by substituting Equation 3.9 in Equation 3.5.T(p+1)(m,n) = Tp(m,n) + αd ·∆t(T p(m+1,n) + Tp(m−1,n) − 2T p(m,n)(∆x)2+T p(m,n+1) + Tp(m,n−1) − 2T p(m,n)(∆y)2+1k((hpt + hpb)d(Tambient − T p(m,n))+Ahϵeσ∆x∆yH∑h=1F(h→m,n)((ϕph)4 − (T p(m,n))4)))(3.16)Equation 3.16 is explicit since the unknown nodal temperatures at time(p + 1) are determined with known temperatures at time p in each timestep. The initial condition at t = 0 when p = 0 for the temperatureof each node on the thermoplastic sheet is equal to ambient temperatureT (x, y, 0) = Tambient. Moreover, since the thermoplastic is being kept in-side a clamp frame and the thermal conductivity of the thermoplastic isvery low, the temperature of all elements around the sheet in touch withthe clamp was considered constant. The boundary condition is ∂T (0,y,t)∂t =∂T (l,y,t)∂t =∂T (x,0,t)∂t =∂T (x,w,t)∂t = 0 where l and w are the length and widthof the sheet, respectively. According to Equation 3.16, since radiation andconvection primarily occur through the upper and lower faces, the rate oftemperature change in each element will be slower for thermoplastic sheetswith greater thickness in the thermoplastic sheet which leads to a greatervolume in elements.3.3 UBC Lab-scale Thermoforming SetupThis study focuses on the problem of controlling the temperature ofzones on the thermoplastic sheet within a thermoforming system, utilizingthe physical apparatus depicted in Figure 3.8. The setup comprises a heaterbank of fifteen ceramic heating elements (FTE-500, Ceramicx Company)placed in a layout of five rows and three columns, as depicted in Figure 3.7.The sheet is divided into zones. Each zone is associated with a sheet areabelow a heater, mainly exposed to radiation by the corresponding heating313.3. UBC LAB-SCALE THERMOFORMING SETUPFigure (3.7) Heaters numbering in laboratory-scale thermoforming heaterbank[HRSM23]. Used under open-access license ©.MDPI.2023.element. Every zone within the system receives radiation from all the heaterspresent, and the overall count of zones corresponds to the quantity of heaters.A mobile frame securely holds the thermoplastic sheet during the heatingstep. After the heating phase, the sheet is moved to the mold utilizingthe mobile frame before being drawn into the mold using vacuum pressure.The system incorporates a FLIR thermal camera (FLIR A38) seamlesslyintegrated into the MATLAB environment using the FLIR thermal SDK.This camera offers a 60 Hz frame rate, 16-bit video streaming, and Powerover Ethernet capability. With dimensions of just 29 × 36 × 59 mm (5.1× 5.5 × 3.9 in). The thermal camera provides temperature feedback data,enabling the controller to send the control values to the Arduino, whichcontrols the relay state. Figure 3.9 illustrates the AC/DC conversion toan Arduino powered at 5V DC from a 220V AC source. The Arduino isconnected to a solid-state relay to control using PWM. This image showstwo phases connected with 4 heaters in each phase.The control box consists of two primary functions: controlling the heaterpower by varying the switching rate of the relay using PWM and commu-nicating with the computational unit via a serial connection between thecomputer and the Arduino. A two-circuit AC source is divided into twodistinct circuits to manage the power requirements of the heating elementswhile preventing overloading of the laboratory power supply. The first cir-cuit powers 8 heaters, each drawing 2.5 amps, and the second circuit powers7 heaters, drawing 2.5 amps per heater. This division ensures efficient uti-lization of the available laboratory power supply, operating at 220V and upto 20A.Each heater is connected to a solid-state relay (Omega - SSRL240AC10)323.3. UBC LAB-SCALE THERMOFORMING SETUPFigure (3.8) Laboratory-scale thermoforming setup employed.Figure (3.9) Electrical components diagram for 8 heaters in two-circuitAC source. Line A: heaters 1, 2, 3, 4. Line B: heaters 12, 13, 14, 15.333.4. EXPERIMENT TO DETERMINE HEATING ELEMENT’S RADIATION PATTERNTable (3.2) System parameters.Parameter Description Values Unitρ sheet density 1380 [ kgm3 ]cp specific heat capacity 1465[Jkg·K]εe effective emissivity 0.95 -k material thermalconductivity0.18[Wm·K]h convection heattransfer coefficient5 [ Wm2·K ]∆z sheet thickness 0.002 [m]Lx Sheet length 0.75 [m]Ly sheet width 0.5 [m]with a capacity of 10A. These relays provide high isolation capability (4000Vrms input to output and 2500 Vrms input/output to ground) to preventelectrical leakage and minimize the risk of short circuits. Additionally, therelays feature low capacitance (8 pF input to output max), ensuring mini-mal signal distortion and accurate transmission of control signals. Operatingwithin a line frequency range of 47 to 63 Hz, these relays accommodate di-verse power supply configurations. Control of the relays is facilitated bya digital pin on the Arduino Atmega 2560 microcontroller. The microcon-troller communicates with a computer via a serial port, executing firmwaredesigned to regulate the relays based on input from MATLAB code simu-lating the RL network. The input values from MATLAB, representing thecontrol inputs in the range of 0 to 255 units, are utilized by the Arduinofirmware to generate PWM signals to control the relays effectively. Table3.2 shows this study’s model parameters.3.4 Experiment to Determine Heating Element’sRadiation PatternThe 2D models of the heating phase in the thermoforming process devel-oped using analytical equations for the view factor assume that the heatingelement is a flat, temperature-maintaining plate. However, in practice, theradiation pattern varies significantly based on the geometry of the heatingelement. Given the trough ceramic heating coil used in the laboratory-scalesetup, an experiment is conducted to determine the exact view factor ma-trix between the trough ceramic heating coil and the sheet elements. Thesaturated surface temperature for a trough ceramic heating element at its343.4. EXPERIMENT TO DETERMINE HEATING ELEMENT’S RADIATION PATTERN(a) Flat CeramicHeating Element.(b) Hollow CeramicHeating Element.(c) Trough CeramicHeating Element.Figure (3.10) Radiation pattern in different types of ceramic heating ele-ments [Cer].(a) t=15 sec. (b) t=30 sec.(c) t=45 sec. (d) t=60 sec.Figure (3.11) Laboratory-scale setup: temperature distribution on thesheet using a trough ceramic heater with 500 W maximum powerconsumption[HRSM23]. Used under open-access license ©.MDPI.2023.maximum electrical power level of 500 W is 803 K. Figure 3.11 depicts theresulting temperature distribution variation on a thermoplastic sheet locatedat a distance of 15 cm over a 60 s interval.353.4. EXPERIMENT TO DETERMINE HEATING ELEMENT’S RADIATION PATTERNFigure (3.12) Reflected and absorbed radiation on opaquesurface[HRSM23]. Used under open-access license ©.MDPI.2023.Since the distance between the sheet and the heating element is smalland there are reflectors behind each heater, we assume that all radiationfrom the center heating coil used in the experiment (Heater number 8 inFigure 3.7) strikes the surface of the sheet. Some thermoplastics, such aspolycarbonate and acrylic, are highly transparent, while PVC and ABS areopaque. Since an opaque acrylic thermoplastic sheet loaded with a specificamount of pigment is used in our experiment, all infrared emissions areconsidered either reflected (ρr) or absorbed (α) by the thermoplastic sheet(Figure 3.12). The fractional relationship between these two variables isgiven by Equation 3.17, which totals 1 or 100 percent. For a black surface,which absorbs all incoming radiation and emits the maximum possible, α = 1and ρr = 0 [Mod03].ρr + α = 1 (3.17)Ignoring the conduction and convection parts in Equation 3.16 and usingonly one heater in the center of the sheet (heater number 8 in Figure 3.7), theproduct of the view factor between the elements of the sheet and the heaterand the absorptivity of the sheet can be calculated at the initial heatingtime using Equation 3.18. The emissivity for the ceramic heater (εh) in thesetup is 0.92 [Cer].αFh→m,n =k∆x∆ydαdAhεeσ((ϕph)4 − (T pm,n)4) × (T p+1m,n − T pm,n∆t) (3.18)Figure 3.13 shows the extracted experimental α×Fh−m,n for each elementof the sheet using Equation 3.18 and the thermal and physical properties ofthe thermoplastic shown in Table 3.2. As can be seen, the elements underthe heater receive a larger portion of the emitted energy.363.4. EXPERIMENT TO DETERMINE HEATING ELEMENT’S RADIATION PATTERN00.51,#F h!m;n#10-54001.5300x(mm)y(mm)200 200100 02468101214#10-6(a) Numerical values of α × Fh→m,n insurface graph.250 500x(mm)0150300y(mm)2468101214#10-6(b) Numerical values of α× Fh→m,n incontour graph.Figure (3.13) Extracted experimental α× Fh→m,n for each element of thesheet[HRSM23]. Used under open-access license ©.MDPI.2023.Considering the summation of the view factors rule and the superpositionrule, if the sheet is divided into Nx×Ny equal elements, the view factor fromthe heater to the sheet is equal to the sum of all view factors from the heaterto the elements of the sheet. Further, it is assumed that all radiation fromthe heating coil strikes the sheet’s surface; Equation 3.19 requires that thesummation of view factors from the heater to the sheet elements equals 1.Fheater→sheet =Nx∑m=1Ny∑n=1Fheater→m,n = 1 (3.19)The total α × Fh→m,n of all elements in Figure 3.13 is 0.55. If theabsorptivity is assumed to be the same in all elements on the sheet, accordingto Equations 3.19 and 3.20, the total product of the view factor between theelements of the sheet and the heater and the absorptivity of the sheet willequal α, indicating that approximately 45 percent of the emitted energy isreflected by the thermoplastic sheet due to its white and smooth surface.Nx∑m=1Ny∑n=1αFheater→m,n = αNx∑m=1Ny∑n=1Fheater→m,n = α (3.20)373.5. PHYSICS-BASED HEATING PHASE MODEL VERIFICATIONTable (3.3) Heater bank condition during the test scenario.#Heater Heating Power [W] Heaters’ Steady StateSurf. Temp. [K]4 200 6038 500 803Figure (3.14) Heating elements 4 and 8 powered by 200W and 500W reachsurface temperatures of 603 K and 803 K, respectively[HRSM23]. Usedunder open-access license ©.MDPI.2023.3.5 Physics-based Heating Phase ModelVerificationThe laboratory-scale setup described in 3.3 validates the proposed heat-ing phase models described in section 3.1. A test scenario with the heaterbank condition specified in Table 3.3 ran for two minutes after the heaterbank reached the steady state temperature. Figure 3.14 shows that onlyheating elements 4 and 8 received 200 W and 500 W of electrical power andreached steady state surface temperatures of 603 and 803 [K], respectively.Figure 3.15 compares the thermoplastic sheet’s actual and predicted tem-perature distribution after 120 seconds. Clearly, the model based on mea-sured view factors (Figure 3.15(b)) provides a better temperature distribu-tion than the traditional analytical model based on flat plate assumption(Figure 3.15(c)). Three cut lines are applied in the x and y directions ofthe sheet, as shown in Figure 3.16(a) and Figure 3.17(a), to compare betterthe temperature variation between the actual heated sheet and the model’soutput. Figures 3.16(b)-3.16(d) and Figures 3.17(b)-3.17(d) compare thevariation in actual and simulated temperatures on the sheet along these cutlines.383.5.PHYSICS-BASEDHEATINGPHASEMODELVERIFICATION(a) laboratory-scale setup. (b) designed model based on experi-mental radiation pattern.(c) designed model based on analyticalview factors.Figure (3.15) Temperature distribution on the sheet[HRSM23]. Used under open-access license ©.MDPI.2023.393.5. PHYSICS-BASED HEATING PHASE MODEL VERIFICATION(a) cut lines locations. (b) cut line 1 temperatures(c) cut line 2 temperatures (d) cut line 3 temperaturesFigure (3.16) Actual and simulated temperature variations along three cutlines in the y direction of the thermoplastic sheet[HRSM23]. Used underopen-access license ©.MDPI.2023.To assess the accuracy of the designed model along these cut lines, theMean Squares Error (MSE) and Root Mean Squares Error (RMSE) werecalculated using Equations 3.21 and 3.22. Table 3.4 lists the MSE and RMSEalong six cut lines in the model using two different methods. Comparing theproposed method’s performance to that of the same model when utilizingthe analytical method [ES93] to calculate the view factors demonstrates thatthe proposed method performs significantly better.MSE =1NN∑i=1(Tactual − Tmodel)2i (3.21)RMSE =√MSE =√√√√ 1NN∑i=1(Tactual − Tmodel)2i (3.22)403.5. PHYSICS-BASED HEATING PHASE MODEL VERIFICATION(a) cut lines locations. (b) cut line 1 temperatures(c) cut line 2 temperatures (d) cut line 3 temperaturesFigure (3.17) Actual and simulated temperature variations along three cutlines in the x direction of the thermoplastic sheet[HRSM23]. Used underopen-access license ©.MDPI.2023.Table (3.4) Comparison of mean square errors along six cut-lines in twomethods.Vertical HorizontalExperimentalRadiationPatternAnalyticalViewFactorsExperimentalRadiationPatternAnalyticalViewFactorsMSE RMSE MSE RMSE MSE RMSE MSE RMSECutline 1 4.3 2.07 22.2 4.71 5.5 2.34 16.2 4.02Cutline 2 10.4 3.22 23.7 4.86 3.8 1.94 40.5 6.36Cutline 3 2.4 1.54 5.7 2.38 3.4 1.84 16.8 4.09413.6. CHAPTER SUMMARYFigure 3.18 depicts the implementation of the same scenario in the COM-SOL software [COMb] to compare the execution speed of the proposedmodel. Given the configuration of computation hardware, it was imprac-tical to implement the actual geometry of heating elements in COMSOL.As a result, a flat plate model of heaters equivalent to the analytical viewfactors was implemented. Despite implementing a simple model of heat-ing elements, the scenario required 480 s to be computed by the COMSOLsoftware. Comparing the 45 s computation time of the proposed model forthis scenario with the same configuration of computation hardware indicatesthat the proposed model is nearly ten times faster.Figure (3.18) Temperature distribution on the sheet after 120 sec, usingCOMSOL[HRSM23]. Used under open-access license ©.MDPI.2023.3.6 Chapter SummaryA 2D physics-based model of the heating phase of the thermoformingprocess is developed in this study to explore the temperature distributionon the thermoplastic sheet. This work was motivated by the fact that ge-ometrically exact 3D models available through commercial FEM softwareare computationally too expensive for model-free control strategies such asreinforcement learning, which requires many “trial and error” episodes orinteractions with the model to train the algorithm. Although conventional2D models of the heating phase of the thermoforming process are computa-tionally efficient, they assume the heating element is a flat plate and employanalytical view factor relations. This assumption creates a discrepancy be-tween the model and the actual output radiation pattern of the heater. The423.6. CHAPTER SUMMARYmethod proposed in this thesis for developing a 2D model of the thermo-forming heating phase for the first time uses view factors determined fromphysical experiments that better reflect the effect on flat heater geometries.The proposed model with experimentally measured view factors reduces theMSE and RMSE four and twice, respectively, compared to conventional 2Dmodels using analytical view factors. Moreover, comparing the computa-tion times of the proposed model and the COMSOL model indicates thatthe developed model is almost ten times faster. As potential next steps, thepresented heat transfer model can be integrated into more extensive processmodeling and optimization platforms, e.g., to migrate thermoforming in-duced wrinkles [RM18]. Non-local modeling [DDM+20] of the forming sheetmay also be worthwhile to account for different material scales within thecomposite structure.43Chapter 4Stability Analysis and PIControl of the Heating Phasein Thermoforming Process4.1 IntroductionThe heating phase of the thermoforming process is crucial as it sig-nificantly impacts the quality and efficiency of the final product. Preciseand time-optimal temperature control ensures the thermoplastic materialachieves the desired formability without defects. PI controllers are widelyused in many industrial applications due to their simplicity and robustness.However, their fixed structure can limit their performance, especially inhandling the complex dynamics of the thermoforming process where systemparameters may vary. This chapter examines the performance of the PI con-troller as a fundamental control approach during the thermoforming heatingphase and analyzes the stability of the closed-loop system. I will evaluatethe PI controller’s performance in a thermoforming heating scenario and dis-cuss its limitations. Additionally, I propose a method to assess the feasibilityof implementing a reference temperature distribution on the thermoplasticsheet, considering the configuration of the thermoforming setup.4.2 Control Objectives for ThermoformingSystemsOur main objectives for controlling the heating phase of thermoformingsystems are efficiency, reliability, and scalability. These criteria are explainedin the following:444.2. CONTROL OBJECTIVES FOR THERMOFORMING SYSTEMS4.2.1 Heating Beyond Glass Transition TemperatureOne of the critical control objectives is the handling of overheating (tem-perature overshoot). From a thermo-mechanical perspective, the appeal ofamorphous thermoplastics lies in their wide forming window, which stemsfrom their rubber-like behavior. This behavior occurs above the glass tran-sition temperature (Tg) and extends over various temperatures. The ther-moformability of these materials hinges on finding a balance between with-standing significant deformations and maintaining a certain level of rigidityto prevent excessive sheet flow during heating. Typically, the forming win-dow ranges from 30 to 60 ◦C above (Tg) for most amorphous thermoplastics[Sch72]. The overshoot should be carefully controlled so the thermoformingprocesses do not exceed the forming window of the material. This constraintis essential to prevent material degradation and ensure product quality. Con-sequently, the admissible set, Y is defined as follows:Y = {Ti|i ∈ {1, ..., 15}, Ti − T di < Tover} (4.1)where Ti and T di are each zone’s average and desired temperature on thethermoplastic sheet, respectively. Tover denotes the allowable overshoot tem-perature determined by the material properties, which is considered equalto zero in this thesis.4.2.2 Operating Time EfficiencyTo enhance production efficiency, the control system must minimize therise time and settling time during the heating phase. A rapid and accurateapproach to the desired temperature setpoint, without excessive overshoot,is essential for increasing throughput.4.2.3 Terminal Control ErrorTo maintain precise temperature control and ensure product quality, thecontrol system must minimize the terminal error at the end of each controlhorizon. Specifically, the absolute difference between the desired setpointtemperature and the actual temperature of each zone should be less than5°C. This constraint can be expressed as:|Ti(tf )− T di | < 5◦C, ∀i ∈ 1, . . . , 15 (4.2)where tf denotes the end of the control horizon, and Ti(tf ) represents theactual temperature of zone i at time tf . Enforcing this terminal error con-straint ensures that the temperature distribution across the thermoplastic454.2. CONTROL OBJECTIVES FOR THERMOFORMING SYSTEMSsheet is within an acceptable tolerance range, preventing material degrada-tion and maintaining consistent product quality.4.2.4 RobustnessIndustrial environments present numerous challenges, including varia-tions in system parameters such as ambient temperature, which can affectthe thermoforming process. Another important process parameter is the gapbetween the heater bank and the thermoplastic sheet. A smaller gap reducesdisturbance from neighboring heaters, resulting in system behavior similarto a Single-Input Single-Output (SISO) configuration for each zone. How-ever, a small gap may form hot spots on the thermoplastic sheet, disruptingthe desired smooth temperature distribution. During the heating phase,high temperatures cause the large sheets to lose their structural integrityand sag downwards under their weight; see Figure 4.1. This sagging leadsto different gaps between heaters in the heat bank and their associated zonesduring the heating phase. The proposed method must be robust to sensitivesystem parameters, maintaining optimal control performance despite exter-nal disturbances. This robustness is crucial for ensuring consistent processoutcomes and minimizing the need for manual adjustments.Figure (4.1) An industrial example of large sheet sagging during heatingphase [CMS].4.2.5 ScalabilityAnother aspect is the control system’s scalability, especially in industrialsettings where thermoforming processes may involve hundreds of heaters.The designed control system must be scaled without significant modifica-tions. This scalability ensures that the control strategy can be efficientlyapplied across different machines and setups, facilitating its adoption in di-verse industrial environments.464.3. ASSESSING FEASIBILITY OF REFERENCE TEMPERATURE PROFILES ON THERMOPLASTIC SHEETS4.3 Assessing Feasibility of ReferenceTemperature Profiles on ThermoplasticSheetsThis section investigates the feasibility of the given reference temperatureprofile created on the thermoplastic sheet by the controller, considering thehardware configuration. Three key parameters affect this feasibility: thedimension of the heaters, their radiation pattern, and the gap between theheater bank and the thermoplastic sheet. Smaller heaters in dimension allowfor more accurate temperature profiles on the sheet; however, they maylack sufficient power, and a large number of heaters in the bank can resultin significant computational costs in the control algorithm. Depending onthe geometry of the heaters, they may exhibit different radiation patterns.Wide radiation patterns lead to a smoother distribution, complicating thecontrol problem as they affect adjacent zones. Conversely, narrower patternsare preferable for control algorithms as they simplify the system behavior,making it more akin to the SISO system. Still, they cause hot spots on thethermoplastic sheet as well.(a) high distance between heaterbank and the sheet.(b) high distance between heaterbank and the sheet.Figure (4.2) The distance separating the heater bank and the sheet andits impact on the interaction between system input and output.By considering the first two parameters as the nature of the heater, thegap between the heater bank and the thermoplastic sheet can be studied asa variable here. As this gap increases, heater radiation covers a wider areaon the thermoplastic sheet. Conversely, as the gap decreases, the radiationmainly concentrates on the heater’s primary zone. This behavior has asimilar effect as the radiation pattern of the heater on the control system.Figure 4.3 shows the steady-state temperature distribution when oneheater is above the sheet’s center according to Figure 4.4. Heaters are lo-474.3. ASSESSING FEASIBILITY OF REFERENCE TEMPERATURE PROFILES ON THERMOPLASTIC SHEETS0.80.60x [m]500.4Temperature[/ C]0.2100y [m]0.20.10 030405060708090100110120(a) d = 5cm.0.80.64060x [m]0.480Temperature[/ C]0.2100y [m]0.20.10 030405060708090100110120(b) d = 15cm.0.80.660x [m]0.480Temperature[/ C]0.2100y [m]0.20.10 030405060708090100110120(c) d = 25cm.Figure (4.3) Steady-state temperature distribution on the thermoplasticsheet at different gaps between heater and sheet484.3. ASSESSING FEASIBILITY OF REFERENCE TEMPERATURE PROFILES ON THERMOPLASTIC SHEETSFigure (4.4) Center zone, longitudinal adjacent 1 zone, longitudinal adja-cent 2 zone, and transverse adjacent zone on the thermoplastic sheetcated at three different distances: 5 cm, 15 cm, and 25 cm above the ther-moplastic sheet. To reach the average steady-state temperature of 100 °Cin the center zone, input powers of 100 W, 230 W, and 380 W are appliedto the heaters, respectively.As depicted in Figure 4.3, the increase in the gap between the heaterand the sheet results in the heater’s radiation having a more pronouncedeffect on the adjacent zones. This investigation was extended to all powerlevels of the heater, ranging from 50 to 500 W. The average steady-statetemperatures of the four zones on the thermoplastic sheet—namely, thecenter zone, longitudinal adjacent 1, longitudinal adjacent 2, and transverseadjacent zone (see Figure 4.4)—in three different gaps between the heaterand the thermoplastic sheet are illustrated in Figure 4.5.Figure 4.5 can be used to assess the reachability of a reference temper-ature distribution considering the gap adjustment between the heater bankand the thermoplastic sheet. The area between the adjacent and center zonecurves represents the feasible area to reach. On the other hand, for adjacentzones, the temperatures under their curves are unreachable.For example, for a reference temperature distribution according to Figure4.6(a), the average temperature for each zone would be as follows:T = [103, 49, 101, 79, 47, 78, 48, 50, 49, 38, 68, 41, 36, 101, 40]This average zone temperature vector is depicted in Figure 4.6(b). Theevaluation will focus solely on zones with higher temperatures. For instance,consider zone number 14, with an average temperature of 101 °C. Using Fig-ure 4.5, the adjacent reference temperatures for three different gap statuses494.3. ASSESSING FEASIBILITY OF REFERENCE TEMPERATURE PROFILES ON THERMOPLASTIC SHEETS0 100 200 300 400 500Power [W]0100200300400Zonesteadystatetemp.[o C]d=5 [cm]center zonelongitudinal adjacent zone 1longitudinal adjacent zone 2transverse adjacent zone(a) Steady-state temperatures on the thermoplasticsheet when d = 5cm.0 100 200 300 400 500Power [W]050100150200Zonesteadystatetemp.[o C]d=15 [cm]center zonelongitudinal adjacent zone 1longitudinal adjacent zone 2transverse adjacent zone(b) Steady-state temperatures on the thermoplasticsheet when d = 15cm.0 100 200 300 400 500Power [W]20406080100120140Zonesteadystatetemp.[o C]d=25 [cm]center zonelongitudinal adjacent zone 1longitudinal adjacent zone 2transverse adjacent zone(c) Steady-state temperatures on the thermoplasticsheet when d = 25cm.Figure (4.5) Reachable temperature range in the adjacent zones.504.3. ASSESSING FEASIBILITY OF REFERENCE TEMPERATURE PROFILES ON THERMOPLASTIC SHEETSbetween the heater and the sheet would be as shown in Table 4.1. Table4.1 also shows the minimum values for adjacent zones that can be reached,considering the gap between the heater bank and the thermoplastic sheet.According to Table 4.1, the given reference temperature is not attainablewhen the gap between the heater bank and sheet is 15 cm and 25 cm.0.80.640x [m]600.480Temperature[/ C]0.2y [m]1000.20.10 0406080100120140(a) temperature distribution on thethermoplastic sheet.0.720x [m]0.36500Temperature[/ C]0.075y [m]1000.150.225 00.3406080100120140(b) average temperature of each zone.Figure (4.6) A sample reference temperature distribution on the thermo-plastic sheet.Table (4.1) Adjacent zones temperatures feasibility.Zone Ref. d = 5[cm] d = 15[cm] d = 25[cm]center 101 101 101 101longitudinal adjacent 1 68 61 91 97longitudinal adjacent 2 50 30 66 84transverse adjacent 40 and 36 29 48 63Figure 4.7 visually illustrates the described feasibility. According to thefigure, when the temperature of the specified zone is regulated by heatersat a particular distance from the sheet, reaching the desired temperature inthe zone results in the adjacent zones exceeding their desired temperatures,particularly when the distance between the heater bank and the sheet is 15cm and 25 cm (see Figures 4.7(b) and 4.7(c), respectively). However, thisreference temperature distribution is achievable when d = 5 cm (see figure4.7(a)).514.3. ASSESSING FEASIBILITY OF REFERENCE TEMPERATURE PROFILES ON THERMOPLASTIC SHEETS0.80.60x (m)500.4Temperature[/ C]0.2100y (m)0.20.10 0406080100120140(a) d = 5 cm.0.80.60x (m)0.450Temperature[/ C]0.2y (m)1000.20.10 0406080100120140(b) d = 15 cm.0.80.640x (m)600.480Temperature[/ C]0.2100y (m)0.20.10 0406080100120140(c) d = 25 cm.Figure (4.7) Reference temperature feasibility evaluation in the differentgaps between the heater bank and the thermoplastic sheet.524.4. OPEN-LOOP CONTROL RESULTS4.4 Open-loop Control ResultsConsider u as the input power vector to the heater bank located 15 cmabove the thermoplastic sheet. This establishes the steady-state tempera-ture distribution on the thermoplastic sheet as depicted in Figure 4.8.u = [197, 118, 138, 69, 53, 40, 155, 17, 78, 0, 0, 0, 226, 154, 153]0.6100x [m]0.4Temperature[/ C]1500.2y [m]0.20.10 080100120140160180(a) 3D.0 0.2 0.4 0.6x [m]00.10.20.3y[m]80100120140160180(b) 2D.Figure (4.8) Temperature distribution from open-loop control (d= 15 cm).In an industrial setting, these input values must be consistently appliedto the heater bank for a predetermined duration to achieve this temperatureprofile on the sheet using an open-loop control method. Figure 4.9 illustratesthe sheet heating process over two cycles when the heater bank is placed15 cm above the sheet using an open-loop control approach. The sheetand heaters initially start at 21 °C ambient temperature. The first cyclerequires approximately 1880 seconds to reduce the temperature error across15 zones to under 5 °C. After 1880 seconds; the sheet is replaced with anew one, initiating the second cycle. As the heaters remain warm from theprevious cycle, the starting conditions are the final surface temperaturesfrom the prior cycle, resulting in a reduced heating time of 1460 seconds forthe second cycle.If the distance between the heater bank and the thermoplastic sheet isreduced, the heating time can be shortened, and the average temperature ofthe zones on the sheet can be obtained more precisely. However, this smallerdistance can also cause non-smooth deviations between the zone tempera-tures. As illustrated in Figure 4.10(a), this is the steady-state temperature534.4. OPEN-LOOP CONTROL RESULTS0 500 1000 1500 2000 2500 3000-100-500Error[/ C]0 500 1000 1500 2000 2500 30000200400Heaterstemp.[/ C]0 500 1000 1500 2000 2500 3000Time [sec]0100200Power[W]cycle #2cycle #1New sheet is replacedFigure (4.9) Cycle-to-cycle open-loop control (d= 15 cm)0.6100x [m]0.4Temperature[/ C]0.2150y [m]0.20.10 080100120140160180(a) 3D.0 0.2 0.4 0.6x [m]00.10.20.3y[m]80100120140160180(b) 2D.Figure (4.10) Temperature distribution from open-loop control (d= 1.5cm).distribution on the sheet when the heater bank is positioned 1.5 cm abovethe sheet, and the input power vector to the heater bank is as follows:u = [100, 76, 67, 100, 72, 65, 98, 70, 63, 99, 74, 65, 102, 81, 68]Figure 4.11 illustrates the sheet heating process over the same processtime in Figure 4.9 when the heater bank is placed 1.5 cm above the sheet544.5. PI CONTROL DESIGN AND STABILITY ANALYSIS0 500 1000 1500 2000 2500 3000-100-500Error[/ C]0 500 1000 1500 2000 2500 30000100200Heaterstemp.[/ C]0 500 1000 1500 2000 2500 3000Time [sec]050100Power[W]Figure (4.11) Cycle-to-cycle open-loop control (d= 1.5 cm)using an open-loop control approach. The sheet and heaters initially startat 21 °C ambient temperature. The first cycle requires approximately 1216seconds to reduce the temperature error across 15 zones to under 5 °C. After1216 seconds; the sheet is replaced with a new one, initiating the secondcycle. As the heaters remain warm from the previous cycle, the startingconditions are the final surface temperatures from the prior cycle, resultingin a reduced heating time of 825 seconds for the next cycles.This duration is still considerably lengthy for heating a sheet. Addition-ally, this method is highly susceptible to environmental and system distur-bances such as sheet sagging. These observations underscore the necessityfor feedback control in this process. In the subsequent chapter, I will thor-oughly examine all critical parameters that influence system dynamics andcompare the performance of the designed controller to that of the open-loopcontrol strategy.4.5 PI Control Design and Stability AnalysisProportional-integral (PI) controllers are fundamental components incontrol engineering, offering a simple yet effective means of system regu-lation. These controllers combine the immediate response of proportionalcontrol with the steady-state error elimination characteristic of integral con-554.5. PI CONTROL DESIGN AND STABILITY ANALYSIStrol. Due to their versatility and efficiency in handling dynamic systems, PIcontrollers are common in various applications. The essence of PI control liesin its dual-action strategy, which addresses the magnitude and duration ofthe error signal. This dual action allows for a more nuanced control strategy,balancing aggressive correction and steady-state accuracy. The mathemati-cal foundation of PI controllers involves the tuning of two parameters: theproportional gain (kP ) and the integral gain (kI). These parameters shapethe controller’s response to disturbances and setpoint changes.The current system can be approximated with a second-order system;that is, we only need a PD or a PI controller to place all poles. Since aPI controller eliminates steady-state error, it was chosen for our applica-tion. The PI controller’s objective in dynamic systems is to minimize errorsbetween the desired and actual system outputs. This is achieved by gener-ating a control signal that adjusts the system input based on the currentand accumulated errors. The design process of a PI controller involves de-riving the system’s equations of motion and applying the PI control law tothese equations. This section investigates the mathematical formulation ofthe PI controller design. PI controllers can be applied to MIMO (MultipleInput, Multiple Output) systems. However, the design and implementationof PI controllers for MIMO systems can be more complex than SISO (Sin-gle Input, Single Output) systems due to interactions between inputs andoutputs. The subsequent analysis provides a step-by-step derivation of thesystem’s response under PI control, highlighting the role of the controller’sparameters in shaping the system dynamics.4.5.1 State Feedback StabilityAs described in [Kha02], the back-stepping approach allows for concep-tualizing the system’s thermal model as a cascade connection of two maincomponents. The first component represents the heat transfer dynamicswithin the sheet, which is influenced by the input ϕ (heater surface tem-perature). The second component pertains to the dynamics of the heatingelements, which are responsible for adjusting the input power u to achievethe desired temperature. We initiate the analysis with the heat transfer dy-namics of the sheet, considering ϕ as the input. Subsequently, we devise thefeedback control ϕ within the heater temperature control loop to regulatethe system’s output to the desired temperature (Td). The designed controlinput ϕ is then regarded as the desired trajectory ϕd for the subsequentcontrol loop. Moving to the dynamics of the heating elements, the inputpower u is tailored through the heater temperature control loop. This loop564.5. PI CONTROL DESIGN AND STABILITY ANALYSISensures that the heating elements effectively track the desired trajectory ϕd.Figure (4.12) a cascaded block diagram of the systemThe heating phase of the thermoforming process is a nonlinear systemthat can be represented as:T˙ = f(T, ϕ), (4.3)where T is the NxNy × 1 vector, including each node’s temperature on thethermoplastic sheet, and ϕ is the vector of heating elements’ average surfacetemperature. Let (Top,ϕop) be our operating point, then by linearizing theequation 4.3 will be as follows after linearization:T˙ = AT +Bϕ, (4.4)where A and B are as follows:A(NxNy×NxNy) =[∂f∂T](Top,ϕop)=∂f1∂T1· · · ∂f1∂TNxNy.... . ....∂fNxNy∂T1· · · ∂fNxNy∂TNxNy ,andB(NxNy×H) =[∂f∂ϕ](Top,ϕop)=∂f1∂ϕ1· · · ∂f1∂ϕH.... . ....∂fNxNy∂ϕ1· · · ∂fNxNy∂ϕH .System’s Output MatrixTo define the output of the system, the sheet is divided into X = Nx/nxcolumns and Y = Ny/ny rows, creating Z = XY equal 2D zones, wherenx and ny are integer numbers and indicate the columns and rows in eachzone. The system’s output is represented by the vector T , which contains574.5. PI CONTROL DESIGN AND STABILITY ANALYSISthe average temperature of all elements within each zone. The relationshipbetween T and the individual temperature within each zone is expressed byEquation (4.5):T = CT, (4.5)where T ∈ RZ is the output vector of the system and C ∈ R(Z×NxNy) isa matrix representing the averaging operation over each zone according toEquation (4.6).C =ZNxNyC ′X×(Xnxny) 0X×(Xnxny) . . . 0X×(Xnxny) 0X×(Xnxny)0X×(Xnxny) C′X×(Xnxny) . . . 0X×(Xnxny) 0X×(Xnxny)....... . .......0X×(Xnxny) 0X×(Xnxny) . . . C′X×(Xnxny) 0X×(Xnxny)0X×(Xnxny) 0X×(Xnxny) . . . 0X×(Xnxny) C′X×(Xnxny) ,(4.6)whereC ′X×(Xnxny) =[C(1) C(2) · · · C(ny) ] ,C(1)X×(Xnx) = C(2)X×(Xnx) = · · · = C(ny)X×(Xnx) = 11×nx 01×nx 01×nx01×nx 11×nx 01×nx01×nx 01×nx 11×nx .In matrix C, components with a value of 1 in each row denote that thecorresponding elements on the sheet are part of the respective zone. Weaim to control the system from T to the desired temperature Td that canhave the following dynamics:T˙d = ATd +Bϕd, (4.7)where ϕd is the input vector relevant to Td and represents the average sur-face temperature of heaters that should radiate after reaching the desiredtemperature distribution to keep the sheet at the same state. Now let ussubtract Equation 4.4 and 4.7, then the error dynamics would be as follows:e˙ = Ae+BΦ, (4.8)where e = T −Td and Φ = ϕ−ϕd. Considering the system limitation, thereare not enough heating elements to control every element on the sheet, sothe system is not fully controllable. Hence, we suggest reducing the orderof the system so that the control strategy fulfills the system demand. Forthis purpose, researchers [WED+17] suggested dividing the system modelinto separated parts and using an average error for each part. Let us divide584.5. PI CONTROL DESIGN AND STABILITY ANALYSISthe sheet into Z rectangular zones, one under each heating element. As aresult, the vector of average error of zones (E) is given by Equation 4.9:E = Ce, (4.9)where C ∈ R(Z×NxNy) according to Equation (4.6). By incorporating thesystem dynamics into E˙, we obtain:E˙ = Ce˙ =⇒ E˙ = C(Ae+BΦ). (4.10)If we apply the full state feedback Φ = −KE, the closed-loop systembecomes as Equation 4.11:E˙ = C(Ae−BKE) =⇒ E˙ = C(A−BKC)e, (4.11)where K is the full state feedback matrix, the first theorem stability of Lya-punov is used to show the system’s stability for the mentioned input.Lyapunov Global Asymptotic Stability Theorem:The Lyapunov Global Asymptotic Stability Theorem is a fundamentalresult of control theory that provides conditions for analyzing the stabilityof dynamic systems. Stability analysis is crucial in understanding the be-havior of dynamical systems, particularly in engineering and applied mathe-matics contexts. Stability analysis based on linearization provides informa-tion about the local stability of the nonlinear system around the operatingpoint. It does not necessarily guarantee stability for the nonlinear systemas a whole. However, if the linearized system is stable, the nonlinear systemmay also exhibit stability near the operating point. In summary, linearizinga nonlinear system around an operating point allows you to leverage Lya-punov stability analysis techniques to assess the system’s stability. Whilelinearization provides insights into local stability, it can be a valuable toolfor understanding the behavior of nonlinear systems and designing controlstrategies.In this section, we introduce the Lyapunov Global Asymptotic StabilityTheorem, which offers insights into the stability properties of equilibriumpoints in dynamical systems. The theorem establishes criteria under whichan equilibrium point is globally asymptotically stable, implying that thesystem’s trajectories converge to the equilibrium point regardless of the ini-tial conditions. The Lyapunov Global Asymptotic Stability Theorem relieson the concept of Lyapunov functions, which are scalar functions definedover the state space of a dynamical system. These functions are candidates594.5. PI CONTROL DESIGN AND STABILITY ANALYSISfor assessing stability, and their properties determine the system’s stabilitycharacteristics.Energy and Dissipation Functions: consider the nonlinear systemx˙ = f(x) and function V : Rn → R. We define V˙ : Rn → R as:V˙ (z) = ∇V (z)T f(z), (4.12)V˙ (z) gives ddtV (x(t)) when z = x(t), x˙ = f(x). We can think of V as ageneralized energy function and −V˙ as the associated generalized dissipationfunction. Suppose there is a function V such that:− V is positive definite− V˙ (z) < 0 for all z ̸= 0, V˙ (0) = 0Then, every trajectory of x˙ = f(x) converges to zero as t→∞ (i.e., thesystem is globally asymptotically stable). V is a positive definite generalizedenergy function, always dissipating energy, except at 0. For the currentheating system, let’s define the energy function asV = ETE, (4.13)where E represents the vector of zones’ average error. Consequently, thederivative of V with respect to time, denoted as V˙ , is given byV˙ = E˙TE + ET E˙, (4.14)substituting (4.11) into (4.14), we deriveV˙ = eT (A−BKC)TCTE + ETC(A−BKC)e, (4.15)we introduce a positive definite matrix Q such thatC(A−BKC) = −QC, (4.16)with Q defined asQ = diag(λ1, λ2, . . . , λZ), (4.17)where λi > 0 for all i. This leads toV˙ = −ET (Q+QT )E = −2ETQE, (4.18)indicating that V˙ is negative definite and the system is asymptotically stable.Hence, as long as the state feedback matrix is chosen to satisfy the Equation4.19 and all components of the matrix Q stay positive, the system is stable.K = BT (A+ CTQC)CT . (4.19)604.5. PI CONTROL DESIGN AND STABILITY ANALYSIS4.5.2 Sheet Temperature Control LoopEquation (4.10) can be expressed as Equation (4.20):E˙ = AE + BΦ, (4.20)where A and B are matrices that represent the order reduced system dy-namics. From Eqs. (4.10) and (4.20), we can deduce that:A = CACT (CCT )−1, (4.21)B = CB, (4.22)where C is a transformation matrix, A is the system matrix, and B is theinput matrix of the initial system before order reduction. A and B arederived by linearizing the system around the operating point (Top, uop) ateach step. The feedback control command ΦPI is defined as:ΦPI = −KE, (4.23)where K and E are given by:K =[kI1 kP1], E =[∫E dt E]Tdefining the integral of the error and the error itself as state variables:X1 =∫E dt, X2 = Ewe can write the closed-loop tracking error state equation as:X˙1 = X2,X˙2 = AX2 + BΦ,which leads to the state-space representation:E˙ =[0Z×Z IZ×Z0Z×Z A]E+[0Z×ZB]ΦPI , (4.24)where 0Z×Z denotes the Z × Z zero matrix, and IZ×Z denotes the Z ×Z identity matrix. Substituting Equation (4.23) into Equation (4.24), weobtain the closed-loop tracking error state equation:E˙ =[0Z×Z IZ×Z−BkI1 A− BkP1]E. (4.25)614.5. PI CONTROL DESIGN AND STABILITY ANALYSISFollowing the approach in [ZH10], to stabilize the error, the feedbackcontroller command ΦPI = −KE should be designed such that:−BkI1 = diag(−ω21, . . . ,−ω2Z) , A− BkP1 = diag (−2ω1ζ1, . . . ,−2ωZζZ)if all ω and ζ are considered equal, the proportional and integral gains canbe determined as:kI1 = ω2B−1,kP1 = B−1 (A+ 2ωζIZ×Z) ,(4.26)4.5.3 PI Control Loop StabilityTo investigate the stability of the PI control loop, it can be demonstratedthat for the evaluated K = [kI1 , kP1 ] around any operating point, the polesof the closed-loop system are located in the negative half-plane, ensuringstability. For the desired temperature distribution shown in Figure 4.8,the best performance of the PI controller is obtained when ω = 13500I15×15and ζ = 10I15×15. Around any operating point, K can be calculated usingEquation 4.26. Since the system is linearized around the operating point ateach time step, according to Equation 4.26, K will vary at each time step.For example, the evaluated values of kI1 and kP1 at time step t = 1secduring the control process are provided in Appendix A. Employing this Kmatrix values, the eigenvalues of Equation 4.25 are computed, as illustratedin Appendix A. It is observed that the real part of all eigenvalues is negative,indicating system stability. The analysis demonstrates that the PI controlloop remains stable around any operating point, as evidenced by the negativereal parts of the eigenvalues.4.5.4 Heater Temperature Control LoopTransitioning to the dynamics of the heating elements, as elaborated inSection 3.1, the input power u is adjusted through a secondary control loop.This loop is instrumental in ensuring that the heating elements adhere to thedesignated trajectory ϕd. Given that the major radiation of heaters is down-ward and not to the side where the other heaters are, treating each heatingelement as a Single-Input, Single-Output (SISO) system within this frame-work is pragmatic. The PI control parameters for this part of the systemwere calibrated based on a model linearized around the operational workingtemperature of 250◦C. This calibration process yielded a proportional gainkP2 = 2 and an integral gain kI2 = 0.002.624.5. PI CONTROL DESIGN AND STABILITY ANALYSIS0 100 200 300 400 500 600 700 800 900 1000-150-100-50050Error[/ C]0 100 200 300 400 500 600 700 800 900 1000Time [sec]0200400600Power[W]Figure (4.13) PI control (d= 1.5 cm)4.5.5 PI Results on the SimulatorTo assess the performance of the designed PI controller, the identicaltemperature distributions obtained from the open-loop control scenario inFigure 4.8 and 4.10 are considered the reference temperature profiles. Fig-ure 4.13 shows the performance of the PI controller during 1000 sec on theheating phase simulator when the heater bank is placed 1.5 cm above thesheet. The first figure in 4.13 displays the temperature error across differentzones of the thermoplastic sheet, highlighting a settling time of 460 secondsto achieve an error margin within ±10◦C. At the control process’s conclu-sion, the zones’ average and maximum error are approximately 1.9 and 3.5°C, respectively. The maximum overshoot observed was approximately 20°C. The second figure in 4.13 shows the heating power during the controlprocess.On the other hand, Figure 4.14 shows the performance of the PI con-troller when the heater bank is placed 15 cm above the sheet. As seen inFigure 4.14 the controller reduced the error of zones within the ±10°C bandafter 850 seconds, with the average and maximum error of zones at the con-trol process’s conclusion being approximately 4 and 8 °C, respectively. Themaximum overshoot observed was approximately 27.5 °C. The second figurein 4.14 shows the heating power during the control process.In the thermoforming industry, where the heating phase is repetitive, the634.5. PI CONTROL DESIGN AND STABILITY ANALYSIS0 100 200 300 400 500 600 700 800 900 1000-150-100-50050Error[/ C]0 100 200 300 400 500 600 700 800 900 1000Time [sec]0200400600Power[W]Figure (4.14) PI control (d= 15 cm)actuators’ initial state vectors in each cycle are derived from the final statevectors of the preceding cycle. This linkage underscores the importance ofanalyzing cycle-to-cycle control. Figure 4.15 illustrates the performance ofthe designed PI in a cycle-to-cycle framework. It shows that the initial sur-face temperature of the heaters in each cycle is set to the final temperaturefrom the previous cycle. Each cycle concludes once the temperature erroracross all zones is reduced to below 5 °C. At this point, the heated sheet isexchanged for a new one. A comparison with an open-loop control approach,depicted in Figure 4.9, highlights that PI slightly improves the heating time.However, the sheet experienced a high overshoot of 27.5 °C in the first and34 °C in the second cycles.The first cycle requires approximately 1405 seconds to reduce the tem-perature error across 15 zones to under 5 °C according to Figure 4.15. After1405 seconds, the sheet is replaced with a new one, initiating the secondcycle. As the heaters remain warm from the previous cycle, the startingconditions are the final surface temperatures from the prior cycle. Therewas an expected reduction in the heating time of the second cycle com-pared to the first. However, this new initial condition of heaters led to agreater overshoot (34 °C) in the second cycle that the PI controller couldnot regulate, resulting in a cycle time of 1425 seconds.I have defined ten key metrics to comprehensively evaluate the PI con-644.5. PI CONTROL DESIGN AND STABILITY ANALYSIS0 500 1000 1500 2000 2500 3000-100-500Error[/ C]0 500 1000 1500 2000 2500 30000500Heaterstemp.[/ C]0 500 1000 1500 2000 2500 3000Time [sec]0500Power[W]Figure (4.15) Cycle-to-cycle PI control (d= 15 cm)troller’s performance, as outlined in Table 4.2. These metrics include factorssuch as the reliance of the control on the plant’s model, the need for heaterbank temperature feedback, offline training time, online computational ef-fort, average error after 1000 seconds, maximum error after 1000 seconds,overshoot, settling time within an error band of ±10 °C, and first and secondcycle time. Table 4.2 presents the results of these metrics specifically for thescenario where the heater bank is positioned 15 cm above the thermoplas-tic sheet. It is worth noting that, due to the controller’s utilization of thelinearized model around the system’s operational point, feedback regardingthe heaters’ temperature is required at each step of the control process.654.5. PI CONTROL DESIGN AND STABILITY ANALYSISTable (4.2) PI performance on the simulator for temperature control ther-moforming heating phase.Metric Open Loop PIReliance on model - All the timeRequires heater bank’s temperature feedback No YesOffline training time [sec] - -Online computational effort [sec] - 90Average error after 1000 sec [°C] 19.5 4Maximum error after 1000 sec [°C] ±25 ±8Overshoot [°C] - 27.5Settling time within an error band of ±10 °C [sec] 1360 850First cycle time 1880 1405Second cycle time 1460 1425Robustness AnalysisThe robustness of the PI controller was assessed by varying the followingsystem parameters: the convection heat transfer coefficient (h), the gapbetween the heater bank and the thermoplastic sheet (d), and the radiationabsorptivity of the thermoplastic sheet (α). As these parameters change,the PI is evaluated for performance metrics, such as maximum overshootand final average error across zones.The convection heat transfer coefficient on the top face of the sheet wasstudied since it has a higher value than the bottom face of a horizontal hotplate. The gap between the heater bank and the thermoplastic sheet is oneof the most important system parameters. Furthermore, during the heatingphase, high temperatures cause the sheet to lose its structural integrity andsag downwards under its weight. This sagging leads to different gaps be-tween heaters in the heat bank and their associated zones during the heatingphase. The radiation absorptivity of various thermoplastic sheet materialssignificantly impacts the controller’s adaptability to different materials. Asa result, it is critical to consider this parameter during the analysis andevaluation stages.Figure 4.16 showcases the average error across zones and the maximumovershoot under variations in the three parameters for the same referencetemperature distribution in Figure 4.8. Figures 4.16(a) and 4.16(c) illustratethe zones’ average error at the end of the control process under different plantparameter values. Also, Figures 4.16(b) and 4.16(d) show the maximumovershoot in the control process under the same conditions. The absorptivitycoefficient of the sheet in the plant was kept at 0.8 in Figures 4.16(a) and4.16(b), and only (h) and (d) were studied. On the other hand, the gap664.5. PI CONTROL DESIGN AND STABILITY ANALYSIS0 5 10 15h [ Wm2"K ]051015averageerror[o C] d=10 cmd=15 cmd=20 cm(a)5 10 15h [ Wm2"K ]102030overshoot[o C] d=10 cmd=15 cmd=20 cm(b)5 10 15h [ Wm2"K ]5101520averageerror[o C] ,=0.5,=0.6,=0.7,=0.8,=0.9(c)5 10 15h [ Wm2"K ]01020overshoot[o C] ,=0.5,=0.6,=0.7,=0.8,=0.9(d)Figure (4.16) Robustness evaluation of the designed controller on the simu-lator. (a) The average error of zones under varying convection heat transfercoefficients and gaps between the heater bank and the sheet. (b) Maxi-mum overshoot under varying convection heat transfer coefficients and gapsbetween the heater bank and the sheet. (c) The average error of zonesunder varying convection heat transfer coefficients and sheet absorptivitycoefficients. (d) Maximum overshoot under varying convection heat transfercoefficients and sheet absorptivity coefficients.674.5. PI CONTROL DESIGN AND STABILITY ANALYSISbetween the heater bank and the thermoplastic sheet was kept at 15 cm inFigures 4.16(c) and 4.16(d) and only (h) and (α) were studied. The collectivefindings from Figures 4.16(a) and 4.16(c) demonstrate that the PI is notrobust for sensitive system parameter variations, where the average error ofzones vary between 2 °C and 18 °C. Additionally, as illustrated in Figures4.16(b) and 4.16(d), within this range of system parameters, the maximumovershoot varies between 0 and 28 °C. Notably, the convection heat transfercoefficient of the environment exhibits a reverse effect on the zone averageerror and the maximum overshoot. While a higher convection heat transfercoefficient of the environment decreases the controller’s accuracy in terms ofthe zone average error, it results in a lower overshoot in the control process.4.5.6 Chapter SummaryGiven the observations from the previous analyses, it becomes evidentthat the PI controller, while effective under certain conditions, exhibits sig-nificant limitations due to its sensitivity to changes in system parameters,such as the distance between the heater bank and the thermoplastic sheet.The necessity to retune the controller’s gains with every alteration in systemconfiguration underscores an inherent inefficiency in its design for applica-tions where system parameters are prone to variation. This inefficiencycomplicates the operational process and limits the controller’s adaptabilityand responsiveness to dynamic system changes. Furthermore, in optimizingthermoforming manufacturing processes, where the goal is to minimize theheating time for each part to increase manufacturing efficiency, the PI con-troller’s inability to swiftly adapt to varying conditions leads to suboptimalcontrol. This extends the heating duration and consequently reduces theoverall efficiency of the manufacturing process. These limitations highlightthe need to explore more sophisticated controller designs that inherently ac-count for system parameter variations, thereby eliminating the constant needfor manual gain adjustments and ensuring consistent performance across abroader range of operating conditions, aligning with efficiency and produc-tivity goals in high-speed manufacturing environments.68Chapter 5Model Predictive Control ofthe Heating Phase inThermoforming Process5.1 IntroductionThe heating phase in the thermoforming process is a complex operationwhere precise temperature control is essential for product quality. ModelPredictive Control (MPC) is a sophisticated control strategy capable of han-dling such complexity by predicting future system states and optimizing con-trol actions. This optimization procedure is schematically shown in Figure5.1 [FPHG15b]. In this section of the thesis, I will investigate the applica-tion of MPC to the heating phase of thermoforming, employing two distinctmodeling approaches to facilitate the controller’s predictive capabilities.(a) Block diagram of MPC-based ther-mal control loop.(b) Function principle of MPC withprediction and control horizons of Np,Nu respectively.Figure (5.1) Model Predictive Control [HS24]. Used with permission from©.Elsevier.2024.The first approach involves linearizing an analytical model at each oper-695.2. ADAPTIVE MPC FOR THE HEATING PHASEational point. This model, grounded in the heating process’s physical laws,allows for predicting the system’s future behavior by incorporating feedbackfrom the heaters’ surface temperatures. This feedback is crucial as it ensuresthat the control actions are based on the system’s current state, allowingreal-time adjustments to maintain the desired temperature profile.In contrast, the second approach adopts a data-driven Nonlinear Auto-Regressive with Exogenous inputs (NARX) model. This model is linearizedaround several specific operational points rather than continuously updatedoperational points, as in the analytical approach. The NARX model isdesigned to capture the system’s dynamics based solely on historical inputpower data to the heaters and the resulting zone temperatures on the sheet.Consequently, this approach eliminates the need for direct feedback fromthe heaters’ surface temperature, relying instead on the correlation betweenthe input power and the sheet temperature to predict system behavior.The core difference between these two methodologies lies in their relianceon feedback and the nature of the models used. The analytical model-basedapproach requires direct feedback from the sheet and heater temperature toadjust the control actions dynamically. In contrast, the data-driven NARXapproach leverages historical data patterns to inform its predictions, thusobviating the need for heater temperature feedback and only relying on thesheet temperature feedback.The following sections clarify each approach’s theoretical foundations,how it is implemented in the MPC framework, and how it affects controlperformance. I will evaluate the benefits and drawbacks of each approachthrough a comparative study, offering insights into their applicability forvarious operating circumstances during the thermoforming heating phase.The final objective is to find the best MPC approach that complies with thethermoforming process’s operational requirements and limitations.5.2 Adaptive MPC for the Heating PhaseWhile MPC has seen significant success in dealing with linear systems,nonlinear systems pose a greater challenge. Adaptive MPC is a specifictype of MPC designed to handle these nonlinear systems. The underlyingprinciple of adaptive MPC is utilizing an adjustable or adaptive model thatchanges over time or operating conditions. In the context of a nonlinearsystem, adaptive MPC typically employs a linearized model of the systemthat is updated. This essentially allows the controller to approximate thenonlinear system locally with a linear model, which can handle the specific705.2. ADAPTIVE MPC FOR THE HEATING PHASEdynamics and behaviors of the system at that instance. Nonetheless, thisprocess can lead to elevated computational complexity as control inputs needto be recalculated continuously. Furthermore, using a linearized model couldlimit the extended prediction horizon for the controller. This limitationarises because these models lose accuracy when deviating from the linearizedoperating point.5.2.1 Physics Based Model LinearizationThe nonlinear system in this study can be represented by Equations(5.1):xt+1 = F(xt, ut), (5.1)where x=[Tϕ]is the state vector and contains the temperature of the ele-ments on the sheet T and the average temperature of heaters (ϕ). u is theinput vector, including the heating power for heaters in the heater bank.F : Dx×Du → RNxNy+H , Dx ⊂ RNxNy+H , Du ⊂ RH , and F is a nonlinearfunctions representing the system dynamics. A linear approximation of thisnonlinear system can be represented by Equations (5.2):∆xt+1 = At∆xt +Bt∆ut, (5.2)Where ∆xt = xt − xop and ∆ut = ut − uop are the deviations of state andinput from the current operating points respectively. At andBt are the state-space matrices obtained by linearizing the system at the current operatingpoint (xop, uop). The dynamics of Tm,n, which is the ith component of vectorT defined in section 3.2, can be approximated by an ordinary differentialequation as:T˙m,n = α(q1 + q2 +1k(q3 + q4)) (5.3)Whereq1 =Tm,n+1+Tm,n−1−2Tm,n(∆x)2, q3 =ht+hbd (Tambient − Tm,n)q2 =Tm+1,n+Tm−1,n−2Tm,n(∆y)2, q4 =Ahεeσ∆x∆y∑Hh=1 Fh−m,n(ϕ4h − T 4m,n),α = kρcp is the thermal diffusivity of the material, ht and hb are the topand bottom face convective heat transfer coefficient [W/m2.K], Fh−m,n is theview factor between the hth heater and the element (m,n) on the sheet (see715.2. ADAPTIVE MPC FOR THE HEATING PHASEsection 3.1). The Jacobian matrices must be computed at each time steparound operating points to acquire the A and B1 matrices in the linearizedstate-space model. This involves calculating the partial derivatives of thefunction F concerning the state variables x(t) and the input variables u(t).Equations (5.4) and (5.5) show the A and B matrices in the linearized state-space model by considering the average time constant (τave) of 210 sec forthe heater various power levels.A =[∂F∂x](xop,uop)=[A11 A12A21 A22](xop,uop)(5.4)whereA11 = α(a4INxNy×NxNy − Aˆ11), Aˆ11 = A˜1 + A˜2 + A˜T2 ,A˜1 =(a7T3op1)1,1(a1)1,2 0 0 . . . 0(a1)2,1(a7T3op2)2,2(a1)2,3 0 . . . 00 (a1)3,2. . .. . . . . . 00 0. . .. . .. . .............. . .. . . (a1)NxNy−1,NxNy0 0 . . . 0 (a1)NxNy ,NxNy−1(a7T3opNxNy)NxNy ,NxNyA˜2 =[[0]Nx(Ny−1)×Nx a2INx(Ny−1)×Nx(Ny−1)[0]Nx×Nx [0]Nx×Nx(Ny−1)],A12 = 4αa5F1−1ϕ3op1 F2−1ϕ3op2· · · FH−1ϕ3opHF1−2ϕ3op1 F2−2ϕ3op2· · · FH−2ϕ3opH....... . ....F1−NxNyϕ3op1 F2−NxNyϕ3op2· · · FH−NxNyϕ3opHA21 =[0]H×NxNy , A22 = −1τaveIH×H , a1 = 1(∆x)2 ,a2 =1(∆y)2, a3 =ht+hbkd , a4 = −2a1 − 2a2 − a3,a5 =Ahεeσk∆x∆y , a6 =∑Hh=1 Fh−m,n, a7 = −4a5a6,B =[∂F∂u](xop,uop)=[B11B21](xop,uop), (5.5)B11 =[0]NxNy×H ,B21 =1τavediag(c5c6ec6(ui)op + c7c8ec8(ui)op)for i ∈ {1, 2, . . . ,H},1c5, c6, c7, and c8 can be obtained from table 3.1725.2. ADAPTIVE MPC FOR THE HEATING PHASE5.2.2 Optimization ProblemThe AMPC finds the optimal control input u(t) by solving the optimiza-tion problem in Equation (5.6) at each time step.minimizeu∑Npt=1ΠETt Et +∑Nc−1t=0 ΩuTt utsubject to xt+1 = Axt +But,Tt = Cxt,Et = T dt − Tt,Tt ∈ Y, ut ∈ U ,(5.6)where Π and Ω penalize the tracking error and rate of the control inputs,U and Y are admissible sets for inputs and outputs, ensuring the systemoperates within safe and efficient bounds. Np and Nc are the predictionand control horizons, respectively. ut is the control input at time step t.A, B, and C are system matrix, input matrix, and output transformationmatrix, respectively, defined in section 5.2.1 and 4.5.1. This optimizationis subject to the linearized system dynamics and constraints on the controlinput and state. In this study, the system utilizes 500 W ceramic heaters,which means the control input u(t) must satisfy 0 ≤ u(t) ≤ 500 W. Afterthe optimal control input is obtained, it is applied to the actual nonlinearsystem defined in 5.3, and the process repeats with the new system state.5.2.3 Simulation ResultFigure 5.3 visually represents the AMPC’s performance throughout a1000-second duration in controlling the desired temperature distribution, asoutlined in Figure 4.8. Considering the sampling time in the lab-scale setup,which is 6 sec, it is also set in the simulator. The optimal performance ofthe AMPC system was observed under specific conditions, specifically whenboth prediction and control horizons (Np and Nu) were configured at 24seconds each. Specifically, Figure 5.2 portrays the resulting temperaturedistribution on the sheet after 1000 seconds. Meanwhile, Figure 5.3 showsthe error variation across the 15 zones on the thermoplastic sheet during thecontrol process. The controller effectively brought the error of zones withinthe ±10 band after 695 seconds, with the average and maximum error ofzones at the control process’s conclusion being approximately 2.1 and ±5.5°C, respectively. The maximum overshoot observed was approximately 5.5735.2. ADAPTIVE MPC FOR THE HEATING PHASE0 0.2 0.4 0.6x [m]00.10.20.3y[m]80100120140160180Figure (5.2) Temperature profile on the sheet on the simulator.0 100 200 300 400 500 600 700 800 900 1000-100-500Error[/ C]Zone: 1Zone: 2Zone: 3Zone: 4Zone: 50 100 200 300 400 500 600 700 800 900 1000-100-500Error[/ C]Zone: 6Zone: 7Zone: 8Zone: 9Zone: 100 100 200 300 400 500 600 700 800 900 1000Time [sec]-100-500Error[/ C]Zone: 11Zone: 12Zone: 13Zone: 14Zone: 15Figure (5.3) Zone errors on the thermoplastic sheet.°C. Also, Figure 5.4 indicates the control inputs applied to the 15 heatingelements throughout the control process.Table 5.1 indicates substantial improvement compared to the perfor-mance achieved with PI control using the same reference temperature pro-file. Specifically, the average error band after the control process has beennearly halved to ±2.1°C. Additionally, AMPC reduced the overshoot by 22°C, achieving an overshoot of 5.5°C. Moreover, the settling time within anerror band of ±10°C is notably improved, with AMPC achieving 695 sec-onds. These advancements signify the effectiveness of AMPC in enhancingtemperature control performance over traditional PI control methodologies.As observed from the results, the computational cost of the AMPC ishigh, and this challenge intensifies as the number of heating elements in-creases, making the optimization problem more complex. Additionally, be-cause an analytical model is linearized at each step for use in MPC, it745.2. ADAPTIVE MPC FOR THE HEATING PHASE0 100 200 300 400 500 600 700 800 900 10000500Power[W] Heater: 1Heater: 2Heater: 3Heater: 4Heater: 50 100 200 300 400 500 600 700 800 900 10000500Power[W] Heater: 6Heater: 7Heater: 8Heater: 9Heater: 100 100 200 300 400 500 600 700 800 900 1000Time [sec]0500Power[W] Heater: 11Heater: 12Heater: 13Heater: 14Heater: 15Figure (5.4) Control input to 15 heating elements.Table (5.1) A comparison between PI and AMPC on thermoforming heat-ing phase simulator with 15 heaters for reference temperature profiles in 4.8.Metric PI AMPCReliance on model All thetimeAll thetimeRequires heater bank’s feedback Yes YesOffline training time [sec] - -Online computational effort [sec] 90 153Average error after 1000 sec [°C] 4 2.1Maximum error after 1000 sec [°C] ±8 ±5.5Overshoot [°C] 27.5 5.5Settling time within an error band of±10 °C [sec]850 695755.3. NARX-MPCrequires feedback on the surface temperature of the heaters at every step.Given hardware limitations in providing feedback from the heater bank andreducing the system’s order to decrease the optimization computational cost,I will develop a data-driven model in the next section. This model will di-rectly map the input heating power data to the output zone temperatures.5.3 NARX-MPCCreating and refining dynamic models are crucial and often the mostchallenging and time-consuming components of MPC implementation inindustrial settings [Ric93, Zhu98, ZPWZ13]. Regularly updating or re-identifying these models is a critical maintenance task within the MPCframework. Additionally, mastering the complexities of MPC theory, fine-tuning control strategies, and assessing performance are intricate processesthat demand a high level of expertise. Integrating Nonlinear Autoregressivewith Exogenous Inputs (NARX) models within MPC frameworks, known asNARX MPC, has further enhanced the capability to manage complex non-linear systems. NARX MPC offers several advantages, including improvedprediction accuracy for nonlinear systems, enhanced disturbance rejection,and better handling of multivariable systems with constraints. These ben-efits lead to more efficient and reliable control systems in industrial ap-plications, from process control in bioethanol fermentation to temperatureregulation in thermoforming and optimization in manufacturing processes[LNK+18, ZFLL24]. In manufacturing, applying NARX-MPC has demon-strated significant potential in stabilizing systems under various conditions,achieving performance that closely approaches the ideal. This utility under-scores the capability of NARX-MPC to optimize manufacturing operations,particularly in scenarios where the dynamics of the system are not fullyknown but where input-output data is readily available [SGM+21]5.3.1 Data-Driven Modeling: Nonlinear AutoregressiveExogenous (NARX)A NARX model extends the traditional linear ARX structure by in-corporating linear and nonlinear components in parallel; see Figure 5.5. Itconsists of two main components: a regressor vector and an output function.The regressor vector contains past input-output data, indicating the model’sorder. The output function involves offset, nonlinear, and linear functions,indicating the model’s structure. Since the output function benefits from aparallel structure, its nonlinear and offset parts are dedicated to reducing765.3. NARX-MPCFigure (5.5) a single-output NARX model block diagram.the modeling error that cannot be captured by the linear part. A NARXmodel with a wavelet nonlinear function can be represented as follows:Xt =[y[t−1:t−n], u[t−1:t−m]]⊤, (5.7a)yt(Xt) = XTt PL+W (Xt) + S(Xt) + ε, (5.7b)where Xt∈ R(n+m)×1 the regressor vector at time step t, ε is the offset, P∈R(n+m)×p is a projection matrix, transforming the regressors into a lower-dimensional space, L∈ Rp×1 is a linear mapping that scales the projectedregressors to the output. W (Xt) represents the sum of dilated and translatedwavelets, which are functions that capture the localized features of the signal.S(Xt) represents the sum of dilated and translated scaling functions, alsoknown as scales, which capture the global trends of the signal. The functionsW (Xt) and S(Xt) together form the wavelet network’s nonlinear component,allowing the NARX model to approximate complex, nonlinear relationshipsbetween the inputs and outputs [Zha97]. To linearize the NARX model,consider the following equation:Xt+1 = AXt + [B1 B2][ytut], (5.8a)yt = F(Xt, ut) (5.8b)where F is a nonlinear function. Xt and ut are state and the input vectorsat time t, respectively. To establish a linear state-space model around theoperating point yop = F(Xop, uop), the first derivatives of the function F withrespect to the state vector X and the input u are represented as follows:FX =∂∂XF(X,u)∣∣∣∣ u=uopX=Xop, Fu =∂∂uF(X,u)∣∣∣∣ u=uopX=Xop, (5.9)775.3. NARX-MPCUsing the derivatives, the changes in the output and state vector from theirvalues at the operating point can be approximated by the following linearequations:∆Xt+1 = A∆Xt + [B1 B2][∆yt∆ut]. (5.10a)∆yt = FX∆Xt +Fu∆ut. (5.10b)By substituting Equation (5.10b) in (5.10a), the linear approximation forthe state vector around the operational point is defined as follows:∆Xt+1 = (A+B1FX)∆Xt + (B1Fu +B2)∆ut (5.11)where ∆Xt = Xt−Xop, ∆ut = ut− uop, ∆yt = yt− yop are the deviation ofthe state vector, input, and output from their value at the operating point,respectively, B1 and B2 are matrices that relate the input and the output tothe state changes. This method is implemented by the “linearize” commandin MATLAB© - System Identification Toolbox.5.3.2 NARX-MPC for Heating Phase in ThermoformingSystemsThe main bottleneck of this method is to achieve a control-orientedmodel. Often, highly accurate and over-complicated models may not be suit-able for control design as the underlying optimization problem becomes ex-pensive and time-consuming, especially for high-dimensional systems. There-fore, we adopt the linearized NARX model, given in (5.10a-5.10b). For aservo-tracking problem2, we define the MPC problem as follows:minimizeu∑Npt=1ΠETt Et +∑Nc−1t=0 ΩuTt utsubject to Xt+1 = AXt + But,yt = CXt +Dut,Et = r − yt,yt ∈ Y, ut ∈ U ,(5.12)where Π and Ω penalize the tracking error and rate of the control inputs,U and Y are admissible sets for inputs and outputs, ensuring the systemoperates within safe and efficient bounds, A = A+B1FX , B = B1Fu +B2,C = FX and D = Fu are the matrices of the control-oriented model from2The reference, rt, is a constant signal over the prediction horizon.785.3. NARX-MPCFigure (5.6) MPC block diagram(5.10a). Np and Nc are the prediction and control horizons, respectively. utis the control input at time step t. The general architecture of the NARX-MPC method is depicted in Figure 5.6.The model inputs encompass the electrical power of the heaters withinthe heater bank, u, and the outputs comprise the average temperature ofeach zone on the thermoplastic sheet, T . Equation (5.13) shows the regressorvector for the NARX model proposed for thermoforming systems.Xt =[ Tt−1:t−nut−1:t−m], (5.13)where Tt−1:t−n ∈ RZ×n and ut−1:t−m ∈ RH×m are the past n and moutput and inputs at discrete-time index t, respectively. Z is the totalnumber of zones on the thermoplastics, and H shows the number of heaters.The system’s input (heater’s power) and output (sheet’s temperature) ut−nand Tt−n at time step t− n takes the following form:Tt−n =Tt−n(1)Tt−n(2)...Tt−n(Z) , ut−n =ut−n(1)ut−n(2)...ut−n(H) . (5.14)The range of values u(1),. . . ,u(H) is from 0 to 500 W, indicating the min-imum and maximum power levels that are applied to each heater, indicatinginput admissible set U in Equation 5.12. A data set including 20,000 sam-ples created the NARX model of the heating phase with 15 heating elements.Considering the setup’s 6-second sampling time, collecting the dataset fromthe experimental setup would take more than 33 hours. Moreover, usingseveral thermoplastic sheets was impractical due to the necessity to span795.3. NARX-MPCa large temperature range and the possibility of melting at higher tem-peratures. To get around this problem, the validated physics-based modelpresented in Chapter 3 was used to gather the required dataset.Figure (5.7) MPC algorithm flowchart in thermoforming process.The developed NARX model is integrated into the MPC framework topredict the zone temperatures on the thermoplastic sheet. As shown inFigure 5.7, this means using the initial prediction T (k+1) as input for sub-sequent predictions. The predicted temperature data is then processed inthe rolling optimization routine, which allows control signals to be generatedfor future and present effective temperature regulation. According to Equa-tion (5.12), the NARX-MPC algorithm iteratively solves an optimizationproblem to minimize the objective function J while considering constraintsand system dynamics. These limitations include the thermoplastic sheet’smaterial properties, particularly its melting temperature, and the heaters’physical limitations, such as their maximum and minimum power output.5.3.3 Data-driven Modeling for Thermoforming SystemsExperiment design3 for nonlinear systems requires a more profound in-tuition and physical understanding compared to linear systems, which pri-marily need to satisfy certain mathematical conditions such as persistentlyexciting (PE) inputs [NN20]. For a thermoforming system, we outline anexperiment to collect data for modeling purposes with the following criteria:I) The input should drive the system across a broad range of sheet temper-atures to expose the non-linearity of the system. II) The system’s outputshould oscillate around various equilibrium points to capture the system’s3Data collection for SysID purposes.805.3. NARX-MPCgain characteristics. III) Given the distinct nonlinear behaviors during ther-moforming systems’ warming and cooling phases, the input signal must bedesigned to explore these phases comprehensively. We employ a pseudo-random binary sequence (PRBS) signal with the upper and lower boundsincreasing every 200 minutes with the amplitude of 100 [W] until reaching amaximum of 500 [W], after which it decreases similarly; see Figure 5.8.0 200 400 600 800 1000 1200 1400 1600 1800 2000Time [min]0100200300400500Power [W]Figure (5.8) Experiment Design for Data Collection: Excitation Inputs.The parameters n and m in the regressor vector related to the number ofpast inputs and outputs are set to 2. This is because the dynamic relation ofheating power to heater temperature and radiation of heaters to the sheettemperature can be defined via two first-order models. The optimizationproblem to achieve the best parameters and number of wavelets is carriedout using MATLAB© - System Identification Toolbox.5.3.4 NARX Model EvaluationThe accuracy of multi-step-ahead predictions can degrade due to theaccumulation of one-step-ahead prediction errors during model recursion.Therefore, to assess the accuracy and robustness of the NARXmodel’s multi-step-ahead predictions, N-step ahead predictions are also evaluated againstthe test set; see Table 5.2. It should be noted that the model’s accuracy forthe N-step prediction is critical as it will be used in the MPC design. Theeffectiveness of this experimental design is evaluated through the model’saccuracy, quantified by the Normalized Root Mean Square Error (NRMSE)and presented in Figure 5.9.NRMSE =√1N∑Ni=1(yi − yˆi)2max(y)−min(y) , (5.15)where N is the number of data points, y is the noise-polluted measurement,and yˆ is the model’s predicted output. The blue and red lines display thedataset and model’s predicted outputs in Fig 5.9 for 1-step prediction. Ad-ditionally, the yellow line depicts the difference between these two values.815.3. NARX-MPCFigure 5.10 shows each zone’s auto-correlation of residuals. Since mostvalues lie in the 99% confidence intervals and there is a big spark on thecenter, we can conclude that the identified model has captured the essentialdynamics, and the residual errors are close to the measurement noise. Sam-ple cross-correlation in Figure 5.11 indicates how the residuals are linearlydependent on lagged inputs or how well the proposed model has exploredthe information in the lagged inputs. It shows that the residuals and laggedinputs have little correlation as their cross-correlation lies in the 99% con-fidence interval. In other words, Figure 5.10 and 5.11 evaluate the modelcapability of capturing linear dynamics and indicate how the remaining erroris similar to noise.Table 5.2 presents the model’s accuracy for each zone relative to the pre-diction horizon for test data, which was 25% of the entire data set. However,we optimized the NARX model for least one-step prediction error; it stillneeds to be evaluated for multi-step prediction. Since MPC relies on themodel’s accuracy given a long prediction horizon. For N = 100, the model’saccuracy is more than 80% on the test dataset. This is sufficient for MPC’sReceding-Horizon Control (RHC) technique4.4MPC with RHC results in an inherently robust implementation but note that the givenframework enables us to implement a Robust and Stochastic MPC to achieve superiorperformance and better constraint satisfaction.825.3.NARX-MPC1-Step Predicted Response Comparison0 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 10 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 20 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 30 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 40 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 50 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 60 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 70 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 80 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 90 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 100 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 110 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 120 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 130 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 140 500 1000 1500 2000Time [min]0500Temp [° C]02040Zone: 15Figure (5.9) 1-step temperature prediction results refer to comparing the temperature prediction and the actualtemperature in the fifteen zones (blue line: dataset, red line: NARX model’s output, yellow line: error).835.3.NARX-MPC00.51Zone1-20 -10 0 10 2000.51Zone2-20 -10 0 10 2000.51Zone3-20 -10 0 10 2000.51Zone4-20 -10 0 10 2000.51Zone5-20 -10 0 10 2000.51Zone6-20 -10 0 10 2000.51Zone7-20 -10 0 10 2000.51Zone8-20 -10 0 10 2000.51Zone9-20 -10 0 10 2000.51Zone10-20 -10 0 10 2000.51Zone11-20 -10 0 10 2000.51Zone12-20 -10 0 10 2000.51Zone13-20 -10 0 10 2000.51Zone14-20 -10 0 10 2000.51Zone15-20 -10 0 10 20Figure (5.10) Sample Auto-correlation of output yj , j ∈ {1, ..., 15} with 99% Confidence Intervals: Test dataset845.3.NARX-MPC-0.100.1Zone1-20 -10 0 10 20-0.100.1Zone2-20 -10 0 10 20-0.100.1Zone3-20 -10 0 10 20-0.100.1Zone4-20 -10 0 10 20-0.100.1Zone5-20 -10 0 10 20-0.100.1Zone6-20 -10 0 10 20-0.100.1Zone7-20 -10 0 10 20-0.100.1Zone8-20 -10 0 10 20-0.100.1Zone9-20 -10 0 10 20-0.100.1Zone10-20 -10 0 10 20-0.100.1Zone11-20 -10 0 10 20-0.100.1Zone12-20 -10 0 10 20-0.100.1Zone13-20 -10 0 10 20-0.100.1Zone14-20 -10 0 10 20-0.100.1Zone15-20 -10 0 10 20Figure (5.11) Sample Cross-correlation of input uj and output yj , j ∈ {1, ..., 15} with 99% Confidence Intervals:Test dataset855.3. NARX-MPCTable (5.2) Accuracy of the proposed NARX model for the test data:NRMSE for each Zone v.s. N-step closed-loop simulation prediction horizon.N 1 10 20 30 40 50 60 70 80 90 100Z: 1 98.6 97.6 96.0 94.4 92.9 91.5 90.1 88.9 87.8 86.7 85.7Z: 2 98.6 97.5 95.7 94.0 92.3 90.7 89.3 87.9 86.7 85.5 84.4Z: 3 98.5 97.4 95.7 93.9 92.2 90.7 89.2 87.9 86.7 85.5 84.4Z: 4 98.6 97.5 95.7 94.0 92.3 90.8 89.4 88.0 86.8 85.6 84.5Z: 5 98.6 97.3 95.4 93.5 91.7 90.0 88.5 87.0 85.6 84.4 83.2Z: 6 98.6 97.3 95.4 93.5 91.7 90.0 88.5 87.0 85.7 84.4 83.2Z: 7 98.6 97.4 95.6 93.8 92.1 90.5 89.0 87.6 86.3 85.2 84.1Z: 8 98.6 97.2 95.2 93.2 91.3 89.6 87.9 86.4 85.0 83.7 82.5Z: 9 98.6 97.3 95.3 93.4 91.5 89.8 88.2 86.7 85.3 84.0 82.8Z: 10 98.6 97.4 95.7 93.9 92.3 90.7 89.2 87.9 86.6 85.4 84.4Z: 11 98.6 97.3 95.3 93.3 91.5 89.7 88.1 86.6 85.2 84.0 82.8Z: 12 98.6 97.3 95.5 93.6 91.9 90.2 88.7 87.3 85.9 84.7 83.6Z: 13 98.5 97.5 95.9 94.2 92.7 91.2 89.9 88.6 87.4 86.3 85.3Z: 14 98.5 97.3 95.5 93.6 91.9 90.2 88.7 87.3 85.9 84.7 83.6Z: 15 98.5 97.4 95.8 94.1 92.5 91.0 89.6 88.3 87.1 86.0 85.05.3.5 NARX-MPC ResultsThis section will first evaluate the performance of the designed NARX-MPC on the simulator. Then, a comparative analysis will be conductedusing the PI and AMPC. Subsequently, the controller’s robustness will beexamined, focusing on three important parameters. Finally, the experimen-tal results obtained from the lab-scaled setup will be presented and comparedwith the outcomes observed in the simulator.Numerical ResultsTo evaluate the controller’s performance on the simulator and conduct acomparative analysis with PI and AMPC, the same reference temperaturedistribution on the thermoplastic sheet was selected according to Figure 4.8.Figure 5.13 shows the error variation across the 15 zones on the thermoplas-tic sheet during the control process. The controller effectively brought theerror of zones within the ±10 band after 560 seconds, with the average andmaximum error of zones at the control process’s conclusion being approxi-mately 0.7 and 1.4 °C, respectively. The maximum overshoot observed wasapproximately 2 °C. Also, Figure 5.14 indicates the control inputs applied865.3. NARX-MPC0 0.2 0.4 0.6x [m]00.10.20.3y[m]80100120140160180Figure (5.12) Temperature profile on the sheet at 1000 sec in the simulator.0 100 200 300 400 500 600 700 800 900 1000-100-500Error[/ C] Zone: 1Zone: 2Zone: 3Zone: 4Zone: 50 100 200 300 400 500 600 700 800 900 1000-100-500Error[/ C] Zone: 6Zone: 7Zone: 8Zone: 9Zone: 100 100 200 300 400 500 600 700 800 900 1000Time [sec]-100-500Error[/ C] Zone: 11Zone: 12Zone: 13Zone: 14Zone: 15Figure (5.13) Zone errors on the thermoplastic sheet0 100 200 300 400 500 600 700 800 900 10000200400Power[W] Heater: 1Heater: 2Heater: 3Heater: 4Heater: 50 100 200 300 400 500 600 700 800 900 10000200400Power[W] Heater: 6Heater: 7Heater: 8Heater: 9Heater: 100 100 200 300 400 500 600 700 800 900 1000Time [sec]0200400Power[W] Heater: 11Heater: 12Heater: 13Heater: 14Heater: 15Figure (5.14) Control input to 15 heating elements during control processto the 15 heating elements throughout the control process.875.3. NARX-MPCTable (5.3) A comparison between PI, AMPC, and NARX-MPC on ther-moforming heating phase simulator with 15 heaters for reference tempera-ture profiles in 4.8.Metric PI AMPC NARXMPCReliance on model All the time All the time All the timeRequires heater bank’s feedback Yes Yes NoOffline training time [sec] - - 900Online computational effort [sec] 90 153 15Average error after 1000 sec [°C] 4 2.1 0.8Maximum error after 1000 sec [°C] ±8 ±5.5 ±1.4Overshoot [°C] 27.5 5.5 3.8Settling time within an error bandof ±10 °C [sec]850 695 569Table 5.3 compares the performance of the designed NARX-MPC withPI and the AMPC across eight criteria. All methods rely on the systemmodel continuously. PI and AMPC require feedback from the surface tem-perature of heaters since they linearize a physics-based system model aroundoperational points. In contrast, NARX-MPC is developed based on input(heating power) and output (zone temperatures) data, eliminating the needfor heaters’ surface temperature feedback during operation. Regarding ac-curacy, the NARX-MPC method outperforms other methods, exhibiting atighter final error band at the end of the control process. Specifically, theerror band achieved with NARX-MPC is ±1.4, almost 4 times narrowerthan AMPC and almost 6 times narrower than PI. Utilizing the NARX-MPC method reduces overshoot, decreasing from 27.5 and 5.5 °C in PI andAMPC, respectively, to 3.8 °C in NARX-MPC. Furthermore, concerning set-tling time within an error band of ±10 °C, the NARX MPC can achieve thiserror band in 569 seconds. In contrast, the settling time for the AMPC andPI are 695 and 850 seconds, respectively. The lower online computationalcost of NARX-MPC is attributed to the simpler model with fewer inputsit uses compared to AMPC. This streamlined yet accurate model resultsin a reduced computational load. Additionally, the results show that thismodel achieves slightly better accuracy, suggesting that the MPC solver isless prone to getting stuck in local minima than when it uses the analyticalmodel.885.3. NARX-MPC0 500 1000 1500 2000 2500 3000 3500-100-500Error[/ C]0 500 1000 1500 2000 2500 3000 35000500Heaterstemp.[/ C]0 500 1000 1500 2000 2500 3000 3500Time [sec]0200400Power[W]Figure (5.15) The cycle-to-cycle control performance of the Narx-MPC inthe heating phase of the thermoforming process.Cycle-to-cycle Control AnalysisFigure 5.15 illustrates the performance of the proposed NARX-MPC ina cycle-to-cycle framework. It shows that the initial surface temperature ofthe heaters in each cycle is set to the final temperature from the previouscycle. Each cycle concludes once the temperature error across all zones isreduced to below 5 °C. At this point, the heated sheet is exchanged for a newone. A comparative analysis with an open-loop and PI control approach,depicted in Figure 4.9 and 4.15, highlights the efficiency of the NARX-MPC.Table 5.4 compares the performance of NARX-MPC with PI and open-loop control approaches in the cycle-to-cycle control framework. While theopen-loop and PI approaches managed only two and three cycles withinthe same timeframe, the NARX-MPC could handle approximately six cy-cles. Notably, the first cycle under NARX-MPC, starting from ambienttemperature, required 703 seconds to complete. Subsequent cycles weremore efficient, each taking 540 seconds, marking a significant improvementin cycle time compared to the open-loop system. This demonstrates theNARX-MPC’s enhanced performance in managing the cycle-to-cycle tran-sitions effectively.895.3. NARX-MPCTable (5.4) A comparison between open-loop, PI, and NARX-MPC incycle-to-cycle control approach.Method 1st cycletime [sec]2nd cycletime [sec]#cyclesOpen-loop 1880 1465 2PI 1405 1425 3NARX-MPC 703 540 6Robustness AnalysisThe robustness of the NARX-MPC controller was assessed by varyingthe following system parameters: the convection heat transfer coefficient (h),the gap between the heater bank and the thermoplastic sheet (d), and theradiation absorptivity of the thermoplastic sheet (α). As these parameterschange, the NARX-MPC (trained with initial values of h = 5,W/m2K,d = 15, cm, and α = 0.8) is evaluated for its performance metrics, such asmaximum overshoot and final average error across zones.The convection heat transfer coefficient on the top face of the sheet wasstudied since it has a higher value than the bottom face of a horizontal hotplate. The gap between the heater bank and the thermoplastic sheet is oneof the most important system parameters. Furthermore, during the heatingphase, high temperatures cause the sheet to lose its structural integrityand sag downwards under its weight. This sagging leads to different gapsbetween heaters in the heat bank and their associated zones during theheating phase. The radiation absorptivity of various thermoplastic sheetmaterials significantly impacts the NARX-MPC controller’s adaptability todifferent materials. As a result, it is critical to consider this parameterduring the analysis and evaluation stages.Figure 5.16 showcases the average error across zones and the maximumovershoot under variations in the three parameters for a scenario with auniform reference temperature distribution at 150 °C. Figures 5.16(a) and5.16(c) illustrate the zones’ average error at the end of the control processunder different plant parameter values. Also, Figures 5.16(b) and 5.16(d)show the maximum overshoot in the control process under the same condi-tions. The absorptivity coefficient of the sheet in the plant was kept at 0.8in Figures 5.16(a) and 5.16(b), and only (h) and (d) were studied. On theother hand, the gap between the heater bank and the thermoplastic sheetwas kept at 15 cm in Figures 5.16(c) and 5.16(d) and only (h) and (α) werestudied. The collective findings from Figures 5.16(a) and 5.16(c) demon-strate the robust performance of the NARX-MPC controller, even among905.3. NARX-MPC5 10 15h [ Wm2"K ]0.511.52Averageerror[o C] d=10 cmd=15 cmd=20 cm Nominal point(a)0 5 10 15h [ Wm2"K ]02468Overshoot[o C]d=10 cmd=15 cmd=20 cmNominal point(b)0 5 10 15h [ Wm2"K ]0.511.52Averageerror[o C] ,=0.5,=0.6,=0.7,=0.8,=0.9Nominal point(c)0 5 10 15h [ Wm2"K ]02468Overshoot[o C],=0.5,=0.6,=0.7,=0.8,=0.9Nominal point(d)Figure (5.16) Robustness evaluation of the designed controller on the sim-ulator (NARX model is trained with initial values of h = 5,W/m2K,d = 15, cm, and α = 0.8). (a) The average error of zones under varyingconvection heat transfer coefficients and gaps between the heater bank andthe sheet. (b) Maximum overshoot under varying convection heat transfercoefficients and gaps between the heater bank and the sheet. (c) The averageerror of zones under varying convection heat transfer coefficients and sheetabsorptivity coefficients. (d) Maximum overshoot under varying convectionheat transfer coefficients and sheet absorptivity coefficients.915.3. NARX-MPCsubstantial parameter variations, consistently maintaining the average er-ror of zones below 2 °C. Additionally, as illustrated in Figures 5.16(b) and5.16(d), within this range of system parameters, the maximum overshootremains below 7 °C. Notably, the convection heat transfer coefficient of theenvironment exhibits a reverse effect on the zone average error and the max-imum overshoot. While a higher convection heat transfer coefficient of theenvironment decreases the controller’s accuracy in terms of the zone averageerror, it results in a lower overshoot in the control process.Experimental ResultsThis section presents the practical application of the designed NARX-MPC in a lab-scale thermoforming setup, as detailed in section 3.3. Theambient temperature during the test was 21 °C. Two feasible reference tem-perature distributions according to Section 4.3, when the heater bank isplaced 15 cm above the thermoplastic sheet, were chosen to evaluate thecontroller’s performance, as shown in Figure 5.17. The first desired tem-perature profile tasks the NARX-MPC with achieving and maintaining auniform temperature of 70 °C across all zones. The second desired tempera-ture profile requires the controller to reach and sustain distinct temperaturesin each zone, as indicated in Figure 5.17(b).(a) uniform desired temperature profile. (b) non-uniform desired temperatureprofile.Figure (5.17) Reference temperature profiles for experimentsAs mentioned in section 3.3, an IR camera was used to measure thesystem output (zones’ average temperature on the thermoplastic sheet) andprovide feedback to the control algorithm. Figure 5.18 illustrates the mea-sured system’s output (temperature profile on the thermoplastic sheet) by IR925.3. NARX-MPC0 0.2 0.4 0.6x (m)00.10.20.30.40.5y(m)5060708090100(a)0 0.2 0.4 0.6x (m)00.10.20.30.40.5y(m)5060708090100(b)0 0.2 0.4 0.6x (m)00.10.20.30.40.5y(m)5060708090100(c)0 0.2 0.4 0.6x (m)00.10.20.30.40.5y(m)5060708090100(d)Figure (5.18) Experimental temperature profile measured by the IR cam-era: (a) t = 250 [sec], (b) t = 300 [sec], (c) t = 350 [sec], (d) t = 400 [sec]camera in multiple time stamps during the control process through NARX-MPC for the uniform reference temperature distribution. Figure 5.19 depictsthe error variation across the 15 zones on the thermoplastic sheet. The fig-ure shows that the controller effectively brought the error of all zones withinthe ±10 °C band after 297 seconds, with the average and maximum errorof all zones at the control process’s conclusion being approximately 0.4 and1.6 °C, respectively. The maximum overshoot observed was approximately4 °C. Additionally, Figure 5.20 showcases the control inputs applied to the15 heating elements throughout the control process.Figure 5.21 illustrates the measured system’s output (temperature pro-file on the thermoplastic sheet) by IR camera in multiple time stamps duringthe control process through NARX-MPC for the non-uniform reference tem-perature distribution. Figure 5.22 depicts the error variation across the 15935.3. NARX-MPC0 100 200 300 400 500 600 700 800 900-40-200Zoneerror[/ C]Zone: 1Zone: 2Zone: 3Zone: 4Zone: 50 100 200 300 400 500 600 700 800 900-40-200Zoneerror[/ C]Zone: 6Zone: 7Zone: 8Zone: 9Zone: 100 100 200 300 400 500 600 700 800 900Time [sec]-40-200Zoneerror[/ C]Zone: 11Zone: 12Zone: 13Zone: 14Zone: 15Figure (5.19) Experimental zone errors on the thermoplastic sheet for auniform reference temperature distribution.0 100 200 300 400 500 600 700 800 9000200400600Power[W] Heater: 1Heater: 2Heater: 3Heater: 4Heater: 50 100 200 300 400 500 600 700 800 9000200400600Power[W] Heater: 6Heater: 7Heater: 8Heater: 9Heater: 100 100 200 300 400 500 600 700 800 900Time [sec]0200400600Power[W] Heater: 11Heater: 12Heater: 13Heater: 14Heater: 15Figure (5.20) Control input to 15 heating elements during the control pro-cess for uniform reference temperature distribution.945.3. NARX-MPCzones on the thermoplastic sheet during the control process. The figureshows that the controller effectively brought the error of all zones withinthe ±10°C band after 308 seconds, with the average and maximum error ofall zones at the control process’s conclusion being approximately 1 and 3.8°C, respectively. The maximum overshoot observed was approximately 5.3°C. Additionally, Figure 5.20 showcases the control inputs applied to the 15heating elements throughout the control process.0 0.2 0.4 0.6x (m)00.10.20.30.40.5y(m)20406080100(a)0 0.2 0.4 0.6x (m)00.10.20.30.40.5y(m)20406080100(b)0 0.2 0.4 0.6x (m)00.10.20.30.40.5y(m)20406080100(c)0 0.2 0.4 0.6x (m)00.10.20.30.40.5y(m)20406080100(d)Figure (5.21) Experimental temperature profile measured by the IR cam-era: (a) t = 125 [sec], (b) t = 250 [sec], (c) t = 375 [sec], (d) t = 500 [sec]Table 5.5 compares the NARX-MPC controller’s performance betweena physics-based simulator and a real-world lab-scale thermoforming setup.The results demonstrate a close alignment in performance across both plat-forms when subjected to the same reference temperature distribution. No-tably, the settling times within the experimental error band of ±10◦C for955.3. NARX-MPC0 100 200 300 400 500 600 700 800 900-60-40-200Zoneerror[/ C]Zone: 1Zone: 2Zone: 3Zone: 4Zone: 50 100 200 300 400 500 600 700 800 900-60-40-200Zoneerror[/ C]Zone: 6Zone: 7Zone: 8Zone: 9Zone: 100 100 200 300 400 500 600 700 800 900Time [sec]-60-40-200Zoneerror[/ C]Zone: 11Zone: 12Zone: 13Zone: 14Zone: 15Figure (5.22) Experimental zone errors on the thermoplastic sheet for anonuniform reference temperature distribution.0 100 200 300 400 500 600 700 800 9000500Power[W] Heater: 1Heater: 2Heater: 3Heater: 4Heater: 50 100 200 300 400 500 600 700 800 9000500Power[W] Heater: 6Heater: 7Heater: 8Heater: 9Heater: 100 100 200 300 400 500 600 700 800 900Time [sec]0500Power[W] Heater: 11Heater: 12Heater: 13Heater: 14Heater: 15Figure (5.23) Control input to 15 heating elements during the control pro-cess.965.4. CHAPTER SUMMARYnonuniform and uniform reference temperature distributions were recordedat 308 seconds and 297 seconds, respectively. These times are significantlyshorter than those observed in the simulator, which recorded settling timesof 380 seconds and 378 seconds for the identical temperature distributions.Additionally, the overshoot in the real-world setup reached approximately5 °C and 4 °C for nonuniform and uniform temperature distributions, re-spectively, surpassing the overshoots recorded in the simulator for the cor-responding temperature profiles. This discrepancy suggests variations insystem parameters between the simulator and the actual experiment. Forexample, such differences could stem from a lower convection heat transfercoefficient in the lab than that assumed in the simulator or a higher radiationabsorptivity of the sheet than the values used in the simulations.Table (5.5) NARX-MPC validation.nonuniform uniformMetric Sim. Exp. Sim. Exp.Zones’ average error [°C] 0.6 1 0.4 0.8Zones’ maximum error [°C] ±2.6 ±3.8 ±1 ±1.6Overshoot [°C] 0.5 5.3 0.75 4Settling time within an errorband of ±10 °C [sec]380 308 378 2975.4 Chapter SummaryThis chapter has shown that using Model Predictive Control (MPC)techniques during the heating stage of thermoforming processes is feasibleand efficient. The study concentrated on two MPC approaches: based onthe analytical model approach and based on the data-driven model (NARX)approach. The MPC based on the analytical model allows the MPC to con-tinuously update its control actions based on the system states by utilizingfeedback from the heaters’ surface and sheet temperatures. On the otherhand, the NARX model employs historical input power data and consequentzone temperatures and offers predictive capabilities that rely only on sheettemperature feedback.A comparative study of their simulator performances shows that NARX-MPC outperforms AMPC in accuracy, overshoot control, and settling time.Furthermore, because the NARX-MPC algorithm utilizes a data-driven model,it does not require feedback from the heater bank temperatures, resultingin a lower computational load and simpler implementation in real-world se-tups. However, it is important to note that creating a new NARX model975.4. CHAPTER SUMMARYfor systems with more heaters requires a larger data set with more samples,leading to increased computational load. This impacts scalability and prac-tical implementation. In the next chapter, I will explore DRL algorithmsas powerful model-free methods to develop a controller that achieves highaccuracy and ensures scalability and ease of implementation.98Chapter 6Reinforcement LearningControl of the Heating Phasein Thermoforming Process6.1 IntroductionReinforcement Learning (RL) [SB18] is an area of research that combinesprinciples from two key domains: optimal control using value functions anddynamic programming with its inspiration of trial-and-error search. Theconcept of RL emerged from these two fields to address the optimal controlproblem in dynamical systems [NLH20b]. RL is a subfield of machine learn-ing that focuses on training agents to make sequential decisions throughinteraction with an environment. RL algorithms learn optimal policies bymaximizing cumulative rewards obtained from the environment [LSR18].Figure 6.1 outlines the general RL framework. The agent learns by re-ceiving environmental observations, rewards for its actions, and adjustingits behavior to maximize its cumulative reward over time. This RL pro-cess involves two fundamental components: Policy: The policy defines thestrategy the agent employs to choose actions in response to observationsfrom the environment. It’s essentially a decision-making function that mapsobservations to actions. Learning Algorithm: The learning algorithmupdates the policy based on the agent’s experiences.6.1.1 Markov Decision Process (MDP)In the conventional RL framework, an agent is modeled as a MarkovDecision Process (MDP), which is typically characterized by a five-tuple ⟨S, A, P, R, γ ⟩. Each element in this tuple plays a crucial role. We considersystems where states are denoted as ’s,’ and the set S encompasses all validstates. At discrete time intervals, represented as ’t,’ an agent observes thecurrent state, denoted as st, and promptly selects an action, at, from apredefined set, A. This action choice triggers a system transition, governed996.1. INTRODUCTIONFigure (6.1) Reinforcement Learning framework.by the probability set P = {p(st+1|st, at)}, resulting in a new state, st+1.As time progresses, a sequence of historical events is generated, representedas h = (s0, a0, s1, a1, . . .). Along with the state transition, the agent receivesa scalar reward, rt, because of its decision, and R = r(st, at); R ∈ R;S×A → R. The discount factor, denoted as γ, ranges between 0 and 1, and itdetermines the weight given to future rewards in the agent’s decision-making.With γ set to 0, the agent prioritizes immediate rewards, emphasizing short-term gains and disregarding long-term consequences. As γ increases, theagent increasingly considers future rewards, indicating a greater focus onthe long-term benefits of its actions. In this framework, given the currentstate st and action at, the next state st+1 is a random variable characterizedby the probability distribution p(st+1|st, at). typically, p, representing thesystem’s complex stochastic dynamics, is considered unknown.6.1.2 Reward FunctionThe reward function, denoted as R, is a critical component in reinforce-ment learning. It maps a state (or a state-action pair) to a numerical reward.The reward function can depend on the current state of the world (st), theaction just taken (at), and the next state of the world (st+1):rt = R(st, at, st+1) (6.1)However, in many cases, the reward function is simplified to depend onlyon the current state (rt = R(st)) or the state-action pair (rt = R(st, at)).This simplification can make the problem more tractable while allowing theagent to learn effective policies. The goal of an RL agent is to maximizesome notion of cumulative reward over a trajectory. This cumulative reward,1006.1. INTRODUCTIONor return (Gt), can be defined in several ways, depending on the problem’sspecifics:Finite-Horizon Undiscounted Return: This is the sum of rewardsobtained within a fixed window of steps. It is defined as:Gt =th∑t=0rt (6.2)where th is the time horizon.Infinite-Horizon Discounted Return: This return considers all fu-ture rewards but discounts them based on how far off they are. The discountfactor, γ, ranges between 0 and 1. The discounted return is defined as:Gt =∞∑t=0γtrt (6.3)The discount factor γ prioritizes immediate rewards over distant ones andensures the sum converges for infinite horizons. The choice between finite-horizon undiscounted return and infinite-horizon discounted return dependson the specific application and the nature of the environment. For example,tasks with a clear endpoint (e.g., games that end after a certain numberof moves) might naturally fit a finite-horizon model, while tasks requiringlong-term planning and consideration of future states might be better suitedfor an infinite-horizon discounted approach [ZMSZ20, XMW20, Sut99].6.1.3 PolicyThe agent’s actions are guided by a ”policy” (π), which is a state-dependent conditional probability distribution: π : S → P(A) over theaction set A. The goal of this policy, π, is to maximize the discountedcumulative reward, which is referred to as the return function:Gt =∞∑i=tγ(i−t)r(si, ai), (6.4)6.1.4 Value FunctionsThe value function is a fundamental concept that quantifies the expectedcumulative reward an agent can obtain from a given state while followinga particular policy. To assist the policy in making informed decisions, two1016.1. INTRODUCTIONtypes of value functions are being used: The state value function and action-value function, denoted as Vπ(s) and Qπ(s, a), respectively, and are definedas follows:State-Value Function (V (s)): The state-value function predicts theexpected cumulative reward an agent can achieve starting from a specificstate s and following a particular policy π. Mathematically, it is defined asthe expected sum of discounted future rewards:Vπ(s) = E[ ∞∑t=0γtr(st)∣∣∣∣∣s0 = s, π](6.5)Where Vπ(s) is the state-value function for state s. γ is the discount factorthat determines the importance of future rewards compared to immediaterewards. r(st) is the reward obtained at time step t when transitioning fromstate st−1 to state st under policy π. The sum is taken over all time steps t.Action-Value Function (Q(s, a)): The action-value function predictsthe expected cumulative reward an agent can achieve starting from a specificstate s, taking action a, and then following a particular policy π. Mathe-matically, it is defined as the expected sum of discounted future rewards:Qπ(s, a) = E[ ∞∑t=0γtr(st, at) | s0 = s, a0 = a, π](6.6)where Qπ(s, a) is the action-value function for state s and action a. r(st, at)is the reward obtained at time step t when transitioning from state st−1to state st and taking action at under policy π. Both the state-valueand action-value functions play a crucial role in RL algorithms. Theyhelp the agent make decisions by estimating the potential rewards of dif-ferent states and actions. By using value functions, agents can select ac-tions that lead to states with higher expected cumulative rewards, ulti-mately guiding them towards optimal decision-making in various environ-ments [YQZ+21, LKHH20, LFL+22, LHP+15].6.1.5 Bellman Equations in Reinforcement LearningThe Bellman equations are fundamental to understanding the dynamicsand optimization in reinforcement learning (RL). They provide a recursivedecomposition for the value functions associated with policies or actions,enabling the computation or approximation of these values through iterativemethods. The Bellman equations can be categorized based on whether they1026.1. INTRODUCTIONapply to on-policy or optimal value functions, each serving a distinct purposein the RL framework.On-Policy Value Functions: The Bellman equations for on-policyvalue functions describe the expected return for following a specific policyfrom a given state or state-action pair. They are defined as follows:− For the state-value function under policy π:V π(s) = Ea∼π,s′∼P [r(s, a) + γV π(s′)], (6.7)where s′ ∼ P (·|s, a) indicates that the next state s′ is sampled accord-ing to the environment’s transition dynamics given the current state sand action a, and a ∼ π(·|s) denotes that action a is chosen accordingto policy π given the current state s.− For the action-value function under policy π:Qπ(s, a) = Es′∼P [r(s, a) + γEa′∼π[Qπ(s′, a′)]], (6.8)where a′ ∼ π(·|s′) signifies that the next action a′ is chosen accordingto policy π given the next state s′.Optimal Value Functions (V ∗π (s), Q∗π(s, a)): The Bellman equationsfor optimal value functions describe the maximum expected return achiev-able from a given state or state-action pair, regardless of the policy beingfollowed. They incorporate a maximization step over actions to reflect theselection of the best possible action:− For the optimal state-value function:V ∗(s) = maxaEs′∼P [r(s, a) + γV ∗(s′)], (6.9)− For the optimal action-value function:Q∗(s, a) = Es′∼P [r(s, a) + γmaxa′Q∗(s′, a′)]. (6.10)The key distinction between the Bellman equations for on-policy andoptimal value functions lies in the maximization (max) operation over ac-tions in the optimal value functions. This maximization reflects the essenceof acting optimally, where the agent selects the action leading to the high-est value, thereby navigating the decision-making process towards the mostfavorable outcomes.1036.1. INTRODUCTION6.1.6 Learning Optimal PoliciesTwo primary approaches are commonly employed to address reinforce-ment learning problems: methods based on value functions, denoted asCritic-only methods, and those centered around policy search, referred toas Actor-only methods. Additionally, the actor-critic hybrid approach inte-grates both value functions and policy search. In the actor-critic framework,the actor and critic components are explicitly represented and learned inde-pendently.Critic-only (Learning based on value functions): Critic-only meth-ods are grounded in determining the optimal value function and deriving anoptimal policy from this foundation. These methods rely on estimating thevalue, or expected return, of being in a given state. The optimal state-valuefunction V ∗(s) is then defined as the maximum value across all possiblepolicies for a given state s:V ∗(s) = maxπV π(s) ∀s ∈ S (6.11)The optimal policy π∗ corresponds to the state-value function V ∗(s),where actions are chosen by maximizing the expected optimal value for thesucceeding state.Actor-only (Policy search): Policy search methods dispense withmaintaining a value function model and directly seek an optimal policyπ∗. This is achievable when the search space is constrained. typically, aclass of policies is parametrized by a real-valued parameters vector θ, andthese parameters are iteratively updated to maximize the expected returnE[R | θ] using either gradient-based or gradient-free optimization [DNP+13].Policy-based reinforcement learning is an optimization problem to find θthat maximizes the expected return. Some approaches, like Hill climbingand Genetic algorithms, do not rely on gradients, but policy gradients arecommonly used.Actor-Critic Method: The Actor-Critic method represents a hybridapproach that integrates both value functions and an explicit representationof the policy. In this approach, the actor, or policy, learns by incorporatingfeedback from the critic or value function. This method balances reducingthe variance of policy gradients and introducing bias from value functionmethods. The actor-critic framework is designed to optimize learning byleveraging the strengths of both policy and value-based methods, resultingin a more robust and effective approach for reinforcement learning.1046.2. DEEP REINFORCEMENT LEARNING6.2 Deep Reinforcement LearningIn recent years, reinforcement learning (RL) has attracted significantinterest from industry and academia due to its success in resolving com-plex problems. A recent groundbreaking achievement in deep learning hashinged upon the efficient training of deep neural networks on extensivedatasets. This advancement has significantly enhanced various domains,such as computer vision and speech recognition. These models are adeptat direct training from raw inputs, employing stochastic gradient descent toupdate the networks’ weights iteratively. This success has paved the way forreinforcement learning (RL) innovations, where a synergy emerges betweenRL algorithms and deep neural networks. This integration enables these al-gorithms to operate directly on raw images and efficiently generate trainingdata through stochastic gradient updates. The intersection of reinforce-ment learning and deep learning has thus opened new avenues for tacklingcomplex problems and learning directly from raw, unprocessed inputs. Theglobal interest in RL from 2007 to 2022 is depicted in Figure 6.2. Interestin RL was fairly steady until 2015 when it significantly increased. The firstfactor contributing to this growth was the 2013 publication of Mnih et al.[MKS+13], in which the authors introduced a deep RL algorithm, a deepQ-learning network (DQN) capable of playing numerous Atari games usingonly image inputs.Figure (6.2) Growth in search results of reinforcement learning on Googlefrom 2007 to 2022 [Rl-].1056.2. DEEP REINFORCEMENT LEARNING6.2.1 Taxonomy of Reinforcement Learning AlgorithmsReinforcement learning algorithms can be classified from various perspec-tives, including model-based vs. model-free, value-based vs. policy-based(or a combination), Monte Carlo vs. temporal-difference, and on-policy vs.off-policy methods, as illustrated in Figure 6.3. Despite their differences,these approaches share crucial characteristics inherited from deep reinforce-ment learning concepts. The control implemented by reinforcement learn-ing acts as a closed-loop system, where the reward represents the system’sfeedback. This reward is delayed, prompting many algorithms to minimizethis delay. Furthermore, deep reinforcement learning algorithms feature se-quential decision-making behavior, where long-term rewards depend on thesequence of actions taken. The credit assignment problem describes thetime-dependent nature of deep reinforcement learning implementations, assome actions show consequences after a delay and many intermediate systemstates.Model-Based andModel-Free: Model-based DRL involves construct-ing an explicit model or representation of the environment. This model isthen utilized to simulate potential scenarios and plan actions accordingly.These approaches aim to enhance sample efficiency and decision-making byincorporating a model. However, accurately modeling complex environmentscan be challenging. In contrast, model-free DRL bypasses the constructionof an explicit model and directly learns optimal actions based on interactionswith the environment. This category includes methods that focus on learn-ing a policy or value function without relying on a predefined understandingof the environment’s dynamics.Value-Based and Policy-Based: Value-based DRL methods concen-trate on estimating the value function, which represents the expected cu-mulative reward for a given state or state-action pair. Algorithms in thiscategory, such as Deep Q Networks (DQN), aim to learn the optimal valuefunction and deduce an optimal policy from it. Policy-based DRL, on theother hand, seeks to learn the optimal policy directly. Instead of estimatingvalue functions, these methods learn a parameterized policy, defining thestrategy for selecting actions. Policy gradient methods fall into this cate-gory, with algorithms like Trust Region Policy Optimization (TRPO) andProximal Policy Optimization (PPO).Monte Carlo and Temporal Difference: Monte Carlo methods esti-mate the value of states or state-action pairs by averaging sampled returnsover multiple trajectories. These approaches update value estimates basedon the observed returns, providing a robust but high-variance learning ap-1066.2. DEEP REINFORCEMENT LEARNINGRL AlgorithmsMDP BanditAction-ValueMethodsGradient BanditMethodsModel-BasedModel-FreePolicy-BasedGradient-freeValue-BasedGradient-basedOff-PolicyTRPO/PPOPolicy GradientACKTRCross-EntropyMethodSAMUELEvolutionStrategyQT-OptDQNQ-LearningACDDPGA2C/A3CSACTD3Double DQNDueling DQNC51SarsaOn-PolicyGiven the ModelLearn the ModelMCTSI2AWorld ModelFigure (6.3) Map of reinforcement learning algorithms (DP (DynamicProgramming), TD (Temporal Difference), MC (Monte Carlos), I2A(Imagination-Augmented Agent) [WRR+17], DQN (Deep Q-Network)[MKS+13], TRPO (Trust Region Policy Optimization) [SLA+15], ACKTR(Actor-Critic using Kronecker-Factored Trust Region), AC (Actor-Critic),A2C (Advantage Actor-Critic), A3C (Asynchronous Advantage ActorCritic) [MBM+16], DDPG (Deep Deterministic Policy Gradient) [LHP+15],TD3 (Twin Delayed DDPG) [FHM18], SAC (Soft Actor-Critic) [HZAL18],PPO (Proximal Policy Optimization) [SWD+17].)TD Methods:TD LearningQ-LearningSARSAMaximum Entropy MethodsPolicy ItterationValue ItterationGeneralized Policy ItterationDirect Policy Search:Policy GradientTrust RegionEvolutionActor-CriticPolicyValueFigure (6.4) Value-based and policy-based methods [ZY20].1076.2. DEEP REINFORCEMENT LEARNINGproach. On the other hand, temporal Difference (TD) learning updates valueestimates by the difference between predicted and observed returns at eachtime step. This approach combines Monte Carlo and dynamic programmingelements, offering a lower-variance alternative suitable for online learningscenarios.UpdateMont-Carlo Methods Temporal-Difference MethodsPolicy GradientTRPO/PPO DDPGActor-CriticDQNQ-LearningFigure (6.5) Monte Carlo methods and temporal-difference methods[ZY20].On-Policy and Off-Policy: On-policy methods involve learning a pol-icy and updating it based on the current policy’s performance. Algorithmslike REINFORCE fall into this category, optimizing the policy directly fromthe experienced trajectories generated by the current policy. Off-policymethods, such as Q-learning, decouple the optimized policy from the datacollection policy. These algorithms maintain a separate exploration policy,allowing them to learn from past experiences and improve the policy inde-pendently.Understanding these categorizations provides a comprehensive overviewof the diverse strategies employed within the deep reinforcement learninglandscape, aiding researchers and practitioners in selecting the most suit-able approaches for their specific applications. In this study, the ProximalPolicy Optimization (PPO) [SWD+17] algorithm is chosen for DRL dueto its demonstrated ability to achieve robust performance with a limitednumber of samples, a crucial factor given the resource-intensive nature ofinteracting with the thermal system. The efficiency in sample utilizationprovided by PPO aligns well with the complexities of the thermal controlproblem at hand1086.2. DEEP REINFORCEMENT LEARNING6.2.2 Proximal Policy Optimization (PPO)PPO has gained popularity due to its effectiveness in handling contin-uous and discrete action spaces and achieving stable policy updates. PPObalances exploration and exploitation, allowing the agent to learn from theenvironment while ensuring policy improvements. To estimate the policyand value function, a PPO agent employs two function approximators; Ac-tor: The actor is responsible for mapping the current observation (state)st to a probability distribution over the action space. Critic: The criticestimates the value function, reflecting the expected cumulative reward anagent can obtain from a particular state st while following a given policy.The value function considers the long-term rewards obtained by followingthe current policy. This estimation helps the agent to evaluate and com-pare different states based on their potential for achieving higher rewards[SWD+17].Figure (6.6) Actor–Critic PPO algorithm process flow chart [HS24]. Usedwith permission from ©.Elsevier.2024.Figure 6.6 represents a general perspective of the PPO algorithm. θ andµ are parameters of the policy (actor) and value function (critic) networks,respectively. In this algorithm, Aˆt is the advantage function described byEquation 6.12.Aˆt = rt + γVπθ(st+1;µ)− Vπθ(st;µ) (6.12)1096.2. DEEP REINFORCEMENT LEARNINGIn PPO, the objective functions differ for the critic (value function) andactor models. The critic, designed to estimate the cumulative long-termreward when adhering to the policy from the state corresponding to thepresent observation, operates by receiving the current observation as inputand generating a single scalar output. The objective function for the criticmodel in PPO aims to minimize the difference between the predicted statevalue Vµ(st) and the target value Vˆtargetµ (st). The objective function LV (µ)is defined as Equation (6.13) [GCCW21]:LV (µ) = Et[∣∣∣Vˆ targetµ (st)− Vµ (st)∣∣∣] , (6.13)Here, LV represents the empirical average over a batch of samples, and| . . . | denotes the absolute difference. Et represents the expectation takenover time steps within a single episode or trajectory. The actor model inPPO optimizes a surrogate objective function that incorporates a probabilityratio and a clipping mechanism to ensure stability during policy updates.The objective function LCLIP(θ) is given by Equation (6.14) [LCYW19]:LCLIP(θ) = Et [min(rt(θ)At, clip(rt(θ), 1− ϵ, 1 + ϵ)At)] (6.14)Here, LCLIP(θ) represents the empirical average over a batch of samples.At is the advantage function, ϵ represents the clipping parameter a hyperpa-rameter controlling the degree of policy update, and rt(θ) =πθ(at|st)πθold (at|st)is thelikelihood ratio that compares the probability of selecting action at underthe current policy πθ with the probability of selecting the same action underthe previous policy πθold . The parameters of the policy πθ are updated us-ing SGD with the gradient ∇LCLIP for the negative of the clipped objectivefunction according to Equation 6.15 [WNW22]. , where we work with thenegative of this objective function to maximize the policy’s performance.The clipped objective function LCLIP avoids making PPO overly greedy infavoring actions with positive advantages and enables it to quickly avoidactions with a negative advantage function from a mini-batch of samples.θ ← arg maxθLCLIP(θ) (6.15)Algorithm 1 provides a detailed explanation of the PPO algorithm.1106.3. DEEP REINFORCEMENT LEARNING FOR THERMOFORMING SYSTEMAlgorithm 1: Proximal Policy Optimization (PPO)Input: Initial parameters for actor and critic networks (θ and µ)and clipping threshold (ϵ)for k = 1, 2, . . . doCollect set of trajectories Dk = {τi} by running policyπk = π(θk) in the environment for T times;Compute reward-to-go Gˆt;Compute advantage estimates, Aˆt based on the current valuefunction Vµk ;Update the policy by optimizing the clipped surrogate objectivefunction (see (13)):θ ← arg maxθLCLIP (θ)Fit value function by regression on mean-squared error:µ← arg maxµLV (µ)whereLV (µ) =1|Dk|T∑τ∈DkT∑t=0(Vµ(st)− Gˆt)2,end6.3 Deep Reinforcement Learning forThermoforming SystemThe design of a Deep Reinforcement Learning (DRL) agent for the heat-ing phase of the thermoforming process involves five main task before thetraining process:− Defining the observation vector− Defining the action vector− Defining the reward function− Selecting an appropriate DRL algorithm− Selecting appropriate neural networks architecture in DRL1116.3. DEEP REINFORCEMENT LEARNING FOR THERMOFORMING SYSTEMThe PPO algorithm is suitable for the heating phase of thermoforming dueto its sample efficiency and ability to learn effective policies with fewer in-teractions with the environment. Additionally, PPO balances explorationand exploitation, essential for optimizing the heating process. The systemneeds to explore different heating strategies to find the optimal solution.Only the feedback from sheet temperatures is considered to simplify thealgorithm and facilitate real-world implementation. Hence, the observationvector comprises the error and derivative of the error for each zone. Theactuators in this system are heating elements with a maximum power con-sumption of 500 [W]. The PPO algorithm can handle continuous and discreteaction spaces. A discrete action space with multiple power levels is chosento align with the nature of ceramic heaters, where their surface temperaturecan be controlled by adjusting between different power level intervals.The reward function should be accurately defined to satisfy the crite-ria mentioned in Section 4.2. One of the critical control objectives is thehandling of overheating, which is essential for preventing material degra-dation and ensuring product quality. To enhance production efficiency, thetrained agent must minimize the rise time and settling time during the heat-ing phase. Furthermore, it must minimize the terminal error at the end ofeach control horizon, ensuring that the absolute difference between the de-sired setpoint temperature and the actual temperature of each zone is lessthan 5°C. The trained agent must be robust to sensitive system parametersand maintain optimal control performance despite external disturbances, asindustrial environments present numerous challenges. The control system’sscalability is also crucial, especially in industrial settings where thermoform-ing processes may involve hundreds of heaters. Hence, the trained agentshould be easily deployable to other settings.Due to the complexity of the plant model, as the number of heaters in-creases in the heater bank, traditional DRLs cannot converge during thetraining process. As illustrated in Table 6.1, the complexity of the actionspace grows exponentially with the number of heaters (H) in our system, par-ticularly when considering our simplified discrete action space, which repre-sents two heater statuses: On and Off. Experience has shown that achievingconvergence with PPO becomes challenging as the number of heaters in-creases, and extended training durations are necessary. For instance, whendealing with 15 heaters, the estimated training time for the DRL agent us-ing the PPO algorithm to manage all heaters effectively amounts to 113days. To mitigate this prolonged training time and address convergencechallenges, particularly in the early stages of training, I introduce two ap-proaches in this thesis: I) a hybrid approach integrating Model Predictive1126.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEMControl (MPC) within the PPO framework, and II) using multiple agentsto handle each heater in the heater bank. In the next sections, these twoapproaches will be investigated in detail, and the results will be comparedto previously designed model-based algorithms.Table (6.1) Estimated DRL agent training time versus the number of theheaters in the heater bank.Numberof HeatersAction Space(2H)ActualTraining TimeEstimatedTraining Time1 2 10 min −4 16 80 min −6 64 − 5.3 hr9 512 − 42.4 hr12 4096 − 14 days15 32768 − 113 days6.4 MPC-Guided DRL Approach inThermoforming SystemThe action at can be either sampled from the policy or provided by theMPC. In the case of applying the output of the MPC, the action can berepresented as Equation (6.16):at = MPC(st) + εt, (6.16)where at is the selected action at time step t, st is the current state, and εtrepresents exploration noise. In RL, exploration refers to taking actions thatmay be suboptimal in the short term but provide new information aboutthe environment, leading to better long-term performance. The MPC con-troller provides an initial action selection, often considered a good startingpoint. However, MPC’s selected action is augmented with exploration noiseto encourage exploration further. Exploration noise to the MPC’s selectedaction introduces stochasticity in the RL training process. This stochastic-ity helps the RL agent explore different actions and learn more about theenvironment, which can be beneficial for finding better policies. Accordingto Figure (6.7), during the RL training process with the augmented actionat, the agent receives feedback from the environment through rewards and1136.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEMupdates its policy based on the observed rewards and actions. The explo-ration noise ensures the agent explores diverse actions, even if the MPCcontroller’s initial action selection is relatively deterministic. The action atis sampled by combining the policy and the MPC actions with a certainprobability p. With probability p, the MPC provides the action at, whereas,with probability 1− p, the action at is drawn from the policy.Figure (6.7) Flowchart of MPC-guided approach for DRL training [HS24].Used with permission from ©.Elsevier.2024.Equation (6.17) shows the probability function used in this study tosample action between MPC and policy output.p(episode) =C1(C2 + eC3(episode−C5))2[C4(1− Rˆ) + Rˆ(C2 + eC3(episode−C5))],(6.17)In the given equation, ”episode” represents the episode number in the train-ing process. The function Rˆ generates a random number between 0 and1 at each episode. By referring to the parameters listed in Table 6.2, wecan determine the probability distribution of using action at during a 3500-episode training process, as illustrated in Figure 6.8. From Figure 6.8, it isevident that initially, the network heavily relies on the actions generated bythe MPC. However, as the training progresses, the network gradually shiftsits dependence on the policy’s output. From episode number 2500 onward,the network solely relies on the policy’s output. This transition indicatesthat the network progressively improves by incorporating more informationfrom its learned policy.1146.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEMTable (6.2) Probability function parameters.Parameters ValuesC1 1C2 1C3 0.0025C4 0.51C5 2500Figure (6.8) Distribution of Action selection probability during trainingprocess[HS24]. Used with permission from ©.Elsevier.2024.This approach accelerates the DRL training process by providing high-quality action inputs from the MPC algorithm, particularly during the initialstages of training when the RL policy is still in its learning phase and facedwith the challenges posed by extensive action and state spaces. The bestperformance of the MPC system was observed when the sample time was setto 5 seconds. By discretizing the MPC’s continuous output into on and offsignals, these signals were then synchronized with the RL agent’s decision-making process, which operates at a faster sample time of 1 second.6.4.1 Observation and Action SpaceA discrete action space was chosen in the PPO approach to simplifythe RL algorithm. Two statuses are considered for each heater: ON with500 W power consumption and OFF, denoted by a binary vector at =[a1, a2, ..., aH ]T , where ai ∈ {0, 500W} and H represents the number of1156.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEMheaters. This binary vector represents the action taken by the agent to con-trol the power consumption of each heater. This decision was influencedby the nature of ceramic heaters, where their surface temperature can beeffectively controlled by adjusting on and off intervals. In the RL statespace, the error and the derivative of the error for the system’s output inEquation (4.5) were considered. The error is the deviation from the desiredtemperature distribution on the sheet that the controller aims to achieve,represented as e = T − Td. Hence, st is a vector with 2Z components,where Z is the number of zones on the sheet. The state at time t is denotedas st = [et, e˙t]T , where et = [e1, e2, ..., eZ ] represents the vector of errorsat zones, and e˙t = [e˙1, e˙2, ..., e˙Z ] represents the vector of derivative of er-rors at zones. The state space consists of these error and derivative values,which the RL algorithm uses to make decisions and improve the controller’sperformance.6.4.2 Reward FunctionThe reward function (R) in this study is composed of three main parts,denoted as R1, R2, and R3 according to Equation (6.18). It guides the RLagent towards achieving the system designer’s objectives.R = R1 +R2 +R3 (6.18)Figure (6.9) Error bands and overheating area[HS24]. Used with permis-sion from ©.Elsevier.2024.R1 (Initial Error Penalty): The first part of the reward function, R1,is designed to account for the initial conditions and initial errors of the zones1166.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEMon the sheet. The error vector is e = [e1, e2, ..., eZ ], where Z is the totalnumber of zones and ei denotes the error in the i-th zone. If the error in anyzone is smaller than −30 ◦C, the RL agent will receive a constant reward of−10000 at each time step.This component of the reward function pushes the agent to raise thesheet’s temperature when the error remains below -30 °C. Mathematically,the first part of the reward function is defined as Equation (6.19):R1(e) ={ −10000, ifei < −30◦C, ∀i = 1, 2, . . . ,Z0 otherwise(6.19)By providing a large negative reward, the agent is motivated to take ac-tions during the initial stages that effectively reduce the errors, leading tothe cessation of negative rewards and progress toward achieving the desiredtemperature distribution on the sheet.R2 (Accuracy Rewards): The second component of the reward func-tion, denoted as R2, is aimed at motivating the agent to maintain errorswithin tighter thresholds. As shown in Equation (6.20) and visually repre-sented in Figure 6.9, the agent’s rewards increase when errors are constrainedwithin smaller bands. As time progresses and errors are consistently lim-ited to these narrower bands, the denominator in Equation (6.20) decreases,leading to higher rewards. The highest possible accuracy reward is attainedwhen errors across all zones approach zero(M=0); at this point, the agentreceives a reward of 105 at each time step.R2(e) =1051 +M , if ei > −30◦C, ∀i = 1, 2, . . . , Z (6.20)where M = max |ei|, ∀i = 1, 2, . . . , ZR3 (Overheat Penalty): The third part of the reward function, R3,is crucial in penalizing the agent whenever a zone’s temperature exceedsthe desired value. Its purpose is to prevent the agent from overheatingany zone, which could lead to undesirable consequences such as melting thethermoplastic sheet. Mathematically, R3(e) can be expressed as Equation(6.21).R3(e) = −100Z∑i=1ei, if ei > 0◦C, ∀i = 1, 2, . . . , Z (6.21)In this equation, the agent receives a negative penalty at each time stepbased on the sum of positive errors in zones where the temperature ex-ceeds the desired value. This discourages the agent from taking actions that1176.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEMcould lead to overheating and emphasizes the importance of maintainingthe desired temperature distribution on the thermoplastic sheet. Moreover,the magnitude of the negative reward increases as more zones experienceoverheating, further deterring excessive heating actions by the agent.0 100 200 300 400 500 600 700 800Time [sec]-1-0.500.511.52Reward#104Figure (6.10) Reward output[HS24]. Used with permission from©.Elsevier.2024.This approach incentivizes the agent to act cautiously and accurately,preventing potential issues like melting the thermoplastic sheet. Ultimately,R3(e) significantly contributes to achieving better control of the thermo-forming process by guiding the agent’s behavior in a responsible and precisemanner. Figure 6.10 displays the output of the designed reward function forthe output error scenario depicted in Figure 6.9, spanning a duration of 800seconds.6.4.3 Agent TrainingThe network structures for the actor and critic models in PPO are de-signed, adhering to the specifications shown in Table 6.3. Each model con-sists of two hidden layers, with 512 neurons in each layer. The connectionsbetween these layers are established using a rectified linear unit (ReLU) acti-vation function, except for the output layers of the actor and critic networks.Considering the nature of ceramic heaters, where their surface temperaturecan be effectively controlled by adjusting on and off intervals, discrete actionspace in PPO was chosen. As a result, a SoftMax activation function waschosen for the actor network’s output layer. The use of SoftMax allows theoutput values to fall within the range of 0 to 1, providing a practical way todetermine the optimal combination of On and Off statuses out of 2H states1186.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEM(H: number of heaters) for the heater bank at each time step. Moreover, alinear activation function was used for the critic’s output layer since it repre-sents a value estimate rather than a probability distribution. This approachensures precise control over the power levels of the heaters, enabling efficientthermal control in the thermoforming process. The detailed structures ofthe implemented actor and critic networks are presented in Table 6.3.Table (6.3) Network Structures of the implemented artificial neural net-works for PPO.Network Neurons of Layers Activation FunctionsActor Network 30; 512; 512; 215 ReLU; ReLU; SoftMaxCritic Network 30; 512; 512; 1 ReLU; ReLU; LinearDuring the RL training process, the desired temperature distribution onthe thermoplastic sheet was established using a random uniform temperaturewithin the range of 50 to 200 °C during each training episode. As previouslydiscussed, the core parameters of the neural network were determined basedon the specifications outlined in Table 6.3. Furthermore, the hyperparam-eters governing the behavior of the PPO algorithm were carefully selectedand are presented in Table 6.4.Table (6.4) Hyperparameters.Hyperparameters Symbol ValuesThe learning rate for actor-network αac 1e− 6The learning rate for critic-network αcr 1e− 6Batch size |B| 128Experience horizon H 512Entropy loss weight β 0.02Number of training epochs per update N 3Clip parameter ϵ 0.2Generalized Advantage Estimation (GAE) factor GAE 0.95Discount factor γ 0.991196.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEMLearning process for a single-heater systemFigure 6.11 illustrates the DRL episode reward during the training pro-cess when the system comprises only one heater controlling the temperatureof a zone. The figure shows that the agent was effectively trained in approx-imately 1000 episodes, equivalent to about 600 seconds. The fluctuation atthe end of the training period is about the different reference temperaturesthat the agent is randomly assigned to control at each episode. With fewervariables and interactions, the learning process becomes more straightfor-ward and converges rapidly.0 200 400 600 800 1000Episode-50510EpisodeReward#106Episod RewardMoving AverageFigure (6.11) DRL episode reward during training process when the sys-tem comprises only one heater controlling the temperature of a specific zoneLearning process for 15-heater systemFigure 6.12 illustrates the episode rewards observed during the agenttraining process using the PPO algorithm for 3,500 episodes. This agent issupposed to learn how to control 15 heaters to regulate sheet temperature toa uniform temperature set point. Since no guidance by MPC is provided, thiscase illustrates the performance of conventional PPO. The depicted rewardsin the graph, as governed by the reward function outlined in section 6.4.2,span from −12×107 to −1×107 and the moving average is around −4×107.While the algorithm converges for systems with up to four heaters (see table6.1), we see a significant lack of overall progress in episode rewards in the caseof 15 heaters, despite an extensive training period spanning 3,500 episodes,and the algorithm fails to reach convergence. This leveling of rewards canbe primarily attributed to the inherent complexity intrinsic to the system.This complexity is further influenced by the wide action and state space1206.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEMinvolved in the process, making exploration difficult in industrial settingswhere the number of heaters can reach 200.0 500 1000 1500 2000 2500 3000 3500Episode number-12-10-8-6-4-202Episodereward#107Episode RewardMoving AverageFigure (6.12) Episode reward during traditional DRL (PPO) training pro-cess for 15 heaters[HS24]. Used with permission from ©.Elsevier.2024.Figure 6.13 provides a visual representation of the episode rewards recordedduring the training process utilizing the PPO algorithm with the integrationof the MPC-guided approach over a duration of 3500 episodes. The under-lying methodology involves the sampler function block depicted in Figure6.7, which operates based on the probability function outlined in Equation(6.17). For the initial 1000 episodes, the agent’s actions are driven primar-ily by MPC’s output. Figure 6.13 shows that the average episode rewardstabilizes around 1× 106 during this phase. Upon passing the 1000-episode,the agent’s decision-making is influenced by a combination of MPC’s outputand the policy’s output, determined by the probability function describedin Equation (6.17). In this transition phase, there is a notable decline inepisode rewards, reflecting the performance of the policy’s output. How-ever, as the training process advances and temporal progression occurs, theepisode rewards exhibit growth. Notably, around episode 2500, when theagent relies only on the policy’s output, the average episode reward expe-riences an increase compared to the former approach in Figure 6.12 wherethe average episode reward was −4 × 107. This trend persists, and by theculmination of the 3500 episodes, the average episode reward attains a valueof −0.04× 107.After successfully training the agent to effectively control the sheet’stemperature to achieve a desired uniform temperature within the range of50 to 200 °C across the entire sheet, the trained agent can tackle the taskof controlling nonuniform temperature profiles within the same temperature1216.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEMFigure (6.13) Episode reward during DRL (PPO) training process withMPC guided approach for 15 heaters[HS24]. Used with permission from©.Elsevier.2024.range. A modest number of additional training episodes without an MPC-guided approach for any specific desired temperature profile on the sheetcan increase the accuracy and reliability of the agent.6.4.4 Simulation ResultsThe previously trained agent underwent an additional 600 episodes of re-training. Figure 6.14 illustrates how the rewards per episode changed duringthis retraining phase. Notably, the performance of the pre-trained agent toachieve the desired temperature profile, as depicted in Figure 4.8, improvedafter just 635 episodes of this retraining, and the average episode rewardincreased from −0.38× 107 to −0.3× 107.1226.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEM0 100 200 300 400 500 600Episode number-5-4.5-4-3.5-3-2.5Episodereward#106Episode RewardMoving AverageFigure (6.14) Episode reward during retraining process[HS24]. Used withpermission from ©.Elsevier.2024.Figure 6.16 demonstrates the DRL’s performance over a 1000-secondduration in achieving the desired temperature distribution outlined in Figure4.8 with specifications mentioned in Table 6.3 and Table 6.4. Since theaction space has been chosen as on and off commands to heaters, the RLagent has trained with a higher sampling time of 1 second compared to theexperimental setup sample time, which is 6 sec.Figure 6.15 visually represents the MPC-guided DRL’s performance through-out a 1000-second duration in controlling the desired temperature distribu-tion, as outlined in Figure 4.8. Figure 6.15 portrays the resulting temper-ature distribution on the sheet after 1000 seconds. Meanwhile, Figure 6.16shows the error variation across the 15 zones on the thermoplastic sheet dur-ing the control process. The controller effectively brought the error of zoneswithin the ±10 band after 765 seconds, with the average and maximum errorof zones at the control process’s conclusion being approximately 1.8 and ±8°C, respectively. The maximum overshoot observed was approximately 11.8°C. Also, Figure 6.17 indicates the control inputs applied to the 15 heatingelements throughout the control process.Table 6.5 compares the performance of PI, AMPC, NARX-MPC, andMPC-guided DRL. As shown, the online computational load of MPC-guidedDRL is very low. Regardless of the capacity of the heater bank in the numberof heating elements and the system’s complexity, this online computationalload remains consistent once the agent is trained. After 3500+600 trainingepisodes, the average and maximum absolute errors at the end of 1000 sec-onds for MPC-guided DRL are 1.8 and ±8 °C, respectively. These valuesare higher than NARX-MPC, with errors of 0.8 and ±1.4 °C. The table also1236.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEM0 0.2 0.4 0.6x [m]00.10.20.3y[m]80100120140160180Figure (6.15) Temperature profile on the sheet at 1000 sec in the simulator.0 100 200 300 400 500 600 700 800 900 1000-100-500Error[/ C] Zone: 1Zone: 2Zone: 3Zone: 4Zone: 50 100 200 300 400 500 600 700 800 900 1000-100-500Error[/ C] Zone: 6Zone: 7Zone: 8Zone: 9Zone: 100 100 200 300 400 500 600 700 800 900 1000Time [sec]-100-500Error[/ C] Zone: 11Zone: 12Zone: 13Zone: 14Zone: 15Figure (6.16) Zone errors on the thermoplastic sheet.0 100 200 300 400 500 600 700 800 900 10000500Power[W] Heater: 1Heater: 2Heater: 3Heater: 4Heater: 50 100 200 300 400 500 600 700 800 900 10000500Power[W] Heater: 6Heater: 7Heater: 8Heater: 9Heater: 100 100 200 300 400 500 600 700 800 900 1000Time [sec]0500Power[W] Heater: 11Heater: 12Heater: 13Heater: 14Heater: 15Figure (6.17) Control input to 15 heating elements.1246.4. MPC-GUIDED DRL APPROACH IN THERMOFORMING SYSTEMreveals that the performance of MPC-guided DRL in handling overshoot issuboptimal, with the sheet experiencing an overshoot of 11.8 °C during thecontrol process. Additionally, the settling time within an error band of ±10°C is 765 seconds, which is higher than the results from MPC but lower thanthose from PI.The results in Table 6.5 demonstrate that although MPC-guided DRLsignificantly reduces the online computational load, improving the agent’saccuracy requires additional training episodes or a more complex neuralnetwork during offline training, which requires more powerful hardware con-figuration. Furthermore, despite its low online computational load, it canface a high offline computational load when dealing with systems with manyheaters. The next chapter will propose a multi-agent DRL approach to ad-dress these challenges.Table (6.5) A comparison between PI, AMPC, NARX-MPC, and MPC-guided DRL on thermoforming heating phase simulator with 15 heaters forreference temperature profiles in 4.8.Metric PI AMPC NARXMPCMPCGuidedDRLReliance on model All thetimeAll thetimeAll thetimeOnly fortrainingRequires heater bank’sfeedbackYes Yes No NoOffline training time [sec] - - 900 60135Online computationaleffort [sec]90 153 15 11Average error after 1000sec [°C]4 2.1 0.8 1.8Maximum error after1000 sec [°C]±8 ±5.5 ±1.4 ±8Overshoot [°C] 27.5 5.5 3.8 11.8Settling time within anerror band of ±10 °C [sec]850 695 569 7651256.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM6.5 Multi-Agent DRL Approach inThermoforming SystemThis section presents the development and application of multi-agentDRL to control radiation thermal systems primarily used during the heatingphase of thermoforming processes involving multiple heaters. The complex-ity of such control systems is characterized by large action and state spaces,where the actions of all actuators collectively influence the system’s out-put. This complexity introduces substantial challenges regarding the com-putational demands for offline exploration and training of learning-basedalgorithms and the online computational costs associated with a real-worldcontroller deployment. The section presents a novel approach to training anadaptive and robust RL agent system that can control a given single heatingelement on the thermoplastic sheet while dynamically considering interac-tive effects from nearby heaters. Results demonstrated that upon deployingthe pre-trained agent for each heating element within the heater bank, thegroup of agents could then regulate the temperature of the thermoplasticsheet to any physically feasible output temperature profile.(a) Single Agent. (b) Multi Agent.Figure (6.18) Single and Multi-Agent RL in heating phase of thermoform-ingFigure 6.19 shows the PPO-based architecture algorithm used in thisstudy. The PPO algorithm employed in this study follows an Actor-Critic(AC) methodology as explained in section 6.2.2. This approach involvesthe concurrent approximation of the policy and value functions, realizedthrough two distinct neural networks: the actor (θ) and the critic (µ). Theactor network determines the appropriate action to take, while the role ofthe critic model is to learn and evaluate whether the action taken by the1266.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMactor model leads to an improved state. The feedback from the critic modelis then utilized to optimize the actor model.Figure (6.19) The RL architecture in radiation thermal control based onPPO with Actor-Critic styleThe neural network configuration, including the number of hidden lay-ers and neurons, was determined through an iterative trial and error, draw-ing upon experimentation and prior experience in similar problem domains.This approach was employed to balance model complexity and computa-tional efficiency, ensuring the neural networks effectively captured the dy-namics of the thermoforming control problem. The network architecturesemployed in the actor and critic models are detailed in Table 1, showcas-ing a configuration of two hidden layers, each containing 512 neurons andconnections established using a rectified linear unit (ReLU) activation func-tion. The actor network’s output layer incorporates a SoftMax activationfunction tailored to the action space definition in section 6.5.1. This choicefacilitates the generation of output values within the range of 0 to 1, offeringa practical means of determining optimal values from the action space ateach time step.1276.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMTable (6.6) Network Structures of the implemented artificial neural net-works for PPO in multi-agent RL approach.Network Neurons of Layers Activation FunctionsActor Network 15; 512; 512; 6 ReLU; ReLU; SoftMaxCritic Network 15; 512; 512; 1 ReLU; ReLU; Linear6.5.1 Observation and Action SpaceOur multi-agent model aims to exhibit robust control over zone temper-atures in the presence of disruptive radiation from neighboring heaters. Toachieve this, I have introduced an observation mechanism that captures thedynamic behavior of the heaters. Furthermore, a set of boundary conditionsfor heaters positioned in zones around the sheet’s periphery has been de-fined to ensure the trained agent’s applicability across the entire heater bank.This approach enhances the adaptability and effectiveness of the trained RLagent in real-world scenarios, as shown in the results section. Figure 3.13illustrates the portion of the absorbed radiation each element receives from aceramic trough heater (FTE-500, Ceramicx Company) located 15 cm abovethe sheet. While each heating element emits radiation across the entire ther-moplastic sheet, according to the specific radiation pattern of the heater inFigure 3.13 and given the dimensions of the zones on the sheet (measur-ing 10 cm x 24 cm for trough ceramic heaters), the heater’s radiation ismost prominent within its own designated zone and adjacent vertical andhorizontal side zones.In constructing the RL observation vector, a set of variables crucial forachieving precise temperature control amidst the influence of neighboringheaters was selected. This comprehensive array includes the error and itsderivative in the central zone, providing vital insights into the system’s cur-rent state and directional trends. The strategic inclusion of input power datafrom the central and adjacent heating elements’ preceding step addressesthe considerable time constant of the heaters, serving as a valuable resourcefor preventing the overshoot in temperature on thermoplastic sheets. Ex-panding our observation set to encompass the temperature and derivativeof temperatures in adjacent zones reflects cues for the algorithm to considertheir impact on the central zone. Collectively, these considerations estab-lish a holistic observational framework, empowering the algorithm to makewell-informed decisions, particularly in scenarios characterized by extendedtime constants. Consequently, the state vector of the agent is structured as1286.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMFigure (6.20) Central and adjacent zones definition; ZCe: central zone,Zu: up neighboring zone, Zd: down neighboring zone, Zr: right neighboringzone, Zl: left neighboring zone.st = [e(1×1); e˙(1×1); T(4×1); T˙(4×1); a(5×1)]. Here, the error (e) is defined as thedeviation from the desired temperature distribution on the sheet that thecontroller aims to achieve in the central zone (Zce). The temperature vec-tor T = [Tu; Tl; Tr; Td] captures the average temperatures in adjacent zones(Zu: up adjacent, Zl: left adjacent, Zr: right adjacent, Zd: down adjacent)according to Figure 5, and the input power vector a = [ace; au; al; ar; ad]represents the power inputs to the heating elements. This comprehensivestate vector equips the RL agent with the necessary information to makeinformed and adaptive decisions in the thermal control process.As illustrated in Figure 6.21, there are three distinct conditions to con-sider regarding zones around the sheet’s periphery. We must account forthese conditions when constructing the state vector for the agents in thesezones. For any missing zones, constant values should be introduced into thestate vector to maintain completeness and consistency, as detailed in Table6.7. This ensures that the agent has a well-defined and comprehensive set ofobservations for all relevant zones, whether within the primary conditionsor under the specified constant conditions.1296.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM(a) Corner zone. (b) Horizontal sidezone.(c) Vertical side zone.Figure (6.21) Zones around the sheet’s peripherySix discrete power statuses characterize action space as follow:ai ∈ {0, 100 W, 200 W, 300 W, 400 W, 500 W}.The decision to employ a discrete approach stems from the intricaciesassociated with continuous spaces, which pose challenges in terms of train-ing complexity. The number of discretized power levels in the action spacebalances control precision and training intricacy. While a higher number ofpower levels enhances control precision, it concurrently escalates the com-plexity of the training process since the RL algorithm needs to explore andlearn the value function for each action, which can lead to combinatorialgrowth in the number of possibilities. This can increase the time and re-sources required for training. In our case, the choice of six power levels alignsstrategically with the hardware’s sample time of 6 seconds. The steady-statesurface temperatures of the heater at these power levels are 21 ◦C, 223.5 ◦C,344.7 ◦C, 418.9 ◦C, 472.4 ◦C, and 530 ◦C, respectively. To reach any temper-ature between two consecutive values, the heater must switch its power levelwith a specific frequency between two related power levels. The maximumerror within the 6-second hardware sample time for the heater to keep anytemperature between two consecutive steady-state surface temperatures are2.8 ◦C, 0.8 ◦C, 0.7 ◦C, 1.6 ◦C, 1.7 ◦C, respectively. This indicates that inthe worst-case scenario, with a 6-second hardware sample time and an ac-tion space characterized by six discrete power statuses, any heater’s surfacetemperature can be obtained with a maximum error of 2.8 ◦C. This dis-cretization effectively captures the system dynamics within the given timeframe, ensuring a balance between control precision and training efficiency.1306.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMTable (6.7) State vector in zones around the sheet’s periphery.Zone T T˙ aCornerUp left [Ta, Ta, Tr, Td] [0, 0, T˙r, T˙d] [0, 0, ar, ad]Upright [Ta, Tl, Ta, Td] [0, T˙l, 0, T˙d] [0, al, 0, ad]Down left [Tu, Ta, Tr, Ta] [T˙u, 0, T˙r, 0] [au, 0, ar, 0]Downright [Tu, Tl, Ta, Ta] [T˙u, T˙l, 0, 0] [au, al, 0, 0]Horizontal sideUp [Ta, Tl, Tr, Td] [0, T˙l, T˙r, T˙d] [0, al, ar, ad]Down [Tu, Tl, Tr, Ta] [T˙u, T˙l, T˙r, 0] [au, al, ar, 0]Vertical sideLeft [Tu, Ta, Tr, Td] [T˙u, 0, T˙r, T˙d] [au, 0, ar, ad]Right [Tu, Tl, Ta, Td] [T˙u, T˙l, 0, T˙d] [au, al, 0, ad]6.5.2 Reward FunctionThe reward function defined in section 6.4.2 is employed here. The firstpart of the reward function is designed to push the agent to elevate thetemperature of the sheet to decrease the rise time and bring the error towithin the ±30°C range. Failing to achieve this objective results in the RLagent receiving a consistent penalty of -10000 at each time step. The secondpart of the reward function encourages the agent to maintain errors withineven stricter thresholds. Conversely, the third and final part penalizes theagent whenever the temperature in a zone exceeds the desired value. Itsprimary function is to discourage the agent from overheating any specificzone.6.5.3 Agent TrainingAn iterative approach was undertaken to optimize the performance ofthe RL model in the context of the thermoforming process in determin-ing the hyperparameters for the PPO algorithm. The values were selectedbased on domain knowledge, experimentation, and fine-tuning through mul-tiple trial-and-error iterations. The chosen hyperparameters were tailoredto strike a balance between exploration and exploitation, ensuring stabilityand convergence during the training process. Consideration was given tothe thermoforming system’s specific characteristics, the control task’s com-plexity, and the dynamics of the physical lab-scale setup. Table 6.8 lists thehyperparameters employed in this study.1316.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMTable (6.8) Hyperparameters.Hyperparameters Symbol ValuesThe learning rate for actor-network αac 1e− 4The Learning rate for critic-network αcr 1e− 4Batch size |B| 128Experience horizon H 512Entropy loss weight β 0.02Number of training epochs per update N 3Clip parameter ϵ 0.2Generalized Advantage Estimation (GAE) factor GAE 0.95Discount factor γ 0.99During the RL agent’s training process to regulate the zone’s tempera-ture on the sheet, accounting for radiation disturbances from neighboringheaters, the heater bank was positioned 15 cm above the sheet. In RL, anepisode is a single, self-contained sequence of interactions between an agentand its environment. It represents a unit of experience during which theagent observes the state of the environment, takes actions, receives feed-back in the form of rewards or penalties, and learns from this experience toimprove its decision-making strategy. In this study, each training episodespanned 1000 seconds, and the target temperature for the designated zonewas randomly selected within the range of 50 to 200 °C for each episode.According to Equation 6.22, random input power levels were also applied toadjacent heating elements during an episode.u ={randi([0, 500], 4, 1) if(rem(⌊time210⌋, 2)= 0)[0]4×1 otherwise(6.22)This training approach allowed the agent to learn the impact of neigh-boring heaters on the central zone during thermal control. Figure 6.22(a)shows the generated random input power to adjacent heaters during anepisode of 1000 sec by the Equation 6.22, and Figure 6.22(b) shows the ad-jacent heaters’ surface temperature variation along the same period. Theserandom inputs are used in each episode during the DRL agent training toaccount for radiation disturbances from neighboring heaters. Figure 6.23shows the zone error, reward, and the input power to the central heaterduring one training episode.1326.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM0 200 400 600 800 1000 1200Time [s]0100200300400500Heatingpower[W]HuHlHrHd(a) random input power to adjacentheaters.0 200 400 600 800 1000 1200 1400Time [s]0100200300400Surfacetemperature[/ C] HuHlHCeHrHd(b) heaters’ surface temperature varia-tion.Figure (6.22) Random input power to adjacent heaters during agent train-ing0 200 400 600 800 1000 1200Time [s]-100-80-60-40-20020Error[/ C](a) zone error during one episode ofDRL agent training.0 200 400 600 800 1000 1200Time [s]-20246810Reward#104(b) reward function during one DRLagent training episode.0 200 400 600 800 1000 1200Time [s]0100200300400500Heatingpower[W](c) heating power during one DRLagent training episode.Figure (6.23) DRL training during only one episode1336.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMThus, a group of these trained agents can be deployed within a heaterbank to control the temperature distribution across the entire sheet collec-tively. Figure 6.24 depicts the episode rewards of the RL agent throughoutthe training process using the reward function defined in section 6.5.2. Thegraph clearly illustrates the convergence of the PPO algorithm in approx-imately 700 episodes where the ’average episode reward’ reaches 7 × 106.Across various training scenarios, it was observed that this value is typicallyreached between 350 and 700 episodes. Notably, fluctuations in episoderewards during the late stages of training are attributed to differences inreference temperatures. Controlling higher reference temperatures requiresadditional time to minimize errors, yielding slightly lower episode rewards.The training process shown in Figure 6.24 is similar to the single heaterresults shown in the 6.4.3.0 200 400 600 800 1000Episode-50510Episodereward#106Episod RewardMoving AverageFigure (6.24) DRL agent’s episode reward during the training process foronly one heater.6.5.4 Simulation ResultsFigure 6.25 exemplifies the performance of the trained RL agent in track-ing a reference temperature over a duration of 1800 seconds. The referencetemperature undergoes two changes during the control period, starting at100 °C, decreasing to 60 °C, and concluding at 120 °C. As depicted in Fig-ure 6.25(a), the trained RL agent closely tracks the reference temperatureand maintains the zone temperature in line with the desired temperature,exhibiting an error of less than 2 °C.1346.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM0 500 1000 1500Time [s]050100Temperature[/ C]Zone Ave. TemperatureRef. Temperature(a) reference temperature tracking on simulator.0 500 1000 1500Time [s]0100200300400500Heatingpower[W](b) heating power.Figure (6.25) RL agent’s performance on the simulatorFigure 6.27 demonstrates the DRL’s performance over a 1000-secondduration in achieving the desired temperature distribution outlined in Fig-ure 4.8 with specifications mentioned in Table 6.6 and Table 6.8. Figure6.26 visually represents the Multi-agent DRL’s performance throughout a1000-second in controlling the desired temperature distribution, as outlinedin Figure 4.8. Figure 6.26 portrays the resulting temperature distributionon the sheet after 1000 seconds. Meanwhile, Figure 6.27 shows the errorvariation across the 15 zones on the thermoplastic sheet during the controlprocess. The controller effectively brought the error of zones within the ±101356.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM0 0.2 0.4 0.6x (m)00.10.20.3y(m)80100120140160180Figure (6.26) Temperature profile on the sheet after 1000 sec in the simu-lator.0 100 200 300 400 500 600 700 800 900 1000-100-500Error[/ C] Zone: 1Zone: 2Zone: 3Zone: 4Zone: 50 100 200 300 400 500 600 700 800 900 1000-100-500Error[/ C] Zone: 6Zone: 7Zone: 8Zone: 9Zone: 100 100 200 300 400 500 600 700 800 900 1000Time [sec]-100-500Error[/ C] Zone: 11Zone: 12Zone: 13Zone: 14Zone: 15Figure (6.27) Zone errors on the thermoplastic sheet.band after 348 seconds, with the average and maximum error of zones atthe control process’s conclusion being approximately 0.3 and ±1 °C, respec-tively. The maximum overshoot observed was approximately 6.5 °C. Also,Figure 6.28 indicates the control inputs applied to the 15 heating elementsthroughout the control process.Observing the system exhibits lower overshoot as the reference temper-atures increase is intriguing. This phenomenon can be attributed to thehotter surface of the thermoplastic sheet, which enhances thermal dissipa-tion through convection and radiation, consequently reducing the overshoot.Figure 6.29 illustrates the control system’s overshoot for various uniform ref-erence temperature distributions on the thermoplastic sheet, ranging from50 to 200 °C. As depicted, the highest overshoot occurs at a reference tem-1366.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM0 100 200 300 400 500 600 700 800 900 10000500Power[W] Heater: 1Heater: 2Heater: 3Heater: 4Heater: 50 100 200 300 400 500 600 700 800 900 10000500Power[W] Heater: 6Heater: 7Heater: 8Heater: 9Heater: 100 100 200 300 400 500 600 700 800 900 1000Time [sec]0500Power[W] Heater: 11Heater: 12Heater: 13Heater: 14Heater: 15Figure (6.28) Control input to 15 heating elements during the control pro-cess.perature of 50 °C, while the overshoot approaches zero for desired referencetemperature distributions exceeding 160 °C.50 100 150 200Uniform ref. temperature [oC]0246Max.oveshoot[o C]Figure (6.29) Maximum overshoot obtained from the simulator in differentuniform reference temperature distributions on the sheet.The performance of the present multi-agent DRL approach is comparedto the PI, adaptive MPC, Narx-MPC, and MPC-guided DRL approacheson the thermoforming heating phase simulator. The evaluation was con-ducted using the identical reference temperature distribution in Figure 4.81376.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMto compare the effectiveness of these approaches. Table 6.9 compares theperformance of these methods on the simulator. As depicted in Table 6.9,the offline training time of the multi-agent DRL approach requires only 545seconds compared to 60135 seconds for the MPC-guided single-agent DRLapproach. Regarding control precision, the multi-agent DRL achieves andmaintains errors within a narrower margin of ±1°C. Furthermore, the av-erage error of all zones at the end of the control process is only 0.3 °C. Itcan be seen that regarding overshoot control, the NARX-MPC outperformsother methods with 3.8 °C. At the same time, the Multi-agent DRL ex-perienced 6.5 °C. Regarding online computational load, Multi-agent DRLperforms as fast as MPC-guided DRL, requiring only 11 sec for the wholecontrol process. Comparing the settling time within an error band of ±10°C, Multi-agent DRL outperforms other methods with 348 sec to bring theerror of zones in this band.1386.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMTable (6.9) A comparison between PI, AMPC, NARX-MPC, MPC-guidedDRL, and Multi-agent DRL on thermoforming heating phase simulator with15 heaters for reference temperature profiles in 4.8.Metric PI AMPC NARXMPCMPCGuidedDRLMulti-agentDRLReliance onmodelAll thetimeAll thetimeAll thetimeOnly fortrainingOnly fortrainingRequires heaterbank’s feedbackYes Yes No No NoOffline trainingtime [sec]- - 900 60135 545Online comput.effort [sec]90 153 15 11 11Ave. error after1000 sec [°C]4 2.1 0.8 1.8 0.3Max. error after1000 sec [°C]±8 ±5.5 ±1.4 ±8 ±1Overshoot [°C] 27.5 5.5 3.8 11.8 6.5Settling timewithin an errorband of ±10 °C[sec]850 695 569 765 3486.5.5 Cycle-to-cycle Control AnalysisFigure 6.30 illustrates the performance of the designed multi-agent DRLin a cycle-to-cycle framework. The initial surface temperature of the heatersin each cycle is set to the final temperature from the previous cycle. Eachcycle concludes once the temperature error across all zones is reduced tobelow 5 °C. At this point, the heated sheet is exchanged for a new one. Acomparative analysis with open-loop control, PI, and NARX-MPC methodsin a cycle-to-cycle control approach highlights the efficiency of the multi-agent DRL.1396.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM0 500 1000 1500 2000 2500 3000 3500-100-500Error[/ C]0 500 1000 1500 2000 2500 3000 35000500Heaterstemp.[/ C]0 500 1000 1500 2000 2500 3000 3500Time [sec]0500Power[W]Figure (6.30) The cycle-to-cycle control performance of the multi-agentDRL in the heating phase of the thermoforming processTable 6.10 compares the performance of Multi-agent DRL with NARX-MPC, PI, and open-loop control approaches in the cycle-to-cycle controlframework. While the open-loop, PI, and NARX-MPC approaches managedonly two, three, and six cycles within the same timeframe, the Multi-agentDRL could handle approximately eleven cycles. Notably, the first cycleunder Multi-agent DRL, starting from ambient temperature, required 409seconds to complete. Subsequent cycles were more efficient, each taking282 seconds, marking a significant improvement in cycle time compared tothe other methods. This demonstrates the Multi-agent DRL’s enhancedperformance in managing the cycle-to-cycle transitions effectively.Table (6.10) A comparison between open-loop, PI, NARX-MPC, andMulti-agent DRL in cycle-to-cycle control approach.Method 1st cycletime2nd cycletime#cyclesOpen-loop 1880 1465 2PI 1405 1425 3NARX-MPC 703 540 6Multi-agent DRL 409 282 111406.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM6.5.6 Robustness AnalysisThe robustness of the proposed multi-agent DRL strategy was validatedby studying the sensitivity to three model parameters: the convection heattransfer coefficient (h), the distance between the heater bank and the sheet(d), and the sheet radiation absorptivity (α). The DRL agent was trainedwith fixed parameter values (h = 5W/(m2K), d = 15 cm, and α = 0.8) andthen evaluated in terms of maximum overshoot during the control processand the final average error of zones at the end of the control process, withconsiderations for different values of these three parameters. The sensitivityanalysis on variations in the convection heat transfer coefficient was specif-ically considered on the sheet’s top face. This choice was made due to thetypically higher values of the convection heat transfer coefficient on the topface than the bottom. The distance between the heater bank and sheet,ranging from 10 cm to 20 cm in industrial settings, introduces distinct ra-diation patterns on the sheet. The radiation absorptivity of the sheet canvary depending on the sheet material, and it is important to investigate theRL agent behavior when facing different types of materials. Suppose theagent’s performance remains favorable under these altered circumstances.In that case, it instills optimism regarding its effectiveness when deployedin an actual setup, potentially eliminating the need for additional trainingin a real-world environment. Figure 6.31 presents the average error acrosszones after the control period and the maximum overshoot when confrontedwith variations in the three parameters for the same scenario as in Figure(A uniform temperature distribution of 150 °C on the sheet).In Figure 6.31(a) and Figure 6.31(b), the heat transfer coefficient wasvaried for different distances between the sheet and the heaters while keep-ing the absorptivity at 0.8. In Figure 6.31(c) and Figure 6.31(d), the heattransfer coefficient was varied for different absorptivity coefficients whilekeeping the distance between the heater bank and the sheet at 15 cm. Fig-ure 6.31(a) and Figure 6.31(c) collectively indicate that despite significantvariations in parameters, the performance of the multi-agent DRL remainsrobust without additional training. It consistently maintains the averageerror of zones on the sheet below 1.5 °C after the control process. Further-more, as illustrated in Figure 6.31(b), an increase in the distance betweenthe heater bank and the sheet results in a larger overshoot. This outcomeis logical, given that the heater’s radiation pattern covers a broader area onthe sheet, potentially impacting other zones as a disturbance. However, acomparison of the maximum overshoot at h = 5 [W/(m2K)], d = 20 [cm] toh = 5 [W/(m2K)], d = 15 [cm] (where RL agent has been trained) shows a1416.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMslight increase from 4 °C to 7 °C.5 10 15h [ Wm2"K ]00.511.5Averageerror[o C]d=10 cmd=15 cmd=20 cmTraining point(a)5 10 15h [ Wm2"K ]0246810Overshoot[o C]d=10 cmd=15 cmd=20 cmTraining point(b)5 10 15h [ Wm2"K ]00.511.5Averageerror[o C],=0.5,=0.6,=0.7,=0.8,=0.9Training point(c)5 10 15h [ Wm2"K ]0246Overshoot[o C],=0.5,=0.6,=0.7,=0.8,=0.9Training point(d)Figure (6.31) Robustness evaluation of the designed controller on the sim-ulator (NARX model is trained with initial values of h = 5,W/m2K,d = 15, cm, and α = 0.8). (a) The average error of zones under varyingconvection heat transfer coefficients and gaps between the heater bank andthe sheet. (b) Maximum overshoot under varying convection heat transfercoefficients and gaps between the heater bank and the sheet. (c) The averageerror of zones under varying convection heat transfer coefficients and sheetabsorptivity coefficients. (d) Maximum overshoot under varying convectionheat transfer coefficients and sheet absorptivity coefficients.Additionally, Figure 6.31 demonstrates that an elevation in the convec-tion heat transfer coefficient (due to a drop in room temperature or exposureof the sheet to wind) leads to a decrease in both maximum overshoot and the1426.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMfinal accuracy of the agent. The outcomes depicted in Figure 6.31 indicatethe trained multi-agent RL system is very robust in handling disturbancesintroduced by changes in system parameters. These results instill confidencein the potential effectiveness of the trained multi-agent DRL when deployedin a real-world lab-scale setup.Table 6.11 compares the robustness of the trained Multi-agent DRLmethod with the NARX-MPC method from the previous chapter. The tableindicates that both methods effectively control overshoot and zone errors inthe presence of disturbances. Both methods perform similarly regarding ac-curacy, with NARX-MPC and Multi-agent DRL maintaining average zoneerrors below 2°C and 1.5°C, respectively. Regarding overshoot, NARX-MPClimits it to under 7°C when system parameters vary, while Multi-agent DRLmaintains it under 9°C.Table (6.11) Robustness comparison between NARX-MPC and Multi-agent DRL for three system parameter variation.Max. ave. error [°C] Max. overshoot [°C]Range NARX-MPCMulti-agentDRLNARX-MPCMulti-agentDRLh[ Wm2K ] 3−15 2 1.5 7 9d [cm] 10−20 1.6 1.5 7 9α 0.5−0.9 2 1.4 7 6.56.5.7 Experimental ResultsThe trained multi-agent DRL agent was deployed on the laboratory-scalethermoforming setup as introduced in Section 3.3. Each heating elementhad its designated agent, with observation parameters defined in Section6.5.1. The identical reference temperature distributions in NARX-MPCexperimental evaluation were assumed to evaluate the agents’ performancewhen working collectively in a real-world setup, as illustrated in Figure5.17. The first desired temperature profile tasks the agents with achievingand maintaining a uniform temperature of 70 °C across all zones. Thesecond desired temperature profile requires the agents to reach and maintaindistinct temperatures in each zone, as indicated in Figure 5.17(b).Figure 6.33 visually represents multi-agent DRL’s performance during750 seconds in achieving the uniform reference temperature distribution of70 °C. Figure 6.32 shows the measured system’s output (temperature profileon the thermoplastic sheet) by IR camera in multiple time stamps duringthe control process through Multi-agent DRL for the uniform reference tem-1436.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEMperature distribution. To assess the performance, Figure 6.33 displays thevariation in error for the 15 zones of the thermoplastic material. Further-more, Figure 6.34 shows the control input commands to the heater bankin the same timeframe. Settling time, overshoots, and final temperatureerrors observed in the experimental phase correspond well with those in thesimulations.0 0.2 0.4 0.6x [m]00.10.20.30.40.5y[m]5060708090100(a)0 0.2 0.4 0.6x [m]00.10.20.30.40.5y[m]5060708090100(b)0 0.2 0.4 0.6x [m]00.10.20.30.40.5y[m]5060708090100(c)0 0.2 0.4 0.6x [m]00.10.20.30.40.5y[m]5060708090100(d)Figure (6.32) Experimental temperature profile measured by the IR cam-era: (a) t = 200 [sec], (b) t = 250 [sec], (c) t = 300 [sec], (d) t = 750 [sec]Figure 6.36 shows multi-agent DRL’s performance during 845 seconds inreaching and maintaining distinct temperatures in each zone, as indicated inFigure 5.17(b). Figure 6.35 illustrates the measured system’s output (tem-perature profile on the thermoplastic sheet) by IR camera in multiple timestamps during the control process through Multi-agent DRL for the nonuni-form reference temperature distribution. To assess the performance, Figure1446.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM6.36 displays the variation in error for the 15 zones of the thermoplasticmaterial. Furthermore, Figure 6.37 shows the control input commands tothe heater bank in the same timeframe.0 100 200 300 400 500 600 700-40-200Error[/ C]Zone: 1Zone: 2Zone: 3Zone: 4Zone: 50 100 200 300 400 500 600 700-40-200Error[/ C]Zone: 6Zone: 7Zone: 8Zone: 9Zone: 100 100 200 300 400 500 600 700Time [sec]-40-200Error[/ C]Zone: 11Zone: 12Zone: 13Zone: 14Zone: 15Figure (6.33) Experimental zone errors on the thermoplastic sheet for auniform reference temperature distribution.0 100 200 300 400 500 600 7000500Power[W] Heater: 1Heater: 2Heater: 3Heater: 4Heater: 50 100 200 300 400 500 600 7000500Power[W] Heater: 6Heater: 7Heater: 8Heater: 9Heater: 100 100 200 300 400 500 600 700Time [sec]0500Power[W] Heater: 11Heater: 12Heater: 13Heater: 14Heater: 15Figure (6.34) Control input to 15 heating elements during the control pro-cess.1456.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM0 0.2 0.4 0.6x [m]00.10.20.30.40.5y[m]20406080100(a)0 0.2 0.4 0.6x [m]00.10.20.30.40.5y[m]20406080100(b)0 0.2 0.4 0.6x [m]00.10.20.30.40.5y[m]20406080100(c)0 0.2 0.4 0.6x [m]00.10.20.30.40.5y[m]20406080100(d)Figure (6.35) Experimental temperature profile measured by the IR cam-era: (a) t = 150 [sec], (b) t = 200 [sec], (c) t = 250 [sec], (d) t = 845 [sec]Table 6.12 compares the multi-agent DRL’s performance between a physics-based simulator and a real-world lab-scale thermoforming setup. The resultsdemonstrate a close alignment in performance across both platforms whensubjected to the identical reference temperature distribution. The multi-agent DRL exhibits effective temperature control for the sheet in scenariosinvolving both uniform and nonuniform reference temperature profiles, witha final maximum error band of approximately ±4 °C. This is well within therequirements of typical industrial thermoforming machines, and it also in-dicates that the multiagent DRL controller trained on a validated softwaresimulator can be deployed without any retraining onto the correspondingthermoforming machine. The real-world implementation of the control sys-tem exhibits an overshoot of approximately 11.5 °C for nonuniform referencetemperature distribution and 6.8 °C for uniform distribution. These values1466.5. MULTI-AGENT DRL APPROACH IN THERMOFORMING SYSTEM0 100 200 300 400 500 600 700 800-60-40-200Error[/ C]Zone: 1Zone: 2Zone: 3Zone: 4Zone: 50 100 200 300 400 500 600 700 800-60-40-200Error[/ C]Zone: 6Zone: 7Zone: 8Zone: 9Zone: 100 100 200 300 400 500 600 700 800Time [sec]-60-40-200Error[/ C]Zone: 11Zone: 12Zone: 13Zone: 14Zone: 15Figure (6.36) Experimental zone errors on the thermoplastic sheet for anonuniform reference temperature distribution.0 100 200 300 400 500 600 700 8000500Power[W] Heater: 1Heater: 2Heater: 3Heater: 4Heater: 50 100 200 300 400 500 600 700 8000500Power[W] Heater: 6Heater: 7Heater: 8Heater: 9Heater: 100 100 200 300 400 500 600 700 800Time [sec]0500Power[W] Heater: 11Heater: 12Heater: 13Heater: 14Heater: 15Figure (6.37) Control input to 15 heating elements during the control pro-cess.1476.6. CHAPTER SUMMARYTable (6.12) Multi-agent DRL validation.nonuniform uniformMetric Sim. Exp. Sim. Exp.Zones’ average error [°C] 0.5 1.6 0.3 1.3Zones’ maximum error [°C] ±1.5 ±4 ±0.5 ±3.8Overshoot [°C] 4.7 11.5 3.5 6.8Settling time within an errorband of ±10 °C [sec]200 236 180 204exceed the corresponding overshoots observed in the simulator by approx-imately 7 and 3 °C, respectively. Notably, the overshoot for the uniformtemperature distribution closely matches that achieved by the NARX-MPCcontroller, which was 4 °C. In nonuniform temperature distributions, themaximum overshoot occurs in zones with the lowest reference temperatures.However, it’s worth noting that even though the overshoot is relatively high,it remains within the forming window of the material. Additionally, as illus-trated in Figure 6.29, there is a correlation between the reference tempera-ture and the magnitude of the overshoot in the control system, with higherreference temperatures resulting in lower overshoot values.6.6 Chapter SummaryThis chapter explores the application of reinforcement learning tech-niques to optimize the heating phase of the thermoforming process. Theintegration of RL into this industrial application demonstrates the poten-tial of machine learning algorithms to enhance efficiency and precision inmanufacturing settings. The chapter begins by introducing the fundamen-tals of RL, emphasizing its capability to optimize decision-making processesthrough interaction with the environment and learning from the outcomes.The chapter underscores the time-consuming and expensive nature of theoffline training process for DRL when dealing with systems with large actionand observation spaces, prompting the need for innovative solutions.Two distinct approaches are proposed to tackle this challenge. Firstly, ahybrid method is introduced, seamlessly integrating MPC with DRL. Thisnovel approach harnesses the predictive capabilities of MPC to accelerate thelearning performance of RL agents, resulting in a remarkable 100-fold im-provement in average episode rewards over traditional DRL methods. Thisapproach is named MPC-guided DRL in this study. Secondly, a multi-agentDRL method is presented, wherein each agent controls a specific heating el-1486.6. CHAPTER SUMMARYement within the thermoforming setup. This approach showcases significantadvantages over conventional single-agent RL, particularly regarding train-ing efficiency and control precision. The offline training process is expedited110 times compared to traditional DRL, leading to an 8 times reductionin the final error margin. Experimental trials on a laboratory-scale setupvalidate the applicability and effectiveness of the multi-agent RL approach,demonstrating robustness in achieving and maintaining desired temperatureprofiles across the thermoplastic sheet.149Chapter 7Conclusions and FutureWorksThis chapter presents the conclusions drawn from the research conductedin this thesis. Section 7.1 discusses the main contributions, Section 7.2summarizes the content in this thesis, and Section 7.3 suggests directionsfor future research.7.1 Research Contributions and ConclusionsThis thesis offers a comprehensive study of the heating phase in thethermoforming process, making the following contributions:1. Development of a Novel 2D Physics-Based Heating Model: Iintroduced a 2D physics-based heating model for thermoforming by in-corporating experimental view factors. This model enhances accuracycompared to traditional 2D models while maintaining their computa-tional efficiency (see Chapter 3).2. Creation of a Data-Driven NARX Heating Model: I devel-oped a data-driven NARX model for the thermoforming heating phase.Based on input-output data, this model simplifies the state space whenlinearizing the system by disregarding the internal states, thereby re-ducing complexity compared to models solely relying on physical prin-ciples (see Section 5.3.1).3. A Novel Evaluation Method for Assessing Feasibility of Ref-erence Temperature Profiles on Thermoplastic Sheets: I intro-duced a method for the first time to assess the feasibility of achieving agiven desired temperature distribution, considering experimental setupconfigurations (see Section 4.3).4. Robust NARX-MPC Controller Implementation: For the firsttime, I designed and successfully implemented a robust NARX-MPC1507.1. RESEARCH CONTRIBUTIONS AND CONCLUSIONScontroller for thermoforming heating phases in an experimental setup.The NARX-MPC can control the sheet temperature with shorter risetimes and lower overshoots than traditional PI controllers (see Table5.3). This can double the production rate (see Section 5.3).5. MPC-guided Deep Reinforcement Learning: Considering thetime-consuming and expensive nature of the offline training processfor DRL when dealing with systems with large action and observationspaces, a hybrid method is introduced, seamlessly integrating MPCwith DRL. This novel approach harnesses the predictive capabilitiesof MPC to accelerate the learning performance of RL agents, resultingin a remarkable 100-fold improvement in average episode rewards overtraditional DRL methods (see Section 6.4).6. Scalable AI-Based Controller Implementation: For the firsttime, I designed and successfully implemented a scalable AI-basedMulti-agent DRL controller for thermoforming heating phases in an ex-perimental setup. I introduced a set of boundary conditions for multi-agent DRL training in thermoforming, resulting in reduced trainingtime, marking a significant advancement with shorter settling timesand lower computational load compared to NARX-MPC (see Tables6.9 and 7.1). This can increase the production rate from 3 cycles (whenusing a traditional PI controller) to 11 cycles (see Table 6.5.5). Fur-thermore, in the centralized DRL approach where traditional single-agent DRL could not converge the training process, I proposed simpli-fying DRL training in thermoforming systems with large action spacesby integrating the MPC into DRL (see Section 6.4) achieving a 100-fold increase in average episode rewards during training compared totraditional single-agent DRL techniques.Two indirect data-driven controllers (NARX-MPC and Multi-agent DRL)were designed and implemented on a real-world lab-scale thermoformingsetup (UBC thermoforming system 3.3). I evaluated their performance us-ing two types of reference temperature profiles (uniform and nonuniform in5.17). Four criteria were considered: average and maximum error of zones atthe end of the control horizon, maximum overshoot, and settling time withinan error band of ±10 °C (see Table 7.1). Furthermore, their performance ina cycle-to-cycle control approach was studied (see table 6.10).The NARX-MPC demonstrated slightly better accuracy on the exper-imental setup than the Multi-agent DRL. However, the Multi-agent DRL1517.2. CONTENT OF THIS THESISeffectively reduced settling time and achieved the desired temperature pro-files more quickly. In terms of scalability and ease of implementation in in-dustrial settings, the Multi-agent DRL offers a more practical solution dueto its adaptability to new configurations without requiring the developmentof a new model, unlike the NARX-MPC.Regarding online computational cost, as the number of heaters in theheater bank increases, the number of samples in the dataset needed to cre-ate the NARXmodel also increases. This leads to a more complex model andhigher computational costs for NARX-MPC. Conversely, the online compu-tational cost for Multi-agent DRL remains consistent regardless of the heaterbank size.Table (7.1) A comparison between NARX-MPC and Multi-agent DRL ex-perimental results for uniform and non-uniform reference temperature dis-tributions.Metric NARX MPC Multi-agent DRLuni. nonuni. uni. nonuni.Average error [°C] 0.8 1 1.3 1.6Maximum error [°C] 1.6 3.8 3.8 4Overshoot [°C] 4 5.3 6.8 11.5Settling time within anerror band of ±10 °C [sec]297 308 204 2367.2 Content of this thesisChapter 1 provides an introduction to the presented thesis. Chapter2 delves into the background of the research presented. It outlines thevarious types of modeling employed in thermoforming, examines existingtemperature-sensing methodologies utilized in thermoforming manufactur-ing, and explores the diverse control strategies previously investigated. Ad-ditionally, it introduces two potential approaches that will be explored inthis thesis: Model Predictive Control and Reinforcement Learning.In Chapter 3, the thermal dynamics of the thermoforming process’s heat-ing phase are modeled in detail. This chapter develops a novel 2D physics-based model to investigate the temperature distribution on the thermoplas-tic sheet during the heating phase of the thermoforming process. Importingthe experimental view factor matrix into the 2D model was suggested and1527.2. CONTENT OF THIS THESISdone for the first time. With this method, Icould create a 2D model withthe low computing effort of a 2D model and the accuracy of a 3D model.The team designed and built a lab-scale thermoforming setup to study pro-cess optimization, and the details of the developed hardware are explainedin this chapter. Integrating an IR camera in the control loop to get andevaluate the thermoplastic sheet temperature feedback in the thermoform-ing systems is introduced and implemented for the first time. This lab-scalesetup was used to design the experiment to produce the experimental viewfactor matrix between the heater and thermoplastic sheets. Also, this lab-scale setup was utilized to verify the performance of the created heatingphase simulator.Chapter 4 of this thesis assesses the system’s stability and explores its dy-namic behavior. The control objectives for thermoforming systems are thendiscussed. A novel method is introduced to assess the feasibility and reach-ability of a reference temperature distribution on the thermoplastic sheet,considering the configuration of the setup. The open-loop heating time ofthe sheet for a particular reference temperature distribution is studied on thesimulator in the cycle-to-cycle. The results of the open-loop approach, whichis used in the industry, are compared to the results of the closed-loop controlin the next chapters to show the importance of automation in this phase.The effectiveness of a proportional-integral (PI) controller is examined asa fundamental control strategy. This chapter explains the shortcomings ofusing only basic control techniques through careful analysis. It emphasizesthe importance of having more complex control strategies adapted to thespecifics of the system in question.Chapter 5 delves into applying Model Predictive Control (MPC) as arobust model-based control technique in the heating phase of the thermo-forming process. The chapter initially focuses on designing an adaptive MPCsystem based on the linearized model derived from the physics-based simu-lator. Even though the results are better than PI on the physics-based sim-ulator, due to the linearized model’s requirement for internal system states(specifically, heater bank temperature feedback) at each operational point,coupled with hardware limitations in implementation and the high compu-tational cost associated with online execution, an alternative approach isexplored. In the second step, a data-driven model of the thermoformingsystem is developed using Nonlinear AutoRegressive with eXogenous inputs(NARX) techniques. This model is then integrated into the MPC frameworkand deployed on the lab-scale setup. The results showcase a substantial en-hancement in MPC performance facilitated by integrating the NARX model.NARX MPC performs better than AMPC in simulations. It is also easier1537.3. FUTURE WORKSto deploy in real-world setups than AMPC. However, every new deploymentto different setups requires the development of a new NARX model.Chapter 6 delves into the realm of Deep Reinforcement Learning (DRL)as a transformative approach for industrial control systems, particularly fo-cusing on its application in the heating phase of the thermoforming process.The chapter underscores the time-consuming and expensive nature of theoffline training process for DRL when dealing with systems with large ac-tion and observation spaces, prompting the need for innovative solutions.Two distinct approaches are proposed to tackle this challenge. Firstly, a hy-brid method is introduced, seamlessly integrating MPC with DRL. Thisapproach harnesses the predictive capabilities of MPC to accelerate thelearning performance of RL agents, resulting in a remarkable 100-fold im-provement in average episode rewards over traditional DRL methods. Thisapproach is named MPC-guided DRL in this study. Secondly, a multi-agentDRL method is presented, wherein each agent controls a specific heatingelement within the thermoforming setup. This approach showcases signifi-cant advantages over conventional single-agent DRL, particularly regardingtraining efficiency and control precision. The offline training process is ex-pedited 110 times compared to MPC-guided DRL, leading to an 8 timesreduction in the final error margin. Experimental trials on a laboratory-scale setup validate the applicability and effectiveness of the multi-agent RLapproach, demonstrating robustness in achieving and maintaining desiredtemperature profiles across the thermoplastic sheet. It is also important tonote that agents are trained for a particular heating element, and differentsets with the same elements can be deployed without additional training.7.3 Future WorksOne avenue for further exploration involves determining the optimal ini-tial surface temperature range for heaters in the thermal control processbased on the desired reference temperature distribution on the thermoplas-tic sheet. Following the completion of the thermoforming process for athermoplastic sheet, and as a new sheet takes its place, the heaters in theheater bank have already attained a specific surface temperature. A crit-ical future study lies in identifying the optimal initial temperature rangefor these heaters, ensuring that the thermal control process time for a cy-cle can be minimized without triggering a substantial overshoot that couldlead to sheet melting. This optimal initial surface temperature range shouldbe reached when the previous sheet forms and cools. This optimization ef-1547.3. FUTURE WORKSfort could significantly enhance the efficiency of the thermoforming industryby streamlining the production process and reducing cycle times for eachmanufactured part.Another critical future study lies in expanding the modeling chapter forthick thermoplastic sheets and considering this parameter in studying therobustness of the designed controllers.Moreover, the research will focus on direct data-driven control techniquessuch as data-enabled predictive control (DeePC) [CLD19] that eliminate thesystem identification step and design the controller directly from the dataset.This may enable us to achieve better performance while simplifying thedesign process.155Bibliography[AAM21] Rym Azzouz, Samir Allaoui, and Raphael Moulart. Com-posite preforming defects: a review and a classification.International Journal of Material Forming, 14(6):1259–1278, 2021. → pages 1[ABD+06] Zahir Albadawi, Benoit Boulet, Robert DiRaddo, PatrickGirard, and Vincent Thomson. Agent-based control forthermoforming processes. IFAC Proceedings Volumes,39:465–470, 2006. → pages 5, 17[Aje04] M. Ajersch. Modeling and Real-Time Control of SheetReheat Phase in Thermoforming. MSc thesis, 2004. →pages 5, 13, 15[AM99] Derek P Atherton and Somanath Majhi. Limitations ofpid controllers. In Proceedings of the 1999 American Con-trol Conference (Cat. No. 99CH36251), volume 6, pages3843–3847. IEEE, 1999. → pages 18[AMM21] KEVIN AMMOURI. Deep Reinforcement Learning forTemperature Control in Buildings and Adversarial At-tacks. MSc thesis, KTH ROYAL INSTITUTE OFTECHNOLOGY, 2021. → pages 20[AORADdJRT19] Roque Alfredo Osornio-Rios, Jose Alfonso Antonino-Daviu, and Rene de Jesus Romero-Troncoso. Recent in-dustrial applications of infrared thermography: A review.IEEE Transactions on Industrial Informatics, 15(2):615–625, 2019. → pages 15[Arn13] DL Arnold. Visor and method of manufacture. 2013,2013. → pages 2156Bibliography[BBM17] Francesco Borrelli, Alberto Bemporad, and ManfredMorari. Predictive Control for Linear and Hybrid Sys-tems. Cambridge University Press, 2017. → pages 17[BBX+21] Nu´ria Banu´s, Imma Boada, Pau Xiberta, Pol Toldra`,and Narc´ıs Bustins. Deep learning for the quality con-trol of thermoforming food packages. Scientific Reports,11:21887, 11 2021. → pages 12[BDGS18] PD Bairagi, AB Darekar, SB Gondker, and RB Sauda-gar. Pharmaceutical packaging materials: a brief review,2018. → pages 2[BDMT98] Bruce Bugbee, Matt Droter, Oscar Monje, and BertrandTanner. Evaluation and modification of commercial infra-red transducers for leaf temperature measurement. Ad-vances in Space Research, 22:1425–1434, 1 1998. → pages15[BHM16] P. Boisse, N. Hamila, and A. Madeo. Modelling the devel-opment of defects during composite reinforcements andprepreg forming. Philosophical Transactions of the RoyalSociety A: Mathematical, Physical and Engineering Sci-ences, 374:20150269, 7 2016. → pages 12[BLAV22] Philip Bean, Roberto A. Lopez-Anido, and Senthil Vel.Integration of material characterization, thermoformingsimulation, and as-formed structural analysis for thermo-plastic composites. Polymers, 14:1877, 5 2022. → pages12[BMBD17] B. Buffel, B. Van Mieghem, A. Van Bael, and F. De-splentere. Optimization of the ir-heating phase in ther-moforming of thermoplastic sheets: Characterization andmodelling. page 020001, 2017. → pages 13, 15[BMDP02] Alberto Bemporad, Manfred Morari, Vivek Dua, and Ef-stratios N. Pistikopoulos. The explicit linear quadraticregulator for constrained systems. Automatica, 38:3–20,1 2002. → pages 5157Bibliography[BNK01] John M. Baker, John M. Noman, and Atsushi Kano. Anew approach to infrared thermometry. Agricultural andForest Meteorology, 108:281–292, 7 2001. → pages 15[BPM+17] Philip Blaser, Florentina Pavlicˇek, Kotaro Mori, JosefMayr, Sascha Weikert, and Konrad Wegener. Adap-tive learning control for thermal error compensation of5-axis machine tools. Journal of Manufacturing Systems,44:302–309, 7 2017. → pages 16[Buf16] Bart Buffel. A numerical study on different heat-ing strategies in thermoforming of thermoplastic sheets.pages 221–226. Belgie¨, 2016. → pages 16[BvMvBD17] B. Buffel, B. van Mieghem, A. van Bael, and F. Desplen-tere. A combined experimental and modelling approachtowards an optimized heating strategy in thermoformingof thermoplastics sheets. International Polymer Process-ing, 32:378–386, 7 2017. → pages 13[BXFS22] Fabio Bonassi, Jing Xie, Marcello Farina, and RiccardoScattolini. An offset-free nonlinear mpc scheme for sys-tems learned by neural narx models. In 2022 IEEE 61stConference on Decision and Control (CDC), pages 2123–2128. IEEE, 2022. → pages 19[CB10] Md. Muminul Islam Chy and Benoit Boulet. A newmethod for estimation and control of temperature pro-file over a sheet in thermoforming process. pages 1–8.IEEE, 10 2010. → pages 17[CB11] Muminul Islam Chy and Benoit Boulet. Development ofan improved mathematical model of the heating phase ofthermoforming process. pages 1–8. IEEE, 10 2011. →pages 14, 15[CB12] M. M. I. Chy and B. Boulet. Iterative learning modelpredictive controller of plastic sheet temperature for athermoforming process. pages 627–633. IEEE, 6 2012.→ pages 19[CBH11] Md Muminul Islam Chy, Benoit Boulet, and AhmadHaidar. A model predictive controller of plastic sheet158Bibliographytemperature for a thermoforming process. pages 4410–4415. IEEE, 6 2011. → pages 5, 16, 17, 19[Cer] Ceramicx company, https://www.ceramicx.com/. →pages xiv, xv, 23, 35, 36[Cha94] BPS Chawla. Pharmaceutical packaging, 1994. → pages2[CHC+21] J Choi, C Han, S Cho, K Kim, J Ahn, D Del Orbe, et al.Customizable, conformal, and stretchable 3d electronicsvia predistorted pattern generation and thermoforming.Sci Adv, 7, 2021. → pages 2[CLD19] Jeremy Coulson, John Lygeros, and Florian Do¨rfler.Data-enabled predictive control: In the shallows ofthe deepc. In 2019 18th European Control Conference(ECC), pages 307–312. IEEE, 2019. → pages 155[CMS] Cms plastic technology company,https://www.scmgroup.com/en/cmsplastic. → pagesxv, 46[Coma] Compuplast na company, https://www.compuplast.ca/.→ pages xiv, 6[COMb] Comsol multiphysics® v. 5.6. www.comsol.com. comsolab, stockholm, sweden. → pages 42[Das22] Aniket Ghosh Dastidar. Modeling and simulation ofthermoforming. PhD thesis, Ecole nationale supe´rieureMines-Te´le´com Lille Douai, 2022. → pages xii, 2[DC02] F. M. Duarte and J. A. Covas. Ir sheet heating in roll fedthermoforming: Part 1 - solving direct and inverse heat-ing problems. Plastics, Rubber and Composites, 31:307–317, 9 2002. → pages 13[DC03] F. M. Duarte and J. A. Covas. Infrared sheet heating inroll fed thermoforming: Part 2 - factors influencing in-verse heating solution. Plastics, Rubber and Composites,32:32–39, 1 2003. → pages 13159Bibliography[DDM+20] Hamid Daghigh, Vahid Daghigh, Abbas Milani, DwayneTannant, Thomas E. Lacy, and J.N. Reddy. Nonlo-cal bending and buckling of agglomerated cnt-reinforcedcomposite nanoplates. Composites Part B: Engineering,183:107716, 2 2020. → pages 43[DGN+22] Giulio D’Emilia, Antonella Gaspari, Emanuela Natale,Antonios G. Stamopoulos, and Antoniomaria Di Ilio.Experimental and numerical analysis of the defects in-duced by the thermoforming process on woven textilethermoplastic composites. Engineering Failure Analysis,135:106093, 5 2022. → pages 12[DHMDK14] Abhishek Dutta, Edward Hartley, Jan Maciejowski, andRobin De Keyser. Certification of a class of industrial pre-dictive controllers without terminal conditions. In 53rdIEEE Conference on Decision and Control, pages 6695–6700. IEEE, 2014. → pages 18[DNP+13] Marc Peter Deisenroth, Gerhard Neumann, Jan Peters,et al. A survey on policy search for robotics. Foundationsand Trends® in Robotics, 2(1–2):1–142, 2013. → pages104[dRTAOR+22] Alejandro del Real Torres, Doru Stefan Andreiana,A´lvaro Ojeda Rolda´n, Alfonso Herna´ndez Bustos, andLuis Enrique Acevedo Galicia. A review of deep reinforce-ment learning approaches for smart manufacturing in in-dustry 4.0 and 5.0 framework. Applied Sciences, 12(23),2022. → pages 21[EA19] Roja Eini and Sherif Abdelwahed. Learning-based modelpredictive control for smart building thermal manage-ment. In 2019 IEEE 16th International Conference onSmart Cities: Improving Quality of Life Using ICT &IoT and AI (HONET-ICT), pages 038–042, 2019. →pages 19[EE14] Ertug˘rul Selc¸uk Erdog˘an and Olcay Eks¸i. Effects ofmanufacturing defects on thermoformed product quality.Usak University Journal of Material Sciences, 3:61–61, 62014. → pages 12160Bibliography[EN07] Fouad Erchiqui and G. Dituba Ngoma. Analyse compar-ative des me´thodes de calcul des facteurs de formes pourdes surfaces a` contours rectilignes. International Journalof Thermal Sciences, 46:284–293, 3 2007. → pages 14, 30[Erc18] F. Erchiqui. Application of genetic and simulated anneal-ing algorithms for optimization of infrared heating stagein thermoforming process. Applied Thermal Engineering,128:1263–1272, 1 2018. → pages 14, 15[Ern05] Aliza Erner. Etude experimentale du thermoformage as-siste´ par poinc¸on d’un me´lange de polystyre`nes. PhDthesis, E´cole Nationale Supe´rieure des Mines de Paris,2005. → pages 3[ES93] J. R. Ehlert and T. F. Smith. View factors for perpendic-ular and parallel rectangular plates. Journal of Thermo-physics and Heat Transfer, 7:173–175, 1 1993. → pagesxiv, 30, 40[ESY07] F. Erchiqui, M. Souli, and R. Ben Yedder. Nonisother-mal finite-element analysis of thermoforming of polyethy-lene terephthalate sheet: Incomplete effect of the formingstage. Polymer Engineering & Science, 47:2129–2144, 122007. → pages 15[FHM18] Scott Fujimoto, Herke Hoof, and David Meger. Address-ing function approximation error in actor-critic methods.In International conference on machine learning, pages1587–1596. PMLR, 2018. → pages xviii, 107[FPHG15a] Michael G. Forbes, Rohit S. Patwardhan, HamzaHamadah, and R. Bhushan Gopaluni. Model predictivecontrol in industry: Challenges and opportunities. IFAC-PapersOnLine, 48:531–538, 2015. → pages 11[FPHG15b] Michael G Forbes, Rohit S Patwardhan, HamzaHamadah, and R Bhushan Gopaluni. Model predictivecontrol in industry: Challenges and opportunities. IFAC-PapersOnLine, 48(8):531–538, 2015. → pages 69[FvHK14] Simone Formentin, Klaske van Heusden, and AlirezaKarimi. A comparison of model-based and data-driven161Bibliographycontroller tuning. International Journal of Adaptive Con-trol and Signal Processing, 28:882–897, 10 2014. → pages16[GAB+05] Guy Gauthier, Mark Ajersch, Benoˆıt Boulet, AmmarHaurani, Patrick Girard, and Robert Diraddo. A newabsorption based model for sheet reheat in thermoform-ing. 5 2005. → pages 15[Gau08] G. Gauthier. Terminal iterative learning for cycle-to-cycle control of industrial processes. PhD thesis, 2008.→ pages 13, 15[GCCW21] Yang Gu, Yuhu Cheng, CL Philip Chen, and XuesongWang. Proximal policy optimization with policy feed-back. IEEE Transactions on Systems, Man, and Cyber-netics: Systems, 52(7):4600–4610, 2021. → pages 110[GLM+04] E Gime´nez, JM Lagaro´n, ML Maspoch, L Cabedo, andJJ Saura. Uniaxial tensile behavior and thermoformingcharacteristics of high barrier evoh-based blends of inter-est in food packaging. Polym Eng Sci, 44:598–608, 2004.→ pages 2[GP21] Swetha Govindaiah and Mikel D Petty. Applying rein-forcement learning to plan manufacturing material han-dling. Discover Artificial Intelligence, 1(1):8, 2021. →pages 20[HARFCLCA18] Jorge Hernandez-Ambato, Jesus Rodriguez-Flores,J Cortes-Llanganate, and F Cabrera-Aguayo. Classiccontrollers design applied to temperature control fora plastic thermoforming machine. In 2018 IEEE 23rdInternational Conference on Emerging Technologies andFactory Automation (ETFA), volume 1, pages 551–558.IEEE, 2018. → pages 17[HBMZ06] Hossein Hosseini, Boris Vasilivich Berdyshev, and Ar-jomand Mehrabani-Zeinabad. Modeling of deformationprocesses in vacuum thermoforming of a pre-stretchedsheet. Polymer-Plastics Technology and Engineering,45:1357–1362, 12 2006. → pages 12162Bibliography[HCLB18] Lianghui Huang, Jian Chen, Zhiyang Liu, and MohamedBecherif. Adaptive thermal control for pemfc systemswith guaranteed performance. International Journal ofHydrogen Energy, 43:11550–11558, 6 2018. → pages 16[HFZ+24] Zongyang Hu, Jiuwen Fang, Ruixiang Zheng, Mian Li,Baosheng Gao, and Lingcan Zhang. Efficient model pre-dictive control of boiler coal combustion based on narxneutral network. Journal of Process Control, 134:103158,2024. → pages 19[HMDS20] John R. Howell, M. Pinar Mengu¨c¸, Kyle Daun, andRobert Siegel. Thermal Radiation Heat Transfer. CRCPress, 12 2020. → pages 29[Hol09] J.P. Holman. Heat transfer, New York. McGraw-Hill,10th edition edition, 2009. → pages 27, 29[HRSM23] Hadi Hosseinionari, Milad Ramezankhani, RudolfSeethaler, and Abbas S Milani. Development of a com-putationally efficient model of the heating phase in ther-moforming process based on the experimental radiationpattern of heaters. Journal of Manufacturing and Mate-rials Processing, 7(1):48, 2023. → pages xiv, xv, 8, 15,26, 29, 30, 32, 35, 36, 37, 38, 39, 40, 41, 42[HS94a] N. Hajji and J. E. Spruiell. Radiation pyrometry on semi-transparent sheets. i: Gray media. Polymer Engineeringand Science, 34:116–121, 1 1994. → pages 14[HS94b] N. Hajji and J. E. Spruiell. Radiation pyrometry on semi-transparent sheets. ii: Media with wavelength dependentabsorption coefficient. Polymer Engineering and Science,34:122–127, 1 1994. → pages 14[HS24] Hadi Hosseinionari and Rudolf Seethaler. The integra-tion of model predictive control and deep reinforcementlearning for efficient thermal control in thermoformingprocesses. Journal of Manufacturing Processes, 115:82–93, 2024. → pages xiv, xvi, xviii, xix, 9, 24, 69, 109, 114,115, 116, 118, 121, 122, 123163Bibliography[HW13] Zhong-Sheng Hou and Zhuo Wang. From model-basedcontrol to data-driven control: Survey, classification andperspective. Information Sciences, 235:3–35, 6 2013. →pages 16[HZAL18] Tuomas Haarnoja, Aurick Zhou, Pieter Abbeel, andSergey Levine. Soft actor-critic: Off-policy maximum en-tropy deep reinforcement learning with a stochastic actor.In International conference on machine learning, pages1861–1870. PMLR, 2018. → pages xviii, 107[IM16] S. M. Ivanova and T. Muneer. Finite-element heat-transfer computations for parallel surfaces with uniformor non-uniform emitting. Journal of Renewable and Sus-tainable Energy, 8:015102, 1 2016. → pages 30[JL23] Howoun Jung and Jay H. Lee. Flexible operation of post-combustion co2 capture process enabled by narx-mpc us-ing neural network. Computers & Chemical Engineering,179:108447, 2023. → pages 19[Kha02] Hassan K Khalil. Control of nonlinear systems. PrenticeHall, New York, NY, 2002. → pages 56[KKPP23] Lotta Kannari, Julia Kantorovitch, Kalevi Piira, andJouko Piippo. Energy cost driven heating control withreinforcement learning. Buildings, 13(2), 2023. → pages20[KMK13] Teck Y Kim, B Min, and Won K Kim. General charac-teristics of packaging materials for food system. ElsevierLtd, 2013. → pages 2[KST13] D Knoeppel, L Sun, and J Tippet. High impactpolystyrene with high gloss and high impact strength.WO 2013/154592 Al, 2013. → pages 2[Kum05] Vijay Kumar. Estimation of absorptivity and heat fluxat the reheat phase of thermoforming process. 2005. →pages 13[LCRMC+18] Wanderson de Oliveira Leite, Juan Carlos Campos Ru-bio, Francisco Mata Cabrera, Angeles Carrasco, and Is-164Bibliographysam Hanafi. Vacuum thermoforming process: an ap-proach to modeling and optimization using artificial neu-ral networks. Polymers, 10(2):143, 2018. → pages 1[LCYW19] Boyi Liu, Qi Cai, Zhuoran Yang, and Zhaoran Wang.Neural proximal/trust region policy optimization attainsglobally optimal policy. arXiv preprint arXiv:1906.10306,2019. → pages 110[LDJ23] Yanrong Li, Juan Du, and Wei Jiang. Reinforcementlearning for process control with application in semicon-ductor manufacturing. IISE Transactions, pages 1–15,2023. → pages 21[LFL+22] Nathan P. Lawrence, Michael G. Forbes, Philip D.Loewen, Daniel G. McClement, Johan U. Backstro¨m,and R. Bhushan Gopaluni. Deep reinforcement learningwith shallow controllers: An experimental application topid tuning. Control Engineering Practice, 121:105046, 42022. → pages 102[LHP+15] T. P. Lillicrap, J. J. Hunt, A. Pritzel, N. Heess, T. Erez,Y. Tassa, D. Silver, and D.Wierstra. Continuous controlwith deep reinforcement learning, 2015. → pages xviii,102, 107[LKHH20] Hyun-Kyo Lim, Ju-Bong Kim, Joo-Seong Heo, and Youn-Hee Han. Federated reinforcement learning for trainingcontrol policies on multiple iot devices. Sensors, 20:1359,3 2020. → pages 102[LM96] J. F. Lappin and P. J. Martin. Sheet temperature mea-surement in thermoforming. pages 878–882, 1996. →pages 15[LMX+10a] Zhen-Zhe Li, Guang Ma, Dong-Ji Xuan, Seoung-YunSeol, and Yun-De Shen. A study on control of heaterpower and heating time for thermoforming. Interna-tional Journal of Precision Engineering and Manufactur-ing, 11:873–878, 12 2010. → pages 13[LMX+10b] Zhen-Zhe Li, Guang Ma, Dong-Ji Xuan, Seoung-YunSeol, and Yun-De Shen. A study on control of heater165Bibliographypower and heating time for thermoforming. Interna-tional Journal of Precision Engineering and Manufactur-ing, 11:873–878, 2010. → pages 17[LNK+18] Won Je Lee, Jonggeol Na, Kyeongsu Kim, Chul-Jin Lee,Younggeun Lee, and Jong Min Lee. Narx modeling forreal-time optimization of air and gas compression systemsin chemical processes. Computers & Chemical Engineer-ing, 115:262–274, 2018. → pages 18, 76[LSR18] Jay H Lee, Joohyun Shin, and Matthew J Realff. Machinelearning: Overview of the recent progresses and implica-tions for the process systems engineering field. Computers& Chemical Engineering, 114:111–121, 2018. → pages 99[LVL19] M. Luzi, M. Vaccarini, and M. Lemma. A tuning method-ology of model predictive control design for energy effi-cient building thermal control. Journal of Building En-gineering, 21:28–36, 1 2019. → pages 16[MA] Peter Milbradt and Wassim Abu Abed. Generalized sta-bilization techniques in computational fluid dynamics.URL: https://www. researchgate. → pages xiv, 25[MA15] Noraini Mohd and Norashid Aziz. Control of bioethanolfermentation process: Narx-based mpc (narx-mpc) ver-sus linear-based mpc (lmpc). Chemical EngineeringTransactions, 45:1297–1302, 2015. → pages 19[MBM+16] Volodymyr Mnih, Adria Puigdomenech Badia, MehdiMirza, Alex Graves, Timothy Lillicrap, Tim Harley,David Silver, and Koray Kavukcuoglu. Asynchronousmethods for deep reinforcement learning. In Interna-tional conference on machine learning, pages 1928–1937.PMLR, 2016. → pages xviii, 107[Mek10] S Mekid. Spatial thermal mapping using thermal infraredcamera and wireless sensors for error compensation viaopen architecture controllers. Proceedings of the Institu-tion of Mechanical Engineers, Part I: Journal of Systemsand Control Engineering, 224:789–798, 11 2010. → pages15166Bibliography[MJ16] Jin Woo Moon and Sung Kwon Jung. Development ofa thermal control algorithm using artificial neural net-work models for improved thermal comfort and energyefficiency in accommodation buildings. Applied ThermalEngineering, 103:1135–1144, 6 2016. → pages 16[MKS+13] Volodymyr Mnih, Koray Kavukcuoglu, David Silver,Alex Graves, Ioannis Antonoglou, Daan Wierstra, andMartin Riedmiller. Playing atari with deep reinforcementlearning. arXiv preprint arXiv:1312.5602, 2013. → pagesxviii, 105, 107[ML99] Manfred Morari and Jay H. Lee. Model predictive con-trol: past, present and future. Computers & ChemicalEngineering, 23:667–682, 5 1999. → pages 5[MLB+] D. G. McClement, N. P. Lawrence, J. U. Backstrom, P. D.Loewen, M. G. Forbes, and R. B. Gopaluni. Meta rein-forcement learning for adaptive control: An offline ap-proach. submitted to Journal of Process Control, 2022.→ pages 11[MLL+21] Daniel G. McClement, Nathan P. Lawrence, Philip D.Loewen, Michael G. Forbes, Johan U. Backstro¨m, andR. Bhushan Gopaluni. A meta-reinforcement learning ap-proach to process control. IFAC-PapersOnLine, 54:685–692, 2021. → pages 11[MM10] R McCool and PJ Martin. The role of process parametersin determining wall thickness distribution in plug-assistedthermoforming. Polym Eng Sci, 50:1923–1934, 2010. →pages 2[Mod03] Michael F. Modest. Radiative Heat Transfer. Elsevier,second edition, 2003. → pages 29, 36[Moo02] B. Moore. In-Cycle Control of the Thermoforming ReheatProcess. MSc thesis, 2002. → pages 13, 15[MRRS00] D.Q. Mayne, J.B. Rawlings, C.V. Rao, and P.O.M.Scokaert. Constrained model predictive control: Stabil-ity and optimality. Automatica, 36:789–814, 6 2000. →pages 5167Bibliography[MSMD98] Serge Monteix, Fabrice Schmidt, Yannick Le Maoult, andR. Diraddo. Profiled infrared radiative heating in blowmoulding and thermoforming. North American Meeting,1998. → pages 13, 14[NB19] Benjamin Neubig and Christian Bonten. Thermoforming– determining temperature profiles across sheet thicknessby pyrometric measurements. page 070029, 2019. →pages 15, 17[NB23] Anastasia Nikolakopoulou and Richard D. Braatz. Poly-nomial narx-based nonlinear model predictive control ofmodular chemical systems. Computers & Chemical En-gineering, 177:108272, 2023. → pages 19[NLH20a] Rui Nian, Jinfeng Liu, and Biao Huang. A review onreinforcement learning: Introduction and applications inindustrial process control. Computers & Chemical Engi-neering, 139:106886, 8 2020. → pages 11[NLH20b] Rui Nian, Jinfeng Liu, and Biao Huang. A review onreinforcement learning: Introduction and applications inindustrial process control. Computers & Chemical Engi-neering, 139:106886, 2020. → pages 99[NN20] Oliver Nelles and Oliver Nelles. Nonlinear dynamic sys-tem identification. Springer, 2020. → pages 80[NPBB+22] Nuria Nievas, Adela Page`s-Bernaus, Francesc Bonada,Llu´ıs Echeverria, Albert Abio, Danillo Lange, and JaumePujante. A reinforcement learning control in hot stamp-ing for cycle time optimization. Materials, 15(14):4825,2022. → pages 20[NSHJ21] L Di Natale, B Svetozarevic, P Heer, and C N Jones.Deep reinforcement learning for room temperature con-trol: a black-box pipeline from data to policies. Journalof Physics: Conference Series, 2042(1):012004, nov 2021.→ pages 20[OF21] Francis Ogoke and Amir Barati Farimani. Thermal con-trol of laser powder bed fusion using deep reinforcement168Bibliographylearning. Additive Manufacturing, 46:102033, 10 2021. →pages 16, 20[PB22] Marcel Panzer and Benedict Bender. Deep reinforcementlearning in production systems: a systematic literaturereview. International Journal of Production Research,60(13):4316–4341, 2022. → pages 21[PNMS21] Jeet P Patil, Vilas Nandedkar, Sushil Mishra, andSandip K Saha. Transient thermal analysis of close pres-sure thermoforming process. Journal of ManufacturingProcesses, 62:513–522, 2021. → pages 17[pre] Thermoformed plastics market - global in-dustry analysis, size, share, growth, trends,regional outlook, and forecast 2022-2030”,https://www.precedenceresearch.com/thermoformed-plastics-market: :text=by→ pages xiv, 4[PRL+21] Bumsoo Park, Alexandra R. Rempel, Alan K.L. Lai,Julianna Chiaramonte, and Sandipan Mishra. Rein-forcement learning for control of passive heating andcooling in buildings∗∗this work was funded by nationalscience foundation cbet-1804218. IFAC-PapersOnLine,54(20):907–912, 2021. Modeling, Estimation and Con-trol Conference MECC 2021. → pages 20[PS02] Rohit S. Patwardhan and Sirish L. Shah. Issues in per-formance diagnostics of model-based controllers. Journalof Process Control, 12:413–427, 4 2002. → pages 16[PTvZ+13] R. H. Poelma, S. Tarashioon, H. W. van Zeijl, S. Gold-bach, J. L. J. Zijl, and G. Q. Zhang. Multi-led pack-age design, fabrication and thermal analysis. Journal ofSemiconductors, 34:054002, 5 2013. → pages xiv, 29[PV10] W Paplham and J Victor. Process for thermoformingacrylic polymer employing foam as a mold and articleformed there from. WO 2010/068787 Al, 2010. → pages2169Bibliography[RCN+21] Milad Ramezankhani, Bryn Crawford, Apurva Narayan,Heinz Voggenreiter, Rudolf Seethaler, and Abbas S. Mi-lani. Making costly manufacturing smart with transferlearning under limited data: A case study on compositesautoclave processing. Journal of Manufacturing Systems,59:345–354, 4 2021. → pages 11[RFLP22] Fre´de´ric Rosin, Pascal Forget, Samir Lamouri, andRobert Pellerin. Enhancing the decision-making processthrough industry 4.0 technologies. Sustainability, 14:461,1 2022. → pages 11[RHSM23] Milad Ramezankhani, Mehrtash Harandi, RudolfSeethaler, and Abbas S. Milani. Smart manufactur-ing under limited and heterogeneous data: a sim-to-realtransfer learning with convolutional variational autoen-coder in thermoforming. International Journal of Com-puter Integrated Manufacturing, pages 1–19, 9 2023. →pages 12[Ric93] Jacques Richalet. Industrial applications of model basedpredictive control. Automatica, 29(5):1251–1274, 1993.→ pages 18, 76[Rl-] Growth in search results of reinforcementlearning on google from 2007 to 2022,https://trends.google.com/trends/explore?date=2010-02-08→ pages xviii, 105[RM18] A. Rashidi and A.S. Milani. Passive control of wrinkles inwoven fabric preforms using a geometrical modificationof blank holders. Composites Part A: Applied Scienceand Manufacturing, 105:300–309, 2 2018. → pages 43[RNSM21] Milad Ramezankhani, Apurva Narayan, RudolfSeethaler, and Abbas S. Milani. An active transferlearning (atl) framework for smart manufacturing withlimited data: Case study on material transfer in com-posites processing. pages 277–282. IEEE, 5 2021. →pages 11[Ros79] SR Rosen. A history of the growth of the thermoformingindustry. Nature, 19:381–382, 1979. → pages 2170Bibliography[SAW+23] Damoon Soudbakhsh, Anuradha M. Annaswamy, YanWang, Steven L. Brunton, Joseph Gaudio, Heather Hus-sain, Draguna Vrabie, Jan Drgona, and Dimitar Filev.Data-driven control: Theory and applications. pages1922–1939. IEEE, 5 2023. → pages 16[SB18] R. Sutton and A. Barto. Reinforcement Learning: AnIntroduction. The MIT Press, Cambridge, 2nd edition,2018. → pages 99[Sch72] Lawrence Raymond Schmidt. Biaxial stretching of heat-softened plastic sheets. University of Colorado at Boulder,1972. → pages 45[Sch03] F Schmidt. Modelling of infrared heating of thermoplas-tic sheet used in thermoforming process. Journal of Ma-terials Processing Technology, 143-144:225–231, 12 2003.→ pages 13[SGM+21] Katrine Seel, Esten I Grøtli, Signe Moe, Jan T Gravdahl,and Kristin Y Pettersen. Neural network-based modelpredictive control with input-to-state stability. In 2021American Control Conference (ACC), pages 3556–3563.IEEE, 2021. → pages 18, 76[SIV22] Seppo Sierla, Heikki Ihasalo, and Valeriy Vyatkin. A re-view of reinforcement learning applications to control ofheating, ventilation and air conditioning systems. Ener-gies, 15(10), 2022. → pages 20[SLA+15] John Schulman, Sergey Levine, Pieter Abbeel, MichaelJordan, and Philipp Moritz. Trust region policy opti-mization. In International conference on machine learn-ing, pages 1889–1897. PMLR, 2015. → pages xviii, 107[SLX+10] Yun-De Shen, Zhen-Zhe Li, Dong-Ji Xuan, Kwang-SuHeo, and Seoung-Yun Seol. Time-dependent optimalheater control using analytic and numerical methods. In-ternational Journal of Precision Engineering and Manu-facturing, 11:77–81, 2 2010. → pages 13171Bibliography[SPG+16] C Schu¨ller, D Panozzo, A Grundho¨fer, H Zimmer,E Sorkine, and O Sorkine-Hornung. Computational ther-moforming. ACM Trans Graph, 35:1–9, 2016. → pages2[SSA+17] MFM Sharif, AA Saad, MZ Abdullah, FC Ani, MYT Ali,MK Abdullah, et al. A study on thermoforming processof stretchable circuit and its performance in manufactur-ing of automotive lighting. In AIP Conf Proc, volume1901, pages 1–5, 2017. → pages 2[SSR23] Sajad Salavatidezfouli, Giovanni Stabile, and GianluigiRozza. Deep reinforcement learning for the heat trans-fer control of pulsating impinging jets. arXiv preprintarXiv:2309.13955, 2023. → pages 21[Sut99] Richard S. Sutton. Open theoretical questions in rein-forcement learning. In European Conference on Compu-tational Learning Theory, 1999. → pages 101[SWD+17] John Schulman, Filip Wolski, Prafulla Dhariwal, AlecRadford, and Oleg Klimov. Proximal policy optimizationalgorithms. arXiv preprint arXiv:1707.06347, 2017. →pages xviii, 107, 108, 109[T-S] T-sim thermoforming company, https://www.t-sim.com/. → pages xiv, 6[Thr96] J. L. Throne. Technology of thermoforming. HanserGardner Publications, 1996. → pages 11, 13[Thr08] James L Throne. Understanding Thermoforming, pagesI–XIII. Hanser, 2008. → pages 2, 11[TMAS17] D Thelen, J Martin, S Allen, and S SA. Thermoformedpolystyrene products. 002:7, 2017. → pages 2[tUN15] the United Nations. The sustainable development goals(sdgs) [online]: https://www.undp.org/, 2015. → pages5[VBHP23] U. Veyna, X. Blasco, J.M. Herrero, and A. Pajares. Pa-rameter uncertainty modeling for multiobjective robust172Bibliographycontrol design. application to a temperature control sys-tem in a proton exchange membrane fuel cell. Engineer-ing Applications of Artificial Intelligence, 119:105758, 32023. → pages 16[VPY+15] J Vanfleteren, B Plovie, Y Yang, J De Smet, R Ver-plancke, F Bossuyt, et al. Free-form 2.5d thermoplasticcircuits using one-time stretchable interconnections. InMater. Res. Soc. Symp. Proc., volume 1798, pages 1–6,2015. → pages 2[WED+17] Jacob Wilson, Scott Everett, Rickey Dubay, Soheil S.Parsa, and Mike Tyler. Spatial predictive control using athermal camera as feedback. Measurement, 109:384–393,10 2017. → pages 13, 15, 16, 58[WHB13] Peng Wang, Nahie`ne Hamila, and Philippe Boisse. Ther-moforming simulation of multilayer composites with con-tinuous fibres and thermoplastic matrix. Composites PartB: Engineering, 52:127–136, 9 2013. → pages 12[Win20] Phil Winder. Reinforcement learning. O’Reilly Media,2020. → pages 20[WKM+22] Simon Wagner, Fabian Kayatz, Manuel Mu¨nsch,Cedric W. Sanjon, Marek Hauptmann, and Antonio Del-gado. Numerical modeling of forming air impact thermo-forming. The International Journal of Advanced Manu-facturing Technology, 120:4917–4933, 6 2022. → pages12[WLPL22] Hao Wang, Xiaobai Li, Mark Phipps, and Bo Li. Numer-ical and experimental study of hot pressing technique forresin-based friction composites. Composites Part A: Ap-plied Science and Manufacturing, 153:106737, 2 2022. →pages 12[WNW22] Hongda Wu, Ali Nasehzadeh, and Ping Wang. A deepreinforcement learning-based caching strategy for iot net-works with transient data. IEEE Transactions on Vehic-ular Technology, 71:13310–13319, 12 2022. → pages 110173Bibliography[WRR+17] The´ophane Weber, Se´bastien Racaniere, David P Re-ichert, Lars Buesing, Arthur Guez, Danilo JimenezRezende, Adria Puigdomenech Badia, Oriol Vinyals,Nicolas Heess, Yujia Li, et al. Imagination-augmentedagents for deep reinforcement learning. arXiv preprintarXiv:1707.06203, 2017. → pages xviii, 107[WWJ+03] M.K. Warby, J.R. Whiteman, W.-G. Jiang, P. Warwick,and T. Wright. Finite element simulation of thermoform-ing processes for polymer sheets. Mathematics and Com-puters in Simulation, 61:209–218, 1 2003. → pages 12[WZZ+19] K Wu, Q Zhou, H Zou, K Leng, Y Zeng, and Z Wu.High precision thermoforming 3d-conformable electron-ics with a phase-changing adhesion interlayer. Microma-chines, 10, 2019. → pages 2[XHB19] Hu Xiong, Nahie`ne Hamila, and Philippe Boisse. Consol-idation modeling during thermoforming of thermoplasticcomposite prepregs. Materials, 12:2853, 9 2019. → pages12[XMW20] Tengyang Xie, Yifei Ma, and Yu-Xiang Wang. Towardsoptimal off-policy evaluation for reinforcement learningwith marginalized importance sampling, 2020. → pages101[Yan08] Shuonan Yang. Cycle-to-cycle control of plastic sheetheating on the aaa thermoforming machine. 2008. →pages 17[YBD02] A. Yousefi, A. Bendada, and R. Diraddo. Improved mod-eling for the reheat phase in thermoforming through anuncertainty treatment of the key parameters. PolymerEngineering & Science, 42:1115–1129, 5 2002. → pages13, 15[YCH+08] Tiebao Yang, Xiang Chen, Henry Hu, Yeou-Li Chu, andPatrick Cheng. A fuzzy pid thermal control systemfor casting dies. Journal of Intelligent Manufacturing,19:375–382, 8 2008. → pages 16174Bibliography[YCH+23] Zhipeng Yuan, Liheng Chen, Hasiaoqier Han, Limin Ren,Shuai Liu, and Renxin Wang. Optimal design of thermalcontrol system for space optical remote sensor based onnsga-ii and opto-mechanical-thermal integration analysis.Case Studies in Thermal Engineering, 43:102813, 3 2023.→ pages 16[YH04] Chyan Yang and Shiu-Wan Hung. Modeling and opti-mization of a plastic thermoforming process. Journal ofReinforced Plastics and Composites, 23:109–121, 1 2004.→ pages 12[YQZ+21] Liang Yu, Shuqi Qin, Meng Zhang, Chao Shen, TaoJiang, and Xiaohong Guan. A review of deep reinforce-ment learning for smart building energy management.IEEE Internet of Things Journal, 8:12046–12063, 8 2021.→ pages 102[ZFLL24] Ye Zhang, Zhuangdong Fang, Changyou Li, and ChengjieLi. Deep-learning-based model predictive control ofan industrial-scale multistate counter-flow paddy dryingprocess. Foods, 13(1), 2024. → pages 18, 76[ZH10] Bing Zhu and Wei Huo. Trajectory linearization controlfor a quadrotor helicopter. pages 34–39. IEEE, 6 2010.→ pages 62[Zha97] Qinghua Zhang. Using wavelet network in nonparamet-ric estimation. IEEE Transactions on Neural networks,8(2):227–236, 1997. → pages 77[Zhu98] Yucai Zhu. Multivariable process identification for mpc:the asymptotic method and its applications. Journal ofProcess Control, 8(2):101–115, 1998. → pages 18, 76[ZMSZ20] Xuezhou Zhang, Yuzhe Ma, Adish Singla, and XiaojinZhu. Adaptive reward-poisoning attacks against rein-forcement learning, 2020. → pages 101[ZPWZ13] Yucai Zhu, Rohit Patwardhan, Stephen B. Wagner, andJun Zhao. Toward a low cost and high performance mpc:The role of system identification. Computers & Chemical175BibliographyEngineering, 51:124–135, 2013. CPC VIII. → pages 18,76[ZY20] Hongming Zhang and Tianyang Yu. Taxonomy of Rein-forcement Learning Algorithms, pages 125–133. SpringerSingapore, Singapore, 2020. → pages xviii, 107, 108176Appendix177Appendix AEvaluation of the PI gainmatrices and closed-loopsystem polesThe PI gain matrices kP1 and kI1 for the sheet temperature control loop(reference temperature in 4.8) evaluated at time step t = 1 second, withω = 13500I15×15 and ζ = 10I15×15.kP1 =[k(1)P1 k(2)P1 k(3)P1]wherek(1)P1 =9369325176 −2663974907 576119859 −14246737088 3590139702−1064248773 9222542943 −2153694883 1236728685 −13763828083123927625 −1286343621 9128105480 −125361172 1557642243−11815299875 2916595308 −559424164 27218661133 −6481656760963377180 −11436934004 2271886276 −2033936391 26152171504−92168826 1232234093 −11477409274 199910500 −26425535387713973638 −1813828394 323231568 −23223162765 5233244278−565759618 7440311266 −1395501368 1498317467 −2220127413545628714 −732849908 7480055788 −126024246 2003866365−3800549188 883263502 −156901200 12636523220 −2801008169270220831 −3657566385 677143485 −784188101 12059864038−21913371 351479680 −3683608513 64600333 −10536602581167704784 −292944010 50273318 −4093077985 952765012−104319495 1136234809 −230523475 295094772 −39229608935244220 −128198551 1132491415 −24254509 382002619178Appendix A. Evaluation of the PI gain matrices and closed-loop system polesk(2)P1 =−712573646 11485655032 −2780092951 526416218 −67507420282838365418 −916688220 11042231373 −2183774450 534989230−13898759069 84085717 −1164661116 11205288690 −499628331236229621 −27590492767 6299089235 −1147695803 17659179012−5039122332 1852526154 −26368237985 4853821955 −117318109526521244660 −163615188 2453177157 −26890201061 102909821−931865457 33943947855 −7772263923 1422773086 −271705659264002038887 −2298623390 32465769574 −5992395737 1814025827−22595999763 207713387 −3039397000 33079606154 −158794497489229653 −23603256588 5319837618 −953262955 27243238265−2132557766 1530760833 −22532877930 4083144862 −203488784812287583020 −130081176 2040940857 −22995325682 199994264−171531379 8414151776 −1980651590 361108380 −12511916829736774175 −627750919 8069747473 −1543624829 1037702811−3982603750 52400253 −806330056 8201719660 −101257139k(3)P1 =1633888623 −308758375 2326356104 −604952501 115147419−6499482187 1275832246 −220663955 2270799329 −477873220682194763 −6571759936 18751790 −271715563 2256898514−4018250431 722811478 −6326649585 1530295786 −28466661116905529033 −3077882907 494614733 −6120763445 1182309068−1557246053 17172097240 −45550709 634731946 −61307257676206319924 −1121559989 10631622307 −2573414565 474942280−26026933176 4756182536 −833136777 10295208994 −19887437672411907747 −26426274843 74271228 −1069583102 10303164937−6492161149 1235667711 −13562057569 3412771815 −66433754626199192335 −5038073839 1156574611 −13173845596 2663685386−2646470454 26524801204 −115037702 1469139339 −131611586083099154756 −602726234 9430891792 −2683570216 578238521−12075089020 2433012547 −1069654772 9296381526 −21629666681319548267 −12191971800 124405187 −1295076902 9175982225179Appendix A. Evaluation of the PI gain matrices and closed-loop system poleskI1 =[k(1)I1 k(2)I1 k(3)I1]wherek(1)I1 =246722 −71963 14926 −389598 102312−28293 248267 −55927 35104 −3872153290 −34877 236153 −3583 45082−321323 81903 −15107 744392 −18430926902 −318360 61339 −57599 735507−2591 34990 −306522 5649 −76258216185 −52505 8999 −644622 151639−16265 213521 −38785 43649 −6341551330 −21447 205761 −3704 59289−108084 25938 −4420 356566 −826257880 −106556 19078 −23300 350344−646 10436 −102794 1934 −3178933489 −8634 1435 −116456 28331−3013 33357 −6505 8793 −114957169 −3791 31865 −728 11587k(2)I1 =−19411 322861 −81700 14761 −18841977211 −26836 319968 −61037 15531−371879 2494 −34820 307682 −145633508 −768976 183514 −31884 489944−136725 54088 −757342 134602 −34075709624 −4818 72917 −732210 3005−25815 940550 −224573 39106 −744592110695 −66307 926696 −164767 52000−613491 5992 −89350 895577 −457813798 −659472 155513 −26596 737693−60018 44955 −648861 113632 −57073339030 −3866 60965 −627677 5599−4868 237894 −58630 10216 −33935820889 −18568 235240 −43462 29322−110778 1582 −24336 226494 −2887180Appendix A. Evaluation of the PI gain matrices and closed-loop system polesk(3)I1 =47524 −8582 64261 −17157 3146−186591 35437 −6166 64110 −1303320176 −179487 535 −7743 61264−116249 20009 −173856 43269 −7733482308 −85090 13866 −171947 32177−45951 466269 −1278 18100 −165575177215 −30647 289011 −72008 12794−733194 129843 −23131 286036 −5359370214 −708912 2072 −30196 275324−182631 33195 −360206 93301 −17497728886 −135465 31403 −357497 70203−75563 703218 −3133 40533 −34374287612 −16327 243567 −70981 14708−336664 65781 −27882 245067 −5515037948 −323776 3237 −34378 233124181Appendix A. Evaluation of the PI gain matrices and closed-loop system polesThe evaluated closed-loop system poles of the PI controller (4.25) at thetime step of t = 1 second for the kP1 and kI1 values calculated above.eigenvalues =−1.43216082668001× 10−5 + 0.00000000000000i−1.43216082667999× 10−5 + 0.00000000000000i−1.43216082668003× 10−5 + 1.35564970415680× 10−19i−1.43216082668003× 10−5 − 1.35564970415680× 10−19i−1.43216082667999× 10−5 + 2.77799367158206× 10−19i−1.43216082667999× 10−5 − 2.77799367158206× 10−19i−1.43216082667998× 10−5 + 0.00000000000000i−1.43216082668002× 10−5 + 1.70082865950005× 10−19i−1.43216082668002× 10−5 − 1.70082865950005× 10−19i−1.43216082668000× 10−5 + 0.00000000000000i−0.00569996410601890 + 0.00000000000000i−0.00569996410601892 + 3.22356480643232× 10−18i−0.00569996410601892− 3.22356480643232× 10−18i−0.00569996410601891 + 0.00000000000000i−0.00569996410601893 + 0.00000000000000i−0.00569996410601893 + 0.00000000000000i−0.00569996410601893 + 0.00000000000000i−0.00569996410601894 + 4.04510585318187× 10−18i−0.00569996410601894− 4.04510585318187× 10−18i−0.00569996410601891 + 0.00000000000000i−0.00569996410601892 + 4.75073589413134× 10−18i−0.00569996410601892− 4.75073589413134× 10−18i−1.43216082667996× 10−5 + 0.00000000000000i−1.43216082668001× 10−5 + 0.00000000000000i−0.00569996410601891 + 0.00000000000000i−1.43216082668001× 10−5 + 0.00000000000000i−0.00569996410601891 + 2.25347189810574× 10−18i−0.00569996410601891− 2.25347189810574× 10−18i−1.43216082668005× 10−5 + 0.00000000000000i−1.43216082668001× 10−5 + 0.00000000000000i182
Thesis/Dissertation
eng
Mechanical Engineering
Vancouver : University of British Columbia Library
Attribution-NonCommercial-NoDerivatives 4.0 International
Graduate
Seethaler, Rudolf
Physics-based simulation and reinforcement learning control of the heating phase in the thermoforming process
Text