Finite Element Modeling of Woven Fabric Composites at Meso-Level under Combined Loading Modes

1. Introduction

Woven fabrics are among the most important materials used in today’s modern industries. Next to their high mechanical properties, they are easy to handle in the dry or pre-impregnated pre-forms, offer good drape-ability and are particularly suited for manufacturing of doubly curved components, membranes, inflatable structures, etc (Cavallaro et al., 2003; 2007). In the dry form, fabrics can be formed into a variety of three-dimensional (3D) shapes and then consolidated with resin via resin transfer molding (RTM) or other manufacturing processes (Boisse et al., 2007). Reliable models capable of predicting the mechanical behaviour of woven fabric materials are not fully developed yet. The biggest challenge in this regard is perhaps the multi-scale nature of the fabric materials. Dry fabrics at macro level are composed of numerous yarns interlaced into each other. The yarns usually have characteristic length in the scale of millimetres and their interaction and behaviour at the fabric level can greatly influence the macro-level material behaviour (Guagliano and Riva, 2001). Yarns themselves are heterogeneous media made of bundles of very thin and long fibers. Figure 1 shows different hierarchical levels in a woven fabric along with their typical dimensions.

Figure 1.

Hierarchical levels in woven fabrics

Among different material scales, meso-level modeling of woven fabrics is known to be a strong tool for predicting their effective mechanical properties at macro-level (Peng and Cao, 2002). It can also be useful for studying their local deformation mechanisms that occur during different manufacturing processes and loading conditions. Of different modeling techniques, the 3D finite element modeling is found to be of great interest to the researchers in the field. However, the multi-scale nature of fabrics makes the applicable numerical procedures different from those of the conventional finite element method. The fact that fabric yarns are heterogeneous media formed by bundles of fibres, and that the loose bounding between fibers in each yarn allows them to slide on each other, makes a considerable distinction in postulating the yarns’ constitutive models as well as the numerical procedures applied to analyze their deformation.

(Kawabata et al., 1973; 1973a; 1973b) presented general theories for modeling woven fabric unit cells using bar and stiffness elements. Their sample model was based on the geometrical simplifications on a unit cell along with some parameters that were determined from experiments. Boisse et al., (1997) used the above model and developed a finite element simulation of a dry fabric forming process. Bi-axial tension tests were used to identify the unknown parameters in the constitutive model of the fabric unit cell. In order to develop more accurate models with more insight towards the local deformation phenomena in fibre yarns, Gasser et al. (2000) developed a 3D finite element model of a unit cell under bi-axial tension. Their results were compared to a set of bi-axial tests and satisfactory agreements were obtained. One of the most important features in their approach was to link the meso-level material model to the micro-level behaviour of yarns. For example, quasi-zero shear modulus and Poisson’s ratio, crushing transverse behaviour of yarns, and the update of direction of material orientation during deformation were taken into account. Later on, the model was extended to simulate the in-plane shear behaviour of dry fabrics (Badel et al., 2007) and an algorithm for implementing a hypoelastic constitutive model was presented by Badel et al. (2008 ; 2009). For implementation of these models in numerical packages, the explicit solver of Abaqus has been frequently used. Recently, Komeili and Milani (2010) used a modified version of the aforementioned algorithm to implement an implicit integrator in Abaqus, which led to an increased accuracy and significantly decreased the simulation run time.

Based on the brief review above, it appears that the meso-level finite element modeling of fabrics has been mostly based on individual axial tension and shear modes. Other researchers have also looked at the homogenization of yarn properties at micro/meso levels, but again under individual deformation modes (Chen et al., 2001; Peng and Cao, 2002). Similarly, at macro-level, Xue et al. (2003), followed by Peng and Cao (2005), developed a constitutive material model for the dry fabric sheets. The model was based on a non-orthogonal local coordinate system whose in-plane axis is coincident with the weft and warp yarns of the fabric. To identify the unknown material parameters, the model was fitted to the experimental data from individual bi-axial and bias-extension tests. Nonetheless, during actual forming processes, a complex combination of the axial and shear deformation modes may be experienced by woven fabrics (Boisse, 2010). Cavallaro et al. (2007) developed a new test fixture with the capability of applying simultaneous axial tension and shear deformation modes to the fabric specimens, which could be advantageously used for a more reliable identification of constitutive modes that are used for simulation of composite forming processes.

The aim of the present work is to first present a general meso-level fabric unit cell model using an implicit integrator in Abaqus. To this end, modifications to the original model developed by Badel et al. (2008) are required. Then, the effect of combined loading on the response of a typical fabric unit cell is studied under different axial-shear combined loading modes. The axial loading is induced through controlled displacement/stretch along the yarns and the shear is applied through controlled rotation on the boundaries of the unit cell (i.e., simulating the picture frame test).

