UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Modeling the creping process in tissue making Pan, Kui 2019

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

Item Metadata

Download

Media
24-ubc_2019_may_pan_kui.pdf [ 3.61MB ]
Metadata
JSON: 24-1.0377758.json
JSON-LD: 24-1.0377758-ld.json
RDF/XML (Pretty): 24-1.0377758-rdf.xml
RDF/JSON: 24-1.0377758-rdf.json
Turtle: 24-1.0377758-turtle.txt
N-Triples: 24-1.0377758-rdf-ntriples.txt
Original Record: 24-1.0377758-source.json
Full Text
24-1.0377758-fulltext.txt
Citation
24-1.0377758.ris

Full Text

MODELING THE CREPING PROCESS IN TISSUE MAKING by  Kui Pan  B.Sc., University of Science and Technology of China, 2011  M.Sc., University of Science and Technology of China, 2014  A DISSERTATION SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY  in THE FACULTY OF GRADUATE AND POSTDOCTORAL STUDIES (Mechanical Engineering)  THE UNIVERSITY OF BRITISH COLUMBIA (Vancouver)  March 2019  © Kui Pan, 2019  ii  The following individuals certify that they have read, and recommend to the Faculty of Graduate and Postdoctoral Studies for acceptance, the dissertation entitled:  MODELING THE CREPING PROCESS IN TISSUE MAKING  submitted by Kui Pan in partial fulfillment of the requirements for the degree of DOCTOR OF PHILOSOPHY in Mechanical Engineering  Examining Committee: A. Srikantha Phani Supervisor Sheldon Green Co-supervisor   Supervisory Committee Member Yusuf Altintas University Examiner Chad Sinclair University Examiner  Additional Supervisory Committee Members: Reza Varizi Supervisory Committee Member Richard Kerekes Supervisory Committee Member iii  Abstract  Creping is a key operation in the manufacturing of low-density tissue paper. In this process, a wet web is pressed and adhered onto a drying cylinder (Yankee) rotating at a high speed, dried on Yankee, and then scraped off by a doctor blade. This controlled and violent interaction between the web moving at a high speed with the stationary blade creates a series of invisible micro-folds, and explodes the thickness of the web through inter-fiber debonding. Various parameters govern creping and finding their optimal combination is currently limited to experience or costly trials. A one-dimensional particle dynamics model is developed to study nonlinear deformations in the creping process, and to understand the underlying mechanisms. Specifically, the web is modeled as a single layer of discrete particles connected by visco-elasto-plastic elements. A mixed-mode discrete cohesive zone model is embedded to represent the adhesive layer. Self-contact of the web is incorporated by a penalty method. First, a systematic parametric study is reported to assess the relative impact of various process parameters on the crepe structure and hence the tissue quality. Then, the model is extended to a multi-layered web to investigate the “sheet explosion”. A phase diagram for the creping regimes is constructed. Next, the effects of inhomogeneities on the creping process are investigated. Three common inhomogeneities are considered separately: the forming fabric pattern; the non-uniform basis weight; and the non-uniform adhesion. Finally, a series of experiments have been conducted on an existing lab-scale creping apparatus to validate the proposed model, and qualitative agreement is observed. The model can serve as a tool to investigate the process-structure-property correlation in tissue making, and the findings in this thesis offer practical guidance to the industry in the choice of forming and creping process parameters. iv  Lay Summary  This thesis concerns the creping process, a key operation in the manufacturing process of low-density tissue paper. Creping consists of pressing and bonding a wet web onto a drying cylinder rotating at high speed, and scraping it off subsequently, in order to enhance the softness, absorbency and bulk of the tissue products. Many parameters govern creping and finding the optimal combination is presently limited to trial and error. This thesis, advances a mathematical model to understand the creping process. The model is able to capture the important features involved in creping, and evaluate the effects of the key governing parameters. The modeling results are validated by conducting experiments on an existing lab-scale creping apparatus. The model provides a useful tool to study the process-structure-property correlation of tissue paper, and the findings in this thesis offer guidance to the industry.   v  Preface  This PhD thesis entitled “Modeling the Creping Process in Tissue Making” is an original work carried out by the author, Kui Pan, under the supervision of Dr. A. Srikantha Phani and Dr. Sheldon Green. The following contributions have been published or are under consideration for publication:  A version of Chapter 2 has been published as: Pan, K., Phani, A.S., and Green, S., 2018. Particle Dynamics Modeling of the Creping Process in Tissue Making. Journal of Manufacturing Science and Engineering, 140(7), p.071003. The author of this thesis was the principal contributor to this publication, and conducted all the analytical work and simulation. Drs. Phani and Green supervised the research and assisted with writing the paper.  A version of Chapter 3 and Chapter 5 has been published as: Pan, K., Das, R., Phani, A.S., and Green, S., 2019. An Elastoplastic Creping Model for Tissue Manufacturing. International Journal of Solids and Structure, Accepted and available online. The author of this thesis was the principal contributor to this publication, and conducted all the analytical work and simulation, and wrote the manuscript. The experimental part was assisted by Ratul Das (M. Sc. student in the same research group). Drs. Phani and Green supervised the research.  A version of Chapter 2 has formed a manuscript: Pan, K., Phani, A.S., and Green, S., 2019. Periodic Folding of a Falling Viscoelastic Sheet. Submitted. The author of this thesis was the principal contributor of this work, and conducted all the analytical work and simulation, and wrote the manuscript. Drs. Phani and Green supervised the research. vi   A version of Chapter 4 is in preparation: Pan, K., Phani, A.S., and Green, S., 2019. Effects of inhomogeneity on the creping process in tissue making.  The author of this thesis was the principal contributor of this work, and conducted all the analytical work and simulation, and wrote the manuscript. Drs. Phani and Green supervised the research. Serval versions of Chapters 2 to 5 have been presented and published in the following conference proceedings: Pan, K., Phani, A.S., and Green, S., 2016. Particle Dynamics Modeling of Buckle Delamination of Thin Film Materials. Proceedings of the 24th International Congress of Theoretical and Applied Mechanics, pp. 2092-2093, Montreal, Canada. Pan, K., Phani, A.S., and Green, S., 2016.  Mechanics of Creping in Tissue Making: Modeling and Experiments. PACWEST, Jasper, Canada. Pan, K., Phani, A.S., and Green, S., 2016.  Particle Dynamics Modeling of Thin Film Folding and Buckle-Delamination. IMECE, Phoenix, USA. Pan, K., Phani, A.S., and Green, S., 2017. Mechanics of Creping Process in Tissue Making: Modeling and Experiments. Tappi Papercon, Minneapolis, USA. The author of this thesis was the principal contributor to these publications, and conducted all the analytical work and simulation. Drs. Phani and Green supervised the research. The experiments in Chapter 5 were conducted by the author and Ratul Das together, with equal contribution. vii  Table of Contents  Abstract ......................................................................................................................................... iii Lay Summary ............................................................................................................................... iv Preface .............................................................................................................................................v Table of Contents ........................................................................................................................ vii List of Tables ................................................................................................................................ xi List of Figures .............................................................................................................................. xii List of Symbols .............................................................................................................................xv Acknowledgements ................................................................................................................... xvii Chapter 1: Introduction ................................................................................................................1 1.1 Background ..................................................................................................................... 1 1.1.1 Overview ..................................................................................................................... 1 1.1.2 Creping Process in Tissue Manufacturing .................................................................. 1 1.1.3 Creping Mechanisms .................................................................................................. 4 1.1.4 Governing Parameters ................................................................................................. 5 1.2 Literature Review............................................................................................................ 7 1.2.1 Experimental Studies on Creping ............................................................................... 7 1.2.2 Creping Models ........................................................................................................... 9 1.2.3 Buckling and Delamination in Thin Film/Substrate Systems ................................... 11 1.3 Motivation and Objectives ............................................................................................ 14 1.4 Outline........................................................................................................................... 15 Chapter 2: A Particle Dynamics Model of the Creping Process .............................................17 viii  2.1 Introduction ................................................................................................................... 17 2.2 Model Description ........................................................................................................ 18 2.2.1 Viscoelastic Deformation of the Web ....................................................................... 19 2.2.2 Interfacial Delamination: Cohesive Zone Models .................................................... 22 2.2.3 Web-Blade Contact ................................................................................................... 25 2.2.4 Self-Contact of the Web............................................................................................ 25 2.2.5 Numerical Implementation ....................................................................................... 27 2.3 Results and Discussion ................................................................................................. 28 2.3.1 Formation of a Single Fold ....................................................................................... 29 2.3.2 Evolution of Creping Force and Delamination Length............................................. 30 2.3.3 Analytical Solution for Maximum Creping Force .................................................... 31 2.3.4 Periodic Folding ........................................................................................................ 34 2.3.5 Creping Angle Effect ................................................................................................ 36 2.3.6 Mixed-Mode Fracture Effect .................................................................................... 37 2.3.7 Evaluation of the Relative Impact of Parameters ..................................................... 38 2.4 Summary and Conclusions ........................................................................................... 40 Chapter 3: Modeling the Effects of Plasticity and Explosive-Bulk in Creping ......................42 3.1 Introduction ................................................................................................................... 42 3.2 Extensions of Particle Dynamics Model for Plasticity ................................................. 42 3.2.1 Constitutive Model for the Uncreped Paper ............................................................. 42 3.2.2 Elasto-Plastic Bending of the Web ........................................................................... 44 3.3 Effects of Plasticity on the Folding Process ................................................................. 47 3.4 Analysis of Process-Structure-Property Relation by Virtual Tensile Test ................... 49 ix  3.5 Stress-Strain Curve of Tissue Paper ............................................................................. 51 3.6 Parametric Study ........................................................................................................... 52 3.6.1 Creping Ratio Effect ................................................................................................. 52 3.6.2 Adhesion Effect ........................................................................................................ 54 3.6.3 Relative Impact of Parameters .................................................................................. 55 3.7 Analysis of the Explosive-Bulk Regime ....................................................................... 56 3.7.1 Multi-Layer Model.................................................................................................... 56 3.7.2 Folding Process in Explosive-Bulk Regime ............................................................. 58 3.7.3 Phase Diagram of Creping Regime........................................................................... 59 3.8 Summary and Conclusions ........................................................................................... 61 Chapter 4: Effects of Inhomogeneity on Creping .....................................................................63 4.1 Introduction: Causes of Inhomogeneity ........................................................................ 63 4.2 Effects of Forming Fabric Pattern................................................................................. 65 4.2.1 Evolution of Web Morphology ................................................................................. 66 4.2.2 Effects of Initial Wavelength and Amplitude ........................................................... 68 4.2.3 Effects of Base Sheet Properties with Initial Pattern Included ................................. 72 4.3 Non-uniform Basis Weight Distribution ....................................................................... 73 4.3.1 Effect of Variation Length ........................................................................................ 73 4.3.2 Basis Weight Ratio Effect ......................................................................................... 76 4.4 Non-Uniform Adhesion between Web and Yankee ..................................................... 77 4.4.1 Effects of Adhesion Variation Length ...................................................................... 78 4.4.2 Creping Force Behavior ............................................................................................ 79 4.5 Summary and Conclusions ........................................................................................... 80 x  Chapter 5: Experimental Validation ..........................................................................................82 5.1 High-Speed Imaging Study of Creping......................................................................... 82 5.1.1 Experimental Setup ................................................................................................... 82 5.1.2 Observation of Creping Process................................................................................ 84 5.1.3 Parametric Studies Based on Creping Rig ................................................................ 85 5.2 Tensile Test of Tissue Paper ......................................................................................... 88 5.3 Validation of Basis Weight Effect by SEM .................................................................. 90 5.4 Summary and Conclusions ........................................................................................... 90 Chapter 6: Summary and Conclusions ......................................................................................92 6.1 Methodology: Particle Dynamics Model ...................................................................... 92 6.2 Summary of Results ...................................................................................................... 94 6.3 Significance of Current Work ....................................................................................... 97 6.4 Limitations and Future Work ........................................................................................ 98 Bibliography ...............................................................................................................................101 Appendices ..................................................................................................................................112 Appendix A Validation of the Elastic Particle Dynamics Model ........................................... 112 Appendix B Validation of the Elasto-Plastic Beam Bending Model ...................................... 115  xi  List of Tables  Table 1.1 Key parameters controlling the creping process. ............................................................ 7 Table 2.1 Impact of the control parameters on creping. ............................................................... 39 Table 3.1 Impact of key parameters on creping predicted by the visco-elasto-plastic model. ..... 56 Table 5.1 Summary of the uniaxial tensile test results of four grades. ......................................... 89  xii  List of Figures  Figure 1.1 A schematic of the dry-creping process with crepe pattern .......................................... 3 Figure 1.2 Stress-strain curves for un-creped and creped paper ..................................................... 4 Figure 1.3 SEM images of cross-section of tissue paper: (a) shaped-bulk (b) explosive-bulk ...... 5 Figure 2.1 (a) A schematic of creping process. (b) A schematic of discrete particle model. ....... 19 Figure 2.2 (a) Elastic bending system (b) Viscous bending system. ............................................ 20 Figure 2.3 Illustration of cohesive zone concept .......................................................................... 23 Figure 2.4 Illustration of  bilinear cohesive zone model .............................................................. 24 Figure 2.5 A schematic of the first scenario of self-contact ......................................................... 26 Figure 2.6 The evolution of the web morphology during creping ................................................ 29 Figure 2.7 The evolution of creping force and the length of delaminated web ............................ 30 Figure 2.8 A schematic of the web being pushed against the blade ............................................. 32 Figure 2.9 Maximum creping force versus (a) creping velocity and (b) Young's modulus ......... 34 Figure 2.10 Evolution of the web shape during the creping process ............................................ 35 Figure 2.11 The evolution of creping force and delamination length under different angle ........ 36 Figure 2.12 The averaged creping force and the creping wavelength under different angle ........ 37 Figure 2.13 Effects of fracture energy on (a) average creping force (b) creping wavelength ...... 38 Figure 3.1 Bilinear elastoplastic material model associated with kinematic hardening ............... 43 Figure 3.2 Local bending systems formed by three consecutive particles ................................... 45 Figure 3.3 Effects of plastic modulus on (a) creping force and (b) delamination length ............. 48 Figure 3.4 Surface pattern of generated tissue paper under different plastic modulus. ................ 48 xiii  Figure 3.5 (a)-(c) The simulated morphology of the web. (d) Top: crepe pattern; bottom: the stretched web after the uniaxial tensile test. ................................................................................. 50 Figure 3.6 Stress-strain curve corresponding to the virtual tensile test in Fig. 3.5(d). ................. 52 Figure 3.7 Effects of creping ratio on (a) crepe pattern and (b) tensile curves ............................. 54 Figure 3.8 Effects of fracture energy on (a) crepe pattern and (b) tensile curves ......................... 55 Figure 3.9 A schematic of the three-layer model .......................................................................... 57 Figure 3.10 (a)-(d) The simulated morphologies of the web by three-layer model. ..................... 59 Figure 3.11 Phase diagram of the creping regime ........................................................................ 60 Figure 4.1 (a) Surface image of tissue paper. (b) Cross-sectional image of tissue paper. ............ 63 Figure 4.2 Schematics of forming fabric (surface view) used in tissue manufacturing. .............. 65 Figure 4.3 Schematic of the initial web pattern due to the forming fabric ................................... 66 Figure 4.4 The evolution of the web morphology with an initial forming fabric pattern. ............ 67 Figure 4.5 Simulated creping pattern under different initial wavelength. .................................... 69 Figure 4.6 FFT of the creping pattern (a):  λ0 λc⁄ = 0.51. (b): λ0 λc⁄ = 1.02. ............................ 69 Figure 4.7 Tensile stress-strain curves of tissue with different initial wavelength ....................... 71 Figure 4.8 Tensile stress-strain curves of tissue with different initial amplitude ......................... 71 Figure 4.9 Tensile stress-strain curves of tissue with different modulus. ..................................... 73 Figure 4.10 (a) Initial flat web. (b)-(d) Crepe pattern with different variation lengths.. .............. 74 Figure 4.11 Evolution of the creping force with variation length l = 4 mm ............................... 75 Figure 4.12 (a) Tensile stress-train curves of tissue paper with different variation lengths. (b) The crepe pattern at different strains during the tensile test with l = 4 mm ....................................... 76 Figure 4.13 Effects of basis weight ratio on (a) Tensile curves and (b) crepe pattern. ................ 77 Figure 4.14 Effects of adhesion variation length on (a) Tensile curves and (b) crepe pattern ..... 79 xiv  Figure 4.15 The evolution of creping force during creping with λa = 200μm. ........................... 80 Figure 5.1 Lab-scale creping rig and the setup of high-speed imaging study. ............................. 83 Figure 5.2 Snapshots of the creping process captured by high-speed camera. ............................. 85 Figure 5.3 Parametric studies based on lab-scale creping rig ....................................................... 87 Figure 5.4 Stress-strain curves: (a) by experiments (b) by simulaton .......................................... 89 Figure 5.5 SEM cross-sectional images of tissue paper with two different basis weights ........... 90  xv  List of Symbols  𝐸1   Elastic modulus of web 𝐸2   Plastic modulus of web ℎ   Thickness of web 𝜌  Density of web 𝑤  Width of web 𝐿0  Initial length of web Δ𝑡  Time step 𝐵𝑊   Basis weight of web 𝛿  Creping angle 𝑉𝑖𝑛  Yankee surface speed 𝑉𝑜𝑢𝑡  Reeling speed 𝐺𝐼𝐶  Mode I fracture energy of adhesives 𝐺𝐼𝐼𝐶  Mode II fracture energy of adhesives 𝜆  Creping wavelength 𝜆0  Initial wavelength due to forming fabric 𝜆𝑐  Natural wavelength 𝐴  Creping amplitude 𝐴0  Initial amplitude due to forming fabric 𝑁  Number of particles 𝜎𝑦  Yield stress of web xvi  𝑘𝑎  Axial stiffness of web 𝑚  Mass of each particle 𝑎0  Initial distance between particles 𝒇𝑎,𝑖  Axial stretching force on particle i 𝒇𝑑,𝑖  Axial damping force on particle i 𝒇𝑏,𝑖  Bending force on particle i 𝒇𝑣,𝑖  Viscous bending force on particle i 𝝎1,𝑖  Angular velocity between particles i and i+1 𝑘𝑏  Bending stiffness of web 𝛾𝑖  Local bending angle 𝑐𝑎  Phenomenological viscous damping coefficient 𝜇  Dynamic viscosity of web 𝐹𝑐  Creping force 𝐹𝑚  Maximum creping force 𝑙𝑑  Delamination length 𝑘𝑏𝑜𝑛𝑑  Inter-fiber bonding stiffness 𝛿𝑖𝑛𝑡𝑐   Critical inter-fiber bonding length 𝜆𝑎  Adhesion variation length 𝑇𝑌  Yankee surface temperature 𝑡𝑎  Adhesives spray time xvii  Acknowledgements  I would like to express my sincere gratitude to my supervisor, Prof. A. Srikantha Phani, and my co-supervisor, Prof. Sheldon Green, for their enormous support and guidance throughout my PhD study. I benefited a lot from their scientific advice and expertise in the course of this research project. Over the past 4 and half years, many insightful discussions and suggestions from them have helped me to overcome the difficulties and come up with new research ideas. This thesis would not have been possible without their great supervision. I would also like to thank my committee members, Prof. Reza Vaziri and Prof. Richard Kerekes, and my external examiner, Prof. Tetsu Uesaka, for their valuable feedback, advice, and comments. This work was supported by Canadian Natural Sciences and Engineering Research Council (NSERC) through the Collaborative Research and Development grant (501467-16) in partnership with FPInnovations (Canada), Kruger Products (Canada), and Solenis (USA). I would like to take this opportunity to thank Xuejun Zou, Hofan Jang, Jimmy Jong, Daniel Ricard, Francois Drolet from FPInnovations, Anna-Karin Ahlman from Kruger Products, Timothy Patterson from Solenis, and other scientists in the industry for providing samples, giving suggestions and access to their facilities to conduct experiments. I would also like to thank Ratul Das, for helping me to conduct the experiments in this project. Last but not the least, special thanks are owed to my wife, my parents and my sister for their endless support, encouragement and unconditional love throughout my years of study. 1  Chapter 1: Introduction 1.1 Background 1.1.1 Overview Tissue is a low-density paper product widely used in our daily life. Annual global demand for tissue is more than 27 million tons and is annually increasing at 4% [1]. Due to the fierce competition in the industry, tissue manufacturers are continuously trying to increase the quality of their products. Compared to other paper products such as printing paper and packaging paper, tissue has a much higher softness, stretchability and absorbency, but lower strength. These properties are attributed to a key operation named creping in the manufacturing process, which breaks some of the fiber-to-fiber bonds and creates microscale folding structures in the paper. Despite its importance to tissue paper making, the published scientific literature on creping is limited. As one researcher puts it: “perhaps creping is considered too lite or just too complex” [2]. For competitive advantage each manufacturer retains the proprietary knowledge developed in-house through painstaking trial and error testing. This thesis aims at understanding the fundamental mechanics involved in creping by theoretical modeling and experiments, and to provide insights into the optimization of the manufacturing process to the industry.  1.1.2 Creping Process in Tissue Manufacturing         Crepe is a word originally borrowed from French, which means curled and represents a type of thin pancake. When used in the paper industry, crepe paper is a paper product patterned with tiny folds or corrugations. Creping generates these microstructures in low-density papers such as tissue and paper towel. There are three types of creping depending on the moisture content of the paper at creping, including wet creping (60-85% solid content), semidry creping (85-93% solid 2  content), and dry creping (93-97% solid content) [3]. In the case of wet creping, the fibers can still form hydrogen bonds after creping during subsequent drying, which makes the sheet less soft. Among these three conditions, dry creping provides the highest strength reduction and softness [3]. All these types of creping involve pressing the sheet when it is wet, which is often called the “conventional” process. Through-air drying (TAD) is another popular process in tissue making, in which the web is dewatered and dried by non-compressive means, thereby avoiding the compaction that occurs during wet-pressing of the web in the conventional process. Vacuum dewatering and through-air drying are the two common water-removal technologies, which result in a product with improved properties such as bulk, CD stretch and absorbency. The biggest disadvantage of TAD is high energy consumption [4].          Figure 1.1(a) shows a schematic of the conventional dry creping process. In this process, a continuous web of wet paper is pressed and adhered onto a heated rotating cylinder (“Yankee Dryer”) whose surface is sprayed with “creping chemicals”, i.e., a combination of adhesive and release agents [5]. After being dried on the surface of the steam-heated Yankee cylinder, the web is scraped off the surface by a doctor blade and then wound onto a reel followed by off-line downstream converting into consumer products. As a result of its violent interaction with the doctor blade, the internal structure of the sheet is changed significantly and a series of micro-folds (“crepes”) are formed (as shown in Figure 1.1(b)). These micro-folds greatly enhance the softness, absorbency and stretchability of the tissue paper [6, 7]. However, the downside of creping is that the strength drops significantly because some of the fiber-to-fiber bonds are completely broken during creping [8]. Figure 1.2 illustrates the stress-strain curves for un-creped and creped paper. While the strengths in both machine direction and cross direction have decreased after creping, the work-to-rupture (toughness) has increased. We can define a creping wavelength 𝜆 and a creping 3  amplitude 𝐴 based on the periodic crepe pattern, as shown in Figure 1.1(c). Although the actual crepe pattern is not exactly periodic due to the inhomogeneity of the web properties, the topological patterns introduced when fibers are formed into a sheet on a fabric and the non-uniform adhesion etc., a dominant creping wavelength may be identified. Statistically, one may interpret creping wavelength as the mean/median of wavelength distribution sampled over different regions to the tissue. Naturally, it can vary with process parameters including but not limited to furnish, creping geometry, blade angle, blade wear, forming, etc.   Figure 1.1 (a) A schematic of the dry-creping process. (b) The web-blade contact area with the definition of creping angle and creping ratio. (c) A schematic of periodic fold pattern in tissue paper with the definition of creping wavelength and amplitude. DoctorBladeCleaningBladeAdhesive Spray NozzleFeltPressureRollWet WebDryerHoodsDryerHoodsReelDiameter   𝑚Surface Speed  1500𝑚 𝑚  𝑇 100  CrepedTissue20 , SC:  20 SC:  40 100 , SC:   5 𝑉𝑖𝑛𝑉𝑜𝑢𝑡(SC: solid content)Creping Ratio:            Creping Angle: 𝛿Take-Off Angle:  Blade Contact Angle:  Blade Bevel Angle:  (a) (b)(c)0 0.2 0.4 0.6 0.8 1-2-1.5-1-0.500.511.52 𝜆2𝐴Wavelength: 𝜆Amplitude: 𝐴𝑉𝑖𝑛𝑉𝑜𝑢𝑡𝛿   Yankee Dryer4   Figure 1.2 Stress-strain curves for un-creped and creped paper.          1.1.3 Creping Mechanisms        Creping is commonly believed to be a periodic delamination and buckling process. As soon as the web hits the blade, compressive stresses develop within the sheet and within the adhesive layer anchoring the sheet to the Yankee. When the stresses in the adhesive layer attain a failure value, interfacial delamination initiates and the web starts debonding from the Yankee. As the incoming web continues to push against the blade, delamination propagates and the length of the separated web increases but not indefinitely. At a critical length, the delaminated segment becomes unstable and buckles under the compressive force exerted by the blade. After buckling, most of the compressive stress is released. This process repeats when the next segment of bonded web impacts the doctor blade. In practice, non-uniform web and adhesive properties make the creping process a quasi-periodic or, a narrow-band stochastic process.         Creping is believed to occur in two distinct regimes [2]. Figure 1.3 shows the cross-sectional images of typical tissue paper in the machine direction, obtained by a scanning electron microscope (SEM). In the shaped-bulk regime, creping involves the buckling of a cluster of fibers (roughly 4 to 6 fibers in the thickness direction of the sheet), as shown in Figure 1.3(a). Typical creping MDCD MDCDStrainStressMD: Machine DirectionCD: Cross DirectionUN-CREPEDCREPEDStretch ( 20 )5  “wavelength” is around 300 𝜇𝑚 . This regime is significant for intermediate grade tissues of moderate softness. Here, the web can be viewed as a homogenized composite material with its own effective properties. However, in the explosive-bulk regime individual fibers separate (but fibers do not break) from each other producing a controlled web explosion, as shown in Figure 1.3(b). This regime is thought to be significant for ultra-soft tissue grade. Individual properties of the fiber, inter-fiber bonding forces, and fiber network architecture are important factors in the explosive-bulk regime. It remains unclear what governs the transition between these two different regimes. Although one can surmise that weaker web cohesion and stronger adhesion to Yankee can explode the sheet, without entering into picking regime where the sheet is torn, or into passing regime where the sheet simply passes beneath the doctor blade.  Figure 1.3 SEM images of the machine direction cross-section of tissue paper: (a) shaped-bulk regime with the definition of h, A and λ (b) explosive-bulk regime. Note that the sheet in (b) has higher specific volume (bulk).         1.1.4 Governing Parameters         Creping is essentially a de-densification process that reduces the density and increases the specific volume of paper. Various parameters govern this complex process. Table 1.1 lists some of the key parameters with their typical industrial values [8, 3, 9]. The properties of the web including the modulus, thickness, basis weight, density and moisture content, are related to the furnish, i.e., the fibers chosen for the papermaking, as well as the pressing and drying operations prior to creping. It is thought that bleached eucalyptus kraft (BEK) can improve the softness and bulk of tissue paper, while northern bleached softwood kraft (NBSK) helps in maintaining higher (b)(a)200 µm 200 µm2𝐴𝜆ℎ6  strength and runnability (sheet does not break) [10]. The adhesive fracture energy (adhesion) is determined by the amount of adhesive chemicals sprayed onto the Yankee surface, temperature, surface chemistry as well as the pressing. If the adhesion is too low, the web can be easily scraped off by the doctor blade, and the desired micro-folds could not be created. On the other hand if the adhesion is too strong, the blade may only “pick out” part of the web that leads to undesirable holes in the tissue paper. The creping angle is defined as the angle formed between the Yankee surface and the blade surface, as shown in Figure 1.1(b). In paper industry, creping angle is also referred to as pocket angle or impact angle. Experiments have shown that if the creping angle is too low, the micro-folds may pile up and form a macro-fold, which decreases the quality of tissue paper [11]. If the creping angle is too large, the web is simply “peeled off” without creating any crepe. Finally, the creping velocity and the creping ratio are controlled by the paper machine. The creping velocity is defined as the surface speed of the Yankee, denoted by 𝑉𝑖𝑛. The creping ratio is defined as (𝑉𝑖𝑛 − 𝑉𝑜𝑢𝑡) 𝑉𝑖𝑛⁄ , where 𝑉𝑜𝑢𝑡 represents the downstream reeling velocity, which is lower compared to 𝑉𝑖𝑛. The creping ratio needs to be set in the right range to ensure sufficient tension in the web to avoid web break.         All the aforementioned parameters can affect the formation of creping microstructure and thus change the final properties of a tissue paper. It is worth mentioning that these parameters are related to other processes during tissue manufacturing, such as sheet formation, drying and pressing. For example, during the formation process, the forming fabric pattern may be embedded onto the web, thus changing its topographic and network properties. During the drying process, changing the drying time and temperature can change the moisture contents of the web and the adhesives, which affects the modulus of the web and the adhesion. In the pressing process, the pressure applied by the roll can also affect the web properties (bonding strength, modulus, and 7  density) and the adhesion between the web and the Yankee. It may not be realistic to include all the potential inter-dependent factors, so in this thesis we mainly focus on studying the effects of the parameters listed in Table 1.1 on creping. Parameters Symbol Typical Value (Range) Units Young’s modulus 𝐸 100 − 1000 MPa Web thickness ℎ 30 − 100 𝜇𝑚 Web density 𝜌 250 − 400 𝑘𝑔 𝑚3⁄  Web basis weight BW 15 − 25 𝑔 𝑚2⁄  Web moisture content MC 3 − 40   Adhesive fracture energy 𝐺𝐼𝐶, 𝐺𝐼𝐼𝐶 20 − 200 𝑁 𝑚⁄  Creping angle 𝛿 70 − 100 Deg. (°) Creping velocity 𝑉𝑖𝑛  00 − 3000 𝑚 𝑚  ⁄  Creping ratio (𝑉𝑖𝑛 − 𝑉𝑜𝑢𝑡) 𝑉𝑖𝑛⁄  20 − 40     N.A. Table 1.1 Key parameters controlling the creping process.  1.2 Literature Review 1.2.1 Experimental Studies on Creping         Performing experiments on full-scale tissue machine running at a high speed is very costly and risky, as any mistake could be catastrophic. Most previous studies are limited to low speed pilot tissue machines and lab-scale creping simulators [9, 11, 12]. The open literature on creping dates back to 1972 when Hollmark performed a high speed imaging study on a pilot tissue machine [11], which has a maximum running speed of 140 𝑚 𝑚  ⁄ . It was found that as the sheet first contacts the doctor-blade, a series of micro-folds are formed and then thus pile up to become a 8  macro-fold. This micro-fold to macro-fold transition only occurs when the creping angle is relatively small (  0°). In modern tissue machines, creping angle is usually set around 90° during operation, thus this transition mechanism may not be relevant to modern commercial operations. Based on their experiments, Hollmark also concluded that larger creping angle, higher adhesion, and a lower basis weight results in smaller creping wavelength and amplitude [11]. Later in 2000, Ramasubramanian et al. built a lab-scale creping rig and studied the effects of various parameters on creping wavelength and creping force [9]. Here creping force is defined as the tangential component of the force exerted by the blade on the web per unit width. It was reported that both creping force and creping wavelength increase as the creping angle decreases. When the adhesion increases, creping force also increases, but creping wavelength decreases. The creping velocity effect was also studied and it showed that creping force increases as the velocity increases. However, the maximum creping velocity was limited to 250 𝑚 𝑚  ⁄ , which is far below the industrial value. Besides, only a finite paper sample was used and the creped sample was not reeled, meaning that the creping ratio effect was not considered. J. Boudreau and C. Barbier developed a laboratory creping device by employing a creping wagon attached to a tensile tester [12]. They observed a decreasing trend in creping force and creping wavelength when the creping angle is increased, which is consistent with the results of Ramasubramanian [9]. A similar discontinuous experiment has been done by Hämäläinen et al., who developed a creping rig by using a long sled as Yankee surface to achieve much higher creping velocity. This allows the drying time of adhesives to be more realistic. Effects of creping angle, speed, dryer temperature and basis weight were shown to agree with results obtained in tissue industry [13]. To mimic a more realistic creping process with creping ratio effect included, Ho et al. built a pilot-scale creping rig which feeds a continuous web from an unwind stand and reels in the crepe paper after the 9  creping section [14]. It was shown that creping ratio, tensile strength, thickness and stretch of generated tissue paper all fall with increasing web tension and temperature.         Previous experiments have shown that creping involves periodic debonding and buckling of the web. The effects of several process parameters such as adhesion, creping angle, basis weight, Yankee surface temperature, and creping velocity have been investigated. However, most of the experiments focus on the process-structure correlation, the structure-property correlation remains unclear. Besides, these experiments are discontinuous (creping ratio is neglected), and limited to relatively low creping speed (< 250𝑚 𝑚  ). The results obtained at low speed need to be validated at the industrial speed range.   1.2.2 Creping Models         Based on the observations in the experiments, it is commonly accepted that the key physics involved in creping are the interfacial delamination between the web and the Yankee surface, and the subsequent buckling of the web. Significant efforts have been made to construct creping models to capture the basic phenomena. Ramasubramanian et al. developed a one-dimensional continuum model by assuming both the web and the adhesive layers are elastic [15]. The failure of the adhesive layer is described by a strength based failure criterion. In their model, the web is quasi-statically pushed against the doctor blade, which gradually leads to the failure of the adhesive layer and the delamination propagation. The critical delamination length corresponding to buckling initiation is recorded and used to interpret the creping wavelength. Based on their model, effects of process parameters such as creping angle, sheet bending stiffness and adhesive strengths on the creping wavelength and creping force are analyzed. The sheet bending stiffness (𝐸𝐼) and adhesive shear strength are found to have the highest impact on the creping wavelength and creping 10  force. The model also predicts that the creping force increases as creping angle increases. However, this trend does not agree with the experimental observations [9, 12]. This is because the model neglects the post-buckling deformation of the web and underestimates the creping force. Gupta constructed a Dynamic Finite Element Model to simulate the creping process in commercial software Abaqus [16]. The web is modeled as a 2-dimensional elastic thin film. The adhesive layer is modeled by a zero-thickness cohesive zone model that captures the initiation and propagation of the delamination. This finite element model successfully simulated the continuous delamination, buckling and post-buckling of paper during creping. The simulation results show that the creping wavelength decreases as the fracture toughness increases, thus resulting in finer crepes. However, the model does not include the creping ratio effect, i.e., the creped sheet is not pulled away. More recently, a discrete element method was applied to model creping [17]. Specifically, the web is treated as a fiber network and each fiber is represented by a series of spherical particles in different geometries (curl, kinks, twist). The network is created by random deposition of the fibers under gravity and the subsequent consolidation under pressure. Inter-fiber bonds are formed during this formation process. Finally, the generated network adheres to the Yankee surface with a given adhesive strength and hits the doctor blade at 2000 𝑚 𝑚  ⁄ . Extensive destruction (“explosion”) of the web is observed and some fibers are broken and completely disconnected from the network. It shows about 30% inter-fiber bonds are broken immediately after the contact. The model is developed based on an open source code named “ESyS-Particle” [18]. To our knowledge, this is by far the only creping model that specifically accounts for the irreversible deformations in the fibers and inter-fiber bonds. However, this fiber-level model requires many input material parameters such as the stiffness and strength of individual fiber and inter-fiber bonds under various deformation modes. Some of these parameters are not accessible from experiments. Based on the 11  simulation results, it was found “many parameters are highly interactive, and a relatively small window exists for successful creping operation”. The creping wavelength and amplitude were not quantified and how exactly the fiber properties change the creping pattern was not investigated. Besides, if one wants to run a continuous simulation with creping ratio included, the computational cost becomes quite expensive. Nevertheless, this 3-D discrete network model provides a useful design tool to investigate the effects of various parameters and improve the creping operation. It is possible to extend this discrete model to account for the inhomogeneity of web properties and the effect of forming fabric pattern.  1.2.3 Buckling and Delamination in Thin Film/Substrate Systems         Since buckling and delamination are the key physics involved in creping, this section reviews the literature related to buckle-delamination in thin film/substrate systems.                Buckling is an instability that often leads to failure, and delamination is essentially interfacial cracking. In compressed thin film/substrate systems, this instability usually involves the simultaneous buckling of the film and the propagation of interfacial delamination [19, 20, 21]. For instance in composite and coating materials, buckle-delamination is a common failure phenomenon that reduces the structure’s reliability and integrity [22, 23]. Recently, buckle-delamination has been utilized to measure the interfacial adhesion based on the buckling profile [24], to enhance the stretchability of flexible electronics [25, 26], and to design nanochannel networks [27]. Buckle-delamination is a coupled nonlinear process that can be affected by many parameters such as the pre-strain, properties of the film, interfacial fracture energy and boundary conditions [28, 29, 30]. Due to the interplay of these parameters, various interesting surface morphologies have been observed including circular blister, straight-sided blister and the well-12  known telephone cord blister [31, 32, 33, 34, 35, 36, 37]. Intensive studies have been made in the past few decades to understand and control the formation of these intriguing patterns due to buckle-delamination. For example, in one-dimensional buckle-delamination the arrested equilibrium delamination length is determined by the balance between the energy release rate at the delamination front with the mixed-mode interfacial fracture energy [21, 38]. When the film is under two-dimensional compression, it is found that the circular delamination front undergoes configurational instability and becomes asymmetric with the increase of the pre-strain [39]. For the complex telephone-cord blister, it has been shown that the formation can be effectively modeled as a secondary-buckling instability [40, 41, 42, 43, 44].         Unlike the case of Euler’s free buckling, buckle-delamination simultaneously involves the nonlinear buckling deformation of the film as well as the interfacial delamination propagation. In the previous analytical studies of quasi-static or dynamic buckle-delamination, Föppl–von Kármán plate theory is commonly used to model the post buckling deformation of the film and Griffith-type fracture criterion is adopted to describe the interface cracking. However, due to the highly nonlinear nature of the Föppl–von Kármán plate equations and the coupling of the mixed-mode fracture energy with the buckling deformation, these studies are restricted to simple geometry and boundary conditions [21, 38, 39, 45], and quasi-static. A series of numerical methods have been developed to account for more complex buckle-delamination problems. Finite element method incorporated with interface cohesive zone model has been a popular approach that can successfully handle a wide range of buckle-delamination problems [46, 47, 48, 49, 50]. Besides, spring network model coupled with strength-based interface failure criteria has been proposed as an effective numerical method to model telephone-cord blisters [44]. More recently, a continuum phase field model based on energy minimization has been developed to handle buckle-delamination problems. 13  It integrates an interfacial cohesive zone model with the Föppl–von Kármán plate theory and uses a kinematic approach to evolve the surface configuration of the film [29, 30, 51].          However, most of the aforementioned grid-based models are not suitable to simulate the buckle-delamination in creping for the following reasons. First, the web undergoes a large post-buckling deformation and self-contact of the web may occur, thus the commonly used Föppl–von Kármán plate theory becomes invalid. Second, creping is a highly dynamic process where inertia and viscosity becomes important. Most of the models assume that buckle-delamination occurs quasi-statically. Third, the boundary condition of the web during creping is very complex. The flat web contacts with a finite-width blade and the subsequent creped web is reeled downstream at lower outgoing velocity, creating a web tension. Finally, paper is essentially a stochastic fiber network and the heterogeneities cannot be conveniently described within a macroscopic continuum mechanics framework. Given these considerations, discrete element method (DEM) seems to be a better approach to model the creping process. In discrete element method the material is treated as an assembly of discrete mass particles, and the interactions among the neighboring particles are defined by the constitutive equations [52, 53]. Particle-based models have advantages to handle discrete and non-homogeneous materials, highly nonlinear dynamic problems, complex geometries, large deformations, friction/contact, arbitrary restrictions, and fracture. Besides, in discrete models we can incorporate other effects such as gravity and adhesion simply by adding the corresponding forces to each particle. In this thesis, particle-based models are developed to simulate and investigate the creping process. First, a one-dimensional model is proposed to capture the bending and stretching of the web in the machine direction. Then, the model is extended to multi-layer model to account for the inter-layer debonding. The torsion and cross direction deformation are neglected, considering that the debonding and buckling mainly occurs along the 14  machine direction during creping. This is also consistent with the fact that, paper properties in machine direction has a more significant change after creping, as shown in Figure 1.2. The details of these models are discussed in Chapters 2 and 3.  1.3 Motivation and Objectives         Based on the literature review, the previous creping experiments are limited to low-speed and lack of creping ratio effect. The mechanism of a continuous high-speed creping process needs to be identified. On the other hand, most of the existing creping models are limited to the shaped-bulk regime in creping. Plasticity of the web, creping ratio, inhomogeneity and explosive-bulk are not investigated. Although these studies have shown some agreement with industrial practice, a systematic comparison and evaluation of the key parameters is still missing. Besides, the correlation between the creping structure (amplitude, wavelength) and the property (softness, stretch) remains unclear. Based on these limitations, the main objectives of this thesis are listed below: 1. Develop theoretical and computational models that include all the essential features of creping, and investigate the underlying mechanisms. 2. Perform systematic parametric studies to evaluate the impacts of various control parameters in creping, and provide guidance to optimize the process. 3. Investigate both the shaped-bulk and explosive-bulk regimes in a unified model, and identify the factors governing different creping regimes. 4. Study the process-structure-property correlation of tissue paper. 5. Perform experiments to reveal the mechanism of a continuous high-speed creping process, and validate the proposed models. 15  1.4 Outline         This dissertation is arranged in six Chapters. A brief description of each Chapter is given as follows: Chapter 1: Introduction. This Chapter has given a background of the creping process in tissue manufacturing. The literature review consists of the experimental and modeling work on creping, and an introduction of buckle-delamination phenomenon in thin film/substrate systems. Based on the limitations of previous studies, the objectives of the thesis have been identified. Chapter 2: A Particle Dynamics Model of the Creping Process. In this Chapter, a one-dimensional particle dynamics model is developed to simulate the shaped-bulk regime in the creping process. The web is treated as an elastic beam and discretized into discrete particles. Parametric studies are performed and an analytical solution of the creping force is derived. Chapter 3: Modeling the Effects of Plasticity and Explosive-Bulk in Creping. In Chapter 3, we extend the elastic model presented in Chapter 2 to include the plasticity of the web and study its effect on creping. A virtual tensile test of the generated tissue paper is performed to understand the structure-property relation. Then, we further extend the one-layer model to be multi-layer in order to simulate the explosive-bulk regime. A phase diagram of creping regime is constructed. Chapter 4: Effects of Inhomogeneity on Creping. In this Chapter, the effects of inhomogeneity on the creping process are investigated. Specifically, three different causes of inhomogeneity are considered, including the non-flat forming fabric pattern embedded within the base sheet, the non-uniform basis weight distribution and the non-uniform adhesion. Chapter 5: Experimental Validation. In Chapter 5, we report the experiments performed to validate the results obtained by modeling. Based on a lab-scale creping rig, we conduct a high-16  speed imaging study to directly observe the creping mechanism. We perform tensile tests on real tissue paper and compare the experimental data with simulation results obtained in Chapter 3.  Chapter 6: Summary and Conclusions. The main findings obtained in this dissertation are summarized in Chapter 6. The novelty and significance of our study are presented. Future work on modeling the creping process is proposed based on the limitations and findings of the research in this thesis. 17  Chapter 2: A Particle Dynamics Model of the Creping Process 2.1 Introduction         A typical creping process involves several coupled mechanical processes: the buckling of the web, the delamination of the interface and the contact between the paper and the blade. Self-contact of the web may also occur when the creping angle is relatively small [6, 8]. In this Chapter, we treat the web as a homogeneous elastic beam and develop a one-dimensional particle dynamics model to account for the above processes. Compared to grid-based method, discrete element method has advantages to handle discrete and non-homogeneous materials, highly nonlinear dynamic problems, complex geometries, large deformations, friction/contact, arbitrary restrictions, and fracture. Discrete element method has been recently applied to model the mechanical behavior of paper materials. The discretization can be performed in two approaches. In the first approach, the paper is treated as a fiber network and each fiber is discretized into a series of discrete particles. Then individual fibers are linked together by inter-fiber bonds to form the network [17, 54, 55, 56]. In this approach, the effects of fiber properties on the mechanical properties of paper such as Young’s modulus, strength and fracture toughness can be explicitly accounted for. In another approach, paper is treated as a continuous beam with its own effective Young’s modulus and bending stiffness [57]. Then the beam is discretized as a series of particles. This approach can be used to model strain distribution, buckling instability and dynamics of paper under various loading conditions. In this thesis, we choose the latter approach to model the creping process due to the following reasons. Before creping, the base sheet is formed by dewatering the fiber suspension on the forming fabric, followed by pressing and drying. The properties of fiber and inter-fiber bond depend on various mechanical and chemical processes, which are extremely difficult to measure and are uncertain. Besides, the creping operation is almost uniform in the cross 18  direction, and the major changes of paper structure occur in the machine direction. Thus, it is reasonable to treat the web as one-dimensional beam during creping. Finally, a 3-D fiber network requires much higher computational cost and more input parameters, which makes it challenging to investigate the effects of all the parameters in a systematic manner. In this Chapter, a one-dimensional particle dynamics model is presented and used to investigate the creping process.  2.2 Model Description         Figure 2.1(a) shows the schematic of a creping process. The Yankee surface is assumed to be flat as the radius of the Yankee ( 3𝑚) is much larger than the web thickness ( 100𝜇𝑚). Since the Yankee is made of cast-iron, its deformation can be neglected during creping. The Yankee is moving with constant velocity 𝑉𝑖𝑛 beneath the doctor blade. The web is initially adhered to the Yankee surface by adhesives and pressing. Figure 2.1(b) represents a discrete model, which models the web as a series of discrete particles connected by viscoelastic elements. Similar discrete particle methods have been successfully applied to model the dynamics of paper web running in the open-draw section [57], as well as in computer graphics to simulate deformable objects [58]. The adhesive layer is modeled by a discrete cohesive zone model, represented by the springs connecting the web and the Yankee surface. The discretization of the web is not unique and alternate discrete particle models exist in the literature [59, 60, 61]. For example, it is possible to prescribe each particle with three degrees of freedom (two translations in the plane and a rotation about an out-of-plane axis) and furnish the link between two adjacent particles with bending moment [54, 62]. In contrast, for the discretization shown in Figure 2.1(b) each particle has two degrees of freedom for planar deformation. The bending effect is accounted for by the additional 19  rotational spring. Nonetheless, discrete cohesive zone modeling scheme proposed here can be embedded in most of particle models.  Figure 2.1 (a) A schematic of the creping process at the web-blade contact area. 𝑭𝒄 represents creping force. (b) A schematic of the corresponding discrete particle model. The particle i is subject to normal and tangential cohesive forces 𝒇𝒄𝒏,𝒊 and 𝒇𝒄𝒕,𝒊 due to the adhesion and 𝒓𝒊 denotes the position vector of particle i.  2.2.1 Viscoelastic Deformation of the Web         The viscoelasticity of the web is described by Kelvin-Voigt model, which has been shown to be suitable for modeling the viscoelasticity of wood products [63]. Phenomenological model using a Wiechert spring-dashpot representation has also been applied to capture the viscoelasticity of paper [64]. Before creping, the web is dried and the solid content of the web is increased from 20-40% to almost 95%. During the rapid shrinkage process, significant tensile strain is developed within the web, which is on the order of 2-3% [65]. Since the failure strain of the base sheet is also in this range, the potential effect of residual stress is significant. When the web hits the doctor blade, large compressive force is applied on the web. Part of the compressive force is used to Yankee(b)0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 2.2x 10-3-2-101234567x 10-4YankeeWeb(a)Adhesive LayerDoctorBlade𝛿𝑛𝛿𝑡𝒇𝑐𝑡,𝑖𝒇𝑐𝑛,𝑖 +1   +1 +2 +2𝑉𝑖𝑛𝐹𝑐  𝒓𝑖20  overcome the initial residual tensile stress, before it compresses the web. Thus, the existence of residual stress potentially increases the creping force. For simplicity, we assume the web is stress-free before it hits the doctor blade. The initial spacing between two adjacent particles is given by 𝑎0 = 𝐿0 (𝑁 − 1)⁄  and the mass of each particle is given by 𝑚 = 𝜌ℎ𝑎0𝑤, where 𝐿0 is the initial length of the web, 𝑁 is the total number of particles, 𝑤 is the width, 𝜌 is the density and ℎ is the thickness. The stiffness of the axial spring is given by 𝑘𝑎 = 𝐸𝑤ℎ 𝑎0⁄  where E is Young’s modulus. After deformation, the axial elastic force acting on particle   is given by: 𝒇𝑎,𝑖 = 𝑘𝑎(|𝒓𝑖 1 − 𝒓𝑖| − 𝑎0)𝒆𝑖,𝑖 1 + 𝑘𝑎(|𝒓𝑖+1 − 𝒓𝑖| − 𝑎0)𝒆𝑖,𝑖+1                            (2.1) where 𝒓𝑖  represents the position vector and 𝒆𝑖,𝑖+1  is the unit vector pointing from particle   to particle  + 1 . |𝒓|  denotes the magnitude of 𝒓 . Throughout this thesis a letter in bold font represents a vector, otherwise it represents a scalar. The particle is also subject to an axial viscous force:  𝒇𝑑,𝑖 = −𝑐𝑎[(?̇?𝑖 − ?̇?𝑖 1) ∙ 𝒆𝑖,𝑖 1]𝒆𝑖,𝑖 1 − 𝑐𝑎[(?̇?𝑖 − ?̇?𝑖+1) ∙ 𝒆𝑖,𝑖+1]𝒆𝑖,𝑖+1                     (2.2) where ?̇?𝑖 represents the velocity of particle   and 𝑐𝑎 represents the phenomenological viscous damping coefficient with 𝑐𝑎 = 𝜇𝑤ℎ 𝑎0⁄ . Here 𝜇 is the dynamic viscosity of the web.  Figure 2.2 (a) A schematic of elastic bending forces in the local bending system formed by particles (i-1)−(i)−(i+1). (b) A schematic of viscous bending forces: 𝒇𝒗𝟏,𝒊 is due to the angular velocity 𝝎𝟏,𝒊, 𝒇𝒗𝟐,𝒊 is due to the angular velocity 𝝎𝟐,𝒊, and 𝒇𝒗,𝒊 is the negative sum of these two.  (a)𝛾𝑖 -1  +1𝒇𝑏1,𝑖 𝒇𝑏2,𝑖𝒇𝑏0,𝑖 𝑖  𝑖 -1   +1𝝎1,𝑖 𝝎2,𝑖𝒇𝑣1,𝑖 𝒇𝑣2,𝑖𝒇𝑣,𝑖(b)21          Since the particles are not endowed with rotational degree of freedom, we could not directly apply the moment. Instead, we introduce bending forces to the particles to achieve the same bending effect. We start with a local bending system formed by three consecutive particles (i-1)−(i)−(i+1), shown in Figure 2.2(a). The system tends to restore to its initial flat state due to the bending effects, thus we can imagine there are two restoring forces 𝒇𝑏1,𝑖 and 𝒇𝑏2,𝑖 acting on the left particle (i-1) and the right particle (i+1). The restoring force is perpendicular to the axial direction and its magnitude is related to the local bending angle 𝛾𝑖, given by: |𝒇𝑏1,𝑖| = |𝒇𝑏2,𝑖| =𝑘𝑏(𝛾0 𝛾 )𝑎02(1+𝜀 )2                                                       (2.3) where 𝑘𝑏 = 𝐸𝑤ℎ3 12⁄  is the bending stiffness, 𝛾0 = 𝜋  is the initial bending angle and 𝜀𝑖 =(|𝒓𝑖+1 − 𝒓𝑖| − 𝑎0) 𝑎0⁄  is the axial strain. The detailed derivation of the bending forces can also be found in [57]. Since the restoring forces are internal forces, there should be a force acting on the central particle i to keep the local bending system in equilibrium, given by 𝒇𝑏0,𝑖 = −𝒇𝑏1,𝑖 − 𝒇𝑏2,𝑖. Note that 𝒇𝑏0,𝑖 is the bending force applied on particle i corresponding to the bending system (i-1)−(i)−(i+1) only. The particle i is also subject to the restoring force 𝒇𝑏2,𝑖 1 corresponding to the bending system (i-2)−(i-1)−(i), and the restoring force 𝒇𝑏1,𝑖+1  corresponding to the bending system (i)−(i+1)−(i+2). Thus the total force acting on particle i due to bending is: 𝒇𝑏,𝑖 = 𝒇𝑏0,𝑖 + 𝒇𝑏2,𝑖 1 + 𝒇𝑏1,𝑖+1.                                                  (2.4) As the sheet buckles, the elastic bending forces tend to restore the sheet to its initial flat state. Concurrently, there are viscous bending forces slowing down this process. The viscous bending force is perpendicular to the axial direction and proportional to the angular velocity but opposite in direction, as shown in Figure 2.2(b). 𝒇𝑣1,𝑖 is due to the angular velocity 𝝎1,𝑖 between particle i and i-1. Similarly, 𝒇𝑣2,𝑖 is due to the angular velocity 𝝎2,𝑖 between particle i and i+1. Note that the 22  angular velocity only depends on the relative motion between two adjacent particles only, thus 𝝎1,𝑖 and 𝝎2,𝑖 are not correlated. The magnitude of the viscous bending force is given by [57]: |𝒇𝑣1,𝑖| = 𝑐𝑎𝜔1,𝑖ℎ212𝑎0,    |𝒇𝑣2,𝑖| = 𝑐𝑎𝜔2,𝑖ℎ212𝑎0.                                         (2.5) As before, to maintain the overall viscous bending force equilibrium the negative sum of 𝒇𝑣1,𝑖 and 𝒇𝑣2,𝑖 is applied on the central particle: 𝒇𝑣,𝑖 = −𝒇𝑣1,𝑖 − 𝒇𝑣2,𝑖. Unlike the elastic bending force, the viscous bending force applied on particle i only depends on its adjacent particles.   2.2.2 Interfacial Delamination: Cohesive Zone Models         In traditional linear elastic fracture mechanics (LEFM), a pre-crack is assumed and based on that, the stress field around the crack front is calculated. Griffith-type fracture criterion is usually used to predict whether the crack propagates. The limitations of LEFM are that the stress solution has singularity at the crack tip and it could not predict the initiation of the crack. To overcome these issues, the concept of cohesive zone was first introduced by Barenblatt [66, 67] and Dugdale [68]. In a cohesive zone model (CZM), fracture occurs in a region in front of the crack tip, called “cohesive zone”. Within the “cohesive zone”, the traction between two interfaces depends on the separation and is limited by its cohesive strength, as shown in Figure 2.3. Different traction-separation laws have been proposed to solve various fracture problems. Tvergaard and Hutchinson applied a Trapezoidal CZM to model crack growth in elastic-plastic materials and peeling of adhesive joints [69]. Needleman introduced a Polynomial and Exponential CZM to model decohesion in particle-matrix [70]. Geubelle and Baylor proposed to use Bilinear CZM to delamination in composite materials [71]. More recently, rate-dependent CZM was proposed to 23  handle the rate dependence of fracture [72]. In order to capture the coupling effect of mode I (opening mode) and mode II (shear mode) fracture, a mixed-mode CZM can be applied [73].  Figure 2.3 Illustration of Cohesive Zone Concept.           Here we adopt a Bilinear CZM (shown as Figure 2.4) to model the delamination in creping, due to its flexibility of changing the strength and toughness independently. The interfacial tractions (𝜎 and 𝜏) are a function of the separations (𝛿𝑛 and 𝛿𝑡). When the separations exceed the critical values (𝛿𝑛𝑐 and 𝛿𝑡𝑐), the damage initiates. Once the maximum values (𝛿𝑛𝑚 and 𝛿𝑡𝑚) are reached, the interfacial tractions drop to zero. The fracture energies for mode I and mode II are: 𝐺𝐼𝐶 =12𝜎𝑐𝛿𝑛𝑚, 𝐺𝐼𝐼𝐶 =12𝜏𝑐𝛿𝑡𝑚.                                                  (2.6) 𝐺𝐼𝐶 and 𝐺𝐼𝐼𝐶 represent the energy required to completely break a unit area in shear direction and normal direction, respectively. A mixed-mode criterion can be adopted to determine the failure of the interface: 𝐺𝐼𝐺𝐼𝑐+𝐺𝐼𝐼𝐺𝐼𝐼𝑐≥ 1                                                              (2.7) Here 𝐺𝐼 and 𝐺𝐼𝐼 are the energy release rates in mode I and mode II [74]: 𝐺𝐼 = {0,         𝛿𝑛 ≤ 𝛿𝑛𝑐  𝛿  𝛿 𝑐𝛿 𝑚 𝛿 𝑐 𝐺𝐼𝐶 , 𝛿𝑛𝑐 < 𝛿𝑛 < 𝛿𝑛𝑚𝐺𝐼𝑐,        𝛿𝑛 ≥ 𝛿𝑛𝑚                                                (2.8) CrackCohesive Zone𝛿Crack Tip24  𝐺𝐼𝐼 = {0,         𝛿𝑡 ≤ 𝛿𝑡𝑐 𝛿  𝛿 𝑐𝛿 𝑚 𝛿 𝑐 𝐺𝐼𝐼𝐶 , 𝛿𝑡𝑐 < 𝛿𝑡 < 𝛿𝑡𝑚𝐺𝐼𝐼𝐶 ,        𝛿𝑡 ≥ 𝛿𝑡𝑚                                                (2.9)  Figure 2.4 Bilinear cohesive zone model: (a) bilinear relation between 𝝈 and 𝜹𝒏. (b) bilinear relation between 𝝉 and 𝜹𝒕 . The shaded area represent the energy release rates 𝑮𝑰  (mode I) and 𝑮𝑰𝑰  (mode II). The total area represent the fracture energies 𝑮𝑰𝒄 and 𝑮𝑰𝑰𝒄.  The cohesive forces applied on each particle, as shown in Figure 2.1(b), are proportional to the normal traction 𝜎 and the shear traction 𝜏: |𝒇𝑐𝑛,𝑖| = 𝜎𝑤𝑎0,  |𝒇𝑐𝑡,𝑖| = 𝜏𝑤𝑎0.                                              (2.10) In creping the adhesive layer is much thinner compared to the web, so it is reasonable to assume the fracture is a “small-scale yielding” process [75]. Thus, the fracture energies (𝐺𝐼𝐶 and 𝐺𝐼𝐼𝐶) are the dominant parameters governing the fracture process of the adhesive layer. We fix 𝜎𝑐 = 𝜏𝑐 =0.01𝐸 and 𝛿𝑛𝑐 = 𝛿𝑡𝑐 = 0.01ℎ in this thesis and only change the fracture energies by altering the maximum separations. The values of cohesive strengths are based on literature [15, 16], while the critical separation lengths are estimated, as they are extremely difficult to measure experimentally. However, since this is a “small-scale yielding” process, the exact values of 𝛿𝑛𝑐 and 𝛿𝑡𝑐 should not matter. In practice, the fracture energies can be varied by tuning the composition of the adhesive polymers.  (b)𝜎𝜎𝑐𝐺𝐼𝐺𝐼𝑐(a)𝛿𝑛𝑐 𝛿𝑛𝑚𝛿𝑛0𝐺𝐼𝐼𝐺𝐼𝐼𝑐𝜏𝜏𝑐0 𝛿𝑡𝑐 𝛿𝑡𝑚𝛿𝑡Damage initiationDamage initiation25  2.2.3 Web-Blade Contact          During creping, the blade is fixed and its deformation is negligible. In our discrete model, the blade is treated as a flat surface and the contact between the web and the doctor blade is described by a simple penalty model. The blade contact force 𝒇𝑏𝑐,𝑖 applied on the particle i is given by:  |𝒇𝑏𝑐,𝑖| = {𝑘𝑐𝑑𝑝,   𝑑𝑝 > 0     0,     𝑑𝑝 = 0                                                         (2.11) where 𝑘𝑐  is the contact stiffness and 𝑑𝑝  is the penetration distance. Since 𝑑𝑝  should be much smaller compared to the thickness of the web,  𝑘𝑐  is chosen three orders larger than the axial stiffness of the web. Based on the Coulomb friction law, the friction force is given by |𝒇𝑏𝑓,𝑖| ≤𝜇𝑓|𝒇𝑏𝑐,𝑖| with 𝜇𝑓 being the friction coefficient. In this thesis 𝜇𝑓 is fixed as 0.3 for simplicity. The same penalty model is employed to prevent the penetration of the particles into the substrate.  2.2.4 Self-Contact of the Web          When the creping angle is small (space is more constrained) or when the creping ratio is high (more material is fed), the micro-folds formed may interact with each other. Under these circumstances self-contact of the web needs to be considered. Similar to the blade contact model, we use a penalty-based model to consider the self-contact effect by introducing self-contact forces to the corresponding particles. There are two scenarios for self-contact in the model. In the first scenario, shown as Figure 2.5(a), the line AB connecting adjacent particles A and B represents a small section of the web. As the web should not interpenetrate, an arbitrary particle C is forbidden to pass through the line AB. To enforce this condition we introduce a self-contact force 𝒇𝑠0 to the particle C normal to the line AB. The magnitude of the self-contact force is related to the distance of particle C to the line AB, given by: 26  |𝒇𝑠0| = {𝑘𝑠[𝑒(𝑑𝑐𝑑𝑠)− 𝑒],     𝑑𝑠 < 𝑑𝑐             0,               𝑑𝑠 ≥ 𝑑𝑐                                               (2.12) where 𝑑𝑠 is the distance of particle C to AB, 𝑑𝑐 is a critical distance and 𝑘𝑠 is a parameter with the dimensions of a force. We found that in our model the simulation results are insensitive to 𝑘𝑠 provided that the simulation converges. Also, the critical distance 𝑑𝑐 has negligible effect on the simulation provided that it is much smaller than the thickness of the sheet. Throughout this paper, we choose 𝑘𝑠 = 𝐸𝑎0𝑤 and 𝑑𝑐 = ℎ 50⁄ . Since the self-contact force is an internal force, we need to apply reacting self-contact forces 𝒇𝑠1 and 𝒇𝑠2 on the particles A and B so that the total force and the total moment caused by the self-contact forces vanish. These two reacting forces are given by: {𝒇𝑠1 = −𝑙2𝑙1+𝑙2𝒇𝑠0𝒇𝑠2 = −𝑙1𝑙1+𝑙2𝒇𝑠0                                                           (2.13)                                                           where 𝑙1 and 𝑙2 are defined in Figure 2.5(a). Note that when particle C approaches segment AB, there is also tangential velocity component. This means that sliding can occur during self-contact of the web. The frictional forces applied on particles A, B and C are related to the self-contact forces, according to the Coulomb friction law.  Figure 2.5 (a) Schematic of the first scenario of self-contact, where 𝒇𝒔𝟎 is the self-contact force on particle C, 𝒇𝒔𝟏 and 𝒇𝒔𝟐 are the reactive self-contact forces. 𝒅𝒄 is the critical distance for self-contact and 𝒅𝒔 is the current distance between particle C and the line AB. (b) Schematic of the second scenario of self-contact, where 𝒇𝒔𝟎 is the self-contact force on particle D and −𝒇𝒔𝟎 is the reactive self-contact force on particle B.  ABC𝑑𝑐𝑑𝑠𝒇𝑠0ABCD𝑑𝑐𝒇𝑠1𝒇𝑠2𝑙1𝑙2𝒇𝑠0−𝒇𝑠0(a) (b)𝑑𝑠27          Figure 2.5(b) shows the second scenario for self-contact. The line AB and line BC formed by three consecutive particles A, B and C represent two adjacent sections of the sheet. To prevent the particle D approaching and passing through the particle B, we introduce a self-contact force 𝒇𝑠0 to particle D: |𝒇𝑠0| = {𝑘𝑠[𝑒(𝑑𝑐𝑑𝑠)− 𝑒],    𝑑𝑠 < 𝑑𝑐           0,               𝑑𝑠 ≥ 𝑑𝑐                                               (2.14) where 𝑑𝑠 represents the distance between particle D and particle B. In this case the reacting force on the particle B is simply −𝒇𝑠0. During the simulation, at each time step for every particle i we check the potential self-contact with its nearby particles. Then we calculate the self-contact forces and the reactive self-contact forces applied on it and sum them to get the total self-contact force 𝒇𝑠,𝑖.  2.2.5 Numerical Implementation         We start the simulation before the web hits the blade. The web is discretized into N particles with the spacing 𝑎0 = 𝐿0 (𝑁 − 1)⁄ . At time 𝑡 = 0, all the particles are perfectly connected to the corresponding substrate nodes. The initial velocity of the particles equals the Yankee surface speed, 𝑉𝑖𝑛. The total force applied on the particle i includes the axial elastic force 𝒇𝑎,𝑖, the axial viscous damping force 𝒇𝑑,𝑖, the bending force 𝒇𝑏,𝑖, the viscous bending force 𝒇𝑣,𝑖, the cohesive forces 𝒇𝑐𝑛,𝑖  and 𝒇𝑐𝑡,𝑖 ,the self-contact force 𝒇𝑠,𝑖 , the contact force with the blade 𝒇𝑏𝑐,𝑖 , and the friction force 𝒇𝑏𝑓,𝑖. We evolve the particle system according to Newton’s equation of motion: 𝒇𝑡𝑜𝑡,𝑖 = 𝒇𝑎,𝑖 + 𝒇𝑑,𝑖 + 𝒇𝑏,𝑖 + 𝒇𝑣,𝑖 + 𝒇𝑐𝑛,𝑖 + 𝒇𝑐𝑡,𝑖 + 𝒇𝑠,𝑖 + 𝒇𝑏𝑐,𝑖 + 𝒇𝑏𝑓,𝑖 = 𝑚?̈?𝑖.             (2.15) A semi-implicit Euler method is adopted for the advancement in time, which converges quickly and shows good stability at low computational cost [54, 76, 77, 78]. Given the particle’s position 28  vector 𝒓𝑖𝑡  and velocity ?̇?𝑖𝑡  at time t, the total force 𝒇𝑡𝑜𝑡,𝑖𝑡  applied on the particle i and the acceleration ?̈?𝑖𝑡 can be calculated. Then the particle’s velocity and position at time 𝑡 + Δ𝑡 can be given by: ?̈?𝑖𝑡 =𝒇   ,  𝑚,                                                               (2.16) ?̇?𝑖𝑡+Δ𝑡 = ?̇?𝑖𝑡 + ?̈?𝑖𝑡Δ𝑡,                                                           (2.17) 𝒓𝑖𝑡+Δ𝑡  = 𝒓𝑖𝑡 + ?̇?𝑖𝑡+Δ𝑡Δ𝑡,                                                      (2.18) where Δ𝑡 is the time step. At each time step, the locations of the substrate nodes are updated according to the Yankee surface velocity. For each simulation, we check the convergence by reducing the time step until the results do not change. The time step used in our creping simulation is typically around 1  𝑠 . The characteristic time of the wave propagation is given by 𝑎0 √𝐸 𝜌⁄⁄  10 𝑠, with 𝑎0 10 𝜇𝑚, 𝐸 100 𝑃𝑎 and 𝜌 300 𝑘𝑔 𝑚3⁄ . The smallest time scale for the contact is related to the self-contact stiffness 𝐸𝑎0𝑤 𝑑𝑐  and given by √𝑚 (𝐸𝑎0𝑤 𝑑𝑐)⁄  1.5  𝑠, where 𝑚 = 𝜌ℎ𝑎0𝑤 is the mass of one particle, and 𝑑𝑐 ℎ 50⁄  is the critical distance for self-contact. Thus, the Courant–Friedrichs–Lewy condition is satisfied [79]. The current discrete model is validated by modeling the buckle-delamination process of a compressed thin film on a rigid substrate. The details can be found in the Appendix A.   2.3 Results and Discussion         We investigate the underlying mechanism of fold formation and the effects of control parameters on the creping wavelength, the amplitude and the creping force. Here the creping force is defined as the tangential force applied on the web per unit width by the doctor blade. This force 29  can be directly measured on an operating tissue machine to monitor the creping process. The symbols and values of the parameters used in this Chapter are based on Table 1.1.  2.3.1 Formation of a Single Fold   Figure 2.6 The evolution of the web morphology during creping: (a) 𝒕 = 𝟑𝝁𝒔, (b) 𝒕 = 𝟖𝝁𝒔, (c) 𝒕 = 𝟏𝟎𝝁𝒔, (d) 𝒕 = 𝟐𝟓𝝁𝒔. (The scale bar only applies for the length of the web and the deflection of the web. To aid in visualization, the thickness of the web and adhesive layer are not scaled)          Figure 2.6 shows the simulated morphologies of the paper web during the creping process. The parameters are chosen as: ℎ = 50 𝜇𝑚 , 𝑤 = 1𝑚 , 𝐿0 = 10 𝑚𝑚 , 𝐸 = 100 𝑃𝑎 , 𝜌 =300 𝑘𝑔 𝑚3⁄ , 𝜇 = 20 𝑘𝑔 𝑚𝑠⁄ ,𝑉𝑖𝑛 = 1200 𝑚 𝑚  ⁄ , 𝛿 = 85°, 𝐺𝐼𝐶 = 𝐺𝐼𝐼𝐶 = 150 𝑁 𝑚, 𝑁 =  01, Δ𝑡 = 1  𝑠. Our simulation indicates there are three typical stages. At time 𝑡 = 3𝜇𝑠 the web first contacts the blade, shown in Figure 2.6(a). As the web continues to be pushed against the blade, interfacial delamination initiates and propagates along the interface, illustrated in Figure 2.6(b). This is the first stage of creping, i.e., delamination propagation. In the second stage, when the delamination length reaches a critical value, the web begins to buckle, shown in Figure 2.6(c). At 1.8 2 2.2 2.4 2.6 2.8x 10-3-20246810x 10-41.8 2 2.2 2.4 2.6 2.8x 10-3-20246810x 10-4(a) (b)(c) (d)100 𝜇𝑚100 𝜇𝑚100 𝜇𝑚 100 𝜇𝑚Stage IDelamination propagationStage IIBuckling initiationStage IIIPost-buckling deformation1.8 2 2.2 2.4 2.6 2.8x 10-3-20246810x 10-4YankeeWebAdhesive LayerDoctorBlade𝑉𝑖𝑛 𝛿1.8 2 2.2 2.4 2.6 2.8x 10-3-20246810x 10-430  the third stage the web undergoes post-buckling deformation until self-contact occurs, shown in Figure 2.6(d). At this point, the formation of one fold is completed.   2.3.2 Evolution of Creping Force and Delamination Length         During the simulation, we have recorded the creping force 𝐹𝑐 applied by the doctor blade and the delamination length 𝑙𝑑. Figure 2.7(a) shows that at the first stage, the creping force 𝐹𝑐 quickly reaches the maximum value 𝐹𝑚, which triggers the delamination to propagate in a purely shear mode (mode II fracture). At the beginning of the second stage, the delaminated web becomes unstable and buckles. Since buckling releases the compressive stress in the web, the creping force drops significantly. The delamination length remains constant during this stage, as shown in Figure 2.7(b). Finally, during the third stage the web undergoes post-buckling deformation, which drives the delamination to propagate again. However, this delamination process is mixed-mode, i.e., the interfacial crack opens in both the normal and tangential directions. It is worth mentioning that currently in practice, only the averaged creping force can be measured because the data acquisition frequency can only reach 10 KHz, which is lower than the creping frequency.  Figure 2.7 The evolution of the creping force (a) and the length of delaminated web (b) during the formation of a single fold. Stage I Stage II Stage III Stage I Stage II Stage III(a) (b)1.8 2 2.2 2.4 2.6 2.8x 10-3-20246810x 10-41.8 2 2.2 2.4 2.6 2.8x 10-3-20246810x 10-41.8 2 2.2 2.4 2.6 2.8x 10-3-20246810x 10-41.8 2 2.2 2.4 2.6 2.8x 10-3-20246810x 10-41.8 2 2.2 2.4 2.6 2.8x 10-3-20246810x 10-41. 2 2.2 2.4 2.6 2.8x 10-3-20246810x 10-431  2.3.3 Analytical Solution for Maximum Creping Force         Creping force is an important parameter that affects the tissue quality. During the first stage of creping, the creping force remains at the maximum value, shown as Figure 2.7(a). The maximum creping force is an important value as it sets the critical length of dynamic buckling, which is related to the final creping wavelength. Here we use energy analysis to derive an analytical solution for the maximum creping force for the first fold and compare it with the simulation. At the first stage, there are mainly two factors contributing to the creping force, i.e., the interfacial debonding and the inertia of the web. Here we consider them separately and then combine them to get the total creping force. Consider a web being pushed towards the blade quasi-statically, shown in Figure 2.8. When the web contacts the blade, the creping force first reaches the maximum 𝐹𝑚. AB represents a section of the web with length 𝑙 before creping, shown in Figure 2.8(a). Next, we consider the moment that the section AB becomes debonded from the Yankee (Figure 2.8(b)). During this process the length of section AB is reduced by ∆𝑙, so the external work equals 𝐹𝑚𝑤∆𝑙. Since the web remains flat during the first stage, the delamination propagates in purely shear mode and there is no bending energy stored in the web. The surface energy has increased by 𝑙𝑤𝐺𝐼𝐼𝐶  and the strain energy stored in the web equals 𝑙𝑤ℎ 𝐸𝜀2 2⁄ . According to energy conservation, we have: 𝐹𝑚𝑤∆𝑙 =  𝑙𝑤𝐺𝐼𝐼𝐶 +𝑙𝑤ℎ𝐸𝜀22                                                      (2.19) where 𝜀 = 𝐹𝑚 𝐸ℎ⁄  is the strain in the section AB and ∆𝑙 = 𝑙𝜀. Substituting the expressions of 𝜀 and ∆𝑙 into Eq. (2.19) we obtain: 𝐹𝑚 = √2𝐸ℎ𝐺𝐼𝐼𝐶                                                           (2.20) 32  This creping force is purely due to the interfacial debonding, and is related to the mode II fracture energy.   Figure 2.8 A schematic of the web being pushed against the blade. (a) When the web first touches the blade. (b) When the section AB becomes delaminated from the substrate. The web is subjected to the maximum creping force 𝑭𝒎 during this process.          Next, we study the effect of the inertia on the maximum creping force. Imagining that the web in Figure 2.8 is free from the substrate and moving with an initial speed 𝑉𝑖𝑛, the web  turns in direction  and slows down  upon contact with the blade,  resulting  in a decrease in  kinetic energy. The upper limit of the contribution of this decrease in kinetic energy to web strain energy can be estimated by considering the  incoming kinetic energy  alone  as shown  below: 𝜌𝑙𝑤ℎ   22=𝑙𝑤ℎ𝐸𝜀22                                                         (2.21) Eq. (2.21) gives the expression of strain as 𝜀 = √𝜌 𝐸⁄ 𝑉𝑖𝑛. Then the maximum creping force is obtained: 𝐹𝑚 = ℎ𝑉𝑖𝑛√𝜌𝐸                                                          (2.22) Combining Eqs. (2.20) and (2.22), we can write the total maximum creping force as: 𝐹𝑚 = √2𝐸ℎ𝐺𝐼𝐼𝐶 +  ℎ𝑉𝑖𝑛√𝜌𝐸                                               (2.23) A BA B𝑙𝑙 ∆𝑙𝐹𝑚𝐹𝑚𝐹𝑚𝐹𝑚(a)(b)33  Here   is a coupling coefficient accounting for the viscous dissipation as well as the interfacial delamination, which can be determined by simulation. This is because in reality, the interfacial delamination and the decrease of kinetic energy are coupled. We need this coupling coefficient to combine these two terms together. Changing the viscosity of the web will result in the change of the coupling coefficient, thus the dissipation is also accounted for. Eq. (2.23) shows that 𝐹𝑚 in the first stage is linear to √𝐸 and is related to 𝐺𝐼𝐼𝐶. The kinetic energy term is quite small compared to the adhesion term for parameters with industrial values (about 10% with 𝑉𝑖𝑛 = 1200𝑚 𝑚  ), thus Eq. (2.23) is a reasonable estimation of the total maximum creping force. Nevertheless, Eq. (2.20) provides the exact solution of the maximum creping force for the quasi-static case.         Figure 2.9(a) plots 𝐹𝑚 as a function of the creping velocity 𝑉𝑖𝑛. It shows that 𝐹𝑚 increases linearly with 𝑉𝑖𝑛 and the simulation results match well with the analytical solution provided that  = 0.7  for the simulation value of 𝜇 = 20 𝑘𝑔 𝑚𝑠⁄ . Under relatively large velocity, the contribution of the inertia to the creping force becomes important. This indicates the lab-scale low speed creping experiments may not be representative of the actual creping process and the creping force at industrial speed. It is worth mentioning that in a previous experimental work [9], the creping force was also reported to increase as the creping speed increases. However, the speed effect is not significant as their maximum speed is limited to 250 𝑚 𝑚  ⁄ . Figure 2.9(b) illustrates the effect of Young’s modulus on the maximum creping force. Larger Young’s modulus results in higher 𝐹𝑚. Here with 𝜇 = 20 𝑘𝑔 𝑚𝑠⁄ ,  = 0. 8 was found to yield the best agreement between the simulations and the theory. In our simulation, it is found that decreasing 𝜇  results in the increase of  . This is because less kinetic energy is dissipated during the impact of the web against the doctor blade. 34   Figure 2.9 (a) The maximum creping force versus the creping velocity. The parameters for the simulation are chosen as: 𝒉 = 𝟓𝟎 𝝁𝒎 , 𝑬 = 𝟏𝟎𝟎 𝑴𝑷𝒂 ,𝝆 = 𝟑𝟎𝟎𝒌𝒈 𝒎𝟑⁄ , 𝝁 = 𝟐𝟎 𝒌𝒈 𝒎𝒔⁄ , 𝑮𝑰𝑪 = 𝑮𝑰𝑰𝑪 = 𝟐𝟎𝟎 𝑵 𝒎 . (b) The maximum creping force versus the Young’s modulus of the web. The parameters are chosen as: 𝒉 = 𝟓𝟎 𝝁𝒎, 𝝆 = 𝟑𝟎𝟎𝒌𝒈 𝒎𝟑⁄ , 𝝁 = 𝟐𝟎 𝒌𝒈 𝒎𝒔⁄ , 𝑽𝒊𝒏 = 𝟓𝒎 𝒔, 𝑮𝑰𝑪 = 𝑮𝑰𝑰𝑪 = 𝟐𝟎𝟎 𝑵 𝒎.  2.3.4 Periodic Folding         In the previous section we have studied the formation of a single fold. This folding process repeats many times during creping. To determine the steady-state creping wavelength and creping amplitude, we need to simulate multiple folds. Once the web is scraped off from the Yankee surface, it is pulled out by a rotating reel at the speed 𝑉𝑜𝑢𝑡. We can define an important parameter called the creping ratio based on the incoming and the outgoing web speeds: (𝑉𝑖𝑛 − 𝑉𝑜𝑢𝑡) 𝑉𝑖𝑛⁄ . To mimic this creping ratio effect, in the simulation we apply a constant speed 𝑉𝑜𝑢𝑡 to the first particle along the blade surface direction once it is debonded from the Yankee. Figure 2.10 shows the evolution of the web morphology when the creping ratio effect is included. The parameters in the simulation are chosen as: ℎ = 50 𝜇𝑚 , 𝑤 = 1𝑚 , 𝐸 = 100 𝑃𝑎 , 𝜌 = 300 𝑘𝑔 𝑚3⁄ , 𝜇 =20 𝑘𝑔 𝑚𝑠⁄ , 𝑉𝑖𝑛 = 20 𝑚 𝑠 , 𝑉𝑜𝑢𝑡 = 12 𝑚 𝑠  (creping ratio is 40 ), 𝛿 = 80° , 𝐺𝐼𝐶 = 𝐺𝐼𝐼𝐶 =100 𝑁 𝑚. Figure 2.10(a) and (b) show the formation of the first fold. The front of the fold slides along the blade surface. Here only the tip of the blade is illustrated here and the length of the edge (b)(a)35  which contacts the paper is set to be 3.5 mm. When the second fold is formed, it contacts the first fold and pushes the first fold out, as shown in Figure 2.10(c). As the web keeps moving towards the blade, the second fold becomes fully developed and the third fold starts to form, shown in Figure 2.10(d). Running the simulation for longer time shows that this folding process repeats, as shown in Figure 2.10(e) and Figure 2.10(f). After the whole simulation is completed, we perform a Fast Fourier Transformation on the final folding pattern to find the dominant creping wavelength and amplitude of the creped tissue paper. In the following sections, we study the effects of various parameters on the creping process.  Figure 2.10 Evolution of the web shape during the creping process. For simplicity, the adhesive layer and the Yankee are not shown. (a) 𝒕 = 𝟏𝟎𝝁𝒔, (b) 𝒕 = 𝟐𝟔𝝁𝒔, (c) 𝒕 = 𝟑𝟔𝝁𝒔, (d) 𝒕 = 𝟔𝟔𝝁𝒔, (e) 𝒕 = 𝟏𝟐𝟎𝝁𝒔, (f) 𝒕 = 𝟐𝟎𝟎𝝁𝒔.  (a) (b) (c)(d) (f)(e)𝑉𝑖𝑛𝑉𝑜𝑢𝑡36  2.3.5 Creping Angle Effect         Figure 2.11 plots the evolution of the creping force and the delamination length with three different creping angles. During the formation of the first fold, the creping forces are almost independent of the creping angle (Figure 2.11(a)). This independence is because initially only the front end of the web is in contact with the blade. However as the creping process continues, other parts of the web begin to contact the blade, so the creping forces increase for all different angles. For smaller creping angle the space is more constrained, so the contact area is larger and the creping force is higher, as is seen in the second and third peaks in Figure 2.11(a). Figure 2.11(b) shows that the delamination propagates faster with smaller creping angles and the delamination length is larger at a given time.   Figure 2.11 The evolution of the creping force (a) and the delamination length (b) under three different creping angles. The parameters for the simulation are chosen as: 𝒉 = 𝟓𝟎 𝝁𝒎, 𝑬 = 𝟏𝟎𝟎 𝑴𝑷𝒂,𝝆 = 𝟑𝟗𝟎𝒌𝒈 𝒎𝟑⁄ , 𝝁 =𝟐𝟎 𝒌𝒈 𝒎𝒔⁄ , 𝑮𝑰𝑪 = 𝑮𝑰𝑰𝑪 = 𝟏𝟎𝟎 𝑵 𝒎, 𝑽𝒊𝒏 = 𝟐𝟎 𝒎 𝒔, 𝑽𝒐𝒖𝒕 = 𝟏𝟐 𝒎 𝒔.          Figure 2.12 plots the averaged creping force and the final creping wavelength versus the creping angle. Both 𝐹?̅? and 𝜆 decrease when the creping angle increases. But when 𝛿 is beyond  0°, 𝐹?̅? and 𝜆 gradually reach a plateau. These results agree very well with the trend reported in previous experimental work [9, 12]. It is worth mentioning that the creping angle data from the (a) (b)37  lab-scale experiments [9] did not have creping ratio effect. So a quantitative comparison is not available.   Figure 2.12 The averaged creping force (a) and the creping wavelength (b) versus the creping angle. The parameters are chosen as the same in Figure 2.11.  2.3.6 Mixed-Mode Fracture Effect         During creping process, a proper adhesion between the web and the Yankee surface is crucial for generating quality tissue products. If the adhesion is too low, the web may simply be peeled off without any crepes formed. However if the adhesion is too high, it can exceed the fiber web cohesion and “picking” (i.e., individual fibers are pulled from the web) occurs, which significantly decreases the softness of the tissue. Figure 2.13(a) and (b) show the effect of the mixed-mode fracture energy on the averaged creping force and the creping wavelength. The three curves represent three different mixed-mode ratios, i.e., the ratio between 𝐺𝐼𝐼𝐶 and 𝐺𝐼𝐶. For each curve, it indicates that increasing fracture energy results in an increase of 𝐹?̅? and a decrease of 𝜆. This is because with larger fracture energy, the creping force required to break the interface becomes larger. According to Euler’s buckling theory, the critical delamination length at which the web buckles is shorter for a larger compressive force, resulting in a smaller creping wavelength. (a) (b)38  Comparing the three curves in each figure, it is apparent that 𝐹?̅? decreases and 𝜆 increases as the mixed-mode ratio is increased. This is because for a fixed 𝐺𝐼𝐼𝐶 increasing mixed-mode ratio means decreasing 𝐺𝐼𝐶, and during the third stage of folding the delamination can propagate easier and longer with a lower 𝐺𝐼𝐶, resulting in a smaller 𝐹?̅? and larger 𝜆. Figure 2.13(a) also indicates that the effects of mixed-mode ratio on 𝐹?̅? is more significant under a larger 𝐺𝐼𝐼𝐶. Our results show that in order to achieve finer crepes, both 𝐺𝐼𝐶 and 𝐺𝐼𝐼𝐶 need to be reasonably large. These results are consistent with the previous experiment [9]. However, in this experiment the fracture energy between the web and the substrate is controlled by varying the concentration of the adhesives. Thus, the effects of mixed-mode ratio was not discussed.    Figure 2.13 The effects of fracture energy on (a) the average creping force and (b) the creping wavelength. The parameters for the simulation are: 𝒉 = 𝟓𝟎 𝝁𝒎 , 𝑬 = 𝟏𝟎𝟎 𝑴𝑷𝒂 , 𝝆 = 𝟑𝟗𝟎𝒌𝒈 𝒎𝟑⁄ , 𝝁 = 𝟐𝟎 𝒌𝒈 𝒎𝒔⁄ , 𝑽𝒊𝒏 =𝟐𝟎 𝒎 𝒔, 𝑽𝒐𝒖𝒕 = 𝟏𝟐 𝒎 𝒔.  2.3.7 Evaluation of the Relative Impact of Parameters          In the previous sections we have studied the effect of the creping angle and fracture energy on creping. All the other parameters can be analyzed in a similar way. The impact of all the control (a) (b)39  parameters on the averaged creping force, creping wavelength and creping amplitude is summarized in Table 2.1.  Parameters ℎ 𝛿 𝐸 (𝑉𝑖𝑛 − 𝑉𝑜𝑢𝑡) 𝑉𝑖𝑛⁄  𝜌 𝐺𝐼𝐶,𝐺𝐼𝐼𝐶 𝑉𝑖𝑛 𝜆 197% -132% 23.0% -20% 11.1% -16.6% -3.64% A 199% -144% 23.4% 44.1% 11.0% -15.7% 0.30% 𝐹?̅? 59.1% -79.9% 3.94% -16.4% 32.5% 27.8% 10.9% Table 2.1 Impact of the control parameters on creping.  The numbers in Table 2.1 represent the relative percentage change of 𝜆, A and 𝐹?̅? when the corresponding parameter is increased from the minimum value to the maximum value shown in Table 1.1. For example, when the Young’s modulus 𝐸 is increased from 100 MPa to 600 MPa, the percentage change in 𝐸 equals 500%. We run the simulations for both values of 𝐸 and found the creping wavelength 𝜆 changes from 0.512 𝑚𝑚 to 1.1 𝑚𝑚. So the percentage change of 𝜆 equals (1.1 − 0.512) 0.512 = 115 ⁄ . Then the final impact number is given by 115 500 = 23 ⁄ . In Table 2.1, the impact of the parameters decreases from left to right. We can see that the thickness of the web ℎ and the creping angle 𝛿 have the highest impact on 𝜆, A and 𝐹?̅?. In previous creping modeling work done by Ramasubramanian et al. [15], the bending stiffness of the web (while keeping the Young’s modulus constant) was found to have the highest impact on the creping wavelength, which is consistent with our findings about the thickness. However, the effect of creping angle on the creping force was not correctly captured in their analytical model as they only consider the initial buckling stage of creping. It was reported that the creping force increases as the creping angle increase, which is against our conclusion (Figure 2.12(a)) and other experimental studies [9, 12]. Table 2.1 also indicates that the Young’s modulus 𝐸, density 𝜌, fracture energy 40  and creping ratio have a similar but lower level of impact. Here mode I and mode II fracture energies are assumed identical. Finally, the Yankee surface speed 𝑉𝑖𝑛 in the range studied in this paper has little effect on the final creping structure, but it has an obvious impact on the creping force. For all but one parameter, varying the parameter causes the creping wavelength 𝜆 and the amplitude A to increase or decrease together. The sole exception is that when the creping ratio is increased, 𝜆 decreases but A increases. For tissue paper, larger amplitude crepes result in higher bulk and smaller wavelength crepes give higher softness. Thus increasing the creping ratio is a good option to increase the bulk and the softness at the same time.  In the real tissue paper machine, the properties of the web (𝐸, 𝜌 and ℎ) are dominated by the selected pulp. The creping angle and the creping ratio can be easily adjusted on the machine. Adhesion can be changed by controlling the volume and composition of the sprayed adhesive chemicals. Our results suggest that the creping angle or the pulp properties should be adjusted to achieve a significant change of the web structure. If smaller changes are required, then the adhesive chemicals and the creping ratio can be modified accordingly.  2.4 Summary and Conclusions         In this Chapter, a one-dimensional particle dynamics model has been presented to study the creping process in shaped-bulk regime. As an initial model, the web is simplified as a viscoelastic beam and modeled as a series of discrete particles connected by viscoelastic elements. The failure of the adhesive layer is described by a mixed-mode discrete cohesive zone model. The self-contact of the web is incorporated based on a penalty model. The creping ratio effect is first incorporated by applying velocity boundary condition to the creped web. The main findings are summarized as follow: 41  1. The formation of micro-fold in creping includes three typical stages: delamination propagation, buckling and post-buckling deformation. We first report that during the third stage, the post-buckling of the web can further trigger the delamination to propagate, which leads to a larger creping wavelength.  2. An analytical solution of the maximum creping force in the first stage is derived based on an energy analysis. This solution matches well with the simulation results. It indicates that the creping force is the sum of a contribution from adhesion and a contribution from web inertia.  3. The effects of mixed-mode fracture energy on creping is first investigated. Increasing the mode-mixity while fixing 𝐺𝐼𝐼𝐶, leads to decrease of the average creping force and increase of the creping wavelength.  4. The relative impact of control parameters on the average creping force, wavelength and amplitude are evaluated by the model. The creping angle and the web thickness have the highest impact on creping. The Young’s modulus, density, fracture energy (𝐺𝐼𝐶  𝑎 𝑑 𝐺𝐼𝐼𝐶) and creping ratio have a similar but lower level of impact. The one-dimensional viscoelastic model is shown to be able to capture the folding mechanism in the shaped-bulk regime of creping process. The predicted effects of the control parameters qualitatively agree with the trend reported in practice. Extensions of the model will be presented in the following Chapters to consider other aspects involved in creping.   42  Chapter 3: Modeling the Effects of Plasticity and Explosive-Bulk in Creping 3.1 Introduction Most of the previous creping models, including the one presented in Chapter 2, assume that the web is purely elastic [15, 16, 80]. However, paper is essentially a visco-elasto-plastic material with an unknown residual strain, that clearly shows nonlinear behavior under loading [81]. The plasticity of paper originates from many factors such as the sliding and rotation of fibers, and the breakage of inter-fiber bonding. In the creping process, these mechanisms all can take place in the web when it hits the doctor blade at a high speed, which result into irreversible deformations. Some fibers are plastically deformed during large buckling deformation and a large amount of inter-fiber bonds are broken. In this Chapter, we extend the particle dynamics model presented in the previous Chapter to incorporate the plasticity of the web. We first use the extended model to investigate the effects of plasticity on creping. Then we perform a virtual tensile test on the creped tissue paper to find its stretch and stiffness. Finally, we explore the explosive-bulk regime by modeling the web as three individual layers with inter-layer bonding.   3.2 Extensions of Particle Dynamics Model for Plasticity         We need to modify the constitutive relation of the web to incorporate the plastic deformation, and change the axial forces and bending forces in the model presented in Chapter 2, accordingly.   3.2.1 Constitutive Model for the Uncreped Paper          Different constitutive models have been developed to account for the inelastic behavior of paper. For example, hyperelastic model based on the strain energy function [82], micromechanics model associated with representative volume elements [83], continuum damage mechanics (CDM) 43  model [81], fiber network models [54, 62, 84, 85] and continuum mechanics models [86, 87]. In this thesis, an elastoplastic material model associated with a kinematic hardening rule (Bauschinger effect) is used to describe the constitutive behavior of the web, shown as Figure 3.1. Elastic unloading is assumed here until the yield strength in the opposite direction is reached. This bilinear constitutive relation was shown to be adequate to capture the fiber characteristics during the unloading of paper [88] and the kinematic hardening rule was proven to apply in fibrous materials by a network simulation [84]. Although the paper properties in tension and compression are generally different. For simplicity, the constitutive relation is assumed to be identical in tension and compression. 𝐸1 and 𝐸2 represent the moduli before and after yielding. 𝜀𝑦 and 𝜎𝑦 are the yield strain and yield stress.   Figure 3.1 Bilinear elastoplastic material model associated with kinematic hardening rule (Bauschinger effect).           In the elastic case, the axial force 𝒇𝑎,𝑖 applied on the particle i only depends on its distance from the adjacent particles, as shown in Eq. (2.1). When considering the plasticity, this axial force depends on both the loading history and the current configuration of the web, thus needs to be modified. The axial strain between two adjacent particles  − 1  and   is given by 𝜀𝑖 1,𝑖 =(|𝒓𝑖 1 − 𝒓𝑖| − 𝑎0) 𝑎0, where 𝒓𝑖 is the position vector of particle  . The corresponding axial stress 𝜎𝜀𝜎𝑦𝜀𝑦2𝜎𝑦𝐸1𝐸244  𝜎𝑖 1,𝑖 is updated at each time step based on the incremental strain ∆𝜀𝑖 1,𝑖  and the stress at the previous time step. Initially the web is assumed stress free without any residual stress. With all the axial stresses updated at each time step, we can write the axial elasto-plastic force applied on the particle   as: 𝒇𝑎,𝑖 = 𝜎𝑖 1,𝑖𝑤ℎ𝒆𝑖,𝑖 1 + 𝜎𝑖,𝑖+1𝑤ℎ𝒆𝑖,𝑖+1                                             (3.1) In this Chapter, we use the same notations as those in Chapter 2 unless specifically mentioned. In the next section, we modify the bending forces in the particle system.  3.2.2 Elasto-Plastic Bending of the Web         In the current discrete particle modeling framework, the bending of the web is modeled by introducing bending forces on the particles, which are related to the bending moment. For elastic bending, the bending moment only depends on the bending rigidity and the current local bending angle formed by three consecutive particles. However, for elasto-plastic bending the bending moment depends on the deformation history of the web since the web can be plastically deformed. In this case, the bending moment can be calculated by integrating the axial stresses along the thickness direction.  Figure 3.2(a)-(d) represent four types of local bending systems formed by three consecutive particles (i-1)−(i)−(i+1), where 𝑖 denotes the elasto-plastic bending moment. In order to find out 𝑖 at any given time, we need to know the axial stress distribution in the thickness direction. We imagine the web to consist of n sublayers and then we keep track of the axial stress in each sublayer, as shown in Figure 3.2(e). It is worth mentioning that these virtual sublayers are only introduced for calculating the local bending moment, and they are not representing the actual fiber layers of the web. The strain in the j-th sublayer contains the axial strain and the bending strain, given by: 45  𝜀𝑗 = 𝜀?̅? −1𝑅 𝑧𝑗                                                               (3.2) where 𝜀?̅? = (𝜀𝑖 1,𝑖 + 𝜀𝑖,𝑖+1) 2⁄  represents the averaged local axial strain, 𝑧𝑗 is the distance from the centerline of the sublayer  to the mid-plane and  1 𝑅𝑖⁄  is the curvature of the cross-section. With a total number of   sublayers, we have: 𝑧𝑗 = (𝑗 −1+𝑛2)ℎ𝑛                                                              (3.3) For a discrete system like Figure 3.2(a)-(d), the magnitude of the curvature is approximately given by: |1𝑅 | =2cos (𝛾 2)𝑎0(1+𝜀 )                                                                (3.4) where 𝛾𝑖 is the local bending angle. If the web bends upward, the curvature is set to be positive (Figure 3.2(a)-(b)), otherwise, it is negative (Figure 3.2(c)-(d)).  Figure 3.2 (a)-(d) Local bending systems formed by consecutive particles (i-1)−(i)−(i+1): (a) Positive moment and positive curvature; (b) Negative moment and positive curvature; (c) Positive moment and negative curvature; (d) Negative moment and negative curvature. (e) n sublayers for calculating the bending moment, where partial yielding is allowed in any given cross-section. (f) Schematic of the viscous bending forces due to angular velocities 𝝎𝟏,𝒊 and 𝝎𝟐,𝒊. (a)(e)𝛾𝑖 -1  +1𝒇𝑏1,𝑖 𝒇𝑏2,𝑖𝒇𝑏0,𝑖(b) 𝑖  𝑖 𝛾𝑖 -1  +1𝒇𝑏1,𝑖 𝒇𝑏2,𝑖𝒇𝑏0,𝑖 𝑖  𝑖𝛾𝑖 -1  +1𝒇𝑏1,𝑖 𝒇𝑏2,𝑖𝒇𝑏0,𝑖 𝑖  𝑖(c)𝛾𝑖 -1  +1𝒇𝑏1,𝑖 𝒇𝑏2,𝑖𝒇𝑏0,𝑖 𝑖  𝑖(d) -1   +1𝝎1,𝑖 𝝎2,𝑖𝒇𝑣1,𝑖 𝒇𝑣2,𝑖𝒇𝑣,𝑖(f)Z𝑧𝑗12  2 2⁄ + 1j − 1 46          At every time step, we calculate the strain in each sublayer based on Eqs. (3.2)-(3.4), and incrementally update the stress according to the constitutive relation described in Figure 3.1. By doing so, partial yielding is allowed for any given cross-section and the plastic effects in beam bending such as “spring back” can be captured. Once the stress in each sublayer is calculated, the total bending moment can be derived as:  𝑖 = ∑ (−𝜎𝑗𝑧𝑗𝑤ℎ𝑛)𝑛𝑗=1                                                         (3.5) where 𝜎𝑗  is corresponding to the stress in the j-th sublayer. Since the stresses in some of the sublayers may enter the plastic zone earlier, this approach can allow for partially yielded cross section. Finally, the bending forces can be calculated and applied on the particles, as shown in Figure 3.2(a)-(d). 𝒇𝑏1,𝑖  and 𝒇𝑏2,𝑖  are restoring forces applied on the left and right particles and perpendicular to the axial direction. The exact direction of the restoring force depends on the signs of curvature and moment, and its magnitude is related to 𝑖, given by: |𝒇𝑏1,𝑖| = |𝒇𝑏2,𝑖| =|𝑴 |𝑎0(1+𝜀 )                                                     (3.6) 𝒇𝑏0,𝑖 = −𝒇𝑏1,𝑖 − 𝒇𝑏2,𝑖  is applied on the central particle to maintain the overall equilibrium of bending forces.          The presented elasto-plastic bending model is validated by simulating the deformation of a cantilever beam and comparing the results with the Finite Element Method. The details are provided in the Appendix B. For all the simulations in this paper, the convergence of the model is checked by increasing the number of total particles 𝑁 and the total number of sublayers  . In general 𝑎0 needs to be smaller than h, and   needs to be at least 20 to ensure the discrete model converges to a stable value. 47          After updating the axial force and bending force, the model can account for the visco-elasto-plastic deformation of the web for both loading and unloading. In the following sections, we present the results obtained based on this extended visco-elasto-plastic model.  3.3 Effects of Plasticity on the Folding Process         In the previous Chapter, we have shown that the formation of individual fold typically consists of three stages, i.e., delamination propagation, bucking initiation and post-buckling deformation (Figure 2.6). It makes sense to see how the mechanisms may change by including the plastic deformation of the web. Figure 3.3 demonstrates the evolution of the creping force and the delamination length during the formation of individual fold. The parameters are chosen as: ℎ =30 𝜇𝑚 , 𝑤 = 1𝑚 , 𝐿0 = 12 𝑚𝑚 ,  𝑁 =  01 , Δ𝑡 = 1  𝑠 , 𝐸1 = 500 𝑃𝑎 , 𝜎𝑦 = 10 𝑃𝑎 , 𝜌 =500 𝑘𝑔 𝑚3⁄ , 𝑉𝑖𝑛 = 1200 𝑚 𝑚  , 𝑉𝑜𝑢𝑡 =   0 𝑚 𝑚   (creping ratio is 20%), 𝛿 =  0° , 𝐺𝐼𝐶 =𝐺𝐼𝐼𝐶 = 50 𝑁 𝑚. For the elastic case, 𝐸2 is set to be equal to 𝐸1. In the first stage, as the creping force increases to its peak value, delamination starts to propagate in mode II (shear mode). It shows that both the maximum creping force and the delamination propagating speed (i.e., the slope in Fig. 7b) decrease as 𝐸2 decreases. This is because for lower value of 𝐸2, the web undergoes more plastic deformation and the overall stiffness of the web has reduced. In the second stage, once the delamination length reaches the critical value the web starts to buckle. Meanwhile the compressive stress is released and the delamination stops propagating. In the third stage, post-buckling of the web induces the delamination to propagate in mixed-mode for the elastic case and for the plastic case when 𝐸2 = 0.5𝐸1. However, when 𝐸2 = 0.2𝐸1 the delamination length remains the same. This is because for a relatively low value of 𝐸2, most of the web is plastically deformed after buckling and the effective bending stiffness is significantly reduced, which decreases the driving 48  force for crack propagation and prevents further increase of the delamination length. Figure 3.4 shows the morphology of the simulated tissue paper under different plastic modulus, while keeping other parameters the same as those used in Figure 3.3. It indicates that the size of the fold decreases as the plastic modulus decreases. This is because if the web has more significant plastic deformation, both 𝐹𝑐  and 𝑙𝑑  are reduced, which leads to smaller creping wavelength and amplitude.  Figure 3.4 also shows that with lower plastic modulus, the folds are less uniform and the peaks of the folds are sharper.   Figure 3.3 Evolution of the (a) creping force and (b) delamination length during the formation of one fold under different values of 𝑬𝟐.   Figure 3.4 Surface pattern of generated tissue paper under different plastic modulus.  (b)(a)𝐸2 = 0.05𝐸1𝐸2 = 0.1𝐸1𝐸2 = 0.2𝐸1𝐿0 = 3𝑚𝑚𝐸2 = 0.05𝐸1𝐸2 = 0.1𝐸1𝐸2 = 0.2𝐸149          Based on the results shown in Figure 3.3 and Figure 3.4, we can conclude that the plastic deformation has shown significant effects on the folding process during creping. By including the plasticity, the crepe pattern becomes more realistic.    3.4 Analysis of Process-Structure-Property Relation by Virtual Tensile Test         During the creping process, significant plastic deformation occurs in the web and periodic crepes are formed. Both of them can potentially affect the macroscopic properties of the tissue paper. So far, we have investigated how various process parameters change the structure of the generated tissue paper, but the correlation between the structure and the property remains to be understood. In this section, we perform a virtual tensile test on the web after it is creped to find out the tensile properties of the tissue. Figure 3.5(a)-(c) are the snapshots of the web’s morphology during the creping process with the following parameters: ℎ = 30 𝜇𝑚 , 𝑤 = 1𝑚 , 𝐿0 = 4 𝑚𝑚,  𝑁 = 201 , Δ𝑡 = 1  𝑠 , 𝐸1 = 1 𝐺𝑃𝑎 ,  𝐸2 = 100 𝑃𝑎 , 𝜎𝑦 = 30 𝑃𝑎 , 𝜌 = 500 𝑘𝑔 𝑚3⁄ , 𝑉𝑖𝑛 =1200 𝑚 𝑚  , 𝑉𝑜𝑢𝑡 =   0 𝑚 𝑚   (creping ratio equals 20%), 𝛿 = 80°, 𝐺𝐼𝐶 = 𝐺𝐼𝐼𝐶 = 25 𝑁 𝑚. To mimic the creping ratio effect, we apply a constant velocity (𝑉𝑜𝑢𝑡) to the first particle in the web once it is detached from the Yankee.  During creping the self-contact of the web occurred as the generated new fold pushes the previous fold out along the blade surface, as shown in Figure 3.5(b)-(c). After the entire web is debonded from the Yankee, we stop applying 𝑉𝑜𝑢𝑡 on the first particle and let the generated tissue paper relax for a sufficient time to reach the equilibrium state. The top curve in Figure 3.5(d) represents the relaxed tissue paper.          Next, we take this simulated tissue paper and perform a displacement-controlled virtual tensile test by the extended particle dynamics model. It is worth noting that the tissue paper is already plastically deformed, and the loading history of the paper is taken into account in the 50  subsequent tensile test. We fix the left end of the paper and slowly pull the right end to perform a quasi-static uniaxial tensile test. As the sheet is stretched, the folds are gradually flattened, shown as the bottom curve in Figure 3.5(d). Since in our model the fracture (i.e., breakage) of the web is not considered, we end the tensile test when the tensile stress in the web reaches 30% of the yield stress 𝜎𝑦, considering the fact that the typical strength reduction of paper is around 70% after creping [8, 3]. The stretch is defined as the engineering strain at the end of the tensile test. It is worth noting that the sheet is not fully flattened even at the end of the tensile test, as shown in Figure 3.5(d). The fiber rotation and sliding occurred in a real tissue paper are implicitly accounted for in our one-dimensional model, through the elasto-plastic constitutive material model.  Figure 3.5 (a)-(c) The simulated morphology of the web when 𝒕 = 𝟏𝟎𝝁𝒔, 𝒕 = 𝟒𝟎𝝁𝒔 and 𝒕 = 𝟏𝟑𝟎𝝁𝒔. The Yankee and the adhesive layer are not shown here. (d) Top: final morphology of the web after the creping is completed; bottom: the stretched web after the uniaxial tensile test. (a)(c) (d)(b)𝑉𝑖𝑛𝑉𝑜𝑢𝑡𝑃𝑃51  3.5 Stress-Strain Curve of Tissue Paper         Figure 3.6 shows the stress-strain curve corresponding to the virtual tensile test in Figure 3.5(d). Before creping, the web follows the bilinear constitutive relation shown in Figure 3.1. After creping, the properties of the web changes completely and the stress-train curve indicates significant nonlinearity, including an inflection point close to 𝜀 = 7 . The initial modulus of the simulated tissue paper in Figure 3.6 is approximately 20 𝑃𝑎, which is much smaller compared to the 𝐸1 = 1 𝐺𝑃𝑎 and 𝐸2 = 100 𝑃𝑎 in the base sheet. This is due to the significant plastic deformation as well as the generated micro-folds during the creping process. The inflection point in the stress-strain curve is likely corresponding to the transition from fiber sliding and reorientation to fiber stretching, in a real tissue sample. Apart from the mild inflection point, in general, the slope of the curves increases with increasing strain. This is because initially the extension is mainly due to the unfolding of the folds. The load required to overcome the bending resistance in order to unfold the sheet is relatively low. As the sheet is gradually flattened, the extension is mainly caused by the stretching of the sheet, and thus the stiffness increases. Our simulated tensile stress-strain curve is very similar to the data reported in the previous literature [3, 81] and our experimental data (as shown later in Chapter 5, Figure 5.4), which suggests that the plasticity and folding mechanism captured in our model are essential aspects governing the macroscopic properties of tissue paper. 52   Figure 3.6 Stress-strain curve corresponding to the virtual tensile test in Figure 3.5(d).  3.6 Parametric Study         The parametric studies conducted by the previous creping models are based on the assumption that the web is purely elastic [15, 16], which may lead to inaccurate conclusions. The creping ratio effect has not been included in these models. Also, only the effects on the crepe structures (creping wavelength and amplitude) have been investigated, and how the parameters exactly affect the properties of tissue remains unclear. In this section, we use the extended elastoplastic model and the virtual tensile test to find out the effects of different control parameters on the crepe structure and the tissue properties such as modulus and stretch.    3.6.1 Creping Ratio Effect         We first investigate the creping ratio effect on tissue paper properties. In practice, the creping ratio is usually fixed for a particular tissue grade, by varying 𝑉𝑜𝑢𝑡 and 𝑉𝑖𝑛 simultaneously. In order to change the creping ratio, one way is changing 𝑉𝑜𝑢𝑡 while keeping 𝑉𝑖𝑛 constant. Figure 3.7(a) 53  shows the simulated tissue paper under three different creping ratios. The top green line represents the initial flat sheet before creping. The other parameters used for Figure 3.7 are: ℎ = 30 𝜇𝑚, 𝑤 =1𝑚 , 𝐿0 = 4 𝑚𝑚 ,  𝑁 = 201 , Δ𝑡 = 1  𝑠 , 𝐸1 = 1 𝐺𝑃𝑎 ,  𝐸2 = 100 𝑃𝑎 , 𝜎𝑦 = 30 𝑃𝑎 , 𝜌 =500 𝑘𝑔 𝑚3⁄ , 𝑉𝑖𝑛 = 1200 𝑚 𝑚  , 𝛿 =  0°  and 𝐺𝐼𝐶 = 𝐺𝐼𝐼𝐶 = 25 𝑁 𝑚 . Figure 3.7(a) indicates that by increasing the creping ratio during the creping process, the creping amplitude is increased but the wavelength of the creped sheet is decreased. This is because for higher creping ratio, the creped sheet is pulled out with lower 𝑉𝑜𝑢𝑡, resulting in a less stretched fold with smaller 𝜆 but larger 𝐴. Figure 3.7(b) plots the corresponding tensile stress-strain curves obtained by performing virtual tensile tests on the generated creped sheet shown in Figure 3.7(a). The tensile stress-strain curve for the uncreped base sheet follows the bilinear constitutive relation described by Figure 3.1, which is not shown here. All the curves in Figure 3.7(b) show significant nonlinearity, including the inflection point. Figure 3.7(b) shows that the stretch (i.e., the final strain) is increased as the creping ratio increases, which agrees with the experiments reported in [14]. This is because for a higher creping ratio, the creped sheet is less stretched during the creping process. Then in the tensile test, the sheet will experience more extension to “unfold” the folds in the tensile test. Also note that the initial stiffness is lower with a higher creping ratio, which potentially increases the bulk softness of the end product, and the caliper (thickness) to basis weight ratio. Also note that the ratio between the creping amplitude and the creping wavelength (𝐴 𝜆⁄ ) increases as the creping ratio increases. Higher 𝐴 𝜆⁄  indicates that there is more material “packed” per unit length in the machine direction, thus the tissue can be stretched longer before it breaks. 54   Figure 3.7 (a) The uncreped base sheet and the simulated creped sheet under different creping ratios. (b) The corresponding tensile stress-strain curve of the creped sheet.   3.6.2 Adhesion Effect         Figure 3.8 illustrates the effects of the interfacial adhesion (fracture energy). The parameters used are: ℎ = 30 𝜇𝑚 , 𝑤 = 1𝑚 , 𝐿0 = 4 𝑚𝑚 ,  𝑁 = 201 , Δ𝑡 = 1  𝑠 , 𝐸1 = 1 𝐺𝑃𝑎 , 𝐸2 =100 𝑃𝑎 , 𝜎𝑦 = 30 𝑃𝑎 , 𝜌 = 500 𝑘𝑔 𝑚3⁄ , 𝑉𝑖𝑛 = 1200 𝑚 𝑚  , 𝑉𝑜𝑢𝑡 =   0 𝑚 𝑚   and 𝛿 = 0°. For a low ratio of 𝐸2 𝐸1⁄ , Figure 3.3 shows that in the third stage of folding the post-buckling deformation does not trigger the delamination to propagate any further. Under this condition, 𝐺𝐼𝑐 is much less important compared to 𝐺𝐼𝐼𝑐. So in Figure 3.8 we simply assume 𝐺𝐼𝑐 = 𝐺𝐼𝐼𝑐 = 𝐺𝑐. Figure 3.8(a) indicates that as the fracture energy increases both 𝐴  and 𝜆  decrease, which is consistent with the predictions by the elastic model in Chapter 2 and the previous experimental results in [9]. Figure 3.8(b) further shows that the stretch decreases and the initial stiffness increases as 𝐺𝑐 increases. Since lower stiffness results in higher softness of the tissue paper, our results suggest that finer crepes (i.e., smaller 𝐴 and 𝜆) does not necessarily mean higher softness. In fact, by calculating the ratio between 𝐴 and 𝜆 in Figure 3.7 and Figure 3.8, we observe that 𝐴 = 4 𝜇𝑚, 𝜆 = 3 5𝜇𝑚𝐴 𝜆⁄ = 0.13(a) (b)𝐴 = 54𝜇𝑚, 𝜆 = 345𝜇𝑚𝐴 𝜆⁄ = 0.1 𝐴 =  3𝜇𝑚, 𝜆 = 310𝜇𝑚𝐴 𝜆⁄ = 0.2055  larger 𝐴 𝜆⁄  corresponds to lower stiffness (higher bulk softness) and higher stretch. This can be understood as follows. A larger 𝐴 𝜆⁄  means the web is more “compressed” and undergoes more plastic deformation during creping, so after creping the web can be stretched longer with lower resistance.  Figure 3.8 (a) The uncreped base sheet and the simulated creped sheet with different fracture energies. (b) The corresponding tensile stress-strain curve of the creped sheet. 𝑮𝑰𝒄 and 𝑮𝑰𝑰𝒄 are assumed to be identical.  3.6.3 Relative Impact of Parameters         We investigate the effects of the other parameters in addition to the creping ratio and fracture energy by using the extended model. Table 3.1 summarizes the impact of the key parameters considered in this paper. The impact factors in the Table indicate the relative percentage change of 𝜆, A, 𝐹?̅? and stretch when the value of a certain parameter is increased from minima to maxima within its typical range. Here the basis weight is used to represent the base sheet property. The thickness and the density are related to the basis weight according to two empirical equations: ℎ =( + 𝐵𝑊)  ⁄  and 𝜌 =  ∙ 𝐵𝑊 ( + 𝐵𝑊)⁄ , where   and   are two constants [89]. In this table 𝐸2 is fixed as 0.1𝐸1 and 𝐺𝐼𝐶 is assumed to be equal to 𝐺𝐼𝐼𝐶. Table 3.1 shows that both the basis weight (a) (b)𝐴 = 54𝜇𝑚, 𝜆 = 345𝜇𝑚𝐴 𝜆⁄ = 0.1 𝐴 = 40𝜇𝑚, 𝜆 = 301𝜇𝑚𝐴 𝜆⁄ = 0.13𝐴 = 35𝜇𝑚, 𝜆 = 2 1𝜇𝑚𝐴 𝜆⁄ = 0.1256  and the creping angle have the highest impact on 𝜆, A and 𝐹?̅?, which is consistent with the trend predicted by the elastic model in Chapter 2 (Table 2.1). However, for the stretch the creping angle and the creping ratio exhibit a higher impact compared to the rest parameters.  Parameters Basis Weight 𝛿 (𝑉𝑖𝑛 − 𝑉𝑜𝑢𝑡) 𝑉𝑖𝑛⁄  𝐸1 𝐺𝐼𝐶, 𝐺𝐼𝐼𝐶 Typical Range 15 25 𝑔 𝑚2 75°  5° 15  30  100-1000 MPa 25 100 𝑁 𝑚 𝜆 85.2% -142% -8.4% 23.1% -13.5% 𝐴 78.8% -140% 50.9% 34.3% -11.9% 𝐹?̅? 31.5% -122% 2.3% 13.4% 69.5% Stretch 53.1% 280% 262% 142% -101% Table 3.1 Impact of key parameters on creping predicted by the visco-elasto-plastic model.  3.7 Analysis of the Explosive-Bulk Regime 3.7.1 Multi-Layer Model         “Sheet explosion” occurred in creping considerably increases the bulk of the tissue paper [8]. The current model assumes the web is a single-layer and hence can only capture the shaped-bulk regime. To explore the explosive-bulk regime, we model the web as three individual layers with the same effective properties, as shown in Figure 3.9. This simplification is justified because a typical tissue paper contains three to five layers of fibers. Both the top and the bottom layers are bonded to the middle layer at the discrete points, but only the bottom layer is bonded to the Yankee. The bonding between two layers is assumed to be elastic and it fails once the separation 𝛿𝑖𝑛𝑡 exceeds a critical distance, 𝛿𝑖𝑛𝑡𝑐 . The initial separation between the layers (𝛿𝑖𝑛𝑡0 ) is set to be equal to h and the bonding stiffness is represented by 𝑘𝑏𝑜𝑛𝑑. The bonding force 𝑓𝑏𝑜𝑛𝑑 is given as: 57  𝑓𝑏𝑜𝑛𝑑 = {𝑘𝑏𝑜𝑛𝑑(𝛿𝑖𝑛𝑡 − 𝛿𝑖𝑛𝑡0 ),      𝛿𝑖𝑛𝑡0 ≤ 𝛿𝑖𝑛𝑡 ≤ 𝛿𝑖𝑛𝑡𝑐               0,                   𝛿𝑖𝑛𝑡 > 𝛿𝑖𝑛𝑡𝑐                                      (3.7) The contact between two adjacent layers is accounted for by the same penalty method described in section 2.2.4. If the bonding stiffness between the layers is strong, they tend to remain bonded to each other as the web buckles. However, for weaker bonding stiffness inter-layer debonding may occur and leads to “sheet explosion”. Here we have assumed that the inter-fiber bonds are uniformly distributed with the same strength. In reality, the inter-fiber bonds are quite stochastic. We leave this aspect for future investigation. Also note that the in-plane tensile properties of each layer is affected only by the plastic deformation occurred during creping. In reality, each individual layer can break, which significantly reduces the strength.  Figure 3.9 Schematic of the three-layer model. The stiffness of inter-layer bond is represented by 𝒌𝒃𝒐𝒏𝒅 and the axial stiffness between two adjacent particles in each layer is represented by 𝒌𝒂.          In theory, the inter-fiber bonding can also be modeled by cohesive zone model, i.e., the bonding force first increases and then gradually decreases to zero as the separation keeps increasing. However, based on the AFM measurements of fiber–fiber contact, the adhesive force drops to zero rapidly once a critical separation is reached, as shown in [62]. So, it is reasonable to 𝑘𝑎𝑘𝑏𝑜𝑛𝑑58  use the failure model in Eq. (3.7). Nevertheless, the inter-fiber bonding behavior is quite significant to the moisture, temperature and other effects. It is possible that in some particular cases, a cohesive zone model is more suitable to describe the inter-fiber bonding.  3.7.2 Folding Process in Explosive-Bulk Regime         Figure 3.10 shows the simulated morphologies of the web at different time steps during creping by the multi-layer model, with the following parameters: ℎ = 25 𝜇𝑚, 𝑤 = 1𝑚 , 𝐿0 =  𝑚𝑚 ,  𝑁 =  01 , Δ𝑡 = 1  𝑠 , 𝐸1 = 0.8 𝐺𝑃𝑎 , 𝐸2 = 80 𝑃𝑎 , 𝜎𝑦 = 24 𝑃𝑎 , 𝜌 = 300 𝑘𝑔 𝑚3⁄ , 𝑉𝑖𝑛 = 1200 𝑚 𝑚  , 𝑉𝑜𝑢𝑡 =  00 𝑚 𝑚  , 𝐺𝐼𝐶 = 𝐺𝐼𝐼𝐶 = 50𝑁 𝑚 , 𝑘𝑏𝑜𝑛𝑑 = 0.1𝑘𝑎 , 𝛿𝑖𝑛𝑡𝑐 = 0.25ℎ and 𝛿 = 85°. Here 𝑘𝑎 = 𝐸1𝑤ℎ 𝑎0 represents the axial stiffness. These parameters are estimated based on the data reported in previous literature [8, 9, 12, 16, 62]. When the web starts to buckle, the debonding between the layers start to initiate at the peak of the fold, as demonstrated in Figure 3.10(a). After that, the separation grows larger during the post-buckling of the web. However, the three layers are still partially bonded, typically at the bottom of each fold where the web “sits” on the blade, as shown in Figure 3.10(b)-(d).  59   Figure 3.10(a)-(d) The simulated morphologies of the web at 𝒕 = 𝟖𝝁𝒔, 𝒕 = 𝟖𝟎𝝁𝒔, 𝒕 = 𝟏𝟐𝟎𝝁𝒔 and 𝒕 = 𝟐𝟏𝟎𝝁𝒔 by three-layer model.   3.7.3 Phase Diagram of Creping Regime         In creping, more “sheet-explosion” is thought to increase the overall softness and bulk of the tissue paper. Premium-grade products like facial tissue need more explosive-bulk, while low-grade products like bathroom tissue can have more shaped-bulk. Thus, it is useful to construct a phase diagram to predict the creping regime based on the control parameters. Figure 3.11 is the phase diagram constructed by using the multi-layer model. We vary the inter-layer bonding stiffness and the web-Yankee interfacial fracture energy, while keeping the other parameters fixed as: ℎ =25 𝜇𝑚 , 𝑤 = 1𝑚 , 𝐿0 = 12 𝑚𝑚 , 𝛿𝑖𝑛𝑡𝑐 = 0.25ℎ , 𝑁 =  01 , Δ𝑡 = 1  𝑠 , 𝐸1 = 2 𝐺𝑃𝑎 , 𝐸2 =200 𝑃𝑎 , 𝜎𝑦 = 40 𝑃𝑎 , 𝜌 = 300 𝑘𝑔 𝑚3⁄ , 𝑉𝑖𝑛 = 1200 𝑚 𝑚  , 𝑉𝑜𝑢𝑡 =   0 𝑚 𝑚   and 𝛿 =85°. Here 𝐺𝐼𝐶  and 𝐺𝐼𝐼𝐶  are assumed equal: 𝐺𝐼𝐶 = 𝐺𝐼𝐼𝐶 = 𝐺𝐶 . 𝑘𝑎 = 𝐸1𝑤ℎ 𝑎0 and 𝐺0 = 𝐸1𝑤ℎ𝑎0 (a) (b)(c) (d)𝑉𝑖𝑛𝑉𝑜𝑢𝑡60  are used to normalized 𝑘𝑏𝑜𝑛𝑑  and 𝐺𝐶 . In practice, the adhesion is related to the amount of adhesives, drying, temperature and pressing. The bonding stiffness is mainly governed by the pressing, drying and “ingredients” used in furnish.   Figure 3.11 Phase diagram of the creping regime. The inter-layer bonding stiffness 𝒌𝒃𝒐𝒏𝒅  and the fracture energy 𝑮𝒄 are varied while keeping the rest parameters fixed.          In the simulation when the number of failed inter-layer bonds reaches 25% of the total bonds, we observe obvious “sheet explosion”. We use this criterion to separate the shaped-bulk regime and the explosive-bulk regime in the phase diagram. Figure 3.11 indicates that the explosive-bulk regime is associated with high fracture energy and low bonding stiffness. This is because when the fracture energy is higher, the folds formed during creping becomes smaller. This causes the curvature at the peak of each fold to increase and the bonding stress between two layers increases as well. Thus, high fracture energy makes the inter-layer bond at the peak of each fold more likely to fail during creping, shown as the inserted figure in Figure 3.11. It is thought that explosive-bulk is required to achieve higher softness tissue products. Thus, it seems plausible that in order to increase the softness, one might increase the amount of adhesive chemicals applied on the Yankee “Shaped bulk”“Explosive bulk”61  surface and reduce the pressure applied by the “pressure roll” on the web at the beginning of creping.  3.8 Summary and Conclusions         In this Chapter, we extend the elastic particle dynamics model developed in Chapter 2 to take account of the plastic deformation of the web as well as the explosive-bulk in the creping process. A bilinear elastoplastic material model associated with a kinematic hardening rule (Bauschinger effect) is used to describe the effective constitutive behavior of the web. The axial force and the bending force are modified to include the plasticity effect. The main findings are summarized as follow: 1. Including the plasticity of the web reduces the creping force and slows down the delamination propagation during the first stage of fold formation. By decreasing the plastic modulus of the web, the creping wavelength and creping amplitude both decrease, and the micro-folds become non-uniform.  2. A virtual tensile test on the generated tissue paper is established to find the correlation between the crepe structure and the property. The virtual tensile stress-strain curves show significant nonlinearity. As the ratio between creping amplitude and creping wavelength increases, the initial stiffness of the tissue paper reduces and the stretch increases. 3. A systematic parametric study is performed based on the one-layer visco-elasto-plastic model. Basis weight and creping angle show relative highest impact on 𝜆, A and 𝐹?̅?. For the stretch, creping angle and creping ratio exhibit a higher impact compared to the rest parameters. 62  4. A multi-layer model is developed to account for the explosive-bulk regime. The constructed phase diagram indicate that explosive-bulk regime is more likely to occur with high adhesion and low bonding strength. Our simulation results indicate that plasticity of the web has significant impact on the creping process, thus it is important to account for that in a realistic creping model. The virtual tensile test enable us to fill the gap between the crepe structure and the tissue properties. By using the extended multi-layer model, the explosive-bulk regime of creping is investigated for the first time. In the next Chapter, we use this visco-elasto-plastic model to explore the effects of inhomogeneity on the creping process. 63  Chapter 4: Effects of Inhomogeneity on Creping  4.1 Introduction: Causes of Inhomogeneity         In Chapter 2 and 3, we have developed a particle dynamics model to simulate the creping process and investigate the effects of various parameters on the creping structure. For the elastic model, once the creping process reaches the steady state the generated micro-folds are uniform with almost identical size, as shown in Figure 2.10. For the plastic model, in the shaped-bulk regime the simulated micro-folds also have similar size, as shown in Figure 3.5. However, in reality the crepe folds in commercial tissue products are quite irregular. Figure 4.1 shows the crepe structure of typical commercial tissue paper. Figure 4.1(a) is the surface view of the tissue, which clearly indicates that the micro-folds are not uniform along the machine direction and not continuous along the cross direction. The typical length of an individual fold in the cross direction is close to the length of a fiber (1 3𝑚𝑚). Figure 4.1(b) shows the SEM cross-sectional image of the same tissue at a different region, which indicates that both the creping wavelength and amplitude are varying from fold to fold.  Figure 4.1 (a) Surface image of tissue paper. (b) Cross-sectional image of tissue paper.          ZDMD1 𝑚𝑚 (b)(a) 100 𝜇𝑚CDMD64          In the tissue manufacturing process, many factors can cause the non-uniformity in the tissue paper. A very common cause is the non-uniform adhesion between the web and the Yankee. The adhesive chemicals are spayed on the surface by a series of spray nozzles below the Yankee, as shown in Figure 1.1. If the flow of the adhesives is not held constant during production, it can cause inhomogeneous adhesion and pressing, which leads to irregular crepe pattern in the tissue paper. The defects on the Yankee surface can also cause the non-uniform adhesion distribution. Thus, the characteristic length for adhesion variation has a wide range. Another common cause is the poor formation. Paper formation refers to the uniformity and distribution of the fibers within a sheet. The main factors affecting formation are the quality of the wood pulp and “ingredients” used, and the speed of the moving fabric. A smaller speed gives the wood fibers more time to line up evenly on the fabric. Poor formation leads to the non-uniform distribution of the basis weight, which is an important parameter governing the size of the crepe based on our simulation results in Chapter 3. The typical length scale of basis weight variation is very wide, starting from the fiber-width scale, in terms of its power spectrum, but a visual formation is normally judged in the length scale of a few millimeters or more [90, 91, 92]. Lastly, the forming fabric pattern is also an important cause of inhomogeneity. In a typical paper machine, the web is created by deposition of the fibers on the forming fabric and the subsequent consolidation and dewatering. During these processes, the pattern of the forming fabric is embedded in the web, meaning that microstructures are already formed before the web enters into the creping section. The typical forming fabric used in tissue making has a variation length of 200-500 𝜇𝑚. How this embedded structure in the base sheet affects the creping process remains unclear.          In this Chapter, we extend the particle dynamics model to investigate the effects of forming fabric pattern, non-uniform adhesion and basis weight distribution on the creping structure and the 65  tissue properties. Although it is possible that all these causes of inhomogeneity can be coupled and arise at the same time during the tissue manufacturing process, here we investigate them separately in order to differentiate the effect of each case.  4.2 Effects of Forming Fabric Pattern         The main functions of forming fabrics in papermaking are to permit water removed from the sheet to flow through the fabric; support, retain, and form the sheet; and to convey the sheet from the headbox to the press section. Previous studies have shown that forming fabric creates embedded sheet structures that can affect the properties of end products [93]. Figure 4.2 is the surface view of a typical forming fabric used in tissue manufacturing process, which shows an intrinsic periodic pattern. The fibers are deposited onto the top surface of the forming fabric to form the sheet. Due to the 3-dimensional geometry of the forming fabric, the sheet is embedded with an uneven surface pattern, i.e., the “wire mark”.   Figure 4.2 A schematics of forming fabric (surface view) used in tissue manufacturing.  1 𝑚𝑚66          In the previous simulations, we modeled the initial web as flat. Here, in order to account for the effects of forming fabric, we assume that the web has a sinusoidal pattern along the machine direction. Although the exact pattern on the web varies depending on the type of fabric used, the sinusoidal form can be used as a general case. Figure 4.3 shows the schematic of an initial web pattern caused by the forming fabric. 𝐴0  and 𝜆0  are used to denote the initial amplitude and wavelength. Note that 𝐴0 is typically much smaller compared to the creping amplitude, but 𝜆0 is comparable to the creping wavelength. The web with this initial sinusoidal pattern is discretized into a series of particles, with the same distance 𝑎0 between two adjacent particles. Since the web is not flat, the initial local angle 𝛾0 is calculated based on the geometry.   Figure 4.3 A schematic of the initial web pattern due to the forming fabric.  4.2.1 Evolution of Web Morphology         Figure 4.4 shows the evolution of the web morphology with an embedded forming fabric pattern. The parameters for this simulation are chosen as: 𝐴0 = 10 𝜇𝑚 , 𝜆0 = 200 𝜇𝑚 , ℎ =35 𝜇𝑚 , 𝑤 = 1 𝑚 , 𝐿0 = 8 𝑚𝑚 ,  𝑁 = 401 , Δ𝑡 = 1  𝑠 , 𝐸1 = 500 𝑃𝑎 , 𝐸2 = 50 𝑃𝑎 , 𝜎𝑦 =10 𝑃𝑎 , 𝜌 = 428 𝑘𝑔 𝑚3⁄ , 𝑉𝑖𝑛 = 1200 𝑚 𝑚  , 𝑉𝑜𝑢𝑡 =  00 𝑚 𝑚  , 𝐺𝐼𝐶 = 𝐺𝐼𝐼𝐶 = 50 𝑁 𝑚 and 𝛿 =  0°. It is worth mentioning that in this section, we assume that adhesion and basis weight are still uniform along the machine direction, despite that the web has embedded forming fabric pattern. We define the creping wavelength corresponding to an initial flat web as the natural wavelength of the system, denoted by 𝜆𝑐. For the parameters used in Figure 4.4 (other than 𝐴0 and 𝜆0), we find out 𝜆𝑐 = 38  𝜇𝑚 by performing a separate creping simulation with 𝐴0 = 0. Figure 5𝐴0𝜆067  4.4(a)-(c) indicate that the folding process is not exactly periodic, and the sizes of the generated folds are significantly different. This is because when the web buckles, the boundary condition of the web varies from fold to fold due to the initial pattern. Even though 𝐴0 is relatively small compared to the typical creping amplitude, it can significantly change the onset of buckling and hence affect the final structure. Figure 4.4(d)-(f) further show that several stacked small folds evolve into one big fold, which makes the final creping pattern even more irregular.  Note that in this simulation, there was no steady-state solution within the time period of simulation performed. Tissue paper with non-uniform creping structure tends to have a lower softness and strength [3], so it makes sense to optimize the design of the forming fabric to achieve better uniformity. Next, we modify the initial embedded forming fabric wavelength and amplitude to investigate their effects on the creping process and the final properties of generated tissue.  Figure 4.4 The evolution of the web morphology with an initial forming fabric pattern.  (a)(e)(d) (f)(c)(b)𝑡 = 10𝜇𝑠𝑡 = 150𝜇𝑠𝑡 = 120𝜇𝑠𝑡 = 110𝜇𝑠𝑡 = 80𝜇𝑠𝑡 = 50𝜇𝑠68  4.2.2 Effects of Initial Wavelength and Amplitude         Figure 4.5 shows the simulated crepe pattern under various initial wavelengths. The parameters are chosen the same as those used in Figure 4.4 except 𝜆0. Here 𝜆0 is normalized by the natural wavelength 𝜆𝑐, which was found to be 38  𝜇𝑚 for this set of parameters. In Figure 4.5 𝜆0 𝜆𝑐⁄  increases from left to right. It illustrates that when 𝜆0 𝜆𝑐⁄  is smaller than 1, the crepe pattern of the generated tissue paper is quite non-uniform and “chaotic”. The size of the fold varies significantly along the web.  Figure 4.6(a) is the Fast Fourier Transform (FFT) of the crepe pattern corresponding to 𝜆0 𝜆𝑐⁄ = 0.51. It clearly shows multiple peaks and a wide spread across the frequency domain, which further verifies the non-uniformity of the structure. However, as the initial wavelength increases and gets close to the natural wavelength, the crepe structure becomes much more uniform and regular, as shown by the three middle plots in Figure 4.5. Figure 4.6(b) plots the FFT of the crepe structure corresponding to 𝜆0 𝜆𝑐⁄ = 1.02, which shows a single peak and narrow spread in the frequency domain. In this case, a dominant creping wavelength can be identified. However, as 𝜆0 continues to increase and deviates from 𝜆𝑐, the crepe pattern becomes non-uniform again with irregular fold size. It is worth mentioning that the results shown in Figure 4.5 is not sensitive to the initial velocity, and the embedded structure is the main factor governing the crepe structure. 69   Figure 4.5 Simulated creping pattern under different initial wavelength.   Figure 4.6 Fast Fourier Transform of the creping pattern. (a):  𝝀𝟎 𝝀𝒄⁄ = 𝟎. 𝟓𝟏. (b): 𝝀𝟎 𝝀𝒄⁄ = 𝟏. 𝟎𝟐. 𝜆0𝜆𝑐= 0.77𝜆0𝜆𝑐= 0.51𝜆0𝜆𝑐= 0. 0𝜆0𝜆𝑐= 1.02𝜆0𝜆𝑐= 1.15𝜆0𝜆𝑐= 1.2 𝜆0𝜆𝑐= 1.541𝑚𝑚(a) (b)70          The results demonstrated by Figure 4.5 and Figure 4.6 can be understood as follows. When 𝜆0 is away from 𝜆𝑐, the buckling process during the fold formation keeps changing, which leads to irregular crepe pattern. However, when 𝜆0  is close to 𝜆𝑐 , the initial pattern and the natural wavelength becomes “in phase” and have similar periodicity. This makes the boundary condition of the web during folding more consistent and results into a more uniform crepe structure. This suggests that in order to produce a more consistent product, the periodicity of the forming fabric should match the natural wavelength.           After the creping simulation, we perform a series of virtual tensile tests on the generated tissue paper. Figure 4.7 illustrates the corresponding tensile stress-strain curves. In general, when 𝜆0 𝜆𝑐⁄  increases, the initial modulus increases and the stretch decreases. This is because for smaller 𝜆0 𝜆𝑐⁄ , the web is more “packed” initially with more material within a unit length. Thus, the web can be stretched easier and longer. However, for the cases when 𝜆0 is close to 𝜆𝑐, the stress-strain curves are quite similar, as shown by the dashed lines in Figure 4.7. This is because in this range, the properties of the tissue paper are dominated by the material itself, instead of the embedded forming fabric pattern. 71   Figure 4.7 Tensile stress-strain curves of the generated tissue paper with different initial wavelength.  Figure 4.8 Tensile stress-strain curves with different initial amplitude.          Next, we fix 𝜆0 while changing 𝐴0 to see its impact on creping, as shown in Figure 4.8. The parameters used here are the same as those in Figure 4.7, except that we vary 𝐴0 but fix  𝜆0 𝜆𝑐⁄ =1.02. The blue curve in Figure 4.8 with 𝐴0 = 0 represents an initial flat web. As 𝐴0 increases, the initial modulus of the generated tissue paper drops and the stretch increases. This is because 72  increasing 𝐴0 results into a larger final creping amplitude 𝐴. But the change of the final creping wavelength 𝜆 is negligible, as 𝜆0  is fixed as 1.02𝜆𝑐 . As we have shown in Chapter 3, as 𝐴 𝜆 increases, the modulus of the tissue paper decreases and the stretch increases.  4.2.3 Effects of Base Sheet Properties with Initial Pattern Included         In the previous section, we have shown that as the initial wavelength approaches the natural wavelength, the creped tissue paper becomes uniform and has better overall properties. For a specific tissue grade, ideally one can change the forming fabric to make sure its pattern matches with the natural wavelength of that particular grade. However, in practice it is not common to frequently change forming fabric as it adds to additional cost and wastes production time. So it makes sense to see the effects of base sheet properties on creping while fixing the forming fabric. Figure 4.9 illustrates the effects of base sheet modulus on the tensile properties of tissue paper. Here the parameters are chosen as: 𝐴0 = 10 𝜇𝑚 , 𝜆0 = 400 𝜇𝑚 , ℎ = 35 𝜇𝑚 , 𝑤 = 1 𝑚 , 𝐿0 =8 𝑚𝑚 ,  𝑁 = 401 , Δ𝑡 = 1  𝑠 , 𝜎𝑦 = 0.02𝐸1 , 𝜌 = 428 𝑘𝑔 𝑚3⁄ , 𝑉𝑖𝑛 = 1200 𝑚 𝑚  , 𝑉𝑜𝑢𝑡 = 00 𝑚 𝑚  , 𝐺𝐼𝐶 = 𝐺𝐼𝐼𝐶 = 50 𝑁 𝑚  and 𝛿 =  0° . Here the elastic modulus 𝐸1and the plastic modulus 𝐸2 are changing simultaneously while fixing 𝐸1 𝐸2 = 10. It indicates that as the base sheet modulus increases, the initial modulus and the stretch of the generated tissue paper both increase. This is consistent with the results obtained with an initial flat sheet (see Table 3.1). 73   Figure 4.9 Tensile stress-strain curves of tissue paper generated with different modulus.  4.3 Non-uniform Basis Weight Distribution 4.3.1 Effect of Variation Length         Paper is essentially a network of fibers interconnected by bonds. During the formation process, if the fibers are not evenly deposited onto the fabric, it causes a non-uniform basis weight distribution in the web. In this section, we introduce a periodic basis weight variation in the web to study its effect on the creping process, as shown in Figure 4.10(a). 𝑙 denotes the variation length of basis weight, which is typically close to the fiber length scale. We assume that the length of low and high basis weight region are the same. The dark area is corresponding to low basis weight region, while the gray color represent the high basis weight region. The density and thickness of the web are also varied based on the local basis weight value, according to the empirical equations ℎ = ( + 𝐵𝑊)  ⁄  and 𝜌 =  ∙ 𝐵𝑊 ( + 𝐵𝑊)⁄ , where   and   are two constants [89]. This allows the web has different bending stiffness at different regions. Figure 4.10(b)-(d) are the crepe structures of simulated tissue paper under different values of variation length, but with the same 74  initial total length. The parameters are chosen as: 𝐿0 = 8 𝑚𝑚 ,  𝑁 = 401 , Δ𝑡 = 1  𝑠 , 𝐸1 =500 𝑃𝑎 , 𝐸2 = 50 𝑃𝑎 , 𝜎𝑦 = 10 𝑃𝑎 , 𝑉𝑖𝑛 = 1200 𝑚 𝑚  , 𝑉𝑜𝑢𝑡 =  00 𝑚 𝑚  , 𝐺𝐼𝐶 =𝐺𝐼𝐼𝐶 = 50 𝑁 𝑚 and 𝛿 =  0°. The web properties are given by 𝐵𝑊1 = 12 𝑔 𝑚2 , ℎ1 = 32𝜇𝑚, 𝜌1 = 375 𝑘𝑔 𝑚3⁄  and 𝐵𝑊2 = 18 𝑔 𝑚2 , ℎ2 = 38 𝜇𝑚, 𝜌2 = 474 𝑘𝑔 𝑚3⁄ . These parameters are in the range of typical industrial values. Here the subscripts 1 and 2 represent the low and high basis weight regions, respectively. Figure 4.10 clearly shows that the folds formed in higher basis weight region are larger. This is mainly because the bending stiffness in high basis weight region is larger. It also shows that when 𝑙 is relatively small, the folding process within one region is more likely affected by the other region, thus creating an overall non-uniform and “chaotic” crepe structure. However, when 𝑙 is sufficiently large, the folds formed in two different regions only interact at the location where two regions are connected. Within each region the folds are uniform, as shown in Figure 4.10(d).   Figure 4.10 (a) Initial flat web with periodic basis weight distribution. 𝒍 represents the variation length. (b)-(d) Crepe structures of simulated tissue paper with different variation lengths.           Figure 4.11 illustrates the evolution of the creping force during the creping process corresponding to 𝑙 = 4𝑚𝑚. It indicates that the creping force in the high basis weight region is Low Basis WeightHighBasis WeightWeb𝑙 𝑙(a)(b)(c)(d)𝑙 = 1𝑚𝑚𝑙 = 2𝑚𝑚𝑙 = 4𝑚𝑚75  higher and the number of peaks is lower, i.e., less folds are formed. This is consistent with our earlier conclusion that the creping wavelength and amplitude increase as the basis weight increases (see Chapter 3, Table 3.1). Figure 4.11 shows that the creping force variation time scale is around 200 𝜇𝑠. Since the data acquisition frequency in practice can reach 10KHz, the averaged creping force within each region can be measured. This may provide some insights for the tissue machine operators to better understand the creping force data captured during production.   Figure 4.11 Evolution of the creping force during creping, with variation length  𝒍 = 𝟒𝒎𝒎.          We take the generated tissue paper in Figure 4.10 and perform virtual tensile tests to investigate the effects on the paper properties. Figure 4.12(a) demonstrates that as the variation length decreases, the initial modulus of the tissue decreases and the stretch increases. Figure 4.12(b) further shows the evolution of the crepe structure during the tensile test of the tissue paper with 𝑙 = 4𝑚𝑚. It indicates that both regions are stretched, but the folds formed within the lower 76  basis weight region (dark region) are stretched longer. This is because in that region there is less material taking the tensile load and the stiffness is lower.    Figure 4.12 (a) Tensile stress-train curves of simulated tissue paper with different variation lengths. (b) The crepe pattern at different strains during the virtual tensile test corresponding to 𝒍 = 𝟒𝒎𝒎.  4.3.2 Basis Weight Ratio Effect         In the previous section, we have studied the effect of basis weight variation length on creping. Here, we fix the variation length and change the basis weight ratio (𝐵𝑊2 𝐵𝑊1⁄ ) to see how it influence the creping process. A high basis weight ratio is typically associated with a poor formation. Figure 4.13(a) illustrates the tensile stress-train curves of tissue paper generated under different basis weight ratios. The parameters used are: 𝐿0 = 8 𝑚𝑚,  𝑁 = 401, Δ𝑡 = 1  𝑠, 𝐸1 =500 𝑃𝑎 , 𝐸2 = 50 𝑃𝑎 , 𝜎𝑦 = 10 𝑃𝑎 , 𝑉𝑖𝑛 = 1200 𝑚 𝑚  , 𝑉𝑜𝑢𝑡 =  00 𝑚 𝑚  , 𝐺𝐼𝐶 =𝐺𝐼𝐼𝐶 = 50 𝑁 𝑚 and 𝛿 =  0°. The variation length for all three cases is fixed as 𝑙 = 2𝑚𝑚, and the average basis weight is 15 𝑔 𝑚2. Figure 4.13(b) represent the generated crepe structures before the virtual tensile test. The sheet is uniform when the basis weight ratio is low. As the ratio increases, the folds generated in different regions become quite different and exhibit a higher 𝜀 = 0.110𝜀 = 0.088𝜀 = 0.0  𝜀 = 0.044𝜀 = 0(a) (b)77  irregularity. The tensile curves show that as the ratio increases, the stretch decreases. This suggests that a more uniform basis weight distribution can provide better stretchability in tissue paper.   Figure 4.13 (a) Tensile stress-strain curves of tissue paper with different basis weight ratios. (b) The corresponding crepe patterns at the beginning of tensile test.  4.4 Non-Uniform Adhesion between Web and Yankee         There are two main causes for non-uniform adhesion between the web and the Yankee surface. First, the amount of adhesive chemicals sprayed by the nozzles are not consistent during the creping process. Second, the forming fabric pattern is embedded within the web. This non-flat pattern can result in uneven contact between the web and the Yankee, which further leads to non-uniform adhesion. Based on the periodic nature of the creping process and the forming fabric pattern, it is reasonable to assume that the adhesion has a sinusoidal distribution along the machine direction: 𝐺𝐼𝐶( ) = 𝐺𝐼𝐼𝐶( ) =𝐺𝑚𝑎𝑥(sin(2𝜋𝑥𝜆𝑎)+1)2                                           (4.1) Here x denotes the coordinate in the machine direction, 𝐺𝑚𝑎𝑥 represents the maximum adhesion and 𝜆𝑎 represents the variation length of the adhesion. Mode I and mode II fracture energies are (a) (b)𝐵𝑊2 𝐵𝑊1 = 1.0⁄𝐵𝑊2 𝐵𝑊1 = 1.5⁄𝐵𝑊2 𝐵𝑊1 = 2.0⁄78  treated the same in this section, for simplicity. Next, we use the one-dimensional particle dynamics model to investigate the effects of non-uniform adhesion on creping.  4.4.1 Effects of Adhesion Variation Length          We assign different values of adhesion variation length (𝜆𝑎) and run separate simulations to see the effects of non-uniform adhesion. Figure 4.14 shows the crepe structures of the simulated tissue paper and the subsequent tensile stress-train curves. The parameters are chosen as: 𝐿0 =8 𝑚𝑚 ,  𝑁 = 401 , Δ𝑡 = 1  𝑠 , 𝐸1 = 500 𝑃𝑎 , 𝐸2 = 50 𝑃𝑎 , 𝜎𝑦 = 10 𝑃𝑎 , 𝑉𝑖𝑛 = 1200 𝑚 𝑚  , 𝑉𝑜𝑢𝑡 =  00 𝑚 𝑚  , 𝐺𝑚𝑎𝑥 = 50 𝑁 𝑚 , 𝛿 =  0°  and 𝐵𝑊 = 15 𝑔 𝑚2 . Figure 4.14(b) clearly indicates that the crepe structure is quite irregular when  𝜆𝑎 is relatively small. The size of the folds varies significantly across the machine direction. However, as 𝜆𝑎 increases, the crepe structure becomes much more uniform, with almost similar folds. This is probably because when 𝜆𝑎 is too small, the adhesion is changing frequently during the formation of individual fold, which makes the buckle-delamination process more “chaotic” and less repeatable. However, as 𝜆𝑎 increases the initial modulus of the tissue paper is increased and the stretch is reduced, as shown in Figure 4.14(a). This is because these two properties are more related to the ratio between the creping amplitude and the creping wavelength. This also suggests that even though the crepe structure becomes more uniform as 𝜆𝑎 increases, the “bulk softness” does not necessarily increase. In fact, the uniformity of the structure is thought to be related to the “surface softness”, as a more uniform texture makes the surface smoother [94].   79   Figure 4.14 (a) Tensile stress-strain curves of tissue paper produced under different adhesion variation length. (b) The crepe structure of the produced tissue paper before the virtual tensile test.  4.4.2 Creping Force Behavior         Figure 4.15 illustrates the evolution of the creping force during the creping process with 𝜆𝑎 =200𝜇𝑚. In the second Chapter we have shown that, when the adhesion is uniform the creping force first reaches the maximum value and then stays constant as the delamination propagates, as shown in Figure 2.7. However, here when the adhesion varies along the interface between the web and the Yankee, the creping force keeps oscillating once it reaches the peak value. Since the creping force can be easily measured in practice, this creping force behavior can be used to check whether the adhesion is consistent during the production of tissue paper.    𝜆𝑎 = 200𝜇𝑚𝜆𝑎 = 300𝜇𝑚𝜆𝑎 = 400𝜇𝑚1𝑚𝑚(a) (b)80   Figure 4.15 The evolution of creping force during the creping process with 𝝀𝒂 = 𝟐𝟎𝟎𝝁𝒎.  4.5 Summary and Conclusions         For the first time, the effects of inhomogeneity on creping and the tissue paper properties are modeled in this Chapter. We use the developed one-dimensional particle dynamics model to investigate these effects. Specifically, three common causes of inhomogeneity are studied separately, including the forming fabric pattern, the non-uniform basis weight distribution due to poor formation and the non-uniform adhesion between the web and the Yankee. The main results are summarized as follows: 1. We assume that the forming fabric imposes a sinusoidal pattern on the base sheet before creping. The effects of both the initial amplitude 𝐴0 and initial wavelength 𝜆0 are studied. It is found when 𝜆0  deviates from the natural wavelength of the system, i.e. 𝜆𝑐 , the simulated crepe structure is significantly non-uniform. When 𝜆0 is close to 𝜆𝑐, the crepe structure becomes uniform with almost identical folds. The virtual tensile tests suggest 81  that as 𝜆0 increases or 𝐴0 decreases, the initial modulus of the tissue paper increases and the stretch decreases.  2. The effect of non-uniform basis weight distribution is studied by assuming that the base sheet has alternative low and high basis weight regions. As the variation length decreases, the initial modulus of the tissue decreases and the stretch increases. As the basis weight ratio increases, the folds generated in different regions become quite different and exhibit a higher irregularity. The tensile tests show that as the ratio increases, the stretch decreases, but the initial modulus almost does not change. 3. The effect of non-uniform adhesion is accounted for by assuming a periodic distribution along the machine direction. The simulation results suggest that, as the adhesion variation length 𝜆𝑎 decreases, the generated tissue paper becomes irregular and non-uniform. Also, the creping force oscillates during the delamination propagation stage, which can be used to monitor the uniformity of adhesion during production. It is worth mentioning that these causes of inhomogeneity can arise simultaneously during creping. For example, the 3D topological structure of the forming fabric can lead to both uneven base sheet and non-uniform adhesion. Here, we have studied them solely for the purpose of differentiating the effect of each cause. The assumptions made about the forming fabric pattern, and the distributions of basis weight and adhesion are reasonable as a starting point. The results presented in this Chapter can hopefully provide some new insights to the industry.  82  Chapter 5: Experimental Validation 5.1 High-Speed Imaging Study of Creping          Several experiments have been performed to directly observe the creping process and understand its underlying mechanisms. For example, Hollmark reported that when the creping angle is small, a series of micro-folds pile up and transform into a macro-fold [11]. However, this observation is based on a low-speed pilot tissue machine where the maximum creping speed is limited to 140 𝑚 𝑚  . This micro-fold to macro-fold transition mechanism has not been verified in any industrial or a high-speed tissue machine. Recently, creping is commonly accepted as a periodic debonding and bulking process. This mechanism was first described by Ramasubramanian et. al., based on their observation on a lab-scale tissue machine with a maximum creping speed of 250 𝑚 𝑚   [9]. To our knowledge, there is no direct observation of creping process in the high-speed range. To qualitatively verify these reported mechanisms, we have performed a high-speed imaging study on an existing lab-scale creping rig, which was developed by Solenis (Delaware, USA). In this section, details of the experiment are given.  5.1.1 Experimental Setup         Figure 5.1 shows an existing lab-scale creping apparatus developed by the industry (Solenis, USA). This device mainly includes five components: a heated Yankee cylinder, adhesive sprayers, a transfer roll, a moisture controlling system and a creping blade. During the experiment, first a small piece of wet paper is attached on the transfer roll with controlled moisture content. Then, adhesive chemicals are spayed onto the surface of the Yankee cylinder. The volume of the adhesives can be varied to reach different levels of adhesion. Next, the wet paper is transferred to the Yankee and bonded on its surface. The Yankee is internally heated and the surface temperature 83  is maintained between 100 − 110 . After the wet paper is transferred, the Yankee starts to rotate and the paper is dried. Finally, when the surface speed of the Yankee reaches the desired speed, the doctor blade is triggered by a step motor to engage with the paper and scrap it off. The surface speed of this creping apparatus can reach up to 1524 𝑚 𝑚  . The angle of the doctor blade, adhesive chemicals, surface speed, surface temperature, moisture content, and base sheet properties can be varied to see the corresponding effects. We use a high-speed camera to capture the creping process and observe the underlying mechanism. A high-intensity line laser is used to illuminate the edge of the paper sample, shown in Figure 5.1. This can provide an edge-on view of the paper during creping.   Figure 5.1 Lab-scale creping rig and the setup of high-speed imaging study.  CameraLine LaserBladeYankeeTransfer Roll84  5.1.2 Observation of Creping Process          Figure 5.2 shows the edge-on view of the creping process captured by the high speed camera. The Yankee surface and the doctor blade are not visible in the images as only the edge of the tissue paper was illuminated by the line laser. The basis weight of the sample used for the experiment is 25 𝑔 𝑚2 and the thickness is  3.  𝜇𝑚. The other parameters calibrated during the experiment are given as: 𝑉𝑖𝑛 = 1120 𝑚 𝑚  , 𝛿 =  5 . For the current creping rig, the sample is not a continuous web and the creping ratio effect is not included. Figure 5.2(a) shows that the flat web is moving towards the blade, during which a compressive stress builds up in the web and delamination starts to propagate. Once the delamination length reaches a critical length, buckling of the debonded web occurs, as shown in Figure 5.2(b). After buckling initiation, the web undergoes post-buckling deformation and self-contact occurs, as demonstrated in Figure 5.2(c). The new fold pushes the previously formed fold out along the blade edge. Finally, the formation of one fold is completed and similar process repeats. The high-speed imaging experiment has confirmed that the periodic debonding and buckling mechanism of creping applies to the high-speed range as well. The micro-fold to macro-fold transition was not observed here, due to the relatively large creping angle. To our knowledge, it is the first time that a creping process of industrial speed was captured by experiments. 85   Figure 5.2 Snapshots of the creping process captured by high-speed camera.  5.1.3 Parametric Studies Based on Creping Rig         We use the creping rig shown in section 5.1.1 to conduct parametric studies. The parameters that can be changed are creping velocity 𝑉𝑖𝑛, temperature of the Yankee surface 𝑇𝑌, creping angle 𝛿, and spraying time of adhesives 𝑡𝑎. Due to the limited access to the creping rig in the industry, each parameter was only varied once. These parameters are assigned with typical values used in a real tissue machine. The size of base sample is 76 mm wide and 508 mm long, with basis weight (a)(d)(c)(b)𝑉𝑖𝑛𝑉𝑜𝑢𝑡BladeBucklingInitiationPost-Buckling One Fold Completed𝑡 = 0𝜇𝑠𝑡 = 20𝜇𝑠𝑡 = 15𝜇𝑠𝑡 = 8𝜇𝑠DelaminationPropagation86  equals  25 𝑔 𝑚2 and thickness equals  3.  𝜇𝑚. Figure 5.3 illustrates the effects of each control parameter. The horizontal axis represents the averaged creping force measured during the creping process. The vertical axis represents the length of the tissue paper after creping, denoted by 𝐿𝑡. Each scatter point in Figure 5.3 represents one individual tested sample. In the legends, the units for 𝑉𝑖𝑛, 𝑇𝑌, 𝛿 and 𝑡𝑎 are 𝑚 𝑚  , , degree and second, respectively. Since the width of the sample is small and the adhesion near the edge of the sample is significantly lower compared to the middle, the overall adhesion is not well controlled, which leads to quite large variation of the data under the same condition.         Figure 5.3(a) clearly indicates that when the creping angle decreased from  5° to 80°, the creping force increased significantly. This agrees with our simulation results shown in Figure 2.12. Note that as 𝛿  decreases, the total length of the generated tissue paper (𝐿𝑡 ) also decreased. However, here 𝐿𝑡  should not be used to interpret the creping wavelength. This is because the shrinkage of paper after creping also depends on the creping amplitude. Here as 𝛿 decreases, the generated tissue paper becomes more difficult to move out due to the reduced space, thus the tissue paper is more compressed and 𝐿𝑡 drops. Figure 5.3(b) demonstrates that as the creping velocity is increased, the creping force is also increased, which agrees with the analytical solution of creping force we derived in Chapter 2. However, to further validate the linear relation between creping force and creping velocity requires more data. Figure 5.3(c) shows that by increasing the spraying time from 8 seconds to 12 seconds, the effective adhesion increases, resulting in increased 𝐹𝑐 and decreased 𝐿𝑡. This is also consistent with our simulation results, as shown in Figure 2.13. Finally, Figure 5.3(d) illustrates the effect of Yankee surface temperature. 𝑇𝑌 = 100  and 𝑇𝑌 = 110  were used, which are the typical values in the real tissue machine. However, the effects on 𝐿𝑡 and 𝐹𝑐 are not conclusive in our experiments. The variation of temperature can cause the change of 87  moisture of the web and the cross-linking process of the adhesive polymers, so the web properties and the adhesion are also subject to change. To fully understand the temperature effect, more specific experiments are needed in a future study.          Since creping is a very complicated process with various control parameters, it is unrealistic to measure all the parameters and aim for a quantitative comparison. However, the experimental results shown in Figure 5.3 qualitatively validates our simulation results.  Figure 5.3 Parametric studies based on lab-scale creping rig. (a) Effect of creping angle. (b) Effect of creping velocity. (c) Effect of spraying time of adhesives. (d) Effect of Yankee surface temperature.  (a)(d)(c)(b)88  5.2 Tensile Test of Tissue Paper         The correlation between the tissue properties and the crepe structures have not been specifically investigated before. In Chapter 3, we have developed a virtual tensile test based on the particle dynamics model to study the correlation. The simulated stress-strain curves show significant nonlinear behavior, as shown in Figure 3.7 and Figure 3.8. To verify the simulation results, in this section we perform real tensile test on commercial tissue paper and compare the experiments with simulation.         The uniaxial tensile tests of tissue paper were performed on Instron Machine (5566) according to the standards ISO 12625. The sample test length is 100 mm with 25.4 mm in width. The extension rate is fixed as (50 ± 2) 𝑚𝑚 𝑚   and the experiments were done in a humidity (50%) and temperature-controlled room (23 °𝐶) to make sure the results are consistent. Four types of samples with different grades were taken directly from a commercial tissue paper machine right after the creping process. We have also changed the extension rate to check the viscous effect. But the difference between different extension rates is negligible.           Figure 5.4(a) shows the uniaxial tensile stress-strain curves along the machine direction. The test for each grade was repeated ten times. All the curves consistently show nonlinear behavior with a mild inflection point where the strain is around 7%, and the slope keeps increasing until the sample breaks. Although there are more mechanisms (fracture, strain localization, heterogeneity…) involved in the experiments, our simulated tensile curve as shown in Figure 5.4(b) (also in Figure 3.7 and Figure 3.8), exhibits the same nonlinearity in Figure 5.4. This indicates that the folding pattern and the plastic deformation of the web captured in our model are the key aspects governing the tensile properties of the tissue paper. The tensile test results are summarized in Table 5.1. ℎ and 𝜆 were measured before the tests. Based on the data, grade 1 and 89  grade 3 have smaller ℎ and 𝜆 compared to grade 2 and grade 4. The initial modulus corresponding to grade 1 and grade 3 are much higher, and the stretch values are slightly lower. It is worth mentioning that we have also changed the extension rate to check the viscous effect, which turned out to be negligible within the extension rates (10 − 100 𝑚𝑚 𝑚  ) we considered.    Figure 5.4 (a) Stress-strain curves of uniaxial tensile tests on commercial tissue paper in machine direction. (b) Typical stress-strain curves by simulation.   Stretch (%) Strength (N/m) Initial Modulus (MPa) BW (gsm) Caliper (µm) 𝜆 (µm) Caliper/BW  (𝑐𝑚3 𝑔) Grade 1 18.32 ± 2.02 97.95 ± 7.05 15.81 ± 1.33 14.3 50 249±45 3.50 Grade 2 23.89 ± 0.88 107.14 ± 9.11 6.13 ± 0.27 20.1 79 306±65 3.93 Grade 3 21.14 ± 1.13 111.76 ± 11.56 10.92 ± 0.47 16.9 58 248±48 3.43 Grade 4 22.72 ± 0.81 85.4 ± 6.07 4.32 ±  0.14 17.9 76 325±61 4.24 Table 5.1 Summary of the uniaxial tensile test results of four grades.  (a) (b)90  5.3 Validation of Basis Weight Effect by SEM         Figure 5.5 shows the cross-sectional images of two commercial tissue samples manufactured under the same condition except the basis weight. In a tissue machine, the basis weight is typically controlled by the amount of fibers deposited on the forming fabric. Figure 5.5 clearly indicates that with higher basis weight, the creping wavelength and creping amplitude are both larger. This finding is consistent with our simulation results reported in Table 3.1 and Figure 4.13. More details about this creping wavelength analysis based on cross-sectional images can be found in [95]. In our model, the basis weight is modified by changing thickness and density of the web. However, in practice the change of basis weight is often related to the modification of headbox consistency, which may affect other properties such as formation, moisture content and residual stress.  Figure 5.5 SEM cross-sectional images of commercial tissue paper with two different basis weights: 𝟏𝟓. 𝟑 𝒈 𝒎𝟐 (top) and 𝟐𝟎. 𝟏 𝒈 𝒎𝟐 (bottom).  5.4 Summary and Conclusions         In this Chapter, we have reported a series of experiments to qualitatively validate the particle dynamics model presented in the previous Chapters. The main findings are summarized as follow: 𝟏𝟓. 𝟑 𝒈 𝒎𝟐500 µm20. 𝟏 𝒈 𝒎𝟐500 µm91  1. A high-speed imaging study has been performed on a lab-scale creping apparatus. A high-speed camera was used to capture the folding process of the web during creping. The images indicate similar periodic debonding and buckling mechanism as our simulations. 2. Effects of creping velocity, temperature of the Yankee surface, creping angle, and spraying time of adhesives have been studied based on the creping rig. The trends of these parameters are consistent with our simulations. 3. Uniaxial tensile tests of four different grades of commercial tissue paper have been performed. The experimental stress-strain curves demonstrate consistent nonlinearity, and qualitatively matches with the results from virtual tensile tests. 4. The effect of basis weight on the crepe structures was investigated, based on observation of the cross-sectional images of commercial tissue paper. The images obtained by SEM show that the creping amplitude and the creping wavelength both increase as basis weight increases, which agrees with the predictions based on particle dynamics model. The creping rig used for parametric study is a smaller version of a real tissue paper machine, with important components such as a heated Yankee cylinder, adhesive sprayers, a transfer roll, a moisture controlling system and a creping blade. Realistic values for the control parameters are used. However, the major disadvantage is that the creping ratio is not included. Nevertheless, the experimental results are still trustworthy and useful for qualitative validation of our model.   92  Chapter 6: Summary and Conclusions         This thesis has studied the creping process in the manufacturing of low-density tissue paper in which a series of micro-folds are created and inter-fiber debonding occurs within the sheet. Creping is commonly understood as a periodic delamination and buckling process based on experimental observations, and significant efforts have been made to model the creping process. However, most of the previous models have only considered the shaped-bulk regime, and neglected the plasticity of the web and the creping ratio effect, which are important in creping. The mechanism behind the explosive-bulk regime (“sheet explosion”), the correlation between the crepe structures and the tissue properties, and the effects of inhomogeneity on creping remains unclear. These aspects are studied in this thesis. The main results, limitations and directions for future work are summarized in this Chapter.  6.1 Methodology: Particle Dynamics Model         We have developed a one-dimensional particle dynamics model as our main approach to systematically investigate the creping process and evaluate the tissue properties. The advantages of particle-based models are discussed at the end of Section 1.2.3. In this model, the web is discretized into a series of particles that are connected by spring-damper systems. The elastoplastic response of the web is modelled by the springs, while the viscosity is accounted for by the dampers. A bilinear elastoplastic material model with associated kinematic hardening rule (Bauschinger effect) is used to describe the effective constitutive relation of the web, and a mixed-mode cohesive zone model is adopted to model the interfacial delamination. It is worth mentioning that the inelastic behavior of the web, such as the fiber rotation and sliding, are indirectly accounted for by the constitutive model. The self-contact of the web as well as the contact between the web and the 93  doctor blade, are modeled based on a penalty method. At each time step, the total force applied on each particle is calculated, and used in Newton’s equation of motion to evolve the particle system. The developed model is validated (see Appendix A and B), and shown to be able to capture the large visco-elasto-plastic deformation of the web, the mixed-mode interfacial fracture and the contact mechanisms. The main parameters investigated in the model include the creping velocity, creping ratio, creping angle, adhesion between the web and the Yankee surface, and the web properties (thickness, density, basis weight and Young’s modulus). Realistic values of these parameters from industrial data and previous literature are used in our model. A virtual tensile test on the simulated tissue paper is implemented to find the correlation between the crepe structure and the tissue properties. The model is also extended to multi-layer to account for the explosive-bulk regime in creping. The model is able to capture the large nonlinear deformations of the web in machine direction and thickness direction. Note that the angular momentum is not precisely conserved, because the rotational degree of freedom of each particle is not taken into account. The torsion and cross direction displacement are not included in the proposed model, since the debonding and buckling mainly occurs along the machine direction during creping. This is also consistent with the fact that the paper properties in machine direction has a more significant change after creping, as shown in Figure 1.2.         To our knowledge, this is the first continuum model for creping that includes both the plasticity of the web and the creping ratio effect, and be able to simulate the shaped-bulk regime and the effects of inhomogeneity. The presented model has provided a useful tool to study the process-structure-property correlation of tissue paper. The main results are summarized in the next section.  94  6.2 Summary of Results         In Chapter 2, we investigate the shaped-bulk regime of creping using a particle dynamics model. The formation of an individual fold is found to include three distinct stages: delamination propagation, buckling initiation, and post-buckling deformation. In the first stage, only the mode II fracture energy (𝐺𝐼𝐼𝐶) is important as the web remains flat. While in the third stage, the post-buckling of the web further triggers the delamination to propagate in mixed-mode (both 𝐺𝐼𝐶 and 𝐺𝐼𝐼𝐶  are important), which leads to a larger creping wavelength. An analytical solution of the maximum creping force in the first stage is derived based on an energy analysis and assuming the web remains elastic. This solution matches well with the simulation results. It indicates that the creping force is the sum of a contribution from adhesion and a contribution from web inertia. The effects of mixed-mode fracture energy on creping is investigated. Increasing the ratio between 𝐺𝐼𝐼𝐶 and 𝐺𝐼𝐶, leads to decrease of the average creping force and increase of the creping wavelength.          In Chapter 3, we extend the model to include the plasticity of the web. A bilinear elastoplastic constitutive relation with associated kinematic hardening rule is used to model the web. The plastic deformation of the web is found to reduce the creping force and the delamination propagation speed, which leads to the decrease of both the creping amplitude and wavelength. We take the simulated creped tissue paper and perform a virtual tensile test to quantify the initial modulus and the stretch. As the ratio between the creping amplitude and the creping wavelength increases, the modulus of the tissue paper reduces and the stretch increases. A systematic parametric study shows that the basis weight and the creping angle have the highest impact on 𝜆, A and 𝐹?̅? , while the creping angle and the creping ratio exhibit a significant impact on the stretch. By using the extended multi-layer model, a phase diagram for creping regimes is constructed, which indicates that the explosive-bulk regime favors high adhesion and low inter-fiber bonding strength. 95          In Chapter 4, the effects of inhomogeneity on the creping process and the tissue paper properties are studied for the first time. Three common causes of inhomogeneity are studied separately, including the forming fabric pattern, the non-uniform basis weight distribution due to poor formation, and the non-uniform adhesion between the web and the Yankee. First, we assume an initial sinusoidal pattern is imposed on the base sheet due to the forming fabric. Our simulations show that the crepe structure is uniform only when the initial wavelength 𝜆0 is close to the natural wavelength of the web, i.e. 𝜆𝑐. It is also found that the initial modulus of the tissue paper increases and the stretch decreases, when 𝜆0 increases or 𝐴0 decreases. Second, the effect of non-uniform basis weight is studied by assuming that the base sheet has alternative low and high basis weight regions, with variation length in the fiber length scale. Decreasing the variation length results in the decrease of the initial modulus and the increase of the stretch of tissue paper. It is also found that, as the basis weight ratio increases the generated microfolds in the web exhibit a higher irregularity. Finally, we consider the non-uniform adhesion by assuming a periodic distribution along the machine direction. The simulate crepe structure becomes irregular and non-uniform as the adhesion variation length 𝜆𝑎 decreases. Also, the creping force is found to oscillate around the peak values at the delamination propagation stage of creping.         In Chapter 5, we have conducted a series of experiments on a lab-scale creping apparatus developed by the industry (Solenis, USA), to qualitatively validate our model. The creping process is first observed at an industrial speed by a high speed camera, which has shown similar periodic debonding and buckling mechanism as our simulations. Effects of creping velocity, creping angle, and adhesion have been investigated, and the results agree with the trends predicted by our model. Tensile tests of four different grades of commercial tissue paper in Machine Direction have been performed. The experimental stress-strain curves have shown similar nonlinearity as our virtual 96  tensile tests, i.e., the stiffness first decreases and then increases until the paper breaks. This implies that the visco-elasto-plastic behavior of the web and the folding mechanism captured by our model, are two essential aspects governing the macroscopic properties of the tissue paper. Finally, SEM cross-sectional images of commercial tissue paper have shown that increasing the basis weight result in larger creping amplitude, creping wavelength and caliper to basis weight ratio, which agrees with the predictions from the model.         The known knowledge of creping and the new findings of our study are summarized below to highlight the important outcomes of this work: Known Creping Mechanisms:  Adhesion 𝐺𝑐: increasing 𝐺𝑐 reduces 𝜆 and 𝛢, increases 𝐹𝑐.  Creping angle 𝛿: increasing 𝛿 reduces 𝜆, 𝛢, and 𝐹𝑐 ; micro-fold to macro-fold transition exists at low 𝛿.  Creping velocity 𝑉𝑖𝑛: increasing 𝑉𝑖𝑛 increases 𝐹𝑐, reduces 𝜆 and 𝛢.   Basis weight 𝐵𝑊: increasing 𝐵𝑊 increases 𝜆 and 𝛢.  Sheet bending stiffness 𝐸𝐼: increasing 𝐸𝐼 increases 𝜆 and 𝐹𝑐.  Yankee surface temperature 𝑇𝑌 : increasing 𝑇𝑌  reduces tensile strength, thickness and stretch. New findings:  Analytical solution for the maximum creping force: 𝐹𝑚 = √2𝐸ℎ𝐺𝐼𝐼𝐶 +  ℎ𝑉𝑖𝑛√𝜌𝐸 (adhesion + inertia), assuming the web is elastic.  Mixed-mode fracture during the third stage of folding (post-buckling deformation).  Elastic modulus 𝐸1: decreasing 𝐸1 reduces 𝜆, 𝛢, 𝐹𝑐, modulus and stretch of tissue. 97   Plastic modulus 𝐸2: decreasing 𝐸2 reduces 𝜆, 𝛢, and 𝐹𝑐, and makes crepes non-uniform.  Creping ratio 𝐶𝑅: increasing 𝐶𝑅 reduces 𝜆 and modulus, but increases 𝛢 and stretch.  𝐵𝑊 and δ have the highest impact on λ, 𝛢, Fc and stretch.  Explosive-bulk favors high adhesion 𝐺𝑐 and low inter-layer bond strength.  Initial wavelength 𝜆0: when 𝜆0 is close to 𝜆𝑐, crepe pattern becomes more uniform.  Initial amplitude 𝐴0: increasing  𝐴0 increases stretch and reduces modulus.  Basis weight variation length 𝑙: decreasing 𝑙 increases stretch and reduces modulus.  6.3 Significance of Current Work         Our study in this thesis provides some new practical guidance to the industry. For example, the parametric study (Table 3.1) suggests that, the creping angle or the furnish (basis weight) should be adjusted to achieve a significant change of the crepe structure. If smaller changes are required, then the adhesive chemicals and the creping ratio can be modified accordingly. It is thought explosive-bulk is required to achieve higher softness tissue products. Based on our constructed phase diagram of creping regime (Figure 3.11), it seems plausible that in order to increase the softness, one might increase the amount of adhesive chemicals applied on the Yankee surface and reduce the pressure applied by the “pressure roll” on the web at the beginning of creping. However, excessive adhesion should be avoided as it can lead to “picking”, chatter and blade wear. The results obtained from virtual tensile test (Figure 3.7 and Figure 3.8) implies that the ratio between the creping amplitude and the creping wavelength governs the initial modulus and stretch of the tissue paper. Thus, it is insufficient to judge the quality (such as softness) of different tissue papers just by comparing their “crepe counts” (1 𝜆⁄ ). The study on the forming fabric pattern (Figure 4.5 and Figure 4.6) indicates that, in order to create a uniform crepe pattern 98  and achieve a higher surface softness, one should design the forming fabric to match the natural wavelength of the targeted tissue grade. When a web has non-uniform basis weight or adhesion distribution, the corresponding creping force exhibits a unique behavior (Figure 4.11 and Figure 4.15). Since the creping force is an accessible data in the tissue machine, our simulation results and analytical solution (Eq. 2.23) can potentially help the industry to better interpret the creping force data in practice.  6.4 Limitations and Future Work         The limitations of the current study and the potential future work are discussed in this section as follows: 1. The fracture of the web is not considered in the developed models. We end the virtual tensile test of the simulated tissue paper once the stress exceeds 30% of the yield stress, considering the fact that the typical strength reduction of paper is around 70% after creping. The virtual tensile test has shown that the modulus of the sheet typically decreases by more than 90% after creping. One possible extension is to implement damage evolution model for the web and predict the fracture. 2. The cohesive zone model describing the interfacial delamination is assumed rate-independent. In reality, the adhesive polymers may exhibit viscous characteristics and the adhesion may depend on the creping velocity. In the future, a rate-dependent cohesive zone model can be used to account for that. Peel tests of adhesives with different peel angles and peel speeds can be conducted to estimate the mixed-mode fracture energies and their rate-dependence.  99  3. Currently in our model we did not consider the residual tensile strain in the base sheet caused during the drying process. In the future, this effect can be incorporated by setting the initial distance between the particles to be 𝑎0(1 + 𝜀𝑟), where  𝜀𝑟 denotes the residual strain and can be estimated by experiments. 4. The effect of three types of inhomogeneity on creping are investigated separately. But in practice, these causes of inhomogeneity can arise simultaneously. In the future, the model can be used to study the coupled effect of forming fabric pattern, non-uniform formation and adhesion. 5. The inter-fiber bonds in the multi-layer model are assumed to be identical and evenly distributed. However, incoming web is essentially a random fiber network, and both the strength and the location of inter-fiber bond vary significantly. This can be accounted for by applying a stochastic distribution of cohesion in the multi-layer model in the future. 6. The current model is restricted to 2D (the thickness-direction and the machine-direction). To account for the effects of cross-direction variation, the current model should be extended to a 3D fiber network model. However, this requires more fiber-level parameters that are difficult to measure, and also significantly increases computational cost. We leave this aspect for future investigation.  7.  Our simulation results are limited to qualitative comparison with the experiments. This is because it is extremely difficult to accurately measure the incoming web properties (density, viscosity, Young’s modulus, bonding strength), and the adhesive properties (strength and fracture energy) right before creping, when the temperature is approximately 100 . However, as a future work it is possible to prepare hand-sheet (with controlled wire marks) and creping adhesives, calibrate their properties by designing tensile tests and peel 100  tests in an environmental chamber, then crepe the sheet and perform tensile test to quantitatively validate the models developed in this work.  101  Bibliography  [1]  IBIS World Industry Report 32229a CA, 2017.  [2]  McConnel, W., 2004. The science of creping. Tissue World Americas, Miami Beach, FL.  [3]  Ramasubramanian, M.K., 2001. Physical and mechanical properties of towel and tissue. In Handbook of Physical Testing of Paper, pp. 683-896. .  [4]  Tysén, A., 2014. Through air drying: The influence of formation and pulp type on non-uniform drying and air flo. Ph.D. thesis, Karlstads university. [5]  Oliver, J.F., 1980. Dry-creping of tissue paper-a review of basic factors. Tappi, 63(12), pp.91-95.  [6]  Raunio, J.P. and Ritala, R., 2012. Simulation of creping pattern in tissue paper. Nordic Pulp and Paper Research Journal, 27(2), p.375.  [7]  Patterson, T., 2013. Evaluating and enhancing tissue softness. In Tappi Papercon Conference, Atlanta, GA, Apr.  [8]  Hollmark, H., 1983. Mechanical properties of tissue. Chap, 11, pp.479-519.  [9]  Ramasubramanian, M.K. and Shmagin, D.L., 2000. An experimental investigation of the creping process in low-density paper manufacturing. Journal of Manufacturing Science and Engineering, 122(3), pp.576-581.  [10]  Zou, X., Jang, H.F., Hu, T., Jong, J., Ricard, D., Drolet, F., Allem, R., 2013, “Softwood kraft (NBSK) for tissue makingan overview,” Tissue 360 Forum, TAPPI PaperCon 2013.  102  [11]  Hollmark, H., 1972. Study of the Creping Process on an Experimental Paper Machine. STFI-meddelande, Serie B, (144).  [12]  Boudreau, J. and Barbier, C., 2014. Laboratory creping equipment. Journal of Adhesion Science and Technology, 28(6), pp.561-572.  [13]  Hämäläinen, P., Hallbäck, N. and Barbier, C., 2016. Development and evaluation of a high-speed creping simulator for tissue. Nordic Pulp & Paper Research Journal, 31(3), pp.448-458.  [14]  Ho, J., Hutton, B., Proctor, J. and Batchelor, W., 2007. Development of a tissue creping test rig. Chemeca, Melbourne, Australia, Sept, pp.23-26.  [15]  Ramasubramanian, M.K., Sun, Z. and Chen, G., 2011. A mechanics of materials model for the creping process. Journal of Manufacturing Science and Engineering, 133(5), p.051011.  [16]  Gupta, S.S., 2013. Study of delamination and buckling of paper during the creping process using Finite Element Method-a cohesive element approach. Ph.D. thesis, North Carolina State University.  [17]  Hossain, S., Bergström, P. and Uesaka, T., 2017. Computational design of fibre network by discrete element method. In Fundamental Research Symposium, Oxford, UK.  [18]  Weatherley, D., Boros, V. and Hancock, W., 2011. ESyS-particle tutorial and user’s guide version 2.1. Earth Systems Science Computational Centre, The University of Queensland.  [19]  Evans, A.G. and Hutchinson, J.W., 1984. On the mechanics of delamination and spalling in compressed films. International Journal of Solids and Structures, 20(5), pp.455-466.  103  [20]  Gioia, G. and Ortiz, M., 1997. Delamination of compressed thin films. Advances in Applied Mechanics, 33(08), pp.119-192.  [21]  Freund, L.B. and Suresh, S., 2004. Thin film materials: stress, defect formation and surface evolution. Cambridge University Press.  [22]  Evans, A.G., Mumm, D.R., Hutchinson, J.W., Meier, G.H. and Pettit, F.S., 2001. Mechanisms controlling the durability of thermal barrier coatings. Progress in Materials Science, 46(5), pp.505-553.  [23]  Bao, G. and Cai, H., 1997. Delamination cracking in functionally graded coating/metal substrate systems. Acta Materialia, 45(3), pp.1055-1066.  [24]  Vella, D., Bico, J., Boudaoud, A., Roman, B. and Reis, P.M., 2009. The macroscopic delamination of thin films from elastic substrates. Proceedings of the National Academy of Sciences, 106(27), pp.10901-10906.  [25]  Rogers, J.A., Someya, T. and Huang, Y., 2010. Materials and mechanics for stretchable electronics. Science, 327(5973), pp.1603-1607.  [26]  Su, Y., Wu, J., Fan, Z., Hwang, K.C., Song, J., Huang, Y. and Rogers, J.A., 2012. Postbuckling analysis and its application to stretchable electronics. Journal of the Mechanics and Physics of Solids, 60(3), pp.487-508.  [27]  Mei, Y., Kiravittaya, S., Harazim, S. and Schmidt, O.G., 2010. Principles and applications of micro and nanoscale wrinkles. Materials Science and Engineering: R: Reports, 70(3-6), pp.209-224.  104  [28]  Ruffini, A., Durinck, J., Colin, J., Coupeau, C. and Grilhé, J., 2012. Effects of sliding on interface delamination during thin film buckling. Scripta Materialia, 67(2), pp.157-160.  [29]  Pan, K., Ni, Y. and He, L., 2013. Effects of interface sliding on the formation of telephone cord buckles. Physical Review E, 88(6), p.062405.  [30]  Ni, Y. and Soh, A.K., 2014. On the growth of buckle-delamination pattern in compressed anisotropic thin films. Acta Materialia, 69, pp.37-46.  [31]  Matuda, N., Baba, S. and Kinbara, A., 1981. Internal stress, Young's modulus and adhesion energy of carbon films on glass substrates. Thin Solid Films, 81(4), pp.301-305.  [32]  Colin, J., Cleymand, F., Coupeau, C. and Grilhé, J., 2000. Worm-like delamination patterns of thin stainless steel films on polycarbonate substrates. Philosophical Magazine A, 80(11), pp.2559-2565.  [33]  Moon, M.W., Chung, J.W., Lee, K.R., Oh, K.H., Wang, R. and Evans, A.G., 2002. An experimental study of the influence of imperfections on the buckling of compressed thin films. Acta Materialia, 50(5), pp.1219-1227.  [34]  Faou, J.Y., Parry, G., Grachev, S. and Barthel, E., 2015. Telephone cord buckles—a relation between wavelength and adhesion. Journal of the Mechanics and Physics of Solids, 75, pp.93-103.  [35]  Ni, Y., Yu, S., Jiang, H. and He, L., 2017. The shape of telephone cord blisters. Nature Communications, 8, p.14138.  105  [36]  Yu, S.J., Xiao, X.F., Chen, M.G., Zhou, H., Chen, J., Si, P.Z. and Jiao, Z.W., 2014. Morphological selections and dynamical evolutions of buckling patterns in SiAlNx films: from straight-sided to telephone cord or bubble structures. Acta Materialia, 64.  [37]  Cordill, M.J., Bahr, D.F., Moody, N.R. and Gerberich, W.W., 2007. Adhesion measurements using telephone cord buckles. Materials Science and Engineering: A, 443(1-2), pp.150-155.  [38]  Hutchinson, J.W. and Suo, Z., 1991. Mixed mode cracking in layered materials. In Advances in Applied Mechanics, 29, pp. 63-191.  [39]  Hutchinson, J.W., Thouless, M.D. and Liniger, E.G., 1992. Growth and configurational stability of circular, buckling-driven film delaminations. Acta Metallurgica et Materialia, 40(2), pp.295-308.  [40]  Audoly, B., 1999. Stability of straight delamination blisters. Physical Review Letters, 83(20), p.4124.  [41]  Parry, G., Cimetière, A., Coupeau, C., Colin, J. and Grilhé, J., 2006. Stability diagram of unilateral buckling patterns of strip-delaminated films. Physical Review E, 74(6), p.066601.  [42]  Moon, M.W., Lee, K.R., Oh, K.H. and Hutchinson, J.W., 2004. Buckle delamination on patterned substrates. Acta Materialia, 52(10), pp.3151-3159.  [43]  Audoly, B., Roman, B. and Pocheau, A., 2002. Secondary buckling patterns of a thin plate under in-plane compression. The European Physical Journal B-Condensed Matter and Complex Systems, 27(1), pp.7-10.  [44]  Jagla, E.A., 2007. Modeling the buckling and delamination of thin films. Physical Review B, 75(8), p.085405.  106  [45]  Lee, Y.J. and Freund, L.B., 1992. Dynamic buckling delamination of a bonded thin film under residual compression. In Advances in Dynamic Systems and Stability. pp.147-160. Springer, Vienna.  [46]  Gruttmann, F. and Pham, V.D., 2008. A finite element model for the analysis of buckling driven delaminations of thin films on rigid substrates. Computational Mechanics, 41(3), pp.361-370.  [47]  Faou, J.Y., Parry, G., Grachev, S. and Barthel, E., 2012. How does adhesion induce the formation of telephone cord buckles?. Physical Review Letters, 108(11), p.116102.  [48]  Mei, H., Gowrishankar, S., Liechti, K.M. and Huang, R., 2010, June. Initiation and propagation of interfacial delamination in integrated thin-film structures. 2010 12th IEEE Intersociety Conference, pp.1-8.  [49]  Nilsson, K.F. and Giannakopoulos, A.E., 1995. A finite element analysis of configurational stability and finite growth of buckling driven delamination. Journal of the Mechanics and Physics of Solids, 43(12), pp.1983-2021.  [50]  Jensen, H.M. and Sheinman, I., 2002. Numerical analysis of buckling-driven delamination. International Journal of Solids and Structures, 39(13-14), pp.3373-3386.  [51]  Pan, K., Ni, Y., He, L. and Huang, R., 2014. Nonlinear analysis of compressed elastic thin films on elastic substrates: From wrinkling to buckle-delamination. International Journal of Solids and Structures, 51(21-22), pp.3715-3726.  [52]  Munjiza, A.A., Knight, E.E. and Rougier, E., 2011. Computational Mechanics of Discontinua. John Wiley & Sons.  107  [53]  Hoover, W.G., Hoover, C.G., Kum, O., Castillo, V.M., Posch, H.A. and Hess, S., 2006. Smooth particle applied mechanics. Advanced Series in Nonlinear Dynamics, 25.  [54]  Persson, J. and Isaksson, P., 2013. A particle‐based method for mechanical analyses of planar fiber‐based materials. International Journal for Numerical Methods in Engineering, 93(11), pp.1216-1234.  [55]  Persson, J. and Isaksson, P., 2014. A mechanical particle model for analyzing rapid deformations and fracture in 3D fiber materials with ability to handle length effects. International Journal of Solids and Structures, 51(11-12), pp.2244-2251.  [56]  Lavrykov, S., Ramarao, B.V., Lindström, S.B. and Singh, K.M., 2012. 3D network simulations of paper structure. Nordic Pulp & Paper Research Journal, 27(2), pp.256-263.  [57]  Edvardsson, S. and Uesaka, T., 2010. System dynamics of the open-draw with web adhesion: particle approach. Journal of Applied Mechanics, 77(2), p.021009.  [58]  Etzmuss, O., Gross, J. and Strasser, W., 2003. Deriving a particle system from continuum mechanics for the animation of deformable objects. IEEE Transactions on Visualization and Computer Graphics, 9(4), pp.538-550.  [59]  Bergou, M., Wardetzky, M., Robinson, S., Audoly, B. and Grinspun, E., 2008. Discrete elastic rods. ACM Transactions on Graphics (TOG), 27(3), p.63.  [60]  Bergou, M., Audoly, B., Vouga, E., Wardetzky, M. and Grinspun, E., 2010, July. Discrete viscous threads. ACM Transactions on Graphics (TOG), 29(4), p.116.  [61]  Krenk, S., 2009. Non-linear modeling and analysis of solids and structures. Cambridge University Press.  108  [62]  Kulachenko, A. and Uesaka, T., 2012. Direct simulations of fiber network deformation and failure. Mechanics of Materials, 51, pp.1-14.  [63]  Marynowski, K. and Kapitaniak, T., 2002. Kelvin–Voigt versus Bürgers internal damping in modeling of axially moving viscoelastic web. International Journal of Non-Linear Mechanics, 37(7), pp.1147-1161.  [64]  Roylance, D., McElroy, P. and McGarry, F., 1980. Viscoelastic properties of paper. Fibre Science and Technology, 13(6), pp.411-421.  [65]  Page, D.H., Tydeman, P.A. and Hunt, M., 1962. The formation and structure of paper. Tech. Sect. Brit. Paper & Board Makers' Assoc., London, p.397.  [66]  Barenblatt, G.I., 1959. The formation of equilibrium cracks during brittle fracture. General ideas and hydrotheses: Axially-symmetric cracks. J. Appl. Math. Mech., 23, pp.622-636.  [67]  Barenblatt, G.I., 1962. The mathematical theory of equilibrium cracks in brittle fracture. In Advances in Applied Mechanics, 7, pp. 55-129.  [68]  Dugdale, D.S., 1960. Yielding of steel sheets containing slits. Journal of the Mechanics and Physics of Solids, 8(2), pp.100-104.  [69]  Tvergaard, V. and Hutchinson, J.W., 1992. The relation between crack growth resistance and fracture process parameters in elastic-plastic solids. Journal of the Mechanics and Physics of Solids, 40(6), pp.1377-1397.  [70]  Needleman, A., 1987. A continuum model for void nucleation by inclusion debonding. Journal of Applied Mechanics, 54(3), pp.525-531.  109  [71]  Geubelle, P.H. and Baylor, J.S., 1998. Impact-induced delamination of composites: a 2D simulation. Composites Part B: Engineering, 29(5), pp.589-602.  [72]  Xu, C., Siegmund, T. and Ramani, K., 2003. Rate-dependent crack growth in adhesives: I. Modeling approach. International Journal of Adhesion and Adhesives, 23(1), pp.9-13.  [73]  Li, S., Thouless, M.D., Waas, A.M., Schroeder, J.A. and Zavattieri, P.D., 2006. Mixed-mode cohesive-zone models for fracture of an adhesively bonded polymer–matrix composite. Engineering Fracture Mechanics, 73(1), pp.64-78.  [74]  Xie, D. and Waas, A.M., 2006. Discrete cohesive zone model for mixed-mode fracture using finite element analysis. Engineering Fracture Mechanics, 73(13), pp.1783-1796.  [75]  Shet, C. and Chandra, N., 2002. Analysis of energy balance when using cohesive zone models to simulate fracture processes. Journal of Engineering Materials and Technology, 124(4), pp.440-450.  [76]  Domokos, G. and Holmes, P., 1993. Euler's problem, Euler's method, and the standard map; or, the discrete charm of buckling. Journal of Nonlinear Science, 3(1), pp.109-151.  [77]  Cromer, A., 1981. Stable solutions using the Euler approximation. American Journal of Physics, 49(5), pp.455-459.  [78]  Press, W.H., Teukolsky, S.A., Vetterling, W.T. and Flannery, B.P., 2007. Numerical recipes 3rd edition: The art of scientific computing. Cambridge University Press.  [79]  De Moura, C.A. and Kubrusly, C.S., 2013. The Courant-Friedrichs-Lewy (CFL) Condition. Communications on Pure and Applied Mathematics, 10(2), pp.363-371. 110  [80]  Chen, G., 2011. A computational mechanics model for the delamination and buckling of paper during the creping process. Ph.D. thesis, North Carolina State University.  [81]  Hägglund, R. and Isaksson, P., 2006. Analysis of localized failure in low-basis-weight paper. International Journal of Solids and Structures, 43(18-19), pp.5581-5592.  [82]  Castro, J. and Ostoja-Starzewski, M., 2003. Elasto-plasticity of paper. International Journal of Plasticity, 19(12), pp.2083-2098.  [83]  Ramasubramanian, M.K. and Perkins, R.W., 1988. Computer simulation of the uniaxial elastic-plastic behavior of paper. Journal of Engineering Materials and Technology, 110(2), pp.117-123.  [84]  Wilbrink, D.V., Beex, L.A.A. and Peerlings, R.H.J., 2013. A discrete network model for bond failure and frictional sliding in fibrous materials. International Journal of Solids and Structures, 50(9), pp.1354-1363.  [85]  Borodulina, S., Kulachenko, A., Nygårds, M. and Galland, S., 2012. Stress-strain curve of paper revisited. Nordic Pulp & Paper Research Journal, 27(2), pp.318-328.  [86]  Xia, Q.S., Boyce, M.C. and Parks, D.M., 2002. A constitutive model for the anisotropic elastic–plastic deformation of paper and paperboard. International Journal of Solids and Structures, 39(15), pp.4053-4071.  [87]  Li, Y., Stapleton, S.E., Reese, S. and Simon, J.W., 2016. Anisotropic elastic-plastic deformation of paper: In-plane model. International Journal of Solids and Structures, 100, pp.286-296.  111  [88]  Ramasubramanian, M.K. and Wang, Y., 2007. A computational micromechanics constitutive model for the unloading behavior of paper. International Journal of Solids and Structures, 44(22-23), pp.7615-7632.  [89]  Yoshihara, H. and Yoshinobu, M., 2014. Off-axis tensile strength and evaluation of the in-plane shear strength of paper. Holzforschung, 68(5), pp.583-590.  [90]  Haglund, L., Norman, B. and Wahren, D., 1974. Mass-distribution in random sheets-theoretical evaluation and comparison with real sheets. Svensk Papperstidning-Nordisk Cellulosa, 77(10), p.362.  [91]  Haglund, L., Norman, B. and Wahren, D., 1971. The Power Spectrum of the Basis Weight of Random Sheets. Svenska Traforskningsinstitutet.  [92]  Sampson, W.W., 2001. The structural characterisation of fibre networks in papermaking processes—a review. In The Science of Papermaking, Transactions of the XIIth Fundamental Research Symposium, Baker CF pp. 1205-1288.  [93]  Archer, S. and Furman, G., 2005. Embedded sheet structures: impact on tissue properties. Tissue World, Nice, France.  [94]  Hollmark, H. and Ampulski, R.S., 2004. Measurement of tissue paper softness: A literature review. Nordic Pulp & Paper Research Journal, 19(3), pp.345-353.  [95]  Das, R., 2019. Experimental studies on creping and its influence on mechanical properties of tissue paper products. M.Sc. thesis, University of British Columbia.     112  Appendices  Appendix A   Validation of the Elastic Particle Dynamics Model  We use the particle dynamics model to simulate the buckle-delamination process of a thin film/substrate structure under residual compression, and compare the simulation results with the existing analytical solution. A compressed thin film of thickness ℎ and length 𝐿0 is bonded onto a rigid substrate with a pre-existing delamination zone of length 2𝑏0, as shown in Fig. A.1(a). When the residual compressive strain 𝜀0 exceeds the Euler strain 𝜀𝐸0 = 𝜋2ℎ2 12𝑏02⁄ , the film within the initial delaminated part will first buckle, as shown in Fig. A.1(b). After that, if the energy release rate at the delamination front exceeds the interfacial fracture energy, buckling will drive delamination to propagate until the energy release rate becomes equal to the fracture energy (Fig. A.1(c)). The maximum buckling deflection is given by 𝜁 = ℎ√(4 3⁄ )(𝜀0 𝜀E⁄ − 1)  where 𝜀𝐸 = 𝜋2ℎ2 12𝑏2⁄  is the Euler-buckling strain [38]. Since the film is under residual compressive strain 𝜀0, the initial spacing between particles is set as  (1 − 𝜀0)𝑎0 and both ends of the film are fixed. First we check the convergence of the model. Fig. A.2 shows that by increasing the particle numbers the maximum deflection converges to the analytical solution. To satisfy the convergence of the model, in general we need to make sure there are at least 5 particles within the cohesive zone length 𝑙𝑐𝑧 𝐸𝐺𝐼𝐶 𝜎𝑐2⁄ . For the parameters chosen here, 𝑁 > 500 is sufficient. Fig. A.3 shows the maximum deflection of the film at equilibrium state under different mixed-mode ratios. The mode I fracture energy is fixed as 𝐺𝐼𝐶 = 0.8𝛤0  where 𝛤0 = 𝐸ℎ𝜀02 2 is the initial strain energy density in the film. Our numerical results agree well with the analytical solution. 113      Figure A.1 Schematic of the buckle-delamination process of thin film/substrate system under residual compression with an initial delamination zone. (a) Initial compressed state. (b) Buckling within pre-existing delamination zone. (c) Buckling driven delamination.   Figure A.2 Convergence check of the model. The normalized maximum deflection of the film converges to a constant value as the number of particles increases. The parameters are chosen as 𝜺𝟎 = 𝟎. 𝟎𝟐, 𝒃𝟎 = 𝟏𝟎𝒉, 𝑳𝟎 =𝟏𝟎𝟎𝒉, 𝑮𝑰𝑪 = 𝟎. 𝟖𝜞𝟎, 𝑮𝑰𝑰𝑪 = 𝟏𝟎𝑮𝑰𝑪, 𝜞𝟎 = 𝑬𝒉𝜺𝟎𝟐 𝟐.  𝜀0(a)(b)(c)𝜁2𝑏02𝑏Thin filmSubstrate2𝑏0114   Figure A.3 The normalized maximum deflection of the film in the quasi-static buckle-delamination process as a function of the mixed-mode interfacial fracture energy. The solid line represents the analytical solution 𝜻 =𝒉√(𝟒 𝟑⁄ )(𝜺𝟎 𝜺𝐄⁄ − 𝟏).               115  Appendix B  Validation of the Elasto-Plastic Beam Bending Model          We use the extended model developed in Chapter 3 to simulate the deformation of an elastoplastic cantilever beam under an end load. Fig. B.1 shows the schematic of the cantilever beam, where 𝑃 represents the end load and 𝛿 represents the vertical deflection. We compare our simulation against the results obtained by Abaqus, a commercial software based on Finite Element Method. A rate-independent bilinear elastoplasticity model was used in Abaqus to describe the material properties. The Abaqus simulation was quasi-static. In our discrete model, the simulation was also performed in a quasi-static manner: for every 1 million time steps (∆t), during the first half million time steps we increase the end deflection δ by ∆δ, and during the second half we fix the end deflection and let the system evolve to the equilibrium state. At the end of each 1 million time steps, we record the end deflection δ and the corresponding end load 𝑃. We repeat this process until the whole simulation is completed. The parameters used in the discrete particle model are given as: ℎ = 1 𝑚𝑚 , 𝑤 = 1𝑚 , 𝐿 = 12.5 𝑚𝑚 , 𝑁 = 201 , 𝐸1 = 10 𝐺𝑃𝑎 , 𝐸2 = 2 𝐺𝑃𝑎 , 𝜎𝑦 =10 𝑃𝑎, 𝜀𝑦 = 0.001, 𝜌 = 3000 𝑘𝑔 𝑚3⁄ ,  = 20, ∆t=2 ns and ∆δ=h /100. Fig. B.2 show our simulation matches very well with the results obtained from Abaqus and they both agree with the analytical solution before yielding occurs. The remaining discrepancy is because the beam model used in Abaqus considers the shear stress, while our discrete model did not.  Figure B.1 Schematic of a cantilever beam under end load.  𝑃 𝛿116   Figure B.2 Load-deflection curve of the cantilever beam. The analytical solution before yielding is given by 𝑷 =𝟑𝑬𝑰𝜹 𝑳𝟑.  

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:
https://iiif.library.ubc.ca/presentation/dsp.24.1-0377758/manifest

Comment

Related Items