Scalable Coupling of Multiscale AEH and PARADYN 

Analyses for Impact Modeling 

by Rama R. Valisetty, Peter W. Chung, and Raju R. Namburu 

ARL-TR-3512 June 2005 

Approved for public release; distribution is unlimited. 
NOTICES 

Disclaimers 

The findings in this report are not to be construed as an official Department of the Army position unless 

so designated by other authorized documents. 

Citation of manufacturer’s or trade names does not constitute an official endorsement or approval of the 

use thereof. 

Destroy this report when it is no longer needed. Do not return it to the originator. 
Army Research Laboratory 

Aberdeen Proving Ground, MD 21005-5067 

ARL-TR-3512 June 2005 

Scalable Coupling of Multiscale AEH and PARADYN 

Analyses for Impact Modeling 

Rama R. Valisetty, Peter W. Chung, and Raju R. Namburu 

Computational and Information Sciences Directorate, ARL 

Approved for public release; distribution is unlimited. 
ii

Report Documentation Page Form Approved 

OMB No. 0704-0188

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 information. Send comments regarding this burden estimate or any other aspect of this collection of 

information, including suggestions for reducing the 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)

June 2005 

2. REPORT TYPE 

Final 

3. DATES COVERED (From - To) 

2003–2004 

5a. CONTRACT NUMBER 

5b. GRANT NUMBER 

4. TITLE AND SUBTITLE 

Scalable Coupling of Multiscale AEH and PARADYN Analyses for Impact 

Modeling 

5c. PROGRAM ELEMENT NUMBER 

5d. PROJECT NUMBER 

5UH7FL 

5e. TASK NUMBER 

6. AUTHOR(S) 

Rama R. Valisetty, Peter W. Chung, and Raju R. Namburu 

5f. WORK UNIT NUMBER 

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

U.S. Army Research Laboratory 

ATTN: AMSRD-ARL-CI-HC 

Aberdeen Proving Ground, MD 21005-5067 

8. PERFORMING ORGANIZATION 

 REPORT NUMBER 

ARL-TR-3512 

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

11. SPONSOR/MONITOR'S REPORT 

 NUMBER(S) 

12. DISTRIBUTION/AVAILABILITY STATEMENT 

Approved for public release; distribution is unlimited. 

13. SUPPLEMENTARY NOTES 

14. ABSTRACT 

This report describes scalable coupling of two stand-alone computer codes for a multiscale impact analysis of composites. An asymptotic expansion 

homogenization (AEH)-based microstructural model available for modeling microstructural aspects of modern armor materials is coupled with 

PARADYN, a parallel explicit Lagrangian finite-element code. The first code enables modeling of a material microstructure and provides material 

response in terms of global structural response at a material integration point. Microstructural codes such as this one are typically used in stand-alone 

form and applied in simple loading situations. The coupling provided here enables a micro/macro type multiscale analysis under generalized 

three-dimensional loading conditions. Three sets of results are presented to demonstrate: (1) the verification of the AEH-PARADYN model 

coupling to PARADYN, (2) the scalability of the coupled model, and (3) an application for modeling the impact response of armor materials. 

In the present work, a consistent AEH numerical formulation was selected that was earlier shown to be scalable on different computing architectures. 

In conjunction with a second-order accurate velocity-based explicit time integration method, the formulation is coupled within the message-passing 

interface scheme of PARADYN. PARADYN is a scalable version of the Lawrence Livermore National Laboratory’s serial DYNA3D explicit 

Lagrangian finite-element code used for obtaining large deformation, elastic/plastic response of structures.

15. SUBJECT TERMS

AEH, PARADYN, HPC, scalable computing, impact modeling 

16. SECURITY CLASSIFICATION OF: 

19a. NAME OF RESPONSIBLE PERSON

Rama R. Valisetty 

a. REPORT 

UNCLASSIFIED 

b. ABSTRACT 

UNCLASSIFIED

c. THIS PAGE 

UNCLASSIFIED

17. LIMITATION 

OF ABSTRACT 

UL 

18. NUMBER 

OF PAGES 

44 

19b. TELEPHONE NUMBER (Include area code) 

410-278-6638 

Standard Form 298 (Rev. 8/98) 

Prescribed by ANSI Std. Z39.18
iii

Contents 

List of Figures iv

List of Tables v

1. Introduction 1

2. AEH Development 3

3. Coupling the Microstructural AEH Model in PARADYN 4

3.1 Guidelines for Using the Model ......................................................................................4 

3.2 Details of the Data Exchange ..........................................................................................5 

4. Results 6

4.1 Verification of the Model Insertion.................................................................................6 

4.2 Validation of the Model ..................................................................................................7 

4.3 Scalability of the Results...............................................................................................11

4.4 An Application of the Model.........................................................................................13 

5. Future Work and Conclusions 25

6. References 26

Appendix. Microstructural Equations 28

Distribution List 35
iv

List of Figures 

Figure 1. Finite-element mesh of a cylindrical wedge under Taylor impact..................................7 

Figure 2. DYNA3D and AEH-PARDYN axial velocity solutions for the finite-element mesh 

of figure 1...................................................................................................................................8 

Figure 3. DYNA3D and AEH-PARDYN axial stress solutions for the finite-element mesh of 

figure 1.......................................................................................................................................9 

Figure 4. DYNA3D and AEH-PARDYN effective plastic strain solutions for the finiteelement mesh of figure 1..........................................................................................................10 

Figure 5. Finite-element model and deformed mesh used for the AEH-PARDYN solution for 

the verification problem...........................................................................................................12 

Figure 6. Finite-element model and deformed mesh used for the DYNA3D solution for the 

verification problem.................................................................................................................13 

Figure 7. DYNA3D and AEH-PARDYN axial velocity solutions for the validation problem....14 

Figure 8. DYNA3D and AEH-PARDYN axial stress solutions for the validation problem........15 

Figure 9. DYNA3D and AEH-PARDYN effective stress solutions for the validation 

problem. ...................................................................................................................................16 

Figure 10. DYNA3D and AEH-PARDYN effective plastic strain solutions for the validation 

problem ....................................................................................................................................17 

Figure 11. Scalability results for a fixed mesh and for a constant processor load........................18 

Figure 12. Finite-element mesh for PC74AL37 laminate and its buffer layer, back plate, and 

flyer..........................................................................................................................................19 

Figure 13. Finite-element mesh for PC74SS37 laminate and its buffer layer, back plate, and 

flyer..........................................................................................................................................19 

Figure 14. Finite-element mesh for PC74GS55 laminate and its buffer layer, back plate, and 

flyer..........................................................................................................................................20 

Figure 15. Experimental result and DYNA3D prediction for axial velocity at the back of 

buffer layer for the PC74AL37 laminate. ................................................................................21 

Figure 16. Experimental result and DYNA3D prediction for axial velocity at the back of 

buffer layer for the PC74SS37 laminate. .................................................................................22 

Figure 17. Experimental result and DYNA3D prediction for axial velocity at the back of 

buffer layer for the PC74GS55 laminate. ................................................................................22 

Figure 18. Experimental result and AEH-PARDYN prediction for axial velocity at the back 

of buffer layer for the PC74AL37 laminate.............................................................................23 

Figure 19. Experimental result and AEH-PARDYN prediction for axial velocity at the back 

of buffer layer for the PC74SS37 laminate..............................................................................24 

