Multiscale statistical progressive failure analysis of unbalanced woven composites

Abstract

Woven composites are widely used across various applications, characterized by the warp and weft tows that intersect perpendicularly to form a weave pattern. Often, the weave architecture is neglected during analysis. This study explicitly models the weave pattern and the undulation of the tows using the NASA Multiscale Analysis Tool, incorporating statistical constituent material properties and measurements of unbalanced composite microstructures. Experimental data and multiscale modeling results were compared for carbon fiber composites manufactured with different tow sizes in the two in-plane directions, resulting in four distinct unbalanced reinforcement patterns. Multiple tensile tests were conducted for each pattern. Using statistical input data, 1000 random model instances were created for each unbalanced composite, predicting distributions of the composite’s Young’s modulus and tensile strength. The predictions demonstrated good agreement with test data, highlighting unique aspects of the modeling approach.

Keywords

Statistical micromechanics effect of undulation progressive failure analysis woven composites multiscale recursive micromechanics

Introduction

Yarn based textiles can be broadly classified into three categories, namely meshed yarn systems (e.g., knits), crossed yarn systems (e.g., weaves and braids) also known as crimp systems, and layered yarn systems also known as non-crimp systems. Crimp refers to the waviness or undulation of the yarn within the textile. Crimping may also refer to the filament waviness in a yarn, which usually provides its unique elasticity and bulkiness in a textile, for example in wool. In this study, undulation is considered at the yarn scale. A yarn is made of many filaments, e.g., a T700SC-24000 carbon fiber yarn from Toray has 24000 carbon filaments. A yarn is also called a tow.

Woven composites consist of a woven fabric (see Figure 1) which is embedded in a suitable resin system. A weave consists of two or more tow system that are perpendicular to each other, called the warp tow and the weft tow. The crossing of these tows creates a weave pattern or a weave texture.

Figure 1.

Woven fabric, with plain weave pattern.

Three common weave patterns are plain, twill and satin weaves.¹ Based on these, various woven structures can be designed. They can be either two-dimensional or three-dimensional fabrics. In a 3D-weave, the yarns not only interlace in the plane but also in the thickness direction, giving the material thickness and a true three-dimensional structure. Two-dimensional weaves are most commonly used, primarily due to their low cost of manufacturing and ease of application. Some examples of applications include wind turbine blades, aircraft wings, fuselages, and airbags. Typical characteristics of a weave are the orthogonal tow direction of the weft and the warp tows as shown in the Figure 1.^1,2

Analysis of woven composite has been the subject of research for many decades and continues to be an active area in the field of composite materials.^3–19 Various aspects have been investigated, including their mechanical behavior under different loading conditions, manufacturing process, and failure mechanisms. A variety of models and software have been developed for this purpose (e.g., TexGen, WiseTex, TexCAD). For simple weave geometries, these software tools provide streamlined solutions.

As indicated by Prodromou et al.,¹² analysis techniques available in the literature can be divided into four major classes: (1) analytical models, mostly based on lamination theories and on the assumption of isostrain or isostress.^{7,8,16,18,20–22} Byström et al. compares some of the analytical models and their applicability to predict the effective elastic properties.²³ Ayranci and Carey extended the analytical models from flat woven fabric to analyze 2D braided tubular composites, accounting for the component curvature²⁴; (2) inclusion models, which are based on mean field theories;²⁵ (3) voxel method or method of cells^26,27; (4) Finite Element Method (FEM) for more complex modeling and analysis using micromechanical failure analysis.^17,28–30

Increasing the modeling efficiency means minimization of computational cost (time and processor power), while maximizing model fidelity. Analytical techniques are more efficient, but falter when trying to capture the complex weave geometry and formulate realistic boundary conditions. Hence, most of the analytical models consider an idealized tow cross-section, weave pattern, fiber volume fraction, fiber orientation, and mean property distribution for the analysis. In reality, variation exists across all parameters.³¹

The limitations of the analytical techniques are typically addressed by employing more rigorous techniques such as the FEM, which allows for the consideration of the actual geometries and complex loading conditions. To capture the effects at different scales, the FE² approach^32–34 can be applied. This method bridges the gap between different component scales, i.e., the micro scale, concerning fiber-matrix interactions, and the macro scale, concerning component behavior, which is not feasible in the traditional continuum mechanics approach. However, this approach may result in significant computational costs. Furthermore, in the majority of the models, including those using micromechanics, variations in parameters such as the weave geometry and other material properties like fiber volume fraction, filament strength, diameter, modulus has not been accounted for. For progressive failure analysis, the distribution of these properties may have significant influence.

In this study, a distribution of the weave parameters and material properties have been considered as inputs, from which one thousand model instances have been analyzed to obtain predicted distributions of strength and stiffness of test coupons. The focus of this work is the analysis of unbalanced plain weave composite under tensile loading. The term unbalanced means that the two orthogonal fiber directions do not have the same number of filaments. This is explained in in the next section. For the modeling and analysis, the NASA Multiscale Analysis Tool (NASMAT), developed at the NASA Glenn Research Center, was employed.^35,36 The multiscale recursive micromechanics (MsRM) framework within NASMAT allows for various micromechanical theories to call each other or themselves recursively. This enables consideration of an arbitrary number of length scales, which is ideal for modeling textile composites. For example, Pineda et al.³⁷ used NASMAT to obtain the stress-strain response for a 3D-woven composite accounting for manufacturing defects.

NASMAT software enables use of a number of leading semi-analytical micromechanics theories, such as Mori-Tanaka, the Generalized Method of Cells (GMC), and the High-Fidelity Generalized Method of Cells (HFGMC),³⁸ at any length scale. Though all the three theories can be used, the HFGMC theory is used in this work. GMC’s first-order displacement field lacks coupling between shear and normal stresses. Bednarcyk and Arnold show that this has a significant impact on the predicted properties of woven composites, but it can be largely overcome by using a two-step homogenization approach.²⁶ GMC/HFGMC as implemented in NASMAT require orthogonal (parallelepiped) subvolumes, making modeling non-orthogonal and curved surfaces more challenging, although parametric HFGMC enables isoparametric subvolumes.³⁹ HFGMC employs a quadratic displacement field, and thus includes normal to shear coupling, resulting in accurate predictions for woven composites.^35,40–42This manuscript is organized as follows. First, the coupon manufacturing and experimental results are presented. Following this, the experimental measurements are used as input for the multiscale modeling using NASMAT. Last, experimental results are compared with the simulation predictions.

