## 1.INTRODUCTION

Composite materials partly behave in a nonlinear fashion, although composite materials generally have been modeled as linear elastic material. The nonlinearity of composite materials can be attributed to inherent material nonlinearity of individual constituents and to micromechanical failures such as fiber or matrix microcracking and interfacial debonding (Yener & Wolcott, 1988; Drucker, 1969). Furthermore, the development of advanced composite material (metal matrix composites) and the complex design of composite structures make it more complicated to model material behavior of composites (Shih & Lee, 1978).

Generally, the literature has assumed linear stress-strain relationships for composite materials.

However, it has also been shown that composite materials partly behave in a nonlinear fashion (Yener & Yi, 1989; Petit & Waddoups, 1969; Hahn & Tsai, 1973). The nonlinearity of composite materials can be attributed to inherent material nonlinearity of the individual constituents.

As mentioned in other paper, the most constitutive models are often used in the practical design field for their simplicity and in spite of their limitations (Sandhus, 1974). Their unavoidable shortcoming for a predictive analysis procedure, however, is that they cannot predict the material behavior of permanent strain accumulation by large deformation. That is the reason why the anisotropic plastic constitutive model is developed in the other paper, which includes the anisotropic yield criterion, a non-proportional hardening rule for changing anisotropic parameters. Hence, highly anisotropic material such as unidirectional composite materials, and their model can be analyzed with account for differential between the tensile and compressive yield strengths (i.e., Bauschinger effect) (Hill, 1948, 1950; Valliappan, 1972).

As mentioned earlier, general response prediction of composite structures becomes possible by developing a realistic and comprehensive analysis procedure for general loading conditions. Such an undertaking, among other considerations, requires very efficient constitutive model which can predict realistically nonlinear material behavior. Hence, an anisotropic plasticity constitutive model for fiber-reinforced composite laminates is implemented into a computer program for a predictive analysis procedure of composites (PACS) (Yener & Yi, 1994).

## 2.A NONLINEAR PREDICTIVE FINITE ELEMENT ANALYSIS PROCEDURE

### 2.1.Characteristics of Predictive Analysis

The proposed predictive nonlinear analysis procedure is presented in here. In a predictive analysis procedure, load is applied to a structure incrementally. At each prescribed step of loading, the desired information describing the behavior of the structure are computed and kept to detect the progressive structural behavior. These behavioral characteristics generally necessitate the utilization of numerical methods. Clearly, because of its lengthy incremental nature, a predictive analysis requires the development of a general purpose computer program.

### 2.2.A Nonlinear Analysis Procedure

In Ref. 15, the finite element governing equations have been formulated, which contain the displacement, strain-displacement matrices and the constitutive matrix of the material. For a formulation to be applicable to general nonlinear problems, however, the constitutive description must be appropriate. The anisotropic elastic-plastic constitutive model, which has been presented in the previous section, is employed. By including the anisotropic elastic-plastic constitutive matrix into the governing equation, the current formulation is applicable to the geometric and material nonlinear problem. The details of finite element formulation and numerical implementation of the total Lagrangian approach is presented in the textbook (Chen & Han, 1988), so is omitted in here.

The nature of the nonlinear equation requires the use of an incremental solution procedure. For each load increment, iterations are generally necessary to satisfy internal equilibrium. When appropriate, the behavior may be assumed to be linear for each iteration within the specified load increment. Successive iterations within each increment are performed to avoid divergence of the solution. The process, here referred to as the modified Newton-Raphson method, is schematically illustrated in Fig. 1. For each load increment, unbalanced load is calculated, which then is used to compute the displacement increment in the proceeding iteration step. This internal iteration process is repeated until the equilibrium is approximated to some acceptable tolerance.

### 2.3.Incremental Elastic-Plastic Constitutive Relationship

This section presents numerical calculations of stress increments. A proper elastic-plastic constitutive relationship integration algorithm is very important in maintaining numerical accuracy, stability, and efficiency for the whole solution procedure (owen & Figueiras, 1983). Several integration techniques for incremental elastic-plastic constitutive relationships are available in the literature (Yener & Yi, 1990, 1992; Tsai & Hahn, 1980). This study uses the stress integration scheme, which is a modified version of the elastic predictor-radial return method (Schreyer, Kulak & Kramer, 1979; Ortiz & Popov, 1985). The elastic predictor-radial return method is first due to Mendelson (1986). During application of an load increment, an element or part of an element can yield. The yield condition is checked to see whether or not plastic deformation has occurred at the Gaussian points. As a result, an element can behave elastically in some parts and plastically in other parts.

If the material behaves plastically for each iteration, the stress should be calculated by an elastic-plastic constitutive relationship. The incremental constitutive relationship for elastic-plastic material is written in matrix form, and then the numerical implementation of the stress computation is followed.

In matrix form, the stress increment {*dσ*} can be expressed in terms of the elastic strain increment {*dϵ ^{e}*}, or the total strain increment {

*dϵ*} as

where [*D ^{e}*] and [

*D*] are the elastic and elastic-plastic stress-strain relationship matrix, respectively. The elastic-plastic constitutive relationship was already given as Eq. 36 in a tensor form, rewritten here in matrix form as

^{ep}Stress calculation is performed at all Gaussian integration points.

Now, the first step of the numerical algorithm is to define a trial stress increment {*Δσ ^{e}*}, assuming an elastic behavior,

Accumulate the trial total stress by adding the assumed elastic stress increment.