Figure 20. Experimental result and AEH-PARDYN prediction for axial velocity at the back 

of buffer layer for the PC74GS55 laminate. ............................................................................24 

Figure A-1. Macro/microanalyses neighborhoods........................................................................28 
v

List of Tables 

Table 1. Material properties............................................................................................................6 

Table 2. Material properties............................................................................................................7 

Table 3. Wall clock times for scalability study.............................................................................18

Table 4. Laminate and impact velocity details. ............................................................................20 

Table 5. Material properties..........................................................................................................21 
vi

INTENTIONALLY LEFT BLANK. 
1

1. Introduction 

Heterogeneous materials, such as layered or dispersed particulate composites, are widely used in 

military armor applications for protection from ballistic impacts or for defeating enemy 

structures and vehicles. For effective use of these materials in short duration impact 

applications, one should be able to quantify the local- or microlevel structural response of the 

constituent phases, understand the different failure modes at the microstructural level, and 

translate this information into the effective constitutive material properties at the global level. 

Spanning a heterogeneous material from its component microphases to its global structural level 

with finite elements to obtain micro- and global structural responses leads to large finite-element 

meshes and computations that are both expensive and intractable. Multiscale analyses are used 

to alleviate this problem and are even being extended to include nanoscale analyses to aid 

nanomaterial development. 

At a basic level, these analyses use one or more microstructural models for describing 

microlevel material features, such as fibers, particles, lamina, voids, etc., and express the 

microlevel response as a local perturbation of global response. Each element, or a selected few 

of the global elements, calls the microstructural model to sample the microlevel structural 

response and build the element level stresses. 

While many theories are available for linking the microstructural response to the global response, 

the theory selected here is capable of modeling elastic/plastic large deformation seen in the 

armor impact applications. It is based on the mathematically rigorous asymptotic expansion 

homogenization (AEH) method. While a historical context to it is presented in the next section, 

its most recent proponents are Chung et al. (1999), who extended its range of application to 

include short transient loading and applied it for predicting the elastic/plastic large deformation 

of composite materials. 

In the AEH approach, the microlevel response is expressed in terms of the global response in a 

strict, mathematically seamless approach. An updated Lagrangian scheme for small strains and 

small rotations is employed to account for large displacements, strains, and rotations over many 

time steps. At each time step, the equilibrium solution from the previous time step is used to 

update the local microstructure. 

Since explicit time integrations take quite a large number of time steps and microlevel 

computations are to be performed at each and every time step and for the selected global finite 

elements, the attendant micro computations can take much longer times for solution. Scalable 

computational approaches are needed to use these computationally intensive methods. With the 

aid of domain decomposition and message-passing interface (MPI), the basic equations of 
2

Chung et al. (1999) were recast for parallel computation and were shown to be scalable for a 

Taylor impact of a laminate with several hard and soft layers by Valisetty (2000). 

While the previous two references used the microstructural AEH as a stand-alone code, the 

current work couples it as a material level response provider to PARADYN, a general purpose 

explicit dynamics code. The objective was to extend the applicability of AEH to large-scale 

multi-material impact analyses. PARADYN is a parallel version of the Lawrence Livermore 

National Laboratory’s serial DYNA3D code. It is a natural choice for the present purpose 

because of its large number of material models, contact algorithms, element types, etc. 

Coupling two stand-alone structural codes, such as the microstructural AEH and PARADYN, is 

best accomplished when the interface computer code is written in a form that is as simple as 

possible. While an independent stand-alone interface is desirable, it has to accomplish two 

objectives: (1) the data exchange between the two codes must be transparent to users, and (2) the 

coupling to the driver code, i.e., PARADYN, must be done in such a way that it does not disturb 

any of the driver code’s special features, such as scalable parallel operations, time integration, 

contact algorithms, etc. 

Since PARADYN is a scalable parallel code, it already has its computational logic and data 

structures set for tackling parallel computations. Even though the code has facilities to add 

user-defined materials and elements, the present microstructural AEH is neither a new material 

nor an element. Since potentially large microstructural element models can be considered with 

the microstructural AEH, requiring large amounts of microelement specific material properties, 

stresses and strains data to be stored and made available for conducting time integration, using 

PARADYN’s data arrays for storing and retrieving such microelement data becomes intractable 

and potentially conflicting to its scalability operations. To alleviate this situation, external 

processor specific files are used for writing a global element’s data and microstructural elements’ 

data. This enables the AEH model coupling without disturbing the PARADYN code’s MPI and 

domain decomposition-based scalability operations. 

To enable access to the microstructure AEH model, a reference is provided to it as a “Material 

Type 3” sub-option in PARADYN. This reference is completely arbitrary and could be placed 

even without getting into any material level subroutines. Global elements tabbed with this 

material option then write the current time step and strain increment data to external files and 

read back the effective mechanical properties and effective stresses. With the aid of a 

user-specified microelement finite-element model and after reading global elements’ data, the 

AEH model conducts microlevel analyses and writes the effective properties and stress data to 

files. The data exchange between the PARADYN and AEH codes occurs through external 

processor specific files. Thus, with this interface, a dual level finite-element modeling capability 

is added to the PARADYN code. 

With regards to writing the computer code of the interface, while many environments are 

available, the one that was considered is the “interdisciplinary computing environment (ICE).” 
3

As described by Clarke and Namburu (2002), this environment uses an “extensible data model 

and format (XDMF)” for describing the nature of data being exchanged between the codes and a 

visual environment for the user to monitor the exchanging data while the codes run 

independently in the background. The xdmf uses extensible markup language (XML), 

hierarchical data format version 5 (HDF5), and network distributed global memory (NDGM). 

The XML language uses PYTHON language script to enable users to describe the types of data 

being exchanged between the independently running codes. The HDF5 data format in which the 

codes exchange data via external files is an effective format not only for doing read/writes but 

also for grouping the data for fast scientific visualization. The NDGM provides data access to 

codes running on processor grids in a scalable manner. 

While the ICE experience of data markup, visualization, storage, and handling will be discussed 

in a future publication, the usefulness of the AEH-PARADYN coupling is presented in the 

current report. When coupled with PARADYN, the microstructural AEH can potentially avail 

the various material models for its constituent materials from PARADYN. Also, in application 

runs, PARADYN’s features (e.g., contact and the ability to remove excessively deformed global 

materials) can also be used. The microstructural AEH enhancement of PARADYN was 

demonstrated through simulation of plate impact tests for three composite laminates made of 

extremely dissimilar materials. 

The outline of this report is as follows. Section 2 presents a brief historical context to the AEH 

development. Section 3 presents the details of the microstructural AEH coupling to PARADYN 

and guidelines for its use. Section 4 presents verification, scalability, and application of the 

coupled codes. The micro- and global structural response coupling and equations are presented 

in the appendix. 

2. AEH Development 

While the rule of mixtures represents an earliest concept in composite micromechanics, the idea 

that stress/strain distributions caused by holes and inclusions in an otherwise uniform medium 

can be used to represent the stress/strain distributions around fibers in laminas is the starting 

point for the development of homogenization theories. Licht Frottement (1987), Lene (1984), 

Moulinec and Suquet (1994), Bensoussan et al. (1978), and Sanchez-Palencia (1980) have 

formalized the mathematics of the homogenization theories for obtaining the effective material 

properties. The works of Hashin (1962), Hashin and Strickman (1963), Hill (1964), Willis 

