NASA/CR-2005-213905 

NIA Report No. 2005-01 

 

August 2005 

Modeling of Unit-Cells with Z-Pins Using 

FLASH: Pre-Processing and Post-Processing 

Ronald Krueger 

National Institute of Aerospace, Hampton, Virginia 
The NASA STI Program Office . . . in Profile

Since its founding, NASA has been dedicated to the 

advancement of aeronautics and space science. The 

NASA Scientific and Technical Information (STI) 

Program Office plays a key part in helping NASA 

maintain this important role. 

The NASA STI Program Office is operated by 

Langley Research Center, the lead center for NASA’s 

scientific and technical information. The NASA STI 

Program Office provides access to the NASA STI 

Database, the largest collection of aeronautical and 

space science STI in the world. The Program Office is 

also NASA’s institutional mechanism for 

disseminating the results of its research and 

development activities. These results are published by 

NASA in the NASA STI Report Series, which 

includes the following report types: 

• TECHNICAL PUBLICATION. Reports of 

completed research or a major significant phase 

of research that present the results of NASA 

programs and include extensive data or 

theoretical analysis. Includes compilations of 

significant scientific and technical data and 

information deemed to be of continuing 

reference value. NASA counterpart of peerreviewed formal professional papers, but having 

less stringent limitations on manuscript length 

and extent of graphic presentations. 

• TECHNICAL MEMORANDUM. Scientific 

and technical findings that are preliminary or of 

specialized interest, e.g., quick release reports, 

working papers, and bibliographies that contain 

minimal annotation. Does not contain extensive 

analysis. 

• CONTRACTOR REPORT. Scientific and 

technical findings by NASA-sponsored 

contractors and grantees. 

• CONFERENCE PUBLICATION. Collected 

papers from scientific and technical 

conferences, symposia, seminars, or other 

meetings sponsored or co-sponsored by NASA. 

• SPECIAL PUBLICATION. Scientific, 

technical, or historical information from NASA 

programs, projects, and missions, often 

concerned with subjects having substantial 

public interest. 

• TECHNICAL TRANSLATION. Englishlanguage translations of foreign scientific and 

technical material pertinent to NASA’s mission. 

Specialized services that complement the STI 

Program Office’s diverse offerings include creating 

custom thesauri, building customized databases, 

organizing and publishing research results ... even 

providing videos. 

For more information about the NASA STI Program 

Office, see the following: 

• Access the NASA STI Program Home Page at 

http://www.sti.nasa.gov

• E-mail your question via the Internet to 

help@sti.nasa.gov 

• Fax your question to the NASA STI Help Desk 

at (301) 621-0134 

• Phone the NASA STI Help Desk at 

(301) 621-0390 

• Write to: 

 NASA STI Help Desk 

 NASA Center for AeroSpace Information 

 7121 Standard Drive 

 Hanover, MD 21076-1320
NASA/CR-2005-213905 

NIA Report No. 2005-01 

National Aeronautics and 

Space Administration

Langley Research Center Prepared for Langley Research Center 

Hampton, Virginia 23681-2199 under Contract NAS1-02117

 

August 2005 

Modeling of Unit-Cells with Z-Pins Using 

FLASH: Pre-Processing and Post-Processing 

Ronald Krueger 

National Institute of Aerospace, Hampton, Virginia 
Available from: 

NASA Center for AeroSpace Information (CASI) National Technical Information Service (NTIS) 

7121 Standard Drive 5285 Port Royal Road 

Hanover, MD 21076-1320 Springfield, VA 22161-2171 

(301) 621-0390 (703) 605-6000 
1

MODELING OF UNIT-CELLS WITH Z-PINS USING FLASH:

PRE-PROCESSING AND POST-PROCESSING

Ronald Krueger1

ABSTRACT

Although the toughening properties of stitches, z-pins and similar structures have

been studied extensively, investigations on the effect of z-pins on the in-plane

properties of laminates are limited. A brief summary on the effect of z-pins on the

in-plane tensile and compressive properties of composite laminates is presented

together with a concise introduction into the finite element code FLASH. The

remainder of the report illustrates the modeling aspect of unit cells with z-pins in

FLASH and focuses on input and output data as well as post-processing of results.

LIST OF SYMBOLS

C length of resin pocket parallel to fiber direction

Ar fraction of the total reinforced area covered by z-pins

AZ cross sectional area of a single z-pin

Dz z-pin diameter

Dz’ height of z-pin plus resin pocket normal to fiber direction

d fiber diameter

E11, EL Stiffness of lamina parallel to fiber direction

E22, ET Stiffness of lamina transverse to fiber direction

ETc Compression stiffness of lamina transverse to fiber direction

G Shear modulus

Gf Fiber shear modulus

GLT, G12 Lamina shear modulus in principal material directions

Gsec Lamina secant shear modulus

Hz vertical spacing between z-pins in unit cell

Lz Horizontal spacing between z-pins in unit cell

n Ramberg-Osgood curve fitting parameter

Nx, Ny Axial force resultant on laminate in X,Y direction

Nxy Shear force resultant on laminate in X-Y plane

rz areal density of z-pins

ux, uy Displacement in X, Y direction

Vf Fiber volume fraction

w kink band width

α Ramberg-Osgood curve fitting parameter

β kink band inclination angle

Δx, Δy distance of fiber misalignment data points in X,Y direction

γe Effective shear strain

γs shear strain

 

1 Senior Staff Scientist, National Institute of Aerospace (NIA), 144 Research Drive, Hampton, VA 23666.
2

γp shear strain

γy Yield strain in shear

λ loading parameter

ν12 Lamina Poisson’s ratio

€

φ Fiber misalignment angle

σc Strength, critical value of stress

σTy Yield strength in tension

σxx, σyy Normal stress in X, Y direction

σxy, σyx Shear stress in X-Y plane

σult Strength of skin/stiffener-flange laminate

σultc Compression strength of skin/stiffener-flange laminate

σ11, σ22 Normal stress in 1, 2 direction

σ12, σ21 Shear stress in 1-2 plane

τe Effective yield strength in shear

τy Yield strength in shear

τxy, τyx Shear stress in X-Y plane

τ12, τ21 Shear stress in 1-2 plane

1. BACKGROUND

One of the most common failure modes for composite structures is delamination [1-3].

The remote loadings applied to composite components are typically resolved into interlaminar

tension and shear stresses at discontinuities that create delaminations under mixed-mode

conditions. In the past research has focussed on stitching to increase the delamination resistance

of composites [4-8]. More recently z-pins2 have been proposed to stitch the material together

through a combination of friction and adhesion [9-12]. Z-pins are pultruded rods of carbon fiber

and epoxy matrix with an approximate diameter of 280 μm or 508 μm. The z-pins are inserted

through the thickness of a laminated composite, which is then autoclaved to cure the epoxy as

normal. This approach to through-thickness reinforcement offers an alternative to stitching, and

can provide much higher areal densities of reinforcement [13].

Examination of the literature shows that the toughening properties of stitches, z-pins and

similar structures have been studied extensively, however, investigations on the effect of z-pins on

the in-plane properties of laminates are limited [13-17]. Steeves examined the effect of z-pins on the

in-plane tensile and compressive properties of composite laminates [13]. Disruption in the

alignment of the fibers in a fiber composite leads to a significant reduction in the in-plane

compressive strength [18]. Since the diameter of the z-pins (~280 μm) is large relative to the

diameter of the fibers (~7 μm) in the composite, the z-pins may cause significant misalignment of

the fibers of the composite. A sketch of the typical distribution of fiber misalignment around a z-pin

is shown in Figure 1. The z-pin is surrounded by a resin-rich pocket on two sides, caused by the

fibers deflecting around the z-pin and leaving a gap in the composite. The maximum misalignment