In a predictive analysis procedure, the stress state has to be monitored at every step. At this point, three possible conditions can be used to calculate the stress increment. Assuming that the stress state from the previous step was in elastic state, the initial yield condition must be checked. If the trial total stress vector does not violate the initial yield criterion, the elastic material behavior is presumed. In the second case, the trial stress increment has passed the initial yield surface, as shown in Fig. 2. Then, the strain increment is divided into two parts, such as the pure elastic response *R{Δϵ}*, and the elastic-plastic response *(1 - R ){Δϵ}*, where *R* is called a scaling factor. Hence, the stress increment is also divided into two parts that can be integrated as

If, as in the third case, the stress state at the end of *n - 1 ^{th}* load step is already in a plastic state and in a plastic loading state, the stress increment is obtained by integrating the stress increment vector, Eq. 4 with

*R = 0*as

Numerical implementations of integrations appearing in Eqs. 4 and 5 are detailed in Ref. 31.

If the stress is an unloading state, an elastic constitutive relation should be used as

where the stress increment vector is negative corresponding to negative strain increment.

### 2.4.Computational Procedure for Predictive Analysis

Based on the preceding discussion, a predictive finite element analysis procedure, is presented through the flow chart of Yener and Yi (1994). Also, it is implemented in the computer program, PACS which stands for Predictive Analysis of Composite Structures. As shown in the flow chart of Yener and Yi (1994), the essential steps required in a predictive solution algorithm for nonlinear problems are presented to provide a physical insight into the viability of the proposed algorithm.

## 3.NUMERICAL ILLUSTRATIVE EXAMPLES

Anisotropic plasticity and a predictive analysis procedure have been presented for analysis of fiber-reinforced composite structures. Also, their models are successfully implemented into the computer program PACS (Yener & Yi, 1994). In order to verify the validity of the problem formulation and the efficiency of numerical implementation, several example problems are analyzed by using the computer program PACS. Numerical results are compared with various experimental and other numerical results.

### 3.1.Elastic-Plastic Analysis of Composite Laminates

In this section, the anisotropic elastic-plastic analysis capability of the computer program is verified with several composite laminates on different boundary and loading conditions. The elastic-plastic material behavior is observed and compared. Also, the behaviors of isotropic and anisotropic materials are compared.

### 3.2.Clamped Square Plate

A clamped square plate under uniformly distributed loading is considered. Only a quarter of the plate has been discretized because of symmetry considerations (Fig. 3). Eight-noded isoparametric elements with four-point (2x2) Gaussian integration in the plane are used, and eight equal thickness layers are taken through thickness direction. The plate has been solved with both isotropic and anisotropic materials; material properties are given in Table 1.

The thickness and span of plate is h = 0.2 m and L = 6.0 m, respectively.

In Fig. 4, vertical displacements at the plate's center are plotted for different load levels for both isotropic and anisotropic materials. The clear difference between isotropic and anisotropic material is observed. Fig. 4 also includes comparison with other numerical results, and agreement with present results is rather close.

### 3.3.Simple Supported Square Plate

This section analyzes the simple supported plate on four sides of the plate with geometry of plate and finite element discretization given in Fig. 5. The purpose of this problem is to see effects of different boundary conditions on the plastic behavior of the plate. In order to compare with the clamped plate, the same geometry as the clamped plate is used with a uniformly loaded case. The comparison is performed only on isotropic material. The development of plastic regions for both clamped and simple supported plates are presented in Fig. 6.

Yielding starts at the corners in the simple supported plate, while starting in the middle of the sides of the clamped plate. Next, the center of both plates yields. Then the plastic regions propagate from the vicinity of the corners and the plate centers. The pattern of the plastic region growth agrees with that of Dodds (1987).

### 3.4.Clamped Quadratic Square Shell

An elastic-plastic analysis of a clamped spherical shell under self-weight condition is undertaken. The geometry of the spherical shell is given in Fig. 7. The coordinates of z-direction are obtained as

where L is the width of shell and c is equal to L/10. Material characteristics are the same as employed in Table 1 of the plate problem. Finite element discretization is same as that of the plate problem in Fig. 5, where nine-nodes Heterosis element is used.

Displacements at center of shell for both isotropic and anisotropic materials are plotted in Fig. 8, which show numerical results from Owen and Figueiras (1983) and good agreements.

The spherical shell is again solved under a concentrated load at center. Vertical displacement at center is plotted against the applied load in Fig. 9.

## 4.CONCLUSIONS

The development and use of an analytical procedure which has capability of predicting the progressive material behavior of structures, named as a predictive analysis procedure in here, is developed for more accurate assessment of structural safety and efficiency of composite structures.A quadratic anisotropic yield criterion in stress space is developed for general use with unidirectional and bidirectional composite lamina. The developed anisotropic work-hardening model allows for a nonproportional change of the yield values so that the subsequent yield surface can be a distorted shape. The predictive finite element analysis procedure is developed for accurate and efficient analysis in predicting progressive behavior of composite structures. As a result, a computer code, PACS (Predictive Analysis of Composite Structures) is developed, which adopts the abovementioned anisotropic plasticity model.

The accuracy and efficiency of the computer code PACS are verified with various benchmark problems. Numerical predictions of the computer program PACS compare very well with available experimental, analytical and other numerical results. Comparisons illustrate the capability of the constitutive model and the computer program PACS. The developed constitutive and predictive analysis model predicts progressive nonlinear behavior from the beginning of loading, in plane and through thickness direction of composite laminates. The capability of PACS to predict nonlinear material behavior of composites can be used as a helpful device in parametric studies for design of fiber-reinforced laminates.