Experimental analysis

To analyze the effect of undulation, experiments were conducted at the coupon scale, focusing on tensile loading. The distribution of filament strength and stiffness parameters, filament diameter, tow geometry in a composite, and fiber volume fraction were determined. This section presents some of the critical aspects of the experimental techniques used in this work.

Coupon preparation

A primary consideration for this analysis is the difference in tow undulation in composites with different unbalanced plain weaves. Figure 2 represents two different such weave configurations. Due to the unbalanced weave pattern, where the warp tows have 24k filaments and the weft tows have 12k filaments, the undulation of the 24k tows is smaller compared to the case where both tows are 24k. The loading direction for tensile testing in this work is always the 24k direction. Further undulation angles are obtained by employing 3k tows and mono-filament configurations while maintaining 24k tows in the primary direction.

Figure 2.

Schematic representation of two weave configurations resulting in a variation of fiber undulation.

A radial braiding machine (RF-144) from Herzog GmbH, Germany, was utilized for the textile manufacturing process. In a braiding machine, two sets of bobbins move in opposite directions and follow sinusoidal paths, resulting in a braid. The manufacturing process is illustrated in Figure 3. Tubes with a diameter of 75 mm were braided using the selected tow combination. This combination is obtained by using fibers with different linear weights (tex), the details of which will be provided in the subsequent section. The number of tows in both the directions remain constant for all variations. The braiding angle was set to 45°, thereby achieving orthogonality between the tows. The braided layer was secured using appropriate techniques and subsequently cut along its length. The layer was opened, forming a 2D textile layer. Epoxy resin (Hexion RIMR426) was infused using vacuum bagging techniques. Following room temperature curing of the composite panels for a minimum of 10 hours, coupons were cut, and loading tabs were mounted using Araldite® 2-component glue with tempering at 100°C for 3 hours. The width of the coupons was approximately 25 mm.

Figure 3.

Coupon manufacturing process.

Panels with four different tow configurations were prepared. In all cases, one set of spools with 24k tows was unaltered, while the spool in the other bobbin was replaced with different tows. The braiding pattern was maintained to produce a plain weave. After curing, coupons were cut from the panel as shown in Figure 3. During tensile testing, the loading direction was consistently aligned with the direction of the 24k tows. Table 1 provides the different tow types and their properties. These properties are obtained from the technical data sheet of the manufacturer. Notice the nomenclature: 24k × 24k means both directions are 24k tows, while 24k × 3k means one direction is 24k and the perpendicular direction is 3k. The 24k × UD configuration employs a thin polyamide mono-filament in the perpendicular direction, hence the designation quasi-unidirectional (QUD). Note that failure of the Trofil PA12 monofilaments was not considered as they are quite compliant and also oriented normal to the applied tensile loading. In all cases, one ply was used for the manufacturing of the coupons. This approach minimizes the effects of fiber misalignment within plies and inter-ply interactions.

Table 1.

Fiber combinations to obtain unbalanced weave.

	Type	Tensile strength (MPa)	Tensile modulus (GPa)	Failure strain (%)	Density (g/cm³)
24k	Toray T700SC 1650 tex, 24,000 filaments	4900	230	2.1	1.8
12k	TohoTenax HTS40 800 tex, 12,000 filaments	4400	240	1.8	1.77
3k	TohoTenax HTS40 200 tex, 3000 filaments	4400	240	1.8	1.77
QUD	Trofil PA12 monofilament $\emptyset$ 0.16 mm	--	1.935	--	1.02
Epoxy	Hexion RIMR426	69	3.3	4	1.14 – 1.18

Measurement of physical parameters

To quantify the physical properties, such as tow area and post-cure fiber volume fraction, computer tomographic (CT) images were captured, as shown in Figure 4. Using a Python code with the computer vision (CV) library, the cross-sectional images were processed to determine the area fraction of the tow. When the code read each image, it represented the image as a list of row vectors, each vector containing pixel data. The length of each row vector corresponded to the number of pixels spanning the image width, and the length of the list of row vectors matched the number of pixels along the image height.

Figure 4.

Computer tomographic image of 24k × 24k and 24k x 12k coupon.

Each pixel was assigned a numerical value representing its brightness level in the greyscale color space. In CT imaging, materials of varying densities are depicted as pixels with different brightness levels. The regions containing both carbon fiber and matrix, the tows, exhibited higher brightness levels. To determine the volume fraction of the tows, the pixels corresponding to brightness levels of carbon fiber/epoxy regions were counted and divided by the total number of pixels in the measured domain. This can be expressed with the following formula:

V_{f} = \frac{N u m b e r o f c a r b o n f i b e r p i x e l s}{T o t a l n u m b e r o f p i x e l s i n t h e i m a g e} = \frac{n_{p i x e l s, f i b e r}}{n_{t o t a l}}

(1)

Furthermore, to capture the geometric variations inherent in the manufacturing of each braided specimen, the area of each tow in each image was computed. This was accomplished using the Python CV library, which identified regions of pixels based on their brightness levels and their connectivity to adjacent pixels of the same brightness. Once these regions were identified, the total area of each region was determined by counting the pixels within those regions.

To do this, the area of a single pixel needed to be determined by:

a_{p i x e l} = (P i x e l h e i g h t \times P i x e l w i d t h) = h_{p i x e l} w_{p i x e l}

(2)

h_{p i x e l} = \frac{I m a g e h e i g h t}{n u m b e r o f p i x e l s a l o n g t h e h e i g h t} = \frac{H_{i m a g e}}{n_{p i x e l s, 1}}

(3)