(1982), and Suquet (1982) can be seen as the engineering counterparts of these theories. 

For a short transient impact loading condition, a detailed knowledge of the material flow in and 

around material microstructural constituents and material interfaces is necessary. This flow is 

three-dimensional and elastic/plastic. Classical rule of mixtures theories based on constant stress
4

or strain assumptions lead to an aggregation of the response in the microstructural details and 

thus are of limited use. 

Recent advances in mathematically rigorous AEH methods enable the coupling of micro- and 

global approaches for both linear and nonlinear structural applications (Fish et al., 1997; Terada 

and Kikuchi, 1995; Ghosh and Moorthy, 1995; Lene, 1986; and Guedes and Nikuchi, 1990). 

Most recently, Chung et al. (1999) demonstrated the applicability of the AEH method for 

heterogeneous media subjected to short transient loading by employing explicit dynamics finiteelement formulations in conjunction with elastic/plastic material response. The investigation of 

Chung et al. (1999) considers the dynamic equation of motion discretized according to a second 

order accurate scheme presented in Namburu (1991). 

Using the AEH approach, the microlevel variables are expressed as direct functions of the global 

variables in a strict, mathematically seamless approach. This method is especially suitable for 

conducting explicit transient elastic/plastic analysis of heterogeneous materials. An updated 

Lagrangian scheme for small strains and small rotations is employed to account for large 

displacements, strains, and rotations over many time steps. 

Because an explicit time integration scheme takes quite a large number of time steps and the 

global and microlevel computations are to be performed at each time step and for each global 

finite element, the attendant micro and global computations can take much longer times for 

solution. Scalable computational approaches are needed to explore these computationally 

intensive methods. In Chung et al. (1999), the microstructural AEH was made parallel in a 

stand-alone explicit dynamics code with results showing good scalability for the AEH’s coupled 

micro/macro equations. These equations are presented in the appendix. They form the basis of 

the microstructural AEH that is now coupled with PARADYN. 

3. Coupling the Microstructural AEH Model in PARADYN 

PARADYN is a general-purpose explicit dynamics code with a capability for defining multiple 

element types and material types. Each element of the global finite-element model that is tabbed 

for detailed microlevel analysis is made to call the microstructural AEH model to sample the 

local microlevel behavior and build the element global stresses from the constituent 

microstresses. 

3.1 Guidelines for Using the Model 

The microstructural AEH sub-option can be invoked by specifying a “Material Type 3” and a 

non-zero number in the fifth data location on the eighth material card. The material option 
5

describes a nonlinear elastic/plastic response with kinematic and/or isotropic hardening. The 

added microcode is not limited to microphase materials of this type, however. 

The microstructure is next defined with the usual nodes, elements, and material properties. 

Although any type of format can be used for this purpose, the one that is selected uses the same 

PARADYN input format. This was done to avail PARADYN’s preprocessors, such as INGRID. 

Although not verified, at present there are no restrictions on the number of microstructural 

elements, nodes, and materials. 

All the information defining the microstructural AEH is to be contained in an input file separate 

from the PARADYN’s global structure-defining input file and is to be made available to the 

PARADYN executable along with a domain decomposition part file. 

3.2 Details of the Data Exchange 

While PARADYN focuses on global computations such as the aggregation of nodal forces, 

solution of the equations of motion, updating global nodal positions, enforcing contact, etc., the 

microstructural AEH code focuses on the microstructural AEH calculations which typically 

involve evaluating the microstructure as a perturbation of global deformation and periodicity. 

For the two codes to perform their respective computations, it is necessary that the relevant data 

be exchanged between them for each time step and for each global element identified for 

microcomputations. In the present work, the data exchange involves two parts as follows: 

1. From PARADYN, at the beginning of each time step and for the identified global elements, 

the following data is written to external processor-specific data files for read access to the 

microstructural AEH model: 

Current time step, 

current global strain increments, 

global strains from previous time step, 

global stresses from previous time step, and 

effective plastic strain from previous time step. 

2. After reading the PARADYN-written data files, the microstructural AEH code computes 

and writes the following data to external processor-specific data files for read access to 

PARADYN: 

Effective density, 

effective bulk modulus, 

effective shear modulus,
6

computed current global stresses for the element, and 

computed current effective plastic strain for the element. 

While PARADYN’s native arrays are used for tracking global stresses, global strains, effective 

properties (density, bulk modulus, and shear strain), and effective plastic strain, no such native 

arrays are used to store the microelement quantities. Instead, the microelement data are written 

to, and read from, processor-specific files during each time step and updated. As opposed to the 

global stresses of the global elements, the stored values of these responses reflect the accurate 

stress distribution within the microelements of the microstructure. During a PARADYN run, or 

after the run is complete, the processor-specific files can be visualized independently to analyze 

the microstructural stress and deformation development for any of the global elements. 

4. Results 

Verification of the coupling algorithm, evaluation of the coupling scalability, and validation of 

the coupled codes were achieved through the following numerical analyses conducted on the 

U.S. Army Research Laboratory’s SGI O2k machines. 

4.1 Verification of the Model Insertion 

For verifying the model coupling, a small cylinder with homogeneous material was considered in 

a Taylor impact, and was analyzed alternatively as a heterogeneous material with DYNA3D and 

as a homogeneous material with alternating layers of two different materials but with same 

identical properties using the microstructural AEH in PARADYN. The idea was that the AEH 

model should default into a homogeneous material solution, and the data exchange can point to 

any errors in the coupling. 

Enforcing quarter symmetry on sides parallel to the axis of the cylinder, a 90° wedge of the 

cylinder was considered with 108 elements in five rows, as shown in figure 1. The length and 

radius of the cylinder was 2.592 mm and 3.2 mm, respectively. The material was considered 

elastic/plastic with isotropic hardening. The material properties are shown in table 1. 

Table 1. Material properties. 

Density 

(kg/m3

) 

Young’s Modulus 

(GPa) 

Poisson’s Ratio Yield Stress 

(GPa) 

Tangent Modulus 

(GPa) 

8930 117 0.35 0.4 0.1 

The velocity of the cylinder was 227 m/s. To affect the deformation, one end of the cylinder was 

assumed to be free of velocity. Results for axial velocity, axial stress, and effective plastic strain 

are considered for three nodal locations, near the free end, center, and at the impacted end, as 

shown in figure 1. They are presented in figures 2–4, respectively, for the heterogeneous 
7

Figure 1. Finite-element mesh of a cylindrical 

wedge under Taylor impact. 

DYNA3D solution and for the homogeneous microstructural AEH-PARADYN solution. Both 

solutions appear to be in agreement with each other, thus verifying the model coupling for 

microstructures with homogeneous materials. 

4.2 Validation of the Model 

In this example, a cylinder similar to the one in the previous example was considered but with 

two materials alternating in a total of 96 layers normal to the cylinder axis. The length and 

radius of the cylinder were 31.1 mm and 3.2 mm, respectively. The material was considered 

elastic/plastic with isotropic hardening. The material properties are shown in table 2. 

Table 2. Material properties. 

Layer Density 

(kg/m3

) 

Young’s Modulus 

(GPa) 

Poisson’s Ratio Yield Stress 

(GPa) 

Tangent Modulus

(GPa) 

1 8930 117e9 0.35 0.4e9 0.1e9 

2 8930 90e9 0.35 0.4e9 0.1e9 

Two solutions were obtained. The first solution was a heterogeneous DYNA3D solution 