of the deflecting fibers occurs in four locations on the four flanks of the resin-rich pocket.

 

2 The generic term z-pin will be used throughout the paper versus the trade mark Z-FiberTM registered by Aztex Inc.
3

Steeves’ experiments and numerical finite element simulations using FLASH showed that

z-pins reduce the in-plane strength of composite laminates. His experimental data from

compressive tests showed that the insertion of a field of z-pins into a laminated fiber composite

reduces the in-plane compressive strength of the composite by 30% or more. He assumed that

the reduction in compressive strength is a consequence of the increased fiber waviness caused by

the z-pins [13]. Compression tests within the scanning electron microscope showed that the zpinned specimens develop microbuckles in the region of greatest fiber waviness around the resinrich pocket, as would be predicted by the microbuckling model of Budiansky and Fleck [19].

The objective of this report is to present a brief introduction into the finite element code

FLASH and illustrate the modeling aspect of unit cells with z-pins in FLASH. The report focuses

on the detailed description of input data for FLASH and output obtained from the analysis as

well as the post-processing of results. The report is intended to supplement existing user

manuals [20, 21]. First, input is described for a finite element model and material data from

reference [13]. Second, the correct load input and boundary conditions for compression, tension,

shear and combined loading is illustrated. Third, for a composite made of carbon/epoxy modeling

for large (0.508 mm diameter) and small (0.28 mm diamter) z-pins is described for 2% and 4% zpin areal density. Input data are determined for unit cells subjected to five different load cases

such as axial compression, combined axial compression and transverse tension, combined axial

compression and shear loading as well as combined axial compression, transverse tension and

shear loading. Additionally, the output obtained from FLASH is explained and separate postprocessing options are discussed. The interpretation of the output with respect to the strength of

z-pin reinforced lamina is discussed in detail in a separate report [22].

2. INTRODUCTION TO FLASH

In order to better assess the influence of critical parameters on lamina compression

strength, Fleck and Shu developed a finite element code called FLASH [20, 21, 23]. The FE code

is based on a 2D general Cosserat couple stress theory that assumes the unidirectional composite

lamina is a homogeneous anisotropic material that carries couple stress as well as classical

Cauchy point stress [24, 25]. The constitutive response is deduced from a unit cell consisting of

a fiber, represented by a linear elastic Timoshenko beam, embedded in a non-linear elastic-plastic

matrix [26]. The fiber diameter, d, is the length scale in the constitutive law that controls fiber

bending resistance. The continuum theory was implemented within a two-dimensional finite

element code that uses 6-noded triangular elements with 3 degrees of freedom at each node (twodisplacements and one rotation corresponding to rotation of the fiber cross section). The finite

element procedure is based upon a Lagrangian formulation of the finite deformation of the

composite and can accommodate both geometric and material non-linearities [23]. The code

models finite deformation using a Newton-Raphson incremental solution procedure with a

modified Riks algorithm in the final stage to handle snap-back behavior associated with fiber

micro-buckling [27]. Boundary loading is piecewise proportional with a loading parameter, λ, for

each loading stage [20, 21, 23].

The FLASH code assumes micro-buckling initiates from an imperfection in the form of

fiber waviness. Inputs include material properties (volume fraction Vf, and stiffness properties
4

EL, ET, ETc, GLT, Gf, all normalized by the shear yield strength, τy) and the Ramberg-Osgood

strain hardening law parameters (α, n). FLASH allows options for input of fiber misalignment

angle due to fiber waviness either as (1) an elliptical patch of waviness, or (2) an arbitrary

distribution of initial fiber waviness through initial misalignment angle, 

€

φ , at the Gauss

integration point for each element [20, 21, 23]. The first option prescribes the elliptical patch

along one edge of the unit cell. Steeves used the second option to input fiber misalignment

distribution obtained from analysis of digital images of specimens taken in a Scanning Electron

Microscope (SEM) [13].

3. GENERIC PRE-PROCESSING

3.1 Unit Cell Geometry

The z-pin reinforced composite is divided into unit cells as shown in Figure 2. The size of

the unit cell depends on the areal density, rz (in %) of the z-pins and the diameter of a single zpin, Dz as shown in Figure 2. The spacing Lz=Hz for a perfect, rectangular z-pin field can be

calculated as

€

LZ = HZ = AZ

Ar

 ,

where AZ is the cross sectional area of a single z-pin

 

€

AZ = π DZ

2

4

and Ar denotes the fraction of the total reinforced area covered by z-pins

€

Ar = rZ

100 .

The length of the resin pocket, C, may be determined from micrographs of the reinforced laminate

as shown in Figure 3.

3.2 Input of Geometrical Points

For the analysis with FLASH only one unit cell is modeled as shown in Figure 4. The resin-rich

pocket and the z-pin are modeled as voids within the material. Considering axial compression (xdirection) as the driving failure load, modeling a void is justified because the matrix strength and

stiffness are low relative to the strength and stiffness of the composite, and the resin-rich pocket

can considered to be significantly cracked early in the simulation. Hence, the simulation models

the unit cell of a composite including fiber misalignment and internal holes [13]. All meshes

generated were composed of six-noded triangular plane-strain elements. The size of the elements
5

was varied to provide the greatest mesh refinement near the resin pocket, and in the region of

greatest fiber misalignment.

First geometric parameters were taken from reference [13]. Since the fiber is completely

surrounded by resin the transverse dimension of the modeled void, D’z, is increased by 0.02 mm

compared to the z-pin diameter, Dz as shown in Figure 3a. The finite element code FLASH

requires the input of geometric data normalized with the fiber diameter of the composite, d as

listed in Table 1. Geometrical data for the areas to be meshed is input using the *GEOMETRICAL

POINTS command in FLASH. The coordinates of the points used to create the finite element

model of the unit cell with z-pin are printed below and their physical location is sketched in

Figure 5.

*GEOMETRICAL POINTS

 26

 1 0.00 0.00

 2 125.00 0.00

 3 250.00 0.00

 4 0.00 103.60

 5 125.00 103.60

 6 250.00 103.60

 7 110.00 105.00

 8 140.00 105.00

 9 0.00 124.00

 10 33.00 124.00

 11 217.00 124.00

 12 250.00 124.00

 13 32.00 125.00

 14 217.99 125.00

 15 0.00 126.00

 16 33.00 126.00

 17 217.00 126.00

 18 250.00 126.00

 19 110.00 145.00

 20 140.00 145.00

 21 0.00 146.40

 22 125.00 146.40

 23 250.00 146.40

 24 0.00 250.00

 25 125.00 250.00

 26 250.00 250.00

The geometric parameters were used to generate the finite element mesh of the unit cell including

the resin-rich pocket as shown in Figure 4.
6

3.3 Input of Mesh Data

Mesh seed data for the areas to be meshed is input using the *MESH SETS command in

FLASH. The input to create 10 mesh sets which form the finite element model of the unit cell

with z-pin (see Figure 4) are printed below. The physical location of the mesh sets is sketched in

Figure 6.

*MESH SETS

 10,6

 1, 1, 4, 2, 2, 5, 8, 10, 0.1, 0.1, 0.1, 0.1

 2, 2, 5, 3, 3, 6, 8, 10, 0.1, 0.1, 10, 10

 4, 4, 9, 5, 7, 10, 8, 10, .25, .25, 0.1, 0.1

 5, 8, 11, 6, 6, 12, 8, 10, .25, .25, 10, 10

 9, 9, 15, 10, 13, 16, 3, 10, 1, 1, 0.1, 0.1

 11, 14, 17, 12, 12, 18, 3, 10, 1, 1, 10, 10

 15, 15, 21, 16, 19, 22, 8, 10, 4, 4, 0.1, 0.1

 17, 20, 22, 18, 18, 23, 8, 10, 4, 4, 10, 10

 21, 21, 24, 22, 22, 25, 8, 10, 10, 10, 0.1, 0.1

 22, 22, 25, 23, 23, 26, 8, 10, 10, 10, 10, 10