w_{p i x e l} = \frac{I m a g e w i d t h}{n u m b e r o f p i x e l s a l o n g t h e w i d t h} = \frac{W_{i m a g e}}{n_{p i x e l s, 2}}

(4)

Thus, to determine the area of a specific region in the image, simply take the sum the areas of each pixel in that region:

a_{r e g i o n} = f_{s}^{2} \sum_{i = 1}^{n} a_{p i x e l, i}

(5)

where,

a_{p i x e l, i}

is the area of the i^th pixel in the defined region and

f_{s}

is the image specific scaling factor, defined as:

f_{s} = \frac{l_{r e f, r e a l, i}}{L_{r e f, i m a g e}}

(6)

where

l_{r e f, r e a l, i}

is the reference length in image space (the representation of 1 mm in image space).

L_{r e f, i m a g e}

is the reference length, measured in physical space (the actual length of line representing 1 mm). Using this scaling factor, it was possible to link

a_{r e g i o n}

to the physically correct length scale.

Thus, each region can be quantified based on the number of pixels it contains. Since many sections were measured, it was possible to generate a mean tow area and determine the variation. This directly influenced the fiber volume fraction at the global level and thus the global response to mechanical loading. The results are consolidated in Table 2. Note that normal distributions were assumed for all measured physical properties. The undulation angles, θ₁ and θ₂, were treated as deterministic, hence only average values are given in Table 2.

Table 2.

Physical properties of the tow measured using image analysis of the computer tomographic images.

	Property	Unit	Mean	Standard deviation
24k × 24k	Volume fraction (V_f) of 24k tows in composite	-	0.312	0.047
	Volume fraction (V_f) of transverse tows in composite	-	0.352	0.052
	24k tow cross-section area	mm²	1.46	0.16
	Transverse tow cross-section area	mm²	1.42	0.27
	24k tow undulation angle (θ₁)	°	11.2	-
	Transverse tow undulation angle (θ₂)	°	11.37	-
	Average V_f from NASMAT of all 1000 instances^a	%	41.5	-
24k × 12k	Volume fraction (V_f) of 24k tows in composite	-	0.488	0.064
	Volume fraction (V_f) of transverse tows in composite	-	0.224	0.043
	24k tow cross-section area	mm²	1.40	0.29
	Transverse tow cross-section area	mm²	0.831	0.063
	24k tow undulation angle (θ₁)	°	8.35	-
	Transverse tow undulation angle (θ₂)	°	9.7	-
	Average V_f from NASMAT of all 1000 instances^a	%	43.4	-
24k × 3k	Volume fraction (V_f) of 24k tows in composite	-	0.621	0.119
	Volume fraction (V_f) of transverse tows in composite	-	0.084	0.014
	24k tow cross-section area	mm²	1.54	0.30
	Transverse tow cross-section area	mm²	0.251	0.056
	24k tow undulation angle (θ₁)	°	3.03	-
	Transverse tow undulation angle (θ₂)	°	5.8	-
	Average V_f from NASMAT of all 1000 instances^a	%	40.0	-
QUD	Volume fraction (V_f) of 24k tows in composite	-	0.632	0.085
	Volume fraction (V_f) of transverse filaments in composite	-	0.009	0.002
	24k tow cross-section area	mm²	1.61	0.21
	Transverse tow cross-section area	mm²	0.020	0
	24k tow undulation angle (θ₁)	°	4.75	-
	Transverse tow undulation angle (θ₂)	°	7.2	-
	Average V_f from NASMAT of all 1000 instances^a	%	35.4	-

^aDetermined purely by geometry measurements plus the known number of filaments in a given tow and the tow diameters. Slight variation (∼0.2%) per model instance was present due to the hexagonal fiber packing discretization used at the micro scale.

Experiments

Experimental analysis was divided into two parts. Firstly, to account for the mechanical property variation at the filament scale, single filament testing on both types of carbon fiber was performed. Following this, the tensile testing results from the unbalance woven composite coupons are presented.

Filament testing was conducted on the FAVIMAT + fiber testing machine from the company Textechno GmbH, Germany. A total of fifty samples were tested for each carbon fiber type (see Table 1). The properties, along with their standard deviations, used for numerical modeling, are detailed in Table 3.

Table 3.

Statistical material properties used in the modeling.

	Property	Unit	Mean^a	Standard deviation
T700SC	Young’s modulus in fiber direction	GPa	220.3	10.6
	Strength	GPa	4.64	1.51
	Diameter	µm	6.87	0.24
HTS40	Young’s modulus in fiber direction	GPa	221.4	4.6
	Strength	GPa	4.26	0.64
	Diameter	µm	6.97	0.21
Epoxy RIMR426^b	Young’s modulus	GPa	3.3	0.096
Epoxy RIMR426^b	Strength	MPa	69	1.9

^a50 samples were tested for each fiber type.

^bEpoxy mean data were taken from the manufacturer technical data sheet and the standard deviations were approximated based on the statistics in Hobbiebrunken et al.⁴³

Subsequently, the manufactured coupons were tensile tested using a Zwick Z100 testing machine, based on the DIN EN ISO 527-4 standard. A displacement rate of 2 mm/min was used for testing. Displacement measurements were obtained via an optical extensometer. The resulting experimental stress-strain curves are illustrated in Figure 5 and Table 4. All measurements are plotted as grey line and the mean for each configuration is plotted in black. The color contour shows the bounds for each weave type. The results exhibit considerable variations in strength and stiffness across different weave types.

Figure 5.

Tensile testing results for all the configurations.

Table 4.

Tensile testing results.

Name	Mean of max. Strength [MPa]	Standard deviation of max. Strength [MPa]	Mean Young’s Modulus [MPa]	Standard deviation of Young’s Modulus [MPa]
24k × 24k	436.34	35.96	26443.15	3805.11
24k × 12k	848.74	107.89	45295.76	2044.94
24k × 3k	1349.18	127.05	66899.81	8085.21
QUD	1252.92	187.82	59780.29	6194.59

