Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Numerical investigation of industrial continuous digesters Fan, Yaoguo 2005

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
831-ubc_2005-104880.pdf [ 11.7MB ]
Metadata
JSON: 831-1.0080759.json
JSON-LD: 831-1.0080759-ld.json
RDF/XML (Pretty): 831-1.0080759-rdf.xml
RDF/JSON: 831-1.0080759-rdf.json
Turtle: 831-1.0080759-turtle.txt
N-Triples: 831-1.0080759-rdf-ntriples.txt
Original Record: 831-1.0080759-source.json
Full Text
831-1.0080759-fulltext.txt
Citation
831-1.0080759.ris

Full Text

NUMERICAL INVESTIGATION OF INDUSTRIAL CONTINUOUS DIGESTERS By Yaoguo Fan B.Sc, Huazhong University of Science and Technology, China, 1994 M.Sc, Huazhong University of Science and Technology, China, 1997 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIRMENT FOR THE DEGREE OF  DOCTOR OF PHILOSOPHY in THE FACULTY OF GRADUATE STUDIES (Mechanical Engineering)  THE UNIVERSITY OF BRITISH COLUMBIA September 2005  © Yaoguo Fan, 2005  Abstract Digesters, either of batch or continuous type, are large chemical cooking vessels used to convert wood chips into fiber solution or pulp, from which paper products are eventually produced. In the continuous digesters, the complex liquor flow, which usually has a number of circulation loops, interacts with the downward moving wood chips through mass, momentum, and energy exchanges and therefore has a major effect on the chip motion and the cooking reactions. A set of comprehensive differential equations governing the two-phase flow, heat transfer, and chemical reaction is developed, with the purpose of computationally investigating the cooking process in continuous digesters. A new cooking reaction model is also derived, based on extensive experimental data. To model the chip phase motion, the chip phase is assumed to be a non-Newtonian fluid with a slip boundary condition on the digester walls. A number of sub-models are examined and used to account for the compressibility of chip column and the inter-phase friction force. Based on the governing equations and mathematical models, a general two-phase three-dimensional CFD code has been developed to calculate the flow fields for both phases, solid fraction, temperature, and the concentrations of chemicals. The developed CFD code also has the capability of solving many other potential multiphase reactive flow problems. Liquid flow in a laboratory digester filled with compacted chips is simulated by the developed code. The computations compare reasonably well with the experimental observations by ERT (Electrical Resistance Tomography). Currently, the code cannot predict the tracer motion. Therefore detailed comparison between CFD simulation and ERT images cannot be made. Simulations are also performed in a full-scale industrial digester. The predictions of yield and kappa number at the digester exit are accurate for various cooking conditions. The distributions of pressure and chip phase volume fraction are discussed, the chip ii  phase velocity profiles at different digester levels are plotted, and the kappa number variation in the digester exit plane is also addressed.  Simulations also show that an  increased production rate results in a higher kappa number and a lower solid pressure at the bottom exit, and therefore the digester becomes more difficult or impossible to operate. The impact of variations of chip type and cold blow flow rate on cooking performance is also investigated computationally. The extensive investigation of the parameter variation impact on the process has shown that the model developed constitutes a reliable tool to predict the complex process in industrial digesters.  iii  Table of Contents Abstract  ii  Table of Contents  iv  List of Tables  vii  List of Figures  viii  List of Symbols  ;  x  Acknowledgements  xiv  Chapter 1 Introduction 1.1 Paper making and kraft pulping 1.2 Continuous digesters 1.2.1 Batch cooking and continuous cooking 1.2.2 Modifications of continuous digesters 1.3 Objective of this thesis Chapter 2 Mathematical Models of the Pulping Process 2.1 Literature review of digester modelling 2.2 Modelling assumptions 2.3 Governing equations 2.3.1 Continuity equation 2.3.2 Momentum equations 2.3.3 Species transport equation 2.3.4 The energy equation 2.4 Sub-models of the liquid phase 2.4.1 Flow resistance or inter-phase forces model 2.4.2 Compressibility equation for the chip mass 2.4.3 Concept of "equivalent diameter" of wood chips 2.5 Special treatment of the chip phase 2.5.1 Non-Newtonian fluid 2.5.2 Slippage on the wall boundary 2.6 Thermal properties of liquor 2.6.1 Density 2.6.2 Viscosity 2.6.3 Specific heat capacity and thermal conductivity  iv  1 1 2 2 4 6 8 8 12 15 15 16 17 19 20 20 23 25 27 27 30 32 32 32 33  Chapter 3 The Kraft Pulping Reaction Model 3.1 Representative reaction models 3.1.1 Single variable model 3.1.2 Three phase model: Gustafson's model 3.1.3 "Two lignin" model 3.2 Motivation for deriving a new reaction model 3.3 Regression method 1. Take full advantage of the available experimental data 2. Select suitable regression variables 3. Avoid using global variables 4. Provide sufficient data to ensure the model accuracy 5. Abandon the "three stage" regression form 3.4 Regression procedures and the pulping model 3.5 Model validation 3.6 Conclusion  35 35 35 36 38 38 39 39 39 40 40 41 41 46 48  Chapter 4 Computational Methods 4.1 Solving Navier-Stokes equations by SIMPLE method 4.1.1 Navier-Stokes equations 4.1.2 Collocated grids 4.1.3 Discretization of the equations by control volume method 4.1.4 Solution procedures for steady-state problems 4.1.5 Curvilinear grids 4.2 Additional treatment of two-phase flow equations 4.2.1 Evaluate the mass flux 4.2.2 Discretization of the inter-phase interaction terms 4.3 Discretization of scalar equations 4.4 Imposing boundary conditions 4.5 The code and its structure 4.6 Code validation 4.6.1 Flow in a cavity 4.6.2 Liquid flow through porous media along a flat plate 4.6.3 Power-law non-Newtonian fluid flowing through a tube  49 49 50 51 52 57 57 60 60 60 61 63 64 66 66 71 73  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester 5.1 Introduction 5.2 The model digester and electrical resistance tomography  76 76 77  v  5.2.1 The laboratory digester 5.2.2 Electrical Resistance Tomography 5.3 Computational method 5.4 Computational and experimental results 5.4.1 Approximation of boundary conditions 5.4.2 Liquid flow in the digester 5.5 Conclusions  77 78 80 82 82 84 89  Chapter 6 Simulation of an Industrial Continuous Digester 90 6.1 Introduction 91 6.1.1 Two-vessel hydraulic system in HSPP 91 6.1.2 Mill records of digester operations 93 6.2 Blended chip species 93 6.3 Computational domain 94 6.4 Boundary conditions 95 6.4.1 Boundary conditions for velocity 95 6.4.2 Boundary values for species equation (Lignin, Carbohydrates and EA)... 102 6.4.3 Boundary conditions for the energy equation 107 6.5 Simulation results 107 6.5.1 Boundary values and model constants 107 6.5.2 The mesh 109 6.5.3 Computational results 110 6.6 Discussion of simulation results 118 6.6.1 Mesh convergence check 118 6.6.2 Predictions of kappa number and yield 118 6.6.2 Variation of production rate 119 6.6.3 Cooking performance affected by the chip type 120 6.6.4 Effect of the cold blow rate 125 6.7 Three dimensional simulation 126 6.8 Summary 132 Chapter 7 Summary and Conclusion  134  Chapter 8 Recommendations for Future Work  136  Bioliography  137  Appendix A Mill record of digester operation at HSPP*  145  Appendix B Average mill records of digester operation at HSPP  151  vi  List of Tables Table 2.1 Flow resistance through packed wood chips  24  Table 2.2 The equivalent diameters of tested wood chips  26  Table 2.3 Coefficients for viscosity equations  33  Table 3.1 Constants fitted to pulping reactions  37  Table 3.2 Cooking conditions of cited experimental data  44  Table 6.1 Average fiber length of different wood species  93  Table 6.2 Specific gravity of some common woods  94  Table 6.3 Selection of recorded data from HSPP digester (09-MAR-2004)  98  Table 6.4 Calculation of solid/liquor phase flow rate through the top inlet  100  Table 6.5 Chemical composition of some soft wood species  104  Table 6.6 Chemical composition of some soft wood in North America  104  Table 6.7 Composition of white liquor in HSPP (g/L)  105  Table 6.8 Data relevant to the modified cooking circulation (09-MAR-2004)  107  Table 6.9 Computed boundary values for the case run on March 9, 2004  108  Table 6.10 Computed kappa # and yield versus mill records  119  Table 6.11 Typical combinations of incoming wood chips in HSPP  121  3  vii  List of Figures Figure 1.1 two-vessel hydraulic continuous cooking system (conventional)  3  Figure 1.2 Two-vessel hydraulic Kamyr EMCC continuous digester  5  Figure 2.1 Chip and free liquor phase in a REV  13  Figure 2.2 Inter-phase forces acting on different phases  16  Figure 2.3 Pressure drop over wood chips [7]  22  Figure 2.4 Solid fraction versus solid pressure [7]  23  Figure 2.5 Flow curves for a Bingham plastic and Newtonian fluids  28  Figure 2.6 Bingham plastic shear stress approximated by the Cross model [27]  29  Figure 2.7 Bingham plastic shear stress approximated by Pougatch et al's model [6].... 30 Figure 2.8 No-slip and slip conditions at the bottom wall  31  Figure 3.1 Effective alkali and sulfide consumptions vs lignin removal (h=5)  42  Figure 3.2 Fitted effective alkali and sulfide profiles  43  Figure 3.3 Lignin profiles by model predictions and experiments  46  Figure 3.4 Lignin profiles by model predictions and experiments  47  Figure 3.5 Model prediction versus experiments  47  Figure 4.1 Colocated grids in 2-D situation  52  Figure 4.2 Two-dimensional curvilinear grids  57  Figure 4.3 Illustration of symmetrical boundary conditions  63  Figure 4.4 Flow chart of the developed CFD code  65  Figure 4.5 Geometry and boundary conditions for 2D cavity flows  66  Figure 4.6 Flow field for 2D lid driven problem (Re=1000)  67  Figure 4.7 Velocity profiles in 3D lid driven problem (Re=100, mesh size: 25x25)  68  Figure 4.8 Velocity profiles in 3D lid driven problem (Re=400, mesh size: 25x25)  68  Figure 4.9 Stream function and temperature contours from Chen et al. (Re=100)  69  Figure 4.10 Stream function and temperature contours by this work (Re=100)  70  Figure 4.11 Streamlines in 2D lid-driven cavity with both side walls inclined at 45°  71  Figure 4.12 Streamlines in 2D lid-driven cavity by Ferziger and Peric's  71  Figure 4.13 Boundary layer flow passing through porous medium along a flat plate  72  Figure 4.14 Velocity profiles by code prediction and analytical solutions  73  Figure 4.15 Velocity distributions for power law fluids  74  vm  Figure 5.1 Configuration of the model digester  78  Figure 5.2 The laboratory digester and the ERT sensor planes  80  Figure 5.3 Flow resistance over chip column, for reference pine chips  81  Figure 5.4 Plane jet: water injected from triple holes, m =0.33kg/s  82  Figure 5.5 Plane jet: water injected from a single port,  83  0  ra =0.33kg/s 0  Figure 5.6 Computational domain, mesh, inlet and outlets  84  Figure 5.7 Flow field given by ERT image and CFD modelling  85  Figure 5.8 Flow field given by ERT image and CFD modelling  86  Figure 6.1 Two-vessel hydraulic digester system in Howe Sound Pulp and Paper  91  Figure 6.2 Illustration of the computational domain of the digester  95  Figure 6.3 Illustration of the liquor flow at the digester bottom  101  Figure 6.4 The extraction/injection ports of the modified cooking circulation  106  Figure 6.5 The structured 2-D mesh shown in reduced size for better viewing  109  Figure 6.6 Computational results of the operating case from HSPP  Ill  Figure 6.7 The liquor phase velocity and the chip phase pressure  112  Figure 6.8 The liquor phase velocity and the chip pressure at the extraction port  113  Figure 6.9 Velocity distribution of the chip phase  114  Figure 6.10 Chip phase velocity profiles at different  115  Figure 6.11 Kappa number variation at the bottom exit of the digester  116  Figure 6.12 Chip phase residence time versus chip locations  117  Figure 6.13 Computed downward liquor phase velocities using different meshes  118  Figure 6.14 Cooking performance affected by variation in production rate  120  Figure 6.15 Cooking performance affected by the chip type  121  Figure 6.16 Cooking affected by equivalent diameter of wood chips  122  Figure 6.17 Chip phase pressure affected by the equivalent diameter  124  Figure 6.18 Chip phase volume fraction affected by the equivalent diameter  124  Figure 6.19 Cooking performance affected by the cold blow rate  125  Figure 6.20 Liquor extraction shifted between screen sectors  .'  126  Figure 6.21 The structured mesh shown in reduced size (original size:30x20xl80)  127  Figure 6.22 Velocity and pressure distributions for the liquor phase  128  Figure 6.23 Liquor phase streamline and the chip phase pressure  129  Figure 6.24 The velocity of the chip phase  130  Figure 6.25 Computed distributions of scalar variables  131  ix  List of Symbols  ^li '  coefficients in discrete equations  ' ''  constants for liquor viscosity empirical constant in modified Ergun equation  B  A,  coefficients of discrete equations A , Ay, A x  surface area of control volume empirical constant  Ck C  z  /"  PJ  C  heat capacity of cooking liquor, kJ/(kg-°C) empirical constant in equation (2.22), MJ/(kg-m ) 3  C  mass fraction of carbohydrates  c  convective flux in x direction  x  coefficient of diffusion rate; diffusivity, m /s  Z)  2  molecular diffusion coefficient, m /s D  dls  hydraulic dispersion, m /s  D  eq  equivalent diameter, meter active energy per mol for lignin, kJ/mol  d ,exp //c//  diffusive flux on e-face approximated by explicit scheme inter-phase friction force, Newton the constant of gravity, 9.8 m/s the hydrosulfide ion concentration, g/L  k  relative reaction rate of lignin  k  constant parameter in cross model species related constants in cooking reactions  L, L  r  mass fraction of lignin, %  X  characteristic slip length, in meter m, n  power index in non-Newtonian models  m  mass transfer rate between chip phase and free liquor phase  m  mass flow rate of the entrapped liquor phase, kg/s  OH  hydroxyl ion concentration  OH'  hydroxyl ion concentration  Pc  pressure in chip phase, Pa  Pf  pressure in free liquor phase, Pa  Pi  prediction pressure, Pa  p  correction pressure, Pa  Q,Q,  heat generated or supplied, source term in i direction  QA, QR  volume flow rate in axi- and radial directions respectively  R  gas constant, kJ/mol/K  n  m  R ,R  constants in flow resistance equation  R,  the generation rate of species i  REV  representative element volume  [OH]  mole concentration of hydroxide, mol/L  [SH]  mole concentration of sulfide, mol/L  th  heating time, s  tc  cooking time, s  T  time, s  T  temperature, °C  T  cooking temperature, °C  t  c  2  velocity of chip phase, m/s velocity of free liquor phase, m/s u]  prediction velocity, m/s  xi  Uj  correction velocity, m/s  Uj  calculated velocity without pressure gradient, m/s  ( *]  interpolated prediction velocity on the east face of a control volume, m/s  Vo  superficial velocity, m/s  V  the volume of a REV, m  mv  V  fm  volume of free liquor, m  tigmr  V  ,  soM ma  3  volume of solid material, m  3  Usiip, w iip  slip velocity at the wall, m/s  Wj  molecular weight of species i, g/mol  Xj  denotes x, y, z when 1=1,2,3  X  dry solids concentration, %  [Xj]  mole concentration of species i, mol/L  Y  yield  Yj  mass fraction of species i  Yj  species mass fraction in cooking liquor  s  Greek letters  a  dispersivity of porous media, meter, equation (2.17)  a  empirical constants in A-S equation (2.26)  a  dimensionless parameter in equation (2.44)  a b  mass fraction of carbohydrates based on oven-dry wood  ctij  mass fraction of lignin based on oven-dry wood  a  over-relaxation parameter in Equation (4.18)  car  g  u  xii  p  empirical constants in A-S equation  p  dimensionless parameter in equation (2.44b)  K  kappa number  V  the kinematic viscosity, m Is molecular conductivity of cooking liquor and water, w/(m-°C)  c  volume fraction of the chip phase  c\  volume fraction of wood material  c2  volume fraction of entrapped liquor  f  volume fraction of the free liquor phase  S  S  £  £  Vf  dynamic viscosity of the free liquor phase, Pa-S dynamic viscosity of the chip phase, Pa-S stress tensor for the chip phase stress tensor for the free liquor phase  Pc  density of the chip phase, kg/m  Pf  density of the free liquor phase, kg/m  3  Viscosity, Pa-S TleflF  equivalent viscosity, Pa-S  A.  thermal conductivity, w/(m-°C)  7  shear rate, s~'  T  shear stress, Pa volume of computation cell, m  A  viscous friction coefficient, Pa-S/ m  Xlll  Acknowledgements I would like to express my gratitude to all those who gave me the possibility of completing this thesis. I would like to express my deep and sincere gratitude to my supervisors, Professor Martha Salcudean and Professor Ian Gartshore, whose gracious support, warm encouragement, and timely guidance helped me all the time in this research work. I am also extremely grateful to both my supervisors for spending significant amounts of time reading and correcting my thesis, and for giving me timely feedback even in some difficult times. I am also deeply indebted to Professor Chad Bennington in the Pulp and Paper Center, Chemical and Biological Department, at the University of British Columbia, for offering the opportunity to simulate the liquid flow in his laboratory digester; this work has served as a good basis for this thesis study. Professor Bennington has also given me useful advice and suggestions for improvement. I would like to acknowledge with much appreciation the help and suggestions given by Professor Steven Rogak and Professor Carl Olliver-Gooch. As my research committee members, both professors have given me much assistance in better representing the transport behavior in the porous medium. I would like to give my special thanks to my master program supervisor, Dr. Jerry Yuan, for his kind support, encouragement, and valuable input to my research. Dr. Yuan has always been of great help, and his knowledge and ideas contributed significantly to my research work. Many thanks to Garry Pageau, Senior Process Engineer in Howe Sound Pulp and Paper, who has generously provided me with the mill records and other useful data; Garry has also offered lots of suggestions and advice, which contributed significantly to  xiv  this thesis work. I am also grateful to Drs. Peter Gorog and Marguerite Gorog from Weyerhaeuser Company, they have given me very useful information and advice. Also, thanks to Konstantin Pougatch for the very helpful discussions of the equations and non-Newtonian fluid model. Thanks to my colleagues, Xiaosi Feng, Suqin Dong, Zheqiong Wang, Mohammand Shariati for their help and friendship, which made my student life more memorable. Finally, I would thank my parents and my wife, Yamei, whose tremendous love, understanding, and encouragement enabled me to complete this work.  xv  Chapter 1 Introduction  1.1 Paper making and kraft pulping In the paper manufacturing, wood timber is transformed into wood chips which are then pulped through one of the following three pulping processes: sulphate or kraft, sulphite and Chemi-Thermo-Mechanical-Pulp (CTMP). The kraft pulping process involves the digesting or cooking of wood chips at elevated temperature and pressure in "white liquor", which is an aqueous solution of sodium sulfide and sodium hydroxide. The white liquor chemically dissolves the lignin that binds the cellulose fibers together. The kraft pulping process is used to produce the majority of pulp worldwide. The aim of kraft pulping is to chemically remove enough lignin so that the fibres are free and to give them some of the required characteristics for paper-making at the lowest possible cost. To evaluate the quality of the produced fiber solution or pulp, a large number of properties are measured or tested, among which are kappa number and yield. The kappa number and yield are important control variables; they are also chief variables in our model predictions; the definitions of these two variables can be found in Michelsen's thesis [1]: "Kappa Number is a term used to define the degree of delignification. It is commonly measured from laboratory tests, typically with one hour intervals. In some modern mills the kappa number is measured on-line in the blow line. The measurement is based on reaction with acidic permanganate solution (KMn0 ), giving the amount of 4  lignin components in a sample from the pulp." "The yield expresses the total residual content of solid material in the pulp. It is defined simply as the sum of the residual content of lignin and carbohydrates. Hence,  1  Chapter 1 Introduction  high yield means that a large part of the wood remains in the pulp, and the consumption of wood is low." Rydholm [2] defined the kappa number (K) mathematically by Equation (1.1): K = (%(a l(a +a \))lc lig  lig  carh  (1.1)  K  The factor c is 0.14 according to Christensen [3]. In most pulping wood, the initial K  mass fraction of lignin and carbohydrates ranges from 0.95 to 0.97; therefore equation (1.1) can be approximated by, /cx0.15 = % «  %  For instance, when wood chips are cooked to a kappa number of 30, the corresponding lignin content is 4.5%. The yield is expressed as: Y = %(a +a ) n  (1.2)  carh  1.2 Continuous digesters 1.2.1 Batch cooking and continuous cooking Kraft pulping occurs in either a batch or continuous digester system. For batch cooking, a mill normally uses several vessels. Wood chips and cooking liquor are added simultaneously to the digester. Then the digester is sealed and raised to a target operating pressure and temperature (approximately 800 kPa, 170°C) for one to two hours; the pulp is then blown into an atmospheric tank, from which it is pumped to the brown stock washing system. For continuous cooking, the most common continuous digester is the vertical down flow type.  2  Chapter 1 Introduction  White liquor  Wash liquor  Figure 1.1 Two-vessel hydraulic continuous cooking system (conventional)  Figure 1.1 shows a conventional two-vessel hydraulic continuous cooking system. The wood chips are preheated in a steaming vessel to remove air and some of the volatile wood constituents before the chips enter the digester. The chips are mixed with the cooking liquor and are fed into an impregnation vessel, where the chips are supposed to be fully impregnated. Chips are then sent to the top of the digester vessel so that they move downward by gravity. The hydraulic pressure in the digester vessel is kept at approximately 1,140 kPa, and the temperature is approximately 160 to 170°C. For typical softwood wood chips, one to one and a half hours of residence time is needed to ensure complete pulping reaction. The reaction is stopped by introducing a low temperature liquor flow (which is called quench flow) in the lower region of the digester, where pulp washing is carried out, normally using a countercurrent flow of the filtrate from the first  3  Chapter 1 Introduction  brown stock washer. The pulp is blown from the bottom of the digester to a tank at atmospheric pressure. A batch pulping system enables a mill to pulp several different grades at once by using different digesters for different pulping conditions or fiber furnishes. The continuous process that produces pulp at a consistent rate with simple operations is more space efficient than the traditional batch process. Continuous cooking has gradually become the dominating pulp producing process. Today, over 65% of the kraft pulp produced in the world comes from continuous digesters, and most new digester installations are now of a continuous type.  1.2.2 Modifications of continuous digesters Until the early 1980s, most pulping processes used a "conventional" cooking mode, but modified cooking practices became available thereafter. Conventional cooking in continuous digesters means adding all the required white liquor within the chip feeding system, thus impregnation and cooking occur entirely in the cocurrent cooking zones. To decrease the bleaching chemical cost and lower the environmental impact, it is important to dissolve the lignin as much as possible in the cooking process. However, conventional cooking in very low lignin concentrations will cause severe cellulose degradation and thus decreases the yield and pulp quality. "Extended Cooking" originally proposed by the Swedish Forest Products Research Institute (STFI) in the later 1970s gained increasing popularity since the late 1980s [4]. With this technique, lower lignin content can be achieved with retained strength and pulp yield. Figure 1.2 shows a typical installation of a two-vessel hydraulic Extended Modified Continuous Cooking (EMCC) digester developed by Kamyr AB.  4  Chapter 1 Introduction  Figure 1.2 Two-vessel hydraulic Kamyr E M C C continuous digester As compared with conventional cooking, the cook starts at a reduced concentration of effective alkali and ends at a higher concentration of effective alkali. The lignin concentration in the chip phase is reduced at the end of the cook. This is accomplished by splitting the alkali charge into different insertion points and by ending the cook in the countercurrent zone where the liquor flows in an opposite direction to the chips. As figure 1.2 shows, the cooking phase in the digester is divided into four different zones. The white liquor charge, comprised of sodium hydroxide and sodium sulfide, is added at four different points in the digester, instead of at just one point as with traditional continuous cooking. Two benefits are gained by adding white liquor to the countercurrent zone. Firstly the initial hydroxide concentration is lowered, which prevents fibers being partially damaged; secondly, some cooking is performed in the countercurrent zone where the dissolved solids are lower and cooking proceeds more effectively.  5  Chapter 1 Introduction  1.3 Objective of this thesis In kraft pulping, digesters play a critical role in shaping pulp and paper qualities. An ideal cooking process leads to reduced consumption of cooking liquor and external heating, thus the recovery load is reduced; it also minimizes the kappa number variation of the blowing pulp and improves pulp washing. However, the physical phenomena in continuous cooking process are quite complex involving chemical reactions and transport phenomena like multiphase flow, mass, momentum and heat transfer. In a commonly used EMCC digester, the counter-current liquor flow combined with a number of liquor circulations retards the motion of the solid chip mass, provides the chips appropriate reaction conditions such as high cooking temperature, enough chemicals and appropriate residence time. The cooking performance is essentially a combined result of many variables. Although industries have achieved great progress in continuous cooking, many underlying physical and chemical mechanisms are still not well understood. Recent studies [5,6] have shown that there is a need to improve further the cooking performance, and therefore efforts are still required to optimize the process and equipment design, operation and control. Some experimental work has been carried out to understand the liquor and chip flow in digesters either in plants or laboratories, but this is still very limited since making measurements in a digester is very difficult. Moreover, these measurements cannot fully represent the process. On the other hand, Computational Fluid Dynamics (CFD) has become a powerful and widely used tool in the investigation of industrial processes. The growing success in CFD technology has sparked an increasing interest in a broad range of industries. In the pulp and paper industry, kraft pulping has also drawn attention of CFD researchers; important CFD work on kraft pulping was carried out by Harkonen [7], Michelsen [1,8], He and Salcudean [5], Pougatch and Salcudean [6,9], Miyanishi and Shimada [10]. The objective of the present work is to investigate computationally the continuous digester in order to improve process understanding and develop a tool that can be used for process optimization and apply the tool to an industrial digester. 6  Chapter 1 Introduction  The thesis presents the research undertaken to achieve this objective as follows: In chapter 2, the two-phase two-fluid models to describe the flow process in digesters are presented; issues related to chip compressibility and the wall stress model are also addressed. A general pulping reaction model is derived in chapter 3. Compared with most existing models, the correlation does not appear in a piece-fitting form, thus no transitional parameters are needed; besides, this model is obtained from extensive experimental data and thus gives more accurate predictions over a wide range of cooking conditions. Based on the governing equations, a general 3-D two-phase computational code is developed in chapter 4. The code uses non-staggered mesh and has a second-order of accuracy. The discretization of the governing equations is described in detail, and various validation cases are also given in this chapter. Chapter 5 gives the simulation of the liquid flow in a laboratory digester where stationary wood chips are added with compaction. The simulation results are then checked using Electrical Resistance Tomography images. In chapter 6, simulations of full-scale industrial digesters are presented. The velocity profiles, pressure distributions etc. are discussed; the impacts on cooking performance caused by production rate, chip type and cold blow rate are also investigated. A brief summary of the thesis is given in chapter 7, important conclusions are also summarized in this chapter. Chapter 8 contains the recommendations for future work.  Chapter 2 Mathematical Models of the Pulping Process  This chapter gives a brief review of previous work in digester modelling. The mathematical model used in the thesis is then described; and the main assumptions on which the model is based are stated.  2.1 Literature review of digester modelling As mentioned in Chapter 1, in the kraft pulping process, strong interactions through mass, momentum and energy exchanges occur between the cooking liquor and the downward moving chip phase; there are also multiple flow circulations and various disturbances from inflows. Awareness of the complexity of the pulping process and its importance to the pulp and paper industry, prompted work on digester modelling over the years lead to two different methods. In the first approach, both liquor and chips are stationary and chip compressibility is usually neglected, whereas diffusion and chemical reactions occurring in the pulping process are of particular interest. One of the earliest kinetic models was developed by Vroom [11], who established an Arrhenius type of relation in the form of (2.1)  %- = -k(T) at  where k(T) = ae-  plT  Here, L is lignin content and T is the temperature, a and P are positive constants. Vroom then further derived the concept of " H " factor,  8  Chapter 2 Mathematical Models of the Pulping Process  H-factor= \ e  { A x m i T )  dt  (2.2)  By the above integration, the H-factor combines temperature and time into a single variable representing the extent of the cooking; it has been widely used in the control schemes in today's pulping industry. In Kerr's model [12], the pulping process consists of five reaction steps: (1) Transport of chemicals to the site of reactions by flow of liquor into a chip (penetration) and molecular diffusion; (2) Adsorption of chemicals on the reaction site; (3) Chemical reactions; (4) Desorption of reaction products; (5) Transport of reaction products away from the site by diffusion. Kerr expresses the delignification rate as a function of effective alkali, lignin content and temperature, which is more complicated than Vroom's model. By virtue of Vroom's method, the connection between lignin content and " H " factor is established through a single equation in which numerous pulping variables are incorporated. In this model, the relationship between effective alkali concentration and pulp lignin content in initial and bulk phases was approximated by two different expressions. This makes the model very complex. Kerr later found that the relationship could be well approximated by a single expression, and this greatly simplifies the model, while the predictions are still accurate [13]. Gustafson [14] provided a good review of early pulping models, and presented a comprehensive model to study the combined diffusion and kinetics within a wood chip during kraft pulping. Assumptions made in his work include fully impregnated wood chips, well-stirred homogeneous bulk phase, constant sulfide concentration and so on. Kinetic models and data were carefully selected, and the model predictions were shown to be remarkably good when compared to some available laboratory cooking experiments.  9  Chapter 2 Mathematical Models of the Pulping Process  The other approach, mainly applied in the continuous cooking process, recognizes the importance of fluid motion on the pulping process performance, and attempts to simulate the overall pulping performance resulting from liquor and chip flows, cooking reactions, heat transfer etc. By assuming plug flow for both chips and alkali solution, Christensen et al. [15] proposed a comprehensive one-dimensional model for continuous digesters. To describe the pulping chemistry, wood is considered as a mixture of four families of compounds: cellulose, hemicellulose, lignin and extractive compounds. In this study, entrance effects, exit effects and cross flows near the extraction ports were assumed to be of minor importance; in addition, it was assumed that temperature and concentration gradients occur only in the vertical direction. A good attempt at simulating the liquor flow two-dimensionally in a digester was made by Harkonen [7]. In his study, a chip is composed of soluble and insoluble materials, and the chip mass is treated as an isotropic porous continuum. Mass and momentum conservation equations are written separately for the solid phase and liquid phase and inertial terms are ignored. Interphase drag is calculated using Ergun's Equation. An energy balance in both phases includes the convective and conductive terms, reaction heat and heat exchange between the chip and liquor phases. A test digester with a diameter of 190 mm and a height of 360 mm was constructed to investigate experimentally the chip compressibility and flow resistance through the chip mass. In 1995, Michelsen [8] proposed a more comprehensive one-dimensional twophase model, where convection, reaction, diffusion and dispersion were all included. In this model, wood chips are assumed to consist of lignin and carbohydrates, with the free liquor consisting of effective alkali, dissolved solids and water. The first order pulping reaction model given by Christensen et al. [15] was adopted. The inter-phase diffusion rate coefficients for the alkali and dissolved solids were adopted from Gustafson et al. [14] and Christensen et al. [15]. In steady state situations, the model gives reasonable predictions to the compaction and pressure of the chip mass, the velocities of free liquor and chip plug, kappa number etc. However there are limitations due to the one-  10  Chapter 2 Mathematical Models of the Pulping Process  dimensional approach as there could be significant radial variations in flow and associated phenomena, especially in modern digesters with large dimensions. The dynamic behavior was also studied for particular variations that occur in Michelsen et al's work. Miyanishi and Shimada [10] simulated the Low-solids cooking process and made some comparison with conventional cooking. The computational results are impressive because even the cellulose degradation has been well captured. The previous work cited has made a valuable contribution to the understanding of the kraft pulping process. A general deficiency of these models lies in their one-dimensional basis, so that many important cooking features, like secondary flow, cross-flow diffusion and non-uniformity cannot be predicted. Harkonen's two-dimensional simulation is based on a potential flow treatment and does not address sufficiently the chemical kinetics. On the other hand, the constantly increasing need for a high quality pulp at low cost always requires better understanding and more precise control. Over the past years, continuous digesters have been greatly enlarged and have incorporated novel technologies, such as E M C C and Low-Solids cooking. The increased process complexity has raised many practical concerns, such as safe operation and uniformity of cooking. Hence, it is important to develop a more comprehensive model capable of simulating the pulping process that occurs in digesters. It would be worthwhile to compute the flow three-dimensionally as in some cases the liquor injection or extraction does not take place axi-symmetrically and that can result in circumferential non-uniformity of the pulp. Preliminary work [5,6] has demonstrated the possibilities for numerical modelling of digesters but resulted in a numerical code that was not robust and exceedingly difficult to converge. The assumption that the chip phase was a Newtonian fluid is a weakness of this model that resulted in some unphysical behavior. The work in this thesis starts with a fresh and different approach based on a collocated grid (as opposed to the staggered grid used in the earlier approachs). The present work is also focused on the simulation of a modern industrial digester (Howe Sound Pulp and Paper M i l l at Port Mellon BC) that would allow extensive testing and  11  Chapter 2 Mathematical Models of the Pulping Process  comparisons. Currently, after significant development of the preliminary numerical scheme, both codes have the capability of performing computations of a digester. The earlier approach has now been developed successfully and a description of its use and application has been accepted for publication [6]. A comparison between sets of numerical results produced by the two different codes for a single application has been made. This will be described in Chapter 4 of the present thesis. In the following sections, the assumptions and differential equations used in the present thesis are described; sub-models of inter-phase forces and chip compressibility are also discussed.  2.2 Modelling assumptions Wood chips used in the typical industrial cooking process are very small, with a thickness ranging from 1 to 10 mm. Practical industrial digesters, on the other hand, have large dimensions; for a digester with a 1,000 ton/day capacity, typical dimensions are 50 meters in height, and 6 meters in diameter. Even a small digester has dimensions which are much bigger than wood chips, so that detailed CFD computations at a scale of a single wood chip are obviously computationally prohibitive. It is therefore necessary to simulate the process from a "macroscopic" point of view. The cooking conditions in this industrial pulping process must be made as uniform as possible to minimize the non-uniformities in the pulp. In practice, digesters are operated carefully to ensure a long impregnation, enhanced diffusion, low gradients in temperature and alkali concentration. This makes a macroscopic approach more acceptable. In this approach, the mesh size is so big that a typical control volume or representative element volume (REV) contains a large number of chips, as schematically shown in Figure 2.1. The CFD calculation will provide the mean quantities such as velocities and the pressure in such an REV.  12  Chapter 2 Mathematical Models o f the Pulping Process  Chip phase includes solid wood chips and liquor entrapped in the pores o f wood chips Liquor phase refers to the liquor that can move freely in the interstitial space o f wood chips  Figure 2.1 Chip and free liquor phase in a R E V As shown in Figure 2.1, wood chips and the entrapped liquor are defined as the chip phase, while the liquor outside the chips is called free liquor phase or liquor phase, for simplicity. To simulate mathematically the flow in a continuous digester, we make the following assumptions: 1. Air is completely removed from the wood chips; as a result, the voids or cavities in the wood are totally filled with the liquid. We further assume that this part of the liquid (bound liquid) will have the same bulk velocity as the wood chips, so that we will describe the bound liquor and the wood chips as the chip phase. The remaining liquor, often referred to as free liquor or the free liquor phase, will interact physically and/or chemically with the chip phase. In Figure 2.1, sj refers to the volume fraction of the free liquor phase, and s refers to c  the volume fraction of the solid/chip phase, which includes two parts, one is solid material s i and the other is entrapped liquor s 2. Specifically, we have, c  C  V, free  liquor  Vchip phase  (2.3)  REV  Vsolid  material REV  Ventrapped V v  13  REV  liquor  (2.4)  Chapter 2 Mathematical Models of the Pulping Process  Immediately, we can get the following relationships: (2.5)  e +e = 1 c  f  (2.6) These relationships are with reference to a control volume. 2. Wood chips have no specific shape. The chip mass is assumed to be isotropic and homogeneous. Individual chips, as we know, are neither isotropic nor homogeneous. However the chip mass, namely a collection of chips, tends to be distributed randomly when compacted in typical pressure range, therefore assuming the chip mass to be isotropic and homogeneous is usually a good approximation. This assumption is supported by the observation that no significant difference was found from a macroscopic point of view, when measuring diffusivities in the three primary directions [14,16]. 3. Solid material in wood chips consists of lignin, carbohydrates, and insoluble material. When lignin and carbohydrates dissolve and enter into the free liquor phase, the corresponding volume will be simultaneously filled by free liquor. 4. The chip mass is compressible, but each individual chip retains its volume throughout the whole cooking process. In other words, when solid pressure (stress) is imposed, the interstitial space in the chip mass will decrease; each single chip may change its shape but not change its volume. This proves to be a valid assumption from pilot plant experiments [8]. 5. The chip mass will be regarded as a moving porous medium.With this assumption, the digestion process can be considered as a continuous two-phase flow, where the free liquor (or unbound liquor) flows through a compressible moving chip bed.  6. The cooking liquor is assumed to be a solution of sodium hydroxide (NaOH) and sodium sulfide (Na2S), with the addition of a certain amount of dissolved wood substance. The main feature of the cooking liquor is that its ingredients (or species concentrations) and density are changing as chemical reactions proceed. The chemical  14  Chapter 2 Mathematical Models of the Pulping Process  kinetics are characterized in terms of effective alkali concentration and sulfite concentration. 7. The flow is laminar. A continuous digester usually involves several liquor circulations, which supply the liquor to the digester typically through central downcomers. These inflows are basically liquor jets, where turbulence tends to occur even at low Reynolds numbers. However the flow in the digester has a very high solid fraction and the occurrence of turbulence over large regions is unlikely. Additionally, studies on fluid flow through porous media could give some idea about the flow pattern in a digester. Michelsen [8] estimated the Reynolds number based on the average grain diameter and found it to be in the range of 0.8-15 for typical digester operations; he then concluded that the flow is mainly in the laminar domain. 8. Pulping reactions are irreversible. Alkali diffusivity depends on temperature and the ingredients of the cooking solution. 9. Chips are in thermal equilibrium with the local cooking liquor, i.e., chips and cooking liquor in the same REV always have the same temperature. Because the temperature in a continuous digester does not change sharply in space, and solid chips only move on the order of millimeters per second, this assumption should be acceptable.  2.3 Governing equations Based on the previous work [5], the macroscopic equations for the two-phase flow are the following. 2.3.1 Continuity equation Mass conservations are expressed as: (2.7)  dt  (2.8)  —mex  15  Chapter 2 Mathematical Models of the Pulping Process  Please note, in this thesis, subscript "f' and "c" refer to free liquor phase and chip phase respectively. In Equation (2.7) and (2.8), m denotes the mass exchange between cx  the chip phase and the free liquor, and will be represented by the chemical kinetics and mass transfer models. 2.3.2 Momentum equations The single-phase momentum equation can be generally written as, f  dp dx,  +  dr dx  N  (2.9) ./ J  where p is the static pressure, r, is the stress tensor (described below), pgj and Fj are the y  gravitational force and other body forces in the i direction, respectively. du, v  dxj  du.  2  du,  M—  dx,. 3  s  (2.10)  a  dx, "  where p is the molecular viscosity and the second term on the right hand side is the effect of volume dilation. Based on the assumptions made in section 2.2, a two fluids two-phase model will be used to predict the flow in the digester. Figure 2.2 schematically shows the inter-phase forces in the x direction for the free liquor and chip phases, where Ff denotes the friction force per unit volume between the two phases, and is considered as a body force.  F,f <:  (a) Free liquor phase  (b) Chip phase  Figure 2.2 Inter-phase forces acting on different phases (within the REV, x direction only)  16  Chapter 2 Mathematical Models of the Pulping Process  From Equation (2.9) and the forces shown in Figure 2.5, the momentum equations for each phase can be written as, J (Pf f fj) £ u  t  + -^T (Pf f fJ fJ £ u  )=.  u  e  V  dx:  -r (Pc c"cj ) + -r— (p e u u ) =  dx  c  dx,  c  CJ  f  ex  c  (2.11)  JJ + cPcgi - F j - s £  £  dt  + EfPfgt + F ,i + m u ,  CJ  dx  JJ  f  c  dx,  (2.12) where F is the inter-phase friction force acting on the fluid per unit volume in the i fj  direction, and is actually the pressure gradient caused by the inter-phase friction force. This will be further discussed in a subsequent section. The stress rl is defined in the same way as in (2.10), and z is defined by Equation (2.14), as the chip phase is assumed c  tj  to be a Bingham plastic fluid. 2 K  dXj  dx,  du  r  hC  dxj  K  |  du  (2.13)  lf  y  dUj^ dx,  2  du  lc  (2.14)  ;  2.3.3 Species transport equation In general, every species equation is determined by chemical kinetics as well as diffusion. Early models in kraft pulping assume that the cooking process is kinetically controlled. This assumption is acceptable for batch cooking. But for continuous cooking, where both chip mass and liquor flow continuously and liquor circulations are often used, ignoring diffusion is no longer reasonable. Apart from this, many of the commercial chips used today are quite large, making the overall reaction rates essentially dominated by both chemical kinetics and diffusion.  17  Chapter 2 Mathematical Models of the Pulping Process  As a matter of fact, efficient diffusion is desirable for industrial digesters because it results in fast mixing of chemical reactants giving high production; high diffusivity also helps to distribute products evenly, ensures complete reaction and therefore decreases non-uniformity of cooking. In the history of the development of continuous digesters, the introduction of a central pipe was considered to be a major advance as it enhances the diffusion and greatly decreases the pulp non-uniformity. Therefore, solving general transport equations for some important species or components is necessary. In the free liquor phase, transport equations for effective alkali or sulfite can be written generally as,  er  1  ''  dx,  y  j  '  1  l j  dx  + R,  (2.15)  ./ J  J  where Yj is the mass fraction of effective alkali or sulfite ions over the total mass of the free liquor; R is the mass rate of creation or depletion of each species by the chemical f  reaction model, which will be described in detail in chapter 3; Dj is the equivalent diffusion coefficient, which is considered as the sum of a contribution from molecular diffusion and from hydraulic dispersion: D = D +D 0  (2.16)  dls  In Equation (2.16), D is the molecular diffusion coefficient, D 0  dis  from hydraulic dispersion. Scheiddegger [17] found that D  dis  is the dispersivity  depends linearly on the  velocity and a characteristic of the porous media referred to as dispersion coefficient; he proposed the following model: D =a\U\  (2.17)  dls  where U is superficial velocity in m/s and a is the dispersion coefficient of the porous media in meter. The average grain diameter is considered as a typical value for the dispersion coefficient of homogeneous sands [18]. However, the dispersion coefficient for wood chips has not been investigated yet. Huang et al. [19] computationally studied  18  Chapter 2 Mathematical Models of the Pulping Process  the diffusion of salt water in soil, using a dispersion coefficient of 0.01 meter. In Huang's computations, the porosity was in the range of 0.35 to 0.62, which is typical for our digester simulations. Since there is no information available in literature for wood chips, we assume in our study a dispersion coefficient of 0.01 meter. However, with the same porosity as the soil, compacted wood chips have different mean grain diameter, thus different mean free path, and therefore the dispersion coefficient is expected to be different as well. The superficial velocity difference in a digester is in the order of 10"  3  m/s. Based on equation (2.17), the diffusivity caused by hydraulic dispersion is expected 5  2  to be of the order of 10" m/s. A sensitivity study in digester modeling has shown that when a is less than 0.1 meter, the variation of a does not influence significantly the distribution of species concentrations. However, large values of a result in noticeable difference in the distribution of species concentrations. The mass fraction Yj in Equation (2.15) can be easily converted to its molar concentration by the following expression: (2.18) Pf  And equivalently, Equation (2.15) can be written as dip  [X,])  f  dt  dx  v  7  /L  , J ;  ;  dx  dXj j  + ^Pf- L  (2.19)  W,  The species equations of lignin and carbohydrates are quite similar to Equation (2.15). Because lignin and carbohydrates do not diffuse freely in the solid phase, i.e., D, =0, the corresponding species equations can be simplified to,  ^(Pc^)+/-U"//)=^ jr  dt  (2-20)  dx ,  2.3.4 The energy equation Chemical reactions occurring in cooking processes are temperature sensitive. In normal cooking conditions, the delignification rate approximately doubles for every 8-  19  Chapter 2 Mathematical Models of the Pulping Process  10°C increase in temperature, and the diffusivity of the reaction agents also strongly depends on temperature [20]. Therefore accurately predicting the temperature distribution in the digester is very important. The heat released in the cooking reactions is not very high, and therefore in normal cooking conditions, the chip phase and free liquor are assumed to have the same temperature within each REV. Thus the following energy equation applies to both phases: —(p c„ dt f  v  y  p j  ,T)+—{pfUf ' cbc, ' v  ,c„ f j  p J  '  — dxj  . dT f  Dp,  +  ^LL  Dt  d{u,T,,) +  ^±JL>  dxj  (dl dC\ Q = C \— + — \dt dt)  +  Q  (2.21)  (2.22)  R  here, C is an empirical constant, with a value of 0.218 MJ/kg/m for dry wood as R  suggested by Gullichsen [21]. On the right hand side of Equation (2.20), the first term is the heat transport by conduction; the second term accounts for the heat converted from flow work; the third term is the heat generated by viscous dissipation; and the last term Q includes the heat of chemical reaction and external heat supply. The second and third terms are negligible for the digester computations.  2.4 Sub-models of the liquid phase 2.4.1 Flow resistance or inter-phase forces model  Flow through a porous medium is very complex, and extensive experimental work has been conducted for different combinations of solids and fluids. To predict the pressure drop resulting from flow through a porous medium, various models have been published in the literature, often in the following two popular forms [22]. A. Ergun and modified Ergun equation [23] form (M-E form)  F .— k  = B\ ( ~ ) A  l  (2.23)  £  where 20  Chapter 2 Mathematical Models of the Pulping Process  dPIdL pV  2  1  0  N  =  M  (2.25)  M  where s is the porosity of the porous medium and D is an appropriate characteristic eq  length for the medium. B. Ahmed and Sunada equation [24,25] form (A-S form) _ pV  dP/dL  0  = a + p—°pV  (2.26a)  p  0  dP  Or:  = apV +fipV 0  (2.26b)  2 0  dL  where dP/dL is the pressure gradient in the direction of the flow; Vo is the superficial velocity and is given by Vo=Q/A, where Q is volumetric flow rate and A is the total cross-sectional area perpendicular to the flow direction; a and P are model parameters to be established empirically. The Modified Ergun (M-E) equation [23] is well known as it includes the quantity porosity, which characterizes the porous medium. Also, it is believed that the M-E equation has a good theoretical basis. This equation does give very good results when the porous medium can be characterized fairly accurately. The A-S equation, on the other hand, does not explicitly include a parameter related to the porous medium. Therefore the parameters a and P must be functions of the medium and accordingly must be empirically established for each different medium. However, The A-S equation has the advantage of being simple and appears to be more accurate when fitted by a large number of experimental data, as concluded by MacDonald [22], especially when the porous medium is not well characterized. It is probably for this reason that most flow resistance models for liquor flow through wood chips are given in the A-S form.  21  Chapter 2 Mathematical Models of the Pulping Process  Harkonen [7] investigated experimentally the flow resistance in the presence of wood chips and established the following correlation: dP dL  f  \  (2.27) f  v  f  s  J  £  For pine Ri=4.6xl0 , R =3.9xl0 . The pressure drop strongly depends on velocity 2  6  2  difference and solid fraction, as plotted in Figure 2.3. 5.0  r  E Q.  Q.  <  CL  2  •D  £>  3 at to  £  CL  2  4  6  8  10  12  14  16  Velocity difference, AV (mm/s)  Figure 2.3 Pressure drop over wood chips [7] When both phases are in motion, such as in the digester, the formula should change accordingly. Michelsen extended the use of Equation (2.27) in his dissertation [8] and gave, dP dL  f  s  v  6  s  2  V-V,  e  n  f  (2.28)  t-v,)  f  where A is an expression of the viscous friction coefficient, and is a function of the porosity of the porous media and the velocity difference between the two phases. Here, both velocities are physical velocities or local velocities instead of superficial velocities which are also commonly used in literature. Equations (2.11) and (2.12) result in, dP_ ~dLJi  A  R ^  22  + Rs 2  c  u. c  •(u -u ) eJ  fJ  (2.29)  Chapter 2 Mathematical Models of the Pulping Process  2.4.2 Compressibility equation for the chip mass As stated before, compressibility is one of the important features of the chip mass. The chips become weaker as lignin is dissolved from the cell walls. Thus compressibility is a function not only of solid pressure, but also of lignin content in the chip mass: the less lignin there is, the softer the wood chip becomes, and the higher the chip compaction is. Harkonen [7] has established the following correlation based on his experiments, and plotted in Figure 2.4. (2.30a)  (2.30b)  Solid pressure (kpa)  Figure 2.4 Solid fraction versus solid pressure [7]  Further work [29, 30, 31] on chip compressibility and flow resistance over wood chips was conducted and correlations of chip compressibility and flow resistance in the form of Equations (2.27) and (2.30a) were obtained. Brief description of the chips tested and the resulting correlations are listed in Table 2.1.  23  Chapter 2 Mathematical M o d e l s o f the Pulping Process  T a b l e 2.1 F l o w resistance through p a c k e d w o o d c h i p s ( U n i t s : P c — P a , d P / L — k P a / m , superficial v e l o c i t y u - mrn/s ) Ref  R a w material  Equations e  Pine, size distribution not specified  = 0.0046  x ( - 0 . 8 3 1 + 0.139ln(x-))  u + 0.0039  f ' 0.663 + (P 110 )° new 4  silvestrus) Reference Chips: Fraction o f accepts chips o f sawmill chips (chips passed through a 13 mm round hole but retained on a 7 mm round hole, no other details) " N e w " Chips: 4 mm thick, 40 mm in length Pinus silvestrus (Scandinavian) Mix 1 4~6mm:22.1%; 2~4mm:44.2% 2~3mm:29.1%; <2mm: 4.6% C h i p bulk Density=165kg/m Mix 2 4~6mm:23.2%; 2~4mm:46.3% 2 - 3 m m : 30.5%; <2mm: 0% C h i p bulk Density =163kg/m Mix 3 4~6mm: 0%; 2~4mm: 60.3% 2~3mm: 39.7%; <2mm: 0% Chips: p =157 k g / m Basic wood: p =401 kg/m 3  3  4  L  Eucalyptus camaldulensis >45mm: 0.6%; 8~45mm: 9.2% 8~13mm:51.4%; <3itim: 1.1% 7~13mm:31.2%; 3~7mm:6.4% Scandinavian birch >45mm: 1.1%; 8~45mm: 4.8% 8~13mm:78.8%;<3mm: 0.5% 7~13mm:12.6%;3~7mm:2.2%; C h i p bulk Density=175 kg/m Basic wood Density =556kg/m 3  (l-s )  (\-£ )  2  = 0.052-  f  2  x ( - 1 . 0 9 2 + 0.161 ln(*r))  39  c  dP  ±L u  x ( - 0 . 7 8 8 + 0.133 ln(*r))  56  = 0.840 + (P 110 )° l—u  f  + 0.0015-  2  ^-u  2  ref  dP  s'T »iix2  c  (l-e )  2  0.082  f  1— u - 0.00011  = 0.604 + (P 110 )° 4  = 0.615+ ( P / 1 0 )  0 7 9  = 0.647 + (P / 1 0 )  0 7 4  4  •v  4  c  dP =  dP L  0  ^  K  2  x ( - 1 . 0 2 1 + 0.205 ln(O)  -w-3.5xl0"  1  u  x (-0.910 +0.181 ln(O)  . 0 2 8 ^ ^ , - 1 . 3 x l O -  = 0.051  2  (l-£f)  x ( - 0 . 9 5 6 + 0.191 ln(/c))  63  c  C  4  4 V  ^ - ^  ,  J  2 W  -  2  'mix!  dP  3  31  5 9  c  3  30  4  c  —  Scandanavian pine, (Pinus  29  = 0.644 + ( P / 1 0 ) °  f  ^u  w-5.7xl0"  = 0.0055  2  -- 0.591 + (P 110 )°  56  x ( - 0 . 6 4 5 + 0.148 ln(/c))  hir 0.630 + (P 110 )° f  64  x ( - 0 . 6 9 7 + 0.151 ln(*-))  4  c  4  c  dP L,„„. dP  0.552-  0.28  3  24  l—u  (I-*/)  + 0.00074-  2  u + 0.0012  ^—u  2  0-ff )  2 f ^ - u 2  Chapter 2 Mathematical Models of the Pulping Process  2.4.3 Concept of "equivalent diameter" of wood chips  The flow resistance through the chip phase could be very big and is usually a dominant factor that affects the chip compaction and hence flow in typical digester operations. Experiments have been conducted to measure the flow resistances over different wood chips and empirical correlations have been proposed as listed in Table 2.1. As can be seen, the correlations are all in the same form, where the flow resistance is linked to the velocity and the porosity of the chip bed; the influence of chip size, however, is not explicitly established in these correlations. Although wood chips have irregular shapes and non-uniform sizes, it is possible to derive an average diameter to characterize the whole chip mass; with this diameter, some popular semi-empirical models that describe the resistance of flow through porous media could also be used. In this section, the "equivalent diameter" of wood chips is proposed and determined so that the modified Ergun equation can also be used to calculate the flow resistance through wood chips. Suppose porous medium A] consists of uniform spherical particles and its flow resistance is given by the modified Ergun equation (2.23). By some mathematical manipulations, the following equation can be obtained, dP  f.iA  dL  D  1-g  V eq  e  )  dL  Dt \ f  f  \ £  (2.31) e  J  e = s and s = 1 - s into (2.31) further results in,  f  pA f  1-e v  Substitutions of V =  dP  pB  +  c  r  D \f  1  J  £  eq  f  \  •V  (2.32)  J  Now, if chip column A has the same flow resistance as the porous medium A], 2  Equation (2.27) and (2.32) will give identical results. Comparing these two expressions, we have,  25  Chapter 2 Mathematical Models of the Pulping Process  juA  (2.33a)  ' ~ D eq 2  R, =  pB_  (2.33b)  Theoretically, when the constants A and B are chosen, R] and R are experimentally 2  determined, both equation (2.33a) and (2.33b) can then be used to determine D . D is eq  e q  regarded as the equivalent diameter of chip column B. Macdonald [22] has found a correlation between A and D which works very well for Dudgeon's [32] and Gupte's eq  [33] data with porosities ranging from 0.366 to 0.64. When D is expressed in meters, eq  the correlation can be written as, A = 75.1 +1.148 xl0 £>  (2.34)  5  eq  Substituting (2.34) into (2.33a) and assuming a liquor viscosity of u=0.001kg/m-s, the equivalent diameter of the chip column is, 114.8 + 7114.8 + 0.3/?, 2  Deq =  (2.35)  2R,  With the calculated D and known constant R , the constant B can be determined eq  2  from equation (2.33b). The equivalent diameters of all the tested wood chips in Table 2.1 are calculated by (2.35) and listed on Table 2.2. Table 2.2 The equivalent diameters of tested wood chips Authors  Wood material  D , mm  Reference  eq  Harkonen  Scandinavian Pine  3.033  7  Wang Wang  Scandinavian Pine (Reference chips) Scandinavian Pine (New chips)  29 29  Lindqvist  Scandinavian Pinus silvestrus, Mix 3  2.735 1.885 1.394  Lindqvist  Scandinavian Pinus silvestrus, Mix 2  1.069  30  Lindqvist  Scandinavian Pinus silvestrus, Mix 1  0.838  30  Lammi  Scandinavian Birch Eucalytus Camaleulensis  0.762  31  0.488  31  Lam mi  26  30  Chapter 2 Mathematical Models of the Pulping Process  As can be observed here, the equivalent diameters differ significantly, from 0.488 to 3.033 mm as seen in the Table; even for chips from the same wood species, the equivalent diameter could be quite different. For instance, for the two types of pine chips tested by Wang [29], the equivalent diameters are 1.885 and 2.735 mm respectively. This indicates that the size distribution of wood chips could significantly affect the equivalent diameter and therefore the flow resistance of wood chips. In Chapter 6, the effect of the equivalent diameter of wood chips on pulping performance is computationally investigated.  2.5 Special treatment of the chip phase 2.5.1 Non-Newtonian fluid The chip mass has been assumed to be a moving porous medium. The motion of the solid chip phase cannot be well described by a Newtonian fluid model, because some onsite observations indicate that the velocity profile of the chip phase is quite flat near the digester periphery. The chip phase is therefore assumed to be a non-Newtonian fluid, whose behavior can be generally described by the Cross model [27], UlZJL = (Ky)  (2.36)  m  where r|o and r|oo refer to the asymptotic values of viscosity at very low and very high shear rates respectively, K is a constant parameter with the dimensions of time and m is a dimensionless constant. With some approximations, the Cross model can be reduced to some other popular viscosity models. For example, for i"|«r|o and r | » r | o o , the Cross model reduces to  27  Chapter 2 Mathematical Models of the Pulping Process  (2.37)  (Ky)"  V  With the definitions of K = 2  V=  K  and n = 1 - m , we further obtain,  7o  (2.38)  (Ky)"  which is the widely used power-law model and n is called the power-law index, K2 is called the 'consistency'. By intuition, the chip phase internal structure should be rigid enough to prevent movement for values of shear stress less than a specific value, and when the stress exceeds this value, the chip phase will exhibit some shear rate like a Newtonian fluid. This characteristic can be best approximated by a Bingham plastic fluid model; the stressstrain relationships of a Bingham fluid and a Newtonian fluid are shown in Figure 2.5.  Newtonian fluid Bingham fluid at  P  re o  shear rate, y  Figure 2.5 Flow curves for a Bingham plastic and Newtonian fluids  The rheological equation for a Bingham plastic may be written T-T =7jy 0  (2.39)  (r>r ) y  28  Chapter 2 Mathematical Models of the Pulping Process  However, an equivalent viscosity is desirable for computational purposes, that is, r = y  (2.40)  Vcff  The Bingham plastic behavior given by equation (2.39) can be approximated by the Cross model, with a specific example given in Figure 2.6. 80 I  0.8 Shear stress (Bingham plastic fluid) Shear Stress (Cross Model) Viscosity (Cross Model)  \ —  in  :\  [Pa  t 60 V 50 r  \  Viscosil  - \ /•  re CL  ,  40 - V  H 0.4  _ ' \ -/ \ 30 4 \  H0.3  in m  2  n  CD  20  0.2  10  Ho.1  .j 0.02  i  i  I  i  i  i  u  0.08  0.04 shear rate, y (s" ) 1  Figure 2.6 Bingham plastic shear stress approximated by the Cross model [27] Olo=100.Pa-S, 7 =l,K=150.0, m=1.05Pa) ra  In some situations, to in Equation (2.39) may also change. To better describe the chip phase motion in a digester, Pougatch et al. [6] assumed that To is a function of both shear rate and chip pressure, with the following form, (2.41)  r =Kp (l-e- ) mf  0  e  The effective viscosity is then expressed as, — = ri + r  (2.42)  -(\-e-'" ) r  r  Similarly, the Bingham plastic behavior can be approximated by Pougatch et al's expression (2.42); a typical example is given in Figure 2.7.  29  Chapter 2 Mathematical Models of the Pulping Process  60 55 50 45  «T  i  40 35  £  30  8 u  25  *  20 15 10 5 0 0.02  0.04  0.06  0.08  0.10  shear rate, y (s ) 1  Figure 2.7 Bingham plastic shear stress approximated by Pougatch et al model [6] (110=1.5^-5^=100.0, K=0.001, Pc=500 Pa) In Chapter 4, a power law fluid flowing through round tubes is simulated and computed velocity profiles are compared with available analytical solutions. 2.5.2 Slippage on the w a l l boundary  The no-slip condition always applies at wall boundaries for Newtonian fluid flow; in non-Newtonian fluid flow, however, relative slippage between the boundary and the fluid may exist, e.g. when squeezing tooth-paste in its tube. In Figure 2.8, the no-slip and slip conditions are illustrated at the fixed bottom wall. In Figure 2.8, the fluid velocity has a finite value U IJ at the bottom wall, which is different from that in 2.8. The slip condition S  p  is typically used for micro scale or "free" surfaces at macroscopic scale [34].  30  Chapter 2 Mathematical M o d e l s o f the Pulping Process  (a) No-slip condition  (b) Slip condition  Figure 2.8 No-slip and slip conditions at the bottom wall  To better represent the motion of the phase in the digester, slip conditions are used at both the downcomer and digester walls. A generalized slip flow model [34] is used here,  'P  wall  sl  v  -L  wall  ~  v  (2.43)  ^ OX  wall  where Aw is the velocity jump at the wall, L is the slip length, w is the velocity of the s  c  chip phase. Thompson and Troian [35] have experimentally observed that the slip length is a function of shear rate in the liquid flow at the solid surface. By analogy, we assume that the slip length for the chip phase flow at the digester wall is a function of chip pressure and has the following form:  a  =i-  (2.44a)  +1  4=4  ^arcfan^J  2 ( 1  (2.44b)  71  where L is the characteristic length scale represented by the radius of the downcomer s  and the digester; P is the chip pressure on the corresponding wall boundary in Pa; A is c  w  the surface area between the chip phase and the wall in m ; a and (3 are dimensionless 2  parameters.  31  Chapter 2 Mathematical Models of the Pulping Process  From Equation (2.43) and (2.44), the slip velocity increases with the velocity gradient at the wall and decreases with the chip pressure at the wall.  2.6 Thermal properties of liquor 2.6.1 Density Liquor density depends on the solid concentration and temperature. A good approximation [21] for the density up to 50% dry solids is the following: (2.45)  p = 997 + 649X 25  p/p =l 25  .008 - 0.2377 /1000 -1.94(7 /1000)  2  (2.46)  where p is density of black liquor, kg/m ; T is the temperature, °C; and X is the dry solids 3  concentration, kg dry solids/kg total suspension. 2.6.2 Viscosity The viscosity of black liquor (mixture of white liquor and dissolved solid material) depends on several factors, particularly the temperature and solids content. The source of the black liquor is also important. The following equations provide an estimate of the kinematic viscosity [21]: (2.47)  \nu = A + BIT  (2.48)  A = -2.4273 + a X + a X + a X 2  x  3  2  3  73 = 6.1347 xlO +b X + b X 7  2  l  2  +b,X  3  where v is the kinematic viscosity, in mm /s; T is the temperature in degrees, K. 2  Constants a,,&, are shown in Table 2.3.  32  (2.49)  Chapter 2 Mathematical Models of the Pulping Process  Table 2.3 Coefficients for viscosity equations [21] Softwood  Hardwood  Tropical  a, a a  9.1578  3.3532  10.482  2  -56.723  3.7654  -54.046  3  72.666  6.  -42.178xl0  b  335.12xl0  h  -349.23xl0  2  -2.4907 -5.442x10  7  7  7  61.933 -40.165xl0  7  21.915xl0  7  300.55xl0  17.042xl0  7  -266.47x10  7  7  7  2.6.3 Specific heat capacity and thermal conductivity The following equations by Masse et al. [28] estimate the black liquor heat capacity, and thermal conductivity. c  = A.l 16(1 - X) + (1.675 + 3.31771000)X + (4.87 - 207 /1000)(1 - X)X  3  p  A = A (l-X)  + aX + bX  2  H20  (2.51)  <3 = 0.3176 + 2.268xl0" r  (2.52)  3  b =-1.394x 10" -3.069xl0" 7 2  (2.50)  (2.53)  3  where C is the black liquor specific heat, kJ/kg • ° C ; T is the temperature of the black p  liquor, °C; and X is the thermal conductivity, W/m°C. 2.6.4 Molecular diffusivity Development of continuous digesters has been focused on improvements in production capacity and pulp quality. The largest continuous digester currently in operation has a capacity of over 2000 tons per day. Obviously, the dimension of the digester increases with its capacity. With large dimensions, solid pressure/stress in the chip mass tends to be high, making the chip mass more compacted, and efficient diffusion becomes a key issue in practice. Therefore it is important to have good estimates of the diffusivity.  33  Chapter 2 Mathematical Models of the Pulping Process  To estimate the molecular diffusivity o f alkali, M c K i b b i n s [16] measured the diffusivity o f sodium i n kraft cooked chips for a variety o f temperatures and presented the following result: D =3Ax\0- T e2  U2  (2.54)  W70/(HT)  0  where Do is the molecular diffusivity, cm /min; T is the temperature, K ; R is the gas 2  constant (cal I mol • K).  Because M c K i b b i n s ' work was done with cooked chips, the  result needs to be corrected with respect to the p H value and lignin content. Taking these into account, Gustafson [14] proposed the following expression for diffusivity: D =5.7xlO" r e" 2  0  l / 2  4 8 7 0 / ( /  " ' ( - 0 . 0 2 Z + 0.13[(9//f )  55  +0.58)  (2.55)  where L is the lignin content in the over-dry wood, %; [OH] is the mole concentration o f hydroxyl in the liquor solution, mol/L.  34  Chapter 3 The Kraft Pulping Reaction Model As introduced in the previous chapter, the whole pulping process consists of many sub-processes, including reaction between lignin, carbohydrates and cooking liquor, diffusion of active components in the cooking liquor (mainly OH" and HS") through wood cells, associated mass and heat transfer and etc. In chapter 2, the governing equations of the pulping process are described. This chapter will focus on discussing the pulping kinetics.  3.1 Representative reaction models Numerous pulping reaction models have been documented in literature. These models fall into the following three categories. 3.1.1 Single variable model A classical kraft pulping reaction model is based on the so-called H-factor, which combines temperature and time into a single variable representing the extent of the cooking. The delignification is assumed to be one single reaction, and the reaction rate k has the Arrhenius form, k = A-e  EI(RT)  (3.1)  From some experimental data and setting a relative reaction rate k as 1 at 100°C, Vroom [11] gives, (3.2)  ln/t = 43.181 - 16113/T And, defines the H-factor = \k{t)dt ={(43.181-16113/7^  35  (3.3)  Chapter 3 The Kraft Pulping Reaction M o d e l  Once the H-factor versus yield relationship is determined for a particular set of pulping conditions, the yield for a given H-factor will be known. This simple model finds widespread use in industry, and the H-factor is an extremely important operation parameter in modern process control. However, the H-factor only provides an accurate fit for data with similar cooking chemical charges and wood compositions. It does not include other important factors, like the concentration of alkali or sulfide. In addition, it does not account for any change of chip species. Nevertheless, a single variable method still has its attraction. Due to its simplicity, it serves widely in the control scheme of modern digesters. There are several other similar methods: Edwards and Norberg [35] developed an extension of the H-factor called the T factor, which combines the effects of time, temperature, initial effective alkali and liquor to wood ratio. Kubes et al. [37] developed a G-factor that works for any alkali pulping method. Paulonis and Krishnagopalan [38] discussed a sophisticated process control system for kraft batch digesters that uses liquid density, electrical conductivity, and other techniques in addition to the usual parameters in order to predict kappa number despite variations in wood chip supply and digester operation.  3.1.2 Three phase model: Gustafson's model  In the second category of approaches, delignification reactions in kraft cooking are divided into three distinct phases: initial delignification, bulk delignification and residual delignification (Gullichsen [21]). The dissolving rates of lignin and carbohydrates are experimentally obtained for each phase. The following method describing the chemical kinetics has been proposed by Gustafson et al. [14]. Initial delignification: (Lignin>22.5%) dL  ^(17.5-8760/7^  dt  (3.4)  36  Chapter 3 The Kraft Pulping Reaction M o d e l  ^ = dt  (3.5)  k [OHf ^ dt U  lcV  J  ic  Bulk delignification: (2.2%<Lignin<22.5%) dL  dt  ioH-y+k e - °  ( 3 5 . 5 - 1 7 2 0 0 / ' /r '  ,  =k  e  K  (29A  >[0H-]*[HS-]  mo  tT  h2  A  (3.6)  L  b\  dC _ dt  dL  (3.7)  dt  bc  Residual delignification: (Lignin<2.2%) dL _ i (19.64-10804/77) Ke ]pH~f L dt 7  rl  dC _ , dt  (3.8)  dL dt  (3.9)  In the above equations, L is the lignin content, C is the carbohydrates content, [OH~] is the hydroxyl ion concentration, and [HS~] is the hydrosulfide ion concentration. The parameters k , k , k , k , k , k a  bx  b2  r/  ic  bc  and k  rc  are species related constants. For Pinus  Silvestris (lignin 27%, carbohydrates 67%, extractives and ash 6%), these constants are listed in Table 3.1. Table 3.1 Constants fitted to pulping reactions [14] Phase Initial Delignification  Equation  Constant  Value  (3.4)  k  1.0  b\  0.15  b2  1.65  (3.6)  k  (3.6)  k  Bulk Reactions  Carbohydrate Dissolution  a  Residual  (3.8)  Ki  2.2  Initial  (3.5)  k  ic  2.53  Bulk  (3.7)  k  bc  0.47  Residual  (3.9)  k  rc  2.19  37  Chapter 3 The Kraft Pulping Reaction M o d e l  3.1.3 "Two lignin" model Smith and Williams [39,40] proposed a comprehensive kinetic model for the conventional kraft process. This model, often referred to as the Purdue model, has been a basis for several other digester studies. In this model, the wood components are high- and low- reactive lignin (also called fast and slow lignin), cellulose and the hemi-cellulose components  galactoglucomannan  and araboxylan. The cooking liquor consists of  sulphide and effective alkali. The reaction rates of all these components are described by first order equations. This model considers only the initial and bulk stages of deliginification. One should also mention that the two classes of lignin have different reaction rates, making the calculations more complicated.  3.2 Motivation for deriving a new reaction model Unfortunately, all existing reaction models have some drawbacks. For instance, computations using the three-stage model given by Gustafson etc. [14] exhibit very good agreement with the experimental data. However, the calculation of alkali reaction rate requires the knowledge of the degradation rate of the acetyl group contained in wood chips and the transition points characterizing the reaction stages may change unexpectedly. Therefore, work is continuing on developing an improved kinetic model. Recognizing the effect of sulfide concentrations on the delignification rate, Aurora Santos et al. [41] established a delignification model where hydrosulfide concentration is included. This model is derived for Eucalyptus globulus; when species or operating conditions change, a so-called "stoichiometric coefficient" needs to be determined. Ling Yang et al. [42] presented a mechanistic pulping model where three key kinetic steps are considered: adsorption of hydroxide and hydrosulfide ions on the lignin; chemical reaction on the solid surface to produce degraded lignin products; and desorption of degradation products from the solid surface. This is physically more acceptable. However, the empirical constants in the model are also dependent on the operating conditions and need to be determined case by case, narrowing its application in some situations.  38  Chapter 3 The Kraft Pulping Reaction Model  Obviously, further CFD work is still needed to develop an easily applicable, more reliable and accurate pulping reaction model. From the computational modeling point of view, a good reaction model should give the reaction rate for each involved species, with only local chemical and physical conditions specified. In continuous cooking, this is more important because both temperature and chemical concentrations may exhibit significant spatial variation. A common drawback of most pulping reaction models is that they all need some modifications to fit pulping processes under different conditions. It is outside the scope of this research to address issues of fundamental chemical kinetics. Instead, we will take advantage of published work and use the available experimental data to derive an improved and easier to apply reaction model by regression analysis. As will be shown, this model has the capability of predicting pulping kinetics in various situations with good accuracy, although it does not give much physical insight.  3.3 Regression method The following considerations have been taken into account in deriving the general pulping model.  1. Take full advantage of the available experimental data So far, all the existing empirical pulping models were based on one set of specific experimental data. This is probably part of the reason why these models need appropriate modifications when used in other situations. An expanded database of experiments can better reflect the physics and lead to a more extensively applicable kinetic pulping model. In view of this, a relatively large number of data were used in the model regression, with cited data including Aurora Santos [41], Aurell et al. [44], Bugajer [45], Matthews [47], Kleinert [48], Wilder and Daleski [49,50] and Chiang et al. [51].  2. Select suitable regression variables Prior to discussion, a few definitions are introduced as follows. Residual lignin content is defined as the mass percentage of undissolved lignin in wood chips per unit mass of wood chips (oven-dry basis) in a macroscopic view. The 39  Chapter 3 The Kraft Pulping Reaction Model  mass fraction of residual lignin is defined as the percentage of residual lignin over the total lignin, which consists of dissolved and undissolved lignin. The residual carbohydroxide and its mass fraction are defined similarly. Delignification is a complex reaction process. As it proceeds, the dissolved solid material will gradually diffuse into the surrounding cooking liquor. The relatively slow diffusion on the other hand will result in a high concentration of dissolved lignin near the chip surface, which retards further removal of lignin from solid chip surfaces. Therefore, even though residual lignin contents are the same, the one with less dissolved lignin will host a faster lignin reaction. From this point, it is more appropriate to express the delignification rate in terms of the mass fraction of residual lignin. Likewise, the mass fraction of the residual carbohydrates has a closer connection to carbohydrate degradation. Therefore, the mass fractions of residual lignin and carbohydrates, instead of the traditionally used residual lignin and carbohydrate contents are chosen as regression variables. Hydrosulfide concentration is another regression variable as its influence on delignification was documented by Bugajer [45]. 3. Avoid using global variables In some pulping models, global variables such as liquor to wood ratio were used. However, chemical reactions always occur locally; therefore, they are dominated by the local temperature and chemical concentrations. In large industrial digesters, both temperature and chemical concentrations vary from place to place; the fixed quantity liquor to wood ratio may not be a good representation. Further, in continuous cooking, liquor to wood ratio is hardly used as a parameter. 4. Provide sufficient data to ensure model accuracy In normal operating conditions, delignification proceeds quite rapidly/efficiently in the bulk stage but slowly in initial and residual stages. This makes the lignin content versus time a very sharp curve in the middle of the range, while quite flat on both ends.  40  Chapter 3 The Kraft Pulping Reaction M o d e l  To characterize this feature in regression, more experimental data are needed, but in practice, only limited measurements are taken. It is therefore very hard to obtain an accurate regression model. Because of this difficulty, all the experimental data to be used in regression are fitted by piece-wise functions, from which a desired number of data are then extracted. This way, a reasonably large number of data are generated for each experiment. Thus the whole database is expanded. This strategy turned out to be very helpful in leading to a higher accuracy regression model. 5. Abandon the "three stage" regression form  An extensive review given by Santos [41] showed that the majority of the pulping models divide the delignification process into three distinct stages. Accordingly, three different expressions are used to calculate the delignification rates. In these "three stage" models, the transitional points are considered to be constant. However experiments indicate that these transitional points may change in different cooking conditions. In this work, the "three stage" idea will not be used in order to eliminate the need for transitional points.  3.4 Regression procedures and the pulping model In typical cooking conditions, Santos [41] has experimentally showed that both the effective alkali and sulfide concentrations appear as only functions of the mass fraction of residual lignin content, as shown in Figure 3.1. If we check some previous pulping models, the carbohydrate degradation rate can also be expressed solely by the rate of lignin removal. This is quite encouraging. It clearly indicates that once the delignification rate is known, the calculation of the concentrations of other important species becomes fairly simple. The accuracy of the model prediction therefore depends to a large degree on the computation of the deliginification rate. These experimental results also imply that  41  Chapter 3 The Kraft Pulping Reaction Model  the mass fraction of residual lignin is a better variable than the residual lignin for characterizing the reaction process as mentioned earlier. Functions of effective alkali and sulfide concentrations related to mass fraction of residual lignin are obtained using Santos's experimental data. They can be well represented by polynomials: [OH] = 0.4429 + 3.262Z -17.19L) +42.3SL -47.49/4 3  r  r  (3.10)  + 20.04 +{OH\ -1-40  (3.11)  [SH] = 0.04765 + 0.15834 -0.3236Z . +0.2616Z? + [SH~\ -0.144 2  where L = LI L , represents the mass fraction of the residual lignin. r  0  The overall F statistics are significant, with F=4.04, p<0.00035 and R-square=0.918 for equation (3.10) and F=8.23, /?<0.00004 and R-square =0.827 for equation (3.11). Calculated concentration profiles are shown in Figure 3.2.  1.20 \-  V x • O * • A O + t>  Exp. Case 1 Exp. Case 2 Exp. Case 3 Exp. Case 4 Exp. Case 5 Exp. Case 6 Exp. Case 7 Exp. Case 8 Exp. Case 13 Exp. Case 14  0.16 • A  4 _l v  x > •  o E I  C* *<&  1  p.  0.20,  o  *A A*  0.08 h  qp -  0.60 fc  0.12 h  Exp. Case 1 Exp. Case 2 Exp. Case 3 Exp. Case 4 Exp. Case 5 Exp. Case 6 Exp. Case 7 Exp. Case 8  x+  0.2  0.0  0.4  _l I I L-J I I0.6 I I I L-  0.2  Lr  I  0.4  I  Lr  (a) Effective Alkali  (b) Sulfide  Figure 3.1 Effective alkali and sulfide consumptions versus lignin removal (liquor to wood ratio h=5) [41]  42  Chapter 3 The Kraft Pulping Reaction M o d e l  1.60 Exp. Case 1 Exp. Case 2 Exp. Case 3 Exp. Case 4 Exp. Case 5 Exp. Case 6 Exp. Case 7 Exp. Case 8 Exp. Case 13 Exp. Case 14 Polynomial fit R =0.918  1.40  1.20  !  0.12  !  E  I  Exp. Case 1 Exp. Case 2 Exp. Case 3 Exp. Case 4 Exp. Case 5 Exp. Case 6 Exp. Case 7 Exp. Case 8 Polynomial fit, R =0.827  I  0.80  0.40  0.6  0.8  1.0  Lr  (a) Effective Alkali  (b) Sulfide  Figure 3.2 Fitted effective alkali and sulfide profiles [41] After developing the correlations between effective alkali, sulfide and the mass fraction of residual lignin, efforts were directed towards obtaining the delignification rate. From theory of chemical kinetics [45], we assume the delignification rate can be given in the following form: dL dt  =  (3.12)  -e ^ [OHf[SHfL c  ir  Or equivalently  In  dL  = C -jC  r  0  dt  + \n[OH\C  y  + \n[SH]-C + ln(4) • C  2  When sufficient data sets of ln provided, the constants C ,  3  dL  r  dt  ~  4  (3.13)  , m[OH], \n[SH] and ln(Z ) are r  i- 0,1,2,3,4 can be determined by multi-variable linear  t  regression. In typical cooking experiments, the lignin content and temperature at different times are measured, but measurements are usually not taken that often. As a remedy, the continuous form of function L (i) in each individual cook can be fitted with the measured r  data, then ln  dL  r  dt  and ln(Z )can be calculated at any desired time, with ln[0//] and r  43  Chapter 3 The Kraft Pulping Reaction Model  ln[57/] being determined similarly. As the delignification rate can be affected substantially by heating conditions, cited experimental data are classified as normal heating and fast heating data. Although there is no specific rule, we suggest cases with heating time less than 80 minutes belong to the fast heating group, while the rest belong to the normal heating group. Therefore, data from Bugajer [45] and Aurora Santos [41] fall into the first category; data from Aurell [44], Mathews [47] and Wilder belong to the latter one, and corresponding conditions are listed on Table 3.2. Accordingly, two different regression expressions are obtained for the fast heating and normal heating conditions respectively. dL  17.19-7201 /7  r  dt dl\ dt  =  _  e  '[O/7]- [5/7]  8 6„-5663/r  ,041  ^10  2075  Z°  (Fast heating condition)  925  373 ^ ^ 0 , 4 8 ^ 0 . 7 6 ,  (3.14) (3.15)  c o n d i t i ( m )  R-squares for these two expressions are 0.927 and 0.965, indicating successful regressions. Table 3.2 Cooking conditions of cited experimental data  3  Reference  Case L  0  Co  %  %  E A as [NaOH] mol/L  Sulfidity as [SH], mol/L  L / W [NaOH] [Na S] t tc min min g/L g/L  7 °C  2  h  C  Aurell 1965  1  27.3  67.7  0.937  0.134  4  32.14  10.44  120  120  170  Aurell, 1965  2  27.3  67.7  1.562  0.223  4  53.57  17.40  120  120  170  Matthews, 1974  27.2  66.5  0.92  0.114  5  32.25  8.86  90  90  170  Wilder, 1965  28.6  71.4  1.50  0.375  200  45.0  24.0  -  205  150  Bugajer, 1984 1  28.77 71.23 1.129  0.161  4  38.71  12.58  65  100  159  Bugajer, 1984 2  28.77 71.23 1.129  0.161  4  38.71  12.58  65  100  168 150,160,  Santos, 1997  C l - 4 22.1  73.6  1.20  0.18  5  40.80  14.04  45  150  Santos, 1997  C9  22.1  73.6  1.8  0.18  5  64.80  14.04  45  150  160  Santos, 1997  C10 22.1  73.6  0.7  0.18  5  20.80  14.04  45  150  160  Santos, 1997  C13 22.1  73.6  1.2  0.43  5  30.79  33.54  45  150  160  Santos, 1997  C14 22.1  73.6  1.2  0.07  5  45.20  5.46  45  150  160  26.8  70.0  1.193  0.211  5  39.31  16.42  100  120  170  C l - 3 28.0  72.0  1.226  0.216  10  40.38  16.87  -  60  170,175,180  Chiang etc., 1988  170,180  b  Kleinert, 1966 b  Note: " 7c, cooking temperature; t heating time; t cooking time. Data not used in the model regression, but will be used to check the model prediction. h  z  b  44  Chapter 3 The Kraft Pulping Reaction Model  Unfortunately, experimental data of carbohydrates degradation are not well documented. The degradation rate of carbohydrates given by Aurora [41] is used directly here: 0.25^  dC  while  L < 0.4  while  L > 0.4  r  r  — - =<  dt  dl  -  ( 3  0.14—^ dt  1 6 )  From Equation (3.10), the reaction rate of effective alkali (EA) can be further derived. Because the mole concentrations of hydroxyl and sulfide appear in the expression of lignin reaction rate (3.14) and (3.15), the EA reaction rate in its mole concentration form is considered here. To be consistent with Santo's expression (3.10), EA is defined by, EA = (NaOH + 0.5Na S) 2  as  g/L  N a 0 H  (3.17)  Note: In North America, EA is defined as (NaOH + 0.5Na S) 2  as  N a l o  g/L.  The mole concentration of EA in the free liquor is then, FA  FA  [EA] = -=^— = — mol/L W 40  (3.18)  Also, the mole concentration of EA is essentially the mole fraction of [OH] with the definitions of (3.17) and (3.18), i.e., [EA] = [OH]  (  3  1  9  )  Manipulations of Equations (3.10) and (3.19) leads to, d\EA\ dt  =  dyOH] dt  =  /  _  2  +  V  _  1  Q  0  Z  3  +  1  0  0  r  Z  4  ) ^ r J  ( 3  21)  dt  Similarly, if we define, Sulf = Sulfidity =  Na S 2  NaOH + Na S 2  ax  (3.22)  NaOH  The derived reaction rate of sulfidity will be, dSu£ dt  =  d[HS] dt  f  +  8  ^  ( 3  ' dt  v  45  2  3  )  Chapter 3 The Kraft Pulping Reaction Model  3.5 Model validation With Equations 3.12, 3.16, 3.21 and 3.23, all of the four species L , C , OH and SH r  r  can be computationally determined. Naturally, the accuracy of the model prediction is a major concern. In view of this, a computational code was developed to simulate the batch cooking experiments. A four-stage Runge-Kutta time advance scheme [43] was used to solve the four differential equations 3.12, 3.16, 3.21 and 3.23. Time steps are reasonably small and vary case by case; measured temperatures are directly used in the computation. As discussed before, the prediction of lignin content is of great importance to the model accuracy. Also, concentration measurements of the other species are not available, and therefore only the lignin content will be checked here. Figure 3.3, 3.4 and 3.5 show the variation of lignin content as time proceeds. A l l lines represent model predictions, while symbols represent experimental measurements, by Aurell et al. [44] and Bugajer [45] in figure 3.3 and Wilder [49,50], Matthews [47] in figure 3.4. As can be seen, fairly good agreement is achieved generally, but the model underestimates Aurell's experiments in the later cooking period. Nevertheless, the biggest predicting error of lignin content on wood is about 11%. 30 •  Aurell case 1  •  Aurell case 2  •  Bugajer case 1  O  Bugajer case 2 Prediction A.1 Prediction A.2 Prediction B.1 Prediction B.2  200  250  :300  time, min  Figure 3.3 Lignin profiles by model predictions and experiments Experiments by Aurell et al. [44] and Bugajer [45] 46  Chapter 3 The Kraft Pulping Reaction Model  35 Mathews Wilder and Daleski Prediction Mathews Prediction Wilder  — TJ O O  c o c  CI)  _!  i  i  i  I 50  i  i  i  i  L 100  150  200  time, min  Figure 3.4 Lignin profiles by model predictions and experiments Experiments by Wilder [49,50] and Matthews [47]  Kleinert170°C Kleinert175°C Kleinert180°C Chiang etal. Prediction L.170 Prediction L.175 Prediction L.180 Prediction Chiang  • o A  •a o o  5  N \ A \  \  c 'E  \ A \  \ A \  N A ^ ,  100  150  time, min  Figure 3.5 Model prediction versus experiments by Kleinert [48] and Chiang et al. [51] In Figure 3.3 and 3.4, however, all the experimental data have been used in the model regression. To further check the quality of the model predictions, some other cases 47  Chapter 3 The Kraft Pulping Reaction Model  not previously used in the model regression were modeled and the simulated results are given in Figure 3.5. Again, the model gives very good prediction to Chiang's experiments, with the error mostly limited to 1%. For experiments done by Kleinert [48], the simulations appear promising in the early cooking stage, but under-estimate the measurements thereafter, with the biggest errors in the range of 2-3% occurring in the bulk cooking period. As listed in Table 1, Kleinert's experiments were conducted with high liquid-to-wood ratio and relatively higher cooking temperature, this is probably the reason why model predictions give relatively larger errors.  3.6 Conclusion Based on extensive experimental data, a kinetic pulping model is obtained exclusively by mathematical regression. In this model, two different expressions are used to calculate the delignification rate respectively for normal heating and fast heating conditions; Santos equation of carbohydrates degradation rate is directly employed because of lack of experimental data.  By introducing the mass fractions of lignin and carbohydrates as regression variables, better utilizing the experimental results and choosing suitable regression form, the resulting pulping model predicts cooking experiments conducted under various conditions by different investigators, without the need to determine some operating conditions related model constants. In addition, the model is not based on the "three cooking stages", thus no transitional points are needed in the computation. Simulations using this model exhibit good agreement with the experimental data, with prediction errors less than 11% for all the tested cases. It is also expected that the model can be used in the simulation of continuous cooking.  48  Chapter 4 Computational Methods In Chapter 2, the governing equations describing the two-phase flow, energy conservation and chemical reactions in the cooking process have been introduced. In this chapter, the general discretization schemes and computational procedures will be presented. As discussed in previous chapters, the solid chip phase is treated as a fluid, thus Navier-Stokes equations apply for both phases. Instead of giving the detailed discretization techniques for each individual equation, we will first give the approach to solve the single phase Navier-Stokes equations, and then address some additional issues in order to model the two phase flow problem. Examples of imposing typical boundary conditions are also provided in the later part of this chapter.  4.1 Solving Navier-Stokes equations by SIMPLE method In compressible flows, the continuity equation can be used to determine the density and the pressure is calculated from an equation of state. Unfortunately, this approach is generally not appropriate for incompressible flows or low Mach number flows, as density remains constant. By adding a time derivative of pressure into the continuity equation, the compressible computation method can be used in incompressible situations. This is called artificial compressibility method which was proposed by Chorin [52]. Usually, this method is used in conjunction with the collocated variable arrangement, and the resulting pressure oscillations can be lowered numerically. However, for fast convergence, this method requires a good compressibility parameter, which varies significantly from case to case. Selecting such a parameter is not always easy, and that proves to be a major deficiency of the artificial compressibility method. Methods based on the SIMPLE algorithm [53] are fairly efficient for solving both steady and unsteady problems, and this kind of method has been used widely in incompressible flow computations. There are  49  Chapter 4 Computational Methods  many variants of this method available; the SIMPLE and SIMPLEC methods have been implemented in the code. In this section, the basic SIMPLE method will be explained. 4.1.1 Navier-Stokes equations We consider a 2-D, steady, laminar, incompressible flows of a Newtonian fluid. In a Cartesian coordinate system, the conservation laws for mass and momentum can be written as: 8(pu)  |  d(pv) d(pw)  dx  |  dy  =  (4.1)  Q  dz  dz d \ d , \ d , x dP — (puu) + (pvu) + — (pwu) = - — dx dy dz dx y dx  dv  rr  (  (  (  dx  dy  KJJ\,  dz„„ yy  r  (  dz  dz  •+  dx  dx  \JJ\,  ^dz „xy d \ d \ d N dP — (puv) + — (pvv) + — (pwv)=- — dy dx dy dz dy  dz,  yx  —  + -  dz„  dy  KJJ\.  yz  dz  dz  P  |  (4.2a)  X  J  dzzy ^  + -  dz  |  + B  (4.2b)  + pB  v  dy dz. dz  + pB  z  (4.2c)  For incompressible flow, we have,  XX  *  r =2p yy  Vr-  = 2//  du o dx  du dv^ xv  5v  X  dy  — Z vz  A*  vx  zv  dy dx  dv dw =u dz dy r*  dw  'du  dw  dz  dx  *:x=*x:=P  ~dz  y  s  (4.3)  j  With these relations, the momentum equation can be arranged into its velocity form: d ~ (puu) + -|- (pvu) + -|- {pwu) = - ^- + p dx dx dy dz dx  du  dv du  dx  dy dx dy  f d , \ d i \ d \ dP du dv — (puv) + — (pvv)+—(pwv)=- — + p + dy • dx dy dx dx dy dz dy (  50  dw ^du dz V dx  dv_  dz j  dv dw  + •dz dz dy  + pB  y  Chapter 4 Computational Methods  dw d \ d , \ d i x dP — (puw) + — (pvw) + —(pww)=- — + p ox dy dz dz dx dx  du  (  —  +  dv  + dz dy dz —  •  dw dy  +  dz V dz J  (4.4) After some mathematical manipulation, the momentum equations can be further arranged into the following uniform form: V -{pvt- pV</>)=-VP + S,  (4.5a)  In Cartesian coordinates, (4.5a) can be written as, ^  d  +—  PV0-JU  dx  dx  x • = x, y, z;  (/> = u, v, w  d0_  dy  dz  pwcj)  - JU  d£ dz  dP_ dx,  (4.5b)  where t  r  du dv' dw + pB M — +• M — +• dx dz dx dx v dxj dy d du^ dv) d_ dw M — +• M — + pB dz dx v y J dy dy v dy j d ( dv^ d_ d ( du^ + pB M — + dz V dz dx ^ dz j dy V dz  x  when 0 = u  y  when 0 = v  z  when <j> = w  f  (4.6)  d  In case of constant viscosity p, the viscous terms in S4, are eliminated by continuity. 4.1.2 Collocated grids While solving the Navier-Stokes equations, one can store the velocity variables and pressure at the same grid points. Such a grid system is called a collocated arrangement. As many of the terms in the transport equations are identical, the number of discretization coefficients is minimized and the programming is simplified. For problems with more variables or a more complex geometry, or when the multigrid approach is used, the advantages of a collocated arrangement become more apparent. However, when used in incompressible flows, the collocated arrangement has difficulty handling the pressurevelocity coupling resulting in pressure oscillations. On the other hand, the staggered arrangement, first introduced by Harlow and Welsh [54], offers a strong coupling 51  +  Chapter 4 Computational Methods  between the velocities and the pressure, and thus successfully solves the incompressible Navier-Stokes equations without oscillations. It is for this reason that the collocated arrangement was seldom used in the early days of CFD development. In the 1980s, improved algorithms were developed to solve the pressure and velocity equations, and the collocated arrangement began to win popularity. Today, it is adopted in many commercial codes, such as FLUENT. In the present code, the collocated variable arrangement is used based on the scheme of Ferziger and Peric [55,56]. In collocated grids, only one set of control volumes is required, the computational points (storage locations for variables) are then located in the center of the control volumes, an example of such arrangement in a 2-dimensional case is given in Figure 4.1.  o n o W o  lii il•  o  °N m  A  w ;p e s ; se  ii  *S  o E o  Figure 4.1 Collocated grids in 2-D situation In the collocated arrangement, careful consideration should be given to the pressure gradient calculation. Ferziger and Peric [55] showed that if pressure derivatives at cell faces are evaluated by interpolation of cell center derivatives, unrealistic pressure distributions could be produced. The main reason is the use of the 2Ax approximation for the first derivatives. To cure this problem, they suggested an effective scheme where pressure derivatives at cell faces are evaluated by central differences and Ax spacing, thus introducing the need for interpolation. 4.1.3 Discretization of the equations by control volume method  52  Chapter 4 Computational Methods  For simplicity, the discretization scheme is described below only for a twodimensional situation, although the code is developed three-dimensionally. Finite Volume Method starts from the conservation equation in integral form, based on equation (4.5a). The integral form of the generic conservation equation can be written as, jpv</> • MS = jfNj s s  • MS + JVPJQ + J S ^ Q n n  (4.7)  Integrating the equations (4.7) over the control volume leads to the following discrete equation: dP dx,  f  4 - 4 + 4-4=  j AV  ^  (4.8)  W =S  AV  4  where the I's represent the average fluxes through the control volume faces in the direction of increasing coordinates, and SAV is the volume integral of the right hand side of Equation (4.5a). The terms appearing in Equation (4.8) have to be evaluated differently for the mass and momentum conservation equations; this will be outlined below. For convenience, we define C = ( UA )  C„ = ( UA I  C =(pVA )  C =(pVA \  x  P  y  X  P  y  yn  C„ = (pUA )  X  v  x  C =(pVA \  y  ys  y  (4.9)  w  (4.10)  Performing integration over the control volume for equations (4.5b) yields the discrete equations. 4.1.3.1 Continuity equation In the continuity equation (4.1), the fluxes are the mass flow rates through the control volume faces, and are evaluated as: I =C e  xe  = {pUA \  I =C  x  n  yn  = (pVA \ y  (4.11)  The discrete continuity equation is then pU A e  x  - UA P  W  X  + pV A - pV A =0 n  y  s  y  4.1.3.2 Momentum equations  53  (4.12)  Chapter 4 Computational Methods  For simplicity, the convection fluxes are approximated with upwind difference scheme (UDS), while the diffusive part is estimated by central difference scheme: I" = r U - ju — (U -U ) V Sx E  I" =c  f  u  A^  P  (U -U ) N  8y  P  =  C U -D (U -U )  (4.13)  =  C U -D (U -U )  (4.14)  xe  yn  e  e  n  E  n  P  N  P  The pressure gradients are also evaluated with central difference approximations. C U -D (U -U )]-[C^U -D (U -U )] xe  e  e  E  P  W  W  P  +  IV  C U -D (U -U )]-[C U -D (U -U )] yn  n  n  N  P  ys  s  s  P  =  s  -A ^P -P ) x  E  (4.15)  w  where cell face velocities are evaluated by upwind scheme, for instance, \U ,  if  C >0  u  if  C <0  P  u„ =  (4.16)  xe  xe  After some rearrangements and applying Successive Over-Relaxation (SOR) scheme [58], the discrete U-momentum equation can be written as: aU p  = Ia„ £/„* + (1 - a )U°  P  6  u  - A , (P  P  x P  - PJ  e  (4.17)  nb  where the index "nb" runs over all neighboring points E, W, N and S; (4.18a)  <*E = A + maxC-C^.O) a  (4.18b)  =D +max{C ,0)  w  w  xw  ®N = n + max(-C ,0)  (4.18c)  a = D + max(C ,0)  (4.18d)  r=-rYu n  (4.18e)  D  >  s  s  a  >  a  nb  where oc is the over-relaxation parameter, having a typical value of 0.7. u  Similarly, the V-momentum equation is, aV p  P  =I a V nb  nh  +{l-a  u  )V" - A , (P - P ) P  y P  n  s  nb  The coefficients in equation 4.19 are similar with 4.17, as listed in formula 4.18.  54  (4.19)  Chapter 4 Computational Methods  4.1.3.3 Pressure-velocity-coupling The U - and V-momentum equations are solved using guessed values for the pressure and mass fluxes. The velocity components U* and V* calculated with these guessed values will in general not satisfy the continuity equation, thus: P  U  - PU A  e A  W  (4.20)  + pV A - pV A = S  X  n  y  s  y  m  The discretized U-momentum equations for nodes P and E can be written as: u;=  — X a U* +(l-a  1  U =  nb  I nb  dp  E  nb  nb  U +U  P  (4.21a)  )U° - A , (P* - P*) =  U +U  (4.21b)  P  X a U* +(l-a nb  (P* - P *) =  )U° - A  u  u  XtP  P  e  X E  P  w  E  e  E  nb  The cell face value is obtained by linear interpolation of the terms in equations 4.2 with the exception of the pressure difference, which is evaluated as in the staggered approach. This yields:  ^ l a  A JP -P ) =U U  2Za U +(l-a )U nb  \ p  1 ^  f  u  j  nb  x  u  E  P e  e+  e  (4.22)  nb  With equation 4.21, we have 1 ^  f  U] \  A  P  ' 1 ^  X^t/;+(i-«„)t/° nb  JE  \  j£  1 ^  r  U  Q  P  1 ^  f  (4.23a)  J  f  A p(P x  \  A  P  Jplnb  T^bKb+a-ccM  So,  KP A  a.  a.  V  f  A  Pw)p  (4.23b)  Jp.  (4.24a)  u:  U = average  K  e  Jp  P  J  dx \oxj  \Sxj  =U +  A~x 4  A  A  K p a  p  E  55  A v  J  f  Sp^ Sx  Chapter 4 Computational Methods  Ax-A ^  'Sp)  r  U„ =U„  x  (Ps-Pp)e=u:-  (Sp  = U +U e  v  P  e  (4.24b)  J  where the overbar indicates linear interpolation. Neglecting the interpolated pressure gradient in equation 4.22 and 4.24b, the velocity correction U is linked to pressure by the so-called velocity correction equation, e  Ax • A V  a  P  Ax • A,  Sp  r  \SxJe  J  \  Sp_  ' 1  P J  Q  \p a  N  j  (4.25)  Substituting equation 4.24b into equation 4.20, and imposing mass conservation, results in a pressure correction equation: a,,P  = Y u  P  a  ^  b - S  p n  (4.26)  m  nb  with the coefficients defined by, a  E  = pA  ' 1^  x,e  \  a  P  w = PA  a  KP Q  J  f 1 >  a = pA y,s  N = P y,n  A  s  \  Q  P  , A  r  J  f  1  (4.27)  J  A  A  (4.28)  S  (4.29)  a =Y^ nb a  p  nb  As can be seen, the SIMPLE algorithm with collocated arrangements is very similar to the staggered approach, but requires some interpolation. After discretization, the differential equations are converted to a system of algebraic equations. For structured grids, the resulting algebraic equations always have very sparse coefficient matrices and they can be efficiently solved by some iteration methods, such as strongly implicit procedure (SIP) proposed by Stone [56] and alternating direction implicit (ADI) method introduced by Hageman and Young [58].  56  Chapter 4 Computational Methods  4.1.4 Solution procedures for steady-state problems 1. Start the calculation of the fields at new iteration step using the latest solution u" and p" as starting estimates for u" and p"^ . +l  2. Assemble and solve the algebraic equation systems for the velocities to obtain u . j  3. Assemble and solve the pressure-correction equation to obtain p . 4. Correct the velocities and pressure to obtain the velocity field w,, which satisfies the continuity equation, and the new pressure p 5. Return to step 1 and repeat, using updated w, and p as improved estimates for u"  +x  and  , until all corrections are negligibly small.  4.1.5 Curvilinear grids In case of curvilinear grids shown in Figure 4.2, the discretization schemes and the approach solving Navier-Stokes equation introduced above still apply, but some modifications are needed when evaluating the fluxes, pressure term, and some source terms. The discretization of corresponding terms in the liquor phase equations is outlined in this section.  Figure 4.2 Two-dimensional curvilinear grids 4.1.5.1 Mass flux in continuity equation On the east face of a control volume, the mass flux is evaluated as:  57  Chapter 4 Computational Methods  (4.30)  m = pV • n dS e  }  e  The mid-point rule approximation of the above integral leads to, m « (pV •n) A = (puA \ e  e  e  (4.31)  + (pvA \  x  y  4.1.5.2 Convectivefluxesin momentum equation Having evaluated the mass flux in continuity equation (4.1), we can further evaluate the convective flux term in the momentum equation (4.5a). On the east face of a control volume, the convective flux integral is, F = \p</>V • MS « (p<t>V • n\A  (4.32)  = mj  c  e  e  e  The computation of F° varies in different discretization schemes. For example, in explicit central difference scheme (CDS), m and <j) are estimated by interpolation with e  e  the values from the last iteration; in implicit upwind difference scheme (UDS), however, Ff is formulated as, c,implicit  max(m 0)u  F.  e9  p  (4.33)  +min(m ,0)u e  E  4.1.5.3 Diffusivefluxesin momentum equation As an illustration, we will again evaluate the diffusive flux integral on the east face of a control volume. By using central difference scheme (CDS), the explicit and implicit diffusive flux integral are calculated as, pd^ucu  j s70.MS*p  =  M  ±A dx  +  d e  x  (4.34)  ±A,  d  dy  x  y  % -nA ={pA) (§ e  Se  V  EP  J  L  58  e  (4.35) ?J-n<t>E-<t>F L EP  J  Chapter 4 Computational Methods  When the line connecting nodes P and E is orthogonal to the cell face, E, • n = 1; if not, this discrepancy will be corrected by deferred correction, which will be described shortly. 4.1.5.4 Approximation of pressure and source terms The pressure can be treated as conservative forces on CV surfaces, applied with Gauss theory, the pressure integral term can be calculated as, F  = - J — dn = $P-ndS*Y c4? n s c  p  > (x,=x,y,z)  P  (4.36)  On the other hand, the pressure can also be regarded as a body force, so we can compute it another way, 8P)  r  F =-\^dQ p  AQ  (4.37)  J Ply n i dx  Generally speaking, the first approach is fully conservative, but the second does not always remain so. The source term is basically a volume integral. The midpoint rule is commonly used to approximate it, which approximately gives second order of accuracy. F '=-{5^Q«^ AQ  (4.38)  s  p  4.1.5.5 Deferred correction scheme As discussed earlier the adoption of collocated grids produces an oscillatory solution that is not real. A remedy for this is the so-called deferred correction scheme with appropriate discretization schemes; in deferred correction schemes, both convective and diffusive fluxes are computed in the following way: F =F  imp  + [F  exp  - F'"' ] p  (4-39)  M  59  Chapter 4 Computational Methods  where 'imp' and 'exp' denote the implicit and the explicit flux approximation, 'old' means value from previous iteration. In Muzaferija's approach [59], all the fluxes are evaluated by CDS scheme except the implicit convective flux, which is computed by UDS scheme. Then the approximation of convective and diffusive flux integrals at the east face of a control volume is written as: F° = max(m ,O)0 + min^,0)^ . + \m </> -max(m ,0)^ e  P  /;  e  e  c  J  + ("4  e  (4.40)  E  old  d^  d£  (  EP  -min(m ,0)(/) ]"  p  (4.41)  \  = ( 4 V</>E-0 I'EP J  4.2 Additional treatment of two-phaseflowequations 4.2.1 Evaluate the mass flux In continuous digesters, the free liquor and chip phases occupy each individual control volume and flow through the same control volume faces. To evaluate the mass fluxes through the control volume faces in the continuity equations (2.7,2.8), the volume fraction of each phase has to be taken into account; for instance, the liquor phase mass flux through the east face of a control volume is approximated as: m « {pjSjV • n)A e  e  = (p e u A ) f  f  f  +  x  (4.42)  (p e v A ) f  f  f  y  4.2.2 Discretization of the inter-phase interaction terms In the presence of wood chips, the liquor phase momentum equation (see equation 2.11) has two extra terms, representing the momentum exchanges caused by the interphase friction force and mass transfer respectively. ' J (pf f fj) £ u  t  +  (pf f f,, fj) £ u  u  d  =f  P  f  +  dr^  £  dx  v  J J  60  + fPfgi E  +  F fJ  +mu n  cj  (4.43)  Chapter 4 Computational Methods  where m is determined by chemical kinetics andF . is given in Equation (2.15), u  ;  s  f  ,  2  4 —  v  + 2 c\ cJ R  4  £  U  ,  (4.44)  -"/,/  The integral of these two terms over the control volume can be evaluated linearly, which gives, \(F/j + m u )dQ ]2  « {F  ci  + m u ) AQ  fJ  n  ci  p  = (S  p  •u  fJ  + S )AQ 0  (4.45)  where 1(  sc 2  V  f  £  1  (  (4.46)  f  £  sc 2  R\  v  2 c "cy  (4.47)  "/,/  c  f  f  £  £  The extra terms in the momentum equations of the solid chip phase are treated similarly.  4.3 Discretization of scalar equations The energy and species equations are the scalar equations in our modeling. All the scalar equations are solved after the flow fields for both phases are obtained. The discretization methods of these scalar equations are similar; as an example, the discretization of the steady state energy equation is outlined below. Referring to equation 2.17, and neglecting the heat generated by work and viscous dissipation, the steady state energy equation becomes, A  f  dx  -(Pf /J ) U  T  =  8X :  dT  \ PJ C  +  Q  (4.48)  .iJ  dx  Using the deferred correction scheme to approximate both left hand side convective and right hand side diffusive term gives,  61  Chapter 4 Computational Methods  p c  _  pcjmp  p i t  _  pdjmp  _|_ j ^ c . e x p  _  _|_ ^ ( / , e x p  (4.49)  pc,impl^  U  _  pdjmp  (4.50)  j°  On the east face of a control volume, fluxes evaluated by explicit and implicit schemes are, F " c  exp  e  =  F d.exp  mX  _  ex  A  ydx j  PJ  C  F°  (4.51)  xp  = max(m ,0)Tp  Jmp  e  E  r \ dT  Ax  EP  r C  p,f  L  e  (4.53)  +mm(m ,0)T  e  djrnp  F  (4.52)  A ey  +  +  Ay EP  D -(T -T ) e  E  P  (4.55)  A  *f  D„ = °P,f  L  (4.54)  E P  -  E P  Similar work is conducted on all other control surfaces, and the following discrete equation is finally obtained.  flpT^X^+S,  (4.56)  nb  where the index "nb" is taken over all neighboring points E, W, N, S. a  E  = D + max(-m c  c  (4.57a)  ,0)  a = D +max(/n ,0)  (4.57b)  a = D + max(-m„,0)  (4.57c)  a  (4.57d)  w  w  N  s  e  n  = D + max(m ,0) s  n  (4.57e) nb  62  Chapter 4 Computational Methods  The approximation to the reaction heat term using midpoint rule yields, (4.58)  AQ V PJ)P C  4.4 Imposing boundary conditions Typical boundary conditions in CFD computations include wall boundary, inflow and outflow boundary conditions, gradient and symmetrical boundary conditions. Although it is usually not hard to implement these boundary conditions, careful consideration is necessary in some cases. As an example, the way to impose symmetrical boundary conditions for velocity equations is given in the following. Suppose the control volume surface OAB in figure 4.3 is the symmetrical plane of the velocity, then the net mass flow through plane OAB should be zero, further, the velocity gradient along the normal direction of plane OAB should be also zero, i.e., (4.59) dV dn ,  (4.60)  =0  b  n  C j  V (u ,v ,w ) h  h  h  h  A  B Figure 4.3 Illustration of symmetrical boundary conditions The norm of plane OAB can be computed by n = OC = OA x OB  63  Chapter 4 Computational Methods  The second equation can be satisfied by the following two expressions: V„ • OA = K • OA  (4.61)  V • OB = V, • OB  (4.62)  h  Where Vb is the velocity at the boundary surface and Vj is the velocity at the center of the boundary cell. The above equations can be written in matrix form as shown below. "OA  y  'OB  "OB  \ oc x  ft.  "OA  0 A  yoc  OC  Z  \  OA  yOA  Z  OB  yOB  Z  0  0  X  X  J  OA  OB  (4.63)  0  Thus the velocity at the boundary can be calculated.  4.5 The code and its structure The code was written in Microsoft C, with the capability of solving either single- or twophase flow, two or three-dimensional problems. It uses collocated mesh with curvilinear grids and has a second order accuracy. As shown in the flow chart on Figure 4.4, the code consists of three parts, namely "Grid", "Compute" and "Plot". The function of "Grid" is to generate user-defined mesh; as the main part of the code, "Compute" reads the pre-generated mesh, satisfies the userrequired boundary conditions and delivers the converged solution; "Plot" functions as a post-processor, which reads the mesh and solution to produce data files recognizable by TecPlot.  64  Chapter 4 Computational Methods  (  Grid  Compute^)  )  Geometry definition Generate mesh  I Physics definition Single/Two phase flow, Newtonian/non-Newtonian fluid  (  Plot  )  Plot definition  H Read mesh file >H  Read mesh file  T  Save mesh file  Set boundary conditions  ( End )  Set iteration parameters and convergence rule Estimate initial solution  H Read solution Generate data files in specific formats  Parameters check Solve flow field Liquid phase N-S equation Solid phase N-S equation  Solve energy & species conservation equations  Solve chip compressibility equation  ( End ) Figure 4.4 Flow chart of the developed CFD code  65  ( End )  Chapter 4 Computational Methods  4.6 Code validation Validation is a very important stage in the process of code development. The main purpose of validation is to ensure that the implementations in the code are correct and free of errors. The following section describes some of the validation work that has been carried out. 4.6.1 Flow in a cavity As classical study cases, flows in a cavity, either lid-driven or buoyancy driven, two or three-dimensional, are well documented [59,61]. The geometry and boundary conditions in a two-dimensional situation are shown schematically in Figure 4.5 (a) and (b). 1 As//////////////////, /  Cold  Adiabatic 'A Adiabatic  /777777777777777777777 Hot  Moving lid  (a) Lid-driven  (b) Buoyancy-driven  Figure 4.5 Geometry and boundary conditions for 2D cavity flows In Figure 4.5(a), the square cavity has unit height and width, the bottom wall moves with unit speed in the positive x direction, the other walls are all at rest and the flow inside the cavity is assumed to be two-dimensional, steady, incompressible, and laminar with constant transport properties. Figure 4.6 shows the streamline and U-velocity along the vertical centerline when Re number equals 1000, the mesh size used in the 66  Chapter 4 Computational Methods  computation is 80x80. In Figure 4.6b, data reproduced from Ferziger and Peric's work [55] using center difference scheme are also plotted, as can be seen, the difference between the two velocity profiles is very small. In Figure 4.6c, the computed U-velocity at the vertical centerline using different mesh sizes are plotted, which indicates the solution is grid-independent. The convergence history of a typical computation given in Figure 4.6d shows that the maximal residual of P, V and V-equations drops below 10"  6  after 700 iterations and then remains almost constant. This is because single precision data are used exclusively in the code.  -0.4 -0.2  0  0.2  0.4  0.6  0.8  1  U  (b) U-velocity along vertical centerline  (a) Streamlines  600  900  1200  1500  iteration times  (d) Convergence history with mesh size of 40x40  (c) Computed U-velocity at vertical centerline using different mesh size  Figure 4.6 Flow field for 2D lid driven problem (Re=1000) 67  Chapter 4 Computational Methods  In Figure 4.7 and 4.8, velocity profiles of lid-driven flow in a three-dimensional cavity are plotted, where Re numbers are 100 and 400 respectively. The results are very close to those given by Ching et al. [59].  u  x  (a) U-velocity at plane: x=y=0.5  (b) V-velocity at plane: y=z=0.5  Figure 4.7 Velocity profiles in 3D lid driven problem (Re=100, mesh size: 25x25)  u  x  (a) U-velocity at plane: x=y=0.5  (b) V-velocity at plane: y=z=0.5  Figure 4.8 Velocity profiles in 3D lid driven problem (Re=400, mesh size: 25x25)  In buoyancy-driven cases, as shown in Figure 4.5b, the bottom wall, at a temperature T=T.O, moves with unit speed in the positive x direction. The other walls are at rest with 68  Chapter 4 Computational Methods  the temperature T=0.0. Chen et al. [60] has simulated this flow using the finite element method, with contour plots of stream function and temperature solutions for Re=100 shown in Figure 4.9. For comparison, the same problem is solved by the present code and corresponding results are provided in Figure 4.10.  (c) Temperature for Re=100, Pr=1.0  (d) Temperature for Re=100, Pr=10  Figure 4.9 Stream function and temperature contours from Chen et al. [60] (Re=100)  69  Chapter 4 Computational Methods  x  x  (c) Temperature for Re=100, Pr=l .0  (d) Temperature for Re=100, Pr=10  Figure 4.10 Stream function and temperature contours by this work (Re=100)  To make sure the code also works correctly for curvilinear grids, the twodimensional flow in an inclined cavity is also simulated. This time, the upper wall moves to the right with a unit speed while all the other walls are at rest. The streamline when the cavity is inclined to 45° and the Re=1000 is given in Figure 4.11; the computational results of the same problem given by Ferziger and Peric's [55] are shown in Figure 4.12. Again, the differences between these two simulations are very small.  70  Chapter 4 Computational Methods  x  Figure 4.11 Streamlines in 2D lid-driven cavity with both side walls inclined at 45° (Re=1000)  Figure 4.12 Streamlines in 2D lid-driven cavity by Ferziger and Peric [55] (with both side walls inclined at 45°, Re=1000) 4.6.2 Liquidflowthrough porous media along aflatplate We have presented the validation of a number of single-phase flow problems with or without heat transfer, occurring in cavities. However, the cooking process is essentially a two-phase flow problem. Although both phases are described by Navier-Stokes equations, there are inter-phase interaction terms contained in the governing equations for both the chip and free liquor phase. It would be very helpful to further validate the code if these terms can be tested. Lauriat and Vafai [63] have obtained the analytical solution to the  71  Chapter 4 Computational Methods  fully developed boundary layer flow when liquid passes through a porous media along a plate. The velocity profile is found to be (4.68)  u = l--^-sec/zV [(v* +c,)/2] /2  c, is the integration constant given by c, = 2cT where u =ulu ;  1/2  y' =  x  (4.69)  sec/T'^a//?)  (4.70)  yI4KIS  a = (l + 2 F - R e ) ; B = 2F-Re /3  (4.71)  Re,  (4.72)  p  p  As an example, when p=1.0, p=2xl0" , « =1.0, the porosity (or void fraction) 3  M  8=0.53, the permeability K=4xl0" , and the inertial coefficient F=0.55, from the above 4  expressions, Re =l,a=2.1, (3=0.3667, c, =7.617, thus it is easy to compute the p  dimensionless velocity u versus the dimensionless height y*. The same problem with identical parameters is also computed with our code, with velocity profiles of a typical case shown in Figure 4.13. As can be seen, the fully developed velocity profiles are obtained downward at a certain distance, where the entrance effect is negligible.  N  0.03  0.07  0.05  0.13  0.09  Figure 4.13 Boundary layer flow passing through porous medium along a flat plate (p=1.0, p=2xl0" , s=0.53, k=4xl0 , Re =l) 2  p  72  Chapter 4 Computational Methods  The code provides very good predictions, as the results compare well with the analytical solutions in a wide range of Re numbers (based on pore size), as exhibited in Figure 4.14. As pore size based Re number increases, the flow becomes closer to the plug flow, and the code prediction is more accurate.  A O •  Re =1 (Analytical) Re =10 (Analytical) Re =100 (Analytical) Re =1 (Computed) Re =10 (Computed) Re =1 (Computed) p  p  p  p  p  'fj  00  p  5-  A/<fr /  2  A/ /  <i>  A/ A /• r* A^.' Oil A ,y A  o /" o / II  0 0.0  0.2  0.4  0.6  0.8  1.0  1.2  U*  Figure 4.14 Velocity profiles by code prediction and analytical solutions  4.6.3 Power-law non-Newtonian fluidflowingthrough a tube Skelland [62] studied the velocity profiles for power law fluid flowing through round tubes and found that under laminar flow conditions, the fully developed dimensionless velocity can be written in the form: u  3» + l  V  n+ l  (4.73)  1  where u is the velocity, V is the average velocity, n is the power law index, R is the radius of the tube. To ensure the non-Newtonian models are properly implemented in the code, the power law fluid flowing through a round tube is simulated and compared with equation 73  Chapter 4 Computational Methods  (4.73). Figure 4.15 shows the comparison of the velocity profiles predicted with the code and the analytical solution (4.73). 1.0 —  0.8  Analytical solution Numerical solution  Analytical solution Numerical solution  0.6 0.4 0.2  0.2 0.0,0.0  0.5  1.0  1.5  2.0  2.5  0.0 0.0  3.0  I • • - •l  0.5  1.0  1.5  2.0  2.5  3.0  u/V  u/V  (b) Power law index n=l/3  (a) Power law index n=l/9  1.0 Analytical solution Numerical solution  Analytical solution Numerical solution  0.0  0.5  1.0  1.5  2.0  2.5  3.0  0.0  0.5  1.0  1.5  2.0  2.5  3.0  u/V  u/V  (d) Power law index n=3.0  (c) Power law index n=1.0  Figure 4.15 Velocity distributions for power law fluids in laminar flow through round tubes (zero-slip between the fluid and the wall, power law index n varies from 1/9 to 3.0) As can be seen in Figure 4.15, over a large range of power law index n (from 1/9 to 3.0), the code can predict the velocity profile reasonably well; the largest error occurs at the tube center when n equals 3.0; also, fluids with low index n approach plug flow over a large portion of the tube diameter.  74  Chapter 4 Computational Methods  For comparison purpose, the code is also applied to an industrial digester (in Kamloops, BC) which was simulated and well documented by Pougatch et al [9]; the simulation results given by the current code are then compared with those in the literature. The comparison shows that there are no essential differences in the velocity distributions of both chip and free liquor phases, and the differences in the chip phase pressure and chip compaction are minimal; however, there are relatively large differences in distributions of temperature and alkali concentration. Discrepancies also exist in the distributions of residual carbohydrates and lignin contents. That is expected, as in the two computational codes, different chemistry models are used, and some of the empirical correlations used are also significantly different. Although extensive validation work has been conducted, on-site experimental data will also be used to assist the validation for the overall cooking process in later chapters.  75  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester In this chapter, the liquid flow in a laboratory digester filled with compacted wood chips is computationally modelled under various conditions. The experimental observations through electrical resistance tomography (ERT) images generally agree with the calculated results, but detailed comparison cannot be performed due to difficulties both in CFD modelling and experimental work.  5.1 Introduction Since the early 1980s continuous cooking has been growing rapidly. Many new cooking techniques have also been emerging. Some of them, including Extended Modified Continuous Cooking (EMCC) and Iso-Thermal Cooking (ITC) [21], have quickly found widespread use in the industry. The cooking performance is essentially a combined result of many variables. Although industries have made great progress in continuous cooking, many underlying physical and chemical mechanisms are still not well understood; continued efforts are still required to further improve the cooking performance, optimize the equipment design and get better control of the operation. The liquid flow in a digester is of significant interest as it affects the motion of chip column, and helps to distribute the chemicals and maintain a high enough cooking temperature. It is expected that a better understanding of the liquor flow will lead to an improved digester process.  76  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  In the past several years, a model digester was built by Prof. Bennington's group (Department of Chemical and Biological Engineering, UBC) and the Electrical Resistance Tomography (ERT) technology was implemented and used to visualize the liquid flow in this digester both with and without wood chips [64,65]. However, the ERT technology can only provide flow information in a limited number of digester crosssections and therefore, it is sometimes very difficult to get a good understanding of the flow structure. Significant progress has been made in the past decade in Computational Fluid Dynamics that allowed it to become a useful tool in achieving a better understanding of industrial processes. In this chapter, the liquid flow in the model digester will be examined numerically with chips in place.  5.2 The model digester and electrical resistance tomography 5.2.1 The laboratory digester Figure 5.1 schematically shows the basic configuration of the laboratory digester [64]. The cylindrical digester has a volume of 40x10" m with an inner diameter of 0.28 3  3  m and a height of 0.75 m. Two water circulation loops can be established in the digester. In the first loop (shown in blue arrows), water pumped from a reservoir moves up through the digester bottom, and flows out via three nozzles located on the top of the digester periphery. In the second loop (shown in red arrows), water is extracted through the circumferential screen and sent back into the digester through the downcomer. The screen is partitioned into sections by eight equally located tubes through which the extracted water flows.  77  Radial Flow  Downcomer 0.75 m  Screen Section  Distribution Plate  Axial Flow  Figure 5.1 Configuration of the model digester To model the flow patterns in real digesters, two fluid circulation loops were created. The first one provides axial upward flow Q , which is injected from the bottom; while A  the second one gives radial flow Q which is injected through the central downcomer. R  The liquid flows out from the extraction screen at the top of the digester. The downcomer can move up and down so as to model different flow situations. Cooked or uncooked wood chips (usually collected from mills) can be added and compacted against the perforated distribution plate near the bottom of the digester. The flow rates Q  A  and QR  can be set at different values. An industrial Electrical Resistance Tomography (ERT) system was installed to generate an image of the flow field. A number of different experiments have already been conducted with this device [64, 65]. 5.2.2 Electrical Resistance Tomography There are numerous methods/ways to make measurements in fluid flows. The usual experimental probes, such as pitot tube, the hot wire and the more advanced laser-doppler  78  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  equipment, can be characterized as point measurements, i.e., they only give information in a single spatial point in the flow field. Although useful, the information is quite limited, especially when flow phenomena are dominated by spatial structures. Newly developed techniques such as particle tracking velocimetry (PTV), particle image velocimetry (PIV) and laser induced fluorescence (LIF) are capable of obtaining spatial information, but they are quite expensive, and usually are used in gaseous and liquid flows. Electrical Resistance Tomography (ERT) belongs to a group of imaging methods collectively termed electrical tomography techniques, which have been extensively reviewed [66,67,68,69,70]. ERT visualizes the flow by imaging the spatial difference in conductivity through the sample volume. In other words, ERT images are generated from measurements of electrical resistance. The ERT technique has greatly progressed since it was developed in the 1980s and offers many benefits, including relative low-cost, easy data collection, high sensitivity, and fast image rates; it can even be used online to monitor and control the operation of existing process equipment with careful design. ERT has been used extensively in various aqueous flows, and recently, it has been employed in studying gas/liquid two-phase flow [70]. A typical ERT system is composed of three main parts: sensors, data acquisition system and image reconstruction system/host computer. As shown in Figure 5.2, in the model digester, there are eight sensor planes equally spaced along the digester height, each sensor plane consists of 16 electrodes evenly configured in the vessel periphery.  79  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  Figure 5.2 The laboratory digester and the ERT sensor planes  The ERT system uses a P2000 data acquisition system to inject current between any pair of electrodes and then measures the potential difference. For the configuration of 16 electrodes,  104 independent measurements  are made. Afterwards, the image  reconstruction by linear backprojection algorithm is conducted on that plane. The ERT generates an image of conductivity distribution in the cross section by measuring the electrical resistance.  In the experiments conducted in the model digester, two techniques were used to image fluid flow. One technique uses controlled pulse of the tracer, the other one uses a step change in water conductivity. The higher conductivity solution can be stored in one of the two storage tanks and introduced into the system as desired.  5.3 Computational method In the laboratory digester, wood chips are stationary; cooking reactions do not occur as cold water instead of hot cooking liquor is used. Therefore, inter-phase mass transfer is not involved, as this is basically a liquid flow through a compacted porous media - chip phase.  80  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  The porosity Sf is measured when wood chips are loaded. Then the water flow is introduced. As an approximation, the porosity is considered to be a constant, although the water flow does have some impact on the chip column compaction. To model the liquid flow with chips in place, equation 2.7 and 2.11 are used, with the term accounting for mass transfer eliminated. Because £ f is treated as a given constant, the compressibility equation is not solved. The inter-phase friction force determined by Wang [29] (for reference pine chips) is used in the modelling. Figure 5.3 gives the flow resistance over pine chips.  Q.  5  Q.  < a  e £ 3  (A in £>  CL  2  4  6  8  10  12  14  16  18  20  Velocity difference AV, mm/s  Figure 5.3 Flow resistance over chip column, determined by Wang [29], for reference pine chips  In this study, purely liquid flow, i.e., without the presence of wood chips is not simulated, as in this case, the flow is in the laminar/turbulent transition regime.  81  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  5.4 Computational and experimental results 5.4.1 Approximation of boundary conditions As already described, there are two liquid loops in the model digester. One is injected through the perforated holes located in the downcomer periphery; the other one is introduced from the bottom of the digester. The injection port on the downcomer consists of many densely distributed small open holes (3mm in diameter); for simplicity, the entire perforated region is considered as a single inlet. To check this assumption, two different cases were computed. In the first case the plane water jets exiting from three slots were modelled with the dimensions and flow rate being typical of those occurring in the downcomer. In the second case an equivalent rectangular geometry was used. The port geometry and velocity distributions are shown in Figure 5.4 and 5.5.  3 mm 1mm  l  0 ls 0  i I i  I  i  j  ii 0.4  i  j  i  i  i  j  ii 0.8  i  I  i  i  i  i  i I 1.2  i  I  i  X, cm  (a) Port geometry  (b) Calculated velocity distribution  Figure 5.4 Plane jet: water injected from triple holes, m =0.33kg/s 0  82  i  i  i  i 1.6  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  11m  0.0  11  1 1  1 1  1 1  1  ' ' 0.4 1  1 s  1  '  1 5  1  ' ' 0.8 1  1  '  1  5  1 5  1  = 1.2  1 5 1  1.6  X, cm  (a) Port geometry  (b) Calculated velocity distribution  Figure 5.5 Plane jet: water injected from a single port, m =0.33kg/s 0  As can be seen in Figure 5.4, the velocity between the slots can quickly increase to the value in the main stream; when the distance is larger than twice the slot height, i.e., x>6mm, the velocity profiles in Figure 5.4 and 5.5 are almost identical. So treating the whole perforated area as a single open port will not cause a large error, especially further downstream from the injection port. Because the downcomer has a dead end at the bottom, the inside water pressure will attain its highest value there. Thus the injection pressure is high in the lower portion and low in the upper portion of the injection port. This leads to higher injection velocities in the lower part of the perforated area. This effect has been observed in the experiments. For the bottom inlet, as a distribution plate is installed, uniform inlet velocity can be assumed. Further, we assume that the outflow through the screen section is uniform. The exit flow from the top of the digester is imposed as an outflow boundary condition (zero gradient condition).  83  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  As an approximation, the flow is assumed to be symmetrical, thus only twodimensional computation is needed to model the liquid flow; the corresponding geometry and mesh are given in Figure 5.6.  Downcomer  Figure 5.6 Computational domain, mesh, inlet and outlets  5.4.2 Liquid flow in the digester Vlaev and Bennington [65] investigated experimentally the liquid flow in the model digester in the presence of wood chips. In these experiments, chips can be cooked or uncooked, compressed or uncompressed, and the flow can be symmetrical or asymmetrical. In the following, both the experimental and computational results performed on the symmetrical cases with uncooked chips are presented. Figure 5.7, 5.8 and 5.9 show the ERT images and modelled velocity fields for three different flow conditions with the ratios of axial to radial flows of 0.5, 1.0 and 2.0 respectively. In all the cases, the volume fraction of chip phase is 0.53; tracer is injected from the bottom of the digester. In all the ERT images, the red color represents the highest concentration of the tracer, while the blue indicates the lowest. From the color images, one can 84  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  qualitatively infer the flow path and the velocity magnitude. However, it is generally hard to obtain the exact velocity field because of the low resolution and two-dimensional nature of the ERT images and therefore detailed comparison with the computations is difficult; moreover, the ERT images may not represent accurately the real flow if the images are obtained too long after the injection, especially when the tracer diffuses quickly.  r, cm  (a) ERT image, after 95 sec  r, cm  (b) Velocity  (c) Streamline  Figure 5.7 Flow field given by ERT image and CFD modelling Q =125xlO" m /s, Q =250 xlO" m /s 6  3  6  A  3  R  From Figure 5.7 (c), all the liquid injected at the digester bottom exits from the screen section. As part of the bottom inflow, the majority of tracer should flow out through the screen section as well, but there might be a small amount of tracer diffused. However, this upflowing tracer will be likely restrained below the liquid injection level from the downcomer, as part of this injection flows down- and outward and exits from the screen section. This is reflected in the ERT image given in Figure 5.7 (a), where the tracer concentration remains high in the lower portion of the digester, while very low in the upper portion above sensor plane 4.  85  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  2  5  8 11 14  r, cm  (a) ERT image, after 229sec  (c) Streamline  (b) Velocity  Figure 5.8 Flow field given by ERT image and CFD modelling Q =125xlfJ m /s, Q =125 xlO" m /s 6  3  6  A  3  R  Figure 5.8 (b) and (c) show the velocity and streamline when a flow rate of 125xl0"  6  m /s is used in both the axial and radial loops. The flow field is similar with the previous case shown in Figure 5.7, but most of the liquid injected from the downcomer flows upward and leaves the digester from the top exit. In Figure 5.8 (a), relatively high tracer concentrations are observed up to the sensor plane 4. Again, the tracer concentrations in the upper sensor planes are low, indicating most liquid injected from the bottom exits from the extraction screen.  50P"''  40|simm  I £ '3  ;  1 0 k ••  2  5  8 11 14  2  r, cm  5  8 11 14  r, cm  (a) ERT image, after 95 sec (b) Velocity (c) Streamline Figure 5.9 Flow field given by ERT image and CFD modeling, with flow rates Q =250xl0" m /s, Q =125 xlO" m /s 6  3  6  A  R  86  3  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  Figure 5.9 shows the flow field for a higher axial flow rate. From the streamline in Figure 5.9 (c), part of the tracer flows up above the extraction screen, moves up continuously near the digester wall, and leaves the digester from the top exit. After a short time (95 seconds), the tracer is detected as high as the top sensor plane as shown in Figure 5.9 (a); in addition, higher tracer concentrations are found near the digester wall. In all the three cases, both experimental ERT images and computational flow fields show that beyond the circulation region, the upward liquid velocity near the digester periphery is higher than that in the center area, indicating that the flow is not uniform, but after some distance away from the screen section, the liquid flow gradually become uniform. Also, as the flow ratio of Q A : Q R increases, the flow pattern in the digester changes as well. In Figure 5.7 and 5.8, all the fluid injected from the bottom exits from the screen section. In Figure 5.9, roughly half of the bottom liquid exits from the extraction screen; the rest of the liquid flows upward above the screen zone near the digester wall and finally exits from the top of the digester. The upward axial velocity in the lower part of the digester increases from case (1) to (3). In case (1), the upward axial velocity is so small that the tracer is experimentally detected mainly in the bottom four plates, as shown in Figure 5.7. In case (2), the tracer can reach the top plate 7, but takes a very long time; the reason is that the bottom flow is constrained to the lower part of the vessel. When compared with (1), the magnitude of the axial velocity in the lower part of the digester increases, and the tracer can be easily transported to plate 3 and 4, as afterwards, it spreads further by diffusion and convection; therefore a moderate tracer concentration is observed, as indicated by the green color shown in Figure 5.8 (a). In case (3), computational results shows that a portion of the bottom flow is pushed to the digester wall and continues to the digester top when it spreads beyond the screen section. On the other hand, the tomographic images in Figure 5.9 (a) indicate that the tracer can rapidly reach the top plate 7. Images in plate 5, 6, and 7 in Figure 5.9 (a) exhibit higher tracer concentrations in the outer layer and this agrees with the computational results in Figure 5.9. 87  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  Similar circulation loops exist in a number of different places in industrial digesters. Knowing the above flow patterns may bring a better understanding that can lead to a better design or operation. Taking wash circulation in an EMCC digester as an example, the flow pattern shown in Figure 5.8 is probably the best in terms of cooking and heating. In the case of Figure 5.7, part of the heated wash liquor and added white liquor are extracted and sent back to the wash liquor circulation loop again; in the flow conditions of Figure 5.9, however, only a part of the wash liquor is extracted and heated; the temperature and white liquor distribution at the bottom of the digester are not uniform, which may undermine the cooking uniformity. The distributions of the liquid dynamic pressure for the three cases are given in Figure 5.10.  r, cm  r, cm  (a)  r, cm  (b)  (c)  Figure 5.10 Calculated liquid dynamic pressure (with chips in place) (a) Q =125xlO" m /s, Q =250xlO" m /s; (b) Q = Q =125xlO" m /s; (c) Q =250xlO" m /s, Q =125xlO" m /s; 6  3  A  6  3  R  6  A  6  3  A  3  R  6  3  R  A reference pressure of 22 Pa is set to the digester top. From Figure 5.10, the liquid pressure is high at the injection places and low at the exits. In addition, the circulation ratios have a significant effect on the pressure distributions: in Figure 5.7, the radial flow is large and thus the pressure at the injection port of the downcomer is high; in Figure  88  Chapter 5 Modelling of Liquid Flow in a Laboratory Digester  5.9, the axial flow is twice as large as the radial one, so the upward velocity is large and therefore the liquid pressure at the digester bottom is high.  5.5 Conclusions The liquid flow in a laboratory digester where chips are filled in is modelled computationally. Calculated results agree reasonably well with the experimental observations obtained by ERT images. The computations have led to the following conclusions: In the presence of compacted wood chips, the momentum of the liquid diffuses very efficiently and the liquid flow tends to be uniform away from the injection and extraction places. Depending on the liquor circulation ratios, different flow patterns can occur in a digester. Appropriate flow rates will be beneficial to uniform cooking. As expected, a higher circulation rate requires a higher injection pressure. Because of the low resolution of the ERT images, detailed comparison between computation and experiments is difficult.  89  Chapter 6 Simulation of an Industrial Continuous Digester As pointed out in previous chapters, continuous cooking is a very complex process. It is also a process hard to control mainly for the following reasons [71]. - In practice, there are usually big variations in chip quality, such as moisture content, size distribution etc; sometimes, the wood species also changes; - It is almost impossible to make online measurements inside a digester because of the presence of moving solid chips and the hostile environment; - There is a long process delay; it can take even 5 hours for an individual chip to move from the chip bin to the digester blow line. To maximize the economical benefits, better knowledge and control of the pulping process are required. CFD technology can help meet this demand in two ways. CFD technology can help to identify some key parameters and assess their influence on the pulping process. This is very useful, as for instance, predicting with reasonable accuracy the impact of upstream variations will allow operators to make adequate adjustments in a timely manner. CFD can also help to achieve a better pulping process by predicting the pulping performance before any desired changes are actually made. This helps to give a good retrofit or design that could significantly decrease the capital cost, and improves the pulping quality. This chapter introduces the computational work performed on an EMCC continuous digester operated by Howe Sound Pulp and Paper Ltd. (HSPP). The cold blow is assumed to be a uniform up-flow through the digester bottom; all boundary conditions are derived from the mill operating data. A typical simulation produces the velocity profiles, pressure distributions and volume fractions for both phases; the temperature, effective alkali concentration and kappa number are also obtained. The impacts of variations in production rate, chip type and cold blow rate on cooking performance are also  90  Chapter 6 Simulation of an Industrial Continuous Digester  computationally evaluated. In this Chapter, a three-dimensional simulation for the continuous digester is also presented to evaluate the non-uniformity of cooking when the liquor extractions are made in a special manner.  6.1 Introduction 6.1.1 Two-vessel hydraulic system in HSPP Howe Sound Pulp and Paper Ltd. (HSPP) produces both mechanical and kraft pulp. The kraft pulp in HSPP is sold directly in the world market, instead of being transformed into different paper products. Currently, HSPP produces 1100 tonnes of kraft pulp each day. The cooking equipment used in HSPP is a two-vessel hydraulic continuous type, which has a separate impregnation vessel. A technique called Extended Modified Continuous Cooking (EMCC) is being used. Chips, Liquor Chips  1.  Chip Bin  Steam Steam  J  High Pressure Feeder  Cold Blow  Brown Stock  Figure 6.1 Two-vessel hydraulic digester system in HSPP  91  Chapter 6 Simulation of an Industrial Continuous Digester  Figure 6.1 is a schematic diagram illustrating this cooking system. In this system, the white liquor charge is divided into four parts. The first charge is delivered to the impregnation vessel. The second one is delivered to the top of the digester where the chips are transferred from the impregnation vessel into the digester. A third charge is introduced to the loop of modified cooking liquor circulation, while the last one is introduced to a place above the washing screen through the central downcomer. In the digester, chips move vertically downward continuously, the spent cooking liquor and washing liquor are extracted through a circumferential screen located in the middle of the digester, thus a co-current flow zone is formed in the upper portion of the digester while a counter-current flow occurs in the lower part of the digester. There are also several heating circulations to maintain the required reaction temperatures. The digester can be divided into four zones, namely cooking zone, MCC (Modified Continuous Cooking) zone, EMCC (Extended Modified Continuous Cooking) zone and washing zone. As compared with conventional cooking, EMCC has several advantages: (1) Because part of the cooking is conducted in the countercurrent zone where the liquor has a much lower solid content, the lignin can be further removed without much cellulose degradation. Therefore EMCC has the capacity to cook chips to very low kappa numbers while retaining pulp strength and yield. It was reported that with EMCC, the kappa number could be lowered to 22-24 for softwood and to 14-16 for birch (hard wood). (2) The alkali is charged into several different places to make the alkali content more uniform during the cook. Consequently, the time for delignification becomes longer, but the non-uniformity of cooking decreases. (3) The digester volume is more effectively utilized as cooking reactions take place in the lower part of the digester as well. (4) This process reduces the need for bleaching chemicals and consequently can make a major contribution to the reduction of emissions of chlorinated organic compounds and some other pollutants in the pulp and paper industry.  92  Chapter 6 Simulation o f an Industrial Continuous Digester  6.1.2 Mill records of digester operations To facilitate our computational study, HSPP has generously provided us with comprehensive data sets on a number of operations of the EMCC digester. It was a serious effort to keep all operations stable for a minimum of 6 hours; the time-averaged values are obtained from the hourly measured data. Appendix A shows the mill records for two typical cases; Appendix B lists the averaged data for a total of 15 cases. In all computations, only the averaged data are used to determine the boundary conditions.  6.2 Blended chip species Wood pulp consisting of longer fibers is stronger than that of short fibers; so many pulp mills blend different wood species, usually two at a time, to ensure that the produced pulp meets the target strength required by customers. Table 6.1 gives some laboratory measurements of the average fiber length of different wood species in British Columbia [72]. Table 6.1 Average fiber length of different wood species [72] Average length, mm  Species BC interior varieties  BC coastal varieties  Lodgepole pine  2.3  Spruce  1.9  Douglas fir  2.0  Hemlock  3.0  True fir  3.3  Douglas fir  3.4  Western red cedar  2.4  Yellow cedar  2.5  For the simulations, the oven-dry wood density is required. For blended wood chips, the oven-dry wood density can be evaluated by the weighted sum of all the species involved, or mathematically,  93  Chapter 6 Simulation o f an Industrial Continuous Digester  PoDw=YX rPoDWj)  C-)  a  6 1  Here, a\ is the mass fraction of the i-th species in the mixture, p , ODW  is the oven dry  wood density of the i-th species. In Howe Sound Pulp and Paper Ltd., four main species are used, namely hemlock, cedar, spruce and Douglas fir. The average specific gravities of these four species are given in Table 6.2. Table 6.2 Specific gravity of some common woods [72] Specific Gravity  Wood species Western hemlock  0.44  Western red cedar [74]  0.31  Engelmann spruce  0.36  Douglas-fir, coast-type  0.51  6.3 Computational domain The schematic of the digester used in HSPP is given in Figure 6.2. The digester has a total height of 51.2 meters and a radius of 4.228 meters of the bottom shell. The digester produces around 1300 tons of air-dry pulp per day, with a target kappa number of 28. The chips form a free surface (which is called digester chip level) at a short distance from the digester top. In the digester operation, a fixed chip level is desired so as to maintain a stable operation. In practice, the chip level changes in a certain range. Above this level, wood chips and cooking liquor fall down by gravity. To avoid the difficulty of modelling the complex flow in the top part of the digester, we start the computation below the chip level and assume that chips and cooking liquor enter the digester with uniform velocities.  94  Chapter 6 Simulation of an Industrial Continuous Digester  While developing the energy equation in Chapter 2, we assumed that the chips and liquor are i n thermal equilibrium, namely the temperatures are identical for both phases. This is a commonly used assumption for liquid flowing through porous media. Then we solve only one energy equation. However, for the digester bottom portion where blowing occurs, the temperatures o f the cold blow liquor and chips are significantly different. In addition, the cold blowing liquor is usually injected upward through a number o f nozzles and that is difficult to simulate in detail. A s an approximation, the C F D computational domain ends 1 meter below the wash liquor screen, and the cold blow is assumed to be a uniform upflow at the bottom plane; this dramatically simplifies the blowing process. In Figure 6.2, the computational domain is shown in the enclosed box with dotted line. Chips & Liquor  •7T>— BCL Trim liquor (not used)  MCL  WL Brown Stock  Cold Blow  Figure 6.2 Illustration o f the computational domain o f the digester (BCL: bulk cooking liquor; B L : black liquor; M C L : modified cooking liquor; W L : wash liquor; Trim liquor: liquor is drawn from chips and circulated to bottom circuit.)  95  Chapter 6 Simulation of an Industrial Continuous Digester  6.4 Boundary conditions 6.4.1 Boundary conditions for velocity  In Chapter 2, the chip phase is defined as wood chips and the entrapped liquor, and liquor phase is defined as free liquor. As the Navier-Stokes equations are solved for both the liquor and chip phases, boundary conditions for the velocity should be provided separately for each phase. For the liquor phase, the bottom plane is considered as an outflow boundary, while all the rest of the ports are treated as inflow boundaries. For the chip phase, the bottom plane is also considered as an outflow boundary and the top plane is considered as an inflow boundary; in addition, the slip condition is used for both the downcomer and the digester walls. For all inflow boundaries, the flow rate of the corresponding phase is needed; however, this is not so obvious, for instance, at the digester top, both solid and liquid flow through the same port, but neither of the flow rates is measured there. To compute the velocity at the boundaries for each phase, the flow rate of each phase and the flow area occupied by each phase have to be determined. In this section, the way of estimating the velocities at the digester top and bottom is presented. 6.4.1.1 The top inlet (1) Compute flow areas  Both solid and liquid are introduced at the digester top inlet. Based on the assumptions made earlier, the chip phase is a porous matrix; and the free liquor flows through the interstitial space of the porous matrix. When the porosity s is determined, the flow area of each phase can be calculated as, (6.2) A =e A c  c  =  (6.3)  (\-e,)A  In Equation (6.1) and (6.2), A is the total area through which the chip phase and free liquor phase flow. The chip phase pressure at the chip level is assumed to be zero, i.e., P =0. Then from the chip column compaction equation given by Harkonen, c  96  Chapter 6 Simulation o f an Industrial Continuous Digester  e = 0.644 + f  (-0.831+ 0.139 In*-)  10 y 4  V  (2.27a)  the volume fraction of liquid is st=0.644, thus the volume fraction of solid is s =Tc  Sf=0.356.  Further, we can calculate the volume fraction of wood material. Suppose the wood material density is /^lSOOkg/m  3  If the blended chips consist of 43.4% cedar and 56.6% spruce, according to Equation (6.1), the oven dry wood density can be calculated as, p  =0.434x310+0.566x360= 338.3 kg/m  0DW  3  (In case of blended chips, the oven-dry wood density can be calculated by equation 6.1). We have, (6-4)  = Poo, c  £  P0DW c £  =  =  Q  0  8  0  3  (  6  5  )  Pp.  =0.276  (6.6)  The volume fractions of the wood material e  pw  and the entrapped liquor sn will be e  used to determine the chip phase density in the subsequent section.  (2) Calculate flow rate for each phase The digester has a complex inflow/outflow configuration and requires careful consideration when determining the velocity at its boundaries. In the two-vessel system (Refer to Figure 6.1), the total inflow via the digester top includes chips (with moisture), white liquor, consumed steam in the chip bin and steaming vessel, cold blow to the 97  Chapter 6 Simulation o f an Industrial Continuous Digester  impregnation vessel, and the bulk cooking circulation. (Some liquor is extracted near the digester top, its temperature is increased and then it is directed to the bottom of the impregnation vessel so as to carry wood chips to the main digester vessel; this is called the bulk cooking circulation. This is shown in Figure 6.1). For instance, selected records on 09-MAR-2004 from HSPP digester operations (Refer to Appendix B ) are listed on Table 6.3. Table 6.3 Selection of recorded data from HSPP digester (09-MAR-2004) Column in  Recorded  Description  Unit  data  Appendix B AE  Weighted average moist  54.26  %  AF  B l o w Consistency C |  10.00  %  AY  Total chips  2428  BDt/d  Digester B l o w F l o w  131.9  L/s  BE  Upper Slice F L  121.6  L/s  BD  B u l k C o o k i n g Circulation  150.4  L/s  0  HPF W L PURGE  3.94  L/s  P  W L F l o w to W A S H  4.69  L/s  Q  W L F l o w to B C  1.74  L/s  R  W L F l o w to M C  4.68  L/s  S  W L F l o w to I V  45.99  L/s  T  C B F l o w to IV  0.41  U  C o l d B l o w to Digester Bottom I  47.78  L/s  V  C o l d B l o w to Digester Bottom 11  10.91  L/s  X  C o l d B l o w to Digester Bottom III  39.54  L/s  BF  M C Extraction  120.1  L/s  W  B l a c k Liquor Extraction  139.5  L/s  Z  B l o w Rate  1267.6  ADt/d  L  b  o w  98  Chapter 6 Simulation of an Industrial Continuous Digester  The solid/liquor phase flow rate through the top inlet can be obtained by the following procedure: (1) Compute the total inflow M j (from the digester top and bottom) n  At steady state flow conditions, the total inflow equals the total outflow, and the total outflow is the sum of the extracted black liquor and the digester blow rate. In Table 6.4, the computation of the total outflow is outlined in rows 1 to 3. (2) Determine the total inflow from the digester top  Mi _to n  P  Wood chips enter the digester from the top; and cooking liquor enters the digester either from the top or the bottom plane. In the bottom, the incoming flow is the cold blow, in the form of three splits. When the cold blow is subtracted from the total inflow, the inflow at the digester top is obtained. This is shown in rows 4 to 5 in Table 6.4. (3) Calculate the flow rate of the solid/liquor phase The inflow of solid chips is known; however, this is not the flow rate of the chip phase. To determine the flow rates of chip and liquor phases, we assume that the chips are fully impregnated. In other words, the interstitial space of wood chips is completely occupied by the cooking liquor. In this situation, the chip phase flow rate will be, 1  j m  c.p  =  p  m  W  +  en  =  m  1050x0.276 = 28.10 1 + 1500x0.0803  Pliq ^ en  pw  m  P £ r  pw  pw  And the chip phase density is, Pcr.  =  P  p  w  S  p  w  +  P  " "  £  c  n  =1151 kg/m  3  Then the flow rate of the incoming liquor phase is, m  fp  = m -m m  cp  = 183.35-95.71=87.64 kg/s  99  =95.71 kg/s  Chapter 6 Simulation of an Industrial Continuous Digester  Table 6.4 Calculation of solid/liquor phase flow rate through the top inlet Row  Description  Vi  Unit  Cited data in Table 6.3  1  Black liquor extracted  147.5  kg/s  px(W)  2  Digester blow rate  140.5  kg/s  px(L)  3  Total mass exiting from the digester M  288.0  kg/s  4  Total cold blow from the bottom M B  104.6  kg/s  px(U+V+X)  5  Total mass inflow via the top M  183.4  kg/s  V3-V4  6  Total solid inflow (chip flow rate)  28.10  kg/s  (AY)  7  Chip phase flow rate via the top M  95.7  kg/s  8  Net liquor phase flow rate via the top M  87.6  kg/s  9  Blow rate M i  13.2  kg/s  10  Total flow area A (for both phases)  37.7  m  OUT  C  b  I N TOP  C P  F P  o w  Note: Suppose pi, =1050, p =1000, PcoidBiow 1065 kg/m =  q  (Z)  2  3  water  (3) Calculate the velocity at boundary The average velocities of solid and liquor phase are given by, m W =— -*-  (6.7)  c  c  Pc c £  A  YYl  W = f  (6.8)  peA f  f  6.4.1.2 Bottom exit The liquor flow at the digester bottom is complex. As illustrated in Figure 6.3, the upflowing cold blow is split into three parts, which enter the digester via different channels. A cap and a distributing device are set immediately above the digester exit to  100  Chapter 6 Simulation o f an Industrial Continuous Digester  distribute the liquor flow. In practice, the liquor flow in this part could be very complex, for instance, a high injection flow-rate from the side tubes could result i n channeling near the digester wall. Consequently, most of the cold blow would move up rapidly through the channeling. In the current stage of our computation, however, this occurrence cannot be captured and assumptions are made to approximate the liquor flow at the digester bottom. A s shown i n Figure 6.3, the cold blow is assumed to be a uniform upflow through the bottom plane A i A ] of the computational domain.  Brown Stock  Cold Blow Figure 6.3 Illustration of the liquor flow at the digester bottom  It is important to note, however, that not all the cold blow moves upward because a portion of the liquor exits along with the chip phase and dilutes the resulting fibers i n the subsequent blowing process. To estimate the net upflow of the free liquor phase, we start with the mass balance of the liquor phase,  mliq,top  1 7 1  liq,bottom + \2 m  m  (6.9)  extract  Equation (6.9) states that the incoming free liquor from the digester top and bottom plane plus the mass transferred from the chip phase equals the extracted black liquor. In  101  Chapter 6 Simulation o f an Industrial Continuous Digester  the previous calculation, the free liquor via the digester top and the extracted black liquor are obtained, as given in Table 6.4 row 1 and 8. From equation (6.9), if the mass transfer rate m is known, the free liquor flow from the digester bottom plane can be calculated. n  In column Z of Appendix B (also see Table 6.3), the digester blow rate is 1267.55 air-dry tones per day, which is equivalent to 13.2kg/s bone-dry pulp; in Table 6.4, the total solid inflow (chip flow rate) is 28.10 kg/s, so the dissolved solid material is 14.9 kg/s, which corresponds to a dissolved volume AQ. If one assumes that this volume is fully replaced by the surrounding cooking liquor, then the net mass transfer from the chip phase into free liquor phase m can be approximated by, n  rh = (p - A,, )AQ = (1500-1065) x i ^ = 4.32 kg/s n  w  (6.10)  Using Equation (6.9), the upflow of the free liquor is then =  - m , ~ ™n = 147.5-87.64-4.32=55.54 kg/s lig lop  (6.11)  6.4.2 Boundary values for species equation (Lignin, Carbohydrates and EA)  In Chapter 3, the proposed pulping reaction model involves four species: lignin, carbohydrates, effective alkali, and sulfidity. However, the sulfidity is usually thought to be a constant through the entire digester and is only measured in the fresh cooking liquor stream. Therefore the measured sulfidity is given as a constant in the digester modelling, instead of solving the conservation equation of sulfidity. Obviously, boundary values are required while solving the rest of the three species conservation equations. 6.4.2.1 Lignin and Carbohydrates  The lignin content is directly related to the kappa number, which is an important factor in determining the compressibility and porosity of wood chips. The final contents  102  Chapter 6 Simulation of an Industrial Continuous Digester  of lignin and carbohydrates also determine the yield of the pulping process. Therefore, it is important to predict accurately the contents of lignin and carbohydrates. In the species equations given in Chapter 2, the mass fractions, instead of the actual contents of the lignin and carbohydrates, are used. For lignin and carbohydrates, the mass fraction is defined as the residual lignin/carbohydrates content in wood chips over the initial lignin/carbohydrates content. Based on the definition, the mass fractions of lignin and carbohydrates at the digester top inlet should both be 1.0 when a single vessel digester is used. In HSPP, however, a two-vessel hydraulic digester system is used. Thus the initial delignification reactions taking place in the impregnation vessel have to be considered. It was estimated that 5-15% lignin is consumed in the impregnation process. In our computations, it is assumed that 10% of lignin and carbohydrates are dissolved since the exact amounts cannot be determined. Therefore, the mass fractions of lignin and carbohydrates at the digester top inlet are both 0.9. Solving the species equations will give the mass fractions of residual lignin and carbohydrates; in order to obtain the actual contents of residual lignin and carbohydrates, initial contents of lignin and carbohydrates should be provided in a specific pulping process. The contents of lignin and carbohydrates in wood species are affected by many variables, including geographical location, soil and weather conditions, and location of the wood within a given tree. Table 6.5 and Table 6.6 give the cellulose and lignin contents of some soft wood species by different authors. From these two Tables, the cellulose contents among different soft wood species vary only about 3%, while the lignin contents differ significantly; for instance, the lignin content of western hemlock is 11.4%) higher than that of eastern hemlock. Both authors have provided the cellulose and lignin contents of Douglas-fir with some differences, which might be a consequence of different testing methods and experimental errors, or a difference in the chosen samples.  103  Chapter 6 Simulation o f an Industrial Continuous Digester  Table 6.5 Chemical composition of some soft wood species [75] (% of oven-dry wood) Cellulose  Lignin  Norway spruce  43  29  White spruce  45  27  Black spruce  44  30  Red spruce  44  28  Douglas-fir  44  32  Northern white cedar  43  31  Eastern hemlock  44  33  Species  Table 6.6 Chemical composition of some soft wood in North America [76] (% of oven-dry wood) Cellulose  Lignin  Western hemlock  42.9  21.6  Western red cedar  45.1  20.5  Engelmann spruce  43.6  26.5  Douglas-fir  42.2  29.6  Species  However, the contents of lignin and carbohydrates in wood chips are not measured in HSPP. Because HSPP has a very stable chip source, variations of chemical compositions in the same species are expected to be small. Here, we simply use the data of Table 6.6 to determine lignin contents since these wood species are from North America. The hemicellulose content is fixed at 27% of oven-dry wood as an averaged value for softwood [76]. As pointed out earlier, blended wood chips are used in HSPP; the initial contents of lignin and carbohydrates are then computed by,  8o  Li  (6.12)  =T ( i- So,,) a  Li  J  104  Chapter 6 Simulation of an Industrial Continuous Digester  Carb =  (a • Carb )  0  t  (6.13)  0t  Here, Lig ,. and Carb , are the initial contents of lignin and carbohydrates in the i-th Q  0  species, and the carbohydrates are considered as the sum of cellulose and hemi-cellulose. 6.4.2.2 Effective Alkali concentration In digester operations, effective alkali (EA) concentration is only measured in the white liquor, so EA concentrations in the digester inlet, outlet and extraction ports need to be computed. In the proposed pulping reaction model, mole concentration of effective alkali is used. It is therefore convenient to solve EA concentration equation in its mole concentration form, i.e., Equation (2.18). In the following, the way to calculate EA mole concentrations at the top inlet and at the MC extraction/injection ports is illustrated. 1. The top inlet The white liquor flows to the high-pressure feeder (O), impregnation vessel (S) and bulk cooking circulation and eventually enters the digester via its top plane, and the total amount is, MWL_TOP=3.94+45.99+1.74=51.67 L/s  (6.14)  The composition of white liquor used in HSPP is given in Table 6.7. Table 6.7 Composition of white liquor in HSPP (g/L) NaOH  69.6  Na S 0  Na S  33.6  EA as NaOH  86.84 103.2  2  2  2  3  Na C0  3  16.8  AA  Na S0 2  4  4.9  TTA  Na S0  3  1.0  Sulphidity as NaOH  2  2  2.0  120 32.60%  As listed in the Table, EA as NaOH is 86.4 g/L, so the mole concentration of EA in the white liquor is,  105  Chapter 6 Simulation o f an Industrial Continuous Digester  EA  86.84  NaOH  In Table  6.4  rows  2.171 mol/L  40  5  and  6,  (6.15)  the total mass coming in via the digester top is  183.35  kg/s, among which the solid inflow is 28.10 kg/s, so the liquid inflow via the digester top is 155.24  kg/s.  From Figure  6.1,  the bulk cooking circulation liquor QBC in the amount of  150.44 L/s  (see column B D on Appendix B ) is also introduced at the digester top, with a residual alkali of  12.79 g/L  (see column C in Appendix B ) , or equivalently,  0.3198  mol/L.  So the mole concentration of EA at the digester top inlet is, [EA]  Top  =  Total EA (mol) Total liquid (liter)  51.67x2.171 + 150.44x0.3198  = 0.499  155.24 /1.05 + 150.44  mol/L  (6.16)  2. M C extraction port/injection ports Figure 6.4 schematically shows the modified cooking circulation loop. A portion of the cooking liquor is extracted from the M C screen. After being mixed with some added white liquor and heated up, it is re-circulated into the digester through the injection port located on the downcomer.  V  <$> M C heater A M C pump „ -  M C screen  White liquor addition  Figure 6.4 The extraction/injection ports of the modified cooking circulation ((g) ~ extraction port; (§) - injection port)  106  Chapter 6 Simulation of an Industrial Continuous Digester  The selection of the relevant data from Appendix B is given on Table 6.8. Table 6.8 Data relevant to the modified cooking circulation (09-MAR-2004) Description  Column in Appendix B  Recorded data  Unit  R  White Liquor Flow to MC, Q L MC  4.68  L/s  BF  Flow rate of MC Extraction, Q  120.1  L/s  BP  MC Residual Alkali, E A c  14.2  g/L  W  MC  M  Note: Suppose pi;=1050 q  The mole concentration of EA at the extraction port is then, 14.16 = 0.354 mol/L 40  FA  [ E A ]  = J^ML W n  (6.17)  ISA  The mole concentration of EA in the white liquor [EA] is given in Equation (6.15); 0  at the injection port of MC circulation loop, we have, [EA] = 3  Q [EA] +Q _ [EA] MC  2  WL  MC  0  QMC+Q, WL MC  120.08x0.354 + 4.68x2.171 120.08 + 4.68  = 0.422 mol/L (6.18)  6.4.3 Boundary conditions for the energy equation  To compute the energy equation, an outflow boundary condition is set at the bottom exit, while the inflow boundary conditions are imposed at all the remaining. The downcomer and digester walls are assumed to be adiabatic, thus zero heat flux condition is used.  107  Chapter 6 Simulation o f an Industrial Continuous Digester  6.5 Simulation results 6.5.1 Boundary values and model constants Because of the symmetrical geometry and operating conditions, the HSPP digester can be well simulated by two-dimensional computations. In Appendix B, 15 sets of data obtained in stable operations of HSPP digester are listed. For each of these cases, the boundary values are determined with the procedures described in Chapter 6, section 4. As an example, Table 6.9 summarizes the computed flow rates of the liquor phase for the digester operation run on March 9, 2004; the temperature and effective alkali concentrations in the corresponding streams are also shown. Table 6.9 Computed boundary values for the case run on March 9, 2004 Liquor flow (kg/s)  Temp (°C)  EA concentration (mol/L)  245.8  157.3  0.4998  L i q u o r upflow v i a bottom port (in)  58.6  151.2  -  Wash liquor v i a downcomer (in)  105.0  160.6  0.254  M C circulation v i a downcomer (in)  131.0  162.6  0.422  Extracted wash liquor (out)  100.1  153.2  0.159  Extracted M C circulation (out)  126.1  159.9  0.354  B l a c k liquor v i a screen (Out)  146.5  159.8  0.260  B u l k cooking liquor circulation (out)  158.0  157.7  0.320  Liquor S t r e a m Liquor flow v i a top inlet (in)  Chip phase inflow: 107.43 kg/s; Production rate: 1267.55 A D T / D a y * ; [  Other conditions and assumptions  ]  Lignin/carbohydrate contents: 23.93%, 71.07%; Assume 10% o f lignin and carbohydrates are dissolved in the impregnation vessel.  [  * Note: The term bone-dry ( B D ) or oven-dry ( O D ) , means pulp from which all water has ]  been removed by heating, but this condition is seldom obtained or sought in the pulp m i l l . A i r - d r y ( A D ) pulp is understood to contain 10 percent moisture in mass; " A D T / D a y " means tons o f airdry pulp per day.  For the chip phase, the computed inflow rate via the top of the digester is 107.43 kg/s. The slip condition is applied on both inner and outer walls with p=0.85. The details of the "slip model" are given in Equation (2.43) and (2.44). To describe the non-Newtonian behavior of the chip phase, the power law model in Equation (2.38) is used with K2 and n equal to 300 Pa • 5 and 0.8 respectively. 1  2  108  Chapter 6 Simulation o f an Industrial Continuous Digester  6.5.2 The mesh As discussed before, two-dimensional computations are adequate for the HSPP digester as no circumferential flow is expected under normal operating conditions. To model the axisymmetrical flow, one can either solve the governing equations in cylindrical coordinates with a plane mesh or solve the equations in Cartesian coordinates with a special mesh. Here, the second approach is used. A structured two-dimensional mesh of 40x200 is generated for a small sector of the digester by the developed code, as shown in Figure 6.5. The angle of the sector can be arbitrary, but should be reasonably small so that the curvature of the digester wall can be well fitted by the mesh cells. Using such a mesh is advantageous, as it eliminates the need to solve equations in the cylindrical coordinates. This approach has also been used in FLUENT, a commercial CFD software. In 6.5, the mesh is shown in reduced size for the viewing purpose. | Chips & Liquor  Figure 6.5 The structured 2-D mesh shown in reduced mesh size for better viewing (The original mesh size is 40x200)  109  Chapter 6 Simulation of an Industrial Continuous Digester  6.5.3 Computational results With the determined boundary conditions and the generated mesh, the computed results are given in Figure 6.6. uiimiMiuuuu •umii milium Miimiiituuui  '5  15fc  12  3  4  r, m  (a)  (b)  r, m  (d)  (e)  110  r, m (C)  Chapter 6 Simulation of an Industrial Continuous Digester  451 EA (g/L)  161.9 161.5 161.1 160.7 160.4 160.0 160.0 159.6 159.3 158.9 158.5 158.1 157.8 157.4 157.0 156.6 156.3 155.9 155.5 155.2 154.8 148.0 132.3 122.2 102.4  'as  12 3 4 r, m  (f)  47.3 46.1 44.8 43.5 42.3 41.0 39.8 38.5 37.3 36.0 34.8 33.5 32.2 31.0 29.7 28.5 27.2 26.0 24.7 23.5 22.2 21.0  'a  r, m  (g)  kappa #  •— 111.0 — 106.4 — 101.7 97.1 9Z5 87.8 83.2 78.6 73.9 69.3 64.7 60.0 55.4 50.8 46.1 41.5 36.8 32.2 27.6  O)  r, m  (h)  Figure 6.6 Computational results for the operating case from HSPP on March 9, 2004 (a) Liquor phase velocity; (b) Liquor phase streamline; (c) Liquor phase dynamic pressure; (d) Chip phase pressure; (e) Volumefractionof the chip phase; (f) Temperature; (g) Effective alkali; (h) kappa number distribution.  6.5.3.1 Distributions of pressure and volume fraction Figure 6.6 (a), (b) and (c) show the distributions of the velocity, streamline and pressure of the liquor phase; (d) gives the pressure distribution in the chip phase; and (e) shows the volume fraction of the chip phase. Figure 6.6 (f) shows the temperature distribution for both phases, as thermal equilibrium is assumed between the chip phase and liquor phase. Figure 6.6 (g) shows the mole concentration of effective alkali (over the total liquid in the unit of mol/L), which is high in the central downcomer area and low near the digester periphery. Figure 6.6 (h) shows the distribution of kappa number, thus lignin content, as these two quantities are directly proportional to each other.  Ill  Chapter 6 Simulation of an Industrial Continuous Digester  As can be seen, there is a relatively larger pressure gradient near the liquor injection and extraction ports in both free liquor and chip phases. This indicates a strong interaction between the liquor and chip phase exists. To see this clearly, the liquor velocity and chip phase pressure near the modified cooking circulation level are plotted in Figure 6.7.  E  a  r, m  r, m  (a) Liquor phase velocity (b) Chip phase pressure Figure 6.7 The liquor phase velocity and the chip phase pressure near the modified cooking circulation level As shown in Figure 6.7 (a), the free liquor has a very high velocity near the injection port due to a small surface area of open nozzles in the downcomer; it spreads quickly in the presence of the compacted wood chips. Figure 6.6 (a) and (b) indicate that a large portion of the injected liquor moves down and outward, being extracted via the modified cooking (MC) screen; the rest of the injected liquor moves upward and is extracted through the black liquor extraction screen. While flowing out through the porous chip phase, the free liquor gradually loses its pressure; the resulting pressure gradient in the liquor phase then raises the chip phase pressure along that direction. Therefore, in the horizontal direction of the injection level, the chip phase pressure is higher on the right end and lower in the injection port, as shown in Figure 6.7 (b). A similar situation happens near the black liquor extraction port, as shown in Figure 6.8. It can also be noted  112  Chapter 6 Simulation of an Industrial Continuous Digester  from Figure 6.8 (a) that the compacted wood chips greatly enhances the momentum transfer between phases, resulting in a very flat velocity profile even in the wall vicinities. 36 DUiimiiiiiiiiuiiiii piimuiuiuiui uiiiinim 35  siiiuuimiiuuiiP niuiiimmiimupi 34 niiimmmumiuw uuuumwuuuup BMimummuuwu  33 F  tumww\ \ \\  32 F 31b• i n / / / ; / / / /  / / /  ////////ML  tmtuin ttttffltmk  ?;  30 F  'ttlllttttttttlJI  29 m m r, m  r, m  (a) Liquor phase velocity  (b) Chip phase pressure  Figure 6.8 The liquor phase velocity and the chip phase pressure near the black liquor extraction port As introduced earlier, wherever a liquid flows across the porous media - the chip phase - there is a pressure drop. It is easy to understand that the pressure of the liquor phase is high in the top and bottom, low in the middle; high near the injection ports, low near the extraction screens. The chip phase moves downward continuously, while the liquor phase flows downward in the upper part, and moves mostly upward in the lower part of the digester, the velocity difference between the liquor and chip phases in the lower part is larger than that in the upper part, and therefore the interphase drag is larger. Also the upflowing liquor flows through a longer path before being extracted. Therefore the liquor phase pressure at the bottom is much higher than at the top of the digester, as shown in Figure 6.6 (c).  113  Chapter 6 Simulation of an Industrial Continuous Digester  The chip phase pressure in the lower part of the digester is decreased significantly by the counter-current liquor flow from the bottom. Nevertheless, the chip phase pressure is still high in the bottom while low at the top due to the body force, as shown in Figure 6.6 (d). Consequently, the volume fraction of the chip phase is also high in the bottom, as given in Figure 6.6 (f). 6.5.3.2 Velocity of the chip phase  Figure 6.9 presents the chip phase velocity in vector form. As can be seen, the chip phase exhibits close to a plug-flow behavior. In the upper part of the digester, the chip phase volume fraction is small, as a result, the cross section area that the chip phase flows through is small, and therefore the downward velocity of the chip phase is large; in the lower part of the digester, the chip phase volume fraction is large so the downward velocity of the chip phase is small. In Figure 6.9 (b) and (c), the detailed views of the chip phase velocity in two different locations are shown. Obviously, the chip phase has a slippage on the digester wall. 1IIIUIIIIIIIUIHI  iiimiiiiiuiaj •imiiiiiiiiu •mmiuuu  24 23  liimiuiiiuiBl iiiuumtmuiB inimiimiwi aiiiiiinuuiKi.i 35IIIIUII m u u i u liiiiiimiiumiw wiiiiiiiiiiiiiiimw iiiiiiiiiiniiiiiiiii'ii iiiiiiiiiiiiiuiiiimii  E 22 •J 21  a>  in1111 iiiiiiiniiil  T  UWIIIIIIIIlllllllllll  aumiiiuiiiiiiiiiiii  20  llllil i mum  19  III 11111 iinmiiiil  18  WllUlllllllllllllllll UlliltlilHHIIllllllll  1 2  3  4  r, m  mum  Ul  wwwm  (b) High chip pressure regime  'a  X 20 milium lllllllli  4|  liniiiiii iiiiiiinii  E  s £  at  3 2  3  1 2  4  r, m  r, m  (a) Overall chip phase velocity  (c) Near the digester bottom  Figure 6.9 Velocity distribution of the chip phase  114  Chapter 6 Simulation o f an Industrial Continuous Digester  In Figure 6.10, the predicted velocity profiles of the chip phase at three different digester levels are shown. As can be seen, the averaged downward velocity is decreased by about 20% from the top portion to the middle part of the digester; this is mainly due to an increased chip phase volume fraction as shown in Figure 6.6 (f). In the lower portion of the digester, the volume fraction of the chip phase further increases, and the mass of the chip phase decreases because of the increased delignification rate. As a result, the averaged downward velocity of the chip phase further decreases, being 0.2mm/sec lower at the bottom than at the middle of the digester.  1 — T o p (Z = 35.1 m) 2 — Middle (Z = 15.8 m) 3 — Bottom (Z = 5.0 m)  -2.5 -  L_i  i  i  A  l 1.0  i  i  i  l  i  i  2.0  '  i  i  i  l 3.0  i  i  i  i  LJ 4.0  r, m  Downcomer Figure 6.10 Chip phase velocity profiles at different digester levels (Z — digester level) In the upper portion of the digester, the chip phase velocities near the outer wall of the digester shell are slightly smaller than those in the center region, as indicated by curve 1 in Figure 6.10; in the middle part of the digester, this trend is augmented (curve 2); but in the lower portion of the digester, the velocity curve 3 becomes flat. This is due to the friction force between the chip phase and the digester wall. As plotted in Figure 6.5 (e),  115  Chapter 6 Simulation o f an Industrial Continuous Digester  the chip phase pressure is relatively high in the upper portion and becomes higher in the middle of the digester. Therefore the friction force at the digester shell retards the chips near the wall, causing the chip phase velocities to decrease in the vicinity of the wall. In the lower portion of the digester, the chip phase pressure is significantly decreased by the upflowing liquor and therefore the friction force at the wall is almost negligible. Further, chips near the digester shell gradually accelerate to the velocities of the chips in the middle entrained by the shear stress caused by the neighboring chips. Near the downcomer, the situation is essentially the same, but the surface area of the downcomer outer wall is much smaller than the digester wall. Consequently, the friction force is small and the velocity difference is minimal.  6.5.3.3  Kappa number variation at the bottom exit  Figure 6.11 presents the prediction of kappa number distribution at the bottom exit of the digester. The kappa number near the downcomer and the digester wall are lower than that in the middle, indicating a non-uniform cook; the total kappa number variation is about 9.  24  I.  i 11 i i 11 * i i 111 i i i 11 i i 111 i i 111 i i 11 11 i 11 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0  il  r, m  Figure 6.11 Kappa number variation at the bottom exit of the digester  116  Chapter 6 Simulation of an Industrial Continuous Digester  The cooking extent is determined by temperature, liquor concentration and residence time. Since all the circulating liquor was injected from the central downcomer, the concentration of effective alkali and temperature are high in the downcomer vicinity and low near the digester wall, as can be seen in Figure 6.6 (d) and (h). However, Figure 6.11 indicates that the chips near.the digester wall have lower kappa numbers and are thus over-cooked. Therefore the dominating factor of the non-uniform cook in this particular case is probably the residence time. To confirm this assumption, the total residence time of the chip phase is calculated based on its flow field and plotted against the radial locations; as shown in Figure 6.12, the chips near the digester wall have much longer residence time which allows further delignification and results in lower kappa number.  3  o £ 0)  E  o o c  0)  •a '53 £  a.  !E o  0.5  1.0  1.5  2.0  2.5  3.0  3.5  4.0  r, m  Figure 6.12 Chip phase residence time versus chip locations  117  Chapter 6 Simulation of an Industrial Continuous Digester  6.6 Discussion of simulation results 6.6.1 Mesh convergence check To ensure a mesh converged solution, three meshes are used for the same cooking case and no obvious difference is found either in the resulting flow fields or the scalar distributions. For simplicity, only the downward velocities of the liquor phase along the vertical line r=2m are shown, as in Figure 6.13. As can be seen, the increased mesh size does not significantly influence the liquor phase velocity, and computations using a mesh of 40x200 should give reasonably good predictions.  Liquor velocity, mm/s  Figure 6.13 Computed downward liquor phase velocities along the vertical line r=2m using different meshes 6.6.2 Predictions of kappa number and yield To check the code's prediction capability, two-dimensional computations have been conducted in 15 available cooking cases, with corresponding cooking conditions listed in Appendix B. In these cases, the production rates range from about 1000 to 1547 tons airdry-pulp per day; other operating parameters change from case to case. Operations in all these cases are fairly stable over a minimum period of 6 hours. In Table 6.10, the computed kappa number and yield, and the corresponding mill records are listed. As can be seen, the prediction errors of kappa number and yield are less than 10%.  118  Chapter 6 Simulation of an Industrial Continuous Digester  Table 6.10 Computed kappa # and yield versus mill records Production Rate (ADT/d)  Chip Meter (RPM)  1  12686  15.06  29.63 /31.40  5.97  46.4/48.2  3.88  2  1313  15.01  27.99 / 29.05  3.79  48.2/50.7  5.19  ->  888  10.40  28.35 /30.86  8.85  46.7/48.7  4.28  4  1317  15.66  27.99/26.74  -4.47  46.4/48.2  3.21  5  1237  15.37  30.99/32.65  5.36  44.4 / 46.4  4.50  6  1223  14.95  29.35 /31.10  5.96  45.1 /47.8  5.99  7  1265  14.53  29.25 /31.31  7.04  48.0/47.2  -1.67  8  1410  14.89  28.26/29.80  5.45  52.2/49.0  -6.13  9  1128  13.79  35.66/33.95  -4.80  45.0/48.9  8.67  10  1100  12.85  40.94/37.36  -8.74  47.1 /49.1  4.25  11  937  11.07  29.42/30.38  3.26  46.4/48.9  5.38  12  1058  11.95  31.97/32.01  0.12  48.6/52.0  7.0  13  1398  14.79  32.27/31.71  -1.74  52.1 /50.3  -3.45  14  1084  13.33  32.29/29.06  -10.0  44.7/43.2  3.36  15  1190  15.38  35.07/33.13  -5.53  42.7/44.0  3.04  No  Yield, («Vo)  kappa #  (Record/Comp.) Error (%) (Record/ Comp) Error (%)  6.6.2 Variation of production rate In daily operations, the production rate of a digester may change as desired; in some pulp mills, digesters are often operated well above their original design capacities. To know how the variations of production rate affect the exit kappa number and yield is helpful for operational personnel. Figure 6.14 shows the simulated results with production rate increasing from 1075 to 1455 tons/day. As plotted, both kappa number and yield increase as production rate increases. In the computed production range, the kappa number increases from 25 to 43, while the yield increases from 47.5% to about 51.7%). The chip phase pressure at digester bottom exit drops from 8.8 to 4.7 kPa, as shown in Figure 6.14 (b).  119  Chapter 6 Simulation of an Industrial Continuous Digester  1100  1200  1300  1100  1400  1200  1300  1400  Production rate, ton/day  Production rate, ton/day  (b) Pc at the bottom exit plane  (a) Yield and kappa number  Figure 6.14 Cooking performance affected by variation in production rate (Except production rates, data from HSPP, March 9, 2004)  6.6.3 Cooking performance affected by the chip type In the coastal area of British Columbia, many pulp mills use different types of wood chips or blended chips to meet the diverse needs of the pulp market. Shifting from one type of chips to another often causes significant variation in the exit kappa number. Without adequate knowledge, it is very difficult to make an appropriate and quick adjustment in order to meet the target kappa number. For different wood chips, the lignin content is different, the chip compressibility and flow resistance over different chips are also different, possibly leading to differences in cooking. Table 6.11 lists typical chip combinations used in HSPP; the corresponding contents of lignin and carbohydrates are computed based on data in Table 6.6.  120  Chapter 6 Simulation of an Industrial Continuous Digester  Table 6.11 Typical combinations of incoming wood chips in HSPP Hemlock Chip combination (%) 60 1 2 50 3 4 40 50 5 30 6  Cedar (%) 40 45  Fir (%)  Spruce (%) 55 50 60  50 70  Ligo  Carbo  21.16 23.80 24.05 24.54 25.60 27.20  70.78 71.28 70.25 70.32 69.55 69.41  Simulated cooking effects are shown in Figure 6.15 when the above six chip combinations are used. In these simulations, the same operating conditions from March 9, 2004 in HSPP are used; Harkonen's correlations [7] are used to calculate the chip compressibility and flow resistance. 40 • A  38  Kappa # Yield(%)  36 34  .a E 3 C  32  re o. a. ns  > '< in  30 28 26 24 22 20  J  I 2  I 3  I 4  I 5  L 6  Chip combinations  Figure 6.15 Cooking performance affected by the chip type (See Table 6.11 for chip combinations) From Table 6.11 and Figure 6.15, the initial lignin content increases from chip combination 1 to 6, and the exit kappa number also increases, from about 24.5 to 33; the yield on the other hand, decreases from about 47.5 to 45. In past years, extensive work [7,29,30,31] has been published on chip compressibility. It is of interest to investigate the effect of chip compressibility on cooking performance.  121  Chapter 6 Simulation of an Industrial Continuous Digester  In Chapter 2, the compressibility of different types of wood chips is listed; the concept of "equivalent diameter" characterizing the chip compressibility is also given. To address the effect of chip compressibility on the cooking performances, simulations have been carried out for different cases, where all the cooking conditions are the same as in the case from March 9, 2004. The six softwood chips listed in Table 2.2 are used with the equivalent diameter ranges from 0.85 to 3.0 mm. As seen in Figure 6.16 (a), when the equivalent diameter increases, both kappa number and yield increase at the beginning and decrease when the equivalent diameter exceeds about 1.4mm; the averaged residence time for wood chips exhibits the opposite trend as plotted in Figure 6.16 (b). Because longer residence time results in more complete cooking, lower kappa number of wood chips, the residence time is probably the immediate cause for kappa number changes. From Figure 6.15 and 6.16, variation in chip combination in HSPP and the equivalent diameter affect the kappa number almost to the same extent, with both leading to kappa number variation ranging from about 24.5 to 33.  •  Kappa #  •  Y i e l d (%)  49 •  •  •  •  2  o >-  H 47  1.0  1.5  2.0  2.5  3.0  46  1.0  1.5  2.0  2.5  3.0  Equivalent diameter, mm  Equivalent diameter, mm  (b) Averaged residence time  (a) kappa number and yield  Figure 6.16 Cooking affected by equivalent diameter of wood chips  To illustrate how the chip residence time is affected by the equivalent diameter of wood chips, the chip phase pressure distributions corresponding to three different  122  Chapter 6 Simulation of an Industrial Continuous Digester  equivalent diameters are shown in Figure 6.17; the chip phase volume fractions are given in Figure 6.18. When the equivalent diameter increases, the flow resistance over wood chips decreases. Consequently, the upflowing liquor decreases the solid pressure less. Therefore with bigger equivalent diameters in cases of Figure 6.17 (b) and (c), the chip phase pressure in the lower part of the digester is higher than that in the case of Figure 6.17(a). The volume fractions of the chip phase shown in Figure 6.18 are determined by the chip phase pressure and chip compressibility. Higher chip phase volume fraction means larger surface area that the chip phase can flow through, thus lower chip phase velocity. Compared with case (b), case (a) has higher chip phase volume fraction in the upper part of the digester, while case (c) has higher chip phase volume fraction in the lower part. As a result, chips in case (c) have the longest residence time, while chips in case (b) have the shortest residence time.  123  Chapter 6 Simulation of an Industrial Continuous Digester  1 2 3 4  1 2 3 4  r, m  r, m  (a) D = 0 . 8 3 8 mm  1 2 3 4  r, m  (b) D = l . 3 9 4 mm  eq  eq  (c) Deq=3.03 mm  Figure 6.17 Chip phase pressure affected by the equivalent diameter  . "3 I  0.546 0.540 0.535 0.529 0.524 0.518 0.513 0.507 0.501 0.496 0.490 0.485 0.479 0.474 0.468 0.463 0.457 0.451 0.446 0.440 0.435 0.429 0.424 0.418  25  20  r, m  (a) D = 0 . 8 3 8 mm eq  • I  0.546 0.540 0.535  W  0.529  I  ^ 0.524 •-' 0.518 0.513 0.507 0.501 0.496 0.490 0.485 0.479 0.474 0.468 0.463 0.457 0.451 0.446 0.440 0.435 0.429 0.424 0.418  ffl  X  0.546 0.540 0.535 0.529 0.524 0.518 0.513 0.507 0.501 0.496 0.490 0.485 0.479 0.474 0.468 0.463 0.457 0.451 0.446 0.440 0.435 0.429 0.424 0.418  20|  r, m (b) Deq=l.394 mm  r, m (c) D = 3 . 0 3 mm e q  Figure 6.18 Chip phase volume fraction affected by the equivalent diameter  124  Chapter 6 Simulation of an Industrial Continuous Digester  6.6.4 Effect of the cold blow rate In HSPP, the cold blow rate is adjusted to maintain the free liquor level and control the pulp consistency (the percentage of fiber in the pulp solution) at the blow line. Figure 6.19 shows the variation of the kappa number and yield when the cold blow rate increases; in the simulations, the black liquor extraction increases accordingly so as to maintain steady state operations. 52  36 Kappa # Yield (%) 34  re a.  51  *X 'S a>  32  re a. o. re  50  E & o .Q +J  2  30  re  a."  49 28  26  J 55  i  L 65  75  J 85  i  l_ 95  48 105  Cold blowrate,kg/s  Cold blowrate,kg/s  (a) Kappa number and yield (b) Chip phase pressure at bottom exit Figure 6.19 Cooking performance affected by the cold blow rate (Except cold blow rate, data from HSPP, March 9, 2004) As plotted in Figure 6.19 (a), both the kappa number and the yield increase with the cold blow rate; the chip phase pressure at the digester bottom decreases considerably, from about 7 to 2 kPa, as shown in Figure 6.19 (b). If the cold blow rate increases further, the chip phase pressure at the digester bottom may approach zero, which means wood chips lose contact with each other. In this situation, the whole chip column moves quite slowly or even stops, leading to the occurrence called "hang-up" that should be prevented.  125  Chapter 6 Simulation of an Industrial Continuous Digester  6.7 Three-dimensional simulation In large digesters, chips may plug the extraction screen when the extraction flow is suddenly changed or the flow rate is too high. As a remedy, the black liquor and modified cooking liquor extraction screens are divided into two parts and switched periodically from one part to the other. The extractions can be switched from either "up-and-down" or "left-and-right" mode. In the first mode, the upper half and the lower half of the screens are used alternatively; in the second mode, the extraction screens are divided into a number of sections and the liquor extractions are shifted in the circumferential direction. As a hypothetic case given in Figure 6.20, both the black liquor and modified cooking liquor extraction screens are divided into 12 equal sectors; the liquor extractions are shifted in the circumferential direction. In this situation, the uneven extraction of liquor may result in different compaction on the chip phase; it may also affect the temperature distribution near the extraction screens, and thus may contribute to a nonuniform cooking. To investigate this effect, a three-dimensional simulation is performed for the digester in HSPP.  (a) Extraction at time T  (b) Extraction at time T+At  Figure 6.20 Liquor extraction shifted between screen sectors In the steady state computations, extractions using half of the screen are modelled; the dynamic behavior when extraction screens are shifted is not considered. Because of the symmetrical condition, the computation only considers one sector enclosed by OA  126  Chapter 6 Simulation of an Industrial Continuous Digester  and OB, but the results in a cross section view are extended to the region enclosed by OA and OC, as shown in Figure 6.20 (b). In Figure 6.21, the structured mesh used in the computation is shown; the liquor flow configuration is also illustrated.  (a) Overall mesh and liquor flow  (b) Top view at black liquor extraction plane  Figure 6.21 The structured mesh shown in reduced size (original size:30x20xl80) (BCL: bulk cooking liquor; BL: black liquor; MCL: modified cooking liquor; WL: wash liquor.) Again, the computation is performed for the operating case on March 9, 2004 from HSPP. Figure 6.22 presents the velocity and pressure distributions for the liquor phase. As can be seen in Figure 6.22 (b) and (c), higher velocity always corresponds to a larger pressure drop; when the circulation liquor is injected from the downcomer at high velocity, the liquor phase pressure drops significantly. Therefore, the total pressure drop in the liquor phase at the modified cooking liquor injection level is much larger than that at the black liquor extraction level.  127  Chapter 6 Simulation of an Industrial Continuous Digester  (c)  (c)  Figure 6.22 Velocity and pressure distributions for the liquor phase (a) Front view in the centered symmetrical plane; (b) cross-section view at the black liquor extraction level (z=31.5m); (c) cross-section view at the modified cooking liquor injection level (z=21.95m).  128  Chapter 6 Simulation of an Industrial Continuous Digester  Figure 6.23 shows the streamline in the liquor phase and the chip phase pressure. In Figure 6.23 (b), the chip phase pressure distribution in the central symmetrical plane is very similar to the distribution obtained in a two-dimensional simulation, as given in Figure 6.6 (d). From the streamline of the free liquor phase and the pressure contours in Figure 6.23 (c) and (d), the liquor velocity contributes greatly to the chip phase pressure distribution. Near the black liquor extraction screen, the extracted liquor tends to push the chip column outward, causing a high chip phase pressure near the extraction port. At the injection level of the modified cooking liquor, the free liquor moves down- and outward, resulting in a continuous rise in the chip phase pressure. The magnitude of the injection velocity is high but decreases quickly; therefore the chip phase pressure increases quickly near the injection port and increases slowly afterwards, as given in Figure 6.23 (d).  (a)  (b)  Figure 6.23 Liquor phase streamline and the chip phase pressure (a), (b) velocity and pressure in the center symmetrical plane; (c), (d) pressure in cross-section view at the black liquor extraction level (z=31.5m) and the modified cooking liquor injection level (z=21.95m).  129  Chapter 6 Simulation of an Industrial Continuous Digester  Figure 6.24 contains the chip phase velocity distributions. In Figure 6.24 (a), the chip phase exhibits a plug behavior as a whole, but near the digester wall, chip motion is retarded because of the friction force from the digester wall. Chip phase velocity profiles at the black liquor extraction level and wash screen level are given in Figure 6.24 (b). The chip phase velocity at the digester wall has a larger gradient at the black liquor extraction level than at the wash screen level, as the chip phase pressure is higher in the former case.  mi. . nil: 35^lllll!IIU!IH. iiiuimnuA llllllllllimf iiiiiiiiiinimii 30 Liilllllliianiii  IIIIIIIIIIIIIIII .IIIIIIIIIIIIIIII  25  1 — black liquor extraction (Z = 31.70 m) 2 — wash screen (Z =2.05 m)  JZ  lilt CT iiiiiiiiiuiiiii l l l l l l l l l l U l ll»\ < 'U X 20 J I I I I M 11 iir  Illllll 111 nun  IIIIIIIIIUIIIII IIIIIIIIIIIIIIII  15)111111111 iiniii  minimi  -4.0  m i n i 11 l i m n  IIIIIIIIIUIIIII IIIIIIIUIIIII 10I IIIIIIII 11 iiniii IIIIIII 11 IIIIIII IIIIIII 11 IIIIIII l l i i i i i i 11 m i n i  IIIIIIIIIIIIIIII IIIIIII 11 iiiiiui IIIIIII i I I I I I I I IIIIIIIJIII.IIII, 12 3 4  r, m  (a)  (b)  Figure 6.24 The velocity of the chip phase (a) velocity in the center symmetrical plane; (b) Axial velocity profiles at two different places. Figure 6.25 presents the distributions of volume fraction of the chip phase, temperature and kappa number. In the center symmetrical plane, the results in Figure 6.25 (1) are all close to those obtained in a two-dimensional computation, as given in Figure 6.6 (f), (g) and (h). The contour plots in the cross section views show that for compaction, temperature and kappa number, the difference in the radial direction is dominant while the difference in the circumferential direction is negligible. This indicates that the nonuniform extraction causes only local variations, and therefore is not a main factor in nonuniform cooking.  130  Chapter 6 Simulation of an Industrial Continuous Digester  Kappa #  <°C)  '5  0.513 0.508 0.503 0.499 0.494 0.489 0.485 0.480 0.475 0.470 0.466 0.461 0.456 0.452 0.447 0.442 0.438 0.433 0.428 0.424 0.419 0.414 0.409 0.405  201  I 20  r, m  (a)  161.9 161.5 161.1 160.7 160.4 160.0 160.0 159.6 159.3 158.9 158.5 158.1 157.8 157.4 157.0 156.6 156.3 155.9 155.5 155.2 154.8 148.0 132.3 122.2 102.4  r, m  r, m  (b) (1) In the center symmetrical plane  (c)  T(°C) Kappa # 96.32 95.81 95.29 94.76 94.27 93.75 93.24 92.72 92.21 91.70 91.18 90.67 90.15  (a)  (b)  (c)  (2) Cross section view at the black liquor extraction level (z=31.5m) T(°C)  1.0  2.0 3.0  Kappa # 30.95 30.05 29.15 28.25 27.35 26.45 25.55 24.65 23.74 22.84 21.94 21.04 20.14  4.0  (b) (3) Cross section view at the bottom exit  (c)  Figure 6.25 Computed distributions of scalar variables (a) chip phase volume fraction; (b) temperature; (c) Kappa number. 131  Chapter 6 Simulation of an Industrial Continuous Digester  6.8 Summary This chapter presents the computational study on an industrial digester in Howe Sound Pulp and Paper Ltd. Firstly, the boundary conditions for both phases are described; boundary values, including inlet flow rates and species concentrations are determined. Two-dimensional simulations are then conducted for 15 data sets for stable operations of the digester in Howe Sound Pulp and Paper Ltd. The predicted kappa number and yield agree reasonably well with the mill records. The computed distributions of pressure, velocity, volume fraction etc. from a typical simulation are shown. The interactions between the free liquor and chip phases are discussed. Further efforts were aimed at the investigation of the influence of different parameters on the pulping process. The following was concluded: 1. Both yield and kappa number increase with production rate. However the chip phase pressure at the digester bottom exit decreases very significantly when the production rate increases. 2. Wood chips differ considerably in lignin contents but not in carbohydrates. When wood species changes, cooking conditions should be changed accordingly. Chips with higher lignin content require longer residence time. 3. Chip compressibility greatly affects the cooking performance. Based on simulations in typical cooking conditions, chips with an equivalent diameter of 1.4 mm have the shortest residence time in the digester and are least cooked. Smaller chips have higher compaction in the upper part, while bigger chips have higher compaction in the lower part of the digester, therefore both have a longer residence time and are more cooked, resulting lower yields and kappa numbers. 4. Increasing the cold blow results in higher kappa number and yield but lower chip phase pressure at the digester bottom. Because of this, a too high cold blow could result in "hanging chip column".  132  Chapter 6 Simulation of an Industrial Continuous Digester  6.8 Summary This chapter presents the computational study on an industrial digester in Howe Sound Pulp and Paper Ltd. Firstly, the boundary conditions for both phases are described; boundary values, including inlet flow rates and species concentrations are determined. Two-dimensional simulations are then conducted for 15 data sets for stable operations of the digester in Howe Sound Pulp and Paper Ltd. The predicted kappa number and yield agree reasonably well with the mill records. The computed distributions of pressure, velocity, volume fraction etc. from a typical simulation are shown. The interactions between the free liquor and chip phases are discussed. Further efforts were aimed at the investigation of the influence of different parameters on the pulping process. The following was concluded: 1. Both yield and kappa number increase with production rate. However the chip phase pressure at the digester bottom exit decreases very significantly when the production rate increases. 2. Wood chips differ considerably in lignin contents but not in carbohydrates. When wood species changes, cooking conditions should be changed accordingly. Chips with higher lignin content require longer residence time. 3. Chip compressibility greatly affects the cooking performance. Based on simulations in typical cooking conditions, chips with an equivalent diameter of 1.4 mm have the shortest residence time in the digester and are least cooked. Smaller chips have higher compaction in the upper part, while bigger chips have higher compaction in the lower part of the digester, therefore both have a longer residence time and are more cooked, resulting lower yields and kappa numbers. 4. Increasing the cold blow results in higher kappa number and yield but lower chip phase pressure at the digester bottom. Because of this, a too high cold blow could result in "hanging chip column".  132  Chapter 6 Simulation of an Industrial Continuous Digester  In addition, a three-dimensional simulation is carried out. The computation shows that the uneven liquor extraction investigated and used in many mills does not cause significant non-uniformity in cooking.  133  Chapter 7 Summary and Conclusion  Chapter 7 Summary and Conclusion In this thesis, a general computational code is developed to simulate the complex cooking process in industrial continuous digesters. In the following, the main contributions and conclusions are listed. A set of comprehensive mathematical equations governing the two-phase flow, heat transfer and chemical reaction in continuous digesters are presented; sub-models accounting for the inter-phase friction, thermal dispersion and chip compressibility are also discussed. In addition, a new, more reliable cooking reaction model is derived based on extensive experimental data. As a main advantage, this model can be used in various cooking conditions without the need to pre-determine some model constants, which in other models are dependent on operating conditions. Unlike most existing pulping models, this model abandons the piece-fitting method, thus parameters characterizing the cooking stages are no longer used, which is another attractive feature of this model. Despite its simple form, the model predictions agree well with a wide range of cooking experimental data. A three-dimensional two-phase CFD code based on a Finite Volume Method is then developed to calculate the flow fields for both phases; the volume fraction of solid phase, the concentrations of chemicals and the residual contents of lignin and carbohydrates are also calculated. The code is based on collocated curvilinear grids and is of second order accuracy; it has the capability of solving many different steady state laminar flow problems with or without the involvement of heat transfer or chemical reactions. The liquid flow in a laboratory digester with chips filled in is simulated using the developed code. Computations indicate that different flow patterns can be formed depending on the circulation rate. Simulations show that the compacted wood chips greatly  134  Chapter 7 Summary and Conclusion  enhance the momentum dispersion of the liquor phase; therefore the liquid flow tends to be uniform away from the injection and extraction places.  Simulations are performed for a full-scale industrial digester for numerous operating conditions. The predictions of yield and kappa number at the digester exit are fairly accurate for various cooking conditions. Further simulations lead to the following conclusions: 1. With increased production rate, yield and kappa number increase; however the chip phase pressure at the digester bottom exit decreases very significantly; wherever the chip phase pressure approaches zero, chips lose contact with each other and "hanging chip column" may occur. 2. Increased cold blow results in higher kappa number and yield, but lowers the chip phase pressure at the digester bottom. 3. Chip compressibility affects significantly the cooking performance. Chips with an equivalent diameter of around 1.4 mm exhibit a minimum residence time and are least cooked, therefore giving the highest yield and kappa number. A three-dimensional simulation is carried out to investigate the circumferential variation of different parameters. The computation shows that uneven liquor extraction investigated and used in some mills does not cause significant non-uniformity in cooking. The developed three-dimensional tool constitutes a reliable and robust tool to investigate industrial digesters. It is a tool that can be used to assist the design and optimization of the pulping process.  135  Chapter 8 Recommendation for Future Work The mathematical model, as well as the computational code developed in this thesis is a useful tool for the investigation and optimization of the pulping process. The tool can be further developed to improve its predictive capability and widen its use. Future work is suggested along the following lines: •  An improved model accounting for the cold blow flow injected into the digester should be developed. The actual geometry is very complex and the simplification used in this thesis allows for good over-all simulation but might miss some local details if the injection is non -uniform or occurs with very high velocity, possibly leading to channeling.  •  The energy equation should be solved separately for the chip phase and the free liquor phase in the lower part of the digester. The temperature difference between the phases could be significant locally and the evaluation of this effect would be useful.  •  Chip compressibility and flow resistance are very important for the predictive capability of the model. Therefore further work is necessary for better definition of these two important parameters, especially when blended wood chips are used.  •  The model should be extended to batch digesters. They are widespread in industry, have been enlarged and incorporate new cooking techniques, and therefore an improved understanding of the process would be a useful contribution.  136  Bibliography  1. Michelsen, F.A., A Dynamic Mechanistic Model and Model-Based Analysis of a Continuous Kamyr Digester, PH.D. Thesis, Dept. of Eng. Cybernetics, The Norwegian Institute of Technology, University of Trondheim, 1995. 2. Rydholm, S., Pulping Processes, Interscience Publishers, New York, 1967. 3. Koch Christensen, P., Treforedlingskjemi, Department of Chemical Engineering, Norwegian Institute of Technology, Trondheim, 1991. 4. Harder, N . "Extended Delignification in Kraft Cooking." Svensk Paperstidning. 81(15):483, 1978. 5. He, P., Salcudean,M., Gartshore I., Bibeau E., Modeling of Kraft Two-Phase Digester Pulping Processes, TAPPI Conference, Anaheim, Sept. 1998. 6. Pougatch, K., Salcudean, M., Gartshore, I., "A numerical model of the reacting multiphase flow in a pulp digester", submitted to Applied Mathematical Modelling, (personal communications) 7. Harkonen, E.J. (1984), A Mathematical Model for Two-Phase Flow, Acta Polytechnica Scandinavia, Mechanical Engineering Series No. 88. 8. Michelsen, F. A., Christensen, R. T., Lunde, G. G., Lundman, G. and Johansson, K., Model Predictive Control of a Continuous Kamyr Digester at SCA-Nordliner Munksund, Sweden. Preprints Control Systems 1992, Dream vs. Reality: Modern Process Control in the Pulp and Paper Industry, Sept.-Oct.1992, Whistler, B.C. Canada, 1992. 9. Pougatch, K., Salcudean, M . , and Gartshore, I., "Mathematical Modelling of Digesters", presented at 52nd Canadian Chemical Engineering Conference, Vancouver, October 20-23, 2002. 137  10. Miyanishi, T. and Shimada, H., Improvement of Pulp Strength and Yield by Computer Simulation of Lo-solids Kraft Cooking, 2001 Tappi Journal Peer Reviewed Paper. 11. Vroom, K.E., The " H " factor: A Means of Expressing Cooking Times and Temperatures as a Single Variable, Pulp and Paper Magazine of Canada, 58(3): 228 (1957). 12. Kerr, A.J., The Kinetics of Kraft Pulping, Appita 24(3): 180-188,1970. 13. Kerr, A.J. and Uprichard, J.M., The Kinetics of Kraft Pulping - Refinement of a Mathematical Model. Appita, 30:48-54, 1976. 14. Gustafson, R., Sleicher,C, Mckean,W., et. al., Theoretical Model of the Kraft Pulping Process, Ind. Eng. Chem. Process Design & Development, 22(l):87-96, 1983. 15. Christensen, T., Albright, L. and Williams,T., TAPPI 1983, Annual Meeting, Atlanta, Georgia, March 2-4, 239(1983). 16. McKibbens, S.W., Application of Diffusion Theory to the Washing of Kraft Cooked Wood Chips. TAPPI Journal, 43(10):801-805, 1960. 17. Scheiddegger, A.E., Statistical hydrodynamics in porous media, J. Appl. Phys. 25(8), 994-1001, 1954. 18. Spitz, K. and Moreno J., A practical guide to groundwater and solute transport modeling, ISBN: 0-471-13687-5, by John Wiley & Sons. Inc., New York, 1996. 19. Huang, H., Li, Z. and Usmani,A.S., Finite element analysis of non-Newtonian flow: theory and software, ISBN: 1852330244. London ; New York : Springer, 1999. 20. Biermann, C.J., Essentials of Pulping and Papermaking, ISBN 0-12-097360-X, Academic Press, 1993. 21. Gullichsen, J., Chemical Pulping, Papermaking Science and Technology, Book 6A, ISBN952-5216-06-03, Fapet Oy, Finland, 1999a. 138  22. Macdonald, I.F., et al., "Flow through Porous Media - the Ergun Equation Revisited", Ind. Eng. Chem. Fundam., Vol.18, No.3, 1979, pi99-208. 23. Ergun, S. "Flow through packed columns", Chem. Eng. Progress., Vol.48, No.2, p89, 1952. 24. N . Ahmed and D. K. Sunada, Non-linear Flow in Porous Media, Journal of Hydr Division, ASCE, Vol.95, No.6, 1969, pi847. 25. N. Ahmed and D. K. Sunada, Non-linear Flow in Porous Media, Journal of Hydr Division, ASCE, Vol.97, 1971, pl233. 26. Barnes, H.A., Hutton J. F. and Walters, K., An Introduction to Rheology, ISBN: 0444871403, Elsevier Science Publishers B.V., 1989. 27. Cross, M.M., Rheology of non-Newtonian fluids: a new flow equation for pseudo-plastic systems, J. Colloid Sci., (20), 417-437. 28. Masse, M.A., Kiran, E., and Fricke, A.L., A Thermodynamical Model of the Heat Capacity of Compositionally Complex Multicomponent Polymer Solutions: Kraft Black Liquor, Chemical Engineering Communications, 50: 81-91 (1987). 29. Wang, Z. and Gullichsen, J., Flow Compressibility and Resistance of Chips Made with a new Chipping Technique, Pro. 1998 TAPPI Pulping Conf, Montreal (Oct25-29), Vol2, p.753-760. 30. Lindqvist, J., The Effects of Flow Conditions in Kraft Cooking of Wood Chips by Various Chip Size Distribution, Masters Thesis, Helsinki University of Technology, 1994. 31. Lammi, L., The Utilization of Eucalyptus Camaldulensis-species with the Superbatch process, Masters Thesis, Helsinki University of Technology, 1996. 32. Dudgeon, C. R., An Experimental Study of Flow of Water through Coarse Granular Media, La Houllie Blanche, France, Vol.21, No.7, p.785, 1966.  139  33. Gupte, A. R., Dissertation, Fluid-Dynamic characteristics of two-phase flow, Karlsruhe, 1970. 34. Nam-Trung Nguyen, Steven T. Wereley, Fundamentals and applications of microfluidics, ISBN: 1580533434, Artech House, Boston, MA, c2002. 35. Thompson, P. A., and Troian, S.M., " A general boundary condition for liquid flow at solid surfaces", Nature, Vol.389, 1997, p. 360-362. 36. Edwards, L. and Norberg, S.E., Alkaline delignification kinetics, A General Model Applied to Oxygen Bleaching and Kraft Pulping, TAPPI J. 56(11): 108-111 1973. 37. Kubes, G.J., Fleming, B.I. and MacLeod, J.M., Viscosities of Unbleached Alkaline Pulps, II. The G factor, J.Wood Chem. Tech. 3(3): 313,1983. 38. Paulonis, M.A. et. al, Adaptive Inferential Control of Kraft Batch Digesters as Based on Pulping Liquor Analysis, TAPPI J. 74(6): 169-175,1991. 39. Smith, C.C. and Williams T.J., A Dynamic Model of a Kamyr Digester and its Application to Pulp Washing and Pulping Reaction Control, Modeling & Control Kraft Production Systems. Proc. ISA Symp. Milwaukee, Oct., 1975. 40. Smith, C.C. and Williams, T.J., "Mathematical Modeling, Simulation, and Control of the Operation of a Kamyr Continuous Digester for the Kraft Process", Report No.64, Purdue Laboratory for Applied Industrial Control, Purdue University, Lafayette, Indiana, 1974. 41. Aurora Santos, Francisco Rodriguez et. al, Kinetic Modeling of Kraft Delignification of Eucalyptus Globulus, Ind. Eng. Chem. Res, 1997, 36, 41144125. 42. Ling Yang, Keng Chuang and Shijie Liu, A mechanistic Model for Kraft Pulping, PAPTAC 90 Annual Meeting, 2004. th  140  43. Hairer E., Lubich, C. and Roche M., The numerical solution of differentialalgebraic systems by Runge-Kutta methods, ISBN: 3540518606, Berlin; New York: Springer-Verlag, cl989. 44. Aurell, R., Harder, N., Kraft pulping of pine. I. Changes in the compounds of the wood residue during the cooking process. Sven. Papperstidn. 1965, 68 (3), 59-68. 45. Kuo, K. K , Principles of Combustion, Chapter 2, pl09, ISBN 0-471-09852-3, John Wiley & Sons, Inc., 1986. 46. Bugajer, S. Kinetics of delignification reactions in the Eucalyptus pulping process (in Portuguese). Ph.D. Thesis, Polytechnic School, University of Sao Paulo, Brazil, 1984. 47. Matthews, C.H., Carbohydrate losses at high temperature in kraft pulping. Sven. Papperstidn. 1974, 77(17), 629-635. 48. Kleinert, T.N., Mechanisms of alkaline delignification. I. The overall reaction patern. Tappi 1966, 49 (2), 53-57 49. Wilder, H.D., Daleski, E.J., Kraft Pulping Kinetics. I. Literature Review and Research Program. Tappi 1964, 47 (5), 270-275. 50. Wilder, H.D., Daleski, E. J. Delignification rate studies. Part II of a series on kraft pulping kinetics. Tappi 1965, 48 (5), 293-297. 51. Chiang, V.L., Puumala, R.J., Takeuchi, H., Eckert, R. E., Comparison of softwood and hardwood kraft pulping. Tappi J.1988, 71 (9), 173-176. 52. Chorin, A.J., A numerical method for solving incompressible viscous flow problem, J. Comput. Phys., (2) 12, 1967. 53. Patankar Suhas V., Numerical Heat Transfer and fluid flow, Hemisphere Pub. Corp.; New York, McGraw-Hill, P37, 1980.  141  54. Harlow, F.H., Welsh, J.E., Numerical Calculation of Time Dependent Viscous Incompressible Flow with Free Surface. Phys. Fluids, (8) 2181-2189,1965. 55. Ferziger, J.H. and Peric, M., Computational methods for fluid dynamics, second Edition, Springer, 1999. ISBN3-540-65373-2. QA911 .F434 1999 c l 56. Peric, M . , et al. Comparison of finite-volume numerical methods with staggered and collocated grids, Computers and Fluids, 16(4):389-403, 1998 57. Stone, H.L., Iterative solution of implicit approximations of multi-dimensional partial differential equations, SIAM Journal of Numerical Analysis, (5) 530-558, 1968. 58. Hageman, L.A. and Young, D.M., Applied iterative methods, Wiley, New York, 1981. 59. Muzaferija, S., Adaptive finite volume method for flow predictions using unstructured meshes and multigrid approach. PhD Thesis, Univeristy of London, 1994. 60. Chen, C.J., Bernatz, R., et al, Finite Analytic Method in Flows and Heat Transfer, ISBN 1-56032-898-3, Taylor and Francis, 2000. 61. Hortmann, M . , Peric, M., Scheurer, G., Finite Volume Multigrid Prediction of laminar natural convection: bench-mark solutions. International Journal of Numerical Methods in Fluids, No.l 1, 189-207, 1990. 62. Skelland, A.H.P, Non-Newtonian Flow and Heat Transfer, John Wiley & Sons, Inc. New York, 1967. 63. Lauriat G. and Vafai K., Forced convective flow and heat transfer through a porous medium exposed to a flat plate or a channel, Convective and Mass Transfer in Porous Media, p289-327, 1991. 64. Vlaev, D.S. and Bennington, C , Using Electrical Resistance Tomography to Image Liquor Flow in a Model Digester, to be published. 142  65. Vlaev, D.S. and Bennington, C , Flow Uniformity in a Model Digester Measured with Electrical Resistance Tomography, 3 World Congress on Industrial Process rd  Tomography, Banff, Canada. 66. Tapp, H.S. and Wilson R.H., Developments in low-cost electrical imaging techniques, Process Control Quality 9 (1-3) (1997) 7-16. 67. Dickin, F., Wang, M., Electrical resistance tomography for process applications, Meas. Sci. Technol. 7 (3) (1996) 247-260. 68. Ma Y., Zheng Z. et al., Application of electrical resistance tomography system to monitor gas/liquid two-phase flow in a horizontal pipe, Flow Measurement and Instrumentation 12, (2001) 259-265 69. Peyton A. J., Yu, Z. Z., Lyon G., AlZeibak S., et al., An overview of electromagnetic inductance tomography: description of three different systems, Meas. Sci. Technol. 7 (3) (1996) 261-271. 70. Reinecke, N . and Mewes, D., Recent developments and industrial/research applications of capacitance tomography, Meas. Sci. Technol. 7 (3) (1996) 233246. 71. Kenneth E. Smith, Process control for pulp and paper mills, Miller Freeman Publications, ISBN: 087930149X, 1983. 72. Graham Vandegriend, pure species LWAFL and coarseness values, ECOLO Talk, No.4, 1999. 73. Wood Handbook: Wood as an Engineering Material, Agriculture Handbook 72, U.S. Department of Agriculture, Forest Service, Forest Products Laboratory, Madison, Wise. Revised, 1987. 74. Robert H. White and Erik V. Nordheim, Charring Rate of Wood for ASTM E 119 Exposure, Fire Techology,Vol.28, No.l, 1992. 75. R. Pettersen, The Chemical Composition of Wood, in R. Rowell, ed. The Chemistry of Solid Wood. American Chemical Society, Washington, D. C , 1984.  143  76. Chen H., Harmon, M. E., Sexton J., and Fasth, B., Fine-root decomposition and N dynamics in coniferous forests of the Pacific Northwest, U.S.A., Can. J. For. Res. 32: 320-331 (2002). 77. Blum L., The Production of Bleached Kraft Pulp, Environmental Defence Fund, 1996.  144  Appendix A Mill record of digester operation at HSPP* B  A Tagname  411AI9129  Engineering Units  KAPPA  C g/L  BC Dig Blowline RESIDUAL ALKALI Kappa  Description  D  E  F  G  H  I  J  K  411AYI9125 411AYI9125 411AYI9125 411AYI9125 411AYI9125 411AYI9125 411CT9036 411DRPM B D E F 411AT9606 A C  H-Fctr  H-Fctr  H-Fctr  H-Fctr  H-Fctr  H-Factor to Wash  H-Factor to Mid EMCC  H-Factor to End MCC  H-Factor to Mid MCC  H-Factor to Extraction  H-Fctr  G/L  H-Factor in WL IV STRENGTH  L 411FI9121  RPM  US  CHIP METER RPM  DIG BLOW FLOW(F12A/B)  9-Mar-2004 15:00  30.1  12.8  1494.0  1255.4  852.8  645.7  355.8  0.0  84.2  15.7  133.8  9-Mar-2004 16:00  29.9  12.8  1513.2  1228.1  826.5  656.0  382.8  0.0  84.4  15.3  134.0  9-Mar-2004 17:00  29.1  12.8  1467.4  1210.9  862.6  703.3  409.7  0.0  84.5  15.4  137.7  9-Mar-2004 18:00  29.7  12.8  1500.6  1264.7  914.8  713.2  375.5  0.0  84.6  13.4  130.9  9-Mar-2004 19:00  29.4  12.8  1565.6  1317.1  883.3  713.3  420.8  0.0  84.8  15.4  129.3  9-Mar-2004 20:00  29.6  12.8  1568.7  1268.1  905.2  662.9  336.6  0.0  84.9  15.2  125.6  Min  29.1  12.8  1467.4  1210.9  826.5  645.7  336.6  0.0  84.2  13.4  125.6  Avg  29.6  12.8  1518.2  1257.4  874.2  682.4  380.2  0.0  84.6  15.1  131.9  Max  30.1  12.8  1568.7  1317.1  914.8  713.3  420.8  0.0  84.9  15.7  137.7  23%  3%  1%  7%  8%  10%  10%  22%  1%  15%  9%  11-Mar-2004 14:00  28.4  12.8  1423.9  1173.2  794.1  625.7  349.9  0.0  88.6  14.8  130.8  11-Mar-2004 15:00  27.9  13.1  1424.3  1168.6  791.5  629.8  362.1  0.0  87.1  14.7  136.5  11-Mar-2004 16:00  28.1  12.3  1405.1  1150.2  805.1  647.9  375.6  0.0  87.1  15.4  120.8  11-Mar-2004 17:00  27.1  12.4  1387.3  1154.0  832.5  645.6  353.7  0.0  87.1  15.3  128.5  11-Mar-2004 18:00  27.7  12.5  1410.2  1181.9  800.9  633.7  373.1  0.0  87.0  14.8  141.7  11-Mar-2004 19:00  28.7  12.3  1421.3  1162.0  823.7  659.4  388.5  0.0  87.0  14.9  142.7  Min  27.1  12.3  1387.3  1150.2  791.5  625.7  349.9  0.0  87.0  14.7  120.8  Avg  28.0  12.6  1412.0  1165.0  808.0  640.3  367.2  0.0  87.3  15.0  133.5  Max  28.7  13.1  1424.3  1181.9  832.5  659.4  388.5  0.0  88.6  15.4  142.7  2%  4%  16%  (max-min)/avg  (max-min)/avg  18%  5%  6%  3%  3%  *HSPP is an abbreviation of "Howe Sound Pulp and Paper". 145  5%  5%  11%  Appendix A Mill record of digester operation at HSPP (continued) A  M  Tagname  N  O  P  Q  R  S  T  u  V  W  411FT9011 411FT9021 411FT9031 411FT9032 411FT9037 411FT9038 411FT9039 411FT9041 411FT9111 411FT9141 411FT9160  Engineering Units  T/H  Description  HI PRESS STEAM FLOW  9-Mar-2004 15:00  40.0  27.4  3.9  5.0  2.5  4.9  48.3  0.4  50.0  11.6  151.0  9-Mar-2004 16:00  36.3  28.1  3.9  4.9  1.3  4.9  47.1  0.4  46.4  10.4  136.6  9-Mar-2004 17:00  31.7  27.1  3.9  4.7  1.5  4.7  47.2  0.4  47.5  10.8  133.2  9-Mar-2004 18:00  32.1  24.0  3.9  4.5  1.7  4.5  40.5  0.4  47.6  11.0  125.2  9-Mar-2004 19:00  36.8  29.7  3.9  4.5  2.1  4.5  47.1  0.4  49.0  11.1  148.5  9-Mar-2004 20:00  33.3  27.8  4.0  4.5  1.4  4.5  45.8  0.4  46.3  10.6  142.3  Min  31.7  24.0  3.9  4.5  1.3  4.5  40.5  0.4  46.3  10.4  125.2  Avg  35.0  27.3  3.9  4.7  1.7  4.7  46.0  0.4  47.8  10.9  139.5  Max  40.0  29.7  4.0  5.0  2.5  4.9  48.3  0.4  50.0  11.6  151.0  23%  23%  21%  1%  12%  10%  17%  9%  8%  11%  18%  11-Mar-2004 14:00  38.7  26.1  3.5  4.7  4.4  4.7  43.7  0.3  50.8  5.7  146.1  11-Mar-2004 15:00  35.7  25.9  3.6  4.8  0.0  4.7  44.3  0.3  52.8  5.8  137.4  11-Mar-2004 16:00  32.9  26.9  3.5  4.0  0.0  4.0  46.6  0.3  47.5  5.6  144.2  11-Mar-2004 17:00  39.2  27.5  3.5  4.0  0.0  4.0  46.3  0.3  50.9  6.0  145.2  11-Mar-2004 18:00  38.6  25.7  3.5  4.5  0.0  4.5  44.7  0.3  54.6  5.8  136.9  11-Mar-2004 19:00  38.7  26.4  3.4  4.8  0.0  4.8  45.0  0.3  55.4  6.4  140.3  Min  32.9  25.7  3.4  4.0  0.0  4.0  43.7  0.3  47.5  5.6  136.9  Avg  37.3  26.4  3.5  4.4  0.8  4.5  45.1  0.3  52.0  5.9  141.7  Max  39.2  27.5  3.6  4.8  4.4  4.8  46.6  0.3  55.4  6.4  146.1  17%  7%  4%  18%  18%  6%  0%  15%  13%  6%  (max-min)/avg  (max-min)/avg  146  18%  T/H  L/S  L/S  F 3EWL LP STEAM F3A HPF FLOW TO WASH FLOW WL PURGE  L/S F3B WL FLOW TO BC  L/S  L/S  L/S  L/S  L/S  L/S  F 3D WL F16AEXTR FLOW TO F3F WL F4 CB F11 CB FOR F 1 4 C B T O SCRN MC FLOW TO IV FLOW TO IV DIG PRES FLOW DIG OD  Appendix A Mill record of digester operation at HSPP (continued) A  X  Tagname  Y  Z  AA  AB  AC  AD  AE  AF  AG  AH  411FT9181 411FYI9025 411FYI9121 411HI9022B 411HT9022B 411LY9102 411LY9752 411MI9022 411NIC9183 411TT9085 411TT9165  Engineering Units  L/S  T/H  ADT/D  FURNSH  N/A  %  CALC STEAM DIGESTER DIGESTER EXIT PRODUCTIO CHIP F18 C B T O FLOW TO DIGESTER BLOW RATE FURNISH N GRADE LEVEL DIG BOT BIN  Description  9-Mar-2004 15:00 9-Mar-2004 16:00 9-Mar-2004 17:00 9-Mar-2004 18:00 9-Mar-2004 19:00 9-Mar-2004 20:00  41.8 38.3 39.0 39.2 40.7 38.2  29.2 30.1 28.6 26.1 30.6 28.5  1286.4 1289.6 1325.2 1259.1 1238.5 1206.6  400.0 400.0 400.0 400.0 400.0 400.0  400.0 400.0 400.0 400.0 400.0 400.0  29.3 29.8 33.2 24.3 27.8 24.0  %  %  % BD  DEG C  DEG C  TRIM LQR IV CHIP WEIGHTED PD18 BLOWTEMPERAT EXTRACTI LEVEL AVE MOIST CONS MV URE ON TEMP  8.8 31.1 27.3 16.5 3.9 19.9  54.0 54.4 54.2 54.0 54.4 54.5  10.0 10.0 10.0 10.0 10.0 10.0  157.9 158.8 158.6 158.3 157.2 156.3  160.3 160.6 161.6 161.2 161.3 159.7  Min  38.2  26.1  1206.6  400.0  400.0  24.0  3.9  54.0  10.0  156.3  159.7  Avg  39.5  28.8  1267.6  400.0  400.0  28.1  17.9  54.3  10.0  157.9  160.8  41.8  30.6  1325.2  400.0  400.0  33.2  31.1  54.5  10.0  158.8  161.6  9%  16%  9%  0%  0%  1%  1%  2%  1%  11-Mar-2004 14:00 11-Mar-2004 15:00 11-Mar-2004 16:00 11-Mar-2004 17:00 11-Mar-2004 18:00 11-Mar-2004 19:00  40.0 41.7 37.1 39.4 42.8 43.5  27.6 27.7 28.3 29.1 27.6 28.3  1299.8 1346.6 1189.8 1264.8 1378.4 1397.7  400.0 400.0 400.0 400.0 400.0 400.0  400.0 400.0 400.0 400.0 400.0 400.0  40.0 39.6 24.8 33.0 42.2 38.3  33.7 35.0 14.0 16.7 33.4 27.4  54.5 54.5 54.5 54.5 54.5 54.5  10.3 10.3 10.2 10.2 10.1 10.2  156.5 157.4 157.5 157.4 157.4 157.7  159.8 159.6 160.2 159.8 159.8 160.2  Min  37.1  27.6  1189.8  400.0  400.0  24.8  14.0  54.5  10.1  156.5  159.6  Avg  40.8  28.1  1312.9  400.0  400.0  36.3  26.7  54.5  10.2  157.3  159.9  Max  43.5  29.1  1397.7  400.0  400.0  42.2  35.0  54.5  10.3  157.7  160.2  16%  5%  16%  0%  0%  0%  2%  1%  0%  Max  (max-min)/avg  (max-min)/avg  147  23%  18%  Appendix A Mill record of digester operation at HSPP (continued) Al  A  AJ  AK  AL  AM  AN  AO  411TT9195 411TT9197 411TT9205 411TT9207 411TT9600 411TT9605 411TT9607  Tagname  Engineering Units  DEG C  DEG C  T19C MC COLD T19H MC TEMP HOT TEMP  Description  DEG C T20C WASH COLD TEMP  DEG C  DEG C  T20H T60C BC WASH COLD HOT TEMP TEMP  DEG C  DEG C  AP  AQ  411TT9700  411TT9705  DEG C  DEG C  T60H BC T60L BC H IV UPPER IV LOWER HOT TEMP TEMP SP TEMPERATURE TEMPERATURE  AR  AS  411TT9710 411TT9711  DEG C  DEG C  TRIM ZONE TEMP  EXTR ZONE 1 TEMP  161.9 161.8 161.7 161.4 161.8 161.4  161.0 161.0 161.2 160.1 161.2 160.7  150.3 153.6 157.4 152.1 149.6 156.2  159.1 159.1 159.2 159.4 158.9 158.9  158.4 159.2 158.0 159.4 156.7 158.5  177.4 177.6 175.4 174.5 174.2 176.0  158.3 159.0 157.4 158.8 156.4 158.3  122.7 125.7 125.3 124.1 123.1 123.3  129.6 129.7 130.0 130.2 130.1 129.8  157.2 157.8 157.9 157.7 158.1 157.3  160.8 160.5 160.6 160.1 160.0 160.3  Min  161.4  160.1  149.6  158.9  156.7  174.2  156.4  122.7  129.6  157.2  160.0  Avg  161.6  160.9  153.2  159.1  158.4  175.9  158.0  124.0  129.9  157.6  160.4  Max  161.9  161.2  157.4  159.4  159.4  177.6  159.0  125.7  130.2  158.1  160.8  9-Mar-2004 9-Mar-2004 9-Mar-2004 9-Mar-2004 9-Mar-2004 9-Mar-2004  15:00 16:00 17:00 18:00 19:00 20:00  (max-min)/avg  23%  0%  1%  5%  0%  2%  2%  2%  2%  0%  1%  0%  11-Mar-2004 11-Mar-2004 11-Mar-2004 11-Mar-2004 11-Mar-2004 11-Mar-2004  14:00 15:00 16:00 17:00 18:00 19:00  161.1 160.7 160.7 160.5 160.6 160.5  160.4 159.9 160.1 160.0 160.0 159.9  147.0 148.6 150.2 142.8 144.0 143.8  158.6 158.2 157.4 157.8 158.1 158.2  157.0 158.0 157.7 157.8 158.2 158.3  176.2 176.1 175.3 176.3 176.7 176.7  157.0 157.9 157.2 157.2 157.8 157.8  123.6 124.6 123.0 123.1 123.1 124.1  129.2 129.8 129.8 129.1 128.7 128.8  156.4 157.0 157.3 157.0 157.1 157.5  160.1 159.6 159.1 159.4 159.7 159.4  Min  160.5  159.9  142.8  157.4  157.0  175.3  157.0  123.0  128.7  156.4  159.1  Avg  160.7  160.1  146.1  158.1  157.9  176.2  157.5  123.6  129.2  157.0  159.5  Max  161.1  160.4  150.2  158.6  158.3  176.7  157.9  124.6  129.8  157.5  160.1  0%  0%  5%  1%  1%  1%  1%  1%  1%  1%  1%  (max-min)/avg  148  18%  Appendix A Mill record of digester operation at HSPP (continued) A  AT  Tagname Engineering Units  AU  AV  AW  AX  AY  AZ  BA  BB  BC  BD  BE  411X19022 411X19022 D C 411FT9601 411FT9611 411TT9712 411TT9713 411TT9714 411TT9715 411UI9160 411WI9022411XI9022A411XI9022B  DEG C  DEG C  Description  EXTR ZONE 2 TEMP  MCC ZONE EXTR ZONE 3 TEMPERAT TEMP URE  DEG C  9-Mar-2004 15:00  161.0  9-Mar-2004 16:00  DEG C WASH ZONE TEMP  %  T/DBD  %  %  %  %  HEMLOCK CEDAR SPF FIR TOTAL CALC DILU CHIPS TO FURNISH FURNISH FURNISH FURNISH CB FACTOR FRACTION FRACTION FRACTION FRACTION  L/S  L/S  F60 BC FLOW  F61A UPPER SLUICE FL  161.9  160.8  0.0  2523.9  0.0  42.6  57.3  0.0  156.3  115.0  160.6  161.1 160.9  161.7  162.0  -0.7  2659.5  0.0  47.3  52.7  0.0  152.8  125.3  9-Mar-2004 17:00  160.8  161.0  162.9  161.0  -0.8  2348.0  0.0  46.4  53.7  0.0  150.6  126.1  9-Mar-2004 18:00  160.7  160.7  161.6  160.0  -0.5  2253.2  0.0  30.5  68.7  0.0  148.2  125.8  9-Mar-2004 19:00  160.5  160.6  161.3  160.7  0.1  2481.0  0.0  47.2  52.7  0.0  147.4  117.8  9-Mar-2004 20:00  160.7  160.6  162.5  160.6  -0.2  2302.7  0.0  46.5  55.1  0.0  147.2  119.4  Min  160.5  160.6  161.3  160.0  -0.8  2253.2  0.0  30.5  52.7  0.0  147.2  115.0  Avg  160.7  160.8  162.0  160.9  -0.4  2428.1  0.0  43.4  56.7  0.0  150.4  121.6  Max  161.0  161.1  162.9  162.0  0.1  2659.5  0.0  47.3  68.7  0.0  156.3  126.1  0%  0%  1%  1%  6%  9%  11-Mar-2004 14:00  160.3  160.3  160.3  160.4  0.1  2304.8  0.0  45.9  53.1  0.0  138.4  121.1  11-Mar-2004 15:00  159.9  159.8  160.1  161.5  0.0  2338.5  0.0  46.5  53.3  0.0  136.5  122.3  11-Mar-2004 16:00  159.2  159.5  160.0  160.8  0.1  2608.2  0.0  47.1  52.9  0.0  135.6  126.0  11-Mar-2004 17:00  159.4  159.6  160.0  159.9  0.1  2558.6  0.0  47.0  53.9  0.0  135.3  120.0  11-Mar-2004 18:00  159.7  159.9  160.5  154.2  0.0  2175.3  0.0  47.0  53.4  0.0  136.4  120.2  11-Mar-2004 19:00  159.7  159.4  159.9  156.0  0.1  2213.6  0.0  47.2  52.8  0.0  136.9  122.3  (max-min)/avg  17%  Min  159.2  159.4  159.9  154.2  0.0  2175.3  0.0  45.9  52.8  0.0  135.3  120.0  Avg  159.7  159.7  160.1  158.8  0.1  2366.5  0.0  46.8  53.2  0.0  136.5  122.0  Max  160.3  160.3  160.5  161.5  0.1  2608.2  0.0  47.2  53.9  0.0  138.4  126.0  1%  1%  0%  5%  2%  5%  (max-min)/avg  149  23%  18%  18%  Appendix A Mill record of digester operation at HSPP (continued) A  BF  Tagname Engineering Units  BG  BH  Bl  BJ  BK  BL  BM  BN  BO  BP  BQ  411FT9191 411FT9201 411FT9614 411TT9125 411TT9245 411NT9126 411IT0024 411PDT9183 411FT9871 411AT9158 411AT9190 421AI9129  DEG C  %  A  CB COOLER TEMP  DIG BLOW CONSIST  OUTLET DEVICE (DIG)  80.7  56.6  10.9  372.2  83.4  0.0  10.3  14.2  32.4  13.7  83.8  55.5  10.9  373.2  81.8  0.1  10.3  14.3  30.9  95.5  15.5  81.4  53.0  10.7  343.0  82.6  0.0  10.5  14.2  31.5  120.1  95.2  16.9  81.9  53.3  10.9  384.0  82.7  0.0  10.5  14.1  32.3  120.1  95.2  18.4  81.2  57.7  10.9  374.5  83.0  0.0  10.2  14.1  31.4  81.2  57.7  10.8  349.1  80.9  0.0  L/S  L/S  Description  F19AMC UPPER FLOW  WASH CIRC FLOW CONT  9-Mar-2004 15:00  120.1  95.7  13.1  9-Mar-2004 16:00  119.9  95.6  9-Mar-2004 17:00  120.2  9-Mar-2004 18:00 9-Mar-2004 19:00  L/S  DEG C  BC Heater Bypass DIF INLET Flow TEMP  KPA  L/S  g/i  g/i  KAPPA  PD18 BLOW BLOW LINE Extraction MC Resid Decker Vat CONSIST DILN Resid Alkali Alkali Kappa  9-Mar-2004 20:00  120.0  94.8  15.5  10.5  14.0  31.3  Min  119.9  94.8  13.1  80.7  53.0  10.7  343.0  80.9  0.0  10.2  14.0  30.9  Avg  120.1  95.3  15.5  81.7  55.6  10.9  366.0  82.4  0.0  10.4  14.2  31.6  Max  120.2  95.7  18.4  83.8  57.7  10.9  384.0  83.4  0.1  10.5  14.3  32.4  0%  1%  4%  8%  2%  11%  3%  3%  2%  5%  11-Mar-2004 14:00 11-Mar-2004 15:00 11-Mar-2004 16:00 11-Mar-2004 17:00 11-Mar-2004 18:00 11-Mar-2004 19:00  120.1 120.2 120.1 119.9 119.9 120.1  89.4 89.2 89.3 89.4 89.1 89.1  46.6 44.1 42.3 44.0 45.5 43.3  82.1 82.2 80.7 83.3 82.2 81.5  58.9 58.1 61.0 58.5 55.6 58.6  11.4 11.1 10.9 11.0 10.8 10.9  410.1 409.4 399.4 410.8 394.6 407.5  93.8 92.3 90.4 90.6 87.1 89.5  0.0 0.0 0.0 0.0 0.6 0.0  11.6 11.4 10.6 10.3 10.3 10.3  13.6 14.6 14.0 13.8 13.7 14.0  28.1 28.0 27.4 27.4 28.7 31.3  Min  119.9  89.1  42.3  80.7  55.6  10.8  394.6  87.1  0.0  10.3  13.6  27.4  Avg  120.1  89.3  44.3  82.0  58.5  11.0  405.3  90.6  0.1  10.7  13.9  28.5  Max  120.2  89.4  46.6  83.3  61.0  11.4  410.8  93.8  0.6  11.6  14.6  31.3  0%  0%  3%  9%  5%  4%  7%  12%  7%  14%  (max-min)'avg  (max-min)/avg  150  23%  18%  Appendix B Average mill records of digester operation at HSPP * In all the cases below, efforts have been made to keep the operations stable for a minimal of 6 hours. A  Al  Tagname Deviation Engineering % Units  B 411AI9129  C  D  E  F  G  H  I  J  411AYI9125 411AYI9125 411AYI9125 411AYI9125 411AYI9125 411AYI9125 E F 411CT9036 D 411AT9606 A B C  KAPPA  g/L  H-Fctr  H-Fctr  H-Fctr  H-Fctr  H-Fctr  H-Fctr  G/L  K  L  411DRPM  411FI9121  RPM  L/S  CHIP METER RPM  DIG BLOW FLOW(F12A  Description  (max-min)/ avg  Dig Blowline Kappa  BC RESIDUAL ALKALI  H-Factor to Wash  H-Factor to Mid EMCC  H-Factor to End MCC  H-Factor to Mid MCC  H-Factor to Extraction  09-Mar-04  23  29.6  12.8  1518.2  1257.4  874.2  682.4  380.2  0.0  84.6  15.1  131.9  11-Mar-04  18  28.0  12.6  1412.0  1165.0  808.0  640.3  367.2  0.0  87.3  15.0  133.5  26-Mar-04  27  28.4  13.5  1517.6  1338.9  975.4  767.1  413.0  0.0  85.9  10.4  88.2  12-Apr-04  16  28.0  14.6  1438.6  1179.5  809.4  633.5  349.8  0.0  84.0  15.7  135.0  06-May-04  16  31.0  13.5  1344.1  1123.9  799.3  640.1  369.1  0.0  89.7  15.4  125.7  02-Jun-04  225  29.3  14.7  1397.5  1130.8  747.9  573.0  302.8  0.0  84.0  15.0  119.0  05-Jun-04  63  29.2  13.8  1293.1  1106.8  807.0  653.6  377.7  0.0  86.7  14.5  124.3  15-Jun-04  24  28.3  13.5  1437.9  1185.8  824.4  651.6  358.5  0.0  87.0  14.9  139.6  23-Jul-04  19  35.7  12.9  1526.9  1251.5  853.5  665.1  367.5  0.0  90.4  13.8  112.9  07-Aug-04  24  40.9  12.8  1596.3  1291.4  868.7  675.4  369.1  0.0  87.6  12.9  108.9  08-Aug-04  41  29.4  13.6  1632.9  1336.6  912.1  712.0  402.4  0.0  90.0  11.1  92.9  09-Aug-04  36  32.0  13.5  1636.6  1317.4  869.0  672.8  362.8  0.0  90.7  12.0  104.8  11-Aug-04  24  32.3  13.5  1444.4  1171.3  783.6  606.3  329.2  0.0  92.1  14.8  140.0  08-Sep-04  28  32.3  13.3  1217.0  1023.4  734.4  595.1  353.8  0.0  89.6  13.3  108.6  13-Sep-04  19  35.1  13.5  1198.6  1012.3  729.2  592.2  349.1  0.0  85.7  15.4  121.5  151  H-Factor in WL IV STRENGTH  IB)  Appendix B Average mill records of digester operation at HSPP (continued) A  M  N  O  P  Q  R  S  T  U  V  W  Tagname  411FT9011  411FT9021  411FT9031  411FT9032  411FT9037  411FT9038  411FT9039  411FT9041  411FT9111  411FT9141  411FT9160  Engineering Units  T/H  T/H  L/S  L/S  L/S  L/S  L/S  L/S  L/S  L/S  L/S  Description  HI PRESS STEAM FLOW  09-Mar-04  35.0  27.3  3.9  4.7  1.7  4.7  46.0  0.4  47.8  10.9  139.5  11-Mar-04  37.3  26.4  3.5  4.4  0.8  4.5  45.1  0.3  52.0  5.9  141.7  26-Mar-04  24.6  18.0  4.0  3.3  1.9  3.3  30.9  0.3  31.9  3.7  109.2  12-Apr-04  38.1  22.7  4.0  4.4  0.2  4.4  48.3  5.0  52.2  11.4  147.4  06-May-04  32.3  20.7  3.9  4.7  0.8  4.6  43.3  0.4  47.3  15.5  133.6  02-Jun-04  27.3  26.2  4.1  3.9  0.0  3.9  48.7  1.2  30.8  7.3  109.4  05-Jun-04  28.2  33.1  4.0  3.8  2.0  3.7  42.3  4.0  36.8  8.2  110.5  15-Jun-04  29.0  21.4  3.7  4.0  0.5  4.0  44.0  0.5  44.6  10.0  115.9  23-Jul-04  29.4  15.9  3.6  3.9  0.8  3.5  41.2  8.0  45.7  10.7  120.2  07-Aug-04  27.4  16.0  4.1  4.6  1.5  4.4  38.7  6.0  42.2  10.3  115.1  08-Aug-04  26.2  13.1  4.3  3.9  2.1  3.9  32.4  6.0  42.1  10.2  116.8  09-Aug-04  26.1  14.4  4.3  4.3  2.9  4.3  34.9  6.0  40.3  9.7  107.7  11-Aug-04  29.1  17.4  4.3  3.8  0.4  4.2  43.4  8.0  47.8  7.0  120.6  08-Sep-04  25.7  20.4  4.2  3.4  0.4  3.5  38.3  9.0  39.6  11.3  110.4  13-Sep-04  30.6  22.8  3.9  4.5  0.9  4.5  44.6  10.0  40.6  9.7  120.2  152  LP STEAM F3A HPF WL PURGE FLOW  F 3D WL F 3EWL F3BWL FLOW TO FLOW TO WASH FLOW TO BC MC  F4 CB FLOW F11 CB FOR F 1 4 C B T O F16AEXTR F3FWL SCRN FLOW DIG OD TO IV DIG PRES FLOW TO IV  Appendix B Average mill records of digester operation at HSPP (continued) AD  AE  411FYI9025 411FYI9121 411HI9022B 411HT9022B 411LY9102  411LY9752  411MI9022  FURNSH  %  %  AA  AB  AC  AF  AG  AH  X  Tagname  411FT9181  Engineering Units  L/S  T/H  ADT/D  F18 C B T O DIG BOT  CALC STEAM FLOW TO BIN  DIGESTER BLOW RATE  09-Mar-04  39.5  28.8  1267.6  400.0  400.0  28.1  17.9  54.3  10.0  157.9  160.8  11-Mar-04  40.8  28.1  1312.9  400.0  400.0  36.3  26.7  54.5  10.2  157.3  159.9  26-Mar-04  29.7  20.9  888.7  400.0  400.0  38.7  29.9  54.5  10.5  153.7  158.3  12-Apr-04  48.3  24.4  1316.6  400.0  400.0  29.7  15.0  52.5  10.2  157.8  160.6  06-May-04  43.2  22.8  1236.6  400.0  400.0  30.4  15.5  50.5  10.2  158.1  160.7  02-Jun-04  29.0  23.9  1223.1  400.0  400.0  28.4  0.7  49.4  10.7  155.6  158.5  05-Jun-04  30.2  22.1  1265.1  400.0  400.0  34.0  10.6  49.5  10.6  157.2  160.7  15-Jun-04  34.6  23.0  1409.7  400.0  400.0  33.0  33.0  49.5  10.5  157.4  160.4  23-Jul-04  37.6  18.9  1127.8  400.0  400.0  29.9  16.9  46.3  10.4  156.6  159.5  07-Aug-04  37.9  18.6  1100.2  400.0  400.0  28.1  16.4  47.2  10.5  155.6  158.7  08-Aug-04  37.7  15.4  937.1  400.0  400.0  32.1  22.8  46.8  10.5  154.8  158.3  09-Aug-04  35.8  16.6  1057.6  400.0  400.0  32.5  27.9  47.1  10.5  154.4  157.5  11-Aug-04  35.6  19.2  1397.7  400.0  400.0  28.3  19.8  47.1  10.4  156.0  159.4  08-Sep-04  32.1  20.1  1083.9  400.0  400.0  30.2  15.7  46.6  10.4  156.4  157.9  13-Sep-04  35.7  24.1  1189.6  400.0  400.0  30.4  12.8  46.5  10.2  158.2  159.7  Description  153  Y  Z  A  N/A  %  DIGESTER PRODUCTIO DIGESTER EXIT FURNISH N GRADE CHIP LEVEL  IV CHIP LEVEL  411NIC9183 411TT9085  % BD  DEG C  411TT9165  DEG C  TRIM LQR WEIGHTED PD18 BLOW TEMPERATU EXTRACTIO RE N TEMP AVE MOIST CONS MV  Appendix B Average mill records of digester operation at HSPP (continued) AR  A  Al  AJ  AK  AL  AM  AN  AO  AP  AQ  Tagname  411TT9195  411TT9197  411TT9205  411TT9207  411TT9600  411TT9605  411TT9607  411TT9700  411TT9705  Engineering Units  DEG C  DEG C  DEG C  DEG C  DEG C  DEG C  DEG C  DEG C  DEG C  DEG C  DEG C  T60L BC H TEMP SP  IV UPPER TEMP  IV LOWER TEMP  TRIM ZONE TEMP  EXTR ZONE 1 TEMP  Description  T60H BC T19H MC T20C WASH T20H WASH T60C BC T19C MC COLD TEMP HOT TEMP COLD TEMP HOT TEMP COLD TEMP HOT TEMP  AS  411TT9710 411TT9711  09-Mar-04  161.6  160.9  153.2  159.1  158.4  175.9  158.0  124.0  129.9  157.6  160.4  11-Mar-04  160.7  160.1  146.1  158.1  157.9  176.2  157.5  123.6  129.2  157.0  159.5  26-Mar-04  159.9  158.1  146.7  154.0  154.2  166.4  153.6  116.8  127.2  153.7  159.0  12-Apr-04  161.4  160.9  143.8  159.1  157.8  170.4  157.3  122.5  129.1  157.2  160.2  06-May-04  159.8  159.4  145.9  157.1  158.2  170.0  157.5  126.2  131.1  157.9  158.6  02-Jun-04  161.9  160.6  157.1  159.0  155.8  174.5 .  155.4  122.3  128.0  155.0  161.7  05-Jun-04  159.5  157.6  149.9  153.9  157.6  171.7  157.0  122.6  128.5  157.4  159.1  15-Jun-04  161.3  159.9  147.4  157.5  157.4  169.0  156.9  124.8  131.6  157.3  161.0  23-Jul-04  160.5  160.5  145.9  159.1  156.6  168.0  156.1  125.8  133.3  156.4  159.2  07-Aug-04  160.4  160.1  144.9  159.0  155.9  164.5  155.4  125.5  133.0  155.3  159.3  08-Aug-04  159.1  158.6  137.8  157.0  154.4  159.8  153.5  126.0  133.5  154.4  157.1  09-Aug-04  159.6  159.5  144.7  158.2  154.9  161.9  154.4  126.8  133.6  153.8  158.6  11-Aug-04  161.3  160.6  145.9  159.0  156.3  167.7  155.7  126.8  132.5  155.8  160.1  08-Sep-04  155.7  156.0  145.9  153.5  156.5  165.5  155.6  125.3  130.5  155.4  154.2  13-Sep-04  158.2  157.9  148.5  155.3  158.0  170.7  157.1  123.9  128.4  157.1  156.9  154  Appendix B Average mill records of digester operation at HSPP (continued) A Tagname  AT  AU  AV  AW  AX  AY  AZ  BA  BB  BC  BD  BE  411TT9712 411TT9713 411TT9714 411TT9715 411U19160 411WI9022 411XI9022A411XI9022B411XI9022C411XI9022D 411FT9601 411FT9611  L/S  L/S  F60 BC FLOW  F61A UPPER SLUICE FL  0.0  150.4  121.6  53.2  0.0  136.5  122.0  47.2  52.8  0.0  136.3  89.3  0.0  47.1  52.8  0.0  192.2  140.1  2409.0  0.0  46.9  52.9  0.0  164.9  155.0  -1.6  2691.1  0.2  45.7  55.1  0.0  130.0  121.9  157.4  -1.3  2527.6  0.0  46.3  53.1  0.2  164.8  127.6  160.6  160.7  -1.0  2717.8  0.0  47.3  53.1  0.0  167.3  140.9  159.3  159.8  159.7  0.6  2053.4  0.0  43.6  56.2  0.0  149.4  140.1  159.6  159.3  159.6  159.4  0.7  1822.6  0.0  49.2  51.1  0.0  161.2  136.4  08-Aug-04  157.5  157.4  157.2  153.5  2.3  1629.0  0.0  47.6  52.9  0.0  165.3  113.1  09-Aug-04  158.7  158.6  158.7  158.2  0.6  1736.0  0.0  48.0  51.5  0.0  166.7  127.3  11-Aug-04  160.5  160.1  160.1  159.9  0.4  2086.2  0.0  47.4  52.6  0.0  144.6  151.3  08-Sep-04  154.4  154.5  154.8  153.0  -0.2  2064.6  0.0  44.9  55.2  0.0  193.2  122.6  13-Sep-04  157.2  157.1  156.7  156.0  -0.8  2361.3  0.0  45.0  54.5  0.0  184.9  137.4  T/DBD  Engineering Units  DEG C  DEG C  DEG C  DEG C  Description  EXTR ZONE 2 TEMP  EXTR ZONE 3 TEMP  MCC ZONE TEMPERAT URE  WASH ZONE TEMP  09-Mar-04  160.7  160.8  162.0  160.9  -0.4  2428.1  0.0  43.4  56.7  11-Mar-04  159.7  159.7  160.1  158.8  0.1  2366.5  0.0  46.8  26-Mar-04  159.5  159.5  159.4  155.8  0.0  1909.0  0.0  12-Apr-04  160.7  160.0  160.0  158.8  0.5  2517.0  06-May-04  159.0  158.6  157.8  159.1  0.7  02-Jun-04  161.8  161.8  161.6  159.7  05-Jun-04  159.3  159.3  158.1  15-Jun-04  160.5  161.3  23-Jul-04  159.4  07-Aug-04  155  %  %  %  %  %  FIR SPF TOTAL HEMLOCK CEDAR CALC DILU CHIPS TO FURNISH FURNISH FURNISH FURNISH CB FRACTION FRACTION FRACTION FRACTION FACTOR  Appendix B Average mill records of digester operation at HSPP (continued) A Tagname  BF  BG  BH  BI  BJ  BK  BL  411FT9191 411FT9201 411FT9614 411TT9125 411TT9245 411NT9126 411IT0024  DEG C  DEG C  A  BN  BO  BP  BQ  411PDT918 411FT9871 411AT9158 411AT9190 421AI9129 3  KPA  L/S  KAPPA  Engineering Units  L/S  L/S  Description  F19AMC UPPER FLOW  WASH CIRC FLOW CONT  09-Mar-04  120.1  95.3  15.5  81.7  55.6  10.9  366.0  82.4  0.0  10.4  14.2  31.6  11-Mar-04  120.1  89.3  44.3  82.0  58.5  11.0  405.3  90.6  0.1  10.7  13.9  28.5  26-Mar-04  119.8  89.7  47.5  82.3  56.7  11.2  408.7  99.0  0.0  10.4  13.3  30.3  12-Apr-04  120.0  96.7  0.1  82.0  61.0  11.0  378.4  87.6  0.0  11.5  13.8  30.1  06-May-04  120.0  95.3  0.0  81.9  63.0  10.8  402.8  90.7  0.0  10.0  15.2  30.2  02-Jun-04  110.1  92.5  41.5  82.0  40.7  10.9  409.4  105.4  9.0  11.3  20.2  29.1  05-Jun-04  110.1  86.3  -0.1  79.7  44.0  10.7  409.4  102.3  8.9  10.3  20.3  29.7  15-Jun-04  110.1  86.7  0.0  74.4  43.0  10.9  409.1  99.0  8.6  10.1  12.7  28.6  23-Jul-04  119.9  88.2  20.9  84.2  64.3  11.3  403.5  95.2  0.0  9.4  12.6  34.1  07-Aug-04  120.2  84.6  0.2  86.3  64.2  11.2  408.8  98.4  0.0  8.6  12.5  39.9  08-Aug-04  120.2  83.3  0.2  85.1  69.7  11.3  408.5  99.0  0.2  9.1  12.6  28.0  09-Aug-04  120.2  83.1  0.0  86.2  63.9  11.3  407.8  98.3  0.3  9.0  12.8  30.9  11-Aug-04  120.1  90.1  21.0  86.0  65.3  10.1  407.4  95.7  10.3  9.7  13.1  31.5  08-Sep-04  120.0  89.4  -0.6  80.6  57.7  11.2  407.5  95.1  0.2  8.7  12.8  32.9  13-Sep-04  120.1  86.2  3.0  82.6  52.6  10.9  392.8  89.5  0.0  8.5  13.3  35.6  156  L/s  BM  %  OUTLET BC Heater CB Bypass DIF INLET COOLER DIG BLOW DEVICE CONSIST (DIG) Flow TEMP TEMP  g/i  g/i  PD18 BLOW BLOW Extraction MC Resid Decker Vat CONSIST LINE DILN Resid Alkali Alkali Kappa  

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0080759/manifest

Comment

Related Items