3.4 Input of Misalignment Data

The finite element code FLASH also calls for the input of the fiber misalignment angle [20,

21]. The misalignment angle is input using the *MISALIGNMENT TYPE and the *MISALIGNMENT

command in FLASH. The input used for all analysis is printed below.

*MISALIGNMENT TYPE

 2

*MISALIGNMENT

 74,74,3.425,3.425

Preliminary analysis with N=1 as input for *MISALIGNMENT TYPE which was supposed

to generate an elliptical patch of waviness was not successful. Post-processing – discussed later -

only showed a small band of waviness along the left edge of the specimen (x=0). The waviness in

the rest of the specimen was zero. Therefore this approach was abandoned and N=2 was chosen

which requires that the discretized misalignment data is provided in a separate input file

(phiex.dat). Input of an arbitrary distribution of the fiber misalignment is possible as by

sampling a micrograph as shown in Figure 7 [13]. However, if this information is not available,

misalignment angles may be chosen to study the misalignment effect on the reduction of in-plane

strength.

To illustrate the input an example is shown in Figure 8 for 10 by 10 sample points on the

unit cell which corresponds to a matrix of 10 by 10 data points in the phiex.dat file. The

distance of one point to the adjacent points was set to Δx=26.0 units in x and Δy=26.0 units in y-
7

direction which extends beyond the unit cell (250 by 250 units) modeled. The distance of points

has to be chosen so that the field of data points is larger than the unit cell. The user is advised to

check the misalignment input after the execution of prep has generated the flash.inp file and

adjust the distance of adjacent points accordingly. If the misalignment data originates from

sampling a micrograph the sample area has to be larger than the model. See reference [21] for

details. Only the file named phiex.dat will be used during the initial step of the analysis when

the file flash.inp. is generated during the execution of prep.

For the current analyses constant fiber misalignment values of 0°, 1°, 2°, 3°, 4°, 5°, 6°,

7°, 8°, 9°, and 10°, were chosen and 11 phiex.dat files were created. Each phiex.dat file

contains a field of 74 by 74 data points which corresponds to a matrix of 74 by 74 sample points

on the unit cell as shown for the example of 5° fiber misalignment.

5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000

5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000

5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000

5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000

5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000

5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000 5.000

5.000 5.000

5.000 5.000 5.000 5.000 5.000 5.000 ...

...

...

The distance of one point to the adjacent points was set to 3.425 units in x and 3.425 units in ydirection to cover the entire unit cell (250 by 250 units).

3.5. Input of Constraints

As mentioned above the composite is divided into unit cells, and only one cell is modeled

using periodic boundary conditions. For these models, periodic boundary conditions are

prescribed by applying compressive loading symmetrically to the specimen, and constraining

three degrees of freedom to prevent rigid body translation and rotation as shown in Figure 9. The

constraints are input using the *TYPE 1 CONSTRAINTS command in FLASH. The input used for

all analysis is printed below.

*TYPE 1 CONSTRAINTS

 3

 1, 1, 2, 21, 1, 1, 0

 9, 9, 15, 4, 4, 2, 0

 24, 24, 25, 21, 1, 1, 0
8

3.6. Input of Load Data

The loading was input using the *STRESS LOADING command in FLASH. The axial

compression stress sketched in Figure 9 is gradually incremented by FLASH until it reaches the

specified limit defined by the user (σxx/τy =-1000). The limit was deliberately chosen to be well

above failure, i.e. formation of kink bands and fiber microbuckling to assure that the analysis

reached the failure point and did not terminate at a chosen lower stress [20]. The input used for

axial compression loading is printed below.

*STRESS LOADING

 10, 1

 1, 1, 4, 1, 1

 -1000, 0, 0, 0, 0

 4, 4, 9, 1, 1

 -1000, 0, 0, 0, 0

 9, 9, 15, 1, 1

 -1000, 0, 0, 0, 0

 15, 15, 21, 1, 1

 -1000, 0, 0, 0, 0

 21, 21, 24, 1, 1

 -1000, 0, 0, 0, 0

 3, 3, 6, 1, 1

 -1000, 0, 0, 0, 0

 6, 6, 12, 1, 1

 -1000, 0, 0, 0, 0

 12, 12, 18, 1, 1

 -1000, 0, 0, 0, 0

 18, 18, 23, 1, 1

 -1000, 0, 0, 0, 0

 23, 23, 26, 1, 1

 -1000, 0, 0, 0, 0

3.7. Requests for Output

During the analysis stresses, strains and total fiber alignment angle is written to the file

flash.out after every N load steps. The data output request is managed using the *OUTPUT

FREQUENCY commands in FLASH as printed below for an output request every 1000 steps.

*OUTPUT FREQUENCY

 1000

During the analysis, for three boundaries (as shown in Figure 9) data are output to the

output file post.d after each load step. For boundary sections 1 and 2 the net force on each

boundary section, resolved along the nominal fiber direction, are given as output. Boundary 3 is a
9

single node and the displacements of this node are given as output. The data output request is

managed using the *OUTPUT LIST command in FLASH as printed below. The boundary sections

selected are shown in Figure 5.

*OUTPUT LIST

 3

 1, 1, 4, 1, 1

 1, 1, 4, 1, 1

 1, 1, 4, 1, 1

Additionally, the value of the global loading parameter λ is written to the output file at

each loading step.

The format of the output files and post-processing of the results is discussed in a

separate section below.

3.8. Input of Material Data

The finite element code FLASH requires the input of material data normalized with the

composite yield shear strength τy as given in Table 2 for a generic graphite epoxy material [13].

Additional material parameters used by FLASH are the constants α=0.429 and n=3 in the

Ramberg-Osgood strain hardening law, and the fiber volume fraction Vf =0.63 of the composite.

The material data is input using the *MATERIAL CONSTANTS and *VFRACTION commands in

FLASH as printed below.

*MATERIAL CONSTANTS

 7

 1200, 175, 74, 41, 0.429, 3.0, 87

*VFRACTION

 0.63

3.9. Input for a Finite Element Model without Z-Pin

For reference purposes a unit cell without z-pin was modeled and analyzed with FLASH.

The mesh with applied load and boundary condition is shown in Figure 10. The geometrical points

and the required input data to generate the mesh have been included in Figure 10.

4. PREPROCESSING FOR COMBINED LOADING

For the initial analyses performed, load and boundary conditions were taken from

reference [13]. Compressive loading was applied symmetrically to the specimen, and three

degrees of freedom were constrained to prevent rigid body translation and rotation as shown in

Figure 9.

In general, z-pin reinforced composite structures are subjected to multi axial loads and
10

input stresses in the unidirectional ply modeled in FLASH need to be calculated. For composite

plates subjected to external compression and combined compression and shear loading input

stresses for FLASH were calculated. Using classical laminate theory the external loads were

resolved into stresses for the individual plies as shown in Figure 11. Calculated stresses in the

critical unidirectional ply consist of axial compression, transverse compression and tension as well

as shear stresses. For these combined loading conditions the above mentioned boundary

conditions and load input had to be modified.

4.1. FE Models of Z-Pin Unit Cells Subjected to Axial and Transverse Loading

The finite element mesh with load and boundary conditions for axial and transverse

compression is shown in Figure 12 with the corresponding input data. Based on results from

classical laminate theory mentioned above, the transverse compression load was chosen to be 2% of

the compression load applied axially. In order to avoid the collapse of the void under transverse

compression points A and B on the periphery of the void could be constrained in transverse

direction as shown in Figure 12. In reality, the boundary condition of zero displacements on the

periphery of the void is not met due to pin compliance. However, the assumed boundary

