Probabilistic invertible neural network for inverse design space exploration and reasoning

Abbreviations

e is the creasing rule deviation and it was estimated as 0, 0.2, 0.4 mm; ZDTT is the z-directional (out-of-plane) tensile test; MD uses machine direction for paper making; t is the thickness of the worksheet (paperboard) and its average is 0.43 mm; e_z = x/t is the elongation in ZDTT; σ is the tensile stress in ZDTT; a₁–a₄, b₁–b₄, and c₁ are the coefficients of Eqs 1–3; f_ZDTT = σ·L_ZDTT is the tensile line force in ZDTT; σ_Y is the yield stress of the specimen in the in-plane MD; U is the distance between the 1st end and 2nd end of the nonlinear spring which is joined to nodes on the peeled layer; K is the stiffness of the nonlinear spring in USPRNG subroutine; V is the feed velocity of the creasing rule in the experiment and it was 0.0167 mm·s⁻¹; τ_{B(in MD)} is the breaking shear strength in the in-plane MD shearing test, estimated as τ_{B(in MD)} = 1.63 MPa; s_t is the shear strength parameter in the glue contact function; L_S is the length of the specimen in the experiment and length was 60 mm; W_S is the width of the specimen in the experiment and width was 15 mm; L_f is the length of the specimen for 180° folding in the simulation model and the length was 20 mm; B is the width of channel die and width of channel die was 1.5 mm; d is the indentation depth of creasing rule; d_as is the permanent depth after scoring by the creasing rule; h is the height difference (step) of rubber from the creasing rule and it was 1.4 mm; r is the radius of the creasing rule and radius of creasing rule was 0.355 mm; b is the thickness of a creasing rule, and it was 0.71 mm; γ = 2d·B⁻¹ = tanδ is the normalized indentation depth (nominal shear strain); θ is the folding angle (°); c_d is the crease deviation and it was estimated as c_d = |e₁–e₂| as illustrated in Figure 5b; ω is the rotation velocity of the fixture for folding (revolution per second (rps), s⁻¹) ( = 2πω rad∙s⁻¹⁾ and was set to 0.2 rps; c_b is the thickness of 1st layer of paperboard, and it was 0.11 mm ( = t_pZDTT); v is the Poisson's ratio of paperboard and it was assumed to be 0.2; μ_b is the friction coefficient between the creasing rule and paperboard and it was assumed to be 0.1; μ_d is the friction coefficient between the grooved plate (die) and paperboard and it was assumed to be 0.1; μ_r is the friction coefficient between the paperboard and rubber fixtures and it was assumed to be 0 (no friction); μ_el is the friction coefficient between each layer (interface layer) of the scoring area and it was assumed to be 0.7; φ is the crease direction and as φ = 90° with respect to the MD; t_f is the thickness of folding profile; h_b is the height of bulged profile; V_r is the feed velocity of the rolling process in the simulation model and was 0.1 mm·s⁻¹; P_p is the pushing plate for feeding the deformable body into the rolling process; μ_sc is the friction coefficient of self-contact for deformable body and it was assumed to be 0.7; μ_rd is the friction coefficient between rollers and paperboard, and it was assumed to be 0.1.

1.   Introduction

Paperboard is popularly used in the packaging industry due to its superior features such as a high strength-to-weight ratio, high surface smoothness, lightweight material, printability, cheap prices, sustainability, and recyclability.

The creasing and folding performance, the relative tensile strength and creasing strength, the in-plane compressive strength, the effect of multiple creasing lines, the optimization of the number of creases, a large-scale forming simulation, the out-of-plane buckles of the panel, the post-buckling, the in-plane tensile forces, the out-of-plane bending moment, the deformed shapes of the boxes, the behavior of corner panels, the failure load of the box, the positioning accuracy of creases, the box resistance concerning vertical compression, the width effect on the crease geometry, and the failure evaluation related to the effect of creases and friction of the paperboard packages have been studied, respectively ^{[1,2,3,4,5,6,7,8,9,10,11,12,13]}. The cushioning effect of the corrugated structure was studied by Gu et al. ^[14]. The relationship between the detaching resistance of interlayers and the bulging process was experimentally studied from the aspect of the behavior of de-laminated layers ^{[15,16,17,18,19,20,21]}. The anisotropic mechanical properties of paperboard and delamination-based folding resistance in a creasing process were analyzed by Stenberg et al. ^[17]. The in-plane and out-of-plane mechanical behavior of paperboard has been experimentally studied by Li et al. ^[22,23]. The creasing process in the real paperboard converting industry is performed under very different in-plane and out-of-plane properties to accomplish a smooth and adequately folded notch that allows forming and folding without surface cracking ^{[16,18,24,25,26]}. The shifting effect of the crease line on the local compressive strength of boxes with offset flaps and the numerical model development of boxes in compression tests were discussed ^[27,28]. SONMEZ et al. ^[29] evaluated the binder ratio in pigment coating medley on the crease ability of coated paperboard, and the stiffness and the cracking values were investigated. Hine ^[24] evaluated the local damage in the out-of-plane paperboard properties to create a good folding line before the folding process.

To make a smart folding of a certain thick paperboard or a laminated structure sheet of dissimilar materials is empirically processed using two stages of pre-creasing and folding. Nagasawa et al. ^[30] explained the mechanism of two stages of folding as shown in Figure 1. Here the damaged zone, the in-plane state, and the bulging by the in-plane compressing state were demonstrated.

Regarding the creasing and folding process, Hine ^[31] investigated the relationship between crease depth and crease width. Huang et al. ^[32] reported that the folding behavior was influenced by the interface strength, while the creasing force was affected by the ply properties. Leminen et al. ^[33] compared the creasing and scoring effects on the performance in the folding process and on the strength of manufactured packages. Carey ^[34] and Carlson et al. ^[35] clarified that the bending behavior of the creased part was comprised of two parallel beams. The forming force optimization, folding resistance, crease dimension, and level of spring back were analyzed in the forming process ^[36]. Nygårds and Sundström ^[37] compared in-plane compression and bending failure of paperboard. To control the in-plane compressive properties in the forming process, a large bending of a creased paperboard was experimentally studied in a range of folding angles 0°–180° ^[38]. However, the misalignment of the creasing rule versus the groove of the counter face plate, the inside bulge and its internal breaking were not discussed, except for the articles of Nagasawa et al. ^[18,20,30]. According to the usefulness of the tapered groove, the rectangular groove counter plate (0° tapered) was compared with that of the 45° tapered groove counter plate, and the bending characteristics of creased bulging were experimentally investigated ^[30].

The deviation error in the creased line processing of paperboard between a groove of the counter plate and a creasing rule is important for making a good product. Therefore, adept's check verification of deviation error is essential, and a long maintenance time is necessary to reduce this misalignment empirically. These processes were observationally conducted by experts in the past.

To make a stable bulging form that was controlled by the rule deviation, the correlation between the rule deviation and crease deviation was experimentally studied by Nagasawa et al. ^[18,20]. The rule deviation and crease deviation were discussed experimentally in a limited range until the right-angle folding. Nagasawa et al. ^[30] reported the deviation at 180° folding. This method (180° folding) is quick and easy to evaluate the deviated quantity c_d, compared to the previous method (measurement at the 90° folding). All of them were discussed with respect to the normalized indentation depth γ = 2d/B = 2t/B(d/t), but numerical or theoretical estimation was not shown. Hence, the relationship between e and c_d is not known theoretically. Also, as an important operation parameter to characterize the bulging deformation of the crease, the permanent indentation depth (scored depth) d_as was not discussed, although the relation between d and d_as was recognized as the spring back ^[30].

The simulation techniques and cohesion zone model of the delamination layers were discussed with respect to the joint strength of the laminated composite structure ^[39]. Beex and Peerlings ^[40] used the friction model and cohesive zone model to describe the delamination behavior of paperboard. Alam et al. ^[41] used a simulation model to describe the relationship between bending stiffness and fold-crack resistance in coated papers. Carlsson et al. ^[42] used the J-integral and finite element method to study the creasing process. The cohesive crack propagation and the in-plane properties of the paperboard have been studied ^{[43,44,45,46]}. Ortiz and Pandolfi ^[47] proposed a delamination model concerned with the combination of a cohesive friction model. The characterized cohesive properties were examined with the delamination patterns. Park et al. ^[48] explained the delamination patterns by using the mixed-mode bilinear cohesive zone model, and the damaged region after creasing was reported ^[35,49]. Biel et al. ^[50] developed the mixed mode cohesive law, the fracture resistance for normal (modes Ⅰ), and out-of-plane shear (modes Ⅲ) interface properties of paperboard were obtained. Confalonieri et al. ^[51] validated the interface cohesive models for mixed-mode Ⅰ (normal) and mode Ⅱ (in-plane shear) delamination for various composite materials with variable mode-ratio. The cohesive models as mentioned above was well predicted for the failure problem of composite materials. Nygårds ^[52] determined the failure criterion during bending and folding for paperboard, the various paperboard structures behavior was conducted during bending through the simulations. Simulations and experiments of paperboard were performed to determine the properties of top, middle, and bottom plies ^{[53,54,55,56,57]}. The collapse of the paperboard was affected by the delamination of the paperboard and buckling of individual plies ^[58]. Sudo et al. ^[59] studied the possibility of non-linear spring elements to explain the fluffing behavior of a de-laminated zone in a creased paperboard, although the folding process was not adequately behave without the shear resistance. Beex et al. ^[60], Nygårds et al. ^[61], Giampieri et al. ^[62], Huang et al. ^[53] and Jina et al. ^[63,64] used numerical models concerning the de-lamination mechanism and bulging deformation.

Nagasawa et al. experimentally investigated the crease deviation for the 180° folding test ^[30], the delamination layer of the inside layer bulged, and its breaking of the folded part was not sufficiently discussed because its purpose was the effects of the tapered groove on the self-centering of creasing position, not the deviation mechanism. The delamination phenomena and its peeling resistance were studied as a simulation model by Jina et al. ^[63]. The out-of-plane normal detaching force per unit length (line force) was named as the fluffing normal resistance, which was derived from the z-directional tensile test (ZDTT). The out-of-plane tensile stress was evaluated as the breaking resistance between the lower and upper layers. Because the paperboard is generally produced as multiple-layered structure, and each layer's bonding strength is relatively/empirically weaker than that of each layer (ply). This kind of detaching resistance is recently known as the cohesive zone model ^[50,64]. As for the delamination of paperboards, the breaking resistance between the upper and lower layers (interlayers) is composed of two parts: the out-of-plane normal resistance and the in-plane shear resistance. The breaking at the interlayers is characterized by those two components.

Jina et al. ^[64] used a combination model consisting of the fluffing normal resistance and the glue-breaking resistance as the shear-yield. The former was measured by the z-directional tensile testing (ZDTT), while the latter was measured as the in-plane shearing stress (based on the frictional shear sliding and the in-plane shear test) for discussing the creasing conditions. The combination model of detaching resistance is necessary to evaluate the creasing deformation of paperboard. The detaching resistance (line force) between the delaminated paper plies was originally studied by Sudo et al. ^[59]. The detaching force was recognized as the non-linear load response, and this resistance was seen during the fluffing state of paper fibers ^[63]. Jina et al. introduced the non-linear spring element (known as a special function in the MSC.MARC code) to describe this node connection ^[63], for each node of interlayers. Although Jina et al. defined its detaching resistance as the fluffing resistance, this kind of non-linear overshooting resistance was recognized as the cohesive zone model ^[50,64]. Moreover, the simulated model based on the fluffing resistance was convenient for defining the real non-linear resistance which was the same as the ZDTT.

Regarding the simulation of creasing, Jina et al. ^[64] studied a combination model of the out-of-plane fluffing resistance plus in-plane shear breaking model (as a special cohesive zone model) considering a shallow folding angle until the right-angle folding ^[63,64]. But Jina et al. did not touch the large angle over the right-angle up to 180° folding, and the behavior of crease deviation in simulation at the 180°'s folding. Therefore, the simulation-based theoretical behavior of deviation of the creased line was considered in this study. Nonetheless, focusing on the influence of the eccentric deviation of the creasing rule on the rectangular groove counter plate (rectangle groove) when varying the permanent depth d_as, the bending characteristics of the creased specimen by the experimental result were compared with that of the numerical simulation model.

In this study, the combination model, consisting of the fluffing normal resistance (ZDTT) and the glue breaking in-plane shear resistance, was implemented and the following items were reviewed: (1) The experimental conditions and its results were reappraised, and crease deviation e versus c_d when changing the indentation depth with respect to d_as ^[30]. (2) The simulation model was considered in the previous studies by Jina et al. ^[63,64]. Therefore, the previous simulation model was applied to a new condition, the cantilever type (CST-J1) folding was changed to the symmetric bending of the creased part, and a press on both sides from the out-of-plane was compared to the experimental result. (3) The rolling press was newly simulated to confirm the folding profile. Simultaneously, the bulging and its internal breaking were numerically analyzed. Here, the deformation characteristics of the bulged profile under the rolling process are discussed through the experimental and simulation results.

2.   Specimens' material properties and experimental scoring (pre-creasing) procedure of specimens

To discuss the simulation behavior of deviated crease, the experimental behavior of a white-coated paperboard is first confirmed ^[30]. A commercially recycled coated paperboard is composed of a pulp fiber structure matrix and a clay-coated layer. The fiber layer comprises multiple plies, and the coated layer is a mixture of ground calcium carbonate, kaolin, and a binder ^[63,64]. A coated paperboard with a thickness t = 0.43 (0.42–0.44) mm and a nominal basis weight of 350 g·m⁻² was selected. The in-plane tensile properties of the coated paperboard in the machine direction (MD) for paper making are shown in Table 1 ^[64].

In the scoring (pre-creasing) process, a specimen was scored using rubber blocks, which had a hardness of 40 HS(A) as illustrated in Figure 2 ^[30]. Figure 3 illustrates a form of the specimen and a scoring direction of cross angle θ = 90° against the MD of paper making process. Figure 4 shows a scored state (pre-creasing) of a specimen using a round-edge creaser knife (a creaser with a radius of r = 0.355 mm, and a thickness of b = 0.71 mm). Here, the creaser knife was set up as an eccentric state of e and a permanent depth d_as was generated from the indented depth d against a specimen. At the beginning step of this experiment, the rule deviation (misalignment) e was cautiously modulated. In the scoring process, the creaser knife was indented to a specimen with a depth d. Here, the expression: tanδ = (2d·B^–1) = γ is the normalized indented depth, which is a nominal shear strain for γ $\ll$ 1) ^[15]. B is the width of the groove and H is the height of the groove. The thickness of the specimen was t = 0.43 mm, the thickness of the creasing rule was b = 0.71 mm, the groove width was empirically chosen as 2t + b ≈ B = 1.5 mm, and indentation depth d was chosen as d < H. After making the scoring process, the paperboard specimen had a permanent depth d_as. The permanent depth d_as was less than the indentation depth d because of the spring back effect ^[30,64]. When the rule deviation is positive e > 0 in Figure 4, the left side clearance of the creasing rule was larger than that of the right side clearance. Here, the left and right side dents were called the wide side and narrow side, respectively.