Interestingly, the 3K configuration demonstrates greater stiffness than the UD configuration, which might initially appear counterintuitive. This phenomenon could be explained by the larger thickness of the monofilament fiber compared to the 3K tow, resulting in a larger undulation angle (see Table 2).

It can be observed from Figure 6 that 24k × 24k and QUD have a notable difference in their final failure mechanisms. Due to the interlacing of the fibers in the 24k × 24k, the crack propagation is seen in these fibers. The rupture is seen to occur at the point of undulation of the fibers in the loading direction. One the other hand, in the QUD case, due to the absence of loadbearing interlacing fibers, individual fibers fail at various positions. Debonding of individual fibers is evident. This phenomenon is akin to what one typically observes in filament wound components.

Figure 6.

Example failure modes in the coupons.

Multiscale micromechanics

The Multiscale Recursive Micromechanics (MsRM) approach, shown schematically in Figure 7, is designed to model materials with multiple heterogeneous microstructures across different length scales.^41,42 As shown, MsRM utilizes repeating unit cell (RUC) micromechanics models embedded within each other to capture the microstructural geometry at any number levels of scale. The effective stiffness matrix for the RUC is passed up the length scales from Level k to Level 0 during homogenization. The subvolume strains are passed down the length scales from Level 0 to Level k and applied to the subscale RUC. Local and global stresses are calculated using the strains and stiffness. For additional details, the reader is referred to Ref. 35. Multiscale modeling of other physics, such as thermal/electrical conductivity, diffusion, and magnetic permeability⁴⁴ can also be conducted with MsRM.⁴⁵

Figure 7.

Schematic detailing the multiscale recursive micromechanics (MsRM) approach whereby separate micromechanics models can be embedded within each other to capture microstructural features of interest across any number of length scales. Each level represents a length scale, with level 0 being the highest, global length scale, and level k being the lowest length scale.

MsRM has been implemented in the comprehensive software suite, NASA Multiscale Analysis Tool (NASMAT).⁴⁶ This software enables the use of a number of leading semi-analytical micromechanics theories, such as Mori-Tanaka, the Generalized Method of Cells (GMC), and the High-Fidelity Generalized Method of Cells (HFGMC),³⁸ at any length scale. In addition, multiple nonlinear damage and viscoplastic constitutive models are available, and NASMAT links with commercially available finite element codes to enable micromechanics analysis at the integration points in a structural finite element model.⁴⁷

Woven and braided composites can be readily analyzed with MsRM within NASMAT as these materials typically exhibit heterogeneous microstructures at a minimum of two length scales. At the higher length scale, the weave pattern geometry can be considered, while, at the lower length scale, the filaments and matrix, which constitute to the woven/braided fiber tows, can be modeled. This approach, as shown in Figure 8, was taken herein to model the unbalanced woven composites. At the highest length scale, the simplest RUC that still captures the key aspects of the weave pattern was employed. It consisted of 7 × 4 × 4 subvolumes representing the tows and the pure matrix regions between the tows. The HFGMC micromechanics theory³⁸ was used at this scale. The sizes and volume fractions of the tows in each direction are captured, as are the undulation angles (θ₁ and θ₂). In addition to these undulation angles, the dimension a₁, a₂, b₁, b₂, and t were sufficient to fully characterize the higher scale woven composite RUC microstructure. At the lower length scale, representing the tows in the woven composite, a hexagonally packed microstructure with 32 × 56 subvolumes, was employed, and, once again, the HFGMC micromechanics theory was used at this scale. Note that, by changing the number of subvolumes occupied by the filaments (shown in blue in Figure 8), the fiber volume fraction within the tows can be altered as desired. The physical property measurements shown in Table 2, along with the known tow filament counts, were sufficient to fully characterize the needed woven composite RUC geometric parameters, including the V_f within each tow. Note that the fiber volume fraction in each of the two orthogonal weave directions is a primary driver of the material response, and the constructed RUC, while coarse, correctly captures this feature. Furthermore, these fiber volume fractions are maintained across all statistical cases that will be considered.

Figure 8.

Example repeating unit cells (RUCs) for the unbalanced woven composite weave pattern (higher length scale) and within the tows (lower length scale).

Predictions of the unbalanced woven composite stress-strain response were conducted using NASMAT by applying incremental strain-controlled uniaxial loading in the direction of the 24k tows to the higher scale composite RUC (while maintaining zero for all the non-loaded stress components). The effective elastic properties result from first homogenizing the tows, and then homogenizing the composite RUC. Note that appropriate coordinate transformations are necessary for the tows in the undulation subvolumes. Conversely, based on the applied strain level, localization provides the strains in the composite subvolumes, which are then localized again within the tows to obtain micro scale strains in the tow subvolumes. To predict progressive damage and failure, the simple subvolume elimination method⁴² was employed. In this approach, the constituents are assigned strength or allowable values (stresses or strains), and then, once a selected failure criterion is met or exceeded in a given subvolume, that subvolume’s stiffness components are immediately reduced to nearly zero. Herein a reduction by a factor of 10⁴ was used to effectively eliminate the failed subvolume’s contribution to the composite stiffness, along with the maximum stress failure criterion. As such, even though the subvolume elimination approach simulates abrupt, brittle failure, because this occurs at the level of the subvolumes within the tows (and between the tows), the overall failure of the composite can be quite progressive. That is, for example, when a single matrix subvolume fails within a tow, this will only slightly affect the homogenized properties of the tow, which then will have only a slight effect on the homogenized properties of the woven composite RUC. As the incrementally applied global composite strain increases, more and more fiber and matrix subvolume fail, eventually resulting in significant loss of global composite stiffness, at which point the composite material may be considered failed. The maximum global stress that is reached in the simulation can then be treated as the predicted composite strength. The maximum stress criterion is thus applied only at the constituent scale, and through homogenization, the effect of the damage is manifested at the higher length scales. No failure criteria are needed for any length scale higher than the constituents.

Results and discussion

Predictions based on mean input values