condition is adopted in order to obtain an estimate for the extreme case of a rigid pin [28].

Nevertheless, analysis without constraints suggested that constraining points A and B for small

transverse compression loads leads to overly conservative results.

The finite element mesh with load and boundary conditions for axial compression and

transverse tension are shown in Figure 13. As for the previous case the transverse load was chosen

to be 2% of the compression load applied axially. For transverse tension it is not required to avoid

the collapse of the void.

4.2. FE Models of Z-Pin Unit Cells Subjected to Axial Compression and Shear Loading

The finite element mesh with load and boundary conditions for axial compression

combined with shear are shown in Figure 14. Based on results from classical laminate theory

mentioned above, the shear load was chosen to be 10% of the compression load applied axially. For

the combined axial compression and shear load case it was decided not to prevent the collapse of

the void as explained in detail in the appendix of reference 23.

The finite element mesh with load and boundary conditions for axial compression and

transverse tension combined with shear are shown in Figure 15. As for the previous case the shear

load was chosen to be 10% of the compression load applied axially and the transverse tension load

corresponded to 2% of the compression load as discussed in the previous section. For this

combined load case it was decided not to prevent the collapse of the void as explained in detail in

the appendix of reference 23.

5. UNIT CELL MODELS FOR 2% AND 5% Z-PIN REINFORCEMENT.

For another carbon fiber reinforced epoxy material the data used as input for the analyses

with FLASH are listed in Table 3. The specimens were reinforced with small z-pins of
11

Dz=0.279 mm (0.011 in) and large z-pins of Dz=0.508 mm (0.02 in) diameter. Three specimen

types were manufactured containing reinforcement fields with rz=2% areal density for the large

z-pin and rz=2% and rz=4% areal density for the small z-pins respectively. Geometric data such

as z-pin spacing, Lz=Hz, as well as width, D’z, and length, C, of the resin pocket as defined in

Figure 3 was measured from micrographs taken of different specimens. The averaged data used as

input for the analysis is summarized in Table 4.

Based on results from classical laminate theory shown in Figure 11 and Table 5, five load

cases were considered as input for analysis of the unit cells: A pure axial compression load case, a

combined axial compression - 2% transverse tension load case, a combined axial compression -

10% shear load case, a combined axial compression - 50% shear load case, and an axial

compression - 2% transverse tension load case combined with 10% shear loading. The input for all

five cases is summarized in Table 6.

The finite element meshes are shown in Figures 16, 17, and 18 for the small pin with 2%

and 4% areal density and the large pin with 2% areal density. As before the resin-rich pocket and

the z-pin are modeled as voids within the material. It is assumed that the fiber is completely

surrounded by resin as shown in Figure 3 and therefore the transverse dimension of the modeled

void, D’z, is increased by 0.02 mm compared to the z-pin diameter, Dz. The normalized dimensions

required as input were included in the figures as well as the input for the *GEOMETRICAL POINTS

command in FLASH. Loads given in Table 6 were applied as discussed in section 4.

6. POSTPROCESSING OF RESULTS

6.1. FLASH Output Files

During the analysis nodal point coordinates, element topology and results are is written

to the file flash.out after every N load steps. The file flash.out contains a header,

1Ccomnew 0.00000E+00 FIBER MICROBUCKLING 3

the nodal point coordinates

2C

 -1 1 0.00000E+00 0.00000E+00

 -1 2 0.00000E+00 0.15646E+02

 -1 3 0.00000E+00 0.31293E+02

 -1 4 0.00000E+00 0.42553E+02

 -1 ...

 -3

the element topology

3C

 -1 1 8

 -2 1 374 372 155 373 154
12

 -1 2 8

 -2 1 3 374 2 156 155

 -1 3 8

 -2 372 710 708 581 709 580

 -1 ...

 -3

the initial fiber misalignment 

€

φ

-3

 100CLCASE1 0.00000E+00 Fiber Misalignment 3 0

 -4 MISALIGN 1 3 0

 -5 PHIBAR 1 1 1 0 0

 -1 1 8 0 0 6

 -2 1 0.40000E+01

 -2 2 0.40000E+01

 -2 3 0.40000E+01

 -2 4 0.40000E+01

 -2 5 0.40000E+01

 -2 6 0.40000E+01

 -1 2 8 0 0 6

 -2 1 0.40000E+01

 -2 2 0.40000E+01

 -2 3 0.40000E+01

 -2 4 0.40000E+01

 -2 ...

 -3

the normalized displacements ux/d and uy/d at all nodal points

100CLCASE1 0.12052E-01 Fiber Microbuckling 3 1000

 -4 DISPLACE 3 1 0

 -5 UX 1 2 1 0 0

 -5 UY 1 2 2 0 0

 -5 U 1 2 0 0 1ALL 

 -1 1 0.14712E+01 -0.32663E+01

 -1 2 0.80530E+00 -0.32449E+01

 -1 3 0.14974E+00 -0.31895E+01

 -1 4-0.31371E+00 -0.31129E+01

 -1 ...

 -3

the normalized stresses σxx/τy, σyy/τy, σxy/τy and σyx/τy, the normalized effective shear stress
13

τe/τy, the strains γs and γp, the fiber misalignment 

€

φ and others [20]

100CLCASE1 0.12052E-01 Fiber Microbuckling 3 1000

-4 STRSSTRN 12 3 0

-5 SIG_XX 1 4 1 1 0

-5 SIG_YY 1 4 2 2 0

-5 SIG_XY 1 4 1 2 0

-5 SIG_YX 1 4 2 1 0

-5 TAU_E 1 1 2 1 0

-5 E_LP 1 1 2 1 0

-5 E_LM 1 1 2 1 0

-5 E_T 1 1 2 1 0

-5 GAMA_S 1 1 2 1 0

-5 KAPPA 1 1 2 1 0

-5 GAMA_P 1 1 2 1 0

-5 PHI 1 1 2 1 0

-1 1 8 0 0 6

-2 1 -0.12328E+02 0.20839E+00 -0.11766E+00 -0.13105E+00 0.14941E+01 -0.30784E-02

-2 1 -0.25050E-02 0.17207E-03 0.20919E-01 -0.57344E-03 0.10263E-01 0.76157E+01

-2 2 -0.11658E+02 0.44093E+00 -0.74204E+00 -0.73252E+00 0.74208E+00 -0.28156E-02

-2 2 -0.25441E-02 0.18100E-02 0.60815E-02 -0.27141E-03 0.13715E-02 0.66904E+01

-2 3 -0.13010E+02 0.17288E-01 -0.74737E+00 -0.74508E+00 0.79505E+00 -0.30924E-02

-2 3 -0.28856E-02 0.63897E-04 0.72372E-02 -0.20685E-03 0.15640E-02 0.67907E+01

-2 4 -0.11966E+02 0.32614E+00 -0.32097E+00 -0.30034E+00 0.12303E+01 -0.29931E-02

-2 4 -0.24544E-02 0.11352E-02 0.14406E-01 -0.53866E-03 0.57257E-02 0.72148E+01

 -2 ...

 -3

During the analysis, data for three boundaries are output to the output file post.d The

file post.d contains contains seven columns. Column 1 lists the load step number, column 2

lists load increase (+1) or decrease (-1), columns 3 and 4 list the normalized displacements ux/d

and uy/d for the selected node on boundary 3, columns 5 and 6 list the net loads on boundary

sections 1 and 2 as shown in Figure 9. Column 7 lists global loading parameter λ for each load

step, where λ is a measure of how close the simulation is to fulfilling the requested stress loading