In this study, the rule deviation was investigated mainly for the narrow side of the groove edge. Figure 4 shows the outline of a creasing rule versus a rectangular-grooved counter face plate (α = 0° tapered). The lower surface of the paperboard appears to be heavily damaged on the narrow side. Figure 4a shows the rule deviation (misalignment) from the centerline. The damage zone on the narrow side is illustrated in Figure 4b. When indenting the deviated creasing knife against the channel die (groove), there is a spring effect on the lateral deflection of the creasing knife.

The difference between the indentation depth d and the permanent depth d_as was confirmed later to compare the scored depth between the experiment and the simulation.

The misalignment of the creasing knife with the groove was set up as e = 0.0, 0.2, and 0.4 mm, using shimming sheets as shown in Figure 2. The normalized indentation depth γ was chosen as 0.6 and 0.8 (d = 0.45 and 0.6 mm, respectively). In the scoring process, the feed velocity of the creasing rule was chosen as V = 0.0167 mm·s⁻¹.

All paperboard specimens were maintained at a temperature of 296 ± 1 K and a relative humidity of 50% ± 1% in a controlled room for 24 h. The creasing experiments were performed in the same room.

In the previous experiments ^[30,64], the folding process of scored zone was performed by hands from a scored state to a 180° folded state. An intermediate folded state of the specimen was observed using a digital microscope.

