MK: 7a.! . :. AI It \ “Nib. si] ti o I .r ... 1.. .. - ....n.. L5,? . lithihflyukfl 11.1.891374 59.)}: .1: V li‘ ’ 1. 9.. ‘ u :0. A .5. {finnmvfish 2!. .. . A a» .- ‘ inf. mar vxltrr ... .K 3.3!. .734 3 .5 t. 1 av. :11! x. )u93911u bl.’il.lu I, {Irv-1|). :ll; 1.1.1.4! )1)... fi . a! 1 ,1] . \ vii-‘12 l 24!}. c,‘i...f!...¢ \1.9 V .t .419; 2 7:73.... . » .lns. . 711...... .y.‘ 1:.~... :41 .i) 7‘... : it: A.» 8&1— 1...)!» a l: .3 . .. V 6.3!.5. . . . _ , .‘ .Efl ; ""3 UBRARY 2%” Michigan State University This is to certify that the dissertation entitled PROCESS OPTIMIZATION FOR MOIST AIR IMPINGEMENT COOKING OF MEAT PATTIES presented by SANGHYUP JEONG has been accepted towards fulfillment of the requirements for the PhD. degree in Biosystems Engineering WM; Cirtléjot‘lPrfifessor’s Signature December 14, 2005 Date MSU is an Affirmative Action/Equal Opportunity Institution —.--.------o---o—c--------I-v-.-o--~--o----.---a—u---o-l-I---‘-c-0-»-I-a-u-t-n----o----c-----. PLACE IN RETURN BOX to remove this checkout from your record. TO AVOID FINES return on or before date due. MAY BE RECALLED with earlier due date if requested. DATE DUE DATE DUE DATE DUE ”Rf-l (a 2914M 2/05 p:/CIRC/DaleDue.indd-p.1 PROCESS OPTIMIZATION FOR MOIST AIR IMPINGEMENT COOKING OF MEAT PATTIES By Sanghyup Jeong A DISSERTATION Submitted to Michigan State University in partial fulfillment of the requirements for the degree of DOCTOR OF PHILOSOPHY Department of Biosystems and Agricultural Engineering 2005 ABSTRACT PROCESS OPTIMIZATION FOR MOIST-AIR IMPINGEMENT COOKING OF MEAT PATTIES By Sanghyup Jeong Process conditions of a moist air impingement cooking system were optimized to achieve maximum yield and to satisfy safety and quality constraints, simultaneously. To accomplish this goal, various strategies were tested by combining different modeling approaches, global optimization algorithms, and parameterization of control profiles. In this study, a finite element model (FEM) predicting yield, Salmonella inactivation, internal color change, and surface color change was considered as an actual experiment with which all the results were compared. Static neural network models (SNNM) and a dynamic neural network model (DNNM) were utilized as potential, faster alternatives to the finite element model. For the global optimization algorithms, genetic algorithms (GA), simulated annealing (SA), and integrated controlled random search in dynamic system (ICRS/DS) algorithms were tested along with the finite element model and alternative models. In addition, piecewise linear interpolation (PLI) and Fourier series (FS) were used for the control profile parameterization. This study was conducted in two different ways. In the first part, overall aspects of this optimization problem and the effectiveness of the various strategies were investigated to identify the best strategy for ideal dynamic control profiles. Secondly, based on prior knowledge, the optimization strategies were applied to several industrially-relevant case studies. The performance of the alternative models (DNNM and SNNM) was fast, general, and robust, with a few exceptions. Even though the accuracy and the power of classification were not as high as the finite element model results, the neural network models showed potential as reliable alternative models. The highest goal (yield) was 73%, which was obtained by using the ICRS algorithm, FEM, and PLI. However, the optimization strategies with alternative models could not find such a high yield; rather, they committed critical classification errors at the later stages of the optimization process. Generally, all the global optimization algorithms showed convergence to an optimal solution, albeit with different convergence speed and goal achievement. Although comprehensive evaluation was impossible, ICRS was observed as the most recommendable algorithm. Single-stage, double-stage, and multi-zone processes were studied by using three different models (FEM, DNNM, and SNNM) and the ICRS algorithm. The maximum yield (67%) was achieved in the double-stage process. The case studies showed that a simple and minor design change of the single-stage oven might improve the performance. In addition, the objective fianction (yield) for the single-stage oven was replaced with a cost function, and the operating conditions for maximum profit were determined, which were different from the results when the objective function was yield. Finally, Monte Carlo simulation showed that all the optimal profiles were highly sensitive to small perturbations, which implied difficulties in the actual application of the optimal solutions, due to unavoidable control errors of a cooking system. To my family iv ACKNOWLEDGEMENTS I would like to express my sincere appreciation to Dr. Bradley P. Marks, a truly inspiring teacher and counselor; for his guidance and friendship during the course of this investigation. Grateful appreciation is also extended to Dr. Alden M. Booren, Dr. Kirk Dolan and Dr. James Steffe for serving on my guidance committee. The author is indebted to the United States Department of Agriculture (USDA) Cooperative State Research, Education, and Extension Service for supporting project 2003-51110-02081 of which this study is a part. Thanks go to our research team members and colleagues of graduate office 5. Special thanks also goes to college group of the New Hope Baptist Church for their prayer and good memories in God. Dr. Sukjung Chang, MD and his wife deserve to have my sincere thanks for their spiritual encouragement all the way through this journey. Of course, my family and parents have been the essence of refreshment for every new challenge in my life. TABLE OF CONTENTS LIST OF TABLES - _ ix LIST OF FIGURES ...... - _ xi NOMENCLATURE- - -- - - -- - - -- - xvi 1 INTRODUCTION - - - - - - - -- -- l 1.1 Background ......................................................................................................... 1 1.1.1 Food Quality and Safety ............................................................................. 1 1.1.2 The Drive for Food Quality and Safety Innovation .................................... 3 1.1.3 Economic Significance of the Meat Industry .............................................. 5 1.1.4 Recent Innovation in Meat Patty Processing .............................................. 6 1 .2 Objectives ........................................................................................................... 8 2 LITERATURE REVIEW - - -- - - - - 9 2.1 Achieving Extrema in Food Processing .............................................................. 9 2.2 Characteristics of Optimization in Food Processing ........................................... 9 2.3 The State of the Art of Optimization Techniques ............................................. 12 2.4 Applications in Food Process Engineering ....................................................... 14 2.4.1 Sterilization of Canned Product ................................................................ 14 2.4.2 Food Dehydration ..................................................................................... 17 2.4.3 Meat Patty Cooking .................................................................................. 19 2.4.4 Various Processing Areas ......................................................................... 20 3 THEORIES- - - - - 23 3.1 Overview of Optimization Theory .................................................................... 23 3. 1 .1 Introduction ............................................................................................... 23 3.1.2 Theory and Methods ................................................................................. 26 3.2 Global Optimization Techniques ...................................................................... 40 3.2.1 Genetic Algorithm (GA) ........................................................................... 40 3.2.2 Simulated Annealing (SA) ........................................................................ 44 3.2.3 Integrated Controlled Random Search for Dynamic Systems (ICRS/D8) 48 3.3 Artificial Neural Network ................................................................................. 53 3 .3. 1 Fundamentals ............................................................................................ 53 3.3.2 Back-Propagation F eed-Forward Neural Network (FF NN) ..................... 55 3.3.3 Generalized Regression Neural Network (GRNN) ................................... 61 4 METHODS AND PROCEDURES _- - 65 4.1 Overall Methods and Procedures ...................................................................... 65 4.2 Impingement Cooking Technology .................................................................. 65 4.3 Integrated Process Model Development ........................................................... 68 vi 4.4 4.5 4.6 4.7 4.8 5.1 5.2 5.3 5.4 5.5 5.6 4.3.1 Basic Finite Element Model ...................................................................... 69 4.3.2 Kinetic Parameters for Microbial Inactivation in a Meat Patty ................ 73 4.3.3 Kinetic Parameters for Internal Color Change of a Meat Patty ................ 74 4.3.4 Kinetic Parameters for Surface Color Change of a Meat Patty ................ 76 4.3.5 Practical Considerations for Integration of Models .................................. 82 Formulating the Optimization Problem ............................................................ 84 Alternative Modeling by Using Artificial Neural Networks ............................ 86 4.5.1 Parameterization of the Control Function ................................................. 86 4.5.2 Training and Validation Data Groups ....................................................... 88 4.5.3 Static Training ........................................................................................... 90 4.5.4 Dynamic Training ..................................................................................... 92 4.5.5 Neural Network Validation ....................................................................... 96 4.5.6 Neural Network Parameter Optimization ................................................. 96 Process Optimization Strategies ....................................................................... 97 4.6.1 Combinations of Techniques .................................................................... 97 4.6.2 GA Based Strategy .................................................................................... 98 4.6.3 SA Based Strategy .................................................................................... 99 4.6.4 ICRS/DSBased Strategy ........................................................................ 101 4.6.5 Benchmark Test for Optimization Algorithms ....................................... 102 Case Studies .................................................................................................... 103 4.7.1 Single-Stage Oven .................................................................................. 104 4.7.2 Double-Stage Oven ................................................................................. 104 4.7.3 Multi—Zone Oven .................................................................................... 105 4.7.4 An Example of Economic-Based Optimization ...................................... 105 Programming and Computation Tools ............................................................ 109 RESULTS AND DISCUSSION - 111 Integrated Process Model Performance .......................................................... 111 Training and Optimizing Neural Networks .................................................... 114 5.2.1 Dynamic Neural Network Model (DNNM) ............................................ 114 5.2.2 Static Neural Network Models (SNNM) ................................................ 116 Validation of the Trained Neural Networks .................................................... 119 5.3.1 The Performance of DNNM ................................................................... 119 5.3.2 The Performance of SNNM .................................................................... 130 Benchmark Test for Optimization Algorithms ............................................... 135 Optimal Conditions from the Various Strategies ............................................ 137 5.5.1 Summary of the Various Optimization Strategies .................................. 137 5.5.2 Optimization Performance of the Alternative Models Coupled with Algorithms .............................................................................................. 142 5.5.3 Efficiencies of Various Strategies ........................................................... 143 5.5.4 Constraints Satisfaction Problem ............................................................ 146 5.5.5 Sensitivity of the Optimal Control Profiles ............................................ 149 Case Studies .................................................................................................... 152 5.6.1 Single-Stage Process ............................................................................... 152 5.6.2 Double-Stage Process ............................................................................. 157 5.6.3 Multi-Zone Process ................................................................................. 160 vii 5.6.4 Economic-Based Optimization ............................................................... 164 6 CONCLUSIONS - 167 6. 1 General Conclusions ....................................................................................... 167 6.2 Suggestions for Future Research .................................................................... 169 APPENDICES - 172 APPENDIX A Tables and Figures of Optimization Process .................................. 173 APPENDIX B Computer Codes ............................................................................. 198 REFERENCES . -- - -- ..... 257 viii LIST OF TABLES Table 2.1 Analysis of basic elements used to optimize sterilization processes of prepackaged conduction-heated foods. ..................................................................... 16 Table 3.1 Unconstrained multivariable optimization methods categorized by how they generate the search direction (Edgar et al., 2001 ;Venkataraman, 2002). ................. 31 Table 4.1 Comparison chart of thermal inactivation parameters of Salmonella cocktail for ground beef (1 :(Murphy et al., 2002); 2:(Juneja et al., 2001); 3:(Murphy et al., 2004)). ....................................................................................................................... 73 Table 4.2 Reference decimal reduction time (Dr) and z-value of internal color change in the various muscles at the reference temperature of 60 °C (Geileskey et al., 1998). 75 Table 4.3 Digitized data of surface time-temperature plot (recipe A) and time-color change plot (recipe C) of Dagerskog and Bengtsson (1974). ................................... 78 Table 4.4 The names and contents of data groups for training and testing neural networks. ................................................................................................................................... 89 Table 4.5 Conditions of each process type to produce training data groups; D1, single- stage (D3), double-stage (D4), and multi-zone (D5) process. .................................. 90 Table 4.6 Optimization strategies (by strategy number), according to their process model, the method of control profile parameterization, and optimization algorithm. .......... 97 Table 5.1 The performance of DNNM for a random data group (D2 of Table 4.4), in terms of RMSE. ...................................................................................................... 119 Table 5.2 Four possible classification categories of the constraints satisfaction ............ 123 Table 5.3 Classification rate of DNNM for random processes ....................................... 125 Table 5.4 The performance of DNNM for single-stage, double-stage, and multi-zone processes in terms of RMSE. .................................................................................. 126 Table 5.5 Classification rate of DNNM for single-stage, double-stage, and multi-zone processes. ................................................................................................................ 127 Table 5.6 The performance (RMSE) of SNNM for random, single-stage, double-stage, and multi-zone processes. ....................................................................................... 131 ix Table 5.7 Classification rate of SNNM for random, single-stage, double-stage, and multi- zone process. ........................................................................................................... 134 Table 5.8 Comparison of the convergence of three different optimization algorithms in terms of RMSE for the Bezier parametric curve fitting problem. .......................... 135 Table 5.9 Brief results of the various optimization strategies. ....................................... 138 Table 5.10 Predicted goal and constraint values of DNNM and SNNM, compared with the result of FEM. ................................................................................................... 142 Table 5.11 Total execution time and the achieved maximum yield for various optimization strategies. ........................................................................................... 144 Table 5.12 Sensitivity of the optimal solutions for F EM_ICRS_PLI and FEM_SA_FS was obtained by using Monte Carlo simulation ...................................................... 150 Table 5.13 The results of the single-stage optimization process by using ICRS-FEM.. 152 Table 5.14 The accuracy of ICRS-DNNM (strategy #26) and ICRS-SNNM (strategy #27) were validated by using FEM results. ..................................................................... 155 Table 5.15 Maximum yield and optimal control variables were found by using ICRS- FEM for double-stage process. ............................................................................... 158 Table 5.16 Maximum yield and optimal control variables were found by using ICRS- FEM and ICRS-DNNM for multi-zone process. .................................................... 162 Table 5.17 Sensitivity of the optimal solutions was obtained by using Monte Carlo simulation for economic-based optimization problem. .......................................... 166 LIST OF FIGURES Figure 1.1 Schematic diagram of a multi-staged, moist-air impingement cooking system. ..................................................................................................................................... 7 Figure 3.1 Locating minimum point of a quadratic function with first derivative information ................................................................................................................ 27 Figure 3.2 Illustration of the operation of captain “A” and captain “B”. ......................... 30 Figure 3.3 Concept of linear programming problem in a two-dimensional case .............. 33 Figure 3.4 Global and local optimum in 2-dimensional case (The picture was adopted from the help manual of MATLAB®) ....................................................................... 35 Figure 3.5 General procedure of a genetic algorithm. ...................................................... 42 Figure 3.6 An example of simple crossover operation in a genetic algorithm. ................ 43 Figure 3.7 Flowchart of modified ICRS algorithm. .......................................................... 52 Figure 3.8 Physical structure of a neuron. ........................................................................ 54 Figure 3.9 Mathematical model of biological neuron ....................................................... 55 Figure 3.10 Topology of back-propagation of feed-forward neural network. .................. 57 Figure 3.11 Architecture of generalized regression neural network (GRNN). ................. 63 Figure 3.12 GRNN response depending on the size of receptive field (a) in two- dimensional case. ...................................................................................................... 64 Figure 4.1 A schematic diagram of overall procedures and methods. .............................. 66 Figure 4.2 An actual and a cross-sectional image of the ISO-IV Jet Stream® Oven of Stein-DSI (FMC FoodTech). .................................................................................... 68 Figure 4.3 Finite element mesh utilized for one quarter of the 2-D cylindrical model product. ..................................................................................................................... 72 Figure 4.4 A sigmoid function (logistic function) of surface temperature (Ts) showing gradual change in a certain factor around a critical value, Tdcw in this case. ............ 81 Figure 4.5 Examples of random time-varying conditions by using Fourier series control profile parameterization. ........................................................................................... 91 xi Figure 4.6 Conceptual diagram of dynamic training paradigm. ....................................... 93 Figure 4.7 The structure of an input/output data pair for dynamic training to identify dynamic characteristics of meat patty cooking. ........................................................ 94 Figure 4.8 Acceptance probability (Boltzmann probability) along with cooling schedule (T 0:25, T N=O.001, N=150, and average increment of accepted objective value=0.3) for SA. ..................................................................................................................... 100 Figure 5.1 An example of the integrated process model performance under constant process condition: (a) Process conditions, yield, center, and surface temperature; (b) Salmonella inactivation; (0) Internal color change; ((1) Surface color change ........ 112 Figure 5.2 An example of the integrated process model performance under dynamic process condition: (a) Process conditions, yield, center, and surface temperature; (b) Salmonella inactivation; (c) Internal color change; ((1) Surface color change ........ 113 Figure 5.3 Optimal “spread” values at the lowest RSME were determined by trying “spread” for each GRNN in DNNM: (a) Yield prediction; (b) Salmonella inactivation, internal color change, and surface color change. ............................... 115 Figure 5.4 Average normalized RMSE of each combination for the number of neurons in each hidden layer of SNNM_R. .............................................................................. 117 Figure 5.5 RMSE (yield, Salmonella inactivation, internal color change, and surface color change) of SNNM_D for different “spread” values. .............................................. 118 Figure 5.6 RMSE of SNNM_M for different “spread” values: (a) Yield prediction; (b) Constraints prediction. ............................................................................................ 1 18 Figure 5.7 Example qualitative performance of the DNNM, with respect to: (a) Yield; (b) Salmonella inactivation; (c) Internal color change; ((1) Surface color change ........ 120 Figure 5.8 The performance of DNNM for random process as a predictor and classifier. (a) Yield; (b) Salmonella inactivation; (c) Internal color change; ((1) Surface color change (some data in NE region was not plotted). ................................................. 121 Figure 5.9 Accepted patties (170 patties) were superimposed on the 1,000 random data. ................................................................................................................................. 124 Figure 5.10 Goal and constraints prediction performance of the DNNM for single-stage process: (a) Yield; (b) Salmonella inactivation; (c) Internal color change; ((1) Surface color change. ........................................................................................................... 128 Figure 5.11 Goal and constraints prediction performance of the SNNM for random Fourier process condition: (a) Yield; (b) Salmonella inactivation; (c) Internal color change; (d) Surface color change ............................................................................ 132 xii Figure 5.12 Examples of objective firnction value (RMSE) for every generation of three different optimization algorithms in the Bezier parametric curve fitting problem. (a) GA; (b) SA; (c) ICRS. ............................................................................................ 136 Figure 5.13 ICRS-F EM optimization process and convergence (strategy #3): (a) PLI parameterization; (b) F S parameterization .............................................................. 139 Figure 5.14 Optimal control profile of temperature, humidity, impingement velocity, and cooking duration found by ICRS-FEM optimization strategy (strategy #3): (a) PLI parameterization; (b) F S parameterization .............................................................. 140 Figure 5.15 Optimal control profile of temperature, humidity, impingement velocity, and cooking duration found by SA-FEM optimization strategy (strategy #6): (a) PLI parameterization; (b) FS parameterization .............................................................. 141 Figure 5.16 Convergence history of each control parameters shows exploration and exploitation features of the optimization algorithms: (a) GA-FEM-PLI (strategy #1); (b) SA-FEM-FS (strategy #6); (c) ICRS-FEM-PLI (strategy #3). ......................... 145 Figure 5.17 A typical search path committing false-pass classification errors (ICRS- DNNM-F S): (3) Yield; (b) Salmonella reduction; (c) Internal color change; ((1) Surface color change. .............................................................................................. 147 Figure 5.18 Directional constraint satisfaction (DCS) algorithm prevents search from drifting toward false-pass region (ICRS-DNNM-FS): (a) Yield; (b) Salmonella reduction; (c) Internal color change; ((1) Surface color change. ............................. 148 Figure 5.19 Histograms of objectives and constraints were plotted from the population of Monte Carlo simulation for the strategy #6: (a) Yield distribution of passed simulations; (b) Salmonella inactivation distribution; (c) Internal color change distribution; ((1) Surface color change distribution. ................................................ 151 Figure 5.20 Optimization processes for single-stage oven: (a) Process for low temperature and high humidity optimal profile; (b) Process for high temperature and low humidity profile. ..................................................................................................... 153 Figure 5.21 Optimization processes and optimal control profiles for double-stage oven obtained by ICRS-FEM strategy: (a) Optimization process; (b) Optimal control profile of trial #1; (c) Optimal control profile of trial #2. ....................................... 159 Figure 5.22 Optimization process and optimal control profile were obtained by using ICRS-FEM and ICRS-DNNM for multi-zone process: (a) Optimization process of strategy #33; (b) Optimal control profiles of strategy #33; (c) Optimization process of strategy #34; (d) Optimal control profiles of strategy #34. ................................ 163 Figure 5.23 Optimization processes and optimal control variables for maximum profit of single-stage oven were found by using ICRS-F EM strategy: (a) Convergence history; (b) Comparison of net profit and patty yield; (c) Various cost history. 165 xiii Figure A.1 The performance of DNNM for double-stage process: (a) Yield; (b) Salmonella inactivation; (c) Internal color change; (d) Surface color change ........ 174 Figure A.2 The performance of the DNNM for multi-zone process: (a) Yield; (b) Salmonella inactivation; (c) Internal color change; ((1) Surface color change ........ 176 Figure A.3 The performance of the SNNM for single-stage process: (a) Yield; (b) Salmonella inactivation; (c) Internal color change; ((1) Surface color change ........ 178 Figure A.4 The performance of the SNNM for double-stage process: (a) Yield; (b) Salmonella inactivation; (0) Internal color change; ((1) Surface color change ........ 180 Figure A5 The performance of the SNNM for multi-zone process: (a) Yield; (b) Salmonella inactivation; (c) Internal color change; ((1) Surface color change ........ 182 Figure A.6 GA-F EM optimization process and convergence: (a) PLI parameterization; (b) F S parameterization ........................................................................................... 184 Figure A.7 Optimal control profile of temperature, humidity, impingement velocity, and cooking duration found by GA-FEM optimization strategy: (a) PLI parameterization; (b) FS parameterization .............................................................. 185 Figure A.8 SA-FEM optimization process and convergence: (a) PLI parameterization; (b) F S parameterization. ............................................................................................... 186 Figure A.9 GA-DNNM optimization process and convergence: (a) PLI parameterization; (b) FS parameterization ........................................................................................... 187 Figure A.10 Optimal control profile of temperature, humidity, impingement velocity, and cooking duration found by GA-DNNM optimization strategy: (a) PLI parameterization; (b) F S parameterization .............................................................. 188 Figure A.11 SA-DNNM optimization process and convergence: (3) PLI parameterization; (b) FS parameterization ........................................................................................... 189 Figure A. 12 Optimal control profile of temperature, humidity, impingement velocity, and cooking duration found by SA-DNNM optimization strategy: (a) PLI parameterization; (b) FS parameterization .............................................................. 190 Figure A.13 ICRS-DNNM optimization process and convergence: (a) PLI parameterization; (b) FS parameterization .............................................................. 191 Figure A. 14 Optimal control profile of temperature, humidity, impingement velocity, and cooking duration found by ICRS-DNNM optimization strategy: (a) PLI parameterization; (b) F S parameterization .............................................................. 192 xiv Figure A. 15 Optimization process and optimal control profiles found by GA-SNNM_R optimization strategy and Fourier series parameterization: (a) Convergence history; (b) Optimal control profile of temperature, humidity, impingement velocity, and cooking duration. .................................................................................................... 193 Figure A.16 Optimization process and optimal control profiles found by SA-SNNM_R optimization strategy and Fourier series parameterization: (a) Convergence history; (b) Optimal control profile of temperature, humidity, impingement velocity, and cooking duration. .................................................................................................... 194 Figure A.17 Optimization process and optimal control profiles found by ICRS-SNNM_R optimization strategy and Fourier series parameterization: (a) Convergence history; (b) Optimal control profile of temperature, humidity, impingement velocity, and cooking duration. .................................................................................................... 195 Figure A.18 Two different optimization processes and optimal control profiles found by single-stage ICRS-FEM optimization strategy: (a) Convergence history of case I (LTHH); (b) Convergence history of case 11 (HTLH); (c) Optimal control profile of temperature, humidity, impingement velocity, and cooking duration of case I; ((1) Optimal control profile of temperature, humidity, impingement velocity, and cooking duration of case 11. .................................................................................... 196 XV NOMENCLATURE Roman Letters: §rm\a5m laseesreeckos 5 as a parameter of a sigmoid activation function color value of redness coefficient of Fourier series kth coefficients of Fourier series color value of yellowness concentration, [kg/m3] surface color change cooling rate mass capacity, [kg/kg] specific heat capacity, [J/(kg °C)] decimal reduction time, [s] capillary diffusivity, [mz/s] fat capillary diffusivity, [mz/s] expected value activation energy, [J/gmol] failure counter in ICRS/D8 nonlinear activation fimction, or probability density function humidity (% moisture by volume), [%MV] output of hidden layer neuron, or enthalpy, [kJ/kg] convective mass transfer coefficient, [m/s] convective heat transfer coefficient, [W/(m2 °C)] internal color change at center objective function, or performance index Boltzmann constant, or inactivation rate constant, [l/s] heuristic parameter of ICRS/DS heuristic parameter of ICRS/DS moisture conductivity, [kgmoism/(m 5)] thermal conductivity, [W/(m °C)] color value of lightness log cycles of microbial destruction a multiplier, or moisture or fat content-decimal dry basis number of control variables final averaged fat content based on non-fat solids, [kg/kg] initial fat content based on non-fat solids, [kg/kg] final average water content based on non-fat solids, [kg/kg] initial water content based on non-fat solids, [kg/kg] neighborhood function, number of discretization, total number of cooling cycle or number of microorganism, [CFU/g] normal to surface heuristic parameter in ICRS/DS xvi SN~<‘<§‘UV D1 gm parameters of control functions, or probability ideal gas constant, [J/(kmol K)] random number, or r-direction in radial coordinates step size temperature, [°C] time, [3] total execution time, [s] dew point temperature, [°C] impingement exit velocity, [m/s] sum of neuron weight output of neuron patty yield, [%] vertical direction in cylindrical coordinates, or microbial z-value, [°C] total color change Greek Letters: geompmfigm—IQ search direction vector of random numbers (Gaussian distribution) error, or small number efficiency of algorithm [% yield/s] learning rate constant time coordinates of control profile parameterization grid minimum interval, or latent heat [kJ/kg] decision variables vector density, [kg/m3] standard deviation vector, or constant for the size of receptive field control coordinates of control profile parameterization grid Bold style (vector or matrix): x€:~’~: output vector of individual neuron target vector control vector, or training vector weight matrix input column vector, or control vector xvii Superscript: HL hidden layer i ith iteration , point, or initial state j jth neuron k kth iteration, or time step, or neuron OL output layer Subscript: 0 initial state c critical value f final state L lower limit 11 number of step r reference state 3 surface U upper limit Abbreviations: ANN artificial neural network BB branch and bound FEM finite element model FFNN feed-forward neural network GA genetic algorithm GO global optimization GRG generalized reduced gradient GRNN generalized regression neural network ICRS/DS integrated controlled random search for dynamic system LP linear programming MILP mixed integer linear programming MINLP mixed integer nonlinear programming MIP mixed integer programming NLP nonlinear programming OMb oxymyoglobin PLI piecewise linear interpolation QP quadratic programming RBF radial basis-function, or radial basis-function neural network RMSE root mean squared error SA simulated annealing SLP successive linear programming SSE sum of squared error xviii 1 INTRODUCTION 1.1 Background 1.1.] Food Quality and Safety Food processing technology has been advanced along with human history. However, the fundamental concept of food processing has not changed considerably, even though many innovative processes and products are being developed. Basically, manufacturing of foods encompasses two types of conversions, physical and chemical. The main purposes of these conversions are preserving and improving the quality of processed food products. One of the predominant unit operations, inducing both physical and chemical conversions, is heat or thermal processing. Thermal processing is transferring heat to a food material to induce desirable results, but this also can cause concurrent undesirable results. Thermal processing increases digestibility (e.g., protein denaturation), reduces enzyme and microorganism activity, and enhances food characteristics (e.g., carbohydrate gelatinization, color development, texture, and flavor changes). However, thermal processing also produces undesirable results, such as loss of heat sensitive nutrients and undesirable color and flavor changes due to over-cooking. Therefore, a compromise must be found between intensity of thermal processing and its various effects on the product (Trystram, 2004). According to the results of a recent survey published in “Trends In The United States - Consumer Attitudes and The Supermarket 1999”(FMI, 1999), the top food selection concerns and the percentages of the shopping public that consider these factors “very important” in their food selection were as follows: (1) Taste, (92% of those interviewed), (2) Nutrition, (70%), (3) Product Safety (70%), (4) Price, (63%), (5) 1 Storability, (42%); (6) Ease of Preparation, (35%); (7) Food Preparation Time (35%) and (8) Product packaging that can be recycled, (29%). Therefore, food manufacturing must address not just production volume, but also a variety of other competing criteria. For the consumer’s food safety concern, the USDA Economic Research Service (ERS) reports that the percent of consumers “completely confident” in the safety of the food supply increased from a low of 72% in 1992-93 to a high of 83% in 1996, with levels declining to 74% in 2000 (ERS, 2002). From the economic viewpoint, food safety is the most critical and nonnegotiable quality factor. In the United States, foodbome diseases have been estimated to cause 6 million to 81 million illnesses and up to 9,000 deaths each year (Mead et al., 1999). For six specific bacterial pathogens, the costs of human illness are estimated to be $9.3-12.9 billion annually (Buzby et al., 1996). In 2000, ERS estimated the annual costsl due to selected foodbome pathogens2 as $6.9 billion (ERS, 2000). These estimated costs are an enormous burden for society and also for food manufacturers. Therefore, it is clear that the major sectors (i.e., consumers, industry, and regulatory agencies) in the food market must collaborate to improve the current situation. ' includes medical costs, productivity losses, and costs of premature deaths 2 Campylobacter (all serotypes), Salmonella (nontyphoidal), E. coli 0157, E. coli non-0157 STEC, and Listeria monocytogenes 1.1.2 The Drive for Food Quality and Safety Innovation Driving forces can be passive or active. Since the 1993 outbreak of E. coli 0157:H7, consumer awareness and demand for food safety has increased (Golan et al., 2004). The passive driving forces often come from foodbome illness outbreaks and recalls3, and they trigger consumer awareness, regulatory changes, or industrial innovation. “The number and size of recalls have increased dramatically over the last decade. During 1993-96, the number of meat and poultry Class I4 recalls averaged about 24 per year and amounted to 1.5 million pounds annually; during 1997-2000, Class I recalls averaged 41 per year and reached 24 million pounds annually (Ollinger and Ballenger, 2003).” These increasing recall cases are not because of loose control, but because the ability to detect pathogens on products has increased dramatically, which can generate more recalls (AMI, 2002). In addition, the ability to track foodbome disease and tie it to a specific food product has evolved into a practical technology (AMI, 2002). Recalls result in bad reputation and catastrophic financial damage to a manufacturer. Therefore, efficient quality assurance has become a critical issue for consumers, manufacturers, and related government organizations. Instrumentation, food safety practices, and lethality criteria are of central importance, with particular emphasis on very high sanitary and hygienic operating standards. Evolving federal regulations, such as 3 A food recall is a voluntary action by a manufacturer or distributor to protect the public from products that may cause health problems or possible death. The purpose of a recall is to remove meat or poultry from commerce when there is reason to believe it may be adulterated (injurious to health or unfit for human consumption) or misbranded (false or misleading labeling and/or packaging) (F SIS, “F SIS Recalls”, USDA, http://www.fsis.usda.gov/Fsis_Recalls/index.asp, March 22, 2005.). 4 Recalls that involve meat or poultry products that could, especially without cooking to safe temperature, cause serious illness or death. 9CFR318.17 (F SIS, 1999 & 2001), change safety regulations from passive to active compliance required of the food industry. Traditionally, regulations have provided a specific endpoint temperature and holding time to achieve target lethality in a meat and poultry product. Thus, the traditional approach discourages the food industry from adopting new technology and voluntary compliance to the regulations. However, the evolving regulations require manufacturers to prove, via scientifically supportable means, that their process or Operating policy achieves a target lethality performance standard. The transition from command-and-control to performance standards allows more freedom of choosing process design and operation policy, but also moves more responsibility to the industry. Contrary to the above passive driving forces, there might be an active driving force that originates from industry. Food manufacturers invest in the development of new methodologies to improve safety and quality of their food product. When industry successfully innovates to produce safe foods, a win-win situation arises, with the innovating firm, consumers, and government all benefiting from improved food safety (Golan et al., 2004). Many attempts have been made to maximize desirable quality and simultaneously minimize undesirable effects by adding ingredients, developing innovative process equipment, improving process conditions, and so forth. Among those approaches, improving process conditions is advantageous, in that it uses existing systems. Therefore, additional capital investment is not necessary to resolve these contradicting factors to achieve both safety requirements and maximize quality and profit. Finding the best operating condition is critical from the perspective of industry, because the need to improve efficiency, reduce energy consumption, increase productivity, and comply with regulations pushes industry to adopt improved safety practices if they can see simultaneous benefits in yield, which translate to profit. 1.1.3 Economic Significance of the Meat Industry In spite of the increased safety concerns, the US. meat and poultry industry contributes significantly to the US. agricultural economy. Total meat and poultry production in 2000 exceeded 80 billion pounds, a 31 percent increase since 1987. The meat and poultry industry is the largest segment of US. agriculture, contributing over $100 billion in annual sales to the GNP (AMI, 2001). The products affected by regulatory changes related to ready-to-eat products account for over $28 billion in annual sales (F SIS, 2001), and consumer trends for ready- to-eat products also suggest continued rapid growth in this category. The size of this market is important as a spur to greater profit in this industry. Accordingly, the industry aims to develop novel products and to increase the efficiency of its production lines. For an example, even a modest 0.5% improvement of cooking yields based on the $28 billion annual sales in this category would give an impact of approximately $140 million increase in annual revenue for ready-to-eat products in the US. Therefore, given the regulatory changes and the economic importance of ready- to-eat (RTE) meat products, there are compelling needs for integrated simulation tools that will allow industry to design and operate processes that meet the lethality performance standards and simultaneously increase quality and profit. Optimization techniques that find the conditions for the best result from a given situation are needed to meet these demands. Process optimization is the most economic approach for industry to maximize yields while ensuring safety and quality factors with existing facilities. Therefore, the information and tools that will enable the industry to design and operate the optimal processes are essential. I. 1.4 Recent Innovation in Meat Patty Processing Impingement cooking technology has been popular in certain segments of the food processing industry, because of its efficiency. Specifically, moist air impingement cooking systems are widely used in the ready-to-eat meat product industry, because the system (Figure 1.1) results in short cooking time and relatively high cooking yield. The moist air impingement cooking systems jets a steam-air mixture through arrays of nozzles onto products, yielding a high heat transfer rate by reducing the thickness of the boundary layer at the surface of the product. Also, at the initial stage of cooking, steam is condensed on the surface of the product, which results in effective transfer of latent heat into the product at low temperature. Therefore, moist air impingement cooking systems are characterized by fast cooking and suppression of moisture loss. Moist air impingement cooking systems involve many control variables, such as cooking duration, air temperature, air moisture content, impingement exit velocity, and impingement geometry (e.g., jet width, spacing, and height). Impingement unit \ .WWL ..n...s...9...9. ? Product \JWLJMQWL IIaII'IIIsIII’IIifi IIROI’IIIQIII’CIIfi U. Conveyc’r/ AAAAAAAAAA AMAAAAMA AAAAAAAAAA belt Oven 3 J W SIDE VIEW db debdb {Product i > END VIEW “ O O'O'O'O'O'C'O'O'CTO'OTO'O'.'Q’.'.'.'.'.'.'.’O'.'O conveyor belt " Figure 1.1 Schematic diagram of a multi-staged, moist-air impingement cooking system. Currently, most oven operators select the oven operating conditions ‘(i.e., the control variables) based on their experience or simple rules of thumb, which are not necessarily proven scientifically. For a single-stage oven, one might operate the oven within sub-optimal conditions. However, if multiple ovens (Figure 1.1) are connected to increase the rate of production with various cooking zones, then the complexity and size of the problem is too large to select optimal conditions with experience. Also, meat cooking under moist air impingement cooking environments involves multiple control Variables, mass transfer coupled with heat transfer, phase transitions, and complex cOndensing boundary conditions, which results in complex combinations of control PrOfiles. Therefore, optimizing the process conditions of moist air impingement cooking systems is a significant challenge and a worthwhile endeavor from the perspective of researchers, food processors, oven manufacturers, and government agencies. 1.2 Objectives The overall goal of this study was to develop an efficient process optimization method, in terms of speed and objective-achievement, and to evaluate the method for maximizing cooking yield, while ensuring the microbial safety and quality of ready-to-eat meat products (ground and formed meat and poultry products) cooked in commercial moist air impingement cooking systems. To achieve the overall goal, the specific objectives were: 1. To add color prediction capabilities to an existing finite element model. 2. To develop an alternative process model for moist air impingement cooking of meat patties by using artificial neural network (ANN) to replace an existing finite element model. 3. To identify the best strategy for process optimization among many combinations of optimization algorithms, process models, and parameterization of control functions. 4. To apply the developed optimization strategies to three case studies: single- stage, double-stage, and multi-zone oven systems. 5. To examine the performance of the optimization strategy for a single-stage oven system, given an economic-based objective function. 2 LITERATURE REVIEW 2.1 Achieving Extrema in Food Processing An ultimate goal of most production activity is to increase profitability at given conditions, and food processing is no exception. Typically, improving the efficiency, reducing process time, and increasing yield and quality are concerns for most food processors and equipment designers. Generally, a typical industry consists of management, process design and equipment specification, and plant operations (Edgar et al., 2001). Because these components are inter-connected, achieving those improvements are not simple tasks. Therefore, depending on the level of complexity and difficulty, the scope of a problem can be the entire enterprise, 3 plant, a process, a single unit operation, a single piece of equipment in that operation, or any intermediate stage between these (Beveridge and Schechter, 1970). Generally, these improving activities are maximizing the capabilities of existing facilities or equipments by changing their conventional operating policies, conditions, numbers, and so on. All these attempts can be described in a single word: “optimization.” A more formal definition would be “the collective process of finding the set of conditions required to achieve the best result from a given situation (Beveridge and Schechter, 1970).” 2-2 Characteristics of Optimization in Food Processing Food processing is unique in that the process involves materials having irregular Shapes, non-homogeneous compositions, and individual variance even in the same material. In addition to the materials, quality factors for product evaluation can be subjective, such as taste, aroma, and flavor. Also, food processing encompasses various techniques, such as frying, baking, boiling, blanching, fermenting, drying, and so on. Therefore, modeling food processing phenomena is very challenging. Models are essential components of modern process systems engineering (i. e., simulation, optimization, and control), and they are usually classified into three categories, which are first-principle models (or white-box), data-driven models (or black- box), and hybrid models (or gray-box) (Banga et al., 2003). Without adequate models, it is impossible to carry out optimization. From the view point of the application of optimization technique, the first-principle models are highly desirable, because the response time of the models is short, and it is convenient to apply various mathematical operations, such as differentiation. However, because first-principle models are difficult to obtain, data-driven models and hybrid models are popular in food processing. Generally, mathematical modeling of a food process requires knowledge of transport phenomena and reaction kinetics. Transport phenomena involve heat, mass, and momentum transfer into a food, and reaction kinetics cover degradation or inactivation of nutritional and organoleptic factors or microbial and enzymatic activity (Oliveira and Oliveira, 1999). Usually, these multi-physical phenomena are expressed as sets of algebraic, partial, and ordinary differential equations in the mathematical model. Due to the lack of theoretical methods for solving those highly complex, nonlinear systems of equations, most of the problems are solved by using numerical techniques, such as the finite difference method and the finite element method. Considering that most 0Immization techniques require numerous iterations of a process model, the numerical, 10 computational requirements of the model are often significant impediments to optimization. One of the popular methods to overcome the above difficulties of modeling is alternative modeling, which can be fast, simple, and robust, with reliable accuracy. Even though alternative models sacrifice some degree of accuracy, compared to numerical models, the overall benefits of the alternative modeling techniques must be considered when evaluating performance in an actual optimization problem. Artificial neural networks (ANN) have emerged as a potential alternative to physical-based models for food process engineering, because of their simple structure, robustness, no requirement of prior knowledge, and adaptive performance (Torrecilla et al., 2004). ANN have been successfully applied to the modeling of food processes (Mittal and Zhang, 2000; Sablani and Shayya, 2001; Chen and Ramaswamy, 2002; Chen and Ramaswamy, 2003; Horiuchi et al., 2004; Torrecilla et al., 2004), property and quality prediction (Berg et al., 1997; Xie and Xiong, 1999; Raptis et al., 2000; Albert et al., 2001; Therdthai and Zhou, 2001; Tominaga et al., 2001; Hussain et al., 2002; Boillereaux et al., 2003; Ganjyal et al., 2003), machine vision and image analysis (Chao et al., 2002; Marique et al., 2003; Diaz et al., 2004), extrusion (Ganjyal and Hanna, 2002), microbial growth and inactivation (Geeraerd et al., 1998; Garcia-Gimeno et al., 2003), high- pressure processes (Torrecilla et al., 2005), and fluid flow (Adhikari and Jindal, 2000; Sablani and Shayya, 2003; Singh and Jindal, 2003). Once an ANN is established, the Computation time of the network is very small with reliable accuracy. ANN is ideal for SyStem identification and replacement of an existing first-principle model. For example, it 11 was shown that ANN were very effective for replacing a finite difference computer simulation of retort processes (Chen and Ramaswamy, 2002). Another difficulty in food processing optimization arises at the characteristics of the response space of a problem. If the response space has just a single unique maximum (or minimum), finding the extremum can be guaranteed. However, food processing models are generally characterized as a system of nonlinear partial differential algebraic equations, which usually exhibit a multimodal nature (Banga et al., 2003). In this situation, an effective and systematic procedure (or algorithm) is essential for optimization. 2.3 The State of the Art of Optimization Techniques In the previous section, characteristics of food processing were discussed by focusing on model related issues. However, to solve practical optimization problems, effective techniques that are capable of consistently finding the best solution to the problems must be available. A popular optimization technique is nonlinear programming (NLP, Section 3.1.2), which is usually using gradient information to decide the search direction (e.g., steepest ascent path). If NLP is applied to highly nonlinear, constrained, and multimodal problems, it usually converges to the “nearest” local solutions, because its search direction and size are determined from its starting point (Edgar et al., 2001). Therefore, NLP cannot guarantee a global solution”, if the starting point is not close enough to the global x 5 l"3f€=l‘s to local maxima or mrnrma in a section of the entire solution space 6 reft'tl's to the highest or the lowest point among other local maxima or minima. 12 optimum. However, the local nature of NLP algorithms can be overcome by global optimization (GO, Section 3.1.2.5) techniques designed to find a global solution. GO algorithms are designed to escape from local solutions and explore promising regions where a global solution may exist. Also, GO can be utilized with ANN in parallel, without any modifications to the process model. In food process engineering, G0 has been used with numerical models, such as finite difference models (F DM), finite element models (FEM), and ANN (Chen and Ramaswamy, 2002). Banga et al. (2001) and Zorrilla et al. (2003) coupled a FDM model of double-sided cooking of meat patties directly with the Integrated Controlled Random Search for Dynamic Systems (ICRS/D8), which is an adaptive stochastic GO algorithm. Response surface methodology (RSM7) is currently the most popular optimization technique in food science, because of its comprehensive method, reasonably high efficiency, visualization, and simplicity, even though RSM is inefficient and cannot be automated in finding the overall optimum (Arteaga et al., 1994). Also, Banga et al. (2003) pointed out some important drawbacks of RSM methods, such as the empirical, local, and stationary nature of these statistical techniques. RSM has uncertainty of model equations, which means the model might not represent a real physical model, because RSM is a statistically designed experimental optimization method. Generally, the method is not considered as a formal optimization technique. However, RSM has been \ 7 RSM uses quantitative data from an appropriate experimental design to determine and simultaneously solve multivariate problems. The equations describe the effect of the test variables on the response, determine interrelationships among test variables, and represent the combined effect of all test variables in the I‘CSponse. This approach enables an experimenter to make efficient exploration of a process or system (PORCIano S. Madamba, “The response surface methodology: an application to optimize dehydration operations of selected agricultural crops”, Lebensm.-Wiss. u.-Technol., v. 35, p. 584) 13 successfully applied to product development and static process identification when the system behavior or the process is unknown, complicated, and static. 2.4 Applications in Food Process Engineering 2. 4.1 Sterilization of Canned Product Optimization techniques have been rigorously applied to sterilization of prepackaged conduction-heated food, such as retorting of canned foods. The basic form of this application is to find the best combination of retort temperature and process time for a constant heating process or the best retort time-temperature history for optimal control, which can simultaneously achieve the required lethality of the target microorganism and maximize quality factors. The first attempt by Teixeria et al. (1969) found the best combination of retort temperature and process time for a constant retort process of conduction-heated foods. The best set of time and temperature was found by plotting thiamine retention against equivalent process conditions producing the same level of lethality. Even though the attempt introduced optimization concepts to food process engineering, the study did not apply formal optimization methods to the problem. Saguy and Karel (1979) applied Pontryagin’s maximum principle (PMP) to maximize thiamine retention in retort process and found a single optimal variable retort temperature profile. The significance of this study was the first application of formal 0Inimization theory to food process engineering. However, because the application of PMP requires quite a modification of the original problem, it might not be suitable for a ComPlex problem. 14 Numerous attempts have been made since the work of Teixeira et al. (1969) to optimize retort operation. However, the essential features of those optimization problems have many similar aspects (Table 2.1). The popular objectives of sterilizing prepackaged conduction-heated foods were maximizing retention of a single nutrient, minimizing quality degradation, and minimizing energy consumption. Most of the prior studies dealt with microbial lethality as an inequality constraint and a two-dimensional numerical model (heat conduction only) as an equality constraint. Although there was similarity in the formulation of the problem among these many studies, the major difference among them was the method used to find the optimal condition or profile. Application of formal optimization methods has been increasing, so that several studies (Saguy and Karel, 1979; Nadkami and Hatton, 1985; Banga et al., 1991; Chalabi et al., 1999; Kleis and Sachs, 2000; Erdogdu, 2002; Erdogdu and Balaban, 2003) found optimal retort solutions, with respect to specific assumptions. The type of Optimal solutions can be categorized into: (a) combinations of process time and constant temperature, (b) piecewise continuous temperature profiles, and (c) on-off type of control profiles. Even though the piecewise continuous temperature profile is the true optimal solution, the solution cannot always be considered as the best, because it is not practically possible to implement a continuous profile in many conventional processes. The biggest advantage in applying process optimization to retorts is the relatively simple process model, as compared with unpackaged food processes. Usually, the retort optimization encompasses a single control variable, temperature, which limits the burden of computation. Also, the process model does not involve mass transfer, which is a very Complex phenomenon in meat cooking process. 15 Table 2.1 Analysis of basic elements used to optimize sterilization processes of prepackaged conduction-heated foods. References 0b 'ective function a Constraints Optimization Optimal ’ IC 1’ EC " method Solution ‘ Teixeira et al. Max. (thiamine 0 969) retention) VAL 2D FDM GS (T, t)opt. Teixeira et al. Max. (thiamine VAL (1975) retention) TL / TU 2D FDM GS “Tim" Saguy and Karel Max. (thiamine VAL 2D FDM PMP 5:133:33 (1979) retention) TL / TU profile (T, t)opt w/ Min. (surface & volume Comparison linear Ohlsson (1980) averaged cook value) Fc 2D FDM (diagrams) come-up- time Chart . . (Unsteady- Comparison Constant aggro et al. :glfinzg) Fc state heat (time- time & p conduction temperature) temperature ) 2D Nadkrani and . . VAL . On-off Ha tton (1985) Max. (nutrient retention) Tl. / Tu numerical PMDP control solution Ban a e t a 1 Max. (nutrient and VAL (199%) ' surface quality retention) Fc 2D FDM ICRS/DS VRT Min. (process time) F Step . c Davis, Swann function w/ Silva et al. (1992) M“ .(surface ““31”” Tend 1D FDM and Campey linear retentlon) tend method come-up- time . 2D Open-loop 815913;” et al. Max. (nutrient retention) Fc Theoretical optimal Eggglrang 2D F DM control . Max. (vitamins FC 22:03:)?“ Sachs retention) Min. ( energy Tl. / Tu 1D FEM SQP f(T,t) consumption) 1‘93! Fc Equidistant Erdogdu (2002) Max. (thiamine) Tl/ Tu mm) “mm“ ramp FDM method . Tend function Max. (thiamine of two FC Complex Piecewise Erdogdu and . . method . Balaban (2003) different objects) — TL/ Tu FDM Wei ghtin continuous multi-objective Tend method g profile IC: inequality constraints; EC: equality constraints; TL: lower temperature limit; TU: upper temperature limit. 2 Max., maximizing; Min., minimizing VAL, volume averaged lethality; Fc, critical lethality C d FDM, finite difference method; FEM, finite element method GS, graphical search; PMP, Pontryagin’s maximum principle; PMDP, Pontryagin’s minimum distributed principle; ICRS/DS; Integrated Controlled Random Search for Dynamic Search; SQP, Sequential quadratic programming e . . VRT, variable retort temperature; opt., optimum 16 Chen and Ramaswamy (2002) developed an ANN model with training and testing data produced by a finite element model, and coupled the ANN model with GA. In that study, the control function (i. e., retort temperature) was parameterized with sine and exponential functions to replace the traditional constant retort process. The coupled ANN-GA model was able to identify the relationship between the operating variables and control function parameters. Even though the ANN -GA found the optimal processing condition of a retort process, the optimal solution could be improved if additional parameters for the control function were used to increased flexibility. Also, the work of Chen and Ramaswamy (2002) was for prepackaged food retorting, a simpler process model than convection cooking of meat patties. 2. 4.2 Food Dehydration Dehydration is a common method to extend the shelf life of foods. This operation is normally removing water in a foodstuff- via evaporation or sublimation (Brennan, 1990). Among many techniques, such as freeze drying, spray drying, super-heated drying, infiared drying, microwave drying, heated air drying of a solid food block is the focus of this review, because of the physical similarities with meat patty cooking under moist air impingement cooking. Using heated air is the typical method of dehydration, in which a food product is placed in contact with a moving stream of heated air. Therefore, drying can be an optimization problem, in which the optimal conditions are sought to maximize retention of nutrients, such as ascorbic acid, while minimizing enzyme or microbial activities (Banga and Singh, 1994). Drying processes are generally driven by evaporation at the surface, which causes Water transport within a foodstuff. F ick’s equation of diffusion is often used to describe l7 water transport in drying. For hot air drying, the mass transfer model must be coupled with a heat transfer model, often assuming a one-dimensional thin slab and an evaporation term to account for the energy balance. These partial differential equations are usually solved by the finite difference method or rarely by existing analytical solutions with empirical equations for model parameters, such as diffusion coefficient, convective heat transfer coefficient, and first-order rate constant for ascorbic acid degradation (Mishkin et al., 1982). Food drying quality parameters are often modeled by using first-order reaction kinetics (Mishkin et al., 1983; Banga and Singh, 1994). Mishkin et al. (1982) used Pontryagin’s maximum principle (PMP) and the complex method to find the optimal air temperature profile maximizing ascorbic acid retention in a model system (a slab composed of water, cellulose, and ascorbic acid) with fixed relative humidity. The complex method was selected, because the method was convenient to use along with any type of process model and constraints without modification. Mishkin et al. (1983) extended their research to a multi-stage drying process. They found optimal stepwise temperature and humidity profiles for three stages by using the complex method. Banga and Singh (1994) set up four different optimization problems for drying of a thin slab of cellulose: (a) maximizing ascorbic acid retention with a constraint on the final moisture content, using air dry bulb temperature as the control variable; (b) minimizing process time with final retention of ascorbic acid, using air dry bulb temperature control; (c) maximizing ascorbic acid retention with final retention of enzyme, using dry bulb temperature and relative humidity control; (d) maximizing energy efficiency with final ascorbic acid retention, with dry bulb temperature control. These 18 problems were solved by using an ICRS/DS algorithm, and they found optimal piecewise linear control profiles in all cases. In spite of many similarities with meat cooking processes under moist air condition, drying processes are different from meat cooking, in that cooking involves more complex mass transfer phenomena (water-fat mixing and fat dripping), microbial inactivation, phase change, and air humidity sufficiently high to cause condensation. 2. 4.3 Meat Patty Cooking The abundance of prior research in retort operation is due to the availability of the process model, which encompasses relatively simple phenomena, such as depletion of certain nutrients coupled with heat transfer. However, the nature of meat patty cooking is an unpackaged process that generally involves various mass transfer phenomena, such as evaporation and dripping of fat and water. In addition, geometry change and phase transition is typical for this process. These phenomena must be coupled with heat transfer and solved. Other difficulties in the modeling of meat patty cooking arise in the heating medium. The heating medium in retort processes is water or steam, which does not interact with the food material in the package and has simple thermal properties. However, in an unpackaged food product, the heating medium interacts with the surface of the food material, which therefore leads to more elaborate boundary conditions. Banga et al. (2001) and Zorrilla et al. (2003) applied the dynamic optimization technique (Section 3.1.2.6) to contact cooking of meat patties, which is considered as the first attempt to optimize a meat patty cooking process. The objective of those studies was to minimize cooking loss of patties, while ensuring inactivation of E. coli 0157:H7 and final product center temperature. One-dimensional coupled heat and mass transfer solved l9 with finite difference method (Pan, 1998) was used as a process model. The above optimization problem was solved by a global optimization algorithm, ICRS/DS (Integrated Controlled Random Search for Dynamic Systems), which found the optimal piecewise step grill surface temperature profile (Banga et al., 2001; Zorrilla et al., 2003). The process model is an important part of an optimization problem, because the model is used to predict the objective value and constraints. The process model developed by Pan, (1998) and used by Banga et al. (2001) and Zorrilla et al. (2003), assumed the patty as a one dimensional infinite slab, which is less accurate than two- dimensional modeling. Even though the model could predict water and fat transfer and microbial log reduction, other important quality factors, such as internal color change and surface color change, were not considered. Also, the nature of the cooking method was contact cooking, which was modeled using simplified, effective boundary conditions, compared with the condensing-convective boundary conditions during meat patty cooking under a moist air environment. In their optimization problem, grill surface temperature was the single decision variable. However, moist-air impingement meat patty cooking involves the additional control variables of air humidity, impingement velocity, and impingement geometry. Even though the performance of ICRS/DS was good enough to locate the global optimum for the contact cooking problem, comparison with other GO methods was not conducted. 2. 4.4 Various Processing Areas In addition to the application for canning, cooking, and drying, optimization techniques have been applied to the other processes, such as ultra-filtration, baking, extrusion, cheese manufacturing, mixing, and so forth. Each application provides some 20 lesson about the implementation of formal optimization techniques to the unique nature of various food processes. Usually, the objectives of wine filtration are to minimize colloid content, maximize color intensity, and maximize flux by varying pore size and the recycle rate of membranes (Gergely et al., 2003). Gergely et al. (2003) expressed the objective firnction as a second-order form of regression functions of the membrane pore size and recycle flow rate, which was a response surface model (RSM). Even though the RSM can solve some Optimization problems, generally the method is not recognized as a formal optimization method (Section 2.4). In optimizing commercial bread baking, the biggest challenge is getting a reliable process model. Therdthai et al. (2002) used four different temperature zones of a commercial oven and baking time to find optimal condition for minimizing loss while controlling the top crust color, side crust color, and average crust color within acceptable ranges. Statistical methods were used to construct a model equation, which was a RSM. However, some researchers have used neural networks coupled with a Genetic Algorithm (GA) for leavening process optimization in a bread-making industrial plant (Fravolini et al., 2003). Fravolini et al. (2003) used a nonlinear system identification method, called NARMA (Nonlinear Autoregressive-Moving Average), to model the leavening process. Extrusion is very complex process. The most common objectives are expansion ratio, shearing strength, and sensory texture, which are functions of temperature, feed moisture, and process variables, such as screw speed, screw compression ratio, feed speed, and die diameter. RSM is the most popular method to model the process and to apply optimization methods, because the results of RSM are analytically differentiable 21 mathematical expressions (Chévez-Jéuregui et al., 2000); Frazier et al., 1983; Iwe et al., 1998; Karwe and Godavarti, 1997; Olkku et al., 1983; Quintero-Ramos et al., 1998; Vainionpaa, 1991). Therefore, optimal condition could be found by applying theoretical optimization techniques directly to the model. In addition to applications for food processes, well-organized dynamic optimization problems also exist in biochemical processes, such as fermentation (Tartakovsky et al., 1995; Banga et al., 1997; Berber et al., 1998; Tsoneva et al., 1998; Faqir, 1998);Lee et al. 1999; Radhakrishnan et al., 1999; Halsall-Whitney et al., 2003; Levisauskas et al., 2003). The nature of biochemical processes is similar to food processes, in that their dynamic behavior is inherently nonlinear. However, most biochemical systems can be modeled by a set of ordinary differential equations, which makes many optimization theories applicable, because the model equations are differentiable. Generally, biochemical processes have been treated as problems of optimal control, which use a special form of the performance function. 22 3 THEORIES 3.1 Overview of Optimization Theory 3.1.] Introduction 3. 1 .1.1 What is Optimization? Our daily life is always full of choices in various activities, such as traveling, shopping, or eating. Decision—making is the cognitive process of selecting a course of action from among multiple alternatives. The purpose of decision-making is to choose the best alternatives fitting with our goals. Even though the purpose of decision-making is clear, making decisions involves many considerations and uncertainties. In the case of a simple problem, we can guess the results of some trials or test the effect of each possible alternative. However, as the number of variables and the interactions between variables increase, these attempts lose the ability of identifying the best solution among the many possible good solutions. For example, we can increase the thickness of insulation to decrease energy loss, but increased insulation thickness increases cost, which is the trade- off. In that case, the problem is to find the thickness of insulation minimizing the total cost, which cannot be solved easily. Therefore, systematic methods and procedures are necessary to solve those problems containing many variables, restrictions, and factors, which compete with each other for the best solution. Optimization can be defined as “the collective process of finding the set of conditions required to achieve the best result from a given situation” (Beveridge and Schechter, 1970). By the help of optimization techniques, we can explore more complex decision-making situations with more certainty and effectiveness. 23 3.1.1.2 Essential Features of Optimization Problems By observing optimization problems, some common features can be found. For example, if a driver wants to travel from city “A” to city “B” in the least time, then the minimum traveling time will be the primary objective for the driver. To accomplish this goal, the driver has to consider the speed limit of each state, weather conditions, physical limits of driver and car, paths, and so on. The basic components of an optimization problem are the objective function, decision variables, constraints, and mathematical model. These components are well explained by Evans (1982). 1. Objective function (performance function, or cost function): This is the quantity to be maximized (or minimized). It is often referred to as the cost or performance function. Whether measured in dollars, efficiency, or other terms, the performance function evaluates alternative solutions to the problem to determine which one is the best. Decision variables: These are the parameters in the process or system that can be adjusted to improve the objective. They are the free or independent variables that must be specified in the traditional case-study approach to problem-solving. Constraints: All optimization problem have constraints on the allowed solutions. These constraints may limit values of the decision variables or of other dependent variables that describe the behavior of the system. Mathematical Model: If the problem is to be solved other than by trial-and-error physical experimentation, we must have a model of the system. The model is the mathematical representation of the system that determines the objective function 24 in terms of the decision variables. It also determines other dependent variables that may be subject to constraints. Now, we can define optimization by using the terms described above. Therefore, optimization is finding a set of conditions maximizing or minimizing the objective function while satisfying the imposed constraints of the problem. 3.1.1.3 Solving Optimization Problems To solve an optimization problem, one must analyze the essential features of problem. This means that the objective and constraints of the problem should be identified first, but not necessarily in the form of mathematical expressions. Then the objective function must be expressed in terms of the system variables by using mathematical expressions. Interrelationships, internal restrictions, and system models must be expressed in mathematical form. Now, the problem is reconstructed according to the essential components of a formal optimization problem. Upon this mathematically redefined optimization problem, an appropriate optimization method and algorithm can be applied to obtain the optimal conditions to achieve the objective. 3 . 1.1.4 Hierarchy of Optimization Problems Optimization can be employed at any level in a plant, ranging from a small piece of equipment to management of a whole company. Someone may need more workers to maximize the production rate. However, from the management perspective, this 0Dtirnization result might not be favorable, because of labor cost. Therefore, the scope of Optimization is important, because one level of optimization does not guarantee the 25 optimality of another level of a system. Typical industrial company optimization encompasses three levels: management, process design and equipment specification, and plant operation (Edgar et al., 2001). The highest level of optimization of a company is the management to maximize net profit. To accomplish the objective, someone would need to model the whole company, including every piece of equipment, which is practically impossible. Therefore, the scope of an optimization problem must be reviewed before proceeding to the next step. 3.1.2 Theory and Methods 3.1.2.1 Basic Concepts of Optimization The basic concepts of optimization can be clearly illustrated by observing an unconstrained one-dimensional case. “Unconstrained” means that there are no equality or inequality constraints, so that the independent or dependent variables can be simulated without any limitations. Let’s take a simple arithmetic example. If we define a quadratic firnction y with an independent variable x and coefficients a, b, and c, then: y=f(x)=ax2+bx+c [3,1] The objective is to find x minimizing or maximizing the function value. One might plot the function to see the actual shape of the curve and find a maximum or minimum point graphically (Figure 3.1). However, there is a mathematical tool to see the feature of the Curve, which is the first derivative 0») information obtained by differentiating the above equation with respect to the variable x. This is shown in Equation [3.2]. y'=§%r—)=2ax+b [3.21 26 f(X) ‘—y=f(x)=ax2+bx+c df(x)/dx “ ‘ 0 -b/2a Figure 3.1 Locating minimum point of a quadratic function with first derivative information. 27 This first derivative information implies that there is an extreme point between the sign changes of the gradient (Equation [3.2]), which is an inflection point. The inflection point (stationary point) of zero gradient can be calculated by setting the right hand side of Equation [3.2] equal to zero. Thus, the stationary point is x= —b/2a. If the x value is inserted to the Equation [3.1], then we can have an extreme value of y = (—b2 + 4ac) /(4a). However, we still do not know whether the value is maximum or minimum if we do not have graphical information. To determine this, the second derivative information of the Equation [3.1] is necessary. The second derivative (y") is: y" = 2a [3.3] If the second derivative is negative, then the extreme value is a maximum or vice versa. Hence, a necessary condition for a minimum or maximum of f(x) is that the first derivative of f(x) becomes zero at x. However, the necessary condition does not tell whether the extremum is the minimum or maximum. So, we need the second derivative information as a sufficiency condition. Now, we can identify a minimum or maximum of f(x) mathematically with the necessary and sufficient conditions. The above case is very simple, but the basic concepts can be extended to the multivariable cases with some help of mathematical techniques. Let f(xl, x2,...,x,,) be an n- dimensional function. To meet the necessary condition, the first derivative of each Variable x1, x2...x,. must be zero as follow. af(x)=6f(x)=,,,=i(x_)=0 [34] 6x1 6x2 a'xn . 28 A set of variables satisfying the set of Equations [3.4] represents an extreme point. The sufficient condition for the extreme point can be checked with the nature of the matrix of second partial derivatives (Hessian matrix) of f(xl, x2,...,x,.) at the point, which is the extension of Equation [3.3] to the multivariable case. If all the objective functions of our problems were quadratic and differentiable, then life would be simple. However, this analytical method has many limitations in real life applications. Therefore, we need a more general approach to solve optimization problems. Before addressing this issue, imagining a movie scene will give us insight for the generalized optimization method. Let’s imagine a situation that two special soldiers (captain “A” and captain “B”) are dropped from an airplane into an enemy region in a night having no moon light (Figure 3.2). Their mission is to find the highest point from sea level and communicate with the nearby resistance. Their landing location is different from each other, so that they have to complete the mission individually. Captain A has the map of the region, compass, and GPS (global positioning system), but captain B lost his equipment because of a tough landing. Captain A saw the map and located the highest point and also his current location with GPS. Then he set the direction toward the highest point and finally got to the point. However, captain B does not have any information or equipment, except a small flashlight in the complete darkness, which means he only can get terrain features a few yards around him. How could he get to the highest point? So, he marked his current location on the ground and looked around with the flashlight. He found that the steepness of terrain was increasing in the northeast direction. So, he kept moving a few yards from his current location to the next location in this direction if there Was an increase of elevation. He repeated this strategy until he could observe no more 29 increase of elevation. Fortunately, he reached the same location where captain A arrived, because there was only one peak in the region. However, if there were several peaks, it could take more time or he might be stranded at a sub-peak. /: I : : 1: i \ ' f I ' j 3 ' 7 ' Sea Level Sea Level Figure 3.2 Illustration of the operation of captain “A” and captain “B”. In the above story, our focus is captain B’s approach, because his situation is very Slmilar to ours, like finding a maximum or minimum within a region that we cannot 30 figure out, which means no analytically differentiable mathematical expressions. Captain B’s approach gives us a clue for setting up a general rule for finding such an extreme point in an unknown region. There are three factors to determine the next point from the current point. Current point (15"), search direction (of), and step size (5") are needed to determine the next point (xw). Also, the relationship can be expressed mathematically as following: xk“ =xk +a" -S" [3.5] The above iterative search procedure stops based on some criteria. If there is no significant improvement, the search stops. There are many methods to perform the search; however, the differences between methods are mainly in how they generate the search direction. A popular method is using the first derivative information. However, function value, and finite difference approximation are also used in lieu of derivatives (Edgar et al., 2001). Table 3.1 shows some methods according to their approaches. Table 3.1 Unconstrained multivariable optimization methods categorized by how they generate the search direction (Edgar et al., 2001;Venkataraman, 2002). Function values only Gradient based @irect Method) (Indirect Method) Random Search Steepest Descent Grid Search Conjugate Gradient Univariate Search Davidon-Fletcher—Powell Method (DF P) Simplex Search Broydon-Fletcher-Goldfarb-Shanno Method (BFGS) Conjugate Search Until now, we discussed handling unconstrained function optimization. However, the real life optimization problem always has equality or inequality constraints. The Constraints have to be handled so that the constrained optimization problem can be Converted into an unconstrained optimization problem. The basic idea is to set the - 31 objective function free from constraints by using mathematical manipulations of the constraints. Let us first look at the handling of equality constraints. Three methods are mainly used for the solution of such problems: direct substitution, constrained variation, and Lagrange multipliers. Direct substitution is to substitute equality constraints to the objective function when the equality constraints are all linear equations. Then, the constraints vanish, and the problem can be handled as an unconstrained case. Constrained optimization subjected to inequality constraints is treated as in the case of equality constrained problem after transformation of inequality constraints to equality ones by introducing a slack variable, which is a buffer between the original inequality constraint and the transformed equality constraint. Therefore, a general approach to solve a constrained optimization problem is to convert the constrained problem into an unconstrained problem. Then, unconstrained multivariable optimization techniques (Table 3.1) can be used to find optimal solution. 3.1.2.2 Linear Programming A linear programming8 (LP) problem is one in which the objective and all of the constraints are linear functions of the decision variables, so that the linear constraints (lines in 2-dimensional case) form boundaries. Therefore, an objective function, represented as a line, is moving through the bounded region to find an optimal set of variables to get an extreme value in a two-dimensional case (Figure 3.3). 8 “The word programming here does not refer to computer programming, but means optimization.” Edgar, T. F ., D. M. Himmelblau and L. S. Lasdon (2001). Optimization of chemical processes. New York, McGraw-Hill. 32 X2 (Constraint) if Maximum \ ’ Constraints Feasible region ' “\\\ \\ f (Objective Function) §\\\ \\ ‘ £7 X1 (Constraint) Constraints Figure 3.3 Concept of linear programming problem in a two-dimensional case. One of the important characteristics of LP is that the extremum of a linear program always occurs at a vertex or corner of the system boundaries, which is the intersection of constraints in the feasible region (Beveridge and Schechter, 1970). The basic idea is a systematic examination of these boundaries, which is converting a set of constraints into a set of equality equations and applying linear algebra and matrix manipulation (Venkataraman, 2002). For instance, the simplex algorithm is designed to explore one intersection after another to the direction of improving the objective function, according to a set of rules, until the best objective function attainable is found (Beveridge and Schechter, 1970). 33 3.1.2.3 Nonlinear Programming Nonlinear programming (NLP) problems contain at least one nonlinear equation of the objective function or constraints (Venkataraman, 2002). No single optimization algorithm can possibly be efficient or even successful in all cases of interest. However, unique techniques are well developed for each specific case. For instance, a nonlinear objective function with linear equality constraints can be handled with the direct substitution method, which solves the objective function explicitly for one variable and eliminates that variable from the problem formulation (Edgar et al., 2001). A nonlinear objective function with linear inequality constraints can be solved by using Kuhn-Tucker conditions9 and Lagrange multipliers"). If an optimization problem has a quadratic objective function and linear inequality or equality constraints, quadratic programming (QP) can be used (Edgar et al., 2001). Another strategy to solve nonlinear optimization problems is to replace all nonlinear functions in the problem with their Taylor series approximations and apply linear programming, which is called successive linear programming (SLP) (Edgar et al., 2001). The most robust and generally accepted nonlinear optimization technique is the generalized reduced gradient (GRG) method, which is also implemented as the “Solver” in the spreadsheet program Microsoft Excel (Edgar et al., 2001). The concept of GRG is to reduce the dimension of the first derivative of the objective firnction by using the first 9 At any local constrained optimum, no (small) allowable change in the problem variables can improve the value of the objective function (Edgar, T. F., D. M. Himmelblau and L. S. Lasdon (2001). thimization of chemical proces_s_e§. New York, McGraw-Hill.). This is a generalization of the Lagrange multiplier method. '0 By introducing an unknown scalar variable to the constraints, a linear combination is formed, which reduces a constrained problem into an unconstrained problem (Venkataraman, P. (2002). Applied optimization with matlab proggamming. New York, John Wiley & Sons). 34 derivative of constraints. All major NLP algorithms are based on estimation of first derivatives of the problem to obtain a solution and to evaluate the optimality conditions (necessary and sufficient conditions) (Edgar et al., 2001). Because NLP is mainly based on the gradient information, the solution has a local nature, which means the solution is a local maximum or minimum, but not necessarily the global one in general cases of nonlinear optimization problems (Figure 3.4). Global maximum Local maximum , fl" l’ I, ’0’: g::: Figure 3.4 Global and local optimum in 2-dimensiona1 case (The picture was adopted from the help manual of MATLAB®). 35 3.1.2.4 Mixed-Integer Programming Another special type of optimization problem is mixed integer programming (MIP). As the name implies, this type of problem includes integer variables that are not continuous. The integer variable is common in plant operation, design, location, and scheduling (Edgar et al., 2001). For example, number of equipment, yes-no decisions, and number of stages are integer variables. If the objective function and constraints are linear in a MIP, then it is called mixed-integer linear programming (MILP). And if the MIP involves nonlinear objective and constraints, then it is called mixed-integer nonlinear programming (MINLP) (Edgar et al., 2001). If the number of integer variables is small, then exhaustive or complete enumeration is possible. However, the effort grows exponentially to examine all possible solutions (V enkataraman, 2002). Therefore, systematic ways are necessary to get to the optimal solution with less effort and time and also to handle large scale problems. One popular solution for this kind of the problem is branch and bound (BB). The basic concept of BB is systematically finding the closest set of discrete variables from the optimum set of continuous variables. To do this, a relaxation technique is used to convert the discrete variables into continuous variables bounded by their maximum and minimum value. BB uses two strategies to find the best set of integer variables. Branching is the efficient way of covering the feasible region by several smaller feasible sub-regions, and bounding is comparing and selecting the sub-region having its upper bound that is less than the lower bound of any other sub-region. 36 3.1.2.5 Global Optimization We have discussed unconstrained, constrained multi-dimensional cases, and some special cases, such as LP, NLP, MILP, and MINLP. One of the underlying assumptions for optimality is that the problem must be convex, which means there is only one minimum or maximum. Thus, if the convex condition is not met, then the above methods can be stranded at a local optimum, not at a global optimum (Figure 3.4). This situation frequently happens when the model has nonlinear equality constraints, such as a nonlinear material balance, nonlinear physical property relations, nonlinear process models, and so on (Edgar et al., 2001). In addition, if a gradient-based method is involved, then the search method can become stranded around the vicinity of local minima or maxima, without guarantee of optimality. However, in spite of this fundamental limitation, there are several effective and practical optimization techniques that explore infinite solution spaces systematically with minimum effort and high possibility of locating the global optimum. This type of optimization is called global optimization (GO). In a G0 problem, a local solution can be considered as a global optimum if the local solution is the best among many local solutions or if multiple search trials from different starting points reach to the same local solution. Widely used GO methods can be classified as deterministic (exact) or stochastic strategies (heuristic) (Banga et al., 2003). Deterministic global optimization methods are based on the systematic gradient- based local search methods following the systematically divided sub-regions of attraction until they meet their termination criteria. This type of method can guarantee global 37 optimum in certain problems, but not in a general GO problem. Branch-and-boundll methods, methods based on interval arithmetic (Kearfott, 1996), and multi-start methods are classified as deterministic GO techniques (Edgar et al., 2001). Even though deterministic GO methods have sound theoretical convergence properties, the associated computational effort increases very rapidly (often exponentially) with the problem size (Banga et al., 2003). Stochastic GO methods do not follow systematically divided regions of attraction. Instead, they set their search direction based on a logic found in natural processes, such as cooling of metals and genetic evolution processes. Many stochastic GO methods can locate the vicinity of global solutions with relatively good efficiency, but the downside is that global optimality cannot be guaranteed (Banga et al., 2003). Scatter search, tabu search, simulated annealing, and genetic and evolutionary methods can be classified as stochastic GO methods (Edgar et al., 2001). Stochastic GO methods are applicable to almost any problem, without modification of the original process model (Edgar et al., 2001). In the case of Genetic and evolutionary GO methods, the method produces a population that is a set of solutions and keeps updating the population with improved solutions according to the rule of biological processes of crossover and mutation. Although stochastic GO cannot guarantee a global optimum, it is widely adapted to solve real life problems, because of its relative simplicity and robustness. " Branch and bound algorithms are a variety of adaptive partition strategies that have been proposed to solve global optimization models. These are based upon partition, sampling, and subsequent lower and upper bounding procedures: these operations are applied iteratively to the collection of active ('candidate') subsets within the feasible set D. Their exhaustive search feature is guaranteed in similar spirit to the analogous integer linear programming methodology (Eric W. Weisstein et al. "Branch and Bound Algorithm." From Math World--A Wolfram Web Resource. http://mathworld.wolfram.com/BranchandBoundAlgorithm.html, March 22, 2005). 38 3.1.2.6 Static and Dynamic Optimization In the previous sections, we have discussed how to solve optimization problems in each special case, especially when all the constraints are algebraic equations. In other words, the problem is not changing with time. For instance, finding the optimum diameter of a pressure vessel does not involve a time factor. However, what if the constraints have time varying equations, such as ordinary differential equations in which a variable is changing with respect to time? This kind of optimization problem can be called “dynamic” optimization, as opposed to the “static” optimization problem. To be precise, this type of problem is generally referred to as an optimal-control problem, which is the area of control. However, some optimization techniques can be employed to solve such problems. There are two general approaches. The first approach is discretization of the control function, which can be understood as dividing the control function into pieces. This method is replacing a differential term by the first order Eulerian difference expression, which is: g z 36(0) - x(tl—l) dt At [3.6] Now the ordinary differential equation can be expressed as a set of algebraic equality constraints, which is the standard form of a constrained nonlinear programming technique. The second method is parameterization of the control function. We need an infinite number of points to express a varying control function, which is not practical. However, if the control filnction is expressed with Fourier series or a simple polynomial 39 form, the infinite solution space can be explored with a finite number of parameters. By applying NLP, the optimal set of parameters can be obtained as the optimal control function, even though there is no guarantee of optimality. 3.2 Global Optimization Techniques 3. 2. 1 Genetic Algorithm (GA) Genetic algorithms (GA) are well-known stochastic global optimization algorithms based on biological evolution theory. Holland (1962) was the first to use the technique, but its use as an optimization tool began in earnest in the late 19805, developed momentum in the mid-19908, and continues to attract serious interest today (Venkataraman, 2002). Genetic algorithms have been successfully applied to a wide range of problems, such as engineering design, scheduling, signal processing, optimal control, transportation, and so on. The biological evolution theory is the combination of Darwin's theory of natural selection and Mendel's theory of genetics. According to the theory, a chromosome in a gene pool is modified by simple rules of genetics, such as crossover and mutation. The biological system having the gene that is the fittest among others to the environment survives and produces the next generation, having individuals with more desirable characteristics. These biological concepts are implemented in genetic algorithms via numerical operations. For example, a vector of design variables is considered as a chromosome in genetic algorithms. The following is a description of genetic algorithm datails; the explanations are largely adopted from the work of Venkataraman (2002). Figure 3.5 shows the general procedures of a genetic algorithm. The initial population is constructed with a certain number of design vectors that take higher rank in 40 terms of a performance index or an objective function value, in the case of an unconstrained problem among randomly produced individuals”. For the constrained problem, these individuals must satisfy the specific constraints. The initial population size usually remains the same throughout the whole generation. However, there are no specific criteria to determine the optimal initial population size (Venkataraman, 2002). Next generation candidates are now produced by using selected parents and genetic operators. A crossover genetic operator exchanges a piece of the chromosome of each parent. A simple crossover exchanges a piece of chromosome in the same location of each parent, and the ratio of the piece is generated randomly at each operation (Figure 3.6). Arithmetic crossover uses a linear combination of parent chromosomes to produce two children. '2 This refers to a vector of design variables. A piece of the chromosome or a portion of a design vector is called as allele. 41 Set up an initial population Evaluation & Parents selection t Perform Genetic Operations Crossover Mutation Immigration Reevaluation of the population (Fitness) Perform competitive selection < New generation. > Check if generation No reached to a target # of generations Yes C Terminate D Figure 3.5 General procedure of a genetic algorithm. 42 ------- * i 3" =iii¥iliifiii : g x5 i i 3’5 i i X5 l l 3’5 g gysi in: gyvé in: E xs i E Y8 g i X3 3 i YB 3 Parents Children Figure 3.6 An example of simple crossover operation in a genetic algorithm. The children are defined as: C1 = XIX-I- le [3.7] C2 = 12X + lilY where, Ii] +12 =1; 11,12 > 0 X, Y: parents C1, C2: children The parameter A, and A; are generated randomly. Mutation selects a design variable randomly and replaces it with a randomly generated value. Generally, these genetic operators, crossover and mutation, narrow a search to a promising region that might have a local extremum (Venkataraman, 2002). However, it is hard for a global optimization algorithm to find a global optimum if it is stranded around a local solution. Another 43 feature of GA is immigration, in which a set of randomly generated, unbiased population is added to the existing population before the selection is made for the next generation. By using this feature, GA can reduce the possibility of being trapped in a local solution and effectively explore an entire solution space. The newly generated population now must be evaluated to identify suitable parents for the next generation. A simple method of selecting parents is ranking or sorting the individuals of the population according to their objective function value. A fraction of the best individuals can be used for reproduction of the next generation. In a strategy called tournament selection, the remaining portion of the population, excluding the fraction of the best individuals, is fed back to the next generation as new immigrants. The above genetic operation and fitness evaluation are performed until the generation number reaches a certain number or the solution does not show any improvement. GA is very useful for handling ill-behaved, discontinuous, and nondifferentiable problems, because the algorithm generates a possible solution group at each iteration, instead of a search direction, and is effective for handling continuous problem (Venkataraman, 2002). 3.2.2 Simulated Annealing (SA) Annealing is a heat treatment technique to alter material characteristic by removing internal stresses and crystal defects. When a material is in molten state, atoms can move freely before they form crystal structure. However, if a molten material is cooled down rapidly, atoms lose chances to form crystal structure, which creates internal stresses and crystal defects. Therefore, the rapidly cooled material is in a higher energy 44 state than the material cooled slowly. The cooling schedule is also called the annealing schedule, which reduces the surrounding temperature from the critical temperature to a sufficiently lower temperature stepoby-step. In each cooling step, the material being treated is allowed to cool in the firmace until it reaches thermal equilibrium with its surrounding temperature. By using the technique, internal energy of a material can be minimized, which means stable crystalline order. Usually, annealing is used to soften a material and make it more ductile, to relieve residual stresses, and to refine the crystal structure (Shigley and Mischke, 1989). The above ideas can be applied to finding the global optimum in an optimization problem, which has become popular for combinatorial optimization” problems. The molten state of the material can be an initial design space in which every combination of variables has equal opportunity to be searched out to find the best combination. The objective is to reduce the internal energy level to the lowest state, which is crystalline order. As the annealing temperature goes down, the mobility of atoms decreases, and the total internal energy level goes down at the same time, which means that the search direction is being fixed to a direction without considering other opportunities. When the temperature reaches the final scheduled temperature, atoms form crystalline structure, which represents the lowest energy level, a global minimum. The above analogy is realized by the basic SA procedure (Floquet et al., 1994) summarized by Edgar et al. (2001) as follows: 0 Choose an initial solution x , an initial temperature T, a lower limit of temperature T Low, and an inner iteration limit L. '3 A method to search for the best possible solution out of very large number of discrete feasible solutions. 45 0 While (DTLow), ('10 O O O Fork= 1, 2,..., L, do Make a random choice of an element x' e N (56) , where N( x) = x + a - S (or: search direction; S: stepsize) Move_value= f (x') - f (x) If move_value S 0 (downhill move), set x = x' If move_value > 0 (uphill move), set x = x' with probability p = exp(- move_value / T) > r (uniformly distributed random number in [0, 1]) (Kirkpatrick et al. , 1983) End inner loop 0 Reduce temperature according to an annealing schedule. An example is new T =C’T 0, where 0 < c < 1 (c is the rate of annealing schedule). 0 End temperature loop Basically, the algorithm is more supportive of a solution that improves the objective function, while permitting adverse solutions to give potential for the algorithm to discover the global optimum and to escape a local optimum (Venkataraman, 2002). Some of the critical parameters are the perturbation method of neighbor selection, transitional probability,14 and the rate of annealing schedule (c). To obtain a neighbor state of X0, a search direction and a stepsize must be determined. Given a neighborhood structure, simulated annealing can be viewed as an algorithm that continuously attempts to transform the current configuration into one of its " Probability to accept a worse solution. 46 neighbors. Practically, the calculation of the search direction and stepsize could be determined by using the traditional one-dimensional case (Venkataraman, 2002). For the transitional probability firnction, Metropolis algorithm15 and Glauber algorithm16 are related to Boltzman probability distribution17 (Edgar et al., 2001). The way in which the temperature is decreased is known as the cooling schedule. The rate of annealing schedule must be appropriate for its application so as not to get trapped in a local minimum due to fast cooling. The annealing rate can be fixed to a value or can vary in each step of cooling (Laarhoven and Aarts, 1987). For the algorithm to be effective, it is recommended that the probability be in the range of 0.5 S p S 0.9 (Venkataraman, 2002). This mechanism is mathematically best described by means of a Markov chain”, a sequence of trial, where the outcome of each trial depends only on the outcome of the previous one (Laarhoven and Aarts, 1987). An example of a Markov chain is a random walk"). A Simulated Annealing program consists of a pair of nested DO-loops. The outer- most loop sets the temperature, and the inner-most loop runs a Metropolis Monte Carlo simulation at that temperature. '5 p = exp (-move_value / T) MP = exp I-move_value / T) / (l + exp (-move_value / T) ‘7 p = -k / T (k: Boltzmann constant; T: annealing temperature) ‘8 A collection of random variables {X,} (where the index t runs through 0, 1,...) having the property that, given the present, the firture is conditionally independent of the past. In other words, P( X, = j l X0 = to, X] = i1, Xr) = i)- 1) = P( X, = jl X,- I = i,-1). (Eric W. Weisstein. "Markov Chain." From Math World--A Wolfram Web Resource. http://mathworld.wolfram.com/MarkovChain.html, May 18, 2005) 19 . . . . . Idea of taking successrve steps in a random direction. 47 Like other discrete optimization techniques, such as Branch-and-Bound and Tabu Search, SA is a popular method to solve combinatorial optimization problems, which seek the best combination out of many possible combinations. Typically, these techniques are being used in production system planning, scheduling, transportation, and logistics. 3.2.3 Integrated Controlled Random Search for Dynamic Systems (I CRS/DS) ICRS/DS is a modification of the ICRS algorithm (Banga and Casares Long, 1987; Casares and Rodriguez, 1989), which was a generalization of the method proposed by Goulcher and Casares Long (1987) for steady-state optimization problems (Banga et al., 1997). As a generic algorithm, ICRS/DS is not as popular as GA or SA, but it has been applied to various problems, such as bioprocesses (Banga 'et al., 1997), wastewater treatment (Banga and Casares Long, 1987), retorting (Banga et al., 1991), drying (Banga and Singh, 1994), and meat patty cooking (Banga et al., 2001), and proved its ability to solve dynamic Optimization problems. Basically, ICRS/D8 uses two strategies: control vector parameterization and a stochastic direct search procedure. The original constrained dynamic optimization problem is transformed into a constrained nonlinear programming (NLP) problem by using a flexible parameterization of the control function. The constrained NLP problem is solved using the stochastic direct search procedure (Banga et al., 1997). Like GA and SA, the search uses search direction and step size to determine the next feasible move. This search mechanism starts with a user-specified feasible control vector and perturbs the vector randomly with a normal probability distribution, which has the control vector as the average and a standard deviation vector. The standard deviation vector is a set of 48 smaller distances between the control values and their upper and lower bounds. If the algorithm cannot find any improvement within a limited number of trials, the standard deviation vector is reduced by a heuristic parameter, which results in smaller search step. The algorithm can be terminated by checking user-specified tolerances for the decision vector and/or the performance index at the end of each iteration (Banga et al., 1997). The detailed procedures of ICRS/DS can be found in the work of Banga et al. (1997). These procedures consist of control profile parameterization and the modified ICRS algorithm. Control profile parameterization: The control function u(t) over t e [to ,tf] is parameterized using N points, (6,- , (0)) (i = 1...N). The value of u(t) at iteration k can be calculated using variable-length piecewise linear interpolation (Equation [3.8]) within the optimization procedure and Q-k _<. t 3 61/11; k k (t)- - 0)- u" (t) = 3124—“: - 91‘) + of [3.8] By using the above parameterization technique, the original optimal control problem is transformed into a 2N (or 2NM if there are M control variables) dimensional nonlinear programming (NLP) problem. For implementation purpose, the parameterized time and control values (0,- , 60,-) are replaced with the decision variable vector 4‘ for any iteration k by using the following rule: 49 5126,", i=1....N k k [3.9] g, =ro,_N, i=N+l...2N Also, the upper and lower bounds for the decision variables are expressed in a simple way: it” = 0r”, A.3-L = 6}; i=1...N U U L L [3.10] .- = 0,-_N. e“.- = 0,-_N; i: N+ l...2N In addition, the parameterization must satisfy following conditions. For time-intervals: ref 3 l9," s a,” ‘ k k . l9,- s cm r1 =1...N, Vk [3.11] L U L 92' s e,- s 9.41 I For control vector limits: (0,-1’ 3 ref s (in , i=1...N, Vk [3.12] For initial and final time-interval limit: L U U 61 =61 =t0, 6N =tf [3.13] The above parameterization method is designed for variable-length intervals. However, if a fixed constant discretization ratio of time is desired, it suffices to take a,” = of + e =[ tf )(i—1)+ e [3.14] where e is a small number. 50 The Modified ICRS Algorithm (Figure 3.7): The algorithm starts with a feasible decision vector, which can be chosen randomly. Subsequently, a standard deviation vector is calculated by multiplying minimum intervals (smaller interval between the decision vector and upper and lower bound vector) by a heuristic parameter (k,). Based on the initial feasible vector, the algorithm generates a next decision vector by using the standard deviation vector (stepsize) and Gaussian distribution (random direction). Once the next decision vector is evaluated as feasible (satisfying all the constraints) and improved (increasing or decreasing objective function value), the above iteration is repeated until a convergence criteria is satisfied. However, if the next generation decision vector is determined as an infeasible or retrogressive (or remains the same), the failure counter (F) is increased by one up to a pre-set value (recommending ne >< total dimension of a problem). Every failure counter increment, another new decision vector is generated with the same standard deviation vector and tested. If the algorithm fails to get an improved decision vector within a pre-set failure counter value, the standard deviation vector (stepsize) is decreased by multiplying another heuristic parameter (k2). This stepsize reduction is continued until the algorithm encounters a feasible and improved decision vector with the frequency of the maximum failure counter value. Once the feasible and improved decision vector is found, the failure counter is set to zero. In this algorithm, the three heuristic parameters (k1, kg, and n.) are critical to be successful in optimization process. Pan (1998) recommended 1/3, 1/2, and 4 as default value for k), kg, and ne, respectively. 51 Initial feasible decision vector 60 Evaluate objective function: 104a") Calculate minimum intervals of decision vector 6*: if? = mink,” — rifle," - g} )j i=1...2N Calculate standard deviation vector: of = rel/if, i=1...2N Generate new decision vector f H: ."+‘ =g," + 0,-1",, i=1...2N Increase failure counter: Evaluate objective fiinction: F=F+1 Me) No @ Yes Yes Reset F =0 of e—kzoi", i=1...2N Yes Convergence criterion? igk+1—g."l k-i-l _ Jki ———'