2. Modeling

A typical glass plain-weave fabric was selected (Figure 2). Because of its simple textile architecture and balanced properties, this type of fabric has found a wide range of applications in the composite industries.

Figure 2.

A typical balanced plain weave fabric (Boisse, 2010).

2.1. Geometry

As mentioned earlier, the meso-level structure of a woven fabric consists of numerous yarns interlaced into each other to construct the whole fabric structure. In order to model such a complex material system (especially if the goal is to find the equivalent/effective material properties at a fabric level) it may be neither necessary nor computationally feasible to consider all individual yarns and their interactions. Instead, a representative volume element (also known as unit cell) may be considered as a sub-model of the whole fabric structure. Based on a given fabric type and the loading mode, different unit cell models have been employed in the literature (Boisse et al., 2006; Peng et al., 2004). Figure 3a shows the unit cell employed in the present study. The geometrical construction of the yarns in the model is based on the sinusoidal curves shown in Figure 3b and defined via Eqs. (1)-(5).

$$y_1\left(x\right)=\frac{h}{2}\left[cos\frac{\pi x}{s}+1\right] \left(0xs\right)$$ E1

$$y_2\left(x\right)=\frac{h}{2}\left[cos\frac{\pi x}{s}-1\right] \left(0xs\right)$$ E2

$$y_3\left(x\right)=-h cos\frac{\pi x}{\beta } \left(0x\frac{w}{2}\right)$$ E3

$$y_4\left(x\right)=-h cos\left[\frac{\pi \left(x-\left(s-\beta \right)\right)}{\beta }\right] \left(s-\frac{w}{2}xs\right)$$ E4

$$\beta =\frac{\pi w}{2 arcsin\left[si{n}^2\left(\frac{\pi w}{4 s}\right)\right]}$$ E5

Figure 3.

a) Schematic of the unit cell; (b) the yarn generating lines (Mcbride and Chen, 1997); For the current model, $w=2.11 mm,$ $h=0.5$ and $S=5.13$ are used. This means that the total length of the unit cell is $2S=0.26 mm.$

2.2. Material

The multi-scale nature of woven fabrics at meso-level means that the material behaviour of yarns is dependent on the attributes of micro-level fibrous structure. This, in turn, justifies using particular constitutive models of yarns with a close attention to the characteristics of fibers and their interactions at a lower material level. First, fiber yarns cannot have considerable tolerance for shear, compression or bending. This is due to the fact that yarns are made of bundles of thousands of very thin fibers which can slide on each other in the dry form. In addition, the high length to diameter ratio of yarns makes it almost impossible to carry compression without buckling. On the other hand, fibers can go under high tension in the axial direction of yarns. Indeed, the latter property is one of the main reasons for a fabric demonstrating superior mechanical properties. In the transverse direction, however, the yarn behaviour is more intricate. At the initiation of loading, there may be noticeable gaps/voids between the fibers in the cross section of yarns, but with increasing the load they vanish and the fibers begin side-to-side contacts (Figure 4). This phenomenon makes the transverse stiffness of yarns non-linear/strain dependant. It is not straightforward to directly measure a yarn’s transverse stiffness during fabric deformation. Consequently, inverse identification methods along with experimental measurements are commonly used

Figure 4.

The X-ray image from cross section of fibrous yarns (top) before loading; (bottom) after loading (Badel et al., 2008)

to arrive at strain-dependant relationships for the yarns’ transverse stiffness parameters (Gasser et al., 2000; Badel et al., 2008 ).

The general form of material properties in the current model are adapted from (Komeili and Milani, 2010) which were extracted by matching the numerical simulations to the experimental measurements by Buet-Gautier and Boisse (2001) under axial tension and by Cao et al. (2008) under shear loading. The properties used for a simultaneous extension-shear are summarized in the following stiffness matrix:

$$\left[C\right]=\left[\begin{array}{cc}\begin{array}{ccc}{E}_{11}& 0& 0\\ 0& E_{22}& 0\\ 0& 0& E_{33}\end{array}& \begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 0& 0& 0\end{array}\\ \begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 0& 0& 0\end{array}& \begin{array}{ccc}G& 0& 0\\ 0& G& 0\\ 0& 0& G\end{array}\end{array}\right]$$ E6