obtained using 10,368 global elements. To model the repeating 48 sets of the two dissimilar 

materials, eight rows of elements were used for each material layer. The second solution was the 

homogeneous microstructural AEH-PARADYN solution obtained with a 1296 element global 

model. The elements are stacked in 96 layers. In contrast to the first solution, a single layer of 

elements was used for each material. All layers are assumed to be exhibiting homogeneous 

composite behavior in the global sense, but with microstructural AEH computed using a 

two-layer microstructural model. 
8

Taylor Test

-400

-200

0

200

400

0 2 4 6

Time (μs)

Axial Velocity (m/s) at Node 

177

AEH_PARADYN

DYNA3D

Taylor Test

-400

-200

0

200

400

0 2 4

Time ( 6

μs)

Axial Velocity (m/s) at Node 

168 

AEH_PARADYN

DYNA3D

Taylor Test

-400

-200

0

200

400

0 2 4 6

Axial Velocity (m/s) at Node 

Time (μs)

159

AEH_PARADYN

DYNA3D

Figure 2. DYNA3D and AEH-PARDYN axial velocity solutions for the finite-element 

mesh of figure 1. 
9

Taylor Test

-5

-3

-1

1

3

0246 Time (μs)

Axial Stress (Gpa) at Node 

177

AEH-PARADYN

DYNA3D

Taylor Test

-7

-5

-3

-1

1

3

0246 Time (μs)

Axial Stress (Gpa) at Node 168

AEH-PARADYN

DYNA3D

Taylor Test

-9

-7

-5

-3

-1

1

3

5

0246 Time (μs)

Axial Stress (Gpa) at Node 

159

AEHPARADYN

DYNA3D

Figure 3. DYNA3D and AEH-PARDYN axial stress solutions for the finite-element 

mesh of figure 1. 
10

Taylor Test

0

0.02

0.04

0.06

0.08

0.1

0246 Time (μs)

Effective Plastic Strain at 

Node 168

PARADYN-AEH

DYNA3D

Taylor Test

0

0.04

0.08

0.12

0.16

0246 Time (μs)

Effective Plastic Strain at 

Node 177

PARADYN-AEH

DYNA3D

Taylor Test

0

0.02

0.04

0.06

0.08

0246 Time, μs

Effective Plstic Strain at Node 

159

AEH-PARADYN

DYNA3D

Figure 4. DYNA3D and AEH-PARDYN effective plastic strain solutions for the finiteelement mesh of figure 1.
11

Enforcing quarter symmetry on sides parallel to the axis of the cylinder, a 90° wedge of the 

cylinder was considered. The velocity imposed on all but one end face of the cylinder is 

227 m/s. Figures 5 and 6 show the two models both before and after the impact. Solutions at 

nodes near top, center, and bottom of the cylinder are presented in figures 7–10, for axial 

velocity, axial stress, effective stress, and effective plastic strain. For all these responses, the 

microstructural AEH was able to track the DYNA3D’s heterogeneous responses well. 

Deviations in the axial velocity and axial stress can be attributed to three factors: (1) these 

differences are near stress-free edges prone to be affected by the transients; (2) the meshes are 

different in their coarseness which affects the time integration in the two solutions; and (3) the 

values are node averaged but the nodes are at different locations in the two models since the 

meshes are different. Similarly, the differences in the plastic strain predictions can be attributed 

to the extreme length of the time duration of the analysis, which makes the differences 

accumulate. 

4.3 Scalability of the Results 

During the development of the AEH-PARADYN coupling interface, due consideration was 

given to the fact that PARDYN is a parallel code while coupling the microstructural AEH. 

Because the PARADYN code has a well tested MPI/OPENMP-based computational logic, the 

coupling was done by keeping a global element’s microstructural AEH computations local to the 

processor on which the global element resides. This means that information such as 

microelement stresses, strains, and mechanical properties needed in microelement computations 

are kept local and not passed into global data structures. Consequently, a global element’s 

microstress, microstrain, and micromechanical properties are needed to be written to, and read 

from, processor-specific files in between time steps. To demonstrate the fact that the penalty 

associated with this approach is small, results from a limited scalability study are presented in 

this section. 

The Taylor impact model of section 4.2 is used for this purpose. Two types of scalability 

studies, one with a fixed mesh for standard scalability, and the other with varying mesh size, but 

with a fixed number of elements per processor, are performed. The nodes, elements, and wall 

clock times from these studies are presented in table 3. In these simulations, the increase in 

number of elements was achieved by adding element rows in the cylindrical axis direction only, 

and not in the hoop and radial directions. 

The wall clock times are also plotted in figure 11. The plot on the left side of this figure shows 

run times for the 122,880 element mesh for the standard scalability study. A linear speedup can 

be discerned from this plot. 

On the other hand, the plot on the right side of the figure shows the run times for the scaled 

scalability study. In this study, the mesh size was doubled with each doubling of the processors. 

The run times did not remain constant as were the number of elements per the processor. The 
12

Figure 5. Finite-element model and deformed mesh used for the AEH-PARDYN solution for the 

verification problem. 

reason could be that the plastic deformation is very concentrated near base and the domain 

decomposition used to distribute the elements over processors is not sensitive to this fact. 

Another factor could be the small processor load (~2000 elements/processor). The scalability 

might improve with a higher processor load (e.g., 20,000 elements/processor or higher). 
13

Figure 6. Finite-element model and deformed mesh used for the DYNA3D solution for the verification 

problem. 

4.4 An Application of the Model 

In this section, the ability of the microstructural model AEH-PARADYN to predict the 

experimentally observed impact behavior of laminates with widely dissimilar and alternating 

layers of materials is demonstrated. For this purpose, impacts of three composite laminates from 

the experimental work of Zhuang (2002) are considered. In all these impacts, the objective was 

to study the effect of composite lay-ups with alternating hard and soft layers on shock wave 

propagation. Capturing the local mechanical behavior of such lay-ups is very difficult, 

especially in the regime of shock. 
14

-25

25

75

125

175

225

275

0 20 40 60 80 Time(μs)

Axial Velocity (m/s) at top 

center node 

AEH-PARADYN

DYNA3D

-25

25

75

125

175

225

275

0 20 40 60 80 Time(μs)

Axial Velocity (m/s) at middle 

center node 

AEH-PARADYN

DYNA3D

-25

25

75

125

175

225

275

0 20 40 60 80 Time(μs)

Axial Velocity (m/s) at bottom 

center node 

AEH-PARADYN

DYNA3D, near

Figure 7. DYNA3D and AEH-PARDYN axial velocity solutions for the validation 

problem.
15

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0 20 40 60 80 Time(μs)

Axial Compressive Stress(Gpa) at 

top Center

AEH-PARADYN

DYNA3D

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

0 20 40 60 80 Time(μs)

Axial Compressive Stress(Gpa) at 

middle Center

AEH-PARADYN

DYNA3D

-8

-6

-4

-2

0

2

4

6

0 20 40 60 80 Time(μs)

Axial Compressive Stress(Gpa) at 

bottom center

AEH-PARADYN

DYNA3D

Figure 8. DYNA3D and AEH-PARDYN axial stress solutions for the validation problem. 
16

0

0.2

0.4

0 20 40 60 80

Time(μs)

Effective Stress(Gpa) at top 

center

AEH-PARADYN

DYNA3D

0

0.2

0.4

0 20 40 60 80 Time(μs)