Using the NASMAT model presented in the previous section, tensile predictions were first conducted wherein the mean values were employed for all input parameters. Figure 9 shows a comparison of these predicted tensile test results with the experimental results for the four weave types. The large drops in the predicted stress strain curves indicate the points at which local failures have accumulated such that, globally, the composite can no longer support any additional load. This has been considered composite failure, and the maximum stress attained in the model prior to this point has been taken as the predicted composite strength.

Figure 9.

NASMAT progressive failure predictions, compared with tensile test data, using mean values of all input data, for the four weave types.

Clearly, the predictions shown in Figure 9 demonstrate reasonably good agreement with the experimental tensile curves. Notably, the lower stiffness and strength of the QUD case compared to the 24k × 3k case is correctly captured by the model. It is also noteworthy that the stiffening exhibited in the test data has not been captured as no such mechanism is present in the current model. However, this may be considered in future work as NASMAT contains large fiber rotation capabilities, which have been applied in the past to dry woven fabrics.⁴⁸ Each simulation required approximately 1 minute to execute on an i7 processor laptop.

Effect of input property variation

To estimate the impact of the constituent material property and woven composite geometry statistics, a script was written to randomly sample the random input variables given in Table 2 and Table 3. NASMAT is quite amenable to scripting because it can be driven by an ASCII input file. For each random sample of the input variables, the script generated a single instance of a NASMAT input file. For each woven composite type, 1000 such instances were generated, along with a batch file containing the 1000 NASMAT execution commands. This batch file was then used to run the 1000 NASMAT cases serially, with each run taking approximately 1 minute (or 16.7 hours for all 1000 cases), again on an i7 processor laptop. Note that execution of these cases could be quickened by running individual cases simultaneously on multiple processors, or by running each batch of 1000 cases (per weave type) on its own processor simultaneously.

Figure 10 shows the 1000 instances of the NASMAT tensile test simulations compared to the test data for each weave type. The model prediction line colors indicate the mode of the final failure. Red lines indicate simulations which failed, due fiber breakage, whereas blue lines indicate that final failure occurred due to complete matrix section failure without failure of the fibers. First, Figure 10 indicates that the multiscale model instances have done a reasonably good job of reproducing the experimental tensile curves. The initial stiffness predictions match better for the QUD and 24k × 3k weaves, whose experimental responses exhibited less stiffening. For the 24k ×12k and 24k × 24k weaves, which exhibited more stiffening experimentally, the model is overpredicting the initial stiffness, but matches the instantaneous stiffness better as the strain is increased and the experimental stiffening occurs. In terms of strength, the ranges of predicted strengths encompass the test data, as they should. It is interesting that out of the 1000 instances per weave type, there are some very low predicted strengths, which are always red curves, indicating final failure due to fiber failure. In fact, it is clear that the simulations that exhibited final failure due to the expected fiber failure mode were consistently the lower strength predictions. This is because these instances were assigned (based on the statistics, see Table 3) low 24k (T700SC) fiber strengths, causing predicted final failure prior to the development of the alternative mode, complete matrix section failure (blue lines). Note that the coefficient of variation of the T700SC fiber was measured to be 30% (see Table 3). Of course, because many geometric parameters and material properties were random, some intermingling between these predicted final failure modes was observed. Finally, recall that the geometric variations have been applied to reproduce the correct known volume fraction of fibers in each direction for each weave. Therefore, none of the variation in the predictions is due to changes in the overall composite fiber volume fraction; it is due only to microstructural and property variation.

Figure 10.

1000 random NASMAT progressive failure prediction instances, compared with tensile test data (solid black lines), for the four weave types. (a) QUD, (b) 24k × 3k, (c) 24k × 12k, (d) 24k × 24k. The model prediction line colors indicate final failure via fiber breakage (red) and complete matrix section failure (blue).

Investigation into the unexpected complete matrix section failure mode indicated that it occurred when the matrix in the undulation subvolumes (see Figure 8) fails completely, which then results in a near zero effective stiffness for the composite in the applied loading direction. In these cases, the model is thus predicting failure of the composite without failure of the fibers. Obviously, this is not realistic, and it is actually an artifact of the scale separation in the multiscale modeling approach coupled with the aforementioned current inability of the model to consider large fiber rotations in woven composites. It can be explained by considering a subvolume from the tow undulation region (see Figure 8), as depicted in Figure 11. The transformation between the local (x₁ - x₂) and global (X - Y) coordinates results in a global Young’s modulus in the Y direction⁴⁹ of,

\frac{1}{E_{Y Y}} = \frac{1}{E_{11}} \cos^{4} θ_{1} + (\frac{1}{G_{12}} - \frac{2 ν_{12}}{E_{11}}) \sin^{2} θ_{1} \cos^{2} θ_{1} + \frac{1}{E_{22}} \sin^{4} θ_{1}

where E₁₁, E₂₂, G₁₂, and ν₁₂ are the local Young’s moduli, shear modulus, and Poisson’s ratio, respectively, and θ₁ is the undulation angle. This equation shows that, if damage at local scale is such that either E₂₂ or G₁₂ approaches zero, then the global Young’s modulus in the Y direction must also approach zero (unless θ₁ also approaches zero). This is true regardless of any reduction in the axial stiffness, E₁₁. Therefore, when the matrix material at the local scale fails completely and is assigned a near-zero stiffness, this undulation subvolume will have near-zero E_YY, even if the there is no damage in the fiber direction. This then precludes the fibers from being loaded any further, even as the global loading continues. In reality, of course, if the matrix in the undulation zone were to fail completely, the fibers would be free to rotate and align with the loading direction, enabling the fibers to be further loaded. As mentioned previously, this type of fiber rotation is available in NASMAT for dry fabrics,⁴⁸ but this capability has not yet been extended for woven composites like those considered here. There is thus a caveat regarding the predictions shown in blue in Figure 10: it is assumed that once the matrix in the undulation regions fails completely, the fibers will quickly rotate into alignment with loading direction and fail. The validity of this assumption will be investigated in future work. Note that there is some evidence of fiber pullout for the 24k × 24k specimen in Figure 6, which would indicate fiber straightening during pullout prior to fiber failure.

Figure 11.