After passing the scoring process by varying the creasing knife deviation e and the indentation depth d with the rectangular-grooved counter face plate (α = 0° tapered), the following details were performed: (1) To evaluate the effects of crease deviation c_d and permanent scored depth in the scored zone before the folding process, the surface dent profile was measured. (2) To prepare a fully folded state of θ = 180°, specimens were folded by hand to make a hinge loose. After the creased part of the specimen was folded up to θ = 180°, the folded hinge was inserted into a rolling device composed of two hard rubber rollers, which had a diameter of 20 mm, as illustrated in Figure 5a. The pressing force of each spring was 14.8 N on both sides of the bearings. Here, the length of the rubber roller was 170 mm. The folded specimen was passed through the rubber rollers, while releasing the opposite sides of folded specimen. Then, the crushed hinge was perfectly formed. The side view of the crushed hinge was recorded by a digital microscope, and the crease deviation was estimated as c_d = |e₁–e₂|, as shown in Figure 5b. Assuming that the scored position moves on the surface of the folded zone with an arc length ⊿, namely, e = Δ, then c_d becomes 2⊿. The total height of a 180° folded specimen is expressed as t_f = 4t − 2t_out (>2t) and e₁ = t_f/2 + ⊿, e₂ = t_f/2–⊿. Here, t_out ($\approx$0.11–0.25 mm) was the thickness of the outside layer pulled in the in-plane. Nagasawa et al. ^[18] reported a crease deviation at the right-angle folding. It was a little different from the definition of c_d. The experiment was performed under the same conditions as the case of the pre-process (scoring process) experiment ^[30]. The side view of the 180° folded state and the relationship between e and c_d was shown later to compare the deformation with the simulation.