Effective Stress(Gpa) at middle 

center

AEH-PARADYN

DYNA3D

0

0.2

0.4

0.6

0 20 40 60 80 Time(μs)

Effective Stress(Gpa) at bottom 

center

AEH-PARADYN

DYNA3D

Figure 9. DYNA3D and AEH-PARDYN effective stress solutions for the validation 

problem.
17

0

0.0025

0.005

0.0075

0.01

0 20 40 60 80 100

Time(μs)

Plastic Strain at top center 

node

AEH-PARADYN

DYNA3D

0

0.1

0.2

0.3

0.4

0.5

0 20 40 60 80 100 Time(μs)

Plastic Strain at middle center node

AEH-PARADYN

DYNA3D

0

0.25

0.5

0.75

1

1.25

1.5

1.75

2

2.25

2.5

2.75

3

0 20 40 60 80 100 Time(μs)

Plastic Strain at bottom center node

AEH-PARADYN

DYNA3D

Figure 10. DYNA3D and AEH-PARDYN effective plastic strain solutions for the 

validation problem
18

Table 3. Wall clock times for scalability study. 

Nodes Elements No. of Processors Wall Clock Time 

9787 8064 4 3072 

18758 16128 8 4051 

35321 31240 16 6046 

71217 64512 32 7492 

137538 127050 64 8944 

140727 122880 8 88369 

140727 122880 16 44376 

140727 122880 32 22303 

140727 122880 64 11353 

0

8

16

24

32

40

48

56

64

72

80

0 8 16 24 32 40 48 56 64 72 80

No. of Processors

Scalability Factor

Actual

Scalability

Ideal

Scalability

Scalability for constant processor load

0

0.5

1

1.5

2

2.5

3

3.5

0 8 16 24 32 40 48 56 64 72

No. of Processors

Scalability Factor

Ideal Scalability

Actual Scalability

Figure 11. Scalability results for a fixed mesh and for a constant processor load. 

In these laminates, a basic unit consisting of a soft layer and a hard layer is repeated a number of 

times. Laminates are then bonded to soft buffer layers and mounted on to back plates. The soft 

layer is a polycarbonate (PC) plastic sheet of 0.74-mm thickness. Different materials and 

thicknesses are used for the hard layers as follows: aluminum (0.37-mm thickness), steel 

(0.37-mm thickness), and glass (0.55-mm thickness). The soft buffer layer has a thickness of 

0.74 mm and is made from PC. The back plate has a thickness of 12.5 mm and is made from 

polymethyl methacrylate plastic. As shown in figures 12–14, the laminates are impacted in a 

direction normal to the lay-up using flyer plates made of PC. The flyer plate has a thickness of 

2.87 mm. The diameters of the flyer plates, laminates, and the back plates are 34, 38.1, and 

38.1-mm, respectively. The details of the three laminate constructions are in table 4. 
19

Figure 12. Finite-element mesh for PC74AL37 

laminate and its buffer layer, back plate, 

and flyer. 

Figure 13. Finite-element mesh for PC74SS37 

laminate and its buffer layer, back 

plate, and flyer. 
20

Figure 14. Finite-element mesh for PC74GS55 

laminate and its buffer layer, back 

plate, and flyer. 

Table 4. Laminate and impact velocity details. 

Laminate Types 

No. of Repeating 

Units 

Laminate Thickness 

Including Buffer Layer

(mm) 

Flyer Plate Velocity 

(m/s) 

PC74AL37 5 6.62 589 

PC74SS37 8 9.97 561 

PC74GS55 7 9.95 563 

As can be seen from the deformed meshes in figures 12–14, different meshes are used for the 

three laminates. The PC74AL37 laminate, including its buffer layer, back plate, and the flyer 

plate, were all meshed with a total of 69899 nodes and 55095 elements. The PC74SS37 laminate 

and buffer, back, and flyer plates were meshed with a total of 108326 nodes and 82827 elements. 

Finally, the PC74GS55 laminate’s impact was modeled with a total of 167549 nodes and 141798 

elements. Because the experiments showed much delamination, and because the interlaminar 

strengths of the laminates are relatively small, i.e., when compared to the impact generated loads, 

laminar interfaces were used and modeled as zero strength slide surfaces. The effect of these 

slide surfaces can be seen as outflowing layers in the deformed meshes in figures 12–14. The 

material properties used for the present simulations are in table 5. 
21

Table 5. Material properties. 

Before conducting the dual level macro/micro analyses with AEH-PARADYN, single level 

DYNA3D global analyses were conducted using the aforementioned finite-element models. 

Using DYNA3D’s elastic/plastic material Model 3, for layer materials, the laminate is treated as 

materially heterogeneous. The results for the axial velocity at the back of the buffer layer are 

compared with the experimental results in figures 15–17. Even though the overall trends of 

loading and unloading of the shock wave are predicted well, there is much discrepancy in the 

predictions of the initial slope and the peak velocity reached. 

0

50

100

150

200

250

300

350

0.E+00 3.E-06 6.E-06

Time, s

Axial velocity, m/s

DYNA3D

Experiment

Figure 15. Experimental result and DYNA3D 

prediction for axial velocity at the back 

of buffer layer for the PC74AL37 

laminate. 

Material Density 

(g/m3) 

Shear Modulus 

(GPa) 

Yield Stress 

(GPa) 

Plastic Modulus 

(GPa) 

Poisson’s Ratio 

PC 1190 0.94 0 1.60 0.37 

PMMA 1180 1.20 0 1.60 0.34 

6061 AL 2710 30.0 0.32 0.69 0.33 

304 SS 7890 77.0 0.33 1.7 0.29 

D 263 glass 2510 30.1 — — 0.208 

Time (s) 

Axial velocity (m/s) 
22

0

50

100

150

200

250

300

0E+00 5E-06 1E-05 Time, s

Axial Velocity, m/s

DYNA3D

Experiment

Figure 16. Experimental result and DYNA3D 

prediction for axial velocity at the 

back of buffer layer for the 

PC74SS37 laminate. 

0

50

100

150

200

250

300

0E+00 4E-06 8E-06

Time, s

Axial Velocity, m/s

DYNA3D

Experiment

Figure 17. Experimental result and DYNA3D 

prediction for axial velocity at the 

back of buffer layer for the 

PC74GS55 laminate. 

The AEH-PARADYN analyses were run assuming that the laminates were all homogeneous in 

the global sense, but heterogeneous at microlevel. Two-layer microstructural AEH models were 

used to bring in the lamination effect of the hard and soft layers. At the microlevel, the layer Axial velocity (m/s) Axial velocity (m/s) 

Time (s) 

Time (s) 
23

constituents are assumed to exhibit elastic/plastic behavior that is similar to the DYNA3D’s 

elastic/plastic Material Type 3. The predictions from these analyses for the axial velocity at the 

back of the buffer layer are compared with the experimental results in figures 18–20. For all the 

three laminates, the predictions compared well, especially in following the initial slope, peak 

velocity, and the post peak unloading. 

This example clearly demonstrates the ability of the AEH microstructural model to accurately 

predict the shock wave propagation in highly dissimilar layered materials, and the capability of 

the included AEH microstructural model for conducting multiscale macro/micro analyses with 

AEH-PARADYN. 

0

50

100

150

200

250

300

0.E+00 4.E-06 8.E-06 Time, s

Axial Velocity, m/s

AEH-PARADYN