discussed in section 3.6.

 0 1 0.00000E+00 0.00000E+00 0.00000E+00 0.00000E+00 0.00000E+00

 1 1 -0.24101E-01 -0.13840E-01 0.10106E+03 0.10106E+03 0.96809E-03

 2 1 -0.37684E-01 -0.22066E-01 0.15836E+03 0.15836E+03 0.15167E-02

 3 1 -0.52874E-01 -0.31723E-01 0.22269E+03 0.22269E+03 0.21325E-02

 4 1 -0.69664E-01 -0.43129E-01 0.29416E+03 0.29416E+03 0.28163E-02

...
14

 996 1 -0.49693E+00 0.21751E+00 0.19896E+04 0.19896E+04 0.18922E-01

 997 1 -0.49626E+00 0.21723E+00 0.19925E+04 0.19925E+04 0.18949E-01

 998 1 -0.49530E+00 0.21693E+00 0.19966E+04 0.19966E+04 0.18988E-01

 999 1 -0.49464E+00 0.21683E+00 0.19995E+04 0.19995E+04 0.19015E-01

1000 1 -0.49373E+00 0.21680E+00 0.20036E+04 0.20036E+04 0.19053E-01

6.2. Result Visualization

The output file flash.out is formatted to be used with the commercial post-processing

software FEMGV3. A translation routine (make_post.f) was developed, which reads

flash.out and creates output data of the model (flash_pat.out) and displacements

(flash.dis) in neutral file format to be used with PATRAN4. The routine can also create output

of the model and result data (flash.dat) to be used with Tecplot5. The routine can easily be

modified to create an output file in a particular format that is compatible with a graphics package

of the user’s choice.

6.2.1 Examples of Deformation Plots

The deformed finite element meshes for the unit cell under axial compression are shown in

Figure 19. At the beginning of the loading a symmetric deformation with respect to the two axis

of symmetry can be observed as shown in Figure 19a. At the edge where the axial compression

stress is applied a non-uniform displacement distribution can be observed which forms an axial

contraction. The increased displacement is caused by a reduced stiffness which is the result of

modeling the resin pocket and z-pin as a void. At the end of the analysis after passing the failure

point and the formation of kink bands the deformation grows to be unsymmetrical and the unit

cell becomes increasingly distorted as shown in Figure 19b.

The deformed mesh for axial and transverse compression is shown in Figure 20. The

transverse compression load was 2% of the compression load applied axially as shown in

Figure 12. Analysis without constraints as shown in Figure 20 suggested that constraining points

A and B shown in Figure 12, was not required for small transverse compression loads. Compared

to the deformation shown in Figure 19a for the pure compression case the axial displacement

shown in Figure 20 is reduced at the top and bottom edges where the transverse compression is

applied and the void remains open.

The deformed mesh for axial compression and transverse tension are shown in Figure 21.

As for the previous case the transverse load was chosen to be 2% of the compression load applied

axially as shown in Figure 13. Compared to the deformation for the pure compression case as

shown in Figure 19a the axial displacement shown in Figure 21 has increased at the top and

bottom edges where the transverse tension stress is applied; also the void opening has increased.

 

3 http://www.femsys.co.uk/

4 http://www.mscsoftware.com/products/quick_prod.cfm

5 http://www.tecplot.com/
15

6.2.2 Examples of Contour Plots

As an example the contour plots of the fiber misalignment 

€

φ , the effective shear stress τe and the

strains γs and γp, are shown in Figures 24-27 for the unit cell subjected to axial and transverse

compression. The contours correspond to the end of the analysis after passing the failure point

where the severe color gradients indicate the formation and location of kink bands.

6.2.3 Examples of Stress-Displacement Plots to Obtain Compression Strength

As an example typical stress versus displacements plots are shown in Figures 28 and 29. Typically

the post.d file is read into a graphing program such as KaleidaGraphTM or ExcelTM. Here the