A single subvolume from the tow undulation region of the woven composite RUC, with the local and global coordinates, along with the undulation angle, indicated.

Figure 12 presents the model and experimental cumulative distribution function (CDF) plots for the strength distributions for all the four weave types. The CDF is a fundamental concept in probability and statistics, representing the probability that a random variable $X$ (tensile strength, in the present case) takes on a value less than or equal to a specific value $x$ (which serves as the abscissa in a CDF plot). Mathematically, the CDF is defined as $F_{X} (x) = P (X \leq x)$ , where $F_{X} (x)$ denotes the CDF of the random variable $X$ . This function provides a comprehensive description of the distribution of $X$ , encompassing all probabilities across the range of possible values. The CDF starts at 0 for the minimum value of the random variable and increases to 1 at the maximum value.

Figure 12.

Comparison of the test and simulation cumulative distribution functions (CDFs) for the four weave types.

Figure 12 indicates that, statistically, the model predictions for the unbalanced woven composite strength match quite well with the experimental data. While the QUD and 24k × 3k predicted CDFs match better with experiments, there is reasonable overlap for the 24k × 12k and 24k × 24k as well. Of particular note in the predicted CDF of the 24k × 24k composite is the jog that occurs just above a 50% probability, while the predicted CDFs for the other composites remain smooth. Examining the probability distribution function (PDF) for this 24k × 24k case (Figure 13), the reason for the jog is clear: the model has predicted a bimodal strength distribution (which was not predicted for the other weave types). Further, this bimodal strength distribution is clearly not due to the two types of failure modes discussed early as, from Figure 10(d), very few of the 24k × 24k simulations exhibited the fiber failure mode, whereas each mode in Figures 12 and 13 encompass close to half of the simulations.

Figure 13.

Simulated probability distribution function (PDF) for the 24k × 24k woven composite.

To examine the source of the bimodal strength distribution predicted for the 24k x 24k composites, a simulation from each of the populations has been plotted in Figure 14. The embedded contour plots show a top view of a damage magnitude, which considers the change in all of the diagonal components of the compliance matrix (squared, summed, square root taken, and normalized) such that 0 indicates no damage and 1 indicates complete damage. Note that these damage contours represent homogenized values through the thickness, not simply the damage on the top of the RUC.

Figure 14.

Comparison of two 24k × 24k woven composite tensile simulations, one of which resulted in a higher strength in the PDF (see Figure 13) and one of which resulted in a lower strength in the PDF. The embedded contour plots show the damage magnitude at the indicated points along the stress-strain curves.

In Figure 14, the red (lower strength) stress-strain curve exhibits some minor damage in the undulation tows oriented transverse to the loading direction prior to final failure at approximately 300 MPa. At this point, the undulation tows oriented in the loading direction fail completely, as do the previously damaged undulation tows oriented transverse to the loading direction. The blue (higher strength) stress-strain curve also exhibits some minor damage in the undulation tows oriented transverse to the loading, but the damage of the undulation tows oriented in the loading direction is delayed to approximately 350 MPa. More importantly, this damage is less severe compared to the lower strength case and it does not couple with additional damage in the previously damaged undulation tows oriented transverse to the loading direction. As such, while there is a slight drop in the stress, the composite survives to a stress of approximately 500 MPa, at which point additional damage has accumulated in the undulation tows oriented transverse to the loading direction. This then leads to complete failure of the undulation tows oriented in the loading direction and the large stress drop indicative of global failure. Thus, the main discriminator between the two failure mode populations shown in Figure 13 for the 24k × 24k composites is the ability to survive damage initiation the undulation tows oriented in the loading direction (which appears to be dependent on many of the random input variables). Interestingly, the final damage pattern is quite similar between the two cases shown in Figure 14, even though the path taken is very different. Finally, note that the damage contours show that no damage is predicted in the four pure matrix regions between the tows. This is because the pure matrix is very low stiffness, and thus these regions experience little stress and contribute little to the overall composite stiffness (although they will eventually fail).

Obviously, the above discussion is applicable only to the 24k × 24k simulation results. There was no indication in the test results of a bimodal strength distribution, nor the discriminator between the populations demonstrated in the predictions. However, only two 24k × 24k tests were performed.

As a final comparison between the model and the test results, Figure 15 shows the mean Young’s modulus and tensile strength. For the simulated results, both the deterministic prediction using mean input properties (from Figure 9) and the mean of the 1000 random instances (from Figure 10) are shown. As discussed, due to the lack of a stiffening mechanism in the model, the Young’s modulus is over predicted. In addition, the predictions of the model using mean input properties match closely with the mean of the 1000 random instances. In contrast, for the tensile strength, there is quite some difference between the predictions of the model using mean input properties and the mean of the 1000 random instances. Further, the means of the 1000 random instances provides better agreement with the mean of the test data in all cases aside from the 24k × 24k.

Figure 15.

Comparison for the woven composite predicted and experimental (a) mean Young’s modulus and (b) mean tensile strength. For the simulated results, both the deterministic prediction using mean input properties and the mean of the 1000 random instances are shown.

Conclusion

Tensile testing and statistical multiscale model have been executed on unbalanced plain weave composites, created using tows of different sizes in the two in-plane directions. The Multiscale Recursive Micromechanics model in NASMAT software was employed to capture progressive damage at the constituent material level within and between the tows. Microstructural measurements and carbon filament testing provided the necessary statistical constitutive data for the multiscale models. Four types of unbalanced weaves were examined, with 1000 random simulations conducted for each to predict the tensile response and failure strength. Model predictions aligned well with limited experimental results, including the cumulative distribution functions of composite strength. Experimental tensile curves exhibited stiffening that the model did not account for, leading to an overestimation of initial elastic Young’s moduli, but better alignment at intermediate strains. Many simulations indicated failure due to complete matrix section failure, rather than fiber failure, attributed to the model’s scale separation and fiber orientations in the tow undulation regions, resulting in near-zero load-bearing capability once the matrix failed. In reality, significant fiber rotation would occurs, aligning fibers with the load and allowing fiber failure—a capability planned for future model updates. For the 24k × 24k weave, a bimodal strength distribution was predicted. Due to the limited number of experiments, it remains uncertain if this reflects reality, but model results explained it through the ability to withstand damage initiation in the loading-oriented undulation tows. The lower strength group failed shortly after this damage initiation, while the higher strength group did not experience immediate damage progression.

Integrating micromechanical models into the design process makes it possible to account for variable properties at constituent level. In most of the design exercises, though a distribution of material data is available, often a mean value is considered. Incorporating this distribution data in to the micromechanical model developed in this study offers a more comprehensive understanding of the failure mechanisms thereby enhancing the accuracy and reliability of the analysis.

Footnotes

Acknowledgments

The first author acknowledges support from the German Federal Ministry for Economic Affairs and Climate Action for funding the project TowpregBraid (KK5055931KU1) as part of the ZIM innovation program. The authors thank the German federal and state government (BMBF & MWK) for funding the Theodore von Kármán Fellowship (GS132) under the Excellence Strategy. B.A. Bednarcyk acknowledges the NASA Space Technologies Mission Directorate Cryotank Technology for Exploration Applications (CTE-A) project.

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: this work was supported by the German Federal Ministry for Economic Affairs and Climate Action; KK5055931KU1 and Theodore von Karman Fellowship; GS132.

ORCID iD

Kumar C Jois

Data availability statement

All relevant data is available in this manuscript.*

References

Gries

Veit

Wulfhorst

. Textile Fertigungsverfahren: Eine einführung. 3rd edition. München, Germany: Hanser, 2019.

Schürmann

. Konstruieren Mit Faser-Kunststoff-Verbunden. 2nd edition. Berlin, Germany: Springer Berlin, Heidelberg, 2008.

Alif

. Mechanical characterization of woven fabric composite materials. Boca Raton, FL: Florida Atlantic University, 1996.

Aliabadi

Guagliano

(eds) Fracture and damage of composites. 1st edition. Billerica, MA: WIT, 2006. Southampton UK.

Dang

Liu

Zhang

, et al. Theoretical prediction for effective properties and progressive failure of textile composites: a generalized multi-scale approach. Acta Mech Sin 2021; 37(8): 1222–12444. DOI: 10.1007/s10409-021-01098-8.

Dixit

Mali

. Modeling techniques for predicting the mechanical properties of woven-fabric textile composites: a review. Mech Compos Mater 2013; 49(1): 1–20. DOI: 10.1007/s11029-013-9316-8.

Ishikawa

Chou

T-W

. Stiffness and strength behaviour of woven fabric composites. J Mater Sci 1982; 17(11): 3211–3220. DOI: 10.1007/BF01203485.

Ishikawa

Chou

T-W

. One-dimensional micromechanical analysis of woven fabric composites. AIAA J 1983; 21(12): 1714–1721. DOI: 10.2514/3.8314.

Kim

Song

J-I

. Geometry effect on mechanical properties of woven fabric composites. J Cent South Univ Technol 2011; 18(6): 1985–1993. DOI: 10.1007/s11771-011-0932-y.

10.

Jiang

Guo

Yue

. Investigation of the three-dimensional micromechanical behavior of woven-fabric composites. Mech Compos Mater 2006; 42(2): 141–150. DOI: 10.1007/s11029-006-0025-4.

11.

Zhou

, et al. Formulation of a unit cell of a reduced size for plain weave textile composites. Comput Mater Sci 2011; 50(5): 1770–1780. DOI: 10.1016/j.commatsci.2011.01.013.

12.

Prodromou

Lomov

Verpoest

. The method of cells and the mechanical properties of textile composites. Compos Struct 2011; 93(4): 1290–1299. DOI: 10.1016/j.compstruct.2010.09.022.

13.

Naik

Ganesh

. Prediction of on-axes elastic properties of plain weave fabric composites. Compos Sci Technol 1992; 45(2): 135–152. DOI: 10.1016/0266-3538(92)90036-3.

14.

Lomov

Huysmans

Luo

, et al. Textile composites: modelling strategies. Compos Appl Sci Manuf 2001; 32(10): 1379–1394. DOI: 10.1016/s1359-835x(01)00038-0.

15.

Tan

Tong

Steven

. Modelling for predicting the mechanical properties of textile composites—a review. Compos Appl Sci Manuf 1997; 28(11): 903–922. DOI: 10.1016/s1359-835x(97)00069-9.

16.

Whitney

Chou

T-W

. Modeling of 3-D angle-interlock textile structural composites. J Compos Mater 1989; 23(9): 890–911. DOI: 10.1177/002199838902300902.

17.

Huang

Zhao

, et al. Progressive failure prediction of woven fabric composites using a multi-scale approach. Int J Damage Mech 2018; 27(1): 97–119. DOI: 10.1177/1056789516663613.

18.

Zhou

Gao

Liu

, et al. An analytical model for the uniaxial tensile modulus of plain-woven fabric composites and its application and experimental validation. Aerospace 2022; 9(1): 26. DOI: 10.3390/aerospace9010026.

19.

Raju

Wang

. Classical laminate theory models for woven fabric composites. J Compos Technol Res 1994; 16(4): 289–303. DOI: 10.1520/CTR10589J.

20.

Wang

Sun

Yang

, et al. An analytical approach for predicting the tensile strength of plain weave fabric composites. Polym Compos 2019; 40(6): 2391–2399. DOI: 10.1002/pc.25099.

21.

Morozov

. The effect of filament-winding mosaic patterns on the strength of thin-walled composite shells. Compos Struct 2006; 76(1-2): 123–129. DOI: 10.1016/j.compstruct.2006.06.018.

22.

Naik

. Analysis of woven and braided fabric-reinforced composites. In: D30 (ed) ASTM symposium on composite materials: testing and design. West Conshohocken, PA: ASTM International, 1994, pp. 239–263.

23.

Byström

Jekabsons

Varna

. An evaluation of different models for prediction of elastic properties of woven composites. Compos B Eng 2000; 31(1): 7–20. DOI: 10.1016/s1359-8368(99)00061-x.