Experiment

Figure 18. Experimental result and AEHPARDYN prediction for axial 

velocity at the back of buffer 

layer for the PC74AL37 

laminate. 

Time (s) 

Axial velocity (m/s) 
24

0.E+00

5.E+01

1.E+02

2.E+02

2.E+02

3.E+02

3.E+02

0.E+00 4.E-06 8.E-06 Time, s

Axial Velocity, m/s

AEH-PARADYN

Experiment

Figure 19. Experimental result and AEHPARDYN prediction for axial 

velocity at the back of buffer layer for 

the PC74SS37 laminate. 

0.E+00

5.E+01

1.E+02

2.E+02

2.E+02

3.E+02

3.E+02

0.E+00 4.E-06 8.E-06 Time, s

Axial Velocity, m/s

AEHPARADYN

Experiment

Figure 20. Experimental result and AEHPARDYN prediction for axial 

velocity at the back of buffer 

layer for the PC74GS55 

laminate.

Time (s) 

Axial velocity (m/s) 

Axial velocity (m/s) 

Time (s) 
25

5. Future Work and Conclusions 

Microstructures are features of the U.S. Army’s advanced structures and armor materials. 

Homogenization theories such as the AEH are being developed to model the effects of these 

microstructures on global material and structural responses. With the present coupling of the 

AEH-based microstructural model in PARADYN, many in-situ material responses can be 

studied under diverse global structural conditions. The sample applications presented in this 

report demonstrate this. Efforts underway include the following extensions to this work: 

• using PARADYN’s diverse material library to model different microphase materials, 

• modifying the AEH model used in this work to include an equation of state, and 

• extending the application range to include complex microstructures and diverse global 

applications.
26

6. References 

Bensoussan, A.; Lions, J. L.; Papanicolaou, G. Asymptotic Analysis for Periodic Structures; 

New York: North Holland Publishing Co., 1978. 

Chung, P. W.; Namburu, R. R.; Tamma, K. K. Three-Dimensional Elasto-Plastic Heterogeneous 

Media Subjected to Short Transient Loads. AIAA Paper-99-1239 1999. 

Clarke, J. A.; Namburu, R. R. A Distributed Computing Environment for Interdisciplinary 

Applications. Concurrency and Computation: Practice and Experience; Cocurrency 

Computa: Pract. Exper. 2002, 14, 1161–1174. 

Fish, J.; Shek, K.; Pandheeradi, M.; Shepard, M. Computational Plasticity for Composite 

Structures based on Mathematical Homogenization: Theory and Practice. Computer 

Methods in Applied Mechanics and Engineering 1997, 148, (1–2), 53–73. 

Ghosh S.; Moorthy, S. Elasto-Plastic Analysis of Arbitrary Heterogeneous Materials With the 

Voronoi Cell Finite-Element Method. Computer Method in Applied Mechanics and 

Engineering 1995, 121, 373–409. 

Guedes, J. M.; Nikuchi, N. Preprocessing and Post-Processing for Materials Based on the 

Homogenization Method With Adaptive Finite-Element Methods. Computer Methods in 

Applied Mechanics and Engineering 1990, 83, (2), 143–198. 

Hashin, Z. The Elastic Modulii of Heterogeneous Materials. ASME J. Appl. Mech. 1962, 29, 

143–150. 

Hashin, Z.; Strickman, S. A. Variational Approach to the Theory of the Elastic Behavior of 

Multiphase Materials. Journal Mech. Phys. Sol. 1963, 11, 127–140. 

Hill, R. Theory of Mechanical Properties of Fiber Toughened Materials. J. Mech. Phys. Sol. 

1964, 12, 99–212. 

Lene, F. Contribution à Lètude des Matériaux Composites et de Leur Endommagement. Thesis, 

Université Paris 6, 1984. 

Lene, F. Damage Constitutive Relations for Composite Materials. Engineering Fracture 

Mechanics 1986, 25, (5/6), 713–728. 

Licht Frottement, C. Viscoplasticité et Homogenization. Thesis, Université Montpellier 2, 1987. 

Moulinec, M.; Suquet, P. A Fast Numerical Method for Computing the Linear and Nonlinear 

Mechanical Properties of Composites. C. R. Acad. Sci. Paris 318, 1994, 1417–1423. 
27

Namburu, R. R. Unified Finite-element Methodology for Thermal-Structural Problems. Thesis, 

The University of Minnesota, July 1991. 

Sanchez-Palencia, E. Non-Homogeneous Media and Vibration Theory – Lecture Notes in 

Physics; Springer Verlag, 1980. 

Suquet, F. Plasticité et Homogenization. Thesis, Université Paris 6, 1982. 

Terada, K.; Kikuchi, N. Nonlinear Homogenization Method for Practical Applications, 

Computational Methods in Micro-mechanics. American Society of Mechanical Engineers, 

Applied Mechanics Division, 1995, 212, 1–16. 

Valisetty, R. R.; Namburu, R. R.; Chung, P. W. Scalable Implementation of Three-Dimensional 

Heterogeneous Media Subjected to Short Transient Loads; ARL-TR-2351; U.S. Army 

Research Laboratory: Aberdeen Proving Ground, MD, 2000. 

Willis, J. Elastic Theory of Composites. In Mechanics of Solids; Oxford: Pergamon Press, 

pp 353–386, 1982. 

Zhuang, S. Shock Wave Propagation in Periodically Layered Composites. Ph.D. Thesis, 

California Institute of Technology, Pasadena, CA, 2002. 
28

Appendix. Microstructural Equations

A.1 Governing Equations 

These equations are presented in several references, most recently in Chung et al.1 Assume a 

three-dimensional body Ω is an assembly of periodic structures containing different materials, as 

shown in figure A-1. Typically, the unit cell is very small and of order ε (where ε is a small 

positive number) compared to the dimensions of the problem domain. Structural response 

quantities such as displacements, velocities, stresses, and strains are assumed to have slow 

variations (macroscopic) from point to point as well as fast (microscopic) variations within a 

small neighborhood ε of a given point x. Let Y be a (periodic) representative part of Ω. 

Composite 

Neighborhood x’ 

Figure A-1. Macro/microanalyses neighborhoods. 

Here we distinguish two scales: the macroscopic scale (x ∈Ω) and microscopic scale 

(y∈Ω). Lagrangian representations of conservation of mass and momentum equations are used. 

The conservation of mass is used to calculate the current density from the initial density. The 

momentum equation, the kinematic relations, and the constitutive model are represented by 

i ji j i , v f εε ε

ρ  −σ = , (A-1) 

.

, ,

1 ( )

2 ij i j j i e vv ε εε = + , (A-2) 

and 

 1

Chung, P. W.; Namburu, R. R.; Tamma, K. K. Three-Dimensional Elasto-Plastic Heterogeneous Media Subjected to Short 

Transient Loads. AIAA paper-99-1239 1999. 
29

(,) ij ij ij SeE

ε

ε σ ⋅

= , (A-3) 

where ρ is the density, i f is the body force vector, E is the internal energy, the superscript ε

denote micro/macro continuum solutions, vi is the velocity vector, σ ij is the stress tensor, and ij e

is the strain tensor. The initial conditions in the domain ε Ω and boundary conditions on the 

surface of the domain ( ) 1 Γ2 Γ = ∂Γ ∪ ∂ are given by 

ε ε

oi i u (t = 0) = u , (A-4) 