3.   Simulation model concerning deviated scored paperboard

To numerically perform the folding process of deviated scored paperboard, a finite element method analysis was carried out, using the authors' developed model. (1) The out-of-plane detaching resistance of the interlayers and the in-plane shear breaking glue resistance of the interlayers were combined to describe the deformation of a multiple-laminated paperboard (8 layers). This combination was developed in a previous study ^[64] for observing a folding angle of 90°. In this study, we considered the simulation of the deep folding up to 180°. (2) In the scoring simulation of the paperboard, the deviated scoring of the round-edge creasing knife was considered. Here, the permanent scored depth was compared with the experimental results. (3) To make the folding profile stable and easy to measure the arc length on the folded surface, a rolling press of the folded zone was carried out in the simulation, and considered in the experiment. To compare the folded and bulged inside profile with the experimental results, the folded state before rolling and after rolling were compared from the aspect of delamination behavior.

A general purpose finite element code, MSC.MARC was employed for simulating the scoring, folding, and rolling process. The updated Lagrange method and a large strain state were employed for analyzing a two-dimensional model (plane strain).

3.1.   Mechanical properties of paperboard

The number of plies (delaminated layers) was assumed to be 8, including the top clay-coated layer in the simulation model. In this simulation model, the number of plies was the same as that reported by Jina et al. ^[64]. The upper layer imposed as the first layer had a thickness of c_b = 0.11 mm ^[63,64], while the other layers were imposed as the interlayers of 2nd–8th layers of 0.046 mm, and the mechanical properties (Young's modulus, Poisson ratio, Yield stress, Shear yield stress, and friction coefficient) and each layer's thickness are shown in Table 2 ^[64]. In the folding deformation, an isotropic elasto-plastic model was assumed by the in-plane tensile testing properties in the MD (Table 1). The deformable-meshed body was considered a full range model because the deformation profile of the folded paperboard was not symmetric with the scored position.

Although this isotropic model is fairly different in terms of the mechanical properties of the thickness direction since the geometrical profile was matched to the real experimental scored depth, the folding resistance was explained by the in-plane stiffness and resistance ^[64]. The material properties and interlayer's friction were based on the previous research ^[64].

3.2.   Out-of-plane (z-directional) detaching resistance of specimens based on z-directional tensile test (ZDTT)

The interface strength of each layer (multiple plies) of coated paperboard was evaluated as the pulling resistance, which was measured by the z-directional tensile test (ZDTT, out-of-plane detaching resistance in the thickness direction). Jina et al. ^[63,64] analyzed the detaching resistance of a weakly bonded layer using the ZDTT. The relationship between the nominal tensile stress and elongation in the thickness direction is indicated in Figure 6. Here, a nonlinear spring model was developed and estimated for describing the detaching resistance of coated paperboard using Eqs 1–3. The first zone of 0 < e_z < e_z1 illustrates the elastic or elasto-plastic behavior before breaking at the weak-bonded layer. A decreasing resistance was observed for the second and third periods of e_z1 < e_z < e_z3. Since it was caused by a fluffing of fibers, the detaching resistance of the interlayer was called a fluffing model. The tensile stress at the first zone of 0 < e_z < e_z1 and the second and third periods of e_z1 < e_z < e_z3 was performed due to two reasons: (1) The tensile state increased to the maximum stress due to the elasto-plastic uniform elongation without drawing and detaching of fibers when the 0 < e_z < e_z1. (2) After reaching the maximum tensile stress, the drawing and detaching phenomena occurred in fiber plies and the tensile stress decreased. At the same time, the curve is continuous, but there is a maximum peak. The peak is caused by the breaking of fibers or layers. The peak response should be discontinuous if the body is not fibrous. But as the body is fibrous, the breaking response was continuous, and residual reducing response occurred when the e_z1 < e_z < e_z3.

A user-defined subroutine of USPRNG (MSC.MARC code) describes the relationship between the nominal tensile stress σ (MPa) and normalized elongation e_z ( = x/t). The stiffness coefficients in Eqs 1–3 are shown in Table 3 ^[64]. This breaking criteria based on the ZDTT was implemented in the simulation model for describing the resistance of detaching layers during the scoring, folding, and rolling process of a creased paperboard.

3.3.   Estimation of in-plane shear strength resistance

The glue breaking in-plane shear yield is performed by the MSC.MARC user's function CONTACT. In the following, this function is called as the glue breaking in-plane shear yield, which is determined by Eq 5. The glued contact is generally detached when Eq 4 is performed by the combination of normal yield stress s_n and shear yield stress s_t. Here, σ_if and τ_if are the contact normal and shear stress, respectively. Generally, the breakage of a plastic body can be defined by the parameters s_n, s_t and m (MSC software, 2010a ^[65]). In the previous study ^[64], Jina et al. assumed to be σ_if = 0 and m = 1 at the delaminated zone but the fluffing model was used. Namely, Eq 5 was applied to the interlayers and also the fluffing model was used. In this study, the same model was considered. Under the delaminated layer, any node is released due to the breaking criteria of Eq 5, the node obeys the rule of frictional CONTACT (MSC.MARC code) between the upper and lower layers of two bodies. The contact stress is calculated using the contact force divided by equivalent areas for shell elements. The glue breaking criteria was fundamentally defined by the in-plane shear stress of the shearing test.

In the 90° folding deformation by Jina et al. ^[64], the in-plane shear-yielding stress was a primary factor to characterize the bulging profile. The bulging profile at the early stage of folding was controlled by a certain level of shear yielding stress s_t. Hence, the in-plane shear-yielding test was investigated for considering the breaking shear yielding stress. Jina et al. explained the upper/lower bound of shear strength in the simulation of the 90° folding ^[64]. The former (lower bound) is the in-plane shear breaking resistance without any compressive interference, while the latter (upper bound) is the in-plane frictional resistance or the shear yielding resistance under a high compressive state. Regarding the first assumption (lower bound) value in the simulation model, the glue model of the detached interlayer was assumed to be s_t = 1.63 MPa ^[64]. However, this value s_t = 1.63 MPa was not matched with the experimental result of folding, while the shearing yield strength (upper bound) was estimated as σ_Y/1.732 = 15.7 MPa from Table 1 and the upper bound model was too strong to detach the interlayer. As a result, the appropriate intermediate value of 1.63 < s_t < 15.7 MPa was empirically detected as s_t = 11 MPa ^[64] when folding up to 90°. The combined delamination model (composed of the fluffing non-linear spring resistance and the in-plane shear breaking glue resistance) was constructed as the breaking criteria of the bonded interfaces at the scored zone for simulating creasing deformation.

In this study, the combined delamination model of the fluffing spring elements (in the normal direction) and the shear glue strength was considered to simulate a new deep folding deformation of creased paperboard up to 180°, and furthermore to pursue a profile of creased part after rolling by using a skin pass rolling die-set.

3.4.   Condition for scoring

The boundary conditions for scoring and the size of the specimen are shown in Figure 7. The bending axis of the crease was across the MD. Since the scoring and folding conditions were not symmetric with the scored-center position, the subdivided mesh model was constructed with a full space of deformable body. The scored zone was considered de-laminated structure in the early stage and its scored profile was adjusted to the actual experimental permanent-scored depth.

Figure 7a empirically defines the detached area classified into two levels: the green zone has the USPRING and in-plane shear yielding, and the yellow zone has the in-plane shear yielding without the USPRING, Figure 7b illustrates a case of e = 0, whereas Figure 7c was a case of e > 0 (the creasing rule relatively moved to the right side from the original position e = 0), respectively. The longitudinal length of the specimen was assumed to be L_f = 20 mm for the full folding model. As a comparison of these two levels of the USPRING, a special case (e = 0, d_as = 0.24 mm) was considered for simulating with the full-ranged USPRING at the yellow and green zones.

The number of divided elements of the deformable body was 10200, while that of the total nodes was 40246. In the simulation, the deformable body was composed of two parts: (i) Each de-laminated layer was perfectly fixed at the left and right outside, and (ii) at the central creased zone, each layer was bonded by the USPRING joints in Eq 3 and the shear glue of Eq 5.

The scoring tools were composed of three parts: a rigid die (grooved counter plate), a creasing rule (round-edge knife), and both sides' fixtures. The width of the groove (channel die) was B = 1.5 mm, while the edge radius of the groove was 0.1 mm. The lower grooved counter plate (die) was fixed, as illustrated in the vertical and horizontal axes (x, y-axis) in Figure 7b, c. The length of the detached zone of the specimen was empirically chosen as L = 2.4 mm (L/B = 1.6) to investigate the creasing deformation.

According to the experimental result ^[30,64], the indentation depth of the creasing knife d and the permanent depth after scoring d_as was fairly different from each other, as shown in Eq 6 and Table 4, due to the spring back effect. Figure 8 shows an example of Eq 6 when e = 0. Therefore, to compare the experimental result and the simulation, the permanent scoring depth based on the experimental result d_as = 0.15–0.3 mm was prepared in the simulation. Here, d_as was calculated using Eq 6, d = 0.45 mm (middle indentation) and 0.6 mm (deep indentation) for e = 0.0, 0.2, and 0.4 mm.

In the simulation, since the spring back was different from that of the experiment, the corresponding indentation was specified as Eq 7. This is caused by the mismatching of anisotropic material properties.

3.5.   Folding condition

To predict the behavior of folding and bulging of the creased part explained in section 2 ^[64], a simulation of the 180° folding was newly investigated and compared with the experimental result. Figure 9a shows the boundary condition and size definition for 180° folding after the creasing rule and rigid die (grooved counter plate) were removed. Herein, the deviation of the creasing rule against the rigid die (grooved counter plate) e was chosen as 0, 0.2, and 0.4 mm. This folding stage was taken over from the scoring stage explained in section 3.4. The rotating fixtures were made of rigid bodies and are illustrated in Figure 9b. The right side fixture rotated from 0 to 90° on the clockwise (−1.57 radians), whereas the left side fixture rotated on the counterclockwise (+1.57 radians) from 0 to 90°. An incremental angle of 0.0185 radians was considered for each step of the folding angle. This stage was performed with 85 steps of total increment for simulating the folding process.

The folding stage and bulging of the creased part are affected by the variation of e and d_as. The simulation model was compared with related experimental results. Also, the bulging profile at the scored zone was compared with the experimental results. In this study, e = 0.1–0.4 mm was considered.

3.6.   Rolling pass (pressing) condition

To reveal the effects of e and d_as on the folding shape of the creased paperboard, the 180° folding model was pressed using a rolling press machine as shown in Figure 10. The rolling process was simulated under the following conditions: if the top layer thickness t_out = 0.11 mm was assumed to be stretched at the outside without crushing, the bulging thickness t_f was ideally 4t-2t_out = 1.5 mm, but it was empirically assumed to be equal to that of the experimental thickness in the simulation as shown in Table 5 because the rollers were assumed to be rigid bodies. The roller gap r_g in Table 5 was used for the simulation. The friction coefficient between the paperboard and the roller μ_rd was assumed to be 0.1, while the friction coefficient of self-contact with the paperboard μ_sc was 0.7 (Table 2). Here, the rollers and the right side virtual tool P_p were assumed to be rigid bodies. The deformable body of the worksheet and the virtual tool P_p were fixed by the GLUE CONTACT function (Figure 10).

The lower roller R_L rotated up to 360° on the clockwise (−6.28 radians), whereas the upper roller R_U rotated on the counterclockwise (+6.28 radians). An incremental angle of 0.0314 radians was performed for each step of the rolling angles. A total rolling angle of 360° was performed with 200 steps for passing through the roller device.

As for the profile of 180° folded zone, the thickness of folded paperboard t_f and the height of bulged profile h_b were measured after passing through the rollers when varying e and d_as. Additionally, the deformation profile of the 180° folded zone without any rolling process was compared with the previous simulation result and experimental result ^[30].

4.   Results and discussion

4.1.   Simulated folding profile of scored zone

Figure 11 shows the initial scored state when choosing d_as = 0.15–0.3 mm and e = 0–0.4 mm from Table 5. The simulated profile of lower side was almost similar to the previous experiment (Figure 12 ^[30]).

In the previous experiment ^[30], the rule deviation e = 0.0, 0.1, 0.2, 0.3, and 0.4 mm under the normalized indentation depth γ = 0.6 (at d_as = 0.15–0.18 mm, intermediate scoring) and 0.8 (at d_as = 0.24–0.3 mm, deep indentation) was carried out and the side views of the fold zone were observed. The folding simulation was performed at θ = 45°, 90°, 135°, and 180°, in cases of (a–c) d_as = 0.15–0.18 mm, and (d–f) d_as = 0.24–0.3 mm as shown in Figure 13, respectively.

It was found that the folding profile was almost symmetric for 45° < θ < 90° and there was a little right side deviation of bulging for 90° < θ < 180° when e ≤ 0.2 mm, whereas the remarkable asymmetric bulging profile occurred when e > 0.2 mm for 45° < θ < 180°, as shown in Figure 13. Since the ideal clearance of (B-b)/2 was 0.395 mm at e = 0 mm, the ideal clearance became 0 mm at e = 0.4 mm. Hence, the worksheet of 0.43 mm appeared to be sheared with zero clearance when e = 0.4 mm.

Concerning the folding shape at the inside (lower side of each picture), the clearance between the bulged surface and the surfaced clamp of both sides is larger, compared to the experimental behavior ^[15,18]. This seemed to be caused by the assumed isotropic properties (such as the yield stress and the Young's modulus) and mismatching of detached resistance of interlayers. Especially for d_as < 0.2 mm (d_as/t < 0.46), the size of the bulging was fairly different between the experiment ^[15,18] and the simulation, whereas it was similar for d_as > 0.2 mm.

To verify the full ranged USRING model at the yellow and green detached area, a special case of d_as = 0.24 mm and e = 0 mm were simulated and shown in Figure 14. Comparing Figure 13d and Figure 14, it was found that the fluffing model was sensitive to changing the bulging shape.

4.2.   Deformation profile after rolling pass

The size of bulging at the 90° folding was not investigated in the experiment ^[30], due to the difficulty of seamless recording by the digital microscope on the CST-J1 for 90° < θ < 180°. Alternatively, the deviated creased specimens were experimentally investigated by folding up to 180° manual operation and pressed by a rubber roller ^[30]. The simulation of 180° folding was considered using the rolling pass of Figure 10. Figure 15a, d, g, j, m, p show sectional views of folded specimens in the experiment after passing through a pair of rollers as shown in Figure 5. These experimental pictures were arranged from the previous study ^[30]. The intermediate indentation cases of d_as = 0.15, 0.17, 0.18 mm, and the deep indentation cases of d_as = 0.24, 0.25, 0.3 mm are shown. Here, the rule deviation was chosen as e = 0.0, 0.2 and 0.4 mm. Comparing these experimental results, the bulging profile was simulated as a 180° folding plus a rolling process. Figure 15c, f, i, l, o, r shows the simulated deformation profile of a folded specimen before the rolling, while Figure 15b, e, h, k, n, q shows the simulated forms after the rolling press. Considering the effect of USPRING fluffing resistance on the yellow-colored detaching zone (Figure 7a), the special case of e = 0 and d_as = 0.24 mm was simulated using the full ranged USPRING, and its result are shown in Figure 16.

The results of the following features were revealed.

(1) The simulated thickness of t_f was adjusted to the experimental thickness when passing through the rolling press. Under this restriction, the simulated height of bulge h_b has a different tendency (a little larger) from that of the experiment. The experimental h_b (after rolling) increased a little (0.9–1.1 mm) with the scored depth d_as, whereas the simulated h_b (before rolling) was almost similar (1.02–1.14 mm), while that (after rolling) was almost constant and larger than 1.3–1.4 mm.

(2) The spatial delamination at the folding zone in the simulation was larger than that in the experiment. The experimental delamination resistance at the folding zone is maintained in a certain level of residual state and then the spatial delamination was quite restricted to forming the bulged profile. The simulation model assumed that the detaching resistance by the USPRING was completely expired at the yellow-colored zone (Figure 7a), where the creasing knife was indented. Therefore, the simulated h_b was almost independent of the scored depth d_as, and it was larger than that of the experiment.

(3) Synthetically, in the developed simulation model, the simulated cases of deep indentation is similar to that of the experimental, while the simulated cases of intermediate indentation are different from that of the experimental, due to a mismatching of the detaching resistance behavior.

(4) Concerning the movement of scored position, the experimental crease deviation was arranged from a previous study ^[30], while the simulated crease deviation in the surface layer was detected in the simulation. Those representative scored positions are marked with a round open in Figure 15.

(5) Using the full-ranged USPRING at the yellow-colored detached zone, the middle layer appeared to be restricted and separated, but the two surface layers were extremely separated outward, due to the assumption of isotropic material properties. The difference of USPRING at the yellow-colored zone didn't contribute much to shape the local crushing of the detached zone.

Although details of the folded shape at the inside corners were far different from that of the experiment due to the mismatching between isotropic and anisotropic material properties and that of detaching resistance in the interlayer, it was found that the folding deformation of the laminated structure including the asymmetric bulging by the creaser deviation, was stably simulated for estimating the crease deviation in the 180° folding model. Regarding the bending stiffness, the peel detaching resistance, in-plane shear yielding and the shear friction are basic factors, while the macroscopic bending is characterized by the in-plane tensile yielding. From the simulation, it was found that the macroscopic isotropic elasto-plastic could not follow up the anisotropic properties of in-plane and out-of-plane stiffness. Using the isotropic properties, as an example when seeing in Figure 15a, c, it was confirmed that the mismatching deformation due to the local crushing occurred at the inner outside (left/right sides). The experiment showed the folding corner was perfectly crushed, but the simulation kept a certain space at the bending inside zone, due to the large stiffness. It was indicated that the local crushing of the orthotropic elasto-plastic model was fairly similar to the experimental result when compared to the isotropic elasto-plastic model.

The effect of anisotropic properties of paperboard on the non-linear buckling was seen in the simulation of flat crush test report by Komiyama et al. ^[66]. The crushing deformation at the locally folded zone occurred when the anisotropic elasto-plastic model was considered in a certain range. Then, the proposed model must be carefully considered for expanding to other problems when using the isotropic elasto-plastic model. Synthetically, the proposed model is available for a certain large folding and large deviated state, but not applicable to a small non-zero deviation.

The simulation model of crease deviation for the 180° folding is also available for multi-layered paperboard due to two reasons: (1) The multi-layered resistance is characterized by the peel detaching resistance and in-plane shearing yield resistance, and the friction of shear was assumed in a certain reasonable range based on the experiments. Each layer's stiffness is basically considered as the elastic-bending mode, but in-plane yielding is based on the in-plane tensile mode. All of these elemental models perform the deformation of multiple layers. (2) In some small deviations (in Jina's work ^[64]), the experiment and simulation were well matched. Regarding this model of multi-layers, interlayer's resistance shear yielding, friction and normal yielding resistance were investigated by varying. And reasonable state was detected for each mode (peeling in out-of-plane, shearing in-plane). Using such a basic data model, the expanded state was simulated: a large folding state and also largely deviated state. If the deviated state is not so different from the zero-deviated state from the aspects of convergence stability, this simulation seems to be available. The results confirmed that a large folding plus a large deviation state is well estimated from the proposed model.

In the intermediate process of bending from θ = 0° up to 180°, since the fixing and rotating virtual device supported the bent sheet with a symmetric condition in the simulation, the inside bulge was almost symmetric except for the eccentric shape of the inside bulge. This supporting model was necessary for stably executing the simulation. When a simple cantilever bending condition was applied to the scored paperboard in the simulation, the folding deformation was quite different (more eccentric bending occurred).

4.3.   Crease deviation of scored position by rule deviation

In the experiment ^[30], the relationship between the crease deviation c_d and the rule deviation e after passing the 180° folding process was investigated. Here, the indentation depth before release was d = 0.3 mm (shallow indentation, d_as = 0.08 mm), 0.45 mm (intermediate indentation, d_as = 0.16–0.19 mm), and 0.6 mm (deep indentation, d_as = 0.24–0.3 mm). The experimental relationship of the crease deviation c_d and the rule deviation e is plotted in Figure 17. It was found that c_d tended to be larger than 2e for e < 0.3 when the indentation depth was small, whereas c_d $\approx$ 2e for e > 0.3 when the indentation depth was deep.

The 180° folding simulation was performed under two conditions. First, the scored worksheet was folded without a rolling press. Second, the folded worksheet was pressed by the rolling. In the simulation, the measurement method of the crease deviation c_d was based on the formula in Figure 5b c_d = |e₁–e₂|. This measurement method was the same as that of the experiment. The simulated crease deviation c_d is plotted in Figure 17 for two conditions: before rolling and after rolling press.

A linear approximation of Eq 8 was applied to the simulation data with respect to e = 0.0, 0.2, and 0.4 mm. Its coefficients are shown in Table 6. Comparing these four cases (intermediate and deep indentation; after rolling and before the rolling process), the relationship of c_d and e based on Eq 8 was almost similar to each other. Namely, ${\text{c}}_{\text{d}}\text{}\text{≈}\text{}\text{2}\text{e}$.

As mentioned in the items (1)(2)(3) of 4.2, when the indentation depth was relatively small, the bulged profile of the simulation was fairly different from that of the experiment because of the actual delamination resistance and anisotropy of mechanical properties. The mismatching between c_d–e relation for the shallow indentation (d_as = 0.08 mm) and intermediate indentation (d_as = 0.16–0.19 mm) occurs in these situations. It was found that the real detaching resistance and anisotropic properties accelerate the off-set relation of c_d > 2e from c_d$\approx$2e for 0 < e < 0.3, whereas the relation of c_d$\approx$2e was kept in the simulated condition.

5.   Conclusions

To verify the folding deformation of creased 0.43 mm thick paperboard, an eccentric pre-crease condition (scored at an offset as a deviation e of the creasing knife) was numerically investigated within the folding angle of 180°, while the experimental deformation of 180° folded paperboard was compared. A combination model of the fluffing resistance based on non-linear spring model in the z-directional (out-of-plane) tensile test, and the in-plane shear glue resistance was reviewed. It was applied to a 180° symmetric folding and a light rolling simulation, using the isotropic elasto-plastic solid properties. Through this study, the following was revealed.

(1) The deviation of the creased position at the 180° folding was sufficiently predictable, compared with experimental behavior.

(2) When e = 0.0, 0.2, 0.4 mm, the simulated 180° folding was roughly classified into two deformation modes: in the range of d_as = 0.15–0.18 mm (intermediate scoring), and 0.24–0.3 mm (deep scoring).

(3) In case of the shallow depth d_as, the bottom profile of the bulged zone was numerically similar to a V shape, while in the case of the deep depth d_as, the bottom profile of the bulged crease zone was numerically similar to a U-shape. This tendency of simulation corresponds to the size of the bulge height of a creased hinge.

(4) When considering the rolling press, the thickness of the folding profile t_f was adjusted to the experimental result. In this simulation, the height of the bulged profile h_b and bulged shape of the folded specimen were different from the experimental result owing that the bulging profile of creased part was characterized by the anisotropic elasto-plastic properties and real detaching resistance of paperboard.

(5) The spring back of the simulation was quite different from the real response of the experiment. This was one of the reasons why the simulation was fairly different from the experiment. Alternatively, using the permanent scored depth d_as, the folding deformation of the creased part was simulated, and its deviation behavior was well predicted, compared with the related experiment.

(6) Synthetically, the crease deviation of c_d was numerically characterized by the creasing knife deviation of e. It was c_d $\approx$ 2e, whereas the experimental non-linearity of c_d > 2e appeared to be caused by the anisotropic properties and the real detaching resistance at the intermediate scoring for e < 0.3.

The apparent modulus and the tensile strength tended to decrease with increased temperature while the stretch tended to increase with increased temperature (Wink ^[67]). Although the temperature and humidity variation seem to affect the bending resistance of the creased zone, the geometrical deviation of the creased zone is basically determined by the original deviated state (asymmetric delamination), not by the bending resistance itself. On the other hand, as the variation of material properties seems to affect the buckling behavior of delaminated layers, those parameters (temperature and water content variance) must be discussed furthermore.

Acknowledgments

This work was supported by the Nagaoka University of Technology (grant number 31055), a fund for developing a core of excellence on the paperboard die cutting technology.

Conflict of interest

The authors declare no conflict of interest in this paper.