MULTI-OBJECTIVE GAUSSIAN PROCESS APPROACH FOR ROBUST OPTIMIZATION AND PREDICTION OF CARBONIZATION PROCESS by Milad Ramezankhani A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF APPLIED SCIENCE in THE COLLEGE OF GRADUATE STUDIES (MECHANICAL ENGINEERING) THE UNIVERSITY OF BRITISH COLUMBIA (Okanagan) August 2017 © Milad Ramezankhani, 2017 ii The following individuals certify that they have read, and recommend to the College of Graduate Studies for acceptance, a thesis/dissertation entitled: MULTI-OBJECTIVE GAUSSIAN PROCESS APPROACH FOR ROBUST OPTIMIZATION AND PREDICTION OF CARBONIZATION PROCESS submitted by Milad Ramezankhani in partial fulfillment of the requirements of the degree of Master of applied science Dr. Abbas Milani and Dr. Rudolf Seethaler Supervisors Dr. Joshua Brinkerhoff Supervisory Committee Member Dr. Zheng Liu Supervisory Committee Member Dr. Gordon Lovegrove University Examiner External Examiner Additional Committee Members include: Supervisory Committee Member Supervisory Committee Member iii Abstract Carbon fibers are high-performance and high-strength reinforcement materials in advanced industries such as aerospace, automotive and energy sector. This class of materials is often derived from polyacrylonitrile (PAN) fiber precursor. The conversion of polyacrylonitrile precursor to carbon fibers is essentially comprised of two major stages; namely, stabilization and carbonization. The carbonization is a very energy-consuming and expensive step due to its high temperature requirements. In order to achieve desirable physical and mechanical properties of carbon fibers through this step, a large amount of energy is required. A cost-effective approach to optimize energy consumption of this process, however, is the use of predictive modelling techniques. In that goal, herein, a Gaussian Process (GP) approach is proposed to firstly, predict multiple mechanical properties of carbon fibers in the presence of manufacturing noise and secondly, optimize them under a minimum energy consumption criterion and a range of process constraints and bounds. The resulting model for each property of carbon fibers consists of two Gaussian Process models. The first model describes the mean value of the property and the second one predicts its standard deviation. The proposed Gaussian Process approach is compared to a traditional regression approach using the measurements obtained from the Carbon Nexus pilot plant in Australia. The Gaussian Process approach clearly proved to be more effective in terms of both prediction accuracy and robustness. Through employing Gaussian Process models, the modulus and tensile strength mean values of carbon fibers along with their standard deviations (STD) were successfully predicted under different process conditions. Squared exponential covariance and linear mean functions were proven to be most suitable in constructing the Gaussian Process in the performed case study. It iv was found that the modulus and tensile strength responses do not have an evident correlation with respect to each other, hence a multi-objective optimization approach was developed to acquire overall optimum process conditions. To estimate the trade-off between fiber material properties under the multi-objective optimization problem, a standard as well as an adaptive weighted sum method was applied under various constraints and bounds, including the energy consumed during carbonization. v Preface A summary of the work done in this thesis has been submitted as an oral presentation in the following conference. My responsibilities included model development, statistical analyses, writing and presenting the article. Dr. Abbas S. Milani and Dr. Rudolf Seethaler co-supervised the project and co-authored the article. • M. Ramezankhani, B. Crawford, H. Khayyam, M. Naebe, R. Seethaler, A.S. Milani (2017) Multi-objective Gaussian process approach for robust optimization and prediction of carbonization process. 20th International Conference on Composite Structures (ICCS20), 4-7, September, Paris, France A journal article is also under submission from Chapters 3 and 4 of the thesis. Experimental data employed in the research has been obtained through collaboration with Dr. Hamid Khayyam and Dr. Minoo Naebe at Carbon Nexus in Australia. vi Table of Contents Abstract ......................................................................................................................................... iii Preface .............................................................................................................................................v Table of Contents ......................................................................................................................... vi List of Tables ..................................................................................................................................x List of Figures ............................................................................................................................... xi List of Symbols ......................................................................................................................... xviii List of Abbreviations ...................................................................................................................xx Acknowledgements .................................................................................................................... xxi Dedication .................................................................................................................................. xxii Chapter 1: Background and Thesis Organization ......................................................................1 Introduction ..................................................................................................................... 1 Motivation and Objective ............................................................................................... 2 Thesis framework............................................................................................................ 2 Chapter 2: Literature Review .......................................................................................................5 Carbon fibers ................................................................................................................... 5 2.1.1 Early work on carbon fibers’ production .................................................................... 5 2.1.2 Steps in carbon fiber manufacturing ........................................................................... 6 2.1.2.1 Polymerization .................................................................................................... 7 2.1.2.2 Stabilization ........................................................................................................ 7 2.1.2.3 Carbonization ...................................................................................................... 8 2.1.2.3.1 Low temperature (LT) carbonization .......................................................... 10 vii 2.1.2.3.2 High temperature (HT) carbonization ......................................................... 11 2.1.2.4 The effect of the process parameters on carbon fibers’ mechanical properties 13 2.1.3 State-of-art on modeling of carbon fiber production process ................................... 14 Gaussian Processes ....................................................................................................... 15 2.2.1 Bayes’ theorem ......................................................................................................... 16 2.2.2 Gaussian Process basics ............................................................................................ 17 2.2.3 Covariance function .................................................................................................. 22 2.2.4 Model Selection ........................................................................................................ 25 2.2.4.1 Negative Marginal Likelihood (NML) ............................................................. 26 2.2.4.2 Cross validation ................................................................................................ 27 2.2.5 Gaussian Process application example ..................................................................... 28 Chapter 3: Model Development ..................................................................................................30 Overview ....................................................................................................................... 30 Case study introduction................................................................................................. 30 3.2.1 Focus ......................................................................................................................... 31 3.2.2 Carbon Nexus institute .............................................................................................. 31 3.2.3 A black-box problem ................................................................................................ 32 Gaussian Process developments ................................................................................... 34 3.3.1 Modulus mean ........................................................................................................... 35 3.3.2 Modulus standard deviation ...................................................................................... 38 3.3.3 Tensile strength mean ............................................................................................... 39 3.3.4 Tensile strength standard deviation (STD) ............................................................... 40 3.3.5 Energy consumption ................................................................................................. 41 viii Verification of the selected Gaussian process (GP) models ......................................... 57 3.4.1 Modulus mean ........................................................................................................... 57 3.4.2 Modulus standard deviation (STD) ........................................................................... 57 3.4.3 Tensile strength mean ............................................................................................... 60 3.4.4 Tensile strength standard deviation (STD) ............................................................... 61 3.4.5 Energy consumption ................................................................................................. 63 3.4.6 Statistical evaluation of the Gaussian process (GP) models ..................................... 64 Summary of findings..................................................................................................... 65 Chapter 4: Multi-objective Optimization of Carbonization Process ......................................66 Overview ....................................................................................................................... 66 Multi-objective optimization formulation..................................................................... 66 4.2.1 Constraints and bounds ............................................................................................. 67 4.2.2 Underlying generic trade-offs in the carbonization process optimization problem .. 68 Customer-driven optimization ...................................................................................... 72 Selected optimization scenarios .................................................................................... 72 4.4.1 Scenario 1: Modulus optimization ............................................................................ 72 4.4.2 Scenario 2: Concurrent optimization of the modulus and tensile strength ............... 74 4.4.2.1 Adaptive weighted sum (AWS) method ........................................................... 78 4.4.3 Scenario 3: concurrent optimization of the modulus and tensile strength ................ 80 4.4.4 Scenario 4: energy consumption optimization .......................................................... 83 Summary of findings..................................................................................................... 85 Chapter 5: Discussions ................................................................................................................90 Overview ....................................................................................................................... 90 ix 5.1.1 How the Gaussian process (GP) would perform compared to a standard regression model..……...…………………………………………………………..90 5.1.2 The effect of stabilization process temperature on the carbon fiber mechanical properties ................................................................................................................. 93 5.1.3 The benefit of using two-zone furnace in the carbonization process ........................ 95 Summary of findings..................................................................................................... 96 Chapter 6: Conclusions and Future Work ................................................................................98 Conclusions ................................................................................................................... 98 Future work recommendations ................................................................................... 102 Bibliography ...............................................................................................................................104 Appendices ..................................................................................................................................111 Appendix A : Detailed information of the carbonization process .......................................... 111 x List of Tables Table 2.1 Gasses evolved from different temperatures during carbonization process [15] ............. 11 Table 3.1 Negative marginal likelihood (NML) values of different covariance and mean functions used for modulus mean prediction .................................................................................... 35 Table 3.2 Negative marginal likelihood (NML) values of combined (composite) covariance and mean functions used for modulus mean prediction .......................................................... 36 Table 3.3 Summary of the modulus mean prediction using the Gaussian process (GP) model ....... 58 Table 3.4 Summary of the modulus standard deviation (STD) predictions using Gaussian process (GP) model ........................................................................................................................ 59 Table 3.5 Summary of the tensile strength mean prediction using Gaussian process (GP) model .. 61 Table 3.6 Summary of the tensile strength standard deviation (STD) prediction using Gaussian process (GP) model ........................................................................................................... 62 Table 3.7 Summary of the carbonization process energy consumption prediction using Gaussian process (GP) model ........................................................................................................... 64 Table 3.8 Statistical measures’ summary of the carbonization process’s Gaussian process (GP) prediction models .............................................................................................................. 65 Table 4.1 Pareto optimal solutions using different weights under scenario 2. ................................. 76 Table 4.2 Pareto optimal solutions of scenario 3 using different weights ........................................ 82 Table 4.3 Summary of the four optimization scenarios for carbonization process .......................... 87 xi List of Figures Figure 1.1 Organization of the thesis........................................................................................... 4 Figure 2.1 Carbon fiber manufacturing process steps. Adopted from [21] ................................. 7 Figure 2.2 Chemical reactions during the stabilization process. Adopted from [15] ................... 8 Figure 2.3 Chemical reactions during the carbonization process [24] ........................................ 9 Figure 2.4 Diagram of nitrogen consumption in carbonization process [25] ............................ 10 Figure 2.5 High temperature furnace in the carbonization process [28] ................................... 12 Figure 2.6 Typical tensile strength (GPa) of PAN-based carbon fibers over maximum heat treatment temperature [32] ...................................................................................... 14 Figure 2.7 A schematic of Gaussian Process regression [43]..................................................... 19 Figure 2.8 (a) Prior and (b) posterior examples in the Gaussian process (GP) modelling, assuming data are true means. Lines represent the mean value, plus signs are the inputs, and the gray are is the confidence interval[10] ............................................ 20 Figure 2.9 Indication of the mean value and the confidence interval in a Gaussian process (GP) [13] .................................................................................................................. 22 Figure 2.10 (a) Three prior functions with Squared Exponential covariance function. (b) Three prior functions with Ornstein-Uhlenbeck covariance function [45] ............. 23 Figure 2.11 Gaussian process (GP) prediction model using squared exponential covariance function with three different characteristic length scales. Adopted from [10] ........ 25 Figure 2.12 The 545 observation data regarding the carbon dioxide concentration in the Hawaii’s atmosphere, along with the 95% confidence interval of the Gaussian process (GP) prediction model shown as hashed region [46] ................................. 29 xii Figure 3.1 Schematic of inputs and outputs for the carbonization process under study ............ 33 Figure 3.2 PRESS values for different modulus mean Gaussian processes (GP) ..................... 37 Figure 3.3 Modulus mean prediction as a function of time using squared exponential (SE) covariance function and Linear mean function at three different furnace temperature conditions. (a) when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1175°C. It is worth mentioning that in the areas where more observed data exists, the confidence interval become narrower, representing a higher level of prediction reliability around those points. ................................................................................ 42 Figure 3.4 Modulus mean prediction as a function of time using periodic covariance function and Constant mean function at three different furnace temperature conditions. (a) when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1475°C.................................................................... 43 Figure 3.5 Negative marginal likelihood values of different Gaussian process (GP) models for modulus standard deviation (STD) prediction ................................................... 44 Figure 3.6 PRESS values of different Gaussian process (GP) models for modulus standard deviation (STD) prediction ...................................................................................... 45 Figure 3.7 Modulus standard deviation (STD) prediction versus time using squared exponential (SE) covariance function and Linear mean function at three different furnace temperature conditions. (a) happens when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) is the modulus prediction at 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) shows the result at 𝑻𝟏=1175°C and 𝑻𝟐=1475°C. .................................................................... 46 xiii Figure 3.8 Negative marginal likelihood values of different Gaussian process GP models for tensile strength mean prediction .............................................................................. 47 Figure 3.9 PRESS values of different Gaussian process (GP) models for tensile strength mean prediction ....................................................................................................... 48 Figure 3.10 Tensile strength mean prediction versus time using squared exponential (SE) covariance function and Linear mean function at three different furnace temperature conditions. (a) happens when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1475°C. ............ 49 Figure 3.11 Tensile strength mean prediction versus time using squared exponential (SE) covariance function and constant mean function at three different furnace temperature conditions. (a) when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1475°C. ............ 50 Figure 3.12 Negative marginal likelihood values of Gaussian process (GP) models for tensile strength standard deviation (STD) prediction .............................................. 51 Figure 3.13 PRESS values of different Gaussian process (GP) models for tensile strength standard deviation (STD) prediction ....................................................................... 52 Figure 3.14 Tensile strength standard deviation (STD) prediction versus time using squared exponential (SE) covariance function and Linear mean function at three different furnace temperature conditions. (a) when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1475°C. ............ 53 Figure 3.15 Negative marginal likelihood values of Gaussian process (GP) models for energy consumption prediction ............................................................................... 54 xiv Figure 3.16 PRESS values of different Gaussian process (GP) models for energy consumption prediction ........................................................................................... 55 Figure 3.17 Energy consumption prediction versus time using squared exponential (SE) covariance function and constant mean function at three different furnace temperature conditions. (a) ...................................................................................... 56 Figure 3.18 Experimental data against prediction values for the modulus mean Gaussian process (GP) model ................................................................................................. 58 Figure 3.19 Experimental data against prediction values for the modulus standard deviation (STD) Gaussian process (GP) model. ...................................................................... 59 Figure 3.20 Experimental data against prediction values of tensile strength mean Gaussian process (GP) model. ................................................................................................ 60 Figure 3.21 Experimental data versus prediction values of tensile strength standard deviation (STD) Gaussian process (GP) model ....................................................... 62 Figure 3.22 Energy consumption prediction values against the experimental data. ................. 63 Figure 4.1 Experimental Pareto graphs of (a) modulus mean versus tensile strength mean, (b) modulus mean versus energy consumption, and (c) tensile strength mean versus energy consumption. Red points represent an experimental data with the highest mean modulus and mean tensile strength along with a low energy consumption in the carbonizationprocess………...……………………………......70 Figure 4.2 Experimental Pareto graphs of (a) modulus STD versus tensile strength STD, (b) modulus STD versus energy consumption, and (c) tensile strength STD versus energy consumption. Red points represent an experimental data with the highest xv mean modulus and mean tensile strength along with a low energy consumption in the carbonization process. ................................................................................... 71 Figure 4.3 The convex region with non-constant curvature in the standard weighted sum method in scenario 2 ................................................................................................ 77 Figure 4.4 The initial sub-optimization of AWS method in scenario 2 .................................... 79 Figure 4.5 AWS method is applied to smaller regions in order to obtain remaining optimal solutions in scenario 2 ............................................................................................. 80 Figure 4.6 Modulus mean against tensile strength mean plot of Pareto optimal solutions by adaptive weighted sum(AWS) method in in scenario 2. Process variables (T1, T2 and t) are also shown at each point, along with the ΔT which is T2-T1. ... 80 Figure 4.7 Modulus against tensile strength plot of Pareto optimal solutions under scenario 3. Values next to each pint indicate the relative weight of the modulus and tensile strength, respectively. ............................................................................ 83 Figure 4.8 Pareto optimal solutions along with experimental data in the performed scenarios 1-4, for mean mechanical properties and the energy consumption. The black circles are the observed data points and the triangles denote the optimum solution of each scenario. ........................................................................................ 88 Figure 4.9 Pareto graphs of Pareto optimal solutions along with experimental data in the performed scenarios 1-4, for standard deviations (STD) of mechanical properties and the energy consumption. The black circles are the observed data points and the triangles denote the optimum solution of each scenario. ................................... 89 Figure 5.1 Comparison of Gaussian process (GP) models versus two classic regression models for modulus mean prediction ...................................................................... 91 xvi Figure 5.2 Comparison of Gaussian process (GP) models versus classic regression models for tensile strength mean prediction ........................................................................ 92 Figure 5.3 Comparison of Gaussian process (GP) models versus classic regression models for modulus STD prediction .................................................................................... 92 Figure 5.4 Comparison of Gaussian process (GP) models versus classic regression models for tensile strength STD prediction ......................................................................... 92 Figure 5.5 Modulus and tensile strength of polyacrylonitrile (PAN) based and mesophase pitch based fibers against processing temperature [29]. .......................................... 93 Figure 5.6 The temperature of carbon fibers versus time in the carbonization process. 𝑻𝟎 is the stabilization process. The dashed dotted lines indicate the slope of the processes. ................................................................................................................. 95 Figure 5.7 Comparison of two carbonization processes with difference furnace designs: the blue line indicates the process with one zone furnace; the green line shows a two-zone furnace carbonization process .................................................................. 96 Figure 6.1 Recommended design chart (n-dimensional Pareto) for multi-objective optimization of mechanical properties of fibers, while the carbonization process energy is constrained and the robustness of fibers performance is taken into account ................................................................................................................... 102 Figure A.1 Modulus versus time plot of the carbonization process, using raw experimental data at T1=1125 °C …………….………………………………………………..110 Figure A.2 Tensile strength versus time plot of the carbonization process, using raw experimental data at T1=1125 °C……………………………………………...…111 xvii Figure A.3 Modulus versus zone-1 temperature plot of the carbonization process, using raw experimental data at t= 162(s)…..……………………….…………………..112 Figure A.4 Tensile strength versus zone-1 temperature plot of the carbonization process, using raw experimental data at=162(s)……………….…………………………113 xviii List of Symbols C Carbon atom N Nitrogen atom 𝑁2 Nitrogen gas Ar Argon gas 𝑝(𝑦1, 𝑦2) Joint probability of 𝑦1 and 𝑦2 𝑝(𝑦𝐴|𝑦𝐵) Conditional probability 𝑦 Response value 𝑥 Input vector 𝑥′ Test input 𝑤 Weight vector E Expected value 𝑓 Function value 𝑓∗ Function value at 𝑋∗ (Gaussian process prediction) 𝜀 Noise 𝜎𝑛2 Noise variance 𝜎𝑓2 Variance of signal ℓ Characteristic length scale 𝛿 Kronecker’s delta 𝑚(𝑥) The mean function of Gaussian process 𝑘(𝑥, 𝑥′) The covariance function of Gaussian process 𝑁 (𝑚(𝑥), 𝑘(𝑥, 𝑥′)) Gaussian distribution with mean function 𝑚 and covariance function 𝑘 xix 𝑉 Variance 𝐾𝑦 Covariance function of the noisy observations 𝑦 𝑌𝑖′ Respective prediction value in PRESS calculation 𝑌𝑖 Omitted observation in PRESS calculation 𝑇0 Stabilization process temperature 𝑇1 Furnace zone-1 temperature in carbonization process 𝑇2 Furnace zone-2 temperature in carbonization process t Carbonization process time 𝐸 Energy consumption of the carbonization process M Modulus mean UTS Ultimate tensile strength mean Mnorm Normalized modulus mean UTSnorm Normalized ultimate tensile strength mean wn Decision variable weight in weighted sum method xx List of Abbreviations CRN Composites Research Network PAN Polyacrylonitrile CF Carbon fibers HCl Hydrogen chloride LT Low temperature carbonization process HT High temperature carbonization process GP Gaussian Process STD Standard deviation NML Negative marginal likelihood PRESS predicted residual error sum of squares RSS residual sum of squares AWS Adaptive weighted sum xxi Acknowledgements I am sincerely grateful for the support of my co-supervisors, Dr. Abbas Milani and Dr. Rudolf Seethaler. They always believed in me and generously guided me throughout the project. Their practical supervision and optimism eased the difficulties during my research. I also want to thank them for giving me the chance to follow my own research procedure. The CRN group at UBC provides a magnificent learning environment for learning and doing research. I truly enjoyed the time spent with graduate students in the CRN laboratory. I want to especially thank Bryn Crawford and Armin Rashidi for their advice and comments. I would also like to thank my friends at CRN, Hossein Montazerian, Ronak Vahed and K M Safat Rashif for the quality times we had together in the past year. Finally, I am indebted to my lovely family for their invaluable patience, love and support. I could not have accomplished the project without their encouragement and assistance. xxii Dedication To my family, Reza, Mojgan and Mohammad, for all of their endless love and support.1 Chapter 1: Background and Thesis Organization Introduction In 1880, Tomas Alva Edison employed carbon fiber as a filament material for his wonder invention, the electric lamp [1]. Edison and his colleagues spent more than one year in their laboratory at Menle Park, New Jersey, USA working on an appropriate element material to be used inside a glass bulb under a high vacuum [2, 3]. It has been said that Edison and the team had tried 1600 different kind of materials such as paper, cloth, boxwood and coconut shells in order to find the proper material. Carbon fibers, known as a lightweight, high-stiffness and high-strength material type, have been developed for more than 45 years. During these years, significant enhancements have been applied to carbon fibers. New precursors improved the mechanical properties of the resulting carbon fibers. Much research and effort has been put into improving the tensile properties of fibers, and at the same time optimizing the costs and energy consumption of carbon fiber production processes [4, 5]. However, the existence of uncertainties in carbon fiber manufacturing can negatively affect the performance and robustness of the final fiber tows, and accordingly the analysis and decision-making steps for designers. The random behavior of hundreds of thousands of filaments in the production line is one of the primary causes of uncertainty [6, 7], which if not accounted for, it can make the design analysis of carbon fiber production complicated, time consuming and costly. Yet, the implementation of a powerful and flexible stochastic process model can facilitate the tasks of prediction and optimization of carbon fiber production for designers. 2 Gaussian process (GP) has been known as a useful probabilistic inference method, which effectively deals with complex data structures in the presence of uncertainty [8]. It is perceived that Gaussian processes provide robust prediction models and avoid over-fitting [8-10]. Moreover, Gaussian processes have been successfully implemented in modelling problems for many engineering cases and proved that they can perform better than classic regression methods [11-13]. Consequently, Gaussian processes may be deemed as an appropriate statistical framework also for the prediction and optimization of carbon fiber production processes under uncertain data structures. Motivation and Objective The main goal of this thesis is to develop a multi-objective GP framework to predict and optimize the tensile properties (namely the Young’s modulus and ultimate tensile strength) of carbon fibers as well as the energy consumption during the carbonization process under data uncertainty (noise) and numerous process bounds including those on workable furnace temperatures and heat treatments’ total time. The carbonization process is one of the most energy-consuming steps in carbon fiber production. The expected practical outcomes of this thesis are to show how Gaussian Processes can cope with the uncertainties that exist in the carbonization process at an industrial level, and how robust and accurate the resulting models are. Moreover, the advantages of using GP models compared to other classic regression methods will be shown. Thesis framework The thesis consists of six chapters. Chapter 2 reviews the literature regarding the different steps in carbon fiber production with the focus on the carbonization process and the resulting carbon fibers’ 3 properties. It also examines the mathematical background and the application of Gaussian Process (GP). Chapter 3 presents the developments of different GP models for the prediction of fibers’ tensile properties and the carbonization process energy consumption. The statistical evaluation methods are applied to the selected GP models in order to assess their accuracy and robustness. Finally, the developed GP models will be utilized in the next chapter for optimization tasks. Chapter 4 starts with introducing the optimization framework to be used for quality improvement of carbonization process, including the formulation, decision variables, constraints and bounds. Subsequently, four distinct design scenarios regarding the multi-objective optimization of the carbonization process are studied using the developed GP models. Moreover, a new application of an adaptive weight sum method is shown for capturing concave design points in the Pareto surface of the multi-objective optimization problem of the carbonization process. Chapter 5 presents discussions regarding the comparison of the GP and classic regression methods and the physical relations between the carbon fibers’ properties and the optimized carbonization process conditions. Lastly, chapter 6 provides a summary of the findings and outlines some potential future work recommendations. Figure 1.2 schematizes the organization of the thesis and how the chapters are related to each other. 4 Chapter 1:Background and thesis organizationChapter 2:Literature reviewChapter 3:GP models’ prediction and result regarding the carbonization process Chapter 4:Optimization model development and investigation of four optimization scenarios for the carbonization processChapter 5:DiscussionsChapter 6:Conclusions and future work Figure 1.1. Organization of the thesis 5 Chapter 2: Literature Review Carbon fibers Carbon fibers (CFs) are high-performance and high-strength reinforcement materials employed in advanced industries such as the aerospace, automotive and energy sectors. Although CFs suffer from the low resistance to oxidizers like hot air and fire, they have become hugely popular in engineering applications because of their exceptional reinforcing capacities and low weight. Some of the main advantages of carbon fibers are as follows [14, 15]: • They are available in many grades and structures. • They possess high modulus and excellent strength values. • They have high thermal stability in the absence of Oxygen gas along with high thermal conductivity. • They resist properly against chemical reactions. However, commercial CF products have some drawbacks, such as: • They are relatively expensive. • They strongly react with air’s oxygen at a temperature over 450°C. • Their properties can highly depend on the processing conditions applied during production. 2.1.1 Early work on carbon fibers’ production This class of materials is often derived from polyacrylonitrile (PAN) fiber precursor, which was primarily introduced by Shindo in 1961 as a material for CF production [16]. Along with the 6 research in Japan, W. Watt and his teammates in the Royal Aircraft Establishment at Farnborough, UK developed a PAN precursor that obtained a carbon fiber with the modulus of 380 GPa when the precursor was heat treated to a temperature of 2500°C [17, 18]. Later on, in 1969, Harwell succeeded in producing 300 meter long continuous fibers for UK commercial purposes. MS Blackie [19] also investigated the structure of carbon fibers and the effect of heat treatments using X-ray diffraction. 2.1.2 Steps in carbon fiber manufacturing The conversion of PAN precursor to carbon fibers is essentially comprised of four major stages; namely, polymerization, stabilization, carbonization and graphitization. Figure 2.1 illustrates the carbon fiber manufacturing process. The fibers normally go through a two-stage carbonization process (low and high temperature steps) after subjecting them to the oxidative stabilization process. The goal of these heat treatments, which operate at different controlled temperatures, is to remove all the atoms except carbon, which eventually results in producing carbon fibers with high modulus and tensile strength [20]. However, some studies showed that the existence of up to 5% nitrogen in the structure of the fibers improves the strength of the resulting product [20]. 7 Polymerization of PAN fibersStretchingStabilization (Air-250°C) LT Carbonization (Nitrogen,1100°C )HT Carbonization (Nitrogen,1800°C )Graphitization(Nitrogen,3000°C ) Surface treatmentsWinding Carbon fibers Figure 2.1. Carbon fiber manufacturing process steps. Adopted from [21] 2.1.2.1 Polymerization Polymerization is a process in which monomers bond to other similar molecules through chemical reaction, creating a polymer -- a giant network molecule. The PAN, which most carbon fibers are derived from, can be produced by the polymerization of acrylonitrile. The four main methods of polymerization, which employ free radical initiation, are: solution polymerization, bulk polymerization, emulsion polymerization and aqueous dispersion polymerization [15]. 2.1.2.2 Stabilization A stabilization pretreatment is essential in order to produce carbon fibers from PAN precursor. This oxidative process heats up the raw materials exposed to air for two hours at a temperature of 250°C [22]. Due to the fact that the stabilization process is counted as an intensely exothermic reaction, it will cause fibers to burn if the heating steps are not followed with precision. Initially, the fibers are white and gradually they change from yellow to brown. The final non-flammable product is black. Figure 2.2 shows the chemical reactions during the stabilization process. The 8 presence of —C==N—C==N— functions as a chromophore in the final structure and is the primary reason for the black color of carbon fibers [15]. Figure 2.2 Chemical reactions during the stabilization process. Adopted from [15] 2.1.2.3 Carbonization In the next stage of carbon fiber production, the previously stabilized precursor PANs, which are now non-flammable and melt-resistant, are subjected to two different phases of the heating processes in an inert atmosphere [23]. In this process, which is called carbonization, fibers go through two heating processes at a high temperature, up to 800-3000°C, and obtain up to 95% carbon content [24]. Figure 2.3 shows the chemical reactions that occur during the two stages of the carbonization process. In the first stage of the carbonization, hydrogen atoms are removed from the fibers’ body through the evolution of 𝐻2 gas. In the next step, the higher temperature of the heat treatment removes the nitrogen from the fibers’ structure. 9 Figure 2.3. Chemical reactions during the carbonization process [24] The inert atmosphere in the carbonization process usually contains Nitrogen (𝑁2) or Argon (Ar) [25]. It is common to use 𝑁2 because of its cheaper price. Due to the higher viscosity and density of Argon, the strength of the resulting carbon fiber increases when the inert atmosphere contains this noble gas. Consequently, in some cases, Ar is preferred despite the fact that it is eight times more costly than 𝑁2. Figure 2.4 shows a simple diagram of nitrogen consumption in the carbonization process. 10 Figure 2.4. Diagram of nitrogen consumption in carbonization process [25] Hydrogen chloride (HCl) vapors are also used for the PAN fibers’ carbonization process. The resulting carbon fibers under the HCl vapors treatment will have a higher yield, as a consequence of the reduction in the amount of cyanide (HCN). The two main drawbacks of using HCI in the carbonization process are that it is pricy, and it can also cause corrosion in the furnaces [26]. Normally, two steps are employed for carbonizing the fibers. The first process is in charge of the low-temperature carbonization and the second one is designed for high-temperature carbonization. 2.1.2.3.1 Low temperature (LT) carbonization The low temperature (LT) furnace in the first stage of the carbonization process heats the stabilized PAN fibers up to 950°C [15]. The LT furnaces usually consume electricity for generating heat. The temperature increases gradually by a rate of 5°C/min in order to lessen the mass transfer, which is mainly because of low fiber formation capacity [27]. Moreover, temperatures higher than 1000°C decompose tars, which cause inadequate modifications in the fibers’ properties [15]. The primary purpose of this stage is to remove tars from the fibers. As previously mentioned, the carbonization process is undertaken in an atmosphere full of 𝑁2. The proper amount of 𝑁2 flow 11 removes evolved tars and other gasses. Table 2.1 shows the evolved gases during the carbonization up to a temperature of around 1000°C. Table 2.1. Gasses evolved from different temperatures during carbonization process [15] 2.1.2.3.2 High temperature (HT) carbonization During the LT carbonization, the fibers’ structure becomes more durable, and consequently the chance of damage due to the high temperature and heating rate is not significant. Hence, in the high-temperature carbonization process, the pre-carbonized fibers go through a furnace, which increases the material temperature up to 1600°C in the presence of 𝑁2 gas. Such an atmosphere reduces the amount of Oxygen in the furnace as much as possible. It prevents the PAN fibers from burning during the HT process. Figure 2.5 illustrates the HT furnace in the carbonization process. 12 It contains two zones with different temperatures. The temperature of the first zone is lower than the second one. The process 𝑁2 valve (PN) provides the required amount of nitrogen gas with a density of 1.25 𝐾𝑔𝑚3 [28]. The tension of the fibers can be modified by changing the speeds of the drivers (DR4 and DR5). The high temperature furnace improves the modulus of the fibers. However, temperatures of more than 1500°C and a heating rate of 20°C per minute or more would have a negative effect on the strength of the resulting product [15]. Furthermore, some studies show that if the fibers are cured with temperatures lower than 250°C in the stabilization process, then they will not be able to tolerate temperatures above 1700°C and the product will be brittle and unsatisfactory. Figure 2.5. High temperature furnace in the carbonization process [28] 𝑇0 𝑻𝟎 Stabilization temperature 13 2.1.2.4 The effect of the process parameters on carbon fibers’ mechanical properties The temperature in the carbonization process plays a crucial role in the characterization of the carbon fibers’ mechanical properties and it should be chosen based on the resulting fibers’ application. In other words, the tensile properties of the carbon fibers vastly depend on the production temperature. By increasing the process temperature, a steady growth in the modulus value can be observed [29]. A temperature of 1000°C produces a low modulus outcome, whereas 1500°C can result in intermediate modulus carbon fibers (type II carbon fiber) [24, 30]. Furthermore, in order for carbon fibers to reach a higher modulus, higher temperatures, 1600-1800°C and up to 3000°C, are essential [30, 31]. The latter process, in which the fibers are heated up to 3000°C, is called graphitization. Trinquecoste et al. [32] also suggested that although high tensile strength can be attained by a temperature around 1500°C, more heating treatment is required to achieve a higher modulus. On the other hand, the relation between the tensile strength and process temperature does not remain linear at all temperatures. Figure 2.6 indicates how the carbonization process temperature changes the carbon fibers’ tensile strength. It illustrates the fact that temperatures greater than 1600°C will decrease the tensile strength. New improvements in the production of the carbon fibers results in the creating new fibers with smaller diameters, which helps the carbon fibers to possess a higher tensile strength [33, 34]. 14 Figure 2.6. Typical tensile strength (GPa) of PAN-based carbon fibers over maximum heat treatment temperature [32] An appropriate amount of stretching during the heat treatment steps can improve the modulus and tensile strength of final carbon fibers. It is more effective if it comes with high-temperature exposure [35, 36]. It was observed that the modulus has a direct relation with the amount of applied tension. The effect of stretching on the tensile strength is not significant, although it still increases the fibers’ strength. Furthermore, high heating rates during the carbonization and graphitization can cause a considerable amount of shrinkage in the fibers. Stretching can prevent fibers from shrinkage in length and loss of the favorable orientation [37]. Tsani et al. [36] suggested that optimum stretching can help to obtain the satisfactory mechanical properties in the resulting carbon fibers. 2.1.3 State-of-art on modeling of carbon fiber production process Since the carbon fiber production process is extremely expensive and time-consuming, the preference of the industrial sectors is to minimize the number of the process runs during data collection. Subsequently, the limited number of experimental data, along with severe uncertainty 15 in the mean data, makes the mathematical predictive modeling of the carbon fiber production process an essential step for process design and optimization purposes. Very recently, different mathematical methods (i.e., non-linear multivariable regression, Levenberg-Marquardt neural network algorithm, thin plate splines, and convex hull technique) were implemented in order to optimize the mean value of the resulting carbon fiber modulus [38]. However, the variability of the carbon fiber mechanical properties were not taken in to account. In this thesis, a Gaussian process (GP) modeling approach for the prediction and subsequently optimization of the carbonization process is implemented. The existence of large noises and potential autocorrelation between data points, which can be emanated e.g. due to the continuous/multi-stage time dependent nature of the process, are the main reasons for utilizing the Gaussian process. The next section reviews the basics of the GP modeling. Gaussian Processes Machine learning is a category of artificial intelligence, which is able to develop statistical models in order to interpolate, extrapolate and learn an appropriate pattern without the interaction of humans. Gaussian processes are considered to be an effective tool in machine learning. In fact, many models that are used in the field of machine learning are specific types of Gaussian Processes [10]. Due to the high-priced processes that are involved in manufacturing products in industrial sectors, there is often insufficient (limited) sampling data for statistical analysis. Furthermore, some processes have severe uncertainties, which causes decision-making and optimization to be very challenging for experts. The Gaussian Processes (GPs) as a powerful probabilistic method can 16 cope with the uncertainty. This reliable method provides robust and accurate stochastic prediction models in the presence of noise [39] based on Baye’s rule. 2.2.1 Bayes’ theorem Prior to the explanation of Gaussian Process fundamentals, it is worthwhile to go over the basics of Bayes’ theorem. This well-known theorem is used in the mathematical analysis of Gaussian Processes [40]. Thomas Bayes was the first scientist who suggested the idea that new evidence at a given time point can formally revise the previous beliefs [41]. Let’s assume that 𝑦1, 𝑦2, … , 𝑦𝑛 are 𝑛 random variable which have a joint probability 𝑝(𝑦1, 𝑦2, … , 𝑦𝑛) [10]. If these random variables are divided in to two disjoint sets 𝑦𝐴 and 𝑦𝐵 (𝑝(𝑦) = 𝑝(𝑦𝐴, 𝑦𝐵)), the marginal probability of 𝑦𝐴 will be: 𝑝(𝑦𝐴) = ∫ 𝑝(𝑦𝐴, 𝑦𝐵)𝑑𝑦𝐵. 2.1 In presence of discrete random variables, the integral would be substituted with sum. The conditional probability function is written as: 𝑝(𝑦𝐴|𝑦𝐵) = 𝑝(𝑦𝐴, 𝑦𝐵)𝑝(𝑦𝐵), 2.2 An equation similar to 2.2 is used to define the conditional probability 𝑝(𝑦𝐵|𝑦𝐴). Using the definition of both 𝑝(𝑦𝐴|𝑦𝐵) and 𝑝(𝑦𝐵|𝑦𝐴), the Bayes’ theorem law can be expressed as: 𝑝(𝑦𝐴|𝑦𝐵) = 𝑝(𝑦𝐴)𝑝(𝑦𝐵|𝑦𝐴)𝑝(𝑦𝐵), 2.3 17 The Bayes’ theorem is used for Bayesian inference models which is also the basis of Gaussian Processes. 2.2.2 Gaussian Process basics Let 𝒚 ∈ ℝ be the response value, 𝒙 ∈ ℝ be the input vector, 𝒘 be the weight vector and 𝒇 represents the function value, the following equations will be the Bayesian analysis of a standard linear regression model with Gaussian noise [10]: 𝒇(𝒙) = 𝒙𝑻𝒘, 𝒚 = 𝒇(𝒙) + 𝜺 2.4 It is assumed that the observed response value 𝒚 is not necessarily equal to the value of 𝒇(𝒙). The difference between these two values are due to the presence of a noise, which is assumed to be an identically distributed Gaussian distribution. The mean and variance of this distribution are zero and 𝜎𝑛2 respectively: 𝜺 ~ 𝑵(𝑶, 𝝈𝒏𝟐) 2.5 Next, a Gaussian Process can be defined with the means of covariance and mean functions: 𝑚(𝑥) = 𝐸[𝑓(𝑥)], 2.6 𝑘(𝑥, 𝑥′) = 𝐸[(𝑓(𝑥) − 𝑚(𝑥))(𝑓(𝑥′) − 𝑚(𝑥′))], 2.7 Where 𝑚(𝑥) is the mean function, 𝑘(𝑥, 𝑥′) is the covariance function and 𝑓(𝑥) is a real process. Using mean and covariance functions, the Gaussian process can be generally expressed as 𝑓(𝑥) ~ 𝑁 (𝑚(𝑥), 𝑘(𝑥, 𝑥′)). 2.8 In the Bayesian analysis, a prior over process parameters should be defined. Basically, a prior tells what the initial beliefs are about the characteristics of the model before observations (new evidence) are taken into account. 18 Let 𝑋∗’s be the points at which prediction should be made with the means of Gaussian Process. In addition, 𝑓∗ is the unknown function value at the point 𝑋∗. Therefore, the prior is obtained using the following Gaussian process 𝑓∗ ~ 𝑁 (0, 𝑘(𝑋∗, 𝑋∗)) 2.9 where mean function is assumed to be zero for simplicity and covariance function is 𝑘(𝑋∗, 𝑋∗) [42]. This equation creates the random Gaussian vector based on what is chosen for the covariance function. Figure 2.8.a shows the plot of the created values using the squared exponential as the covariance function (covariance functions will be studied in more detail in section 2.2.3). Clearly, generated lines are completely arbitrary. The reason is that the observations have not been included in the model yet. It only represents the prior and expresses the initial beliefs about the model. The gray area in the plot shows the 95% confidence region. In other words, it links to the standard deviation error bars. The upper bound is obtained by adding two times the standard deviation to the mean value at each point. Subtraction of two times the standard deviation from the mean value gives the lower bound. As previously discussed, a posterior distribution plays the primary role in the Bayesian linear model inference. The posterior is the outcome after the observations are added to the model. Usually, the primary objective is not drawing the prior using some random Gaussian vector. However, plotting the posterior using the observations for developing a proper model is the point of interest. Figure 2.7, illustrates the Gaussian Process function and how the test output, 𝑓∗, is obtained using training points. The blue nodes are the latent variables and are going to be predicted 19 based on the observed input data, which are shown by red color. It is worth mentioning that the prediction 𝑦𝑖∗ are only dependent on their relative latent variables. Figure 2.7 A schematic of Gaussian Process regression [43] Using Bayes’ theorem, it is perceived that: 𝑃𝑜𝑠𝑡𝑒𝑟𝑖𝑜𝑟 =𝐿𝑖𝑘𝑒𝑙𝑖ℎ𝑜𝑜𝑑 × 𝑝𝑟𝑖𝑜𝑟𝑚𝑎𝑟𝑔𝑖𝑛𝑎𝑙 𝑙𝑖𝑘𝑒𝑙𝑖ℎ𝑜𝑜𝑑 𝑝(𝑓, 𝑓∗|𝑦) = 𝑝(𝑓, 𝑓∗) × 𝑝(𝑦|𝑓)𝑝(𝑦) 2.10 Equation 2.10 is utilized in order to reach posterior with the means of prior and observations. 𝑝(𝑓, 𝑓∗∗) is a joint distribution of the training set outputs, 𝑓, and the test set outputs, 𝑓∗, which can be defined as: [𝑓𝑓∗] ~ 𝑁 (0, [𝐾(𝑋, 𝑋) 𝐾(𝑋, 𝑋∗)𝐾(𝑋∗, 𝑋) 𝐾(𝑋∗, 𝑋∗)]) 2.11 20 𝐾(𝑋, 𝑋∗) is a 𝑛 × 𝑛∗ matrix where 𝑛 and 𝑛∗ are the number of observations and test points respectively. This covariance matrix contains all pairs of training and test points. In more realistic cases the observations contain noise. Therefore, the equation 2.11 can be modified to: [𝑦𝑓∗] ~ 𝑁 (0, [𝐾(𝑋, 𝑋) + 𝜎𝑛2𝐼 𝐾(𝑋, 𝑋∗)𝐾(𝑋∗, 𝑋) 𝐾(𝑋∗, 𝑋∗)]) 2.12 where 𝑦 is the observed response value, and 𝜎𝑛2 is the variance of the Gaussian noise 𝜀. For the 𝐾(𝑋, 𝑋) matrix, only observed data points are included, and similarly for 𝐾(𝑋∗, 𝑋∗) test points are evaluated. Figure 2.8.b, indicates three different functions generated from the posterior. In fact, the colored lines are the prior functions that take the training points into account [10]. The black plus signs are the observed data points. From the plots it is obvious that in the areas far from data points the confidence interval become wider (𝑥 = 5). On the other hand, at point like 𝑥 = 0 which is close to observations, the confidence interval becomes very narrow. It indicates that the prediction in this area is much more reliable. Figure 2.8 (a) Prior and (b) posterior examples in the Gaussian process (GP) modelling, assuming data are true means. Lines represent the mean value, plus signs are the inputs, and the gray are is the confidence interval[10] 21 The expected value (mean) of function values 𝑓∗ and its variance are two factors that thoroughly characterize the GP prediction models. Following [10], the mean and the variance can be shown to take the final form of: 𝑓?̅? = 𝐾(𝑋∗, 𝑋)(𝐾(𝑋, 𝑋) + 𝜎𝑛2𝐼)−1𝑦 2.13 𝑉[𝑓∗] = 𝐾(𝑋∗, 𝑋∗) − (𝐾(𝑋, 𝑋) + 𝜎𝑛2𝐼)−1(𝐾(𝑋, 𝑋) + 𝜎𝑛2𝐼)−1 2.14 The application of the mean and variance is exemplified in figure 2.9. The black dotted line in figure 2.9.b shows the mean value at each point, and it was calculated using equation 2.13. The upper and lower bounds of the confidence interval, which is shown by a normal distribution for each point (figure 2.9.a) is calculated using equation 2.14. It is worth noticing that the variance (estimated prediction error by GP) in equation 2.14 does not have any dependency on observations 𝑦, and it does not require to be constant at all data points (as opposed to the standard regression which identifies a single fit error/variance for all points). It is only affected by the covariance function selection and its properties (i.e., model hyper-parameters within 𝑘(𝑋, 𝑋)). Also notice that the predictions 𝑦∗ does not directly depends on function values 𝑓∗ (see also figure 2.7). 22 Figure 2.9. Indication of the mean value and the confidence interval in a Gaussian process (GP) [13] 2.2.3 Covariance function As previously mentioned, Gaussian Processes are described only by their mean and covariance functions [44]. In predictions with the means of Gaussian Process, the covariance function plays a vital role. In fact, the covariance functions govern the properties of the GP prediction model. From an unsupervised learning point of view, the ‘similarity’ is a fundamental concept underlying the prediction of test points [10]. Namely, it implies that for the points whose 𝑥’s are close to each other, it is more probable that their target outputs 𝑦 have similar values. In the Gaussian Process, the covariance function is responsible for specifying this similarity and closeness of the data points. If the mean function is assumed to be zero for simplicity (i.e., for a system with a stationary data nature; i.e., when average of true means remains zero over time/under different manufacturing conditions), the whole trend of the GP prediction model can be defined by the covariance function. For example, the Squared Exponential (SE) covariance function provides smoothness for the 23 process; however, the Ornstein-Uhlenbeck function is suitable for the cases that the observations are highly random with arbitrary characteristic lengths [45]. A periodic covariance function can generate a periodic trend in the prediction model. Figure 2.10 show how the covariance function selection modifies characteristics of the prior in the Gaussian Process, which then naturally affects the final predictions. Figure 2.10 (a) Three prior functions with Squared Exponential covariance function. (b) Three prior functions with Ornstein-Uhlenbeck covariance function [45] Numerous covariance functions are available in the literature which differ from each other mainly in their structure and their free (hyper)parameters. Two general types of covariance functions are the stationary and the dot product types. A covariance function defined based on Euclidian distance between points, (𝑥 − 𝑥′), is called a stationary covariance function. On the other hand, dot product covariance functions are functions defined based on (𝑥. 𝑥′) [10, 45]. Squared exponential, Periodic and Rational quadratic are some examples of stationary covariance functions. As employed widely in the literature and past case studies, here we further discuss the squared exponential (SE) function as given by [10]: 24 𝑘(𝑥, 𝑥′) = 𝜎𝑓2 exp (−|𝑥 − 𝑥′|22ℓ2) + 𝜎𝑛2𝛿 2.15 Where ℓ is a characteristic length scale, 𝜎𝑓2 is the variance of the signal, 𝜎𝑛2 is the noise variance and 𝛿 is Kronecker delta. ℓ, 𝜎𝑓2 and 𝜎𝑛2 are free parameters to be identified via evidence data. Due to the fact that these parameters are related to a non-parametric model, they are called hyper-parameters. These hyper-parameters characterize the trend and behavior of the covariance function and therefore, strongly affect the GP prediction models. Figure 2.11 indicates how the hyper-parameters modify the prediction models. Three values of ℓ are chosen for the Squared exponential covariance function in order to predict a sample of training data. The prediction results for three different characteristic length scales are shown. Figure 2.11.a indicates the GP prediction with optimized value of ℓ. Clearly, the prediction mean line appropriately passes through the observations, and the size of the confidence intervals are fairly reasonable. In figure 2.11.b however, an inadequate value of ℓ generates too wide confidence intervals. Figure 2.11.c also shows that choosing a large value of ℓ causes the prediction model not able to fit the observations. Consequently, in order to develop a good GP prediction model, the hyper-parameters should be optimized precisely. It is worth adding that the Squared exponential covariance function is infinitely differentiable [10]. This special property is the main reason of the exceptional smoothness of the GP prediction models with this covariance function. 25 Figure 2.11. Gaussian process (GP) prediction model using squared exponential covariance function with three different characteristic length scales. Adopted from [10] 2.2.4 Model Selection In order to appropriately employ the Gaussian Process for prediction and optimization tasks, some important decisions needed to be made. First, the type of the covariance and mean functions should be chosen. Then, after choosing the proper covariance and mean functions, hyper-parameters of each function should be optimized in a way that the GP prediction model becomes reliable and robust. 26 Several ways exist for selecting the suitable covariance and mean functions. One way is to choose the functions based on the trends of the observations. For example, if the observation values vary rapidly and randomly, the Ornstein-Uhlenbeck covariance function can be a good choice since this covariance function is able to deal with severe fluctuations (figure 2.10.b). On the other hand, for a smooth training dataset, the Squared exponential covariance function can be used for developing the prediction model. For the cases that the number of the observations are limited, this selection approach is not practical, because the trend of the observations cannot be detected clearly. Consequently, the covariance functions should be chosen based on the results of statistical evaluation measures. The negative marginal likelihood and cross validation methods are two powerful tools that are used [10] in order to select the covariance and mean functions and their corresponding hyper-parameters for the GP prediction model. 2.2.4.1 Negative Marginal Likelihood (NML) Hyper-parameters of the mean and covariance functions shape the behavior of the prediction model [12]. The value of hyper-parameters can be identified by maximizing the marginal likelihood, which is the probability of the input data, given the model [10]. The marginal likelihood equation is written as: log(𝑦|𝑋, 𝜃) = −12𝑦𝑇𝐾𝑦−1𝑦 − 12log|𝐾𝑦| −𝑛2log 2𝜋 , 2.16 where 𝜃 are the hyper-parameters of the covariance function and 𝐾𝑦 is the covariance function of the noisy observations 𝑦 and is given by: 𝐾𝑦 = 𝐾𝑓 + 𝜎𝑛2𝐼, 2.17 where 𝐾𝑓 is the covariance function of the observations without noise. 27 Equation 2.16 has three components. The −12𝑦𝑇𝐾𝑦−1𝑦 component is the only part that includes observations 𝑦 and is called data-fit. The 12log|𝐾𝑦| component is called complexity penalty and it only depends on the covariance function. The 𝑛2log 2𝜋 component is a constant where 𝑛 indicates the number of training points [46]. In the statistical analysis, a desired GP prediction model is the one that gives the lowest value of the negative marginal likelihood (highest marginal likelihood value). 2.2.4.2 Cross validation There are different cross validation methods (exhaustive and non-exhaustive cross validations [47]) that can be used in order to evaluate the accuracy and reliability of the developed prediction models. In practice, a proper cross validation method is chosen based on the sample size and other statistical characteristics of a given dataset [48]. The main goal of cross validation is to detect and prevent overfitting and to understand how well the developed model predicts a new unknown data. The predicted residual error sum of squares (PRESS) statistic is a cross validation tool, which is widely employed to assess predictability of regression models [33]. For cases with small number of experimental data (similar to the proposed case study in Chapter 3), the PRESS method is deemed as an appropriate option due to its reasonable computational cost [33]. The goal in the PRESS method is to determine how well the prediction model fits to a sample of the experimental data that is not used in developing the model. In other words, after the prediction model is developed, each observation is removed from the experimental data sequentially and then the previously developed model does the prediction in the absence of the removed observation. Then the predicted value for the omitted observation is measured. Finally, the PRESS value is calculated 28 by adding the squares of the differences between the values of the respective predictions 𝑌𝑖′ and omitted observations 𝑌𝑖: 𝑃𝑅𝐸𝑆𝑆 = ∑(𝑌𝑖 − 𝑌𝑖′)2𝑛𝑖=1 2.18 2.2.5 Gaussian Process application example This section provides an example of the GP application in an engineering case study and shows how GP accomplishes the task of model development and prediction under noise. The goal of the case study in [46] was to provide a reliable statistical model that predicted the concentration of the carbon dioxide in the atmosphere [10, 49]. Figure 2.12 shows the observations and the resulting prediction of the 𝐶𝑂2 concentration. The blue points are the monthly averages of the carbon dioxide in the atmosphere of Hawaii, from 1958 to 2003, and the gray area is the GP prediction’s confidence interval. In order to model the carbon dioxide concentration against time, the trend of the observations needed to be studied. A long term rising behavior, periodic seasonal variation and a small amount of fluctuation can be observed. Consequently, the covariance function needs to contain all of these characteristics. The covariance function chosen for this case study was a summation (composite) of three basic covariance functions. Namely, SE, periodic and rational quadratic covariance functions. SE was employed to characterize the long smooth rising trend, the periodic covariance function modeled the seasonal changes, and the rational quadratic had a responsibility of modelling the fluctuation. The ensuing composite GP, as shown in figure 2.12, accurately forecasted the carbon dioxide concentration in the atmosphere of the Hawaii. 29 Figure 2.12. The 545 observation data regarding the carbon dioxide concentration in the Hawaii’s atmosphere, along with the 95% confidence interval of the Gaussian process (GP) prediction model shown as hashed region [46] 30 Chapter 3: Model Development Overview This chapter provides the analysis methods and results regarding the GP predictive models that are developed and will be applied in the next chapters to optimize the carbonization process of carbon fiber production. Namely, section 3.2 defines the basis of the proposed case study at the carbon Nexus institute, where the carbonization process data were undertaken from. Next, section 3.3 presents different GP prediction models developed, along with their statistical evaluation results (using Negative Marginal Likelihood/NML and PRES measures). These GP models are then compared in order to select the top ranked model for the prediction of three process outputs: carbon fiber tensile properties including (1) the Young’s modulus and (2) the ultimate tensile strength, as well as (3) the energy consumption during the carbonization. Chapter 3.4 provides the results and the Pareto graphs regarding the carbonization process prediction. The top candidate GP models are those that can predict both the mean values of carbon fibers’ mechanical properties and their standard deviations, as well as the process energy consumption. Furthermore, the statistical measurements (RSS and R-squared) are employed to evaluate how well the GP model predicts the carbonization process conditions. Finally, section 3.5 summarizes the chapter’s main conclusions. Case study introduction As described in chapter 2, the production of carbon fibers is very energy-consuming and expensive due to the high-temperature processes involved. In addition, obtaining the desired carbon fiber mechanical properties is a challenging task since these characteristics do not have any linear correlation [50]. The carbon fiber production processes also encompass many uncertainties and 31 randomness. This noise is mainly due to the fact that fibers consist of large amounts of filaments i.e., 1200 filaments in each tow, and each filament may behave differently. Consequently, in order to successfully deal with all of these special conditions, a thorough analysis is required using accurate and reliable stochastic models. 3.2.1 Focus The focus of the current research is on the carbonization stage of carbon fiber production (see also figures 2.1. and 2.5). The carbonization, as described in chapter 2, is a high temperature process where the pre-stabilized fibers are heated up to 1800°C in a nitrogen rich atmosphere. The objective is to determine how a Gaussian Process framework can predict and optimize the mechanical properties of carbon fibers as well as the process’s energy consumption under different processing parameters and bounds. 3.2.2 Carbon Nexus institute The carbonization processes that were investigated in this thesis, completed at the Carbon Nexus institute at Deakin University in Melbourne, Australia. This institute has extensive composite research facilities that are able to produce commercial grade carbon fibers by means of its low and high temperature carbonization furnaces and other related equipment. The multi-zoned carbonization furnace provides a low temperature between 1100°C and 1175°C, and a high temperature of 1400-1475°C. The schematic figure of the furnace is shown in figure 2.5. These temperatures can be modified by altering the manageable furnaces’ heat elements. The pure nitrogen gas, which exists in the inert atmosphere of the furnaces, is responsible for removing the gasses evolved in the furnace during the carbonization process. The Process 𝑁2 valve (PN) 32 provides enough nitrogen gas to be employed in the furnace. The ratio of nitrogen supplement can vary between 13 and 24 liters per min (𝑙𝑚𝑖𝑛) [28]. The adjustable speed of the drivers (DR4 and DR5) allows applying the required amount of tension to the fibers. The drivers’ speed varies from 25.4 to 28.9 meters per hour (𝑚ℎ). Table A.1 [appendix A] includes some important calculated parameters of the HT furnace. The experimental data of the carbonization process including the furnace temperatures (namely for zone 1 and zone 2 in fig. 2.5), the process time, along with the ensuing tensile properties of the carbon fibers and the process energy consumption (i.e. as output or response variables of the system) are provided in table A.2 and figures A.1-A.4 [appendix A]. 3.2.3 A black-box problem Due to the large amount of uncertainty that exists in carbon fiber production, it is almost impossible to determine the relation between inputs and outputs of the carbonization process explicitly. Therefore, any prediction or optimization investigation falls into the black-box problems category. A black-box is a system in which no information and explicit mathematical understanding of the inner function of the system is available [51]. The only attainable knowledge is the characteristics of the system evidenced through input and output data. One of the primary benefits of the Gaussian Process’s application is that it can appropriately deal with black-box problems. This is because of the fact that GP does not require any supervised modeling of input-output relations. It automatically establishes the relations based on the statistical nature of measured data (e.g., their distances and autocorrelations) in the presence of noise. 33 In order to formally define different aspects of the current optimization problem, a schematic of inputs and outputs for the carbonization process is illustrated in figure 3.1. The inputs consist of some controlled process parameters and some uncertainty or noise. On the other hand, the outputs include the mechanical properties of interest from the resulting carbon fibers as well as the process’s energy consumption. The controlled system parameters are the furnace temperatures of zone 1 and zone 2 and the amount of time that the fibers spend in the furnace. The elastic modulus and tensile strength of fibers and their variabilities, along with the energy consumption magnitude are defined as output. Fibers used in the production line consist of a large number of tows (up to 600) and each tow contains 6000-24000 filaments. The random behavior of filaments and tows, oven heat transfer non-uniformity, and the variation in microstructure of fiber precursors can cause noise in the process response. Consequently, 50 random repeats were undertaken for each process condition in order to capture the variability of the outputs (modulus and tensile strength). Unlike the modulus and tensile strength, the value of energy consumption for each carbonization process does not have any STD. This is expected since for a specific set of carbonization process conditions (i.e. when the inputs -- the two furnace temperatures and time—are fixed) there should be one particular value for energy consumed for the entire system. Carbonization Process Temperature of zone 1°C Temperature of zone 2°CTime (s) NoiseModulus (GPa), Modulus STD (GPa) Tensile Strength (GPa), Tensile strength STD (GPa)Energy Consumption (kWh)Figure 3.1 Schematic of inputs and outputs for the carbonization process under study 34 Gaussian Process developments As discussed in section 2.3.4, the initial step in the model development of the Gaussian Process is to select suitable covariance and mean functions. To this end, the negative marginal likelihood (NML) and PRESS (cross validation method) are utilized to evaluate the quality of candidate GP models. In the next step, as a means to precisely select the best combination of covariance and mean functions, the two-dimensional prediction graphs of the top candidates are sketched. This is practically beneficial for engineering judgement, since the statistical analyses (NML and PRESS) do not take the physical behavior of the process variables into account. In other words, they only examine the model in a mathematical perspective. Therefore, researchers and industrial engineers are responsible for assessing the suggested models qualitatively and determining if they mimic the nature of process under question. For instance, a statistical analysis may suggest that the periodic function for GP is the best option for the carbon fiber’s modulus prediction. However, it is clear that the modulus does not have any periodic trend based on Appendix A data. Therefore, the periodic function should be ignored although mathematically it is among the top candidates for the modulus model selection. In the following sections, the procedure used for evaluating and choosing the best covariance and mean functions will be explained. First, the modulus mean and its standard deviation prediction will be investigated. Then, the applied approach for the prediction of the tensile strength mean and standard deviation will be presented. The prediction of the energy consumption during the carbonization process is the last part in this chapter. All the statistical analyses, including calculating negative marginal likelihood, GP perdition, PRESS evaluation and single and multi-35 objective optimizations, were conducted in MATLAB®. The GPML toolbox was used for the implementation of the Gaussian Process models [52]. 3.3.1 Modulus mean In order to obtain an initial understanding of what covariance and mean functions are suitable for the prediction of carbon fiber modulus mean values, the NML of some probable functions are calculated. It is worth mentioning that one unique feature of the GP approach is that new covariance and mean functions can be generated by combining known covariance and mean functions [53]. Tables 3.1 and 3.2 indicate the NML values for selected combinations of covariance and mean functions. Table 3.2 indicates that the Squared Exponential and periodic covariance functions give a better outcome compared to the other functions. Moreover, the Linear and Constant functions are the top choices for the mean function. Table 3.1. Negative marginal likelihood (NML) values of different covariance and mean functions used for modulus mean prediction Covariance function Mean function Negative marginal likelihood Squared Exponential (SE) Linear 66.05 Squared Exponential Constant 67.30 Squared Exponential Zero 76.53 Periodic Linear 64.37 Periodic Constant 66.49 Periodic Zero 70.85 Polynomial Linear 68.66 Polynomial Constant 70.32 Polynomial Zero 69.34 Linear Linear 69.55 Linear Constant 73.61 Linear Zero 84.07 36 Table 3.2. Negative marginal likelihood (NML) values of combined (composite) covariance and mean functions used for modulus mean prediction Covariance function Mean function Negative marginal likelihood SE+Perodic Linear 63.688 SE+Perodic Constant 66.84 SE+Linear Linear 69.73 SE+Linear Constant 67.30 SE Linear+Constant 67.10 Periodic Linear+Constant 67.57 Linear Linear+Constant 71.18 In the next step, the top candidates are evaluated with the means of the cross-validation method (PRESS). Figure 3.2 indicates the PRESS values for the selected functions. It shows that among all of these function combinations, the Squared exponential covariance function and the linear mean function give the best outcome. The periodic covariance function and the constant mean function ranked the second best model. Moreover, combined covariance and mean functions do not seem to be superior choices in the current case study. They only make the model more complex and more expensive to analyse while not performing as well as single functions. In order to select the final model reliabley, however, visual inspection of predictions is also required. Two-dimensional graphs of the modulus mean prediction are plotted. Figure 3.3 shows the plots of the modulus versus time at three different furnace temperatures. The black line in the graphs is the mean value of the prediction, the black stars are the observation points and the gray area is the 95% confidence interval of the predictions. Figure 3.4 also depicts the modulus mean prediction using the periodic covariance function and the constant mean function (second best model in cross-validation). It shows that the error bars are wider than the previous GP model and the prediction mean line and its confidence intervals have a periodic trend that is not desirable due to the fact that no periodic behavior would be physically meaningful for the carbon fiber modulus. 37 To summarize, it was shown that the GP prediction model using the Squared exponential covariance function and the Linear mean function gave the best outcome in modulus mean prediction. The cross-validation results of this model is the best (smallest PRESS) among other models per fig. 3.2. Moreover, in two-dimensional (qualitative assessment) graphs, the prediction of mean using this model showed a reasonable trend and appropriately passed through the input data points. Therefore, the final choice for covariance and mean functions of modulus mean prediction are the Squared Exponential and the Linear, respectively, and will be used in subsequent chapters. Figure 3.2. PRESS values for different modulus mean Gaussian processes (GP) 6.33E+038.12E+03 8.16E+03 8.24E+038.72E+031.01E+040.00E+002.00E+034.00E+036.00E+038.00E+031.00E+041.20E+04SE Covariancefunction andLinear meanfunctionPeriodicCovariancefunction andConstant meanfunctionSE+PeriodicCovariancefunction andLinear meanfunctionSE Covariancefunction andConstant meanfunctionPeriodicCovariancefunction andLinear meanfunctionSE+PeriodicCovariancefunction andConstant meanfunctionPRESS 38 3.3.2 Modulus standard deviation The measured tensile properties of fibers (training data in Table A.2 and/or fig A.1-4) of the carbonization process include both the mean and standard deviation (STD). The STDs under each process condition are unequal and have been clearly influenced by different heat treatment times and furnace temperatures. Statistically, a state like this, in which the variability of random variables are not equal, is called heteroscedasticity [54, 55]. However, in the standard GP, the model assumes that standard deviations of actual data are alike (homoscedasticity). Hence, for heteroscedastic data, two separate GP models should be developed. The first model tackles the mean prediction as explained in section 3.3.1. In the second GP model, the inputs are the same as the first GP model, however, the output response is now the standard deviation. In this section, the suitable covariance and mean functions for the prediction of the modulus STD will be investigated. Figure 3.5 indicates the values of the negative marginal likelihoods (NML) for the different covariance and mean functions. It is clear that the linear covariance function is not a good choice for modulus STD prediction because it gives a large value of NML. The SE and Periodic covariance functions, however, perform reasonably well and have a better NML result than the linear covariance function. Moreover, whenever the constant mean function is chosen, NML values become larger than the cases in which the mean function is linear. Figure 3.6 shows the summary of the cross-validation results of the GP models for modulus STD, which are developed using top covariance and mean functions. PRESS also verifies that the SE covariance function and the linear mean function can be an appropriate choice (second lowest) for the prediction of modulus STD. Although the GP model with the periodic covariance function and 39 the constant mean function has the lowest (best) PRESS value, it is not chosen for modulus STD prediction due to the fact that no apparent periodic behavior can be qualitatively observed in the observations trend. Figure 3.7 illustrates the two-dimensional plots of modulus STD prediction. The mean line of the prediction has a very smooth trend, and satisfyingly fits the experimental data. 3.3.3 Tensile strength mean In the prediction of the carbon fiber modulus, mean and standard deviation were modeled separately by two GPs due to heteroscedasticity. The same procedure is applied here to find best GP models to predict the tensile strength response. The values of the negative marginal likelihood from several GP models using different covariance and mean functions are shown in figure 3.8. Among three covariance functions, polynomial gives the highest and SE gives the lowest NML values. Moreover, the constant mean function along with the SE covariance function creates the best combination for obtaining the lowest NML. Figure 3.9 also demonstrates the cross validation results for the tensile strength mean. It shows that the SE covariance function and the linear mean function give the most reliable outcome. Both NML and PRESS results indicate that the SE would be a good choice for the covariance function. However, choosing the suitable mean function is more complex in this case, because NML and PRESS suggest different mean functions. The constant mean function provides the best outcome in NML and the Linear mean function gives the superior result in PRESS. Hence, the ultimate step for choosing the best mean function would be to visually observe the two-dimensional prediction plots. Figures 3.10 and 3.11 indicate the prediction plots of the tensile strength mean with different 40 mean functions. SE is the covariance function for both figures as it has been chosen as the best function for the GP prediction model. From these figures, it is clear that the model with the linear mean function fits the data more appropriately than the model with the constant mean function. This can be seen in figure 3.11.c where the predictions at input data are assocaited with some confidence interval (i.e. prediction uncertainity), whereas the mean line in figure 3.10.c successfully passes through all the experimental data with almost zero confidence interval. Consequently, the final GP model for the prediction of the tensile strength mean utilizes SE covariance function and linear mean function. As mentioned in section 2.3.3, the covariance function is the factor that describes the main trend of the prediction model and mean function usually has a minor effect. The predictions indicated in figures 3.10 and 3.11 have the same general trend due to the fact that the SE covariance function is selected for both models; however, their confidence intervals are clearly impacted by the choice of mean function. 3.3.4 Tensile strength standard deviation (STD) The summary of NML and PRESS of different GP models for the prediction of the tensile strength STD are illustrated in figures 3.12 and 3.13. The GP model using Linear mean function and SE covariance function comes with the lowest NML and PRESS. Therefore, this model is chosen for the prediction of tensile strength standard deviation. Figures 3.14 shows the prediction plots of tensile strength STD at three different furnace temperature conditions. 41 3.3.5 Energy consumption As explained in section 3.2.3, it is not required to develop two separate GP models for the energy prediction, as the latter response has no STD. Figure 3.15 and 3.16 summarize the NML and PRESS of different types of GP models for the prediction of energy consumption. The observations show that the carbonization process energy consumption values are very close to each other and have a constant trend. Moreover, the NML and PRESS results indicate that the Constant mean function provide much better result compared to other mean functions. Therefore, the SE covariance function and Constant mean function are chosen for the energy consumption prediction in the carbonization process. Figure 3.17 indicates the two-dimensional graphs of the energy consumption GP prediction at different furnace temperature conditions. 42 (a) (b) (c) Figure 3.3. Modulus mean prediction as a function of time using squared exponential (SE) covariance function and Linear mean function at three different furnace temperature conditions. (a) when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1175°C. It is worth mentioning that in the areas where more observed data exists, the confidence interval become narrower, representing a higher level of prediction reliability around those points. 43 (a) (b) (c) Figure 3.4. Modulus mean prediction as a function of time using periodic covariance function and Constant mean function at three different furnace temperature conditions. (a) when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1475°C. 44 Figure 3.5. Negative marginal likelihood values of different Gaussian process (GP) models for modulus standard deviation (STD) prediction 58.6441539 58.74195244 60.0527765671.317263583.64409514103.6652999020406080100120SE covariancefunction andlinear meanfunctionPeriodiccovariancefunction andconstant meanfunctionPeriodiccovariancefunction andlinear meanfunctionSE covariancefunction andconstant meanfunctionLinearcovariancefunction andconstant meanfunctionLinearcovariancefunction andLinear meanfunctionNML45 Figure 3.6. PRESS values of different Gaussian process (GP) models for modulus standard deviation (STD) prediction 1537.531626.831784.932022.752323.272636.880.00500.001000.001500.002000.002500.003000.00Periodiccovariancefunction andconstant meanfunctionSE covariancefunction andlinear meanfunctionPeriodiccovariancefunction andlinear meanfunctionSE covariancefunction andconstant meanfunctionLinearcovariancefunction andlinear meanfunctionLinearcovariancefunction andconstant meanfunctionPRESS46 (a) (b) (c) Figure 3.7. Modulus standard deviation (STD) prediction versus time using squared exponential (SE) covariance function and Linear mean function at three different furnace temperature conditions. (a) happens when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) is the modulus prediction at 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) shows the result at 𝑻𝟏=1175°C and 𝑻𝟐=1475°C. (GPa) (GPa) (GPa) 47 Figure 3.8. Negative marginal likelihood values of different Gaussian process GP models for tensile strength mean prediction 2.242.6752.793.0123.26 3.31100.511.522.533.5SE covariancefunction andconstant meanfunctionSE covariancefunction andlinear meanfunctionPeriodiccovariancefunction andconstant meanfunctionPeriodiccovariancefunction andlinear meanfunctionPolynomialcovariancefunction andlinear meanfunctionPolynomialcovariancefunction andconstant meanfunctionNML48 Figure 3.9. PRESS values of different Gaussian process (GP) models for tensile strength mean prediction 1.21.39 1.482.44.600.511.522.533.544.55SE covariancefunction and linearmean functionSE covariancefunction andconstant meanfunctionPeriodic covariancefunction andconstant meanfunctionPolynomialcovariance functionand linear meanfunctionPolynomialcovariance functionand constant meanfunctionPRESS49 (a) (b) (c) Figure 3.10. Tensile strength mean prediction versus time using squared exponential (SE) covariance function and Linear mean function at three different furnace temperature conditions. (a) happens when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1475°C. 50 (a) (b) (c) Figure 3.11. Tensile strength mean prediction versus time using squared exponential (SE) covariance function and constant mean function at three different furnace temperature conditions. (a) when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1475°C. 51 Figure 3.12. Negative marginal likelihood values of Gaussian process (GP) models for tensile strength standard deviation (STD) prediction 0.7090.7130.7170.7230.7590.7720.6700.6800.6900.7000.7100.7200.7300.7400.7500.7600.7700.780SE covariancefunction andconstant meanfunctionSE covariancefunction andlinear meanfunctionPolynomialcovariancefunction andlinear meanfunctionPeriodiccovariancefunction andconstant meanfunctionPeriodiccovariancefunction andlinear meanfunctionPolynomialcovariancefunction andconstant meanfunctionNML52 Figure 3.13. PRESS values of different Gaussian process (GP) models for tensile strength standard deviation (STD) prediction 0.60.650.680.7 0.70.720.540.560.580.60.620.640.660.680.70.720.74SE covariancefunction andlinear meanfunctionPeriodiccovariancefunction andLinear meanfunctionPeriodiccovariancefunction andconstant meanfunctionSE covariancefunction andconstant meanfunctionPolynomialcovariancefunction andlinear meanfunctionPolynomialcovariancefunction andconstant meanfunctionPRESS53 (a) (b) (c) Figure 3.14. Tensile strength standard deviation (STD) prediction versus time using squared exponential (SE) covariance function and Linear mean function at three different furnace temperature conditions. (a) when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1475°C. 54 Figure 3.15. Negative marginal likelihood values of Gaussian process (GP) models for energy consumption prediction 2.918.35 9.281523370510152025303540SE covariancefunction andconstant meanfunctionSE covariancefunction andlinear meanfunctionPeriodiccovariancefunction andconstant meanfunctionLinearcovariancefunction andconstant meanfunctionLinearcovariancefunction andlinear meanfunctionPeriodiccovariancefunction andlinear meanfunctionNML55 Figure 3.16. PRESS values of different Gaussian process (GP) models for energy consumption prediction 0.011 0.129 0.2512.655.510.066024681012Periodiccovariancefunction andconstant meanfunctionSE covariancefunction andconstant meanfunctionSE covariancefunction andlinear meanfunctionLinearcovariancefunction andconstant meanfunctionLinearcovariancefunction andlinear meanfunctionPeriodiccovariancefunction andlinear meanfunctionPRESS56 (a) (b) (c) Figure 3.17. Energy consumption prediction versus time using squared exponential (SE) covariance function and constant mean function at three different furnace temperature conditions. (a) when 𝑻𝟏=1100°C and 𝑻𝟐=1400°C. (b) when 𝑻𝟏=1150°C and 𝑻𝟐=1450°C and (c) when 𝑻𝟏=1175°C and 𝑻𝟐=1475°C57 Verification of the selected Gaussian process (GP) models In order to closely verify the performance of the developed GP models above, in this section the prediction values and experimental data (details in Appendix A) are compared for each model. 3.4.1 Modulus mean As explained in section 3.3.1, the SE covariance and Linear mean functions were selected for the modulus mean GP prediction. Figure 3.18 represents the experimental data versus the prediction values of this model. The closeness of the points to the 45-degree line shows how well the prediction model has performed. Table 3.3 also provides more detailed information about the goodness of fit of the modulus mean prediction model. 3.4.2 Modulus standard deviation (STD) The summary of the prediction results of the developed SE covariance-Linear mean GP model for modulus STD is provided in figure 3.19. The closeness to ideal prediction line proves how successful the GP has been in the modeling of modulus STD. Table 3.4 also provides detailed information and relative error values. 58 Figure 3.18. Experimental data against prediction values for the modulus mean Gaussian process (GP) model Table 3.3. Summary of the modulus mean prediction using the Gaussian process (GP) model Temperature of zone 1 (°C) Temperature of zone 2 (°C) Time (s) Experimental data (GPa) Prediction values (GPa) Relative error 1100 1400 156 228.8 228.813 0.00006 1100 1400 167 227.7 227.698 0.00001 1100 1400 162 227.5 227.463 0.00016 1100 1400 151 234.6 234.494 0.00045 1125 1425 156 228.3 228.299 0.00001 1125 1425 167 229.2 229.183 0.00007 1125 1425 162 226.7 226.677 0.00010 1125 1425 151 231.9 231.807 0.00040 1150 1450 156 219.8 219.807 0.00003 1150 1450 167 221.8 221.795 0.00002 1150 1450 162 217.6 217.592 0.00004 1150 1450 151 224.3 224.223 0.00034 1175 1475 156 197.7 197.510 0.00096 1175 1475 167 154.2 154.316 0.00076 1175 1475 162 157.5 157.686 0.00118 1175 1475 151 163.02 163.241 0.00136 1100, 1400, 1561100, 1400, 1671100, 1400, 1621100, 1400, 1511125, 1425, 1561125, 1425, 1671125, 1425, 1621125, 1425, 1511150, 1450, 1561150, 1450, 1671150, 1450, 1621150, 1450, 1511175, 1475, 1561175, 1475, 1671175, 1475, 1621175, 1475, 1510501001502002500 50 100 150 200 250Prediction values(GPa)Experimental data (GPa) Temperature 1(°C ), Temperature 2(°C ), time(s)59 Figure 3.19. Experimental data against prediction values for the modulus standard deviation (STD) Gaussian process (GP) model. Table 3.4. Summary of the modulus standard deviation (STD) predictions using Gaussian process (GP) model Temperature of zone 1 (°C) Temperature of zone 2 (°C) Time(s) Experimental data (GPa) Prediction values (GPa) Relative error 1100 1400 156 18.47 18.43 0.002 1100 1400 167 9.288 9.48 0.020 1100 1400 162 9.89 10.02 0.013 1100 1400 151 13.67 13.64 0.002 1125 1425 156 14.48 14.47 0.001 1125 1425 167 34.03 33.90 0.004 1125 1425 162 11.66 11.75 0.008 1125 1425 151 15.9 15.82 0.005 1150 1450 156 7.797 7.85 0.007 1150 1450 167 10.11 10.25 0.014 1150 1450 162 32.99 32.80 0.006 1150 1450 151 10.78 10.74 0.003 1175 1475 156 35.32 35.02 0.008 1175 1475 167 9.5 9.63 0.013 1175 1475 162 9.13 9.21 0.009 1175 1475 151 5.39 5.40 0.002 1100, 1400, 1561100, 1400, 1671100, 1400, 1621100, 1400, 1511125, 1425, 1561125, 1425, 1671125, 1425, 1621125, 1425, 1511150, 1450, 1561150, 1450, 1671150, 1450, 1621150, 1450, 1511175, 1475, 1561175, 1475, 1671175, 1475, 1621175, 1475, 15105101520253035400 5 10 15 20 25 30 35 40Predicted modulus STDExperimental data Temperature1(°C), Temperature2(°C), time(s)60 3.4.3 Tensile strength mean From section 3.3.3, it was observed that the best GP model for prediction of tensile strength mean also uses a SE covariance function and a linear mean function. The outcome of the predictions using this GP is summarized in figure 3.20 and table 3.5. Figure 3.20. Experimental data against prediction values of tensile strength mean Gaussian process (GP) model. 1100, 1400, 156 1100, 1400, 1671100, 1400, 1621100, 1400, 1511125, 1425, 1561125, 1425, 1671125, 1425, 162 1125, 1425, 1511150, 1450, 1561150, 1450, 1671150, 1450, 1621150, 1450, 1511175, 1475, 1561175, 1475, 1671175, 1475, 1621175, 1475, 15100.511.522.533.540 0.5 1 1.5 2 2.5 3 3.5 4Predicted tensile strength mean (GPa)Experimental data (GPa)Temperature1(°C), Temperature2(°C), time(s)61 Table 3.5. Summary of the tensile strength mean prediction using Gaussian process (GP) model Temperature of zone 1 (°C) Temperature of zone 2 (°C) Time(s) Experimental data (GPa) Prediction values (GPa) Relative error 1100 1400 156 3.016 3.016 0.00012 1100 1400 167 3.201 3.201 0.00003 1100 1400 162 3.275 3.274 0.00016 1100 1400 151 3.458 3.457 0.00029 1125 1425 156 2.925 2.925 0.00008 1125 1425 167 3.185 3.185 0.00007 1125 1425 162 3.031 3.031 0.00002 1125 1425 151 3.225 3.224 0.00019 1150 1450 156 2.945 2.945 0.00008 1150 1450 167 3.13 3.130 0.00002 1150 1450 162 3.11 3.110 0.00007 1150 1450 151 3.18 3.179 0.00016 1175 1475 156 2.78 2.780 0.00014 1175 1475 167 2.11 2.111 0.00050 1175 1475 162 2.35 2.351 0.00023 1175 1475 151 2.209 2.210 0.00050 3.4.4 Tensile strength standard deviation (STD) The prediction summary of tensile strength STD is shown in figure 3.21 and table 3.6, again proving GP’s accuracy and reliability as a predictive modelling framework for the current case study. The results, similar to the earlier GPs, are based on the SE covariance function and the linear mean function. 62 Figure 3.21. Experimental data versus prediction values of tensile strength standard deviation (STD) Gaussian process (GP) model Table 3.6. Summary of the tensile strength standard deviation (STD) prediction using Gaussian process (GP) model Temperature of zone 1(°C) Temperature of zone 2 (°C) Time(s) Experimental data (GPa) Prediction values (GPa) Relative error 1100 1400 156 1.06 1.059 0.00090 1100 1400 167 0.68 0.680 0.00047 1100 1400 162 0.77 0.770 0.00003 1100 1400 151 0.76 0.760 0.00017 1125 1425 156 1 0.999 0.00093 1125 1425 167 0.81 0.810 0.00007 1125 1425 162 0.81 0.810 0.00013 1125 1425 151 0.59 0.590 0.00043 1150 1450 156 0.63 0.630 0.00006 1150 1450 167 0.87 0.870 0.00039 1150 1450 162 0.65 0.650 0.00033 1150 1450 151 0.63 0.630 0.00017 1175 1475 156 0.68 0.680 0.00052 1175 1475 167 0.36 0.361 0.00242 1175 1475 162 0.41 0.411 0.00141 1175 1475 151 0.34 0.341 0.00202 1100, 1400, 1561100, 1400, 1671100, 1400, 1621100, 1400, 1511125, 1425, 1561125, 1425, 1671125, 1425, 1621125, 1425, 1511150, 1450, 1561150, 1450, 1671150, 1450, 1621150, 1450, 1511175, 1475, 1561175, 1475, 1671175, 1475, 1621175, 1475, 15100.20.40.60.811.20 0.2 0.4 0.6 0.8 1 1.2PREDICTEDTENSILESTRENGTHSTDExperimental data Temperature1(°C) , Temperature2(°C), time(s)63 3.4.5 Energy consumption Section 3.3.5 showed that the SE covariance function with a constant mean function provides the best GP model for the prediction of carbonization process energy consumption; which was expected due to the fact that the energy consumption experimental values showed a steady behavior over the time variable (see fig. 3.17). Figure 3.22 and table 3.7 illustrate the results of the energy consumption predictions using the corresponding GP model. Figure 3.22. Energy consumption prediction values against the experimental data. 1100, 1400, 1561100, 1400, 1671100, 1400, 1621100, 1400, 1511125, 1425, 1561125, 1425, 1671125, 1425, 1621125, 1425, 1511150, 1450, 1561150, 1450, 1671150, 1450, 1621150, 1450, 1511175, 1475, 1561175, 1475, 1671175, 1475, 1621175, 1475, 15111.81212.212.412.612.81311.8 12 12.2 12.4 12.6 12.8 13Predicted energy consumption (kwh)Experimental data (kwh)Temperature1(°C), Temperature2(°C), time(s)64 Table 3.7. Summary of the carbonization process energy consumption prediction using Gaussian process (GP) model Temperature of zone 1 (°C) Temperature of zone 2 (°C) Time(s) Experimental data (GPa) Prediction values (GPa) Relative error 1100 1400 156 11.87 11.886 0.0014 1100 1400 167 11.95 11.955 0.0004 1100 1400 162 11.92 11.919 0.0001 1100 1400 151 11.86 11.869 0.0007 1125 1425 156 12.19 12.319 0.0106 1125 1425 167 12.44 12.397 0.0035 1125 1425 162 12.4 12.360 0.0032 1125 1425 151 12.37 12.289 0.0065 1150 1450 156 12.49 12.497 0.0006 1150 1450 167 12.58 12.577 0.0002 1150 1450 162 12.51 12.541 0.0025 1150 1450 151 12.47 12.462 0.0007 1175 1475 156 12.7 12.706 0.0005 1175 1475 167 12.81 12.788 0.0018 1175 1475 162 12.77 12.753 0.0013 1175 1475 151 12.65 12.663 0.0010 3.4.6 Statistical evaluation of the Gaussian process (GP) models The overall performances of selected GP models for predicting the carbonization process parameters are assessed in this section, using statistical measures of the residual sum of square (RSS) and R-squared (or coefficient of multiple determination). The RSS determines how close the prediction and observations are [56]. A smaller value for RSS and a higher R-squared indicate a closer fit between the model and data points [57]. Table 3.8 summarizes the RSS and R-squared values of the GP prediction models, clearly showing a notable performance of each model predictability. 65 Table 3.8. Statistical measures’ summary of the carbonization process’s Gaussian process (GP) prediction models GP prediction model RSS (GPa) R-squared Modulus mean 0.161 0.99 Modulus STD 0.258 0.99 Tensile strength mean 0.00005 1.00 Tensile strength STD 0.0000038 1.00 Energy consumption 0.028 0.98 Summary of findings This chapter investigated several covariance and mean functions in order to find the best GP models towards the prediction of the carbonization process response variables. The SE covariance function is chosen for the GP models due to its strong smooth behavior under all responses. It was shown that regardless of how the mean function is chosen, the GP model with the SE covariance function was among top ranked candidates, as judged by the NML and PRESS values. Moreover, it was discussed that although the periodic covariance function results in a reliable NML and PRESS, it is not a physically meaningful choice for the carbonization process prediction, because there is no physical evidence of a periodic trend in the variation of responses with process control variables. Furthermore, this chapter introduced the fibers’ Young’s modulus, tensile strength, their STDs, and the process’ energy consumption as five important decision variables for a multi-objective optimization of the carbonization process. Comparison of predicted and measured values indicates that the selected GP models have accomplished the prediction task extremely well. In particular, the RSS values smaller than 0.25 and R-squared values all above 99% confirmed that the prediction models closely fit the experimental data, with only three hyper-parameters. It is interesting that the GP models were particularly able to predict quality variables such as modulus STD and tensile strength STD despite their very fluctuating (noisy) nature. 66 Chapter 4: Multi-objective Optimization of Carbonization Process Overview This chapter provides the analysis and results regarding the multi-objective optimization of the carbonization process. Sections 4.2 and 4.3 present the choice of the optimization method and determination of the process constraints and bounds. Four different process design scenarios will be then pursued in section 4.4 and for each scenario, the multi-objective optimization solution will be discussed. Namely, section 4.4.1 will investigate a ε-constraint multi-objective optimization task with the primary goal of maximizing the carbon fibers’ modulus mean. Sections 4.4.2 and 4.4.3 will study two multi-objective optimization scenarios with different constraints, with a focus on discovering the optimal solutions for both the modulus and tensile strength means using the weighted sum method. Section 4.4.4 will investigate a solution for ε-constraint minimum of the energy consumption of the carbonization process. Section 4.5 will summarize the results and findings of the performed carbonization process optimization scenarios. Multi-objective optimization formulation Every practical optimization problem falls into one of the two general categories, single- or multi-objective. In single-objective optimization, the goal is to find one solution that usually minimizes or maximizes the objective in the presence of some constraints/process bounds. In a multi-objective optimization, however, there are more than one objective to be satisfed and they conflict one another. The latter causes that there would be more than one optimal solution (each solution would favor one objective more or less than the others while overall being a best compromise). A common practical strategy in solving multi-objective problems is to create the so-called Pareto solutions (tradeoff curves) using multi-criteria decision making [56], and accordingly, the decision 67 maker/designer can select the final optimal solution based on his/her preference, additional practical constraints and bounds, etc. 4.2.1 Constraints and bounds Defining the decision variables, objectives, constraints, and bounds is a fundamental step in any optimization procedure [58]. After completion of model definition, an algorithm for solving optimization is needed for finding the Pareto optimal solutions. The objective of this section is to define the requirements and formulation of the carbonization process optimization. As described in section 3.3, the Gaussian Processes were used in order to develop black-box models for the prediction of mechanical properties (elastic modulus and the tensile strength) of resulting carbon fibers in the carbonization process, the variabilities of these properties (standard deviations), and the energy consumption. The five developed GP models (sections 3.3.1 to 3.3.5) will be used as part of the optimization tasks. Based on the latter five objectives (system responses), identifying the constraints are highly dependent on the method that is to be chosen for the multi-objective optimization. For example, if the weighted sum method [57] is employed to aggregate the means, the standard deviations of modulus and tensile strength may be considered as the constraints, similar to [59, 60]. In another method, called the ε-constraints method, one of the primary decision variables (objectives) is selected to be optimized and the other variables are all converted to constraints [61, 62]. The ε-constraints method can be mathematically written as: 68 minimize 𝑓𝑙(𝑥) subject to 𝑓𝑗(𝑥) < 𝜀𝑗 for all 𝑗 = 1, … , 𝑘, 𝑗 ≠ 𝑙 4.1 Where 𝑙 ∈ {1, … , 𝑘} and 𝜀𝑗 is the upper bound for the converted objectives (𝑗 ≠ 𝑙). Each 𝑓 (decision variable) can be e.g. one of the five system responses identified in the carbonization process. Note that in hybrid modeling cases, the main objective 𝑓𝑙 can itself be e.g. a weighted sum of two of the responses. The bounds for the optimization should consist of the limitations that exist directly in the design of input space of the problem (here, the range of the low-temperature and high-temperature furnace values and the range of heat treatment time based on raw data in Appendix A): • Zone 1’s temperature interval is from 1100(°C) to 1175(°C): 1100 < 𝑇1(°C) < 1175 • Zone 2’s temperature interval is from 1400(°C) to 1475(°C): 1400 < 𝑇2(°C) < 1475 • The amount of time which fibers spend in the furnaces is between 151(s) and 167 (s): 151 < 𝑡(𝑠) < 167. 4.2.2 Underlying generic trade-offs in the carbonization process optimization problem As explained in section 4.2, the multi-objective optimization is essentially utilized whenever the decision variables conflict with each other and there is no distinct correlation between them. One way to visually evaluate these trade-offs is to create Pareto graphs of decision variables. Figure 4.1 depicts the generic (unconstrained) Pareto graphs of the modulus, tensile strength, and energy consumption in the current study, whose values were directly derived from the experimental data of the carbonization process (Table A.2). No linear correlation can be seen between any two of the 69 decision variables, and interestingly, not higher consumed energy necessarily means achieving higher fiber modules to tensile strength. Overall, it can be observed in Figure 4.1 that there are no apparent conflicts among the three criteria trends (i.e., high modulus and high tensile strength are generally achievable at low energy consumption levels). In particular, the red point in the graphs represents an experimental data point that has had the highest modulus and tensile strength along with a low energy consumption at the same time. The above observations then bring up the question whether a multi-objective optimization is really necessary for the carbonization process. To answer, Figure 4.2 shows the Pareto fronts of modulus STD, tensile strength STD and energy consumption. The figures now show that the red point, which successfully had satisfied the mean values of mechanical properties, does not have outstanding STDs. In addition, the figure shows that lower energy consumption is associated with higher STD in fibers mechanical properties, which is a clear practical conflict for a decision maker. Consequently, if all the decision factors are to be taken into account in the carbon fiber manufacturing quality improvement (including both average performance and robustness of fibers), the multi-objective optimization will be an essential step for the designers. 70 (a) (b) (c) Figure 4.1. Experimental Pareto graphs of (a) modulus mean versus tensile strength mean, (b) modulus mean versus energy consumption, and (c) tensile strength mean versus energy consumption. Red points represent an experimental data with the highest mean modulus and mean tensile strength along with a low energy consumption in the carbonization process. 00.511.522.533.540 50 100 150 200 250Tensilestrength (GPa)Modulus (GPa)11.81212.212.412.612.8130 50 100 150 200 250Energy Consumption (kWh)Modulus (GPa)11.81212.212.412.612.8130 1 2 3 4Energy Consumption (kWh)Tensile strength (GPa)71 (a) (b) (c) Figure 4.2. Experimental Pareto graphs of (a) modulus STD versus tensile strength STD, (b) modulus STD versus energy consumption, and (c) tensile strength STD versus energy consumption. Red points represent an experimental data with the highest mean modulus and mean tensile strength along with a low energy consumption in the carbonization process. 00.20.40.60.811.20 10 20 30 40Tensile Strength STD (GPa)Modulus STD (GPa)11.81212.212.412.612.8130 10 20 30 40Energy (kWh)Modulus STD (GPa)11.81212.212.412.612.8130 0.5 1 1.5Energy (kWh)Tensile Strength STD (GPa)72 Customer-driven optimization The manufacturing-based optimization problems should be specified based on their practical objectives/applications. The objectives are usually linked to the customer needs. For example, a customer may require carbon fibers with a small variability in their modulus values so that it can meet design factors. Due to carbon fiber’s specific utilization in a strucutre, the same customer may also want the mean modulus values to be as high as possible. In this situation, the optimization of modulus mean can be pursued by setting the modulus STD as a constraint and applying its upper and lower bounds based on the customer’s need. Selected optimization scenarios In this section, four different optimization scenarios have been developed and solved. These scenarios may be viewed as potential case studies that can occur for carbon fiber industrial sectors based on customers’ need. The weighted sum method and the ε-constraints method are the two methods used in solving the multi-objective scenarios. It is of note that a major drawback of the weighted sum method is that for nonconvex cases and convex problems with non-constant curvature, some of the Pareto optimal solutions can not be captured [63]. However, for the latter cases, an adaptive weight sum method may be implemented, as developed in [64] and will be exemplified under scenario 2 here. 4.4.1 Scenario 1: Modulus optimization Due to high-temperature operation demand, the carbonization process is a very energy-consuming step in carbon fiber manufacturing. Consequently, the manufacturers are looking for a way to obtain suitable carbon fibers’ mechanical properties along with the lowest energy consumption. In 73 the first scenario the modulus of the carbon fiber is maximized in the presence of constraints and bounds. The ε-constraints method is used for this optimization. In this case, the energy consumption, which is a decision variable, converts to a constraint as follows: Objective: maximizing mean modulus M, Constraints: • The energy consumption must be less than 12.4 kW: 𝐸(𝑘𝑊ℎ) < 12.4 • The standard deviation of the modulus cannot be more than 20 GPa: 𝑆𝑇𝐷𝑚𝑜𝑑𝑢𝑙𝑢𝑠 < 20 Bounds: • 1100 < 𝑇1(°C) < 1175 • 1400 < 𝑇2(°C) < 1475 • 151 < 𝑡(𝑠) < 167 Based on these constraints and bounds, and solving the above optimization in Matlab using ε-constraints algorithm, the maximum mean modulus occurs when: • Temperature of zone 1 is 1100(°C), • Temperature of zone 2 is 1452(°C), and • Time is 167(s). Under these conditions, the mean modulus value is 236.14 (GPa). The amount of energy consumed during this carbonization process is 12.39 (kWh) and the standard deviation of the modulus is 16.66 GPa. Both energy consumption and standard deviation values meet the restrictions initially imposed in the problem. 74 4.4.2 Scenario 2: Concurrent optimization of the modulus and tensile strength High modulus and good strength are the primary attributes of carbon fibers. However, reaching the desired values of modulus and strength is expensive and time-consuming, because many carbon fiber manufacturing tests are required to find the best properties. Moreover, the experiments are influenced by uncertainty factors. Operator skills and a large number of filaments within each produced tow are the example sources of the uncertainty. One approach to find the optimum manufacturing conditions for producing the carbon fibers with both high modulus and high tensile strength means is using statistical models in order to optimize the process. The objective of scenario two is to use of the weighted sum method. In this method, the weights of normalized decision variables (normalized modulus and normalized tensile strength means) are chosen based on the importance of each factor and the preference of the decision maker [63]. Normalized values can be found by dividing the modulus and tensile strength values by their respective maxima. The STD of the mechanical properties and the energy consumption are deemed as secondary decision variables in this scenario and were converted to constraints, as follows: Objective: maximize w1Mnorm+w2UTSnorm, Constraints: • The standard deviation of the modulus to be less than 10 GPa: 𝑆𝑇𝐷𝑚𝑜𝑑𝑢𝑙𝑢𝑠 < 10 • The standard deviation of the tensile strength to be less than 0.9 GPa: 𝑆𝑇𝐷𝑠𝑡𝑟𝑒𝑛𝑔𝑡ℎ < 0.9 • The consumed energy to be less than 12.45 (kWh): 𝐸 (𝑘𝑊ℎ) < 12.45 Relative weights (preference values) w1 and w2 can be different for each designer. Bounds: • 1100 < 𝑇1(°C) < 1175 • 1400 < 𝑇2(°C) < 1475 75 • 151 < 𝑡(𝑠) < 167 Various Pareto optimal solutions can be obtained by varying the weights. Optimized modulus and tensile strength means are shown in table 4.1. 76 Table 4.1. Pareto optimal solutions using different weights under scenario 2. Modulus Weight Tensile Strength Weight Temperature of zone 1 (°C) Temperature of zone 2 (°C) Time (s) Modulus (GPa) Tensile Strength (GPa) Modulus STD (GPa) Tensile Strength STD (GPa) Energy consumption (kWh) 0.1 0.9 1124 1453 156 216.80 3.082 9.26 0.99 12.43 0.2 0.8 1126 1455 157 217.78 3.081 9.99 0.96 12.45 0.3 0.7 1124 1454 157 218.15 3.079 9.96 0.90 12.44 0.4 0.6 1120 1457 158 221.42 3.056 9.97 0.39 12.44 0.5 0.5 1118 1457 158 222.09 3.049 9.79 0.25 12.43 0.6 0.4 1114 1454 158 222.54 3.041 9.89 0.20 12.41 0.7 0.3 1114 1454 158 222.54 3.041 9.89 0.20 12.41 0.8 0.2 1113 1453 158 222.58 3.04 9.95 0.21 12.41 0.9 0.1 1113 1453 158 222.58 3.04 9.95 0.21 12.41 77 Figure 4.3 shows the Pareto optimal solutions based on table 4.1. A large gap is observed in the mid region of decision space between the two optimal points with the modulus weights of 30% and 40%. As described in section 4.4, in convex problems with non-constant curvature, the weighted sum method is incapable of predicting all the Pareto optimal solutions regardless of how the weights are chosen [65]. The black dashed line indicates an imaginary trend of the Pareto optimal points in the concave area, which can be mathematically identified using a method as follows. Figure 4.3. The convex region with non-constant curvature in the standard weighted sum method in scenario 2 10 , 90 20 , 8030 , 7040 , 6050 , 5060 , 40 70 , 3090 , 1080 , 203.0353.043.0453.053.0553.063.0653.073.0753.083.085216 217 218 219 220 221 222 223Tensile strength (GPa)Modulus (GPa)Modulus weight, Strength weight78 4.4.2.1 Adaptive weighted sum (AWS) method As it was shown in figure 4.3, the initial result of the standard weighted sum method gave a rough estimation of the optimal decision space solutions, containing a gap between two separate sets of the Pareto points. This gap is in essence due to the non-constant curvature of the convex Pareto where some optimal points have not been detected. The normal boundary intersection (NBI) and the adaptive weighted sum (AWS) are the two methods frequently used in order to capture such Pareto problems [66]. In normal boundary intersection method, the utopia points are found by the means of single objective optimizations. Using utopia line between anchor points the hidden Pareto optimal solutions are found. The AWS method, however, addresses the problem by applying new inequality constraints and providing a new sub-region for the method [64, 65]. Since the current optimization problem is a bi-objective problem, two constraints are required. Namely, these new constraints will be applied to the modulus mean and the tensile strength mean variables. Figure 4.4 shows the new sub-region provided by the initially applied constraints. In the next step, a sub-optimization is carried out in order to explore the sub-convex area and find the undetected points. The blue point in figure 4.4 is the first detected Pareto optimal solution in the convex region. Next, more sub-regions can be defined by using new constraints in order to find the remaining points. Figure 4.5 shows the updated sub-region and the resulting optimal points; for more mathematical details of the ASW method, the reader is refereed to [65]. Finally, figure 4.6 reveals all the detected Pareto optimal solutions and their trend, from which the decision maker can make a final decision to choose one of the points as optimum condition of the fiber carbonization process under scenario 2. The numbers adjacent to each data point indicate zone-1 temperature, zone-2 temperature and the heat treatment time, respectively. 79 Figure 4.4. The initial sub-optimization of AWS method in scenario 2 3.0353.043.0453.053.0553.063.0653.073.0753.083.085216 217 218 219 220 221 222 223Tensile strength (GPa)Modulus (GPa)3.0353.043.0453.053.0553.063.0653.073.0753.083.085216 217 218 219 220 221 222 223Tensile strength (GPa)Modulus (GPa)80 Figure 4.5. AWS method is applied to smaller regions in order to obtain remaining optimal solutions in scenario 2 Figure 4.6. Modulus mean against tensile strength mean plot of Pareto optimal solutions by adaptive weighted sum(AWS) method in in scenario 2. Process variables (T1, T2 and t) are also shown at each point, along with the ΔT which is T2-T1. 4.4.3 Scenario 3: concurrent optimization of the modulus and tensile strength The standard deviations of the modulus and tensile strength in scenario two were considered to have a fixed upper bound. However, in some design cases, the goal is that the standard deviations be less than a specific portion of the modulus and tensile strength means. Mimicking the latter case, the objective and bounds of scenario 3 is the same as scenario 2, however the constraints are different as follows: 3.0353.043.0453.053.0553.063.0653.073.0753.083.085216 217 218 219 220 221 222 223Tensile Strength (GPa)Modulus (GPa)1125°C1456°C157 s331 °C1122°C1455°C157 s333 °C 1120°C1456°C157 s336 °C1124°C1453°C156 s329 °C1126°C 1455°C 157 s 329°C 1124°C 1454°C 157 s 330°C 1120°C 1457°C 158 s 337 1118°C 1457°C 158 s 339 °C 1113°C 1453°C 158 s 340 °C 1114°C 1454°C 158 s 340°C 81 Objective: maximize w1Mnorm+w2UTSnorm Constraints: • The standard deviation of modulus to be less than 5% of the modulus mean : 𝑆𝑇𝐷𝑚𝑜𝑑𝑢𝑙𝑢𝑠 < 0.05𝑀 • The standard deviation of the tensile strength is less than 5% of the tensile strength mean: 𝑆𝑇𝐷𝑠𝑡𝑟𝑒𝑛𝑔𝑡ℎ < 0.05 𝑈𝑇𝑆 • The consumed energy to be less than 12.45 (kWh): 𝐸 (𝑘𝑊ℎ) < 12.45 Bounds: • 1100 < 𝑇1(°C) < 1175 • 1400 < 𝑇2(°C) < 1475 • 151 < 𝑡(𝑠) < 167 The Pareto optimal solutions for different modulus and strength values are shown in table 4.2 and figure 4.7. 82 Table 4.2. Pareto optimal solutions of scenario 3 using different weights Modulus Weight Tensile Strength Weight Temperature of Zone 1(°C) Temperature of Zone 2(°C) Time (s) Modulus (GPa) Tensile Strength (GPa) Modulus STD (GPa) Tensile Strength STD (GPa) Energy Consumption (kWh) 0.1 0.9 1120 1460 159 223.36 3.045 10.36 0.14 12.45 0.2 0.8 1118 1459 160 224.79 3.042 11.20 0.08 12.45 0.3 0.7 1118 1459 160 224.79 3.042 11.20 0.08 12.45 0.4 0.6 1118 1459 160 224.79 3.042 11.20 0.08 12.45 0.5 0.5 1116 1458 160 225.16 3.038 11.17 0.02 12.44 0.6 0.4 1115 1457 160 225.20 3.037 11.23 0.03 12.43 0.7 0.3 1115 1457 160 225.20 3.037 11.23 0.03 12.43 0.8 0.2 1115 1457 160 225.20 3.037 11.23 0.03 12.43 0.9 0.1 1115 1457 160 225.20 3.037 11.23 0.03 12.43 83 Figure 4.7. Modulus against tensile strength plot of Pareto optimal solutions under scenario 3. Values next to each pint indicate the relative weight of the modulus and tensile strength, respectively. 4.4.4 Scenario 4: energy consumption optimization As outlined in section 2.1.2.3, the carbonization process is extremely energy consuming due to its high temperature operating furnaces. Additionally, even small differences in the reported energy values for one run can make a large difference in practice as the production line is often a non-stop process; taking place several hours a day at industrial scale. Therefore, it is crucial for manufacturers to reduce the energy consumption in their facilities in order to reach higher financial and environmental benefits. They try to reach the desired range of carbon fibers’ modulus and tensile strength by consuming the lowest amount of energy, while ensuring the performance robustness of the produced fibers. Scenario 4 is a multi-objective optimization in which the energy consumption is minimized under constraints on the other objectives and bounds. Following the ε-10 , 9020 , 8030 , 7040 , 6050 , 5060 , 4070 , 3090 , 1080 , 203.0373.0383.0393.0403.0413.0423.0433.0443.0453.046223.00 223.50 224.00 224.50 225.00 225.50Tensile strength (GPa)Modulus (GPa)modulus weight, tensile strength weight84 constraints method, the remaining decision variables (modulus and tensile strength) are converted to constraints, including their standard deviations. The summary of the ensuing optimization model is as follows: Objective: Minimize E, Constraints: • The modulus of resulting carbon fiber to be between 190 GPa and 210 GPa: 190 < 𝑀 (𝐺𝑝𝑎) < 210 • The tensile strength of resulting carbon fiber to be between 1 GPa and 5 GPa: 1 < 𝑈𝑇𝑆 (𝐺𝑃𝑎) < 5 • The standard deviation of the modulus is less than 5% of the modulus mean : 𝑆𝑇𝐷𝑚𝑜𝑑𝑢𝑙𝑢𝑠 < 0.05𝑀 • The standard deviation of the tensile strength is less than 5% of the tensile strength mean: 𝑆𝑇𝐷𝑠𝑡𝑟𝑒𝑛𝑔𝑡ℎ < 0.05 𝑈𝑇𝑆 Bounds: • 1100 < 𝑇1(°C) < 1175 • 1400 < 𝑇2(°C) < 1475 • 151 < 𝑡(𝑠) < 167 Upon solving this optimization problem, the minimum amount of energy consumption occurred when: • Temperature of zone 1 is 1130 (°C), • Temperature of zone 2 is 1469 (°C), and • Time is 151 (s). 85 Under these conditions, the energy consumption is 12.44 (kWh) during the carbonization process. The modulus and tensile strength of the resulting carbon fiber are 214.192 GPa and 2.99 GPa. The modulus STD is 2.96 GPa and the standard deviation value of tensile strength is 0.13 GPa, all within the predefined limits. The reason that the optimized value of energy consumption is higher than the lowest energy consumption values of observed data (11.86 kWh) is the presence of constraints and bounds. In particular, requiring narrow ranges for modulus and tensile STDs is deemed a big obstacle for the carbonization process to reach lower energy consumption. Summary of findings In this chapter, the weighted sum and ε-constraints methods were discussed as the potential multi-objective optimization methods for fiber carbonization problem. Through the created experimental Pareto graphs, it was observed that inherently there are conflicts between the energy consumption objective and STD of fiber mechanical properties; which in essence suggest that the robustness of the fiber carbonization process is highly sensitive to the furnace temperatures and/or the heat treatment time, with a complex data trend. Furthermore, different optimization scenarios for the carbonization process were studied. These multi-objective optimization scenarios developed using the weighted sum method and the ε-constraints method, indicated that the constraints modification can vastly affect the optimal solutions of the primary objectives, specially for the energy minimization. Table 4.3 summarizes the four performed optimization scenarios. For the scenarios in which the weighted sum method is employed, only the Pareto optimal solution of equal weightings (modulus and tensile strength weights of 50%) is shown. 86 Figures 4.8 and 4.9 also depict Pareto optimal solutions of all scenarios along with the experimental data. The optimized solutions are indicated by triangles and the experimental data are shown by circles. Figure 4.9 indicates that the optimal solutions for modulus and tensile strength STDs are among the lowest standard deviation values of the experimental data. This is because of the employment of the STD constraints in the optimization scenarios. Furthermore, even though these constraints may act as obstacle for reaching the maximum modulus and tensile strength, as well as the lowest energy consumption, the results in figure 4.8 show that the Pareto optimal solutions possess satisfactory mean values of carbon fiber mechanical properties and the process energy consumption. For instance, the optimal solution of mean modulus in figure 4.8 is are among the highest modulus values found in the carbonization process, while satisfying multiple other constraints per scenario 1. 87 Table 4.3. Summary of the four optimization scenarios for carbonization process Scenario (main objective) Constraints Temperature of zone 1 (°C) Temperature of zone 2 (°C) Time (s) Modulus (GPa) Tensile Strength (GPa) Modulus STD (GPa) Tensile Strength STD (GPa) Energy (kWh) 1: Modulus optimization 1100 1452 167 236.14 2.85 16.66 0.61 12.39 2: Concurrent optimization of modulus and tensile strength 1114 1454 158 222.53 3.04 9.88 0.2 12.41 3: Concurrent optimization of modulus and tensile strength 1115 1457 160 225.20 3.037 11.23 0.029 12.43 4: Energy consumption optimization 1130 1469 151 214.19 2.99 2.96 0.13 12.44 𝐸(𝑘𝑊ℎ) < 12.4 𝑆𝑇𝐷𝑚𝑜𝑑𝑢𝑙𝑢𝑠 < 20 𝑆𝑇𝐷𝑚𝑜𝑑𝑢𝑙𝑢𝑠<10 𝑆𝑇𝐷𝑠𝑡𝑟𝑒𝑛𝑔𝑡ℎ < 0.9 𝐸(𝑘𝑊ℎ) < 12.45 𝑆𝑇𝐷𝑚𝑜𝑑𝑢𝑙𝑢𝑠 < 0.05𝑚𝑚𝑜𝑑𝑢𝑙𝑢𝑠 𝑆𝑇𝐷𝑠𝑡𝑟𝑒𝑛𝑔𝑡ℎ < 0.05𝑚𝑠𝑡𝑟𝑒𝑛𝑔𝑡ℎ 𝐸 (𝑘𝑊ℎ) < 12.45 190 < 𝑀 (𝐺𝑝𝑎) < 210 1 < 𝑈𝑇𝑆 (𝐺𝑃𝑎) < 5 𝑆𝑇𝐷𝑚𝑜𝑑𝑢𝑙𝑢𝑠 < 0.05𝑚𝑚𝑜𝑑𝑢𝑙𝑢𝑠 𝑆𝑇𝐷𝑠𝑡𝑟𝑒𝑛𝑔𝑡ℎ < 0.05𝑚𝑠𝑡𝑟𝑒𝑛𝑔𝑡ℎ 88 (a) (b) (c) Figure 4.8. Pareto optimal solutions along with experimental data in the performed scenarios 1-4, for mean mechanical properties and the energy consumption. The black circles are the observed data points and the triangles denote the optimum solution of each scenario. 00.511.522.533.540 50 100 150 200 250Tensile Strength (GPa)Modulus (GPa)123411.81212.212.412.612.8130 1 2 3 4Energy (kWh)Modulus (GPa)3411.81212.212.412.612.8130 50 100 150 200 250Energy (kWh)Modulus (GPa)12 3 4 1 2 89 (a) (b) (c) Figure 4.9. Pareto graphs of Pareto optimal solutions along with experimental data in the performed scenarios 1-4, for standard deviations (STD) of mechanical properties and the energy consumption. The black circles are the observed data points and the triangles denote the optimum solution of each scenario. 00.20.40.60.811.20 10 20 30 40Tensile Strength STD (GPa)Modulus STD (GPa)123411.81212.212.412.612.8130 10 20 30 40Energy (kWh)Modulus STD (GPa)11.81212.212.412.612.8130 0.5 1 1.5Energy (kWh)Tensile Strength STD (GPa)1 1 2 2 3 3 4 4 90 Chapter 5: Discussions Overview In this chapter, three different aspects of findings of this research regarding the carbonization process and the performed Gaussian Processes will be discussed. First, the GPs will be compared to standard methods such as least square regression. Next, the effect of the stabilization temperatures on the properties of the resulting carbon fibers will be scrutinized from process science point of view. Finally, the advantages of using a furnace with two zones will be highlighted. 5.1.1 How the Gaussian process (GP) would perform compared to a standard regression model Cross validation analysis was performed for both the GPs presented in Chapter 3 and two classic regression models (namely, a linear regression and a 2nd order polynomial). Figures 5.1-5.4 show the PRESS value of cross validations for the prediction of modulus, tensile strength, modulus STD and tensile strength STD. The modulus and ultimate tensile strength PRESS of the GP models are shown in black, and grey indicates the classic regression models’ PRESS. While at least one GP model outperforms the regression PRESS values for each mean response, a huge gap can be observed for the STD responses between several GPs’ and the classic regressions’ cross validation values. This confirms the fact that GPs work much better than other regression method specially in the presence of highly fluctuating data. The existence of hyper-parameters (free parameters) in the body of the covariance and mean functions in GP greatly increases the accuracy and reliability of the GPs’ model for such data. As discussed in section 2.3.4, the hyper-parameters specify the 91 characteristics of the covariance and mean functions, and they can be optimized based on the input and output properties. Figure 5.1. Comparison of Gaussian process (GP) models versus two classic regression models for modulus mean prediction 63358120 8160 824087201011086906600020004000600080001000012000SECovariancefunction andLinear meanfunctionPeriodicCovariancefunction andConstantmeanfunctionSE+PeriodicCovariancefunction andLinear meanfunctionSECovariancefunction andConstantmeanfunctionPeriodicCovariancefunction andLinear meanfunctionSE+PeriodicCovariancefunction andConstantmeanfunctionLinearregressionPolynomialregresionPRESS 1.21.39 1.482.44.61.861.600.511.522.533.544.55SE covariancefunction andlinear meanfunctionSE covariancefunction andconstant meanfunctionPeriodiccovariancefunction andconstant meanfunctionPolynomialcovariancefunction andlinear meanfunctionPolynomialcovariancefunction andconstant meanfunctionLinearregressionPolynomialregresionPRESS92 Figure 5.2. Comparison of Gaussian process (GP) models versus classic regression models for tensile strength mean prediction Figure 5.3. Comparison of Gaussian process (GP) models versus classic regression models for modulus STD prediction Figure 5.4. Comparison of Gaussian process (GP) models versus classic regression models for tensile strength STD prediction 15.4 16.317.820.223.226.439.531.00.05.010.015.020.025.030.035.040.045.0Periodiccovariancefunction andconstantmeanfunctionSEcovariancefunction andlinear meanfunctionPeriodiccovariancefunction andlinear meanfunctionSEcovariancefunction andconstantmeanfunctionLinearcovariancefunction andlinear meanfunctionLinearcovariancefunction andconstantmeanfunctionLinearregressionPolynomialregresionPRESS0.6 0.650.68 0.7 0.7 0.722.50.8800.511.522.53SE covariancefunction andlinear meanfunctionPeriodiccovariancefunction andLinear meanfunctionPeriodiccovariancefunction andconstant meanfunctionSE covariancefunction andconstant meanfunctionPolynomialcovariancefunction andlinear meanfunctionPolynomialcovariancefunction andconstant meanfunctionLinearregressionPolynomialregresionPRESS93 5.1.2 The effect of stabilization process temperature on the carbon fiber mechanical properties In section 2.2.1, it was discussed that the furnace temperature during the carbonization process can closely affect the mechanical properties of the resulting carbon fibers. Early research on the carbon fibers has revealed that the fibers’ modulus has a strong correlation with the process temperature [15]. Figure 5.5 reveals the comparison of mechanical properties between a PAN based carbon fiber and a mesophase pitch fiber [29]. It illustrates that the modulus has a steady growth as temperature increases; however, the tensile strength of PAN-based fibers starts decreasing at the temperature of about 1575°C. Figure 5.5. Modulus and tensile strength of polyacrylonitrile (PAN) based and mesophase pitch based fibers against processing temperature [29]. 94 The experimental data of the current carbonization process, however, state an opposite behavior, such that the carbon fibers’ mean modulus has decreased as the furnace temperatures increase (see fig. A.3). In other words, the highest modulus occurs when the zone-1 temperature is 1100°C (lowest value). This contrary behavior is believed to have originated from the gap between the stabilization and carbonization process furnace temperatures. Schematically, figure 5.6 depicts how the difference between these temperatures can negatively affect the final mechanical properties of the carbon fibers. It shows the carbon fibers’ temperature versus the time the fibers spend in the furnace. 𝑇0 is the stabilization process temperature, 𝑇1 and 𝑇2 are the temperatures of zone 1 and zone 2 in the carbonization process, respectively. Figure 5.6.a indicates the carbonization process with the lowest possible furnace temperatures (i.e. T1=1100°C and T2=1400°C). Similarly, figure 6.6.b depicts the same process procedure with the highest possible temperatures (T1=1175°C and T2=1475°C). It is assumed that the stabilization process temperature is constant and does not change for all the carbonization processes. The larger difference between 𝑇0 and 𝑇1 in the process (b) can cause a more sudden change in the carbon fibers’ behavior. This can be understood by observing the temperature profile slopes, which are drawn at the beginning of each process. The sharp heat rate in the process (b) would cause a sudden shock to the structure of the carbon fibers, resulting in unsatisfying product properties. 95 T0T1=1100T2=1400 T0T1=1175 T2=1475Time (s) Time (s)∆T2∆T1∆T3∆T4 (a) (b) Figure 5.6. The temperature of carbon fibers versus time in the carbonization process. 𝑻𝟎 is the stabilization process. The dashed dotted lines indicate the slope of the processes. Additionally, the above hypothesis may be confirmed by recalling the optimization results of the carbonization process. For instance, in the modulus optimization (i.e., scenario 1 in section 4.4.1), the optimized value of the mean modulus occurred when the furnace zone-1 temperature was at its lowest value (1100°C). Moreover, 𝑇2 was obtained at 1452°C, which is one of the highest temperatures of zone-2. In other words, the closer the 𝑇1 to 𝑇0, and the higher 𝑇2 with a sufficient time for fibers to spend within the furnace (i.e., for a full heat treatment at 𝑇2), the better the mechanical properties will be. Similar indications could be observed under scenario 2 (fig 4.6), where the Pareto curves is highly proportional to 𝑇2 − 𝑇1 , and modulus is maximized when 𝑇1is low. 5.1.3 The benefit of using two-zone furnace in the carbonization process Another important factor that can affect the tensile characteristics of the carbon fibers in the carbonization process is the heating rate during the process. Theoretically, it is essential to keep the heating rate in the desired slow levels, however it is a laborious task since the fibers’ 96 temperature must increase notably in a small amount of time to ensure fast production and hence, cost/energy saving. One solution is to use two separate controlled zones in the carbonization process furnace. This gives the fibers the chance to experience a roughly constant heating rate during the process. To elucidate, the green line in Figure 5.7 indicates how the employment of two zones can make the process closer to a constant heating rate. The blue line shows a process in which a one zone furnace is used. T0T1T2Time (s)∆T2∆T1Constant heating rate Zone 2 Zone 1 Figure 5.7. Comparison of two carbonization processes with difference furnace designs: the blue line indicates the process with one zone furnace; the green line shows a two-zone furnace carbonization process Summary of findings In this chapter the GP models regarding the prediction of the carbonization process were compared with classical regression models (namely, a linear regression and a 2nd order polynomial). The results indicated that GP performs well specially in the prediction of the fiber property standard deviations where high fluctuations existed in the input (measured) data. Furthermore, the investigation of the stabilization process revealed the fact that the huge difference between 97 carbonization and stabilization process temperatures can cause harmful effects on the resulting carbon fibers structures, which result in poor mechanical properties (modulus and tensile strength). Finally, it was discussed that the implementation of the two-zone furnaces instead of one zone can greatly help the carbonization process to be performed with slow heat transfer rates. It is believed the best properties are found when (i) the carbonization and stabilization process temperatures are close, (ii) the heat rate is rather slow and controlled (constant) in the furnace (i.e., by using multiple heating zones), and (iii) the maximum furnace temperature is at the higher bounds. 98 Chapter 6: Conclusions and Future Work Conclusions GPs are an effective and robust tool for capturing the behavior of complex manufacturing processes under noisy responses [67]. Conversely, more classic prediction modelling approaches, e.g. standard regression, would have a weaker performance compared to GPs when the data has a highly non-smooth nature, as was also noticed in the current case study. Additionally, one important benefit of GPs is that unlike regression methods, it is not necessary to provide a predefined explicit relation between inputs and outputs. The covariance and mean functions are the factors that are required to be specified by analysis and GPs establish the relations in an unsupervised framework based on the properties of the prior functions and the given experimental evidence (data) [68]. Finally, GPs favorably predict a fit error at each prediction point separately, as opposed to the standard regression which gives a single overall fit error; or e.g. neural networks which normally predict only one output characteristic (e.g. the response mean) with no error estimation. To select the best GP model with carbonization process response variables, evaluation tools such as NML and PRESS were used. The overall performance of the selected models assessed by statistical measures including RSS and R-squared, to measure the closeness of the prediction model to the experimental data. The model assessment results clearly confirmed that GPs could successfully model and predict the carbonization process response variables, with goodness of fit values over 99% under highly fluctuating data for some of the responses such as modulus and tensile strength STDs. In the optimization stage of the research, common multi-objective 99 optimization methods (namely, the weighted sum method and ε-constraints method) were adapted and used under different process design scenarios, based on potentially different needs of carbon fiber customers. The specific conclusions driven from the study may be: • It was shown that, unlike classic regression models, for the cases that the standard deviation of the training data is unequal (heteroscedasticity), appropriate prediction models could be built using GP. To accomplish this task, two separate GP models were developed. The first model predicts the mean value of the subject variable, and the second model predicts the standard deviation. • The squared exponential covariance function was ranked best for all the GP predictions of the carbonization process response variables; the SE covariance function provided the lowest NML and PRESS values. Smoothness is deemed to be the main feature of the GP models with the SE covariance function. • Advantages of using GPs in the prediction of the carbon fiber mechanical properties and the energy consumption responses under uncertainty, became evident through evaluating the results of the cross validations. In all prediction cases, the selected GP models performed better than classic regression methods. In particular, the prediction of the modulus and tensile strength STDS was where the GP worked significantly better than classic models. 100 • Although the covariance function plays the primary role in the GP model development and it has the most significant effect on the trend and behavior of the posterior, the selection of a good mean function is also important for building a complete GP model. Analysis of the NML and PRESS confirmed that the linear mean function is the preferred choice for the prediction of the fibres modulus, tensile strength and their standard deviations. For the prediction of energy consumption, the constant mean function was chosen since it provided the least NML and PRESS. • The standard weighted sum method (WSM) was employed for the Pareto optimization of the carbon fibers’ modulus and tensile strength. However, this method does not capture all of the Pareto efficient solutions, regardless of how the weights are perturbed. Hence, the application of the adaptive WSM method was shown and proved effective in capturing concave regions of the carbon fiber optimization Pareto. • From the optimization results, the maximum modulus value is reachable when zone-1 temperature is low and zone-2 temperature is high. This confirms the hypothesis that, in order to prevent any probable shock in the carbon fibers’ microstructure, the difference between the stabilization process temperature and the temperature of the furnace zone-1 should be kept as close as possible. Furthermore, the reason that the maximum modulus occurs when zone-2 temperature is at its high levels is higher level of heat revetment, providing sufficient time is given to fibers to spend in the furnace at reach a steady state. 101 • It was concluded from the comparison of the two multi-objective optimization scenarios 2 and 3 that the constraints play a crucial role for obtaining Pareto optimal solutions. It was shown that applying design-driven modifications in the constraints changes the optimal solutions of the optimization scenarios with the same objectives. • It was found that the modulus and tensile strength responses do not have an evident correlation with respect to each other, and with the energy consumption response, in the obtained generic Pareto fronts. In addition, the STDs of mechanical properties showed a clear conflict with energy minimization goal (the less energy consumed, the less robust the fibers’ mechanical properties will be). Hence, a multi-objective optimization approach is required to acquire overall optimum process conditions. • A summary of the results of the performed multi-objective optimization scenario 2 is shown in figure 6.1, as it is believed it would be most practical scenario for carbon fiber customers to optimize design of their composite structures; the hidden mid region of the Pareto has not been included as it showed sharp changes in the design space in fig. 4.6 and hence would not be preferable by manufacturers in practice. The bubble chart in figure 6.1 indicates that there is a nonlinear tradeoff between the fiber modulus and tensile strength means. Thus, the optimized results are highly depended on the designer’s preference regarding the criteria weighing factors. In addition, the marginal rate of substitution between the two fiber mechanical properties is not linear (i.e., the designer should give much higher weight than 50% to the tensile strength in order to increase it notably). The energy consumption has a direct positive relation with the maximization of tensile strength, while maximizing the modulus can be achieved with relatively lower energies. In addition, 102 the modulus STD and the tensile strength STD bars indicate very disproportionate decision space for designers (e.g., a higher mean value of modulus does not necessarily mean a higher STD from the carbonization process). Figure 6.1 is deemed to be a useful general tradeoff graph for the carbonization process that may be used by process engineers and structural designers for fast review and recommendations. Figure 6.1. Recommended design chart (n-dimensional Pareto) for multi-objective optimization of mechanical properties of fibers, while the carbonization process energy is constrained and the robustness of fibers performance is taken into account Future work recommendations Although through this thesis it was shown that GP prediction methods can successfully model and predict the carbonization process variables during carbon fibers manufacturing, it cannot be denied 1124°C 1453°C 156 s ΔT = 329 °C 1126°C 1455°C 157 s ΔT = 329 °C 1124°C 1454°C 157 s ΔT = 330 °C 1113°C 1453°C 158 s ΔT = 340 °C 1114°C 1454°C 158 s ΔT = 340 °C 1118°C 1457°C 158 s ΔT = 339 °C 1120°C 1457°C 158 s ΔT = 337 °C 103 there are also other powerful machine learning methods (e.g., neural network, support vector regression/SVR, thin plate splines and convex hull technique) that can adequately characterize such complex processes. Moreover, GP itself has some disadvantages, which need to be considered. The large number of available covariance and mean functions for utilizing in GP’s structure makes the model selection often a difficult task. In addition, using Gaussian process in cases with large sample sizes can increase the calculation cost exponentially. Hence, investigating the ability of other black-box learning methods for modeling and predicting the carbonization process will be the next step of this study. Moreover, here the multi-objective optimization results suggested that the best fiber properties may be realized when (i) the carbonization and stabilization process temperatures are close, (ii) the heat rate is rather slow and controlled (constant) in the furnace (i.e., by using multiple heating zones), and (iii) the maximum furnace temperature is at the higher bounds. Follow-up validation tests can be performed to validate these practical suggestions. Finally, multi-physics finite element models may be established to avoid the costly experimental trials of the carbon manufacturing process. 104 Bibliography [1] R. Friedel, P.B. Israel, Edison's electric light: the art of invention, JHU Press2010. [2] M. Josephson, The Invention of the Electric Light, Scientific American 201 (1959) 98-114. [3] W.B. Carlson, M.E. Gorman, Understanding invention as a cognitive process: The case of Thomas Edison and early motion pictures, 1888-91, Social Studies of Science 20(3) (1990) 387-430. [4] K. Shimazaki, M. Hirai, Studies on preparation of polyacrylonitrile based activated carbon fiber, Nippon Kagaku Kaishi 7 (1992) 739-744. [5] P. Demmer, M. Fowler, A.A. Marino, Use of carbon fibers in the reconstruction of knee ligaments, Clinical orthopaedics and related research 271 (1991) 225-232. [6] Park, Soo-Jin. Carbon fibers. Netherlands, Dordrecht: Springer, 2015. [7] L. Manocha, O. Bahl, G. Jain, Length changes in PAN fibres during their pyrolysis to carbon fibres, Macromolecular Materials and Engineering 67(1) (1978) 11-29. [8] M.K. Titsias, Variational Learning of Inducing Variables in Sparse Gaussian Processes, AISTATS, 2009, pp. 567-574. [9] J. Hensman, N. Fusi, N.D. Lawrence, Gaussian processes for big data, arXiv preprint arXiv:1309.6835 (2013). [10] Rasmussen, Carl Edward, and Christopher KI Williams. Gaussian processes for machine learning. Vol. 1. Cambridge: MIT press, 2006. [11] R.B. Gramacy, B. Haaland, Speeding up neighborhood search in local Gaussian process prediction, Technometrics 58(3) (2016) 294-303. 105 [12] N. Corrado, M. Gherlone, C. Surace, J. Hensman, N. Durrande, Damage localisation in delaminated composite plates using a Gaussian process approach, Meccanica 50(10) (2015) 2537-2546. [13] D. Petelin, A. Grancharova, J. Kocijan, Evolving Gaussian process models for prediction of ozone concentration in the air, Simulation modelling practice and theory 33 (2013) 68-80. [14] J.D. Whitcomb, Composite materials: testing and design (eighth conference), ASTM International1988. [15] P. Morgan, Carbon fibers and their composites, CRC press2005. [16] A. Shindo, Report No. 317, Government Industrial Research Institute, Osaka, Japan (1961). [17] R. Moreton, W. Watt, W. Johnson, Carbon fibres of high strength and high breaking strain, Nature 213(5077) (1967) 690-691. [18] D. Clark, N. Wadsworth, W. Watt, Handbook of composites 1, strong fibers, Elsevier Science Publisher BV, 1985. [19] M.S. Blackie, D.J. Poynton, Cellulose fibres for cement reinforcement, Google Patents, 1987. [20] L. Peebles, Carbon fibres: Structure and mechanical properties, International materials reviews 39(2) (1994) 75-92. [21] S.-J. Park, G.-Y. Heo, Precursors and manufacturing of carbon fibers, Carbon Fibers, Springer2015, pp. 31-66. [22] A. Gupta, I. Harrison, New aspects in the oxidative stabilization of PAN-based carbon fibers, Carbon 34(11) (1996) 1427-1445. [23] J. Bromley, Carbon fibres, their composites and applications, Proc. Int. Carbon Fibre Conf., London, 1971, p. 3. [24] T.H. Ko, The influence of pyrolysis on physical properties and microstructure of modified PAN fibers during carbonization, Journal of applied polymer science 43(3) (1991) 589-600. 106 [25] M.S.A. Rahaman, A.F. Ismail, A. Mustafa, A review of heat treatment on polyacrylonitrile fiber, Polymer Degradation and Stability 92(8) (2007) 1421-1432. [26] M. Kaburagi, Y. Bin, D. Zhu, C. Xu, M. Matsuo, Small angle X-ray scattering from voids within fibers during the stabilization and carbonization stages, Carbon 41(5) (2003) 915-926. [27] E. Fitzer, W. Frohs, M. Heine, Optimization of stabilization and carbonization treatment of PAN fibres and structural characterization of the resulting carbon fibres, Carbon 24(4) (1986) 387-395. [28] H. Khayyam, M. Naebe, A. Bab-Hadiashar, F. Jamshidi, Q. Li, S. Atkiss, D. Buckmaster, B. Fox, Stochastic optimization models for energy management in carbonization process of carbon fiber production, Applied Energy 158 (2015) 643-655. [29] T. Matsumoto, Mesophase pitch and its carbon fibers, Pure and applied chemistry 57(11) (1985) 1553-1562. [30] J. Mittal, R. Mathur, O. Bahl, Post spinning modification of PAN fibres—a review, Carbon 35(12) (1997) 1713-1721. [31] N. Yusof, A.F. Ismail, Post spinning and pyrolysis processes of polyacrylonitrile (PAN)-based carbon fiber and activated carbon fiber: A review, Journal of Analytical and Applied Pyrolysis 93 (2012) 1-13. [32] M. Trinquecoste, J. Carlier, A. Derré, P. Delhaes, P. Chadeyron, High temperature thermal and mechanical properties of high tensile carbon single filaments, Carbon 34(7) (1996) 923-929. [33] E. Zussman, X. Chen, W. Ding, L. Calabri, D. Dikin, J. Quintana, R. Ruoff, Mechanical and structural characterization of electrospun PAN-derived carbon nanofibers, Carbon 43(10) (2005) 2175-2185. 107 [34] M. Dobb, H. Guo, D. Johnson, C. Par, Structure-compressional property relations in carbon fibres, Carbon 33(11) (1995) 1553-1559. [35] S. Ozbek, D. Isaac, Strain-induced density changes in PAN-based carbon fibres, Carbon 38(14) (2000) 2007-2016. [36] J.S. Tsai, C.H. Lin, The effect of the side chain of acrylate comonomers on the orientation, pore‐size distribution, and properties of polyacrylonitrile precursor and resulting carbon fiber, Journal of applied polymer science 42(11) (1991) 3039-3044. [37] Z. Wu, D. Pan, X. Fan, B. Qian, Thermal shrinkage behavior of preoxidized polyacrylonitrile fibers during carbonization, Journal of applied polymer science 33(8) (1987) 2877-2884. [38] H. Khayyam, S.M. Fakhrhoseini, J.S. Church, A.S. Milani, A. Bab-Hadiashar, R. Jazar, M. Naebe, Predictive Modelling and Optimization of Carbon Fiber Mechanical Properties through High Temperature Furnace, Applied Thermal Engineering (2017). [39] J. Yuan, K. Wang, T. Yu, M. Fang, Reliable multi-objective optimization of high-speed WEDM process based on Gaussian process regression, International Journal of Machine Tools and Manufacture 48(1) (2008) 47-60. [40] Vapnik, Vladimir Naumovich, and Vlamimir Vapnik. Statistical learning theory. Vol. 1. New York: Wiley, 1998. [41] A. Jøsang, Generalising Bayes’ Theorem in Subjective Logic, International Conference on Multisensor Fusion and Integration for Intelligent Systems (MFI 2016), 2016. [42] J. Quiñonero-Candela, C.E. Rasmussen, A unifying view of sparse approximate Gaussian process regression, Journal of Machine Learning Research 6(Dec) (2005) 1939-1959. [43] C.E. Rasmussen, Advances in Gaussian processes, Advances in Neural Information Processing Systems 19 (2006). 108 [44] C.J. Paciorek, Nonstationary Gaussian processes for regression and spatial modelling, Citeseer, 2003. [45] D. Barber, Bayesian reasoning and machine learning, Cambridge University Press2012. [46] R. Grbić, D. Kurtagić, D. Slišković, Stream water temperature prediction based on Gaussian process regression, Expert systems with applications 40(18) (2013) 7407-7414. [47] Y. Bengio, Y. Grandvalet, No unbiased estimator of the variance of k-fold cross-validation, Journal of machine learning research 5(Sep) (2004) 1089-1105. [48] R. Kohavi, A study of cross-validation and bootstrap for accuracy estimation and model selection, Ijcai, Stanford, CA, 1995, pp. 1137-1145. [49] C.D. Keeling, T.P. Whorf, Atmospheric CO2 concentrations derived from flask air samples at sites in the SIO network in Trends: A Compendium of Data on Global Change, Carbon Dioxide Information Analysis Center, Oak Ridge National Laboratory, US Department of Energy, Oak Ridge, Tennessee, USA (2004). [50] X. Qin, Y. Lu, H. Xiao, Y. Wen, T. Yu, A comparison of the effect of graphitization on microstructures and properties of polyacrylonitrile and mesophase pitch-based carbon fibers, Carbon 50(12) (2012) 4459-4469. [51] E. Cauer, W. Mathis, R. Pauli, Life and work of wilhelm cauer (1900 1945), Proc. 14th Int. Symp. Mathematical Theory of Networks and Systems, MTNS, 2000, pp. 1-10. [52] C.E. Rasmussen, H. Nickisch, Gaussian processes for machine learning (GPML) toolbox, Journal of Machine Learning Research 11(Nov) (2010) 3011-3015. [53] R. Urtasun, D.J. Fleet, P. Fua, 3D people tracking with Gaussian process dynamical models, Computer Vision and Pattern Recognition, 2006 IEEE Computer Society Conference on, IEEE, 2006, pp. 238-245. 109 [54] T.S. Breusch, A.R. Pagan, A simple test for heteroscedasticity and random coefficient variation, Econometrica: Journal of the Econometric Society (1979) 1287-1294. [55] R.F. Engle, Autoregressive conditional heteroscedasticity with estimates of the variance of United Kingdom inflation, Econometrica: Journal of the Econometric Society (1982) 987-1007. [56] J. Morgan, J. Tatar, Calculation of the residual sum of squares for all possible regressions, Technometrics 14(2) (1972) 317-325. [57] A.C. Cameron, F.A. Windmeijer, An R-squared measure of goodness of fit for some common nonlinear regression models, Journal of Econometrics 77(2) (1997) 329-342. [58] R.L. Rardin, Optimization in operations research, Prentice Hall Upper Saddle River, NJ1998. [59] L. Zadeh, Optimality and non-scalar-valued performance criteria, IEEE transactions on Automatic Control 8(1) (1963) 59-60. [60] S. Gass, T. Saaty, The computational algorithm for the parametric objective function, Naval Research Logistics (NRL) 2(1‐2) (1955) 39-45. [61] V. Chankong, Y. Haimes, Optimization-based methods for multiobjective decision-making-an overview, Large Scale Systems In Information And Decision Technologies 5(1) (1983) 1-33. [62] Y.Y. Haimes, L. Ladson, D.A. Wismer, Bicriterion formulation of problems of integrated system identification and system optimization, IEEE Transactions on Systems Man and Cybernetics (3) (1971) 296-&. [63] Y. Censor, Pareto optimality in multiobjective problems, Applied Mathematics & Optimization 4(1) (1977) 41-59. 110 [64] I.Y. Kim, O. de Weck, Adaptive weighted sum method for multiobjective optimization, 10th AIAA/ISSMO multidisciplinary analysis and optimization conference, Albany, New York, 2004. [65] I.Y. Kim, O.L. de Weck, Adaptive weighted-sum method for bi-objective optimization: Pareto front generation, Structural and multidisciplinary optimization 29(2) (2005) 149-158. [66] I. Das, J.E. Dennis, Normal-boundary intersection: A new method for generating the Pareto surface in nonlinear multicriteria optimization problems, SIAM Journal on Optimization 8(3) (1998) 631-657. [67] Y. Gal, M. van der Wilk, C.E. Rasmussen, Distributed variational inference in sparse Gaussian process regression and latent variable models, Advances in Neural Information Processing Systems, 2014, pp. 3257-3265. [68] F.M. Gray, M. Schmidt, Thermal building modelling using Gaussian processes, Energy and Buildings 119 (2016) 119-128. 111 Appendices Appendix A : Detailed information of the carbonization process Table A.1 Important measurements of the carbonization process furnace [28] 112 Table A.2. Raw experimental data from the carbonization process (each row/test was repeated 50 times) Zone one's temperature (°C) Zone two's temperature (°C) Time (s) Modulus (Gpa) Modulus STD (GPa) Tensile Strength (GPa) Tensile Strength STD (GPa) Energy Consumption (kWh) 1100 1400 156 228.8 18.47 3.016 1.06 11.87 1100 1400 167 227.7 9.288 3.201 0.68 11.95 1100 1400 162 227.5 9.89 3.275 0.77 11.92 1100 1400 151 234.6 13.67 3.458 0.76 11.86 1125 1425 156 228.3 14.48 2.925 1 12.19 1125 1425 167 229.2 34.03 3.185 0.81 12.44 1125 1425 162 226.7 11.66 3.031 0.81 12.4 1125 1425 151 231.9 15.9 3.225 0.59 12.37 1150 1450 156 219.8 7.797 2.945 0.63 12.49 1150 1450 167 221.8 10.11 3.13 0.87 12.58 1150 1450 162 217.6 32.99 3.11 0.65 12.51 1150 1450 151 224.3 10.78 3.18 0.63 12.47 1175 1475 156 197.7 35.32 2.78 0.68 12.7 1175 1475 167 154.2 9.5 2.11 0.36 12.81 1175 1475 162 157.5 9.13 2.35 0.41 12.77 1175 1475 151 163.02 5.39 2.209 0.34 12.65 113 Figure A.1. Modulus versus time plot of the carbonization process, using raw experimental data at T1=1125 °C. 050100150200250300350400150 152 154 156 158 160 162 164 166 168Modulus (GPa)Time (s)114 Figure A.2. Tensile strength versus time plot of the carbonization process, using raw experimental data at T1=1125 °C. 0123456150 152 154 156 158 160 162 164 166 168Tensile Strength (GPa)Time (s)115 Figure A.3. Modulus versus zone-1 temperature plot of the carbonization process, using raw experimental data at t= 162(s). 0501001502002503001090 1100 1110 1120 1130 1140 1150 1160 1170 1180Modulus (GPa)Zone-1 temperature (°C )116 Figure A.4. Tensile strength versus zone-1 temperature plot of the carbonization process, using raw experimental data at t=162(s). 00.511.522.533.544.551090 1100 1110 1120 1130 1140 1150 1160 1170 1180Tensile strength (GPa)Zone-1 temperature (°C )
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Multi-objective Gaussian process approach for robust...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Multi-objective Gaussian process approach for robust optimization and prediction of carbonization process Ramezankhani, Milad 2017
pdf
Page Metadata
Item Metadata
Title | Multi-objective Gaussian process approach for robust optimization and prediction of carbonization process |
Creator |
Ramezankhani, Milad |
Publisher | University of British Columbia |
Date Issued | 2017 |
Description | Carbon fibers are high-performance and high-strength reinforcement materials in advanced industries such as aerospace, automotive and energy sector. This class of materials is often derived from polyacrylonitrile (PAN) fiber precursor. The conversion of polyacrylonitrile precursor to carbon fibers is essentially comprised of two major stages; namely, stabilization and carbonization. The carbonization is a very energy-consuming and expensive step due to its high temperature requirements. In order to achieve desirable physical and mechanical properties of carbon fibers through this step, a large amount of energy is required. A cost-effective approach to optimize energy consumption of this process, however, is the use of predictive modelling techniques. In that goal, herein, a Gaussian Process (GP) approach is proposed to firstly, predict multiple mechanical properties of carbon fibers in the presence of manufacturing noise and secondly, optimize them under a minimum energy consumption criterion and a range of process constraints and bounds. The resulting model for each property of carbon fibers consists of two Gaussian Process models. The first model describes the mean value of the property and the second one predicts its standard deviation. The proposed Gaussian Process approach is compared to a traditional regression approach using the measurements obtained from the Carbon Nexus pilot plant in Australia. The Gaussian Process approach clearly proved to be more effective in terms of both prediction accuracy and robustness. Through employing Gaussian Process models, the modulus and tensile strength mean values of carbon fibers along with their standard deviations (STD) were successfully predicted under different process conditions. Squared exponential covariance and linear mean functions were proven to be most suitable in constructing the Gaussian Process in the performed case study. It was found that the modulus and tensile strength responses do not have an evident correlation with respect to each other, hence a multi-objective optimization approach was developed to acquire overall optimum process conditions. To estimate the trade-off between fiber material properties under the multi-objective optimization problem, a standard as well as an adaptive weighted sum method was applied under various constraints and bounds, including the energy consumed during carbonization. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2017-09-05 |
Provider | Vancouver : University of British Columbia Library |
Rights | Attribution-NonCommercial-NoDerivatives 4.0 International |
DOI | 10.14288/1.0355395 |
URI | http://hdl.handle.net/2429/62961 |
Degree |
Master of Applied Science - MASc |
Program |
Mechanical Engineering |
Affiliation |
Applied Science, Faculty of Engineering, School of (Okanagan) |
Degree Grantor | University of British Columbia |
GraduationDate | 2017-09 |
Campus |
UBCO |
Scholarly Level | Graduate |
Rights URI | http://creativecommons.org/licenses/by-nc-nd/4.0/ |
AggregatedSourceRepository | DSpace |
Download
- Media
- 24-ubc_2017_september_ramezankhani_milad.pdf [ 2.96MB ]
- Metadata
- JSON: 24-1.0355395.json
- JSON-LD: 24-1.0355395-ld.json
- RDF/XML (Pretty): 24-1.0355395-rdf.xml
- RDF/JSON: 24-1.0355395-rdf.json
- Turtle: 24-1.0355395-turtle.txt
- N-Triples: 24-1.0355395-rdf-ntriples.txt
- Original Record: 24-1.0355395-source.json
- Full Text
- 24-1.0355395-fulltext.txt
- Citation
- 24-1.0355395.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.24.1-0355395/manifest