( 0) i i o vt v ε ε = = , (A-5) 

1 on o

i i u u ε ε = ∂Γ , (A-6) 

and 

2 on ji i i n T ε σ = ∂Γ , (A-7) 

where ui is the displacement vector and ni is a surface normal. 

A.2 Time Integration 

A second-order accurate Lax-Wendroff-based explicit time integration procedure Tamma and 

Namburu2

 is employed for the conservation equation or equation of motion. In this approach, 

the dependent variable velocity is first discretized in time using a second order Taylor series 

expansion. +1 +1/ 2

= + ∆ n n n

i i i v v tv ε ε ε

 . (A-8) 

Taking appropriate time derivatives of equation 1 and substituting for the velocity and 

acceleration gives 

v v ( f ) t n i i ji j i

n n n

− = + ∆ + + + 1/ 2 ,

1 1 / 2

ρ ρ σ ρ ε ε ε . (A-9) 

The stress increment is related to a strain increment through an appropriate constitutive equation. 

The elasto-plastic constitutive equation at midpoint time increment is defined by 

( ) 1/ 2 1 / 2 1 / 2 1 / 2 + = + ∆ + − − n

ij ijkl ij klmn mn M C t n n n ε ε ε σ σ ε , (A-10) 

where Mijkl and Cijkl are radial return in plasticity tensor and elasticity tensor. 

Stresses and strain rates at midpoint time can be evaluated using the predicted displacements as 

shown by

 2

Tamma, K. K.; Namburu, R. R. A Robust Self Starting Explicit Computational Methodology for Structural Dynamic 

Applications: Architecture and Representations. International Journal for Numerical Methods in Engineering 1990, 29, 1441–

1454. 
30

1 / 2 1 / 2

2

+ ∆ +

= + n n n

i i i v

t

u u ε ε ε . (A-11) 

A.3 Asymptotic Expansion and Elasto/plastic Constitutive Equations 

The homogenization method is based on the asymptotic expansion of the primary variables 

together with a unit cell approach for a heterogeneous structure. Assume a 3-D body Ω is an 

assembly of periodic structures. Typically, the unit cell is very small, of order ε (ε is a small 

positive number) compared to the dimensions of the problem domain. Typically, two-scale 

asymptotic expansion can be employed to approximate the displacement, velocity, strain, or 

stress fields. 

2

1 2 (, ) (, ) (, ) (, ) , ii i i o u xy u xy u xy u xy ε =+ + + ε ε " (A-12) 

2

1 2 (, ) (, ) (, ) (, ) ii i i o v xy v xy v xy v xy ε =+ + + ε ε ", (A-13) 

2

1 2 (, ) (, ) (, ) (, ) ii i i o e xy e xy e xy e xy ε =+ + + ε ε ", (A-14) 

and 

2

1 2 (, ) (, ) (, ) (, ) , ii i i o xy xy xy xy ε σ σ εσ ε σ =+ + +" (A-15) 

where (, ) i u xy ε , ( , ) i v xy ε , ( , ) i e xy ε , and ( , ) i x y ε σ are Y-periodic functions. 

The asymptotic expansion homogenization (AEH) approach for nonlinear applications is based 

on the instantaneously linearized assumption for the constitutive model.3

 The fundamental 

assumption of the AEH approach, therefore, is that the true solution in the ε space is decomposed 

into a macro space x and a micro space y. The basic assumption here is that multiple scales exist 

only in the spatial variables and that no such scaling exists for the time variable. Further reading 

in space and time asymptotic expansion approaches can be found in Bensoussan.4 The general 

approach to heterogeneous problems is to separate and draw clear distinction between the microand macrolevel equilibrium equations regardless, per the earlier linearized assumption, of 

material nonlinearity. This is accomplished by asymptotically expanding the primary variables 

where the asymptotic scale is approximated to the second order. 

2

1 2 (, ) (, ) (, ) (, ) ii i i o v xy v xy v xy v xy ε =+ + + ε ε ". (A-16) 

Spatial gradients in ε-space are taken with respect to the x-coordinate system. The 

corresponding gradient for Y-periodic functions (where micro and macro are now 

distinguishable) is given by the chain rule where the scaling is defined by y = x/ε; hence, for any 

Y-periodic function φ, the earlier gradients are replaced by 

 3Terada, K.; Kikuchi, N. Nonlinear Homogenization Method for Practical Applications, Computational Methods in Micromechanics. American Society of Mechanical Engineers, Applied Mechanics Division, 1995, 212, 1–16. 

4

Bensoussan, A.; Lions, J. L.; Papanicolaou, G. Asymptotic Analysis for Periodic Structures; New York: North Holland 

Publishing Co., 1978. 
31

1

ii i

d

dx x y

φ φ φ

ε

∂ ∂

= +

∂ ∂ . (A-17) 

Then, using equation A-1 in equation A-7 while considering equation A-2 gives the rate of strain 

tensor defined by 

1 12 [ ( ) ( )] [ ( ) ( )] x y xy

ij ij o ij ij ij ev ev ev ev ε

ε ε  =++ ++", (A-18) 

where the symmetrized gradient tensors, e

x

ij and e

y

ij are defined by 

1 () ( ) 2

x j i

ij

j i

e

x x

φ φ φ ∂ ∂

= +

∂ ∂ (A-19) 

and 

1 () ( ) 2

y j i

ij

j i

e

y y

φ φ φ ∂ ∂

= +

∂ ∂ . (A-20) 

The relevant expressions for gradients and velocities are first substituted into the governing 

equations of motion. Next, the micro and macro equations are identified by selecting the 

appropriate coefficients to the scaling factors ε which must each be identically zero. Substituting 

the asymptotically-expanded velocities, equation A-1, in the strain rate, equation A-3, and then 

using the constitutive equation, equation A-16, in the equations of motion, equation A-15, yields 

a set of equations dependent on powers of ε. To satisfy the equations of motion, each term 

associated with each of the powers of ε must approach zero identically. This leads to a set of 

equations associated with the microscopic and macroscopic equations of motion. The first 

equation, associated with the powers of −2 ε , is given by 

1

2 * : 0 k

n

o

ijkl

j l

v

C

y y

ε

−

− ∂ ∂ = ∂ ∂ , (A-21) 

where Cijkl = MijpqCpqkl

* . Equation A-21 states that n-1

oi v is a function only of x. Hence, 

derivatives of y are zero. Using this inference, the equation associated with the powers of −1 ε is 

derived as 

( ). 1: 1/ 2

1

* 1

1

1 * −

− −

−

∂

∂

∆ = − ∂

∂

∂

∂

+

∂

∂

∂

∂ n

ijkl kl

l j

n

ijkl

l j

n

o

ijkl

j

M

y t y

v

C

x y

v

C

y

k k ε σ (A-22) 

Equation A-22 relates the perturbative velocity field, v1, to the macroscopic velocity field, v0. 

The equation relating these two quantities will become the important micro/macroequation, 

providing the direct link to relate microscopic to the macroscopic velocities. It provides the 

corrections which account for the shape of the interface separating various phases and the timedependent plastic softening effect due to material nonlinearities. For a homogeneous material 

problem, where gradients over y of l

n

ijkl o C v x i ∂ ∂ − / * 1 are zero, and in which plastic effects have no y
32

dependence, the solution to equation A-22 1

1

n−

i v is identically zero. In the special case where no 