$$E_{11}\left({\epsilon }_{11}\right)=\{\begin{array}{c}100\text{ MPa }{\epsilon }_{11}1.0\times {10}^{-3}\\ 5\text{ GPa }1.0\times {10}^{-3}\le {\epsilon }_{11}1.6\times {10}^{-3}\\ 50\text{ GPa }1.6\times {10}^{-3}\le {\epsilon }_{11}\end{array}$$ E7

$$E_{tt}=2.5\times {10}^5\left|{\epsilon }_{tt}^{}\right|{\epsilon }_{11}^2+3.0\text{ MPa}$$ E8

$\left[C\right]$ is the stiffness matrix, $E_{11}$ and ${\epsilon }_{11}$ are the axial stiffness and strains; $E_{tt}$ , ${\epsilon }_{\text{tt}}$ , $tt=\left\{\mathrm{22,33}\right\}$ are the transverse stiffness and strains, respectively. $G$ is the shear modulus of the yarns which for dry fabrics should be small compared to the axial and transverse stiffness values. Here $G=60 MPa$ as been selected merely for numerical stability purposes (Gasser et al., 2000); although the shear modulus is at the same order of magnitude as the other two stiffness values in the beginning of loading, it becomes less significant as $E_{11}$ and $E_{22}$ increases with the loading magnitude.

The material model of Eqs (6)-(8) was implemented in the Abaqus finite element software via a UMAT (implicit) user-defined subroutine. In doing so, however, it was noted that the large difference between the stiffness in the yarn axial direction compared to the transverse and shear stiffness values highlights the extreme importance of applying proper material orientation updates during loading steps. The point is that the material properties should be defined in a frame which is rotating with the fiber direction in the yarns. On the other hand, conventional methods in the finite element codes use other (e.g., Green & Naghdi, 1965; Jaumann, 1911) methods for updating the material orientation under large deformation. The problem can be handled with user-defined material subroutines. Subsequently, two approaches may be implemented to ensure that the material properties during stress updates is based on the frame attached to the fibers: (1) Either the stiffness matrix defined along the fiber direction can be transformed to the current working frame of the finite element software, or (2) the stress in the working frame of the software can be transformed to the frame of the fiber and transformed back to the working frame after applying the stress updates in the fiber frame. The details of each method are available in (Badel et al., 2008) and (Komeili and Milani, 2010); the former reference employed an explicit and the latter reference an implicit integrator.

2.3. Periodic boundary conditions

A single isolated unit cell cannot be considered as a good representative of the whole fabric structure unless the effect of adjacent cells is taken into account. In other words, suitable kinematic (or dynamic) conditions should be applied on the perimeter of the unit cell where it is attached to the adjacent cells. These conditions are often called periodic boundary conditions. They are very similar (though different) to symmetric boundary conditions. A thorough discussion on their mathematical details and implementation under individual loading modes is given in (Badel et al., 2007).

The method that has been used in this study is based on the periodic boundary conditions reported in (Peng and Cao, 2002). According to their work, the side surfaces of yarns should remain plane and normal to the unit cell mid surface during deformation. More details of the latter kinematic conditions on unit cells are also given in (Komeili and Milani 2010).

2.4. Loading boundary conditions

There is a variety of test setups used for the axial tension and shear testing of woven fabrics (Buet-Gautier and Boisse, 2001; Cao et al., 2008). On the other hand, experimental setups for the combined loading modes are new and limited. First, it should be defined how a combined loading mode is exerted on a fabric specimen. For example, having a bi-axial load on a fabric where the axial loads does not rotate with the rotation of the yarns and stays parallel to its original direction during deformation, even after the shear load is applied, may be considered a special case of combined loading. As another example, one may consider a combined loading condition where the direction of the axial load rotates and realigns along the yarn direction. For a practical analysis of fabrics, the latter case of stretching in the yarn direction is more important than the former case of stretching yarns along a (fixed) off-axis direction (Boisse 2010). A new test setup capable of applying combined loading in the form of shear and biaxial stretching along the yarns (Figure 5) has been developed in (Cavallaro et al., 2007).

In order to simulate the unit cell of the fabric under such combined loading in the aforementioned Abaqus model, a set of kinematic couplings were applied around the unit cell to satisfy the periodic boundary conditions. Namely, the shear loading has been applied

Figure 5.

Experimental fixture for applying combined shear and axial tension on fabrics (Cavallaro et al., 2007).

via rotation on one of the yarn sides and the rest of unit cell boundaries have linked to follow this movement through a periodic boundary condition. For the axial tension, connector elements between the two corners of each side yarn have been used (they can be seen as solid lines around the unit cell in Figure 6). The connector elements are chosen from the Abaqus library and provide an axial degree-of-freedom between their reference nodes. The axial distance between the nodes can be changed to apply/simulate stretching on the yarns. The reference points are not part of the yarns geometry, but they are kinematically connected to the nodes on the cross sectional surfaces of yarns (i.e., the side surfaces of the unit cell) to implement the periodic and loading boundary conditions. Moreover, there are four reference points on the mid-points of the side lines to impose the kinematic conditions on the middle yarns. The latter reference points are also connected to the corner points by kinematic constraints. Figure 6 shows the aforementioned conditions schematically. Eventually, the material resistance to deformation in the form of reaction moment from the rotation boundary condition and the normal force from the axial connector elements are calculated and reported in the post processing of simulations. They can then be used in the normalized form and compared with experimental results.

Figure 6.

The loading boundary conditions used on the unit cell to model the deformation under a combined loading mode; Circles show the location of reference points.

3. A preliminary validation

In order to validate the model with the existing data in the literature, it is compared to two basic cases where the unit cell is under pure bi-axial tension and shear (Komeili and Milani, 2010). Figure 7 shows the results of these comparisons. In the same figure, a set of actual picture frame test data, collected at the Hong Kong University of Science and Technology (HKUST), is replicated from (Cao et al., 2008). For the axial mode, however, data with the same unit cell geometrical parameters was not available. The differences between the resultant forces and moments in each mode can be related to the type of the unit cell used, shear stiffness of yarns, the method of applying boundary/loading conditions, and other details of the two finite element models in controlling their convergence (e.g. hourglass stiffness, mesh size, etc). In addition, one may redo the inverse identification of the yarn model using the current model. However, as the main goal of this chapter is to highlight the relative effect of combined loading on the mechanical characterization of woven fabrics (i.e., compared to the individual deformation modes), the current model and material properties are used without a loss of generality of the approach.

Figure 7.

A validation of the current model under (a) pure bi-axial and (b) shear mode.

4. The effect of combined loading

In this section the effect of combined loading on the response of the material is analysed, when compared to those obtained from the individual biaxial and shear modes under the same loading magnitude. Figure 8 shows the effect of combined loading on the reaction force in the bi-axial tension and the reaction moment under shear loading. The amount of normalized reaction moment while the fabric is under combined loading has increased up to four times. It has also caused ~12% higher axial reaction force under an identical stretching magnitude.

Figure 8.

The effect of combined loading on the reaction force and moment when compared to the individual (a) shear and (b) biaxial modes. The difference between curves in each graph indicates the presence of additional local deformation phenomena/ interactions between shear and axial modes under combined loading.

The obtained numerical results from bi-axial loading well agree with what has been suggested through experimental measurements in the literature. Namely, Boisse, et al. (2001) and Buet-Gautier and Boisse (2001) argued that the effect of shear strain on the axial behaviour of plain fabrics is not considerable. In other words, it may be concluded that the small effect of shear deformation on the axial behaviour (~12%) can be considered as an inherent material noise in the experimental data. On the other hand, Cavallaro et al. (2007) reported that having the yarns under pretension in axial direction can greatly affect the subsequent shear behaviour of the fabrics, which is in fact the case from the simulation results in Figure 8.

After assessing the effect of combined loading on the basic normal and shear response of the fabric, another important notion may be studied. The question is, “Does the sequence of loading steps affect the response too?” In other words, if the axial loading is applied first, followed by the shear loading, or vice versa, are the resultant reaction force and moments the same as those when the two loadings are applied simultaneously?

To study the latter effect, let us define a normalized loading parameter $\alpha$ . It ranges from 0 to 1, where 0 refers to the initiation of loading and 1 represents the end of loading. For example, during a simultaneous/combined loading:

$$\theta \left(\alpha \right)=\alpha {\theta }_{max}; \epsilon \left(\alpha \right)=\alpha {\epsilon }_{max}$$ E9

where $\theta$ and $\epsilon$ , are the shear angle and axial strain in each step of loading and ${\theta }_{max}$ and ${\epsilon }_{max}$ , are the corresponding maximum values. Similarly, for the shear loading followed by the axial loading at $\alpha$ = $\frac{1}{2}$ we have:

$$\theta \left(\alpha \right)=\left(R\left(\alpha \right)-R\left(\alpha -\frac{1}{2}\right)\right) {\theta }_{max}; \epsilon \left(\alpha \right)=R\left(\alpha -\frac{1}{2}\right) {\epsilon }_{max}$$ E10

where,

$R\left(x\right)=\{\begin{array}{c}2 x x\ge 0\\ 0 x0\end{array}$

For the opposite case where the axial loading is followed by the shear loading, we may write:

$$\theta \left(\alpha \right)= R\left(\alpha -\frac{1}{2}\right){\theta }_{max}; \epsilon \left(\alpha \right)=\left(R\left(\alpha \right)-R\left(\alpha -\frac{1}{2}\right)\right) {\epsilon }_{max}$$ E11

Results of the new simulations are presented in Figure 9. It can be clearly seen that for the shear response, the sequence of the loading affects the resultant moment up to four times. However, the axial response is still less sensitive to the effect of deformation from the shear mode and the loading sequence. The results also indicate that if the shear deformation is applied to the specimen first, the shear reaction moment is decreased substantially. Moreover, during the step that the pure shear is applied, there seems to be a small reaction force in the form of tension. This is perhaps due to the fact that during shearing, the sliding of yarns on each other and their replacement in the fabric affect their waviness/crimp. In turn, the crimp interchange would induce a small axial stretch in some regions of yarns, especially if they are constrained at their ends (like in the picture frame test).

Figure 9.

The effect of loading sequence on the response of (a) shear and (b) axial deformation; the second loading step is applied after the dashed line for Axial+Shear and Shear+Axial cases.

The results in Figure 9 can also be linked to the constitutive model of yarns. Recalling Eq. (8), the transverse stiffness is a function of yarns’ axial and transverse strains (crushing formula). The bi-axial stretching induces axial strain in the yarns, which leads to an increase in the yarns’ transverse stiffness. In turn, the normal contact forces at the yarns cross over regions are increased, leading to a higher contribution from friction to the total reaction force. However, the opposite effect is not true. Under the bi-axial mode, the shear deformation (before locking point) does not induce considerable axial and transverse stretches in the yarns. Previously, using a sensitivity analysis under individual deformation modes, it was also reported by Komeili and Milani (2010) that the effect of transverse stiffness on the fabric response in the shear mode is considerable whereas it is ignorable in the bi-axial mode.

5. Summary

A numerical finite element model of a plain weave fabric unit cell at meso-level is developed. The model is capable of simulating specimens under simultaneous axial loading along the yarn directions and the fabric shearing. It can be a useful tool for predicting the meso-level local deformation phenomena in woven fabrics under complex loading conditions, as well as for developing equivalent material models at macro-level for fast simulation of fabric forming processes. Two fundamental deformation modes (shear and equi-biaxial stretching) are applied through two separate kinematics boundary conditions to facilitate extracting the contributions from each mode on the total resultant force and moment.

The analysis on the effect of combined loading has been conducted in two ways. First, the force and moment response of the unit cell under a predefined combined loading with a specific shear angle and axial strain is compared to those of the pure shear and axial modes. It was of interest to see if there is any interaction effect between the fundamental axial and shear deformation mechanisms when a combined loading is applied. Results showed that this interaction in fact exists and it has a dramatic effect on the ensuing reaction moment response (shear rigidity), but it is less important for the axial reaction force. Second, the effect of applying combined loading in two sequential steps was scrutinized. Again, the shear deformation response showed high sensitivity to the sequence of loading if it is applied before the axial deformation. Moreover, it was noted that during shear deformation there is a small tension reaction force, even though no stretching is applied to the yarns. This is perhaps due to the crimp interchanges along with the imposed boundary conditions on the end surfaces of yarns.

In summary, the above mentioned results show a high level of nonlinear interactions between the material response in the axial tension and shear modes. This can be directly related to the geometrical nonlinearities that exist in woven fabrics at meso-level and the effect of crimp interchanges during loading. After each stage of loading, the rearrangement of yarns in the fabric and their interactions should occur before yarns can go through further stretching/shearing. Under the combined loading, the crimp changes due to each loading mode can affect the reaction from the other mode. If loads are applied in sequence (e.g., shear followed by biaxial tension), the crimp changes in each step can affect the global response due to the effect from the previous loading step. Considerably different magnitudes of the shear moment were found between two cases where the shear and bi-axial deformations are applied at the same time and where the shear is applied after the axial loading. This observation clearly showed the higher sensitivity of the shear response to the crimp interchanges. On the contrary, because the axial reaction forces are more related to the stretching in the yarns, the shear deformation has minor influence on their axial force magnitudes. The effect of axial tension on increasing the transverse stiffness of yarns is deemed to be the main reason for the presence of interactions between the axial tension and shear deformation under combined loading modes. Further experimental and/or numerical studies are needed to scrutinize and validate the reported effects.

Acknowledgments

The authors would like to acknowledge financial support from the Natural Sciences and Engineering Research Council (NSERC) of Canada.