24.

Ayranci

Carey

. Predicting the longitudinal elastic modulus of braided tubular composites using a curved unit-cell geometry. Compos B Eng 2010; 41(3): 229–235. DOI: 10.1016/j.compositesb.2009.10.006.

25.

Huysmans

Verpoest

van Houtte

. A poly-inclusion approach for the elastic modelling of knitted fabric composites. Acta Mater 1998; 46(9): 3003–3013. DOI: 10.1016/S1359-6454(98)00021-4.

26.

Bednarcyk

Arnold

. Micromechanics-based modeling of woven polymer matrix composites. AIAA J 2003; 41(9): 1788–1796. DOI: 10.2514/2.7297.

27.

Hang

Cui

Liu

, et al. Micro/meso-scale damage analysis of a 2.5D woven composite including fiber undulation and in-situ effect. Compos Struct 2021; 256(12): 113067. DOI: 10.1016/j.compstruct.2020.113067.

28.

Hayat

Lei

Ali

. Prediction of elastic behavior of woven fabric reinforced plastics composites using two-step homogenization. In: 2018 15th international bhurban conference on applied sciences and technology (IBCAST), Islamabad, Pakistan, 9–13 January 2018, IEEE, 2018.

29.

Liu

, et al. A novel modeling method for the mechanical behavior of 3D woven fabrics considering yarn distortion. Compos Sci Technol 2022; 230(8): 109691. DOI: 10.1016/j.compscitech.2022.109691.

30.

Zhang

Harding

. A numerical micromechanics analysis of the mechanical properties of a plain weave composite. Comput Struct 1990; 36(5): 839–844. DOI: 10.1016/0045-7949(90)90154-T.

31.

Talreja

. Failure analysis of composite materials with manufacturing defects. 1st edition. Boca Raton, FL: CRC Press, 2023.

32.

Zhang

, et al. A hierarchical multiscale model for the elastic-plastic damage behavior of 3D braided composites at high temperature. Compos Sci Technol 2020; 196(12): 108230. DOI: 10.1016/j.compscitech.2020.108230.

33.

Carvelli

Poggi

. A homogenization procedure for the numerical analysis of woven fabric composites. Compos Appl Sci Manuf 2001; 32(10): 1425–1432. DOI: 10.1016/S1359-835X(01)00041-0.

34.

Tikarrouchine

Benaarbia

Chatzigeorgiou

, et al. Non-linear FE2 multiscale simulation of damage, micro and macroscopic strains in polyamide 66-woven composite structures: analysis and experimental validation. Compos Struct 2021; 255(8): 112926. DOI: 10.1016/j.compstruct.2020.112926.

35.

Pineda

Bednarcyk

Ricks

, et al. Efficient multiscale recursive micromechanics of composites for engineering applications. JMC 2021; 19(4): 77–105. DOI: 10.1615/IntJMultCompEng.2021039732.

36.

Aboudi

Arnold

Bednarcyk

. Micromechanics of composite materials: a generalized multiscale analysis approach. Oxford, UK: Butterworth-Heinemann, 2013.

37.

Pineda

Bednarcyk

Ricks

, et al. Multiscale failure analysis of a 3D woven composite containing manufacturing induced voids and disbonds. Compos Appl Sci Manuf 2022; 156: 106844. DOI: 10.1016/j.compositesa.2022.106844.

38.

Aboudi

Arnold

Bednarcyk

. Practical micromechanics of composite materials. Oxford, UK: Butterworth-Heinemann, 2021.

39.

Haj-Ali

Aboudi

. Formulation of the high-fidelity generalized method of cells with arbitrary cell geometry for refined micromechanics and damage in composites. Int J Solid Struct 2010; 47(25): 3447–3461. DOI: 10.1016/j.ijsolstr.2010.08.022.

40.

Bednarcyk

Ricks

Pineda

, et al. Thermal conductivity of 3D woven composite thermal protection system materials via multiscale recursive micromechanics. AIAA SCITECH 2022 forum. 2022. DOI: 10.2514/6.2022-2280.

41.

Ricks

Pineda

Bednarcyk

, et al. Multiscale progressive failure analysis of 3D woven composites. Polymers 2022; 14(20): 4340. DOI: 10.3390/polym14204340.

42.

Pineda

Bednarcyk

Ricks

, et al. Multiscale failure analysis of a 3D woven composite containing manufacturing induced voids and disbonds. Compos Appl Sci Manuf 2022; 156(12): 106844. DOI: 10.1016/j.compositesa.2022.106844.

43.

Hobbiebrunken

Fiedler

Hojo

, et al. Experimental determination of the true epoxy resin strength using micro-scaled specimens. Compos Appl Sci Manuf 2007; 38(3): 814–818. DOI: 10.1016/j.compositesa.2006.08.006.

44.

Bednarcyk

Aboudi

Arnold

. Micromechanics of composite materials governed by vector constitutive laws. Int J Solid Struct 2017; 110-111: 137–151. DOI: 10.1016/j.ijsolstr.2017.01.033.

45.

Bednarcyk

Ricks

Pineda

, et al. Multiscale recursive micromechanics of three-dimensional woven composite thermal protection materials thermal conductivities. AIAA J 2022; 60(12): 6506–6519. DOI: 10.2514/1.J061920.

46.

NASA . NASA multiscale analysis tool (NASMAT): LEW-20244-1 2024. https://software.nasa.gov/software/LEW-20244-1.

47.

Pineda

Ricks

Bednarcyk

, et al. Benchmarking and performance of the NASA multiscale analysis tool. In: AIAA Scitech 2021 forum. Reston, VA: American Institute pf Aeronautics and Astronautics, 2021.

48.

Hearley

Pineda

Bednarcyk

, et al. Predicting fiber breakage failure mode of plain weave fabrics with multiscale recursive micromechanics. In: 37 ASC annual technical conference, Tucson, AZ, 19–21 September 2022, American Institute of Aeronautics and Astronautics, 2022.

49.

Jones

. Mechanics of composite materials. Boca Raton, FL: CRC Press, 2018.