plastic yielding occurs, the problem degenerates to a transient elastic homogenization problem. 

Finally, the equation associated with the powers of ε

o

 is written as 

1/2

1 1

0* * 1

1 1

* * 2

1 1/2

2

: 

1 .

k k

k k

n

i

n n

o

ijkl ijkl

j lj l

n n

o

ijkl ijkl

j lj l

n n

o ijkl kl i

j

v v

C C

x xx y

v v

C C

y xy y

vM f

t tx t

ε

ε

ρ ρ σ −

− −

− −

+ +

∂ ∂ ∂ ∂

+

∂ ∂∂ ∂

∂ ∂ ∂ ∂

+ +

∂ ∂∂ ∂

∂ = ∆− −

∆ ∆∂ ∆

 (A-23) 

Much like equation A-22 provided the basis for the microscale problem, equation A-23 provides 

the basis for the macroscale problem. By taking the volume average of equation A-23, 

substituting for 1 v using the solution of equation A-22 and noting the periodicity, the final 

governing equation of motion is given by 

+1 +1/ 2 +1/ 2 + + ∆

∂

∂ ∆ = ∆ n

i

c

ij

n

ij

j

n

o t f

x

v t i ρ σ σ ρ , (A-24) 

where the “corrector stress,” c σ ij , is defined by 

l

n

ijkl

c

ij y

v

C k

∂

∂ =

−1

* 1

σ , (A-25) 

and where quantities in brackets denote 

1 ( ) ( ) Y

dY

V

= ∫ . (A-26) 

It is of interest to note that equation A-24, for heterogeneous conditions, is similar in form to 

equation A-6, the homogeneous equation, making the integration of a micromechanical problem 

into a macroanalysis tractable and straightforward. 

A.4 Constitutive Equation: Elasto/Plastic 

The stress tensor is split into two parts. The first part, Sij, is deviatoric stress, which is related to 

material strength, and the second part is pressure, P. To simplify notation, the superscript ε has 

been omitted. 

σ δ ij ij ij = − + P S (A-27) 

and 

1 – 3

σ ii = P . (A-28)
33

The first term in equation A-27 accounts for volumetric changes and is typically evaluated from 

equations of state in explicit dynamics. The second term, the stress deviator, is related to the 

deformation of the material and is defined by a constitutive model, and in particular, the time rate 

of change of the deviator is evaluated from the strain rates. 

Sij = G ij +∆ij 2 . (A-29) ε

This model can be defined with any number of local nodes, elements, and materials, but since 

this model is called repeatedly for each time step, for each element for which microstructural 

response is deemed important, a complex microstructure model can overwhelm the overall 

computations, therefore keeping the model simple is important. ∆ij is the correction for rigid 

body rotation. 

An isotropic hardening model with a rate-dependent Von-Mises yield condition is employed in 

the present analysis. A consistency condition ensures that the stress state remains on the yield 

surface at the start and end of a time-step. Such a condition is given by 

n n 1 1 R Q S ij ij+ + = . (A-30) 

The variable R is the rate-dependent radius of the yield surface or the apparent yield stress, Qij is 

the vector specifying the normal direction to the yield surface, and ij S is the deviatoric 

component of the stress. The stresses are understood to be the co-rotated second-order tensor 

according to the Jaumann definition of the co-rotational derivative. The Jaumann derivative of 

the Cauchy stress, 

.

σ j

ij, is related to the material time derivative, 

.

σ ij, by 

ij ik kj ik kj

J σij = σ +σ ω +ω σ , (A-31) 

where ωij is the rotation tensor, the standard skew-symmetric component of the velocity gradient. 

The radius of the yield surface is defined by 

1

1 p

o p

Y

e R H

D

σ ε

⎡ ⎤ ⎛ ⎞ =+ + ′ ⎢ ⎥ ⎜ ⎟ ⎢ ⎥ ⎝ ⎠ ⎣ ⎦

, (A-32) 

where σ

y

o is the static yield stress of the material, and D and p are the so-called fluidity 

parameters. H' is the hardening parameter, and − p

ε is the effective plastic strain. The effective 

rate of deformation, −

ε , is defined by 

ij ij e e e 3

2 = , (A-33) 

where eij is the deviatoric component of the total strain, εij. The normality condition, specified 

by the normal vector Qij, is given by 

ij

ij

mn mn

S

Q

S S = . (A-34)
34

For a homogeneous material (which is applicable in the present derivation since the constitutive 

equations are applied at the microlevel where the material is homogeneous within each phase), 

the incremental relationships for the apparent yield stress and the deviatoric stress are defined as 

1 2

3

n n RRH λ + = + ∆′ (A-35) 

and 

1 1 2 T n n

ij ij ij S S GQλ + + = −∆ , (A-36) 

where T n

ij S +1 is the deviatoric trial stress, ∆λ is a scalar quantity representing the magnitude of 

the radial return correction due to plastic yielding, and G is the elastic shear modulus. The 

deviatoric trial stress is defined by 

S S C e t n

ijkl kl

n

ij

n

ij

T

= + ∆ +  1 , (A-37) 

where Cijkl is the tensor containing elastic properties. 

Finally, substituting equation A-37 and A-35 into equation A-30 and solving for the radial return 

correction parameter gives 

3/2 3 1 1

(3 ') 2

T T n n

ij ij n SS R

G H

λ + + ⎛ ⎞ ∆= − ⎜ ⎟ + ⎝ ⎠

. (A-38) 

Using equation A-38 in equation A-36 and employing the standard assumption that a material is 

elastic in its dilatational behavior and plastic only in shear gives the total stress increment in 

equation A-16. The constitutive equation A-16 can now be employed in the derivation of the 

micro and macro governing equations. 
NO. OF NO. OF 

COPIES ORGANIZATION COPIES ORGANIZATION

35

 1 DEFENSE TECHNICAL 

 (PDF INFORMATION CTR 

 ONLY) DTIC OCA 

 8725 JOHN J KINGMAN RD 

 STE 0944 

 FORT BELVOIR VA 22060-6218 

 1 US ARMY RSRCH DEV & 

 ENGRG CMD 

 SYSTEMS OF SYSTEMS 

 INTEGRATION 

 AMSRD SS T 

 6000 6TH ST STE 100 

 FORT BELVOIR VA 22060-5608 

 1 INST FOR ADVNCD TCHNLGY 

 THE UNIV OF TEXAS 

 AT AUSTIN 

 3925 W BRAKER LN STE 400 

 AUSTIN TX 78759-5316 

 1 US MILITARY ACADEMY 

 MATH SCI CTR EXCELLENCE 

 MADN MATH 

 THAYER HALL 

 WEST POINT NY 10996-1786 

 1 DIRECTOR 

 US ARMY RESEARCH LAB 

 IMNE ALC IMS 

 2800 POWDER MILL RD 

 ADELPHI MD 20783-1197 

 3 DIRECTOR 

 US ARMY RESEARCH LAB 

 AMSRD ARL CI OK TL 

 2800 POWDER MILL RD 

 ADELPHI MD 20783-1197 

 3 DIRECTOR 

 US ARMY RESEARCH LAB 

 AMSRD ARL CS IS T 

 2800 POWDER MILL RD 

 ADELPHI MD 20783-1197 

ABERDEEN PROVING GROUND

 1 DIR USARL 

 AMSRD ARL CI OK TP (BLDG 4600) 
36

INTENTIONALLY LEFT BLANK.