stress is obtained by multiplying the global loading parameter λ from column (#7) of the output

file post.d with the matrix yield stress τy from the Ramberg-Osgood strain hardening law that is

used as input and the requested stress loading (-1000) described in section 3.6. The compression

stress σxx= λ τy 1000 is plotted versus the normalized displacement in x-direction (ux/d) for the

selected node on boundary 3 as shown in Figure 9. The normalized displacement (ux/d) is listed in

column (#3) of the output file post.d. As shown in Figures 28 and 29, the plot typically indicates a

clear first maximum. The exact maximum value, which is the predicted compression strength, is best

found in the data rather than the plot.

ACKNOWLEDGEMENTS.

This research was supported and partially funded by Sikorsky Aircraft and the Aviation

Applied Technology Directorate under Technology Investment Agreement No. DAAH10-02-2-

0001 as part of the Survivable – Affordable - Repairable – Airframe - Program (SARAP).

The author gratefully acknowledge Prof. Noman A. Fleck of University of Cambridge,

United Kingdom and Dr. Craig A. Steeves of Princeton University for providing the FLASH

finite element code and insight into its use and application.

REFERENCES

[1] T.K. O'Brien, Fracture Mechanics of Composite Delamination, in ASM Handbook, Volume

21, Composites: ASM International, 2001, pp. 241-245.

[2] V.V. Bolotin, Delaminations in composite structures: its origin, buckling, growth and

stability, Composites Part B: Engineering, vol. 27B, pp. 129--145, 1996.

[3] T.E. Tay, Characterization and analysis of delamination fracture in composites - An overview

of developments from 1990 to 2001, Applied Mechanics Reviews, vol. 56, pp. 1-32, 2003.

[4] L.K. Jain, K.A. Dransfield, and Y.-W. Mai, On the effect of stitching in CFRPs - II: mode II

delamination toughness, Composites Science and Techonology, vol. 58, pp. 829-837, 1998.

[5] K.A. Dransfield, L.K. Jain, and Y.-W. Mai, On the effect of stitching in CFRPs - I: mode I

delamination toughness, Composites Science and Technology, vol. 58, pp. 815-827, 1998.

[6] D. Shu and Y.-W. Mai., Effect of stitching on interlaminar delamination extension in

composite laminates, Composites Science and Technology, vol. 49, pp. 165-171, 1993.

[7] L.K. Jain and Y.-W. Mai., Determination of mode II delamination toughness in stitched

laminated composites, Composites Science and Techonology, vol. 55, pp. 241-253, 1995.
16

[8] B.N. Cox, Constitutive model for a fibre tow bridging a delamination crack., Mechanics of

Composite Structures and Materials, vol. 6, pp. 117-138, 1999.

[9] D. Cartie and I. Partridge., Z-pinned composite laminates: Improvements in delamination

resistance, presented at 5th International Conference on Deformation and Fracture of

Composites, 1999.

[10] G. Freitas, T. Fusco, T. Campbell, J. Harris, and R. S., Z-fiber technology and products for

enhancing composite design., presented at AGARD Conference Proceedings 590:

Bolted/Bonded Joints in Polymeric Composites, 1996.

[11] B. Budiansky and J.C. Amazigo, Toughening by aligned, frictionally constrained fibers,

Journal of Mechanics and Physics of Solids, vol. 37, pp. 93-109, 1989.

[12] J.D. Barrett, The mechanics of z-fiber reinforcement, Composite Structures, vol. 36, pp. 23-32,

1996.

[13] C.A. Steeves, Mechanics of Failure in Composite Structures, in Department of Engineering.

Cambridge, United Kingdom: University of Cambridge, 2001.

[14] R. Massabo and B.N. Cox, Unusual Characteristics of Mixed-Mode Delamination Fracture in

the Presence of Large-Scale Bridging, Mechanics of Composite Materials and Structures,

8:61–80, 2001, vol. 8, pp. 61-80, 2001.

[15] K.L. Rugg, B.N. Cox, and R. Massabò, Mixed mode delamination of polymer composite

laminates reinforced through the thickness by Z-fibers., Composites A: applied science and

manufacturing, vol. 33, pp. 177–190, 2002.

[16] M. He and B.N. Cox, Crack bridging by through-thickness reinforcement in delaminating

curved structures, Composites Part A: Applied Science and Manufacturing, vol. 29, pp. 377-

393, 1998.

[17] R. Massabo and B.N. Cox, Concepts for bridged Mode ii delamination cracks, Journal of the

Mechanics and Physics of Solids, vol. 47, pp. 1265-1300, 1997.

[18] N.A. Fleck, Compressive failure in fiber composites, Advances in Applied Mechanics, vol. 33,

pp. 43-117, 1997.

[19] B. Budiansky and N.A. Fleck., Compressive failure of fibre composites, Journal of Mechanics

and Physics of Solids, vol. 41, pp. 183-211, 1993.

[20] J.Y. Shu and N.A. Fleck., User's manual for finite element code for fibre microbuckling.

Technical Report CUED/C-MATS/TR.224, University of Cambridge Department of

Engineering. 1995.

[21] D. Liu and N.A. Fleck., User's manual II for finite element code FLASH for fibre

microbuckling,Technical Report CUED/C-MICROMECH/TR.29, Cambridge University,

Engineering Department 1999.

[22] T.K. O’Brien, R. Krueger, and J. Schaff, Influence of compression and shear on the strength

of composite laminates with Z-pinned reinforcement, NASA Langley Research Center,

Hampton, VA, NASA/TM-2005-XXXXX, 2005.

[23] N.A. Fleck and J.Y. Shu, Microbuckle initiation in fibre composites: a finite element study,

Journal of Mechanics and Physics of Solids, vol. 43, pp. 1887-1918, 1995.

[24] E. Cosserat and F. Cosserat, Theorie des Corps Deformables. Paris.: Herman et fils, 1909.

[25] W.T. Koiter, Couple stresses in the theory of elasticity, I and II., Proc. Ned. Akad. Wet., vol.

67, pp. 17-44, 1964.

[26] S.P. Timoshenko and J.N. Goodier, Theory of Elasticity, 3rd ed: McGraw-Hill, 1970.

[27] M.A. Crisfield, Non-linear Finite Element Analysis of Solids and Structures, vol. 1, chapter 9:

John Wiley & Sons, U.K., 1991.

[28] N.A. Fleck, D. Liu, and J.Y. Shu, Microbuckle initiation from a hole and from the free edge of

a fibre composite, International Journal of Solids and Structures, vol. 37, pp. 2757-2775,

2000.
17

LIST OF TABLES

Table 1.

Dimensions of Unit Cell

Graphite/Epoxy UD Prepreg with 2% Z-Pins [13]

from normalized with df

DZ 0.28 mm

D’Z 0.3 mm 42.8

HZ 1.75 mm 250

LZ 1.75 mm 250

C 1.29 mm 184

d 7 μm 1
18

Table 2.

Material Properties for Initial Preliminary Analysis

Graphite/Epoxy Prepreg (A) [13]

from normalized with τxy

E11 117 GPa 1083

E22 (tension) 9.0 GPa 83.33

E22 (compression) 9.5 GPa 87.96

G12 4.8 GPa 44.44

Gf 22 GPa 203.7

τy 108 MPa 1

d 7 μm -

Vf 0.55 -

α 3/7 -

n 3 -

Generic Graphite/Epoxy Prepreg (B) [13]

normalized with τxy

E11 - 1200

E22 (tension) - 74

E22 (compression) - 87

G12 - 41

Gf - 175

τy - -

d 7 μm -

Vf 0.63 -

α 0.429 -

n 3 -
19

Table 3.

Carbon/Epoxy Material Properties

Graphite/Epoxy UD Prepreg

from normalized with τxy

E11 161 GPa 4248

E22 (tension) 11.4 GPa 301

E22 (compression) 12.8 GPa 338

G12 5.17 GPa 136

Gf 22 GPa 580

τy 37 MPa 1

d 5.1 μm -

Vf 0.59 -

α 0.00923 -

n 8.54 -
20

Table 4.

Dimensions of Unit Cells

Graphite/Epoxy UD Prepreg with 2% large diameter Z-Pins

from normalized with d

DZ 0.508 mm -

D’Z 0.528 mm 103.53

HZ 3.175 mm 622.55

LZ 3.175 mm 622.55

C 2.1844 mm 428.31

Graphite/Epoxy UD Prepreg with 4% small diameter Z-Pins

from normalized with d

DZ 0.28 mm -

D’Z 0.3 mm 58.8

HZ 1.2446 mm 244

LZ 1.2446 mm 244

C 0.868 mm 170.2

Graphite/Epoxy UD Prepreg with 2% small diameter Z-Pins

from normalized with d

DZ 0.28 mm -

D’Z 0.3 mm 58.8

HZ 1.7526 mm 343.65

LZ 1.7526 mm 343.65

C 0.868 mm 170.2
21

Table 5.

Load Cases for Laminate Analysis and Results

 external load Nx=-1000 lbs/in

σ11, psi σ22, psi τ12, psi

[0/90]s -1898 -29.5 (~2% σ11) 0

[0/±45]s -2304 +48.9 (~2% σ11) 0

[0/45/-45/90]s -2597 -3.3 (~0.2% σ11) 0

external load Nx=-1000 lbs/in, Nxy=-1000 lbs/in

σ11, psi σ22, psi τ12, psi

[0/90]s -1898 -29.5 -1000 (~50% σ11)

[0/±45]s -2304 +48.9 (~2% σ11) -222 (~10% σ11)

[0/45/-45/90]s -2597 -3.3 -275 (~10% σ11)

Table 6.

FLASH Input for Load Cases Used for Strength Reduction Analysis

axial

compression

compression/

2%transverse

tension

compression

10% shear

compression

50% shear

compression

2% tension

10% shear

σ11/τy -1000 -1000 -1000 -1000 -1000

σ22/τy - +20 - - +20

τ12/τy - - 100 500 100

τ21/τy - - 100 500 100
x

y

area of maximum 

waviness on flanks of 

resin-rich pocket

z-pin

direction of longitudinal fibers

Figure 1: Typical shape of local area around z-pin 

resin rich region

22
unit cell

 Lz

 Hz

 Dz

 C

fiber orientation

C size of resin pocket

Dz pin diameter

D’z pin + resin

Hz vertical spacing

Lz horizontal spacing

d fiber diameter

 D’z 

Figure 2: Detail of perfect z-pin field with typical dimensions and unit cell

Figure 3: Micrograph of individual z-pin with surrounding laminate [13]

0.1 mm

z-pin

resin pocket

C 

23
Figure 4: Typical finite element mesh of unit cell with normalized dimensions and

initial misalignment due to z-pin insertion

x

y

 D’z Hz / f / f 

 C / f 

 Lz / f 

24
Figure 5: Outline of unit cell and void with geometry points

15

21

x

y

1 2 3

4

5

6

7 8

9

10 11

12 13 14

16 17

18

19 20

22

23

24 25 26

*GEOMETRICAL POINTS

 26

 1 0.00 0.00

 2 125.00 0.00

 3 250.00 0.00

 4 0.00 103.60

 5 125.00 103.60

 6 250.00 103.60

 7 110.00 105.00

 8 140.00 105.00

 9 0.00 124.00

 10 33.00 124.00

 11 217.00 124.00

 12 250.00 124.00

 13 32.00 125.00

 14 217.99 125.00

 ...

 ... 

 26 250.00 250.00

25
x

y

Figure 6: Outline of unit cell and void with geometry points and mesh sets

1 2 3

4

5

6

7 8

9

10 11

12

15

16 17

18

19 20

21

22

23

24 25 26

*MESH SETS

 10,6

 1, 1, 4, 2, 2, 5, 8, 10, 0.1, 0.1, 0.1, 0.1

 2, 2, 5, 3, 3, 6, 8, 10, 0.1, 0.1, 10, 10

 4, 4, 9, 5, 7, 10, 8, 10, .25, .25, 0.1, 0.1

 5, 8, 11, 6, 6, 12, 8, 10, .25, .25, 10, 10

 9, 9, 15, 10, 13, 16, 3, 10, 1, 1, 0.1, 0.1

 11, 14, 17, 12, 12, 18, 3, 10, 1, 1, 10, 10

 15, 15, 21, 16, 19, 22, 8, 10, 4, 4, 0.1, 0.1

 17, 20, 22, 18, 18, 23, 8, 10, 4, 4, 10, 10

 21, 21, 24, 22, 22, 25, 8, 10, 10, 10, 0.1, 0.1

 22, 22, 25, 23, 23, 26, 8, 10, 10, 10, 10, 10

13 14

26
(b). Contour map of fiber alignment angle [13]

Figure 7: Measurement of fiber misalignment angle

(a). Electron micrograph of individual z-pin with surrounding laminate [13]

0.1 mm

extent of contour map

27
Figure 8: Outline of unit cell with geometry points and fiber misalignment data points

1 2 3

4

5

6

7 8

9

10 11

12

13 14

15

16 17

18

19 20

21

22

23

24 25 26

*MISALIGNMENT TYPE

 2

*MISALIGNMENT

 10,10,26.0,26.0

x

y

x

y

*GEOMETRICAL POINTS

 26

 1 0.00 0.00

 2 125.00 0.00

 3 250.00 0.00

 4 0.00 103.60

 5 125.00 103.60

 6 250.00 103.60

 7 110.00 105.00

 ...

 ... 

 26 250.00 250.00

28
x

y

uy=0.0

ux=0.0

ux=0.0

Figure 9: Finite element model subjected to axial compression loading 

xx

*STRESS LOADING

 10, 1

 1, 1, 4, 1, 1

 -1000, 0, 0, 0, 0

 4, 4, 9, 1, 1

 -1000, 0, 0, 0, 0

 9, 9, 15, 1, 1

 -1000, 0, 0, 0, 0

 15, 15, 21, 1, 1

 -1000, 0, 0, 0, 0

 21, 21, 24, 1, 1

 -1000, 0, 0, 0, 0

 3, 3, 6, 1, 1

 -1000, 0, 0, 0, 0

 6, 6, 12, 1, 1

 -1000, 0, 0, 0, 0

 12, 12, 18, 1, 1

 -1000, 0, 0, 0, 0

 18, 18, 23, 1, 1

 -1000, 0, 0, 0, 0

 23, 23, 26, 1, 1

 -1000, 0, 0, 0, 0

*TYPE 1 CONSTRAINTS

 3

 1, 1, 2, 21, 1, 1, 0

 9, 9, 15, 4, 4, 2, 0

 24, 24, 25, 21, 1, 1, 0

xx

*OUTPUT LIST

 3

 1, 1, 4, 1, 1

 1, 1, 4, 1, 1

 1, 1, 4, 1, 17

section 1,2

section 3

29
Figure 10: Mesh with axial compression load and boundary contions for UD Ply without Z-Pin

x

y

ux=0.0

ux=0.0

uy=0.0

xx xx

30
Figure 11: Calculation of ply stresses using classical laminate theory 

y

x

z

laminate critical 0 ply

classical laminate Nx theory (CLT)

Nx

Ny

Ny

Nxy

Nyx xy = 12

yy = 22

xx = 11

yx = 21

two load cases 

Nx=1000 lbs/ in

Nx=1000 lbs/ in, Nxy= 1000 lbs/ in

three lamiantes

[0/90]s

[0/±45]s

[0/45/-45/90]s

four critical cases identified

xx (compression) only

xx (compression), yy (tension) = 2% xx

xx (compression), xy = yx =10% xx

xx (compression), xy = yx =10% xx, , yy (tension) 

= 2% xx

2

1

xx = 11

31
x

y

uy=0.0

ux=0.0

ux=0.0

Figure 12: Finite element model subjected to axial and transverse compression loading 

xx

yy

*STRESS LOADING

 14, 1

 1, 1, 4, 1, 1

 -1000, 0, 0, 0, 0

 4, 4, 9, 1, 1

 -1000, 0, 0, 0, 0

 9, 9, 15, 1, 1

 -1000, 0, 0, 0, 0

 ...

 ...

 23, 23, 26, 1, 1

 -1000, 0, 0, 0, 0

 1, 1, 2, 1, 1

 0, -20, 0, 0, 0

 2, 2, 3, 1, 1

 0, -20, 0, 0, 0

 24,24,25, 1, 1 

 0, -20, 0, 0, 0

 25,25,26, 1, 1

 0, -20, 0, 0, 0

*TYPE 1 CONSTRAINTS

 3

 1, 1, 2, 21, 1, 1, 0

 9, 9, 15, 4, 4, 2, 0

 24, 24, 25, 21, 1, 1, 0

 5, 8, 11, 1, 17, 2, 0

 16, 19, 22, 17, 1, 2, 0

yy

xx

A

B

32
x

y

uy=0.0

ux=0.0

ux=0.0

Figure 13: Finite element model subjected to axial compression and transverse tension loading 

*TYPE 1 CONSTRAINTS

 3

 1, 1, 2, 21, 1, 1, 0

 9, 9, 15, 4, 4, 2, 0

 24, 24, 25, 21, 1, 1, 0

*STRESS LOADING

 14, 1

 1, 1, 4, 1, 1

 -1000, 0, 0, 0, 0

 4, 4, 9, 1, 1

 -1000, 0, 0, 0, 0

 9, 9, 15, 1, 1

 -1000, 0, 0, 0, 0

 15, 15, 21, 1, 1

 -1000, 0, 0, 0, 0

 ...

 ...

 23, 23, 26, 1, 1

 -1000, 0, 0, 0, 0

 1, 1, 2, 1, 1

 0, 20, 0, 0, 0

 2, 2, 3, 1, 1

 0, 20, 0, 0, 0

 24,24,25, 1, 1 

 0, 20, 0, 0, 0

 25,25,26, 1, 1

 0, 20, 0, 0, 0

xx xx

yy

yy

33
x

y

uy=0.0

ux=0.0

ux=0.0

Figure 14: Finite element model subjected to axial compression and shear loading 

xy

yx

xx xx

*STRESS LOADING

 14, 1 

 1, 1, 4, 1, 1

 -1000, 0, -100, 0, 0

 4, 4, 9, 1, 1

 -1000, 0, -100, 0, 0

 9, 9, 15, 1, 1

 -1000, 0, -100, 0, 0

 15, 15, 21, 1, 1

 -1000, 0, -100, 0, 0

 ...

 ... 

 23, 23, 26, 1, 1

 -1000, 0, -100, 0, 0

 1, 1, 2, 1, 1

 0, 0, 0, -100, 0

 2, 2, 3, 1, 1

 0, 0, 0, -100, 0

 24,24,25, 1, 1 

 0, 0, 0, -100, 0

 25,25,26, 1, 1

 0, 0, 0, -100, 0

*TYPE 1 CONSTRAINTS

 3

 1, 1, 2, 21, 1, 1, 0

 9, 9, 15, 4, 4, 2, 0

 24, 24, 25, 21, 1, 1, 0

34
x

y

uy=0.0

ux=0.0

ux=0.0

Figure 15: Finite element model subjected to compression, transverse tension and shear loading 

xy

yx

xx xx

yy

yy

*STRESS LOADING

 14, 1 

 1, 1, 4, 1, 1

 -1000, 0, -100, 0, 0

 4, 4, 9, 1, 1

 -1000, 0, -100, 0, 0

 9, 9, 15, 1, 1

 -1000, 0, -100, 0, 0

 15, 15, 21, 1, 1

 -1000, 0, -100, 0, 0

 ...

 ... 

 23, 23, 26, 1, 1

 -1000, 0, -100, 0, 0

 1, 1, 2, 1, 1

 0, 20, 0, -100, 0

 2, 2, 3, 1, 1

 0, 20, 0, -100, 0

 24,24,25, 1, 1 

 0, 20, 0, -100, 0

 25,25,26, 1, 1

 0, 20, 0, -100, 0

*TYPE 1 CONSTRAINTS

 3

 1, 1, 2, 21, 1, 1, 0

 9, 9, 15, 4, 4, 2, 0

 24, 24, 25, 21, 1, 1, 0

35
x

y

Figure 16:Finite element mesh of unit cell for small z-pin with 2% areal density 

laminate with

2% carbon z-pins

C = 0.868 mm 

Dz = 0.280 mm

D’z = 0.3 mm

Hz = 1.7526 mm

Lz = 1.7526mm

d = 0.0051 mm

C /d = 170.2 

D’z /d = 58.8

Hz /d = 343.65

Lz /d = 343.65

*GEOMETRICAL POINTS

 26

 1 0.00 0.00

 2 171.85 0.00

 3 343.7 0.00

 4 0.00 142.45

 5 171.85 142.45

 6 343.7 142.45

 7 158.25 144.35

 8 185.75 144.35

 9 0.00 170.85

 10 86.75 170.85

 11 256.95 170.85

 12 343.7 170.85

 13 86.05 171.85

 14 257.65 171.85

 15 0.00 172.85

 16 86.75 172.85

 17 256.95 172.85

 18 343.7 172.85

 19 158.25 199.35

 20 187.75 199.35

 21 0.00 201.25

 22 171.85 201.25

 23 343.7 201.25

 24 0.00 343.7

 25 171.85 343.7

 26 343.7 343.7

*MISALIGNMENT TYPE

 2

*MISALIGNMENT

 74,74,4.7,4.7

36
Figure 17:Finite element mesh of unit cell for small z-pin with 4% areal density 

laminate with

4% carbon z-pins

C = 0.868 mm 

Dz = 0.280 mm

D’z = 0.3 mm

Hz = 1.2446 mm

Lz = 1.2446 mm

d = 0.0051 mm

C / d = 170.2 

D’z / d = 58.8

Hz / d = 244

Lz / d = 244

*MISALIGNMENT TYPE

 2

*MISALIGNMENT

 74,74,3.4,3.4

x

y

*GEOMETRICAL POINTS

 26

 1 0.00 0.00

 2 122.00 0.00

 3 244.00 0.00

 4 0.00 92.60

 5 122.00 92.60

 6 244.00 92.60

 7 108.10 94.50

 8 135.90 94.50

 9 0.00 121.00

 10 36.90 121.00

 11 207.10 121.00

 12 244.00 121.00

 13 36.20 122.00

 14 207.80 122.00

 15 0.00 123.00

 16 36.90 123.00

 17 207.10 123.00

 18 244.00 123.00

 19 108.10 149.50

 20 135.90 149.50

 21 0.00 151.40

 22 122.00 151.40

 23 244.00 151.40

 24 0.00 244.00

 25 122.00 244.00

 26 244.00 244.00

37
Figure 18:Finite element mesh of unit cell for large z-pin with 2% areal density 

laminate with

2% carbon z-pins

C = 2.1844 mm 

Dz = 0.508 mm

D’z = 0.528 mm

Hz = 3.175 mm

Lz = 3.175 mm

d = 0.0051 mm

C / d = 428.31 

D’z / d = 103.53

Hz / d = 622.55

Lz / d = 622.55

*MISALIGNMENT TYPE

 2

*MISALIGNMENT

 74,74,8.6,8.6

x

y

*GEOMETRICAL POINTS

 26

 1 0.00 0.00

 2 311.275 0.00

 3 622.55 0.00

 4 0.00 259.51

 5 311.275 259.51

 6 622.55 259.51

 7 276.275 262.91

 8 346.275 261.91

 9 0.00 310.275

 10 97.12 310.275

 11 525.43 310.275

 12 622.55 310.275

 13 96.12 311.275

 14 526.43 311.275

 15 0.00 312.275

 16 97.12 312.275

 17 525.43 312.275

 18 622.55 312.275

 19 276.275 359.64

 20 346.275 359.64

 21 0.00 363.04

 22 311.275 363.04

 23 622.55 363.04

 24 0.00 622.55

 25 311.275 622.55

 26 622.55 622.55

38
x

y

a. Deformed mesh after 10 load increments

x

y

undeformed

outline

undeformed

outline

Figure 19: Analysis of unit cell of UD Ply with Z-Pin subjected to axial compression

b. Deformed mesh after 1000 load increments

39
undeformed

outline

undeformed

outline

Figure 21: Deformed finite element model of unit cell subjected to

 axial compression and transverse tension loading

Figure 20: Deformed finite element model of unit cell subjected to

axial and transverse compression loading

40
Figure 23: Deformed finite element model of unit cell subjected to

axial compression, transverse tension and shear loading

undeformed

outline

undeformed

outline

Figure 22: Deformed finite element model of unit cell subjected to

axial compression and shear loading

41
Figure 24: Contour plot of fiber misalignment angle φ

Figure 25: Contour plot of shear stress τe

42
Figure 26: Contour plot of shear strain γp

Figure 27: Contour plot of shear strain γs

43
0

200

400

600

800

1000

1200

0 0.2 0.4 0.6 0.8 1 1.2

Figure 29. Stress-displacement plot

displacement ux/ fiber diameter d

stress, 

MPa

maximum

0

100

200

300

400

500

600

700

800

0 5 10 15 20

stress, 

MPa

displacement ux/ fiber diameter d

Figure 28. Stress-displacement plot

maximum

44
REPORT DOCUMENTATION PAGE Form Approved 

OMB No. 0704-0188 

The public reporting burden for this collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing

data sources, gathering and maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or 

any other aspect of this collection of information, including suggestions for reducing this burden, to Department of Defense, Washington Headquarters Services, Directorate 

for Information Operations and Reports (0704-0188), 1215 Jefferson Davis Highway, Suite 1204, Arlington, VA 22202-4302. Respondents should be aware that 

notwithstanding any other provision of law, no person shall be subject to any penalty for failing to comply with a collection of information if it does not display a currently 

valid OMB control number. 

PLEASE DO NOT RETURN YOUR FORM TO THE ABOVE ADDRESS. 

1. REPORT DATE (DD-MM-YYYY) 2. REPORT TYPE 3. DATES COVERED (From - To) 

4. TITLE AND SUBTITLE 5a. CONTRACT NUMBER 

5b. GRANT NUMBER 

5c. PROGRAM ELEMENT NUMBER 

6. AUTHOR(S) 5d. PROJECT NUMBER 

5e. TASK NUMBER 

5f. WORK UNIT NUMBER 

7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES) 8. PERFORMING ORGANIZATION 

REPORT NUMBER 

9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES) 10. SPONSORING/MONITOR'S ACRONYM(S) 

11. SPONSORING/MONITORING

REPORT NUMBER 

12. DISTRIBUTION/AVAILABILITY STATEMENT 

13. SUPPLEMENTARY NOTES 

14. ABSTRACT 

15. SUBJECT TERMS 

16. SECURITY CLASSIFICATION OF: 17. LIMITATION OF 

ABSTRACT 

18. NUMBER 

OF 

PAGES 

19b. NAME OF RESPONSIBLE PERSON 

a. REPORT b. ABSTRACT c. THIS PAGE 

19b. TELEPHONE NUMBER (Include area code) 

Standard Form 298 (Rev. 8-98) Prescribed by ANSI Std. Z39-18