University of Southern California Dissertations and Theses

Multi-scale modeling of functionally graded materials (FGMs) using finite element methods

(USC Thesis Other)

Multi-scale modeling of functionally graded materials (FGMs) using finite element methods

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content MULTI-SCALE MODELING OF FUNCTIONALLY GRADED MATERIALS
(FGMs) USING FINITE ELEMENT METHODS
by
Richard Sangwon Rhee
_______________________________________________________________
A Dissertation Presentation to the
FACULTY OF THE GRADUATE SCHOOL
UNIVERSITY OF SOUTHERN CALIFORNIA
In Partial Fulfillment of the
Requirements for the Degree
DOCTOR OF PHILOSOPHY
(CIVIL ENGINEERING)
December 2007
Copyright 2007 Richard Sangwon Rhee
ii
Dedication
This dissertation is dedicated to:
My parents who gave me everything and now they are in heaven;
M. Kim, who provides tremendous support;
My two hope, Ryan and Kyle;
Terry, who sacrifices during my education pursuits;
My elder sisters and brother-in-law, Jin Y. Sohn;
And J. Jackman, who gave me wisdom and encouragement.
iii
Acknowledgements
I would like to express my sincere thanks and gratitude to my advisor, Prof. L.
C. Wellford, for giving me a chance to pursuit the Ph.D. program at Civil Engineering
department and during academic years, he always shared his time in order to provide
the information which I need. Also, I would like to thank to my dissertation
committee members, Prof. V. Lee, Prof. S. Masri, Prof. J. Anderson and Prof. Shiflett
for their valuable time and suggestions.
Sincere thanks to Kem Schankereli and Keith Myers, who encouraged me from
the beginning of my graduation study. Also special thanks to my friends, Dr. H. J.
Kwon and Dr. S. Kulkani, for their guidance. Last but not least, I would like to thank
to my sworn brothers, Prof. S. J. Lee, Prof. K. D. Park and I. K. Jung, who always
make me to believe my dream.
iv
Table of Contents
Dedication ii
Acknowledgements iii
List of Tables vii
List of Figures ix
Abstract xiv
Chapter 1: Introduction 1
Chapter 2: Problems Involving a Micro Structure
2.1 Review of Homogenization Theory 9
2.2 Analytical Solution of the Homogenization in 2-D 18
2.3 Macro Homogenized Solution 28
2.4 Variation of Cell Solution
kl
i
in Microcell Geometries 30
2.4.1 A Typical Shape of Microcell Geometry 30
2.4.2 Computation of the Cell Solution
kl
i
in Soft and
Hard Microcell Structure Cases 31
2.5 Verification of the Homogenized Elastic Constants 37
2.5.1 Soft and Hard Isotropic Composite Material with
Variation 37
2.5.2 A Rectangular Hole in the Material 40
Chapter 3: Nonperiodic Homogenization (NPH) Method
3.1 Mathematical Theory of NPH 43
3.2 Evaluation of the NPH Strain Tensors 47
3.3 Explicit Parameterization of the Cell Geometry over the
Element 53
3.4 Coordinate Transformation 55
3.5 NPH Element Stiffness Matrix 58
3.6 NPH Element Stress Calculation 62
3.7 Typical Allowable NPH Microcell Geometries 64
3.8 Data Specification points for hx
2
& hy
2
Matrix Values in
the Nonperiodic Microstructure 66

v
Chapter 4: Numerical Examples of NPH 1-D Cases
4.1 Introduction 68
4.2 Case 1 – Comparison between the NPH and the Homogenization
Solutions 69
4.2.1 Geometry Modeling 69
4.2.2 Local and Global Deformation Analysis 70
4.3 Case 2 – Descending Low Density Microcell Structure 72
4.3.1 Geometry Modeling 73
4.3.2 Local and Global Deformation Analysis 73
4.4 Case 3 – Descending High Density Microcell Structure 75
4.4.1 Geometry Modeling 75
4.4.2 Local and Global Deformation Analysis 76
4.5 Case 4 – Descending and Ascending Microcell Structure 79
4.5.1 Geometry Modeling 79
4.5.2 Local and Global Deformation Analysis 80
4.6 Case 5 – Descending Microcell Structure with a Sudden Jump 82
4.6.1 Geometry Modeling 82
4.6.2 Local and Global Deformation Analysis 83
4.7 Case 6 – Rapidly Varying Descending, Ascending and
Descending Microcell Structures 85
4.7.1 Geometry Modeling 85
4.7.2 Local and Global Deformation Analysis 86
Chapter 5: Numerical Examples of NPH 2-D Cases
5.1 Introduction 88
5.2 Case 7 – Periodic Microstructure 88
5.2.1 Geometry Modeling 88
5.2.2 Local and Global Deformation Analysis 90
5.3 Case 8 – Descending Horizontal Fiber Strips in One Direction 92
5.3.1 Geometry Modeling 92
5.3.2 Local and Global Deformation Analysis 94
5.4 Case 9 – Descending and Ascending FGMs with Square Fibers 96
5.4.1 Geometry Modeling 96
5.4.2 Local and Global Deformation Analysis 98
5.5 Case 10 – Descending and Symmetric Matrix Structure 101
5.5.1 Geometry Modeling 101
5.5.2 Local and Global Deformation Analysis 104
5.5.3 Local Stress Analysis 106
5.6 Case 11– Descending Horizontal and Vertical Fiber Strips 110
5.6.1 Geometry Modeling 110
5.6.2 Local and Global Deformation Analysis 113
5.6.3 Local Stress Analysis 115
vi
Chapter 6: Summary and Conclusions
6.1 Summary 118
6.2 Conclusions 121
Bibliography 123
Appendix A
Nonperiodic Homogenization (NPH) Cell Solution 126
vii
List of Tables
Table 2.1 Comparison of the Homogenized Elastic Constants: Soft and Hard
Material 39
Table 2.2 Elastic Tensors after Homogenization in Rectangular Hole 42
Table 4.1 Deformation Results of NPH and Homogenization Methods 71
Table 4.2 Deformation Results for the Descending Low Density Microcell
Structure Case 74
Table 4.3 The Matrix Size Values for the Descending and Ascending
Microcell Structure: Fiber size hx
1
= 0.0262 76
Table 4.4 The Deformation Results of the High Density Microcell Structure
with 10 percent Reduction 77
Table 4.5 The Matrix Size Values for the Descending and Ascending
Microcell Structure: Fiber size hx
1
= 0.0262 80
Table 4.6 Global Deformation Values in High Density Microcell Structure
Case with 10 percent Descending and 16 percent Ascending
Matrix Sizes 80
Table 4.7 The Matrix Size Values for the Descending Microcell Structure
with a Sudden Jump: Fiber size hx
1
= 0.0194 83
Table 4.8 Table 4.8 Deformation Results for the High Density Descending
Microcell Structure with a Sudden Jump 84
Table 4.9 Table 4.9 The Matrix Size Values for the Rapidly Varying
Descending, Ascending and Descending Structure: Fiber size hx
1
= 0.070 86
Table 4.10 Deformation Results of the High Density Microcell with Rapidly
Varying Descending, Ascending and Descending Microcell
Structure 86
Table 5.1 Comparison of the Deformation between Homogenized Solution
and NPH Solution @ X
1
=1.0, 1.5 & 2.0 91
Table 5.2 Matrix Size Values for the Descending Horizontal Fiber Strips 93
viii
Table 5.3 Summary of Model Sizes for Case 8 94
Table 5.4 Matrix Sizes in Descending and Ascending Square Fibers 97
Table 5.5 Summary of Model Sizes for Case 9 98
Table 5.6 Matrix Sizes in Symmetric and Descending Case 103
Table 5.7 Summary of Model Sizes for Case 10 104
Table 5.8 Von-Mises Stress in Descending and Symmetric Case 109
Table 5.9 Matrix Sizes in Descending Horizontal and Vertical Strips Case 112
Table 5.10 Summary of Model Sizes for Case 11 113
Table 5.11 Von-Mises Stress in Horizontal and Vertical Strips 116
ix
List of Figures
Figure 1.1 Fiber-Reinforced Polymer Bridge 2
Figure 1.2 Reinforcement types of Metal Matrix Composites (MMCs) or
Ceramics Matrix Composites (CMCs); a) Matrix with Fibers b)
whiskers c) particulates 2
Figure 1.3 Different types of the Functionally Graded Materials (FGMs): a)
Continuously Graded Microstructure; b) Discretely Graded
Microstructure with fiber and matrix configuration; c) Multi-
phase Graded Microstructures. 2
Figure 1.4 Replacement scheme used in the layer-wise Homogenization
model in Continuously Graded Microstructure 4
Figure 1.5 Equivalent continuum with layers wise Homogenization in
Discretely Graded Microstructure with fiber and matrix
configuration (MMCs / CMCs types) 6
Figure 2.1 Periodic Microstructure 10
Figure 2.2 Associated Macrostructure 10
Figure 2.3 Base cell of the Composite Microstructure 10
Figure 2.4 Microcell Structure: Model 1 – Cell Size 1.0 x 1.0; Model 2 –
Cell size 0.5 x 1.0 30
Figure 2.5 Microcell Solution for Model 1 – Cell size 1.0 x 1.0 31
Figure 2.6 Microcell Solution for Model 2 – Cell size 0.5 x 1.0 32
Figure 2.7 Superimpose the Cell Solution Pictures of Figure 2.5 and
Figure 2.6 32
Figure 2.8 Microcell Solution for Model 1 – Cell size 1.0 x 1.0 33
Figure 2.9 Microcell Solution for Model 2 – Cell size 0.5 x 1.0 33
Figure 2.10 Superimpose the Cell Solution Pictures of Figure 2.8 and
Figure 2.9 34
x
Figure 2.11 Microcell Solution for Model 1 – Cell Size 1.0 x 1.0 35
Figure 2.12 Microcell Solution for Model 2 – Cell Size 0.5 x 1.0 35
Figure 2.13 Superimpose the Cell Solution Pictures of Figure 2.11 and
Figure 2.12 36
Figure 2.14 Material Properties in a Single Microcell Geometry: (E
soft
= 10
and E
hard
= 1000) 37
Figure 2.15 A Single Microcell Geometry Model Defined by Bendoe and
Kikuchi [3]: 16 x 16 Uniform Square Four-node Isoparametric
Elements 38
Figure 2.16 A Single Microcell Geometry Model Defined by HOMOG.f90:
8 x 8 Uniform Square Nine-node Lagrange Quadratic Elements 39
Figure 2.17 Material Properties in a Single Microcell Geometry with a Hole:
(E
1111
= E
2222
= 30 and E
1122
= E
1212
= 10) 40
Figure 2.18 A Single Microcell Geometry Model Defined by Bendoe and
Kikuchi [3]: 20 x 20 Uniform Square Four-node Isoparametric
Element 41
Figure 2.19 A Single Microcell Geometry Model Defined by HOMOG.f90:
10 x 10 Uniform Square Nine-node Lagrange Quadratic
Elements 41
Figure 3.1 Nodal Degree-of-Freedom in 1-D Element 45
Figure 3.2 Nodal Displacement at the Mesh Element and Microstructure:
i.e.) E
1
= 10 and E
2
= 1000 63
Figure 3.3 A Microcell Structure in 2-D (Shaded areas Represent Fibers,
hx
1
and hy
1
) 64
Figure 3.4 Variation of 2-D Microcell Structures with Fiber (shaded area):
(a) Vertical Direction; (b) Center; (c) Outer Edges; (d)
Distributed in the Center and Four Corners 65
Figure 3.5 A Total of Three hx
2
Data Collection Points per One Finite
Element: Black Dots indicate Data Collection Points and Circles
with Cross indicate the Location of Gauss Points. 66
xi
Figure 3.6 A Total of Nine Data Collection Points per Finite Element: Black
Dots indicate Data Collection Points and Circles with Crosses
indicate the Location of Gauss Points. 67
Figure 4.1 A Periodic Structure with Microcells 69
Figure 4.2 A Single Macro Element for the NPH Method 70
Figure 4.3 Detailed View of the Cell 1 for the Exact Solutions: hx
2
= 0.375
and hx
1
= 0.250 70
Figure 4.4 Three Different hx
2
Data Collection Methods; Black Dots
indicate Data Collection Points and Circles with Cross indicate
the Location of Gauss Points. 72
Figure 4.5 Nonperiodic Low Density Microcell Structure 73
Figure 4.6 Deformation Results for the Descending Low Density Microcell
with Three Different Data Collection Methods 74
Figure 4.7 Descending Microcell Structure Case with 10 percent
Descending 75
Figure 4.8 Deformation Results of Descending Microcell Structure in High
Density Microcell with Three Different Data Collection Methods 78
Figure 4.9 Detailed View of the Figure 4.5 (Red Box Region) 78
Figure 4.10 High Density Microcell Structure Case with 10 percent
Descending and 16 percent Ascending Matrix Sizes 79
Figure 4.11 Deformation Results for High Density Descending and
Ascending Microcell Structure 81
Figure 4.12 Detailed View of the Figure 4.5 (Red Box Region) 81
Figure 4.13 Descending Micro Structure with a Sudden Jump 82
Figure 4.14 Deformation Results of the High Density Descending Micro
Structure with a Sudden Jump 84
Figure 4.15 High Density Microcells with Rapidly Varying Descending,
Ascending and Descending Structure 85
xii
Figure 4.16 Deformation Results of the High Density Microcell with Rapidly
Varying Descending, Ascending and Descending Structure 87
Figure 5.1 Periodic Microcell Geometry and Boundary Conditions 89
Figure 5.2 Homogenization and NPH Deformation Results for the Periodic
Microstructure at X
1
= 1.0, X
1
= 1.5 and X
1
= 2.0 90
Figure 5.3 Descending Horizontal Fiber Strips Microcell Structure: (a)
Descending Horizontal Fiber Strips in X
2
Direction, (b) 16
Macro Elements for NPH 92
Figure 5.4 Detailed View of the Microcell Strip: Fiber Value hy
1
= 0.030
(Constant) 93
Figure 5.5 Conventional FEA Mesh and Boundary Conditions 94
Figure 5.6 Local Deformation Results: the NPH Solution and the FEA
Solution at X
1
= 1.995 95
Figure 5.7 Descending and Ascending Matrix with Square Fibers 96
Figure 5.8 Detailed View of the Square Fibers: Fiber Value hx
1
& hy
1
=
0.030 (Constant) 97
Figure 5.9 Conventional FEA Mesh Sizes and Boundary Conditions 99
Figure 5.10 Detailed View of the FEA Mesh Sizes and Boundary Conditions
(Red Box Region) 99
Figure 5.11 Deformation Results Using the Conventional FEA Method 99
Figure 5.12 Deformation Results of the Conventional FEA Method and the
NPH Method at X
1
= 1.9427 100
Figure 5.13 Descending Matrix Microstructures in X
1
Direction and
Symmetric condition in X
2
Direction with Descending
Microstructures 101
Figure 5.14 Detailed View of Figure 5.13 (Red Box Region) 102
Figure 5.15 Symmetric and Descending Structure: Deformation
Displacement at X
1
= 1.995 105
xiii
Figure 5.16 Detailed View of Symmetric and Descending Structure (Red
Box Region) 105
Figure 5.17 Deformation Results of the FEA Method (ABAQUS) 106
Figure 5.18 Von-Mises Stress Results from FEA (ABAQUS) 107
Figure 5.19 Detailed View of the High Stress Region (Red Box Region):
Max. Stress = 5.811 107
Figure 5.20 Von-Mises Stress Results by NPH method (Red Box Region in
Figure 5.19): Max. Stress Value = 7.473 108
Figure 5.21 Von-Mises Stress Results by Commercially Available FEA
(ABAQUS): Max. Stress Value = 5.777 109
Figure 5.22 Simultaneously Reducing Matrix Values in X
1
& X
2
Directions 110
Figure 5.23 Detailed View of the Microcell Structure 111
Figure 5.24 Deformation Results of Descending Horizontal and Vertical
Strips 114
Figure 5.25 FEA (ABAQUS) Deformation Result 114
Figure 5.26 Von-Mises Stress Results from FEA (ABAQUS) 115
Figure 5.27 Detailed View of Figure 5.26 (Red Box Region): Maximum
Stress Value is 11.35 116
Figure 5.28 Von-Mises Stress Results by the NPH Method at Integration
Points: Maximum Stress = 13.86 117
Figure 5.29 Von-Mises Stress Results by FEA (ABAQUS) at Integration
Points: Maximum Stress = 12.11 117
xiv
Abstract
Functionally Graded Materials (FGMs) have a gradual material variation from
one material character to another throughout the structure. Applications of these types
of materials have significant advantages in civil and mechanical systems including
thermal systems. Analyzing the FGMs at the microstructure level with the
conventional Finite Element Method (FEM) takes enormous pre-processing and
computational time due to the complex material characteristics at the microstructure
level. Essentially, the model contains too many degrees of freedom to be solved
economically.
The homogenization method has been successfully applied to solve periodic
microstructure problems. However, the development of analysis procedures for
structures with nonperiodic material or cell geometry, as occurs in graded materials,
has turned out to be a significant challenge.
A new method is developed which accurately models the nonperiodic
microstructure in FGMs. This method allows the efficient solution of nonperiodic
problems without requiring the simplification of the original models. The
performance of the developed theory is verified through the solution of appropriate
nonperiodic problems associated with graded materials. In the nonperiodic 1-D cases,
the global displacement U(x) was obtained and compared with the exact solution. At
the same time, the proposed data collection point method was investigated. In the
nonperiodic 2-D cases, the global displacement U(x) and the microstructural level
xv
displacements were computed. In the program, the Von-Mises Stress computation
process was included to evaluate the local stress values at the microstructure level and
the results were compared with very fine scale finite element calculations.
The performance of the developed nonperiodic homogenized (NPH) algorithm
indicates that it is a promising tool for estimating the FGMs characteristics in loaded
conditions. The method can be applied to estimate the global and local displacements
in nonperiodic geometries which contain continuously decreasing and/or increasing
microstructures.
1
Chapter 1
Introduction
Functionally Graded Materials (FGMs) are composite materials which have a
gradual variation of the volume fractions from one material to the other. The material
property changes are usually in one direction with two different materials being used.
The microstructure is organized in a nonperiodic manner. The concept of FGMs
originated in Japan during the space plane project in 1984. The Japanese scientists
were developing a thermal barrier material which could withstand a surface
temperature gradient of approximately 1000 K across a 10 mm cross section.
Traditional thermal barrier coatings (TBCs) have been applied with Ni-based alloys as
the oxidation resistant bond coat and a heat resistant ZrO
2
ceramic top coating.
However, conventional plasma sprayed TBCs have a problem of low durability during
thermal cycling and poor bond strength. The uniqueness of the FGMs is the ability to
produce a new composite material with a gradual composition variation from heat
resistant ceramics to fracture resistant metals. Applications of these types of FGMs
have significant advantages in civil and mechanical systems including thermal systems
(e.g. rocket heat shields, heat exchanger tubes, thermoelectric generators, wear-
resistant linings, diesel and turbine engines, etc.).
Fiber-reinforced polymer (FRP) is another application of FGMs for reinforcing
concrete materials as shown in Figure 1.1. The FRP materials improve the corrosion
resistance of the steel and enhance the life cycle of the material strength.
2
Figure 1.1 Fiber-Reinforced Polymer Bridge
Figure 1.2 Reinforcement types of Metal Matrix Composites (MMCs) or
Ceramics Matrix Composites (CMCs); a) Matrix with Fibers b)
whiskers c) particulates
Figure 1.3 Different types of the Functionally Graded Materials (FGMs):
a) Continuously Graded Microstructure; b) Discretely Graded
Microstructure with fiber and matrix configuration; c) Multi-phase
Graded Microstructures.
(a) (b) (c)
(a) (b) (c)
3
FRP has been used in several bridge decks recently constructed in North America; The
Morristown Bridge, in Vermont, has an entire concrete deck slab constructed using
glass FRP (GFRP) bars. Also, there are laminate types of FRP which are contained in
an arrangement of unidirectional fibers or woven fiber fabrics embedded in a thin
layer of light polymer matrix materials - polyester, Epoxy or Nylon, etc.
Other major benefits of the FGMs are in the design and manufacturing of
Metal Matrix Composites (MMCs) and Ceramic Matrix Composites (CMCs). The
associated manufacturing processes provide the best material properties for
composites of metal and ceramics by, for examples, removing the brittleness of
ceramics and making a strong metal lighter and stiffer. The proportions of the matrix
alloy (the metal) and the reinforcement material (the ceramic), as well as shape and
location of reinforcement can be varied form place to place in the structure to achieve
particular desired properties. For example, ceramic reinforcements in the form of
fibers, whiskers or particulates can be introduced into the metal in a varying density
pattern as shown in Figure 1.2. However, the use of this type of material requires an
explicit understanding of the material behavior at each location and over all length of
scales.
Many more applications of the FGMs can be found in the conference papers
and technical journals including those related to solar energy conversion devices [20].
Most of the FGMs microstructures are fabricated in three major types: continuously
graded microstructure, discretely graded microstructure, and multi-phase graded
microstructure as pictured in Figure 1.3.
4
Analyzing the FGMs at the microstructure level with conventional Finite Element
Method (FEM) takes enormous pre-process and computational time. Essentially, the
model contains too many degrees of freedom to solve economically. Therefore, many
researchers have tried to analyze the FGMs using various methods, which increase the
accuracy of numerical solution and enhance the prediction of local stress
concentrations.
One of the common methods is to divide the FGMs domain into multi-layers in
the direction of the material gradation and to apply the traditional homogenization
method within each layer [27, 31] as shown in Figure 1.4. These layer-wise averaged
models are based on the self-consistent method [15]. In order to minimize the errors
in the layer-wise homogeneous model, Vemaganti and Deshmukh [34] used the
adaptive approach to model the FGMs. Some of the researchers, such as Sandra and
Lambros [29] varied the material property matrix in the finite element with a simple
boundary condition. Kim and Paulino [19] used an exponential variation of material
elastic properties to model nonhomogeneous, isotropic and orthotropic, materials.
Figure 1.4 Replacement scheme used in the layer-wise Homogenization
model in Continuously Graded Microstructure
5
Higher order theory was proposed by Aboudi, et. al. [1, 25]. This theory
allows the thermo-inelastic analysis of materials with spatially varying microstructure
based on volumetric averaging of the various field quantities and satisfaction of the
field equations. The use of the finite-element discretization approach utilizing
rectangular cells was compared with the averaging estimation methods in linear,
modified rules of mixtures (average Young’s modulus) and the Wakashima-
Tsukamoto estimate presented by Cho and Ha [7]. Also, Berezovski A. et al. [4] used
the linear rule of mixtures to define the material properties in his study of wave
propagation in FGMs. On the other hand, a discrete micromechanics approach, using
planar hexagonal cells of equal size, was developed by Ghosh et al. [10].
At the case of the MMCs or CMCs, layers-wise process also applied to replace
the heterogeneous microstructure properties to an equivalent continuum with a set of
macroscopic properties as shown in Figure 1.5. However, in this case no coupling
exists between local and global responses and these properties are calculated without
taking into account of the influence of the adjacent variable micro structural details
explicitly.
Modeling of material with microstructure has been carried out using classical
asymptotic homogeneous methods. In the classic methods, the microstructure is
periodic and is associated with a microstructure cell. However, the material
distribution of the FGMs, because of the definition of micro structural properties, is
arbitrary. The micro structure is not periodic in the length scale of the macrostructure.
6
Figure 1.5 Equivalent continuum with layers wise Homogenization in
Discretely Graded Microstructure with fiber and matrix configuration
(MMCs / CMCs types)
This leaves a major difficulty because very few analysis methods exist to solve a
structural problem which has a nonperiodic microstructure.
Therefore, in this proposal, a new theory of coupling the micro-macro
structural models is developed. The new method is capable of dealing with the
nonperiodic microstructure of the FGMs. A nonperiodic homogenization (NPH)
theory is created to analyze the FGMs models. This model will handle microstructure
with fiber/matrix combinations. The nonperiodic homogenization algorithm has been
developed to solve the problems with nonperiodic, arbitrarily spaced inclusions or
continuous fibers composite materials. The developed algorithm links the new
microstructural model with conventional Finite Element Method (FEM) discredited
technique.
Chapter 2 reviews the fundamental homogenization theory and defines the
homogenized elastic constant. Based on the microstructure cell solution and
7
homogenized elastic constant, 2-D cases with microstructure were solved, and the
solution was compare to analytical solutions. These results will be used for creating
the HOMOG algorithm. The algorithm was verified with the results of two different
examples, which were demonstrated from M. P. Bendsoe and N. Kikuchi [3].
Chapter 3 presents the nonperiodic homogenization theory and mathematical
formulations which are transferred from the Cartesian coordinate system to the natural
coordinate system. The details of conversion processes are shown in Appendix B.1.
Based on the mathematical formulations, the finite element stiffness matrix has been
defined and used in the creation of the nonperiodic homogenization algorithm. Many
variations of generic microcell structure were shown to fix frequently used FGMs
geometries.
Chapter 4 presents the results of 1-D cases of the nonperiodic homogenization
program. Most of the verifications are conducted in 1-D cases in order to compare
with the available analytical solutions. The model problems and associated boundary
conditions are described. The verification cases are: Case 1 – Comparison between
the NPH and the Homogenization Solution, Case 2 – Descending Low Density
Microcell Structure, Case 3 – Descending High Density Microcell Structure, Case 4 –
Descending and Ascending Microcell Structure, Case 5 – Descending Microcell
Structure with a Sudden Jump and Case 6 – Rapidly Varying Descending, Ascending
and Descending Microcell Structures.
Chapter 5 presents the 2-D verification cases for FGMs. The verification cases
are Case 7 – Periodic Microstructure, Case 8 – Descending Horizontal Fiber Strips in
8
One Direction, Case 9 – Descending and Ascending FGMs with Square Fibers, Case
10 – Descending and Symmetric Matrix Structure and Case 11 – Descending
Horizontal and Vertical Fiber Strips. Von-Mises stresses are presented for the Case 10
and Case 11.
Finally, Appendix A presents the nonperiodic homogenization coordinate
transformation from the Cartesian coordinate system to the natural coordinate system.
9
Chapter 2
Problems Involving a Microstructure
2.1 Review of Homogenization Theory
The following analysis represents the homogenization theory for the
periodic microstructural case. This is the starting point for the analysis of the
nonperiodic microstructure.
Consider the 2-D case. A solid body is made of two different materials and
has a base cell geometry of order (very small positive number) in size as shown
in Figure 2.1. Suppose the material properties in the base cell vary rapidly from
point to point producing heterogeneity. Thus, it is reasonable to assume that all
quantities have two explicit dependences. One is on the macroscopic level
x coordinate, and the other one is on the microscopic level coordinate / x . If g is
a general function, g = g( x , / x ) and / x can be replaced with ) / ( x y = . Due
to the periodic nature of the microstructure, the dependence of a function on the
micro-coordinate y is periodic. The quantity / 1 can be thought of as a
magnification factor, which enlarges the dimensions of a base cell, / x y = . Let
be an open subset of
2
with a smooth, boundary as described in Figure 2.2.
The Figure 2.2 depicts the associated macrostructure.
10
Figure 2.1 Periodic Microstructure Figure 2.2 Associated Macrostructure
Figure 2.3 Base cell of the Composite Microstructure
t
f

x
Y
y
S

Y
11
Let Y be a rectangular region in two dimensional spaces defined by
Y = ) , 0 (
0
1
y x ) , 0 (
0
2
y , (2.1)
Let be an open subset of Y with boundary
S = (2.2)
and let
/ Y = , (2.3)
where is the solid part of the cell with a material property E
1
, denotes the
closure of and Y represents the base cell of the composite microstructure. The
base cell properties vary inside Y , and the set represents a material property E
2
inside Y . Define now,

=
, if E
, if E
) (
2
1

y
y
y (2.4)
and extend to
2
by periodicity (i.e. it repeats the base cell in all two
direction). The superscript is the characteristic inhomogeneity dimension.
Find

V u , such that

+ + =

t
dS v p d v t d v f d
x
v
x
u
E
i i i i i i
j
i
l
k
ijkl

V v (2.5)

12
Here, it is assumed that the stress-strain and strain-displacement relations are

kl ijkl ij
e E = , (2.6)

+

=
k
l
l
k
kl
x
u
dx
u
e

2
1
, (2.7)
and that the elastic constants have the following properties:

klij ijlk jikl ijkl
E E E E = = = , (2.8)
Since the body forces f , tractions p and the elastic constants very with the small
cells of the composite (i.e., they are functions of both / and x y x = ). The general
functional form in the two scale approximation is,

/ ), , ( ) ( x y y x x = = , (2.9)
The displacement solution

u take this form; that is,

/ ), , ( ) ( x y y x x u = = u , (2.10)
Using a two scale asymptotic expansion,

/ ..., ) , ( ) , ( ) , ( ) (
2 2 1 0
x y y x u y x u y x u x u = + + + = (2.11)
where,
) , ( y x u
j
is defined in ) , ( × y x ,
) , ( y x u y
j
is Y- periodic.
Fact 1: The derivative of a periodic function is also periodic.
Fact 2: The integral of the derivative of a function over the period is zero.
Fact 3: If ) , ( y x = and y depends on x, then
13
x
y
y x dx
d

+

=

, / x y =
y x dx
d

+

=

1
(2.12)
Equation (2.5) becomes,

+

+

=

j
i
j
i
l
k
l
k
ijkl
j
i
l
k
ijkl
y
v
x
v
y
u
x
u
E d
x
v
x
u
E

1 1
0 0
()] +

+

+

+ d
y
v
x
v
y
u
x
u
j
i
j
i
l
k
l
k
...
1 1
2
1 1

(2.13)

+

+

+

=

j
i
l
k
j
i
l
k
j
i
l
k
j
i
l
k
ijkl
y
v
y
u
x
v
y
u
y
v
x
u
x
v
x
u
E
0
2
0 0 0
1 1 1

#
#
$
%
+

+

+

+

+ d
y
v
y
u
x
v
y
u
y
v
x
u
x
v
x
u
j
i
l
k
j
i
l
k
j
i
l
k
j
i
l
k
...
1
1 1 1 1

(2.14)
Now, rearrange for the terms in
2
1

,

1
and then,
&

&

#
#
$
%

+

+

+

=

j
i
l
k
j
i
l
k
l
k
j
i
l
k
ijkl
x
v
y
u
y
v
y
u
x
u
y
v
y
u
E
0 1 0 0
2
1 1

()
&
'
&
(
)
+
#
#
$
%

+

+

+

+ d
y
v
y
u
x
u
x
v
y
u
x
u
j
i
l
k
l
k
j
i
l
k
l
k
...
2 1 1 0
(2.15)
14
Equation (2.14) becomes,

&
&
'
&
&
(
)
&
&

&
&

+
#
#
$
%

+

+

+

+
#
#
$
%

+

+

+

d
y
v
y
u
x
u
x
v
y
u
x
u
x
v
y
u
y
v
y
u
x
u
y
v
y
u
E
j
i
l
k
l
k
j
i
l
k
l
k
j
i
l
k
j
i
l
k
l
k
j
i
l
k
ijkl

(...)
1 1
2 1 1 0
0 1 0 0
2

+ + =

t
S
i i i i i i
dS v p d v t d v f (2.16)
The functions are smooth enough so that the limit when 0 of all integrals exist
then, equation (2.16) holds if the terms of the same power of are equal to zero.
Therefore,

x
0
2
0
1

=

V d
y
v
y
u
E
j
i
l
k
ijkl
v (2.17)

x
0 1 0
1

=
#
#
$
%

+

+

V dS v p d
x
v
y
u
y
v
y
u
x
u
E
S
i i
j
i
l
k
j
i
l
k
l
k
ijkl
v (2.18)

#
#
$
%

+

+

+

d
y
v
y
u
x
u
x
v
y
u
x
u
E
j
i
l
k
l
k
j
i
l
k
l
k
ijkl

2 1 1 0

x

+ =

V d v t d v f
t
i i i i
v (2.19)
Fact 4: ()

*

*
+

,
1
lim
0
dYd
Y
d y
x
()

*

*
+

S S
dsd
Y
d ,
1
lim
0
y
x
15
v is an arbitrary function , thus we choose ) (y v + = then integrating by parts,
applying the divergence thermo to the integral in . Equations (2.17) thru (2.19)
becomes,

x
0
0
1

=

V d dY
y
v
y
u
E
Y
j
i
l
k
ijkl
v (2.20)

x
0 0
0
1

=
&
'
&
(
)
&

&

+
#
#
$
%

,

V d dS v n
y
u
E dY v
y
u
E
y Y
i j
S l
k
ijkl i
l
k
ijkl
i
v (2.21)
From the equation (2.21) becomes,
, y , 0
0
=

,
l
k
ijkl
i
y
u
E
y
(2.22)
S n
y
u
E
j
l
k
ijkl
on 0
0
=

(2.23)
Thus, this indicates that the first term of the expansion equation is depends only on
x , which is the macroscopic scale.
) ( ) , (
0 0
x u y x u = (2.24)
Now inset equation (2.24), ), (x u u
0 0
= into equation (2.18) and multiplying by
and using Fact (4). Since, ) (y v + = which means it’s derivative will be zero,
0 =

x
+
.

x
0 1 0
1

=
#
#
$
%

+

+

V dS v p d
x
v
y
u
y
v
y
u
x
u
E
S
i i
j
i
l
k
j
i
l
k
l
k
ijkl
v (2.25)
16
1
1 0

=
#
#
$
%

+

d dS v p
Y
d
y
v
y
u
x
u
E
S
i i
l
i
l
k
l
k
ijkl

(2.26)
( ) ( ) 1 1
1 0

=

+

d dS v p
Y
dYd
y
v
y
u
x
u
E
Y
S
i i
j
i
l
k
l
k
ijkl

y x
(2.27)
From the equation (2.27),
( ) ( )

=

+

S
i i
j
i
l
k
l
k
ijkl
dS v p dY
y
v
y
u
x
u
E

y x
1 0
is linear with respect to
0
u and
i
p .
Thus, if 0 =
i
p , then,
( ) ( ) ( )

=

,

dY
y
v
y
u
E dY
y
v
x
u
E
j
i
l
k
ijkl
j
i
l
k
ijkl
y y x
1 0
(2. 28)
or,
( ) ( ) ( )

, =

dY
y
v
x
u
E dY
y
v
y
u
E
j
i
l
k
ijkl
j
i
l
k
ijkl
y x y
0 1
(2.29)
Let
l
k kl
p p
x
u
u

, =
) (
) , (
0
1
x
y x (2.30)
and
l
k
q
kl
p
q
p
x
u
dy y
u

, =

) (
) , (
0
1
x
y x
substitute to equation (2.29) then,
( ) ( ) ( )

, =

,

dY
y
v
x
u
E dY
y
v
x
u
dy
E
j
i
l
k
ijkl
j
i
l
k
q
kl
p
ijkl
y x y x
y x
0 0
) (
) , (
(2.31)
17
Since, 1
) (
0
=

l
k
x
u x
, the equation (2.31) becomes,
( ) ( )

=

dY
y
v
E dY
y
v
dy
E
j
i
ijkl
j
i
q
kl
p
ijkl
y y
y x ) , (

V v (2.32)

, =

dY
y
E E
Y
E
m
kl
p
ijpm ijkl
H
ijkl
1
) (x (2.33)
H
ijkl
E defined by equation (2.33) represents the homogenized elastic constant.
18
2.2 Analytical Solution of the Homogenization in 2-D
For the 2-D problems, it would be sufficient to solve them for the cases;
Case A: k = l = 1
Case B: k = l = 2
Case C: k = 1, l = 2 (or k = 2, l = 1)
Case A: k = l = 1
Expand the equation (2.32) and remove terms with zero coefficients,
dY
y
v
y
E
y
E
y
E
y
E
y
v
y
E
y
E
y
E
y
E
y
v
y
E
y
E
y
E
y
E
y
v
y
E
y
E
y
E
y
E

&
&
&
&
&
&
'
&
&
&
&
&
&
(
)
&
&
&
&
&
&

&
&
&
&
&
&

+

+

+

+

+

+

+

+

+

+

+

+

+

+

+

2
2
2
11
2
2222
1
11
2
2221
2
11
1
2212
1
11
1
2211
1
2
2
11
2
2122
1
11
2
2121
2
11
1
2112
1
11
1
2111
2
1
2
11
2
1222
1
11
2
1221
2
11
1
1212
1
11
1
1211
1
1
2
11
2
1122
1
11
2
1121
2
11
1
1112
1
11
1
1111
) (
) (
) (
) (
y
y
y
y

+

+

=

dY
y
v
E
y
v
E
y
v
E
y
v
E
2
2
2211
1
2
2111
2
1
1211
1
1
1111
) ( ) (
) ( ) (
y y
y y
(2.34)
After the simplify the equation (2.34),
19
dY
y
v
y
v
y y
E
y
v
y
E
y
E
y
v
y
E
y
E

&
&
'
&
&
(
)
&
&

&
&

+

+

+

+

+

+

1
2
2
1
1
11
2
2
11
1
1212
2
2
2
11
2
2222
1
11
1
2211
1
1
2
11
2
1122
1
11
1
1111
) ( ) (
) ( ) (
y y
y y

+

=

dY
y
v
E
y
v
E
2
2
2211
1
1
1111
) ( ) ( y y
(2.35)
Define the cell solution as follows:
,
11
i i
= i = 1, 2 (2.36)
Define the components of the material stiffness matrix as follows:
1111 11
E D = (2.37)
2211 1122 21 12
E E D D = = = (2.38)
2222 22
E D = (2.39)
2121 2112 1221 1212 33
E E E E D = = = = (2.40)
Then, the material stiffness matrix is shown:
#
#
#
$
%

=
33
22 12
12 11
0 0
0
0
D
D D
D D
D (2.41)
And
{}
&
'
&
(
)
&

&

=
0
21
11
1
D
D
d ,{}
&
'
&
(
)
&

&

=
0
22
12
2
D
D
d and{}
&
'
&
(
)
&

&

=
33
3
0
0
D
d (2.42)
20
The equation (2.35) can be written as below:
[] []
'
(
)

#
#
#
#
#
#
#
$
%

#
#
#
#
$
%

2
1
1 2
2
1
1 2
2 1
2 1
0
0
0
0
y y
y
y
D
y y
y y
v v

[] dY D
D
y y
y y
v v
&
'
&
(
)
&

&

#
#
#
#
$
%

=

0
0
0
12
11
1 2
2 1
2 1

(2.43)
Further more the arbitrary function v and cell solution can be expressed with its
nodal shape function:
() [ ] {}
1 18 18 2
2
1
X
N
i X
N
v x
v
v
0 =
'
(
)

= v , N = 1, 9 (2.44)
() [ ] {}
1 18 18 2
2
1
X
N
i X
N
x =
'
(
)

= 0 , N = 1, 9 (2.45)
[] [][]
18 2 2 3 18 2
1 2
1
1
) ( ) ( 0
0
X
N
X X
N
x L x
y y
y
y
0 0 =
#
#
#
#
#
#
#
$
%

= B , N= 1, 9 (2.46)
21
Then, equation (2.43) can be formed as shown below:
[] [ ][ ]{} [] dY D
D
B dY B D B
T T T T
&
'
&
(
)
&

&

=

0
12
11

v v (2.47)
[ ] [ ][ ] { } [ ] { }dY d B dY B D B
T T
1

=

(2.48)
or
[ ] { } [ ] f K = for the cell solution (2.49)
where
[ ] [ ] [ ] [ ] B D B K
T
= ,
[ ] [ ] { }dY d B f
T
1

=

After obtained the cell solution from the equation (2.49), the homogenized elastic
constant defines by using equation (2.33).
dY
y
E
y
E
y
E
y
E E
Y
E
H

,

,

,

, =

2
11
2
1122
1
11
2
1121
2
11
1
1112
1
11
1
1111 1111 1111
1
(2.50)
where
0
1211
=
H
E
0
2111
=
H
E
22
Then equation (2.50)
H
E
1111
becomes as shown below:
dY
y
E
y
E E
Y
E D
H

,

, = =

2
11
2
1122
1
11
1
1111 1111 1111 11
1
dY
y
D
y
D D
Y

,

, =
2
2
12
1
1
11 11
1
[]

&
&
&
'
&
&
&
(
)
&
&
&

&
&
&

+

, =

dY
y y
y
y
D D D
Y
0
1
1
2
2
1
2
2
1
1
12 11 11
(2.51)
Then, equation (2.51) becomes:
[] () ( )

, =

1 dY d D
Y
D
1
T
1 11 11
(2.52)
where
( ) 1 =
&
&
&
'
&
&
&
(
)
&
&
&

&
&
&

+

1
2
2
1
2
2
1
1
y y
y
y
For the case of
2211
E from equation (2.33) becomes:
dY
y
E
y
E E
Y
E D
H

,

, = =

2
11
2
2222
1
11
1
2211 2211 2211 21
1
23
dY
y
D
y
D D
Y

,

, =
2
2
22
1
1
21 21
1
(2.53)
or
dY
y
D
y
D D
Y

,

, =
2
2
22
1
1
12 12
1
(2.54)
[] () [] ( )

, =

1 dY d D
Y
1
T
2 12
(2.55)
Therefore, from the equations (2.52) and (2.55),
[] () ( )

, =

1 dY d D
Y
D
H
1
T
1 11 11
(2. 56)
[] () ( )

, =

1 dY d D
Y
D
H
1
T
2 12 21
(2. 57)
where
( ) 1 =
&
&
&
'
&
&
&
(
)
&
&
&

&
&
&

+

1
2
2
1
2
2
1
1
y y
y
y
Equations (2.49), (2.56) and (2.57) will be used for developing the finite element
homogenization algorithm.
24
Case B: k = l = 2
Expand the equation (2.32) and removing terms with zero coefficients,
dY
y
v
y
v
y
E
y
E
y
v
y
E
y
E
y
v
y
E
y
E

&
&
'
&
&
(
)
&
&

&
&

+

+

+

+

+

+

1
2
2
1
1
22
2
1221
2
22
1
1212
2
2
2
22
2
2222
1
22
1
2211
1
1
2
22
2
1122
1
22
1
1111
) ( ) (
) ( ) (
y y
y y

+

=

dY
y
v
E
y
v
E
2
2
2222
2
1
1122
) ( ) ( y y
(2.58)
Define the cell solution as follows:
,
22
i i
* = i = 1, 2 (2.59)
Then, equation (2.58) can be form as below:
[ ] [ ][ ] { } [ ] { }dY d B dY B D B
T T
2

= *

(2.60)
or
[ ] { } [ ] f K = * for the cell solution (2.61)
where
[ ] [ ] [ ] [ ] B D B K
T
= ,
[ ] [ ] { }dY d B f
T
2

=

After obtained the cell solution from the equation (2.61), the homogenized elastic
constant defines by using equation (2.33).
25
dY
y
E
y
E E
Y
E
H

,

, =

2
22
2
1122
1
22
1
1111 1122 1122
1
(2.62)
where
0
1222
=
H
E
0
2122
=
H
E
dY
y
E
y
E E
Y
E
H

,

, =

2
22
2
2222
1
22
1
2211 2222 2222
1
(2.63)
Apply same process as it shown from previous case. Then equations (2.62) and
(2.63) become:
[] () ( ) dY d D
Y
E D
H H

* , = =

1
1
T
1 12 1122 12
(2.64)
[] () ( ) dY d D
Y
E D
H H

* , = =

1
1
T
2 22 2222 22
(2.65)
where
( ) * 1 =
&
&
&
'
&
&
&
(
)
&
&
&

&
&
&

*
+

*

*

*
1
2
2
1
2
2
1
1
y y
y
y
Equations (2.61), (2.64) and (2.65) will be used for developing the finite element
homogenization algorithm.
26
Case C: k = 1, l = 2 (or k = 2, l = 1)
dY
y
v
y
v
y
E
y
E
y
v
y
E
y
E
y
v
y
E
y
E

&
&
'
&
&
(
)
&
&

&
&

+

+

+

+

+

+

1
2
2
1
1
12
2
1221
2
12
1
1212
2
2
2
12
2
2222
1
12
1
2211
1
1
2
12
2
1122
1
12
1
1111
) ( ) (
) ( ) (
y y
y y

+

=

dY
y
v
y
v
E
1
2
2
1
1212
) ( ) ( y y
(2.66)
Define the cell solutions as follows:
,
12
i i
= i = 1, 2 (2.67)
Then, equation (2.66) can be form as below:
[ ] [ ][ ] { } [ ] { }dY d B dY B D B
T T
3

=

(2.68)
or
[ ] { } [ ] f K = for the cell solution (2.69)
where
[ ] [ ] [ ] [ ] B D B K
T
= ,
[ ] [ ] { }dY d B f
T
3

=

After obtained the cell solution from the equation (2.69), the homogenized elastic
constant defines by using equation (2.33).
0
1121 1112
= =
H H
E E
dY
y
E
y
E E
Y
E E
H H

,

, = =

1
12
2
1221
2
12
1
1212 1212 1221 1212
1
(2.70)
27
dY
y
E
y
E E
Y
E E
H H

,

, = =

1
12
2
2121
2
12
1
2112 2112 2121 2112
1
(2.71)
0
2221 2212
= =
H H
E E
From the results of equations (2.70) and (2.71), it shows that
H H H H
E E E E
2121 2112 1221 1212
= = = (2.72)
Then, equations (2.70) and (2.71) become,
[] () ( ) dY d D
Y
D
H

, =

1
1
T
3 33 33
(2.73)
where
( ) 1 =
&
&
&
'
&
&
&
(
)
&
&
&

&
&
&

+

1
2
2
1
2
2
1
1
y y
y
y
Equations (2.69) and (2.73) will be used for developing the finite element
homogenization algorithm.
Finally, the homogenized material stiffness matrix
H
D is defined and
shown below. It will use in conventional finite element analysis process:
[]
#
#
#
$
%

=
H
H H
H H
H
D
D D
D D
33
22 12
12 11
0 0
0
0
D (2.74)
28
2.3 Macro Homogenized Solution
The conventional stiffness matrix for 2-D plane stress element is shown below:
[ ] [ ] [ ][ ]

= tdxdy k
T
B E B (2.75)
where
[] [] [][]
18 2 2 3 18 2
1 2
1
1
) ( ) ( 0
0
X
N
X X
N
x L x
y y
y
y
0 0 =
#
#
#
#
#
#
#
$
%

= B , N = 1, 9
[ ] = E Material stiffness matrix
t = Material thickness
Convert the equation (2.75) from Cartesian coordinate system dxdy to
Isoparametric coordinate systems 2 3d d . The sides of the element limitation are at
1 ± = 3 and at 1 ± = 2 .
[ ] [] [ ][] [ ] [ ][ ] 2 3d d t tdxdy k
T T
J B E B B E B
1
1
1
1

,,
= = (2.76)
where J is called the Jacobian matrix which is defined as below:
#
#
$
%

=
5 5
5 5
i
N
i i
N
i
i
N
i i
N
i
y x
y x
2 2
3 3
0 0
0 0
, ,
, ,
J (2.77)
29
Finally, replace the stiffness matrix [E] with [D
H
]. Then the equation (2.76)
becomes:
[] [ ][][] 2 3d d t k
H T H
J B D B
1
1
1
1

,,
= (2.78)
30
2.4 Variation of Cell Solution
kl
i

in Microcell Geometries
2.4.1 A Typical Shape of Microcell Geometry
Based on the equations (2.32) and (2.33), homogenization program was
created and tested. Two cell geometries, as pictured in Figure 2.4 were defined.
The microcell structures have soft and hard materials, E
soft
and E
hard
. In the test
cases, it was assumed that E
soft
= 10, E
hard
= 1000 and Poisson ratio v = 0.3. The
shaded area represents E
hard
which is a fiber and the brightened area represents E
soft
which is a matrix. The microcell models are pictured in Figure 2.4. The fiber
width hx
1
and hy
1
= 0.250 for both Model 1 and 2. The matrix lengths in Model 1
were hx
2
= hy
2
= 0.375 and Model 2 were hx
2
= 0.125 and hy
2
= 0.375. The cell
solution
kl
i
was computed for the two cell geometry cases.
Figure 2.4 Microcell Structure: Model 1 - Cell Size 1.0 x 1.0; Model 2 - Cell
size 0.5 x 1.0
Y
1
Y
2
hx
2
hy
2
hy
1
hx
1
3
1
3
2
Model 1
Y
1
Y
2
hx
2
hy
2
hy
1
hx
1
3
1
3
2
Model 2
31
2.4.2 Computation of the Cell Solution (
kl
i

) in Soft and Hard Microcell
Structure Cases
Case A: k=l=1
Figure 2.5 and Figure 2.6 showed the results of the microcell solution for
Model 1 and Model 2 after the homogenization program was ran. Model 1 has the
cell size as 1.0 x 1.0 and Model 2 has 0.5 x 1.0.
Micro Cell Solution
(Case A: k = l = 1)
0.00
0.25
0.50
0.75
1.00
0.00 0.25 0.50 0.75 1.00
Y1
Y2
Model 1 - 1.0 x 1.0
Figure 2.5 Microcell Solution for Model 1 – Cell size 1.0 x 1.0
32
Micro Cell Solution
(Case A: k = l = 1)
0.00
0.25
0.50
0.75
1.00
0.00 0.25 0.50 0.75 1.00
Y1
Y2
Model 2 - 0.5 x 1.0
Figure 2.6 Microcell Solution for Model 2 – Cell size 0.5 x 1.0
The pictures of Figure 2.5 and Figure 2.6 were superimposed in order to
demonstrate that the microcell solutions are different if the microstructural shapes
are different (see Figure 2.7).
Micro Cell Solution
(Case A: k = l = 1)
0.00
0.25
0.50
0.75
1.00
0.00 0.25 0.50 0.75 1.00
Y1
Y2
Model 1 - 1.0 x 1.0
Model 2 - 0.5 x 1.0
Figure 2.7 Superimpose the Cell Solution Pictures of Figure 2.5 and Figure 2.6
33
Case B: k = l =2
Figure 2.8 and Figure 2.9 showed the results of the microcell solution for
Model 1 and Model 2 after the homogenization program was ran.
Micro Cell Solution
(Case B: k = l = 2)
0.00
0.25
0.50
0.75
1.00
0.00 0.25 0.50 0.75 1.00
Y1
Y2
Model 1 - 1.0 x 1.0
Figure 2.8 Microcell Solution for Model 1 – Cell size 1.0 x 1.0
Micro Cell Solution
(Case B: k = l = 2)
0.00
0.25
0.50
0.75
1.00
0.00 0.25 0.50 0.75 1.00
Y1
Y2
Model 2 - 0.5 x 1.0
Figure 2.9 Microcell Solution for Model 2 – Cell size 0.5 x 1.0
34
Again, the pictures of Figure 2.8 and Figure 2.9 were superimposed in
order to demonstrate that the microcell solutions are different if the microstructural
shapes are different (see Figure 2.10)
Micro Cell Solution
(Case B: k = l = 2)
0.00
0.25
0.50
0.75
1.00
0.00 0.25 0.50 0.75 1.00
Y1
Y2
Model 1 - 1.0 x 1.0
Model 2 - 0.5 x 1.0
Figure 2.10 Superimpose the Cell Solution Pictures of Figure 2.8 and Figure 2.9
Case C: k = 1, l = 2 (or k = 2, l = 1)
Figure 2.11 and Figure 2.12 showed the results of the microcell solution for
Model 1 and Model 2 after the homogenization program was ran.
35
Micro Cell Solution
(Case C: k = 1, l = 2 or k = 2, l = 1)
-0.20
0.00
0.20
0.40
0.60
0.80
1.00
1.20
-0.20 0.00 0.20 0.40 0.60 0.80 1.00 1.20
Y1
Y2
Model 1 - 1.0 x 1.0
Figure 2.11 Microcell Solution for Model 1 – Cell Size 1.0 x 1.0
Micro Cell Solution
(Case C: k = 1, l = 2 or k = 2, l = 1)
-0.20
0.00
0.20
0.40
0.60
0.80
1.00
1.20
-0.20 0.00 0.20 0.40 0.60 0.80 1.00 1.20
Y1
Y2
Model 2 - 0.5 x 1.0
Figure 2.12 Microcell Solution for Model 2 – Cell Size 0.5 x 1.0
36
Finally, the pictures of Figure 2.11 and Figure 2.12 were superimposed in
order to demonstrate that the microcell solutions are different if the microstructural
shapes are different (see Figure 2.13).
Micro Cell Solution
Case C: k = 1, l = 2 or k = 2, l = 1)
-0.20
0.00
0.20
0.40
0.60
0.80
1.00
1.20
-0.20 0.00 0.20 0.40 0.60 0.80 1.00 1.20
Y1
Y2
Model 1 - 1.0 x 1.0
Model 2 - 0.5 x 1.0
Figure 2.13 Superimpose the Cell Solution Pictures of Figure 2.11 and Figure 2.12
37
2.5 Verification of the Homogenized Elastic Constants
2.5.1 Soft and Hard Isotropic Composite Material with Variation
A finite element cell solution program HOMOG.f90, was created to obtain
approximate cell solutions. The homogenized elastic constant,
H
ijkl
E , was
computed and compared with the results from Bendoe and Kikuchi [3]. For the
first case, the cell had the soft and hard material properties with isotropic and plane
stress. And the soft material’s Young’s modulus E
soft
= 10 and E
hard
= 1000 with
the same Poisson ratio 3 . 0 = + . The structure is pictured in Figure 2.14.
Figure 2.14 Material Properties in a Single Microcell Geometry:
(E
soft
= 10 and E
hard
= 1000)
E
soft
E
hard
E
soft
E
hard
E
soft
E
hard
E
soft
E
hard
E
soft
38
Bendoe and Kikuchi [3] calculated the homogenized elastic constant with a
16x16 mesh of elements on the cell picture in Figure 2.15. They also used an
adaptive method to obtain refined the results.
HOMOG.f90 used an 8 x 8 mesh of elements of the cell geometry which is
shown in Figure 2.16 for calculating the homogenized elastic constant. A nine-
node Lagrange quadratic element and sixteen Gauss points were used to increase
accuracy. The comparison results are shown in Table 2.1.
Figure 2.15 A Single Microcell Geometry Model Defined by
Bendoe and Kikuchi [3]: 16 x 16 Uniform Square Four-node
Isoparametric Elements
39
Figure 2.16 A Single Microcell Geometry Model Defined by
HOMOG.f90: 8 x 8 Uniform Square Nine-node Lagrange
Quadratic Elements
Table 2.1 Comparison of the Homogenized Elastic Constants: Soft and Hard
Material
Mesh or Nodes
H
E
1111
H
E
1122
H
E
2222
H
E
1212
16 x 16
[3]
149.80 71.61 149.80 87.12
1
st
Adapt
[3]
127.12 62.91 127.12 75.90
2
nd
Adapt
[3]
125.79 62.62 125.79 75.28
HOMOG.f90
(9-node)
136.55 68.56 136.55 81.07
40
2.5.2 A Rectangular Hole in the Material
For the second case verification, Bendoe and Kikuchi [3] used rectangular
microcell geometries with a rectangular hole in the middle. The structure is
pictured in Figure 2.17 and has material properties which are characterized by
1111
E = 30
2222
= E and 10
1212 1122
= = E E . In this case, the hole was not
approximated by a very soft material. The material properties of the microcell
only consider the solid area.
Figure 2.17 Material Properties in a Single Microcell Geometry
with a Hole: (E
1111
= E
2222
= 30 and E
1122
= E
1212
= 10)
In Figure 2.18, a 20 x 20 mesh size was used for the cell geometry and
three different adaptation mesh methods were generated to obtain refined the
results of homogenized elastic constant.
41
HOMOG.f90 used a 10 x 10 mesh of elements of the cell geometry which
is shown in Figure 2.19 for calculating the homogenized elastic constant. As same
as the first case, a nine-node Lagrange quadratic element and sixteen Gauss points
were used to increase accuracy. The comparison results are shown in Table 2.2.
Figure 2.18 A Single Microcell Geometry Model Defined by
Bendoe and Kikuchi [3]: 20 x 20 Uniform Square Four-node
Isoparametric Element
Figure 2.19 A Single Microcell Geometry Model Defined by
HOMOG.f90: 10 x 10 Uniform Square Nine-node Lagrange
Quadratic Elements
42
Table 2.2 Elastic Tensors after Homogenization in Rectangular Hole
Mesh or Nodes
H
E
1111
H
E
1122
H
E
2222
H
E
1212
16 x 16
[3]
13.015 3.241 17.552 2.785
1
st
Adapt
[3]
12.910 3.178 17.473 2.714
2
nd
Adapt
[3]
12.865 3.146 17.437 2.683
3
rd
Adapt
[3]
12.844 3.131 17.421 2.668
HOMOG.f90
(9-node)
12.924 3.198 17.487 2.708
43
Chapter 3
Nonperiodic Homogenization (NPH) Method
3.1 Mathematical Theory of NPH
Assume that the microscopic and the macroscopic have a relationship
based on equations (2.24) and (2.30). Then an equation can be written as below:
) (
) (
) , ( ) ( ) , (
) (
x
x
y x x y x i
l
H
k kl
i i i
C
X
u
u u

+ = (3.1)
Macro Structure Micro Structure
Solution Solution
where, ) , ( y x
kl
i
= Microcell solution (i = 1, 2)
) (x i C = Correction coefficient related to the micro structure
) (x
H
k
u = Homogenized displacement solution (k = 1, 2)
l
X = Macro coordinate system (i = 1, 2)
Also, ) (x
i
u can be expressed as function of the homogenized displacement with
correction coefficient ) (x
i
C and ) (x i C through the expression below:
) ( ) ( ) ( ) (
) (
x x x x i
i
H
i i
C C u u + = (3.2)
Then, insert in the equation (3.2) into equation (3.1) then can be rewritten as
below:
44
) (
) (
) , ( ) ( ) ( ) ( ) , (
) ( ) (
x
x
y x x x x y x i
l
H
k kl
i
i
i
H
i i
C
X
u
C C u u

+ + = (3.3)
Now, the nodeless form of the equation can be derived. Let at
N
x be a nodal
position in the macro element. Let
c
y be a single location in the microcell.
Evaluating equation (3.3) at
N
x and
c
y produces the follow equations:
) (
) (
) , ( ) ( ) ( ) ( | ) , (
) ( ) ( N
i
l
N
H
k
C N
kl
i N
i
N i N
H
i i
C
X
u
C C u u
C
N
x
x
y x x x x y x
y y
x x

+ + =
=
=
(3.4)
Therefore, displacement components
i
u at nod N can be formed as shown below:
N
i
l
N
H
k
C N
kl
i
N
i
N
i N
H
i
N
i
C
X
u
C C u u
) (
) , ( ) (
) ( ) (

+ + =
x
y x x (3.5)
Then, correction coefficient component
N
i C , can be defined as a function of the
displacement associates with the cell solution as follows:
N
i
l
N
H
k
C N
kl
i
N
i N
H
i
N
i
N
i C
X
u
C u u C
) (
) , ( ) (
) ( ) (

, , =
x
y x x (3.6)
Now introduce equation (3.6) into the equation (3.4) and rearrange the resulting
expression in terms of the correction coefficients ) (x
i
C and ) (x i C .
) (
) (
) , ( ) ( ) ( ) , (
) ( ) (
x
x
y x x x y x i
l
H
k kl
i i
H
i i
C
X
u
C u u

+ =
]
) (
) , ( ) ( )[ (
) ( ) (
N
i
l
N
H
k
C N
kl
i
N
i N
H
i
N
i
N
C
X
u
C u u

, , +
x
y x x x 0 (3.7)
45
where
N
i
N
i
C C ) ( ) ( x x 0 =
N
i
N
i C C ) ( ) ( x x 0 =
N
i
N
i
u u ) ( ) ( x x 0 = (3.8)
Figure 3.1 indicates the nodal degree-of-freedom for the 1-D element case. Each
nodal point has three degree-of-freedom
i
C , i C and
i
u .
Figure 3.1 Nodal Degree of Freedom in 1-D Element
Introduce equation (3.8) into equation (3.7) to obtain
N
i
N
l
H
k kl
i
N
i
N H
i i
C
X
u
C u u ) (
) (
) , ( ) ( ) ( ) , (
) ( ) (
x
x
y x x x y x 0 0

+ =
#
$
%

, , +
N
i
l
N
H
k
C N
kl
i
N
i N
H
i
N
i
N
C
X
u
C u u
) (
) , ( ) ( ) (
) ( ) (
x
y x x x 0 (3.9)
Finally, the displacement vector ) , ( y x
i
u can be formed as shown below:
{ }
N
i
N
N
H
i
N H
i i
C u u u ) ( ) ( ) ( ) ( ) , (
) ( ) (
x x x x y x 0 0 , =
N
i
N
l
N
H
k
C N
kl
i
N
l
H
k kl
i
C
X
u
X
u
) (
) (
) , ( ) (
) (
) , (
) ( ) (
'
(
)

,

+ x
x
y x x
x
y x 0 0
C
1
, C
1,
u
1
1 2
C
2
, C
2,
u
2
46
N
i
N
u ) (x 0 + (3.10)
or
N
i
N
N
i
N N
i
N
i
u C C u ) ( ) ( ) ( ) , ( x x x y x 0 7 1 + + = (3.11)
where, ) (x
N
1 and ) (x
N
7 are nodeless shape functions:
) ( ) ( ) ( ) (
) ( ) (
x x x x
N
N
H
i
N H
i
N
u u 0 0 1 , =
N
7 = ) (
) (
) , ( ) (
) (
) , (
) ( ) (
x
x
y x x
x
y x
N
l
N
H
k
C N
kl
i
N
l
H
k kl
i
X
u
X
u
0 0

,

47
3.2 Evaluation of the NPH Strain Tensors
The two dimensional strains in Cartesian coordinates are defined as follows:
1
1
1
1
1
1
11
) , ( ) , (
X
Y
Y
u
X
u
e

+

=
y x y x
(3.12)
2
2
2
2
2
2
22
) , ( ) , (
X
Y
Y
u
X
u
e

+

=
y x y x
(3.13)
) , ( ) , ( ) , ( ) , (
1
1
1
2
1
2
2
2
2
1
2
1
12
X
Y
Y
u
X
u
X
Y
Y
u
X
u
e

+

+

+

=
y x y x y x y x
(3.14)
Equation (3.11) through equation (3.13) can be rearranged into matrix format as
follows:
&
&
&
'
&
&
&
(
)
&
&
&

&
&
&

+

+

+

+

+

=
&
'
&
(
)
&

&

1
1
1
2
1
2
2
2
2
1
2
1
2
2
2
2
2
2
1
1
1
1
1
1
12
22
11
) , ( ) , ( ) , ( ) , (
) , ( ) , (
) , ( ) , (
X
Y
Y
u
X
u
X
Y
Y
u
X
u
X
Y
Y
u
X
u
X
Y
Y
u
X
u
e
e
e
y x y x y x y x
y x y x
y x y x
#
#
#
#
#
#
#
#
#
#
$
%

&
&
&
&
'
&
&
&
&
(
)
&
&
&
&

&
&
&
&

+
&
&
&
&
'
&
&
&
&
(
)
&
&
&
&

&
&
&
&

#
#
#
$
%

=
2 2
2
1 1
2
2 2
1
1 1
1
2
2
1
2
2
1
1
1
1 ) , (
1 ) , (
1 ) , (
1 ) , (
) , (
) , (
) , (
) , (
x
0 1 1 0
1 0 0 0
0 0 0 1

Y
u
Y
u
Y
u
Y
u
X
u
X
u
X
u
X
u
y x
y x
y x
y x
y x
y x
y x
y x
(3.15)
where,
2
2
2 1
1
1
1
and
1
X
Y
X
Y

=

=

(3.16)
48
Now, substitute equation (3.10) into equation (3.15). Some of the terms
will become zero after differentiation. For example,
1
) ( 1
) (
X
u
N
H

x
is a constant
relative to variations in X
1
and X
2
. Set value of cell solution at
c
y to zero (i.e.
0 ) , (
) ( 1
=
c N
kl
y x ). Rearrange the equation with respect to common coefficient
terms, then the equation (3.15) becomes as shown below:

&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

,
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

#
#
#
$
%

=
&
'
&
(
)
&

&

N
N
N
H
N
N
N
H
N
N
N
H
N
N
N
H
N
N
H
N
N
H
N
N
H
N
N
H
C
X
u
C
X
u
C
X
u
C
X
u
C
X
u
C
X
u
C
X
u
C
X
u
e
e
e
2
2
) ( 2
2
1
) ( 2
1
2
) ( 1
1
1
) ( 1
2
2
2
2
1
2
1
2
1
1
1
1
12
22
11
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
0 1 1 0
1 0 0 0
0 0 0 1
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
0
0
0
0
0
0
0
0
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

,

,

,

,

+
N N
N
H
N N
H
N N
N
H
N N
H
N N
N
H
N N
H
N N
N
H
N N
H
C
X
u
C
X
u
C
X
u
C
X
u
C
X
u
C
X
u
C
X
u
C
X
u
2
2
) ( 2
2
2
2
2
1
) ( 2
2
1
2
1
2
) ( 1
1
2
1
1
1
) ( 1
1
1
1
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
0 0
0 0
0 0
0 0
49
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

,

,

,

,

+
N
N
l
H
k
c N
kl N
N
l
H
k kl
N
N
l
H
k
c N
kl N
N
l
H
k kl
N
N
l
H
k
c N
kl N
N
l
H
k kl
N
N
l
H
k
c N
kl N
N
l
H
k kl
C
X X
u
C
X X
u
C
X X
u
C
X X
u
C
X X
u
C
X X
u
C
X X
u
C
X X
u
2
2
) ( 2 2
2
2
2
1
) ( 2 2
1
2
1
2
) ( 1 1
2
1
1
1
) ( 1 1
1
1
) ( ) (
) , (
) ( ) (
) , (
) ( ) (
) , (
) ( ) (
) , (
) ( ) (
) , (
) ( ) (
) , (
) ( ) (
) , (
) ( ) (
) , (
x x
y x
x x
y x
x x
y x
x x
y x
x x
y x
x x
y x
x x
y x
x x
y x
0

0

0

0

0

0

0

0

&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

+
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

+
N N
l
H
k kl
N N
l
H
k kl
N N
l
H
k kl
N N
l
H
k kl
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
C
X X
u
C
X X
u
C
X X
u
C
X X
u
C
X
u
X
C
X
u
X
C
X
u
X
C
X
u
X
2
2
2
2
1
2
1
2
1
1
1
1
2
2
2
2
1
2
1
2
1
1
1
1
) (
) (
) , (
) (
) (
) , (
) (
) (
) , (
) (
) (
) , (
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
x
x
y x
x
x
y x
x
x
y x
x
x
y x
x
x y x
x
x y x
x
x y x
x
x y x
0
0
0
0
0

0

0

0

) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
) (
) (
) (
) (
2
2
2
2
2
1
1
2
1
2
2
2
1
2
1
1
1
1
1
1
2
2
2
1
1
2
1
1
#
#
#
#
#
#
#
#
#
#
$
%
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

+
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

+
X
Y
C
X
u
Y
X
Y
C
X
u
Y
X
Y
C
X
u
Y
X
Y
C
X
u
Y
u
X
u
X
u
X
u
X
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
N
N
N
N
N
N
N
N
x
x y x
x
x y x
x
x y x
x
x y x
x
x
x
x
0

0

0

0

0
0
0
0
(3.17)
Let[]
#
#
#
$
%

=
0 1 1 0
1 0 0 0
0 0 0 1
A .
50
Then the equation (3.17) can be written as shown below:
[ ][ ] [ ][ ] [ ] [ ][ ] [ ][ ]
8 7 6 5 4 3 2 1
12
22
11
x x x x x x x x B B B B B B B B A
e
e
e
=
&
'
&
(
)
&

&

(3.18)
where
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

=
N
N
H
N
N
H
N
N
H
N
N
H
C
X
u
C
X
u
C
X
u
C
X
u
B
2
2
2
2
1
2
1
2
1
1
1
1
1
) (
) (
) (
) (
) (
) (
) (
) (
x
x
x
x
x
x
x
x
0
0
0
0
(3.19)
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

, =
N
N
N
H
N
N
N
H
N
N
N
H
N
N
N
H
C
X
u
C
X
u
C
X
u
C
X
u
B
2
2
) ( 2
2
1
) ( 2
1
2
) ( 1
1
1
) ( 1
2
) (
) (
) (
) (
) (
) (
) (
) (
x
x
x
x
x
x
x
x
0
0
0
0
(3.20)
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

=
N N
H
N N
H
N N
H
N N
H
C
X
u
C
X
u
C
X
u
C
X
u
B
2
2
2
2
1
2
1
2
1
1
1
1
3
) (
) (
) (
) (
) (
) (
) (
) (
x
x
x
x
x
x
x
x
0
0
0
0
(3.21)
51
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

=
N
N
l
H
k kl
N
N
l
H
k kl
N
N
l
H
k kl
N
N
l
H
k kl
C
X X
u
C
X X
u
C
X X
u
C
X X
u
B
2
2
2
2
1
2
1
2
1
1
1
1
4
) ( ) (
) , (
) ( ) (
) , (
) ( ) (
) , (
) ( ) (
) , (
x x
y x
x x
y x
x x
y x
x x
y x
0

0

0

0

(3.22)
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

=
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
C
X
u
X
C
X
u
X
C
X
u
X
C
X
u
X
B
2
2
2
2
1
2
1
2
1
1
1
1
5
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
x
x y x
x
x y x
x
x y x
x
x y x
0

0

0

0

(3.23)
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

=
N N
l
H
k kl
N N
l
H
k kl
N N
l
H
k kl
N N
l
H
k kl
C
X X
u
C
X X
u
C
X X
u
C
X X
u
B
2
2
2
2
2
1
2
2
1
2
2
1
1
1
2
1
6
) (
) (
) , (
) (
) (
) , (
) (
) (
) , (
) (
) (
) , (
x
x
y x
x
x
y x
x
x
y x
x
x
y x
0
0
0
0
(3.24)
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

=
N
N
N
N
N
N
N
N
u
X
u
X
u
X
u
X
B
2
2
2
1
1
2
1
1
7
) (
) (
) (
) (
x
x
x
x
0
0
0
0
(3.25)
52
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

=
2
2
2
2
2
1
1
2
1
2
2
2
1
2
1
1
1
1
1
1
8
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
X
Y
C
X
u
Y
X
Y
C
X
u
Y
X
Y
C
X
u
Y
X
Y
C
X
u
Y
B
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
x
x y x
x
x y x
x
x y x
x
x y x
0

0

0

0

(3.26)

53
3.3 Explicit Parameterization of the Cell Geometry over the Element
The equation (3.23) includes the term
l
H
k
i
kl
X
u
X

) ( ) , (
1
x y x
, which is related
to the rate change in the microcell solution with respect to
j
X . Changes in
j
X are
dependent on changes in the parameters, ) , (
2
y x hx and ) , (
2
y x hy , controlling the
cell size. Essentially ) , (
2
y x hx and ) , (
2
y x hy define an explicit parameterization of
the cell size in terms of macro element positions. The microcell
approximation ) , ( y x
kl
i
, is shown below:
) ( ) , ( ) , ( x y x y x
N kl
i
N kl
i
8 = (3.27)
Then evaluating equation (3.27) at all Gauss point ij, it can be seen that:
ij ij
ij
x
N kl
i
x
N
x
kl
i
= =
=
=
x x
x
x y x y x ) ( ) , ( ) , ( 8
) ( ) (
ij N kl
i
ij N
x y 8 =
) (y
ij kl
i
= (3.28)
Thus,
[] ) ( ) (
) ( ) , (
1 1 ij N kl
i
ij N
i
x
i
ij kl
x
i
kl
X X X
ij ij
x y
y y x
x x
8

=

=

= =
(3.29)
Or
i
ij N kl
ij N ij N kl
i
i
ij N
x
i
kl
X X X
ij

+

=

=
) (
) ( ) (
) ( ) , (
1 1
x
y x
y y x
x

8
8
54
) (
) (
) (
) ( ) (
) (
) (
2
2
2
2
x
x
x
y x
x
y
N kl
i
i
ij
ij
ij N
i
ij
ij
ij N
X
hy
hy X
hx
hx

8 8
#
$
%

+

=
#
$
%

+

+
i
ij
ij
N kl
i
ij
ij
N kl
ij N
X
hy
hy X
hx
hx
) (
) (
) ( ) (
) (
) (
) (
2
2
1 2
2
1
x
x
x x
x
x
y

8 (3.30)
See the Appendix B.1 (equations B.14 through B.31) for detail processing
in the equation (3.23). The rate of changes in the cell solution with respect to
) (
) (
2
1
ij
N kl
hx x
x

and
) (
) (
2
1
ij
N kl
hy x
x

can be approximated using finite differences. For
example, if k = l = 1,
2
) (
11
1 ) (
11
1
2
11
1 2 2 2
) ( ) (
) (
) (
hx hx
hx
N
hx hx
N
N
9
,
=

9 + x x
x x
x
x

(3.31)
2
) (
11
1 ) (
11
1
2
11
1 2 2 2
) ( ) (
) (
) (
hy hy
hy
N
hy hy
N
N
9
,
=

9 + x x
x x
x
x

(3.32)
55
3.4 Coordinate Transformation
The equation (3.19) through equation (3.26) will be converted from the
Cartesian coordinate system to the natural coordinate system (see Appendix B.1
for the detail presentation of the transformation process). Both the macro
coordinate x and the micro coordinate y must be transformed. Thus the typical
function ) (x f in macro coordinate system and ) (y g in micro coordinate system
are changed to ) (3 f and ) (3 g respectively where 3 is the natural coordinate for
the macro system and 3 is the natural coordinate for the micro system. After each
transformation has been completed, the following equations can be formed:
Transforming equation (3.19) - (3.26) gives,
[ ] [ ] [ ] { }
1 18 18 4 4 4 1
) ( ) (
x
N
i x X
H
N
C DJBX u B 3 3 = (3.33)

[ ] [ ] [ ] { }
1 18 18 18
N
18 4 2
) (
x
N
i X
H
i x
C u DJBX B 3 = (3.34)

[ ] [ ] [ ] [ ] { }
1 18 18 2 2 18
N
18 4 3
) ( ) (
X
N
i x
N
X
H
i x
C u DJBX B 3 0 3 = (3.35)

[] [ ] [][ ]{} [] { }
1 18 18 4
4 4
1 18
3 18
11
18 2 4
) ) ( ) ( (
x
N
i
X
x
x
H
i
x
N
i x
C DJB u DJBX A B
ij
#
$
%

#
$
%

= 3 3 0
(3.36)
56
[]
[]
[]
[] {} {}
1 18
18 2
3 36
1
2
2
N
1
36 4
3 36
1
2
2
N
1
36 4
3 18
18 4
1
2
2
1
3 18
18 4
1
2
2
1
5
) ( )] ( [
) (
) (
) (
) , (
) (
) (
) (
) , (
) (
) (
) (
) (
) (
) (
x
N
i
x
N H
i
X
ij
ij
ij
kl
X
ij N
X
ij
ij
ij
kl
X
ij N
X
ij
N kl
i
X
ij
ij
ij
X
ij
N kl
i
X
ij
ij
ij
C u DJBX A
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
B
#
$
%

&
&
&
&
&
&
&
&
'
&
&
&
&
&
&
&
&
(
)
&
&
&
&
&
&
&
&

&
&
&
&
&
&
&
&

#
#
$
%

+
#
#
$
%

+
#
$
%

#
#
$
%

+
#
$
%

#
#
$
%

= 3 0 3
3
3
3
2 3 8
3
3
3
2 3 8
3
3
3
8
3
3
3
8
(3.37)

[]
[] [][] () {}
[] [][] () {}
{}
1 18 18 2
2 4
1 18
18 2
2 4
2
22
18 2
2 4
2
21
4 4
1
3 18
11
18 2
1 18
18 2
2 4
2
12
18 2
2 4
2
11
4 4
1
3 18
11
18 2
6
) (
) (
) (
) (
) (
) (
) (
X
N
i
X
X
X
H
i
x
N
X
i j
x
N
X
j i
X
x
N
i x
X
H
i
x
N
X
i j
x
N
X
j i
X
x
N
i x
C
u J A
u J A
B
ij
ij
3 0
3 0
3 3
3 0
3 3
3 0
3 0
3 3
3 0
3 3
3 0
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
$
%

#
$
%

+
#
$
%

#
$
%

#
$
%

+
#
$
%

#
$
%

=
,
,
(3.38)

[ ] [ ] { }
1 18 18 4 7 X
N
i X
u DJBX B = (3.39)

[] [ ][] ()[][ ] {} {} { }
1 18 18 2
1 4
1 18 18 4 3 18 1 18 4 8
) ( ) (
X
N
i
X
X
i
i
X
H
i X X
ij N kl
X
C
X
Y
u DJBX A DJBY B 3 0
#
#
$
%

'
(
)

= x
(3.40)
57
Introducing equation (3.33) - (3.40), the equation (3.18) can be evaluated.
Arranging the equation (3.18) in blocks corresponding to the correction
coefficients (C and C ) and displacement u, the following equation is obtain:
[] {}
1 54 54 3
12
22
11
x
x x
UC P
e
e
e
=
&
'
&
(
)
&

&

(3.41)
where,
[ ] [ ] [ ] [ ] [ ]
54 3 54 3
coeff. of term | coeff. of term | coeff. of term
x x
P u C C =
(3.42)
{} [ ]
T
N N
N N
N N
x
u u u u C C C C C C C C UC
2 1
1
2
1
1
2 1
1
2
1
1
2 1
1
2
1
1 1 54
... | ... | ... = (3.43)

58
3.5 NPH Element Stiffness Matrix
The potential energy in a linear elastic body is defined as follows:
{}[]{} { }[]{} {}{ }

+ , = :
V
T T T
dV D D
0 0 P
2
1
e e e e e
{ } { } { } { } { } { }

, , ,
VS
T T
F dS dV U u F u (3.44)
where
{e}
&
'
&
(
)
&

&

=
12
22
11
e
e
e
is the strain field
[ ] D = Material property matrix
{ }
0
e and { }
0
= Initial strain and initial stress
{}
'
(
)

=
2
1
u
u
u , Displacement field
{ } F = Body force
{ } = surface traction
{ } U = Displacement vector at points where concentration forces are
applied
{ } F = Concentrated loads applied to particular degrees-of-freedom
59
From the equation (3.44), the multi-scale variational model for the potential energy
in a linear elastic body can be obtained by an averaging process,
{}[]{} {}[]{}{}{ }

'
(
)

+ , = :
VA
T T T
dAdV D D
A
0 0
P
2
1 1
e e e e e
{ } { } { } { } { } { }

, , ,
VS
T T
F dS dV U u F u (3.45)
where, A is the microscopic cell.
() dA
A
Area

1
is the “average over” the microstructure cell of the strain energy
density at the particular macro position in x. Substituting equation (3.41) into
equation (3.45) and assuming that there is no initial strain, initial stress, body force
or surface traction, the following equation can be obtained:
{}[] [ ][ ]{} { } {}{ }

, = :
VA
n
ij
n
T
ij
n
T
n
F dAdV UC P D P UC
A
U
1
2
1
P (3.46)
{}[]{}{}{} F U UC UC , = : K
T
2
1
P (3.47)
where,
[] [] [ ][ ] { } []
5

=
= =
elem of
i
i
A
ij
n
T
ij
n
V
k dAdV P D P
A
K
#
1
1
(3.48)

{}
0 =

:
UC
p
yields (3.49)
[ ]{ } { } F UC = K (3.50)
60
Then, by using the Gauss Quadrature method, the equation (3.48) can be formed at
the particular ij in microcell structure. Thus, the element stiffness matrix is shown:
[] [] [][] { }dAdV P D P
A
k
A
ij
n Vp
T
ij
n
V
elem of
i
i
1
#
1

5
=
=
[] [][] { }
ij
N
i
N
j
A
ij
n Vp
T
ij
n j i
DetJac t P D P
A
WT WT
x
y
;
'
(
)

=
55

== 11

1
{}
ij
N
i
N
j
ij
j i
DetJac t Micro WT WT
x
y
; =
55
== 11
(3.51)
where, [] [ ][] { }

=
A
ij
n Vp
T
ij
n
ij
P D P
A
Micro
1
[][][] ()

=
ij
y
ij
x
L
L
ij
Vp
T
ij
y
dy dy thk P D P
A
00
2 1
1
[][][] ()

,,
=
1
1
1
1
2 1
1
3 3 d d thk DetJac P D P
A
ij
Vp
T
ij
y
[][][] () DetJac thk P D P WT WT
A
ij
Vp
T
M
ii
M
jj
jj ii
ij
y
x
y
*
1
11
55
==
= (3.52)
where, WT = Gauss Quadrature weight factor
DetJac = Determinant of Jacobian matrix
thk = Cell solution thickness
t = Macro structure thickness
Vp
D = Variable material property in the micro cell
61
Also,
ij
Micro can be determined while the microcell structure changes
depending on the ( ) 3
2
hx and ( ) 3
2
hy values at the each Gauss point. Finally, the
NPH stiffness matrix can be formed as shown below:
[] {}
ij
N
i
N
j
ij
j i
elem of
i
i
DetJac t Micro WT WT k
x
y
; =
55 5
== = 11
#
1
(3.53)
where
[][][] () DetJac thk P D P WT WT
A
Micro
ij
Vp
T
M
ii
M
jj
jj ii
ij
y
ij
x
y
*
1
11
55
==
=
62
3.6 NPH Element Stress Calculation
The stress and strain relationship in a linear elastic body with an initial
strain and a stress term is defined as follows:
0 0
e e E + , = ) ( (3.54)
where,
{ }
0
e = Initial strain
{ }
0
= Initial stress
Since there is no initial strain and stress are presented, the equation (3.54) can be
written in terms of the elastic constant and strain. The strain term was defined in
the equation (3.41) which is a product of the [P] and nodal displacement, {UC}.
Thus, the equation (3.54) can be written as below:
Ee = (3.55)
{ } [ ] [ ] { }
elemt x x x
UC P E
54 3 3 3 1 3
= (3.56)
where
[ ] [ ] [ ] [ ] [ ]
54 3 18 3 18 3 18 3 54 3
| |
x X X X x
P u C C =
63
At the microcell structure level shown in Figure 3.2, the stresses were
computed based on the material properties at the location of the micro structural
with the matrix [P
i
]. Therefore, the equation 3.56 can be written as shown below:
{ } [ ] [ ] { }
elemt x x x
UC P E
54 3 1 3 3 1 1 3 16 of 1
= (3.57)
{ } [ ] [ ] { }
elemt x x x
UC P E
54 3 2 3 3 2 1 3 16 of 2
= (3.58)

Figure 3.2 Nodal Displacement at the Mesh Element and
Microstructure: i.e.) E
1
= 10 and E
2
= 1000
U
y7
U
x7
U
y6
U
x6
U
y4
U
x4
U
y8
U
x8
U
y9
U
x9
U
y2
U
y3
U
y1
U
y5
U
x5
E1 E1
E1 E1
E2
E2
E2
E2
E2
E2
E2 E2
E2 E2
E2E2
U
x1
U
x2
U
x3
64
3.7 Typical Allowable NPH Microcell Geometries
A genetic microcell structure in 2-D is shown in Figure 3.3. In the research
proposal, the fiber size hx
1
and hy
1
are constant values and matrix sizes hx
2
and hy
2
are variable sizes in the micro coordinate y (Y
1
, Y
2
). Therefore, overall microcell
structure will vary depending on the location of the integration points in the macro
mesh.
Figure 3.3 A Microcell Structure in 2-D (Shaded areas Represent
Fibers, hx
1
and hy
1
)
Various microcell structures can be created in order to fit the best
characteristic structure of the Functionally Graded Materials. For example, in
Y
1
Y
2
hx
2
(3)
hy
2
(3)
hy
1
hx
1
3
1
3
2
65
order to create the dot fiber micro structure as shown in Figure 3.4 (b), change the
vertical and horizontal fiber material properties from the genetic microcell
structure pictured in Figure 3.3 to the matrix material properties except the center
location. Various types of the microcell structure for the Functionally Graded
Materials can be generated as shown in Figure 3.4.
Figure 3.4 Variation of 2-D Microcell Structures with Fiber (shaded area): (a)
Vertical Direction; (b) Center; (c) Outer Edges; (d) Distributed in the Center and
Four Corners.
(a) (b) (c) (d)
66
3.8 Data Collection Points for hx
2
and hy
2
Matrix Values in the
Nonperiodic Microstructure
In order to characterize the nonperiodic microstructure pattern of the
Functionally Graded Materials, the data collection point method was used. The
meaning of the “data collection point” is to identify the variation of hx
2
matrix
values in 1-D or the variation of hx
2
and hy
2
matrix values in 2-D cases throughout
the macro structure. The structure can have two or more data collection points per
finite element. Figure 3.5 shows an example of the three data collection points in
1-D finite element case. Then, these data collection points were interpolated at the
Gauss location. In order to increase the accuracy, a total of four integration points
per finite element were used. These data collection points could be obtained
manually or automatically and placed in the data input file for the NPH
computation.
Figure 3.5 A Total of Three hx
2
Data Collection Points per One Finite
Element: Black Dots indicate Data Collection Points and Circles with Cross
indicate the Location of Gauss Points.
hx
2
67
Figure 3.6 shows the location of the data collection points of hx
2
and hy
2
in 2-D
case. A total of nine data collection sets were obtained at the near of the element
nodal points. As in the 1-D case, these data collection sets were interpolated at the
Gauss integration points.
Figure 3.6 A Total of Nine Data Collection Points per Finite Element: Black
Dots indicate Data Collection Points and Circles with Crosses indicate the
Location of Gauss Points.
4
1 2
3
5
6
7
8 9
Y
1
Y
2
Pt 6 (hx
2
, hy
2
)
68
Chapter 4
Numerical Examples of NPH 1-D Cases
4.1 Introduction
The intent of this chapter is to test the performance of the nonperiodic
analysis procedure on one dimension problem for which as exact solution exists.
Six independent cases are considers including:
Case 1 – Comparison between the NPH and the Homogenization Solution
Case 2 – Descending Low Density Microcell Structures
Case 3 – Descending High Density Microcell Structure
Case 4 – Descending and Ascending Microcell Structure
Case 5 – Descending Microcell Structure with a Sudden Jump
Case 6 – Rapidly Varying Descending, Ascending and Descending
Microcell Structures.
69
4.2 Case 1 – Comparison between the NPH and the Homogenization
Solutions
4.2.1 Geometry Modeling
A structure which contains periodic microcells was created in order to
compare the deformation results from the NPH and the homogenization methods.
The cell structure possesses typical matrix and fiber configurations. The test
configuration is pictured in Figure 4.1. In the case study, the length of the matrix
component hx
2
and the length of the fiber component hx
1
were defined as 0.375
and 0.25, respectively. The Young’s modulus of the matrix and fiber were
assumed to be E
matrix
= 10 and E
fiber
= 1000. A concentrated force, F was applied
on the free end, and its value was F = 3.0. The total length of the structure was 2.0,
and the cross section area of the structure was also 2.0. In the NPH algorithm, a
single macro finite element was used for the numerical computations. Figure 4.2
showing the single macro element in the macro model for the NPH method and
Figure 4.3 showing detailed view of the microcell structure, cell 1, for exact
solutions.
Figure 4.1 A Periodic Structure with Microcells
hx
2
hx
2
hx
1
hx
2
hx
2
hx
1
F = 3
Cell 1 Cell 2
70
Figure 4.2 A Single Macro Element for the NPH Method
Figure 4.3 Detailed View of the Cell 1 for the Exact
Solutions: hx
2
= 0.375 and hx
1
= 0.250
4.2.2 Local and Global Deformation Analysis
The deformation results were obtained at the three different locations: X =
0.0, X = 1.0 and X = 2.0 which are the node point locations for three nodes in the
macro-model. As is shown in Table 4.1, the homogenization and the NPH results
were identical. This is due to the fact that the nonperiodic terms in the NPH
algorithm are zero for this case. The result confirms that the NPH algorithm
produces the exact homogenization solution for periodic structures. Therefore, the
NPH method works correctly for the periodic cases.
X
Y
hx
2
hx
2
hx
1
71
Table 4.1 Deformation Results of NPH and Homogenization Methods
Location in X Homogenization NPH
0.0 0.0000E+00 0.0000E+00
1.0 0.1129E+00 0.1129E+00
2.0 0.2258E+00 0.2258E+00
72
4.3 Case 2 – Descending Low Density Microcell Structure
A Typical nonperiodic low density microcell structure was created to
investigate the accuracy of the NPH algorithm. Three different methods were used
as follows:
Method 1 - mapping the geometry with a single finite element and using
two hx
2
data collection points per finite element, picture in Figure 4.4(a).
Method 2 - mapping the geometry with a single finite element mesh and
using three hx
2
data collection points per finite element, picture in Figure
4.4(b).
Method 3 - mapping with the geometry with four finite elements and
using three hx
2
data collection points per finite element, picture in Figure
4.4(c).
Figure 4.4 Three Different hx
2
Data Collection Methods; Black Dots indicate Data
Collection Points and Circles with Cross indicate the Location of Gauss Points.
hx
2
hx
2
(a) Method 1- 2 Collection Pts.
with 1 Finite Element
(b) Method 2 – 3 Collection Pts.
with 1 Finite Element
(c) Method 3 – 3 Collection Pts.
with 4 Finite Elements
hx
2
73
4.3.1 Geometry Modeling
Consider the structure picture in Figure 4.5. The configuration of the
microcell 1 was the same as cell 1 in Figure 4.3. The length of the matrix hx
2
is
0.375 and the length of the fiber hx
1
is 0.25. Then, the microcell matrix hx
2
values,
for cell 2 and 3, were decreased as compare to cell 1. The values of hx
2
for cell 2
and for cell 3 are 0.225 and 0.025, respectively. The fiber length was not changed
and took the value as 0.25. All boundary conditions were as presented in section
4.2.1.
Figure 4.5 Nonperiodic Low Density Microcell Structure
4.3.2 Local and Global Deformation Analysis
The deformation results, comparing the solution using method 1, 2 and 3,
are displayed in Table 4.2. The results indicate that when more data collection
points were employed in the finite elements, the accuracy of the solution was
increased. However, adding more finite elements while collecting the using the
same number of data points did not significantly improve the results. This
F
Cell 1 Cell 2 Cell 3
hx
2
hx
2
hx
2
74
occurred because the density of the microcell structures in the geometry was too
coarse. Figure 4.6 shows that the NPH results reasonably follow the exact solution
when the number of the finite element is increased.
Table 4.2 Deformation Results for the Descending Low Density Microcell
Structure Case
Computational Methods X = 1.0 X = 2.0
Percent of
the Exact Solution
at x = 2.0
Exact Solution 0.11288 0.18613 -
Method 1 – 2 Collection
points w/ 1 Finite Element
0.10662 0.16960 91.1
Method 2 – 3 Collection
points w/ 1 Finite Element
0.11027 0.17505 94.1
Method 3 – 3 Collection
points w/ 4 Finite Elements
0.14588 0.17257 92.7
Descending Low Density Microcells with
Three Different Data Collection Methods
0.00
0.04
0.08
0.12
0.16
0.20
0.00 0.50 1.00 1.50 2.00
Distance
Deformation
Exact Soln
2 Collect. Pts. w/ 1 Finite Ele.
3 Collect. Pts. w/ 1 Finite Ele.
3 Collect. Pts. w/ 4 Finite Eles.
Figure 4.6 Deformation Results for the Descending Low Density Microcell with
Three Different Data Collection Methods
75
4.4 Case 3 – Descending High Density Microcell Structure
4.4.1 Geometry Modeling
From an application perspective the high density microstructure is very
important. Thus, a nonperiodic high density microcell structure was analyzed to
investigate the accuracy of the NPH algorithm for these problems. As it is shown
in the Figure 4.7, the numbers of microcells were added in the geometry. Three
different methods which were described in section 4.3 were tested: Method 1,
Method 2 and Method 3.
Figure 4.7 Descending Microcell Structure Case with 10 percent
Descending
The geometry contains a total of eighteen microcells. The configuration of
the microcell 1 was set as followed; the length of the matrix hx
2
is 0.090 and the
length of the fiber hx
1
is 0.0262. The length of the matrix value was reduced
continuously by 10 percent from cell to cell, starting from the cell 1. Thus, the
length of the matrix hx
2
at the cell 18 is 0.015. The matrix size values for the
F
Cell 1 Cell 18
76
eighteen cells are presented in Table 4.3. The fiber dimension was kept constant,
which was 0.0262. All boundary conditions were as presented in section 4.2.1.
Table 4.3 The Matrix Size Values for the Descending and Ascending Microcell
Structure: Fiber size hx
1
= 0.0262
Cell
No.
Matrix Size
hx
2
Cell
No.
Matrix Size
hx
2
Cell
No.
Matrix Size
hx
2
1 0.0900 7 0.0478 13 0.0254
2 0.0810 8 0.0430 14 0.0229
3 0.0729 9 0.0387 15 0.0206
4 0.0656 10 0.0349 16 0.0185
5 0.0590 11 0.0314 17 0.0167
6 0.0531 12 0.0282 18 0.0150
4.4.2 Local and Global Deformation Analysis
The calculated deformation values for method 1, 2 and 3 are presented in
Table 4.4. It can be seen that as compared to Case 2 in section 4.3, more
consistent convergence trends are evident, and this is due to the increase in the
number of microcells in the model. For example, the Method 3 used the most
refined model in terms of numbers of elements and data collection points and
produced the best accuracy. It is clear that collecting more matrix size values hx
2
,
in each finite element, helped to increase the accuracy. And adding more finite
77
elements in the geometry also produces better results. In Figure 4.8 and 4.9, the
displacement patterns for the three methods are pictured.
Table 4.4 The Deformation Results of the High Density Microcell Structure with
10 percent Reduction
Computational Methods x = 1.0 x = 2.0
Percent of
The Exact Solution
at x = 2.0
Exact Solution 0.12675 0.23018 -
Method 1 – 2 Collection
points w/ 1 Finite Element
0.13675 0.22448 97.5
Method 2 – 3 Collection
points w/ 1 Finite Element
0.13695 0.22474 97.6
Method 3 – 3 Collection
points w/ 4 Finite Elements
0.12722 0.22766 98.9
78
Descending Microcell Structure in High Density
Microcell
0.00
0.05
0.10
0.15
0.20
0.25
0.00 0.50 1.00 1.50 2.00
Distance
Deformation
Exact Soln
2 Collect. Pts. w/ 1 Finite Ele.
3 Collect. Pts. w/ 1 Finite Ele.
3 Collect. Pts. w/ 4 Finite Eles.
Figure 4.8 Deformation Results of Descending Microcell Structure in High Density
Microcell with Three Different Data Collection Methods
Detail View of the Descending Microcell
Structure in High Density Microcell
0.12
0.14
0.16
0.18
0.20
0.22
0.24
1.00 1.20 1.40 1.60 1.80 2.00
Distance
Deformation
Exact Soln
2 Collect. Pts. w/ 1 Finite Ele.
3 Collect. Pts. w/ 1 Finite Ele.
3 Collect. Pts. w/ 4 Finite Eles.
Figure 4.9 Detailed View of the Figure 4.5 (Red Box Region)
79
4.5 Case 4 – Descending and Ascending Microcell Structure
4.5.1 Geometry Modeling
The NPH algorithm was tested for structure with continuously varying
microcells with the cell size both descending and ascending in terms of matrix size
hx
2
. A test case was created as pictured in Figure 4.10. In this case again the cell
size is slowly varying, while ascending and descending, in terms of the matrix size
hx
2.
The geometry contains a total of fifteen microcells. Four macro finite
elements were used with three hx
2
data collection points per finite element.
The first microcell in the geometry was identical to Case 3; the length of
the matrix hx
2
is 0.090 and the length of the fiber hx
1
is 0.0262. Then, 10 percent
reduction was made from the first microcell matrix value hx
2,
until it reached
approximately 3/4 of the geometry. After that, the matrix value was increased by
116 percent continuously until it reached the end of the structure. Thus, the matrix
value hx
2,
at the end of the structure became 0.0578. The fiber dimension was kept
constant value which was 0.0262. The matrix size values for the fifteen cells are
presented in Table 4.5. All boundary conditions were as presented in section 4.2.1.
Figure 4.10 High Density Microcell Structure Case with 10 percent
Descending and 16 percent Ascending Matrix Sizes
F
Cell 1 Cell 15
80
Table 4.5 The Matrix Size Values for the Descending and Ascending Microcell
Structure: Fiber size hx
1
= 0.0262
Cell
No.
Matrix Size
hx
2
Cell
No.
Matrix Size
hx
2
Cell
No.
Matrix Size
hx
2
1 0.0900 6 0.0531 11 0.0314
2 0.0810 7 0.0478 12 0.0366
3 0.0729 8 0.0430 13 0.0426
4 0.0656 9 0.0387 14 0.0496
5 0.0590 10 0.0349 15 0.0578
4.5.2 Local and Global Deformation Analysis
Table 4.6 compares the NPH deformation results at global coordinates X =
1.0 and X = 2.0 with exact solution and the homogenization solution. The
homogenization solution was added in the table in order to show the non-linearity
of the NPH solution. Even though only four macro finite elements were used in
the NPH model, the deformation results are extremely accurate. In Figure 4.11
and 4.12, the deformation results for NPH are compared to the exact solution.
Table 4.6 Global Deformation Values in High Density Microcell Structure Case
with 10 percent Descending and 16 percent Ascending Matrix Sizes
Computational Methods x = 1.0 x = 2.0
Percent of
The Exact Solution
at x = 2.0
Exact Solution 0.12675 0.24183 -
3 samples w/ 4 Finites 0.12726 0.23948 99.0
Homogenization Soln 0.13116 0.26233 108.5
81
High Density Microstructure Case with 10% descending and
16% Ascending Matrix Sizes
0.00
0.05
0.10
0.15
0.20
0.25
0.30
0.00 0.50 1.00 1.50 2.00
Distance
Deformation
Homogenization Soln
Exact Soln
3 Collect. Pts. w/ 4 Finite Eles
Figure 4.11 Deformation Results for High Density Descending and Ascending
Microcell Structure
Detail View of the High Density Microstructure Case with
10% descending and 16% Ascending Matrix Sizes
0.15
0.17
0.19
0.21
0.23
0.25
0.27
1.50 1.60 1.70 1.80 1.90 2.00
Distance
Deformation
Homogenization Soln
Exact Soln
3 Collect. Pts. w/4 Finite Eles
Figure 4.12 Detailed View of the Figure 4.5 (Red Box Region)
82
4.6 Case 5 – Descending Microcell Structure with a Sudden Jump
4.6.1 Geometry Modeling
One of the FGMs characteristics is that the geometry has a gradual
variation in space. However, in this case, a problem with a sudden jump of
microcell structure geometry was created and investigated. The microstructural
model is pictured in Figure 4.13. The geometry involved a total of fifteen
microcells, and the model incorporated four macro finite elements and used three
hx
2
data collection points per finite element.
In the developed model, the microcell matrix size hx
2
was decreased by 15
percent from cell to cell, starting from the left hand side and a sudden jump
(approximately 5 times more than its neighbor matrix size) was generated at the
middle of structure. After that, the matrix value hx
2
was decreased again by 34
percent from cell to cell. The cell matrix size values are presented in Table 4.7.
All other boundary conditions were as same as Case 1 in section 4.2.1.
Figure 4.13 Descending Micro Structure with a Sudden Jump
F
Cell 1
Cell 15
83
Table 4.7 The Matrix Size Values for the Descending Microcell Structure with a
Sudden Jump: Fiber size hx
1
= 0.0194
Cell
No.
Matrix Size
hx
2
Cell
No.
Matrix Size
hx
2
Cell
No.
Matrix Size
hx
2
1 0.0900 6 0.0399 11 0.0655
2 0.0765 7 0.0339 12 0.0433
3 0.0650 8 0.0289 13 0.0286
4 0.0553 9 0.1500 14 0.0189
5 0.0470 10 0.0992 15 0.0125
4.6.2 Local and Global Deformation Analysis
Table 4.8 presents the displacement results for the case of geometry with a
sudden jump. The comparisons, between the exact solution and the NPH solution
were made at four different locations: X = 0.5, 1.0, 1.5 and 2.0. In Figure 4.14, the
deformation patterns for the NPH method are compared to the exact solution. The
exact solution clearly showed that there are two connected but independent curves:
one from before the sudden jump and the other one after the sudden jump. Figure
4.14 shows that the NPH solution does not accurately follow the exact solution for
the case with sudden jumps. The reason is that the NPH algorithm uses the
homogenized displacement solution
H
i
u , for correcting its deformation values.
Thus, if there is a sudden jump at the micro structural level, the mathematical
correction cannot be accurately determined.
84
Table 4.8 Deformation Results for the High Density Descending Microcell
Structure with a Sudden Jump
Computational Methods X = 0.5 X = 1.0 X = 1.5 X = 2.0
Exact Solution 0.06954 0.13118 0.21585 0.25681
NPH Solution
(3 samples w/ 4 Finites)
0.06895 0.11546 0.19124 0.23861
Descending Microcell Structure
with a Sudden Jump
0.00
0.05
0.10
0.15
0.20
0.25
0.30
0.00 0.50 1.00 1.50 2.00
Distance
Deformation
Exact Soln
3 Collect. Pts. w/ 4 Finite Eles.
Figure 4.14 Deformation Results of the High Density Descending Micro Structure
with a Sudden Jump
85
4.7 Case 6 – Rapidly Varying Descending, Ascending and
Descending Microcell Structures
4.7.1 Geometry Modeling
The tests in this section were designed to explore the ramification of the
sudden jump result. Consider a case in which the matrix size value hx
2
is rapidly
varying, but does not have a sudden jump. Thus, the matrix value hx
2
was
continuously decreased, increased and decreased in a smooth but rapidly varying
way. This generates a smooth transition among the microcells. The geometry
includes a total of fifteen microcells and is pictured in Figure 4.15. The model
incorporates four finite elements with three hx
2
data collection points per finite
element.
The microcell was reduced by 30 percent continuously from the first
microcell until it reached the middle of the structure. After that it was increased
by 130 percent then decreased 25 percent until the end of the structure was reached.
The associated matrix size hx
2
parameters are presented in Table 4.9. All
boundary conditions were as presented in section 4.2.1.
Figure 4.15 High Density Microcells with Rapidly Varying Descending,
Ascending and Descending Structure
F
Cell 1
Cell 15
86
Table 4.9 The Matrix Size Values for the Rapidly Varying Descending, Ascending
and Descending Structure: Fiber size hx
1
= 0.070
Cell
No.
Matrix Size
hx
2
Cell
No.
Matrix Size
hx
2
Cell
No.
Matrix Size
hx
2
1 0.0900 6 0.0151 11 0.0187
2 0.0630 7 0.0197 12 0.0140
3 0.0441 8 0.0256 13 0.0105
4 0.0309 9 0.0332 14 0.0079
5 0.0216 10 0.0249 15 0.0059
4.7.2 Local and Global Deformation Analysis
Table 4.10 shows the deformation results at four different locations: X =
0.49, 0.95, 1.43 and 1.90. Figure 4.16 compares the NPH solution and the exact
solution. These results indicate that the NPH method can accurately follows the
exact solution for the case in which the microstructure is rapidly varying. Again,
the use of three data collection point method, in NPH algorithm, is validated.
Table 4.10 Deformation Results of the High Density Microcell with Rapidly
Varying Descending, Ascending and Descending Microcell Structure
Computational Methods x = 0.49 x = 0.95 x = 1.43 x = 1.90
Exact Solution 0.05273 0.08004 0.11148 0.12911
NPH Solution
(3 samples w/ 4 Finites)
0.04838 0.07785 0.10779 0.12733
87
Microcell with Rapidly Varying
Descending, Ascending & Descending
0.00
0.02
0.04
0.06
0.08
0.10
0.12
0.14
0.0 0.4 0.8 1.2 1.6 2.0
Distance
Deformation
Exact Soln
3 Collect. Pts. w/ 4 Finite Eles.
Figure 4.16 Deformation Results of the High Density Microcell with Rapidly
Varying Descending, Ascending and Descending Structure
88
Chapter 5
Numerical Examples of NPH 2-D Cases
5.1 Introduction
After carefully tested the performance of the NPH algorithm in
nonperiodic 1-D cases, further studies were conducted in two dimensional
cases. In this chapter, not only the global and local deformation values were
determined, but also local Von-Mises stresses were computed. These results,
which were computed by the NPH method, were compared with the results of
commercially available Finite Element Analysis (FEA) software. Five
independent 2-D FGMs cases were considered and they are:
Case 7 – Periodic Microstructure
Case 8 – Descending Horizontal Fiber Strips in One Direction
Case 9 – Descending and Ascending FGMs with Square Fibers
Case 10 – Descending and Symmetric Matrix Structure
Case 11 – Descending Horizontal and Vertical Fiber Strips
5.2 Case 7 – Periodic Microstructure
5.2.1 Geometry Modeling
In order to validate the performance of the NPH method in 2-D cases, a
structure with periodic microcells was created. The deformation results from the
89
NPH method and the homogenization method were then compared. The
microstructural model is pictured in Figure 5.1. A total of sixteen macro elements
were used to mesh the entire geometry in the NPH method. Nine data collection
points were utilized per finite element.
The lengths of the matrix component, hx
2
and hy
2
, are 0.085 and the lengths
of the fiber component hx
1
and hy
1
are 0.030 in both X
1
and X
2
directions. The
geometry includes a total of 100 periodic microcells. Thus, the size of the
geometry was 2.0 by 2.0 with the thickness 1.0. The Young’s modulus of the
matrix E
1
and fiber E
2
were assumed to be 10 and 1000, respectively. The Poisson
ratio for both the matrix and the fiber was assumed to be 0.3. The distribution
force F was applied at the free end of the structure in X
1
direction and its value
was 3.0.
Figure 5.1 Periodic Microcell Geometry and Boundary Conditions
(a) Periodic Microcell (b) 16 Macro Elements in NPH
F
X
1
X
2
Fiber Strip
90
5.2.2 Local and Global Deformation Analysis
As it is shown in Figure 5.2, the deformation results were obtained at the
three different locations: X
1
= 1.0, 1.5 and 2.0. The homogenization and the NPH
deformation values are presented in Table 5.1. The results indicate that there is no
difference between the NPH method and the Homogenized method at the nodal
points. Therefore, the NPH algorithm works correctly for the 2-D periodic case.
Deformation Results for the Periodic Microstructure
0.0
0.5
1.0
1.5
2.0
2.5
0.000 0.004 0.008 0.012 0.016 0.020
Deformation
Distance
NPH at X1=2.0
Homog. at X1=2.0
NPH at X1=1.5
Homog. at X1=1.5
NPH at X1=1.0
Homog. at X1=1.0
Figure 5.2 Homogenization and NPH Deformation Results for the Periodic
Microstructure at X
1
= 1.0, X
1
= 1.5 and X
1
= 2.0
91
Table 5.1 Comparison of the Deformation between Homogenized Solution and NPH
Solution @ X
1
=1.0, 1.5 & 2.0
Deformation at X
1
= 1.0
Location at X
2
HOMOG NPH
2.00 8.917E-03 8.916E-03
1.75 8.906E-03 8.905E-03
1.50 8.899E-03 8.897E-03
1.25 8.897E-03 8.894E-03
1.00 8.896E-03 8.893E-03
0.75 8.897E-03 8.894E-03
0.50 8.899E-03 8.897E-03
0.25 8.906E-03 8.905E-03
0.00 8.917E-03 8.916E-03
Deformation at X
1
= 1.5
Location at X
2
HOMOG NPH
2.00 1.337E-02 1.337E-02
1.75 1.336E-02 1.336E-02
1.50 1.336E-02 1.335E-02
1.25 1.335E-02 1.335E-02
1.00 1.335E-02 1.335E-02
0.75 1.335E-02 1.335E-02
0.50 1.336E-02 1.335E-02
0.25 1.336E-02 1.336E-02
0.00 1.337E-02 1.337E-02
Deformation at X
1
= 2.0
Location at X
2
HOMOG NPH
2.00 1.782E-02 1.782E-02
1.75 1.782E-02 1.781E-02
1.50 1.781E-02 1.781E-02
1.25 1.781E-02 1.781E-02
1.00 1.781E-02 1.780E-02
0.75 1.781E-02 1.781E-02
0.50 1.781E-02 1.781E-02
0.25 1.782E-02 1.781E-02
0.00 1.782E-02 1.782E-02
92
5.3 Case 8 – Descending Horizontal Fiber Strips in One Direction
5.3.1 Geometry Modeling
The high density microcell strip structure was created to evaluate the
performance of the NPH algorithm in 2-D. In the structure, a total of twenty
groups of cell strip were used. The matrix values were decreased in only X
2
direction and the structure model is pictured in Figure 5.3. In the NPH method, a
total of sixteen macro elements were used to compute the global deformation
values.
As is shown in Figure 5.4, the matrix size hy
2
value is 0.0850 at the cell
strip 1 and there are imaginary lines which separate the cell groups. The matrix
size was linearly reduced by 11 percent starting from the cell strip 1 to cell strip 20.
Thus, the last cell matrix value became 0.0093. The matrix size values for the 20
cells are presented in Table 5.2. The material properties and the boundary
conditions were as presented in section 5.2.1.
Figure 5.3 Descending Horizontal Fiber Strips Microcell Structure: (a) Descending
Horizontal Fiber Strips in X
2
Direction, (b) 16 Macro Elements for NPH
hy 2
X 1
(a) Horizontal Fiber
S i
(b) 16 Macro Elements for NPH
X 2
F
93
Figure 5.4 Detailed View of the Microcell Strips: Fiber Value
hy
1
= 0.030 (Constant)
Table 5.2 Matrix Size Values for the Descending Horizontal Fiber Strips
Cell No.
Matrix Size
hy
2
Cell No.
Matrix Size
hy
2
1 0.0850 11 0.0265
2 0.0757 12 0.0236
3 0.0673 13 0.0210
4 0.0599 14 0.0187
5 0.0533 15 0.0166
6 0.0475 16 0.0148
7 0.0422 17 0.0132
8 0.0376 18 0.0117
9 0.0335 19 0.0104
10 0.0298 20 0.0093
hy
2
= 0.0850
hy
2
= 0.0850
Cell
Strip 1
hy
2
= 0.0757
hy
2
= 0.0757
Imaginary Lines
Cell
Strip 2
Fiber strip
94
5.3.2 Local and Global Deformation Analysis
Using the conventional FEA method, models for a high density structure
are required to have massive numbers of macro elements. Figure 5.5 illustrates the
characteristic of the macro mesh. Table 5.3 compares the model sizes associated
with the conventional FEA method and the NPH method.
Table 5.3 Summary of Model Sizes for Case 8
FEA (ABAQUS)
Method
NPH
Method
Total Macro Elements 17,956 16
No. of Integration Points
per Element
4 (2X2) 16 (4X4)
Element Type Four (4) node Lagrange Nine (9) node Lagrange
Nodal DOF/ Total DOF 2/18,225 6/486
Figure 5.5 Conventional FEA Mesh and Boundary Conditions
95
The local and global deformation results for the FEA and the NPH
approaches were compared and pictured in Figure 5.6. In the matrix region, the
distributed force was directly acting on the edge of the matrix. It should be noted
that in Figure 5.6 the NPH displacement at the macro-nodes are presented. This is
because the nature of the materials involved a microstructural behavior is
associated cell solution which only varied in the X
2
direction. This is a
characteristic of the homogenization method.
Descending Horizontal Fiber Strips Case
at X= 1.995
0
0.5
1
1.5
2
2.5
0.000 0.005 0.010 0.015 0.020 0.025 0.030 0.035 0.040 0.045
Deformation Value
X2 - location
FEA Solution
NPH Global Def. U(x)
Figure 5.6 Local Deformation Results: the NPH Solution and the FEA
Solution at X
1
= 1.995
96
5.4 Case 9 – Descending and Ascending FGMs with Square Fibers
5.4.1 Geometry Modeling
The NPH algorithm was tested for the descending and ascending matrix
case with square fibers. The test model is shown in Figure 5.7 (a). A total of the
225 microcell groups were presented in the model. The microcells were reduced
by 10 percent in both X
1
and X
2
directions until the cells reached approximately
3/4 of the length. Then the microcells were increased by 117 percent until they
reached the end of the structure.
To understand the details of the modeling, Figure 5.8 is included. In
Figure 5.8, it can be seen that there are imaginary lines which separate the cell
groups. The initial matrix dimension started with 0.090 for both hx
2
and hy
2
and
the last matrix element value was 0.0578. For example, cell 2 has the matrix
values hx
2
= 0.081, hy
2
= 0.090 and the fiber value hx
1
= hy
1
= 0.030. The matrix
values for the X
1
and X
2
direction are presented in Table 5.4. The material
properties and the boundary conditions were as presented in section 5.2.1.
Figure 5.7 Descending and Ascending Matrix with Square Fibers
(a) Square Fibers (b) 16 Macro Elements for NPH
Descending Ascending
X 1
X 2
Cell 1
Cell 225
97
Figure 5.8 Detailed View of the Square Fibers: Fiber Value hx
1
&
hy
1
= 0.030 (Constant)
Table 5.4 Matrix Sizes in Descending and Ascending Square Fibers
Cell No.
(X
1
- dir.)
Matrix Size
hx
2
/ hy
2
Cell No.
(X
2
- dir.)
Matrix Size
hx
2
/ hy
2
1 0.0900 / 0.0900 1 0.0900 / 0.0900
2 0.0810/ 0.0900 16 0.0900 / 0.0810
3 0.0729/ 0.0900 31 0.0900 / 0.0729
4 0.0656/ 0.0900 46 0.0900 / 0.0656
5 0.0590/ 0.0900 61 0.0900 / 0.0590
6 0.0531/ 0.0900 76 0.0900 / 0.0531
7 0.0478/ 0.0900 91 0.0900 / 0.0478
8 0.0430/ 0.0900 106 0.0900 / 0.0430
9 0.0387/ 0.0900 121 0.0900 / 0.0387
10 0.0349/ 0.0900 136 0.0900 / 0.0349
11 0.0314/ 0.0900 151 0.0900 / 0.0314
12 0.0366/ 0.0900 166 0.0900 / 0.0366
13 0.0426/ 0.0900 181 0.0900 / 0.0426
14 0.0496/ 0.0900 196 0.0900 / 0.0496
15 0.0578/ 0.0900 211 0.0900 / 0.0578
hx 2 =
0.090
hx 2 =
0.090
Imaginary
Lines
hy 2 = 0.090
hy 2 = 0.090
Cell 1
Cell 2 Cell 3
Cell 16
Cell 17
Cell 18
hx 2 =
0.081
hx 2 =
0.081
X
1
hy 2 = 0.081
hy 2 = 0.081
98
5.4.2 Local and Global Deformation Analysis
In the conventional FEA program, the model required a total of 24,235
macro elements in order to obtain accurate results. In the same model, the NPH
method used 16 macro elements to compute the local and global deformation
values. Table 5.5 compares the model sizes for the FEA method and the NPH
method. Also, the FEA mesh configuration is shown in Figure 5.9 and Figure 5.10.
The deformation results of the conventional FEA method are pictured in Figure
5.11.
Table 5.5 Summary of Model Sizes for Case 9
FEA (ABAQUS)
Method
NPH
Method
Total Macro Elements 24,235 16
No. of Integration Points
per Element
4 (2X2) 16 (4X4)
Element Type Four (4) node Lagrange Nine (9) node Lagrange
Nodal DOF/ Total DOF 2/49,084 6/486
99
Figure 5.9 Conventional FEA Mesh Sizes and Boundary Conditions
Figure 5.10 Detailed View of the FEA Mesh Sizes and Boundary
Conditions (Red Box Region)
Figure 5.11 Deformation Results Using the Conventional FEA Method
100
As shown in Figure 5.12, the deformation results for the conventional FEA
method and the NPH method are compared at X
1
= 1.9427. The location is set
because the cell solutions in the microstructure were computed adjacent to the
square fiber in the NPH method. The NPH local and global deformation profiles
accurately matched the FEA results. Both the FEA method and the NPH method
produced less deformation in the area, which had more square fibers.
Descending and Ascending FGMs with Square Fibers
at X1 = 1.9427
0 . 0
0 . 5
1 . 0
1 . 5
2 . 0
2 . 5
0 . 2 6 0 0 . 2 6 5 0 . 2 7 0 0 . 2 7 5 0 . 2 8 0
Deformation
Distance
FEM Fine Mesh
NPH Global Def. U
NPH Local Def.
Figure 5.12 Deformation Results of the Conventional FEA Method and
the NPH Method at X
1
= 1.9427
101
5.5 Case 10 – Descending and Symmetric Matrix Structure
5.5.1 Geometry Modeling
The high density nonperiodic descending structure was created to evaluate
the NPH method. In this case, fiber strips were included in both in the X
1
direction and the X
2
direction. As is shown in Figure 5.13, the vertical strips were
gradually reduced until they reached the end of the structure. On the other hand,
the horizontal strips also were reduced, but the reduction started from the center of
the structure. Thus, the model has a symmetric condition about the horizontal
centerline of the model. The material properties and the boundary conditions were
as presented in section 5.2.1.
Figure 5.13 Descending Matrix Microstructures in X
1
Direction and
Symmetric condition in X
2
Direction with Descending Microstructures
(a) Symmetric FGMs (b) 16 Macro Elements for NPH
X 1
X 2
Decreasing
Decreasing
102
The detailed view of the structure is pictured in Figure 5.14. The
microcells were constructed in the “Fiber-Matrix-Fiber” pattern in the X
2
direction.
Using the imaginary lines between the vertical and the horizontal fibers, a total of
560 microcells were generated in the structure. Each microcell structure has a pair
of the matrix sizes hx
2
and hy
2
: for example, the cell 3 has hx
2
= 0.0673 and hy
2
=
0.0093. The fiber size hx
1
and hy
1
were constant as 0.030. The matrix sizes hx
2
and hy
2
in the cells are presented in Table 5.6. The material properties and the
boundary conditions were as presented in section 5.2.1.
Figure 5.14 Detailed View of Figure 5.13 (Red Box Region)
hy
2
=
2 x 0.0104
hx
2
=
0.0850
hx
2
=
0.0850
hx
2
=
0.0757
hx
2
=
0.0757
hx
2
=
0.0673
hx
2
=
0.0673
hy
2
=
2 x 0.0093
hy
2
=
2 x 0.0117
Cell 1 Cell 2 Cell 3
Cell 21
Cell 41
Imaginary Lines
X
1
X
2
Fibers
103
Table 5.6 Matrix Sizes in Symmetric and Descending Case
Cell No.
(X
1
- dir.)
Matrix Size
hx
2
/ hy
2
Cell No.
(X
2
- dir.)
Matrix Size
hx
2
/ hy
2
1 0.0850 / 0.0093 1 0.0850 / 0.0093
2 0.0757 / 0.0093 21 0.0850 / 0.0104
3 0.0673 / 0.0093 41 0.0850 / 0.0117
4 0.0599 / 0.0093 61 0.0850 / 0.0132
5 0.0533 / 0.0093 81 0.0850 / 0.0148
6 0.0475 / 0.0093 101 0.0850 / 0.0166
7 0.0422 / 0.0093 121 0.0850 / 0.0187
8 0.0376 / 0.0093 141 0.0850 / 0.0210
9 0.0335 / 0.0093 161 0.0850 / 0.0236
10 0.0298 / 0.0093 181 0.0850 / 0.0265
11 0.0265 / 0.0093 201 0.0850 / 0.0298
12 0.0236 / 0.0093 221 0.0850 / 0.0335
13 0.0210 / 0.0093 241 0.0850 / 0.0376
14 0.0187 / 0.0093 261 0.0850 / 0.0422
15 0.0166 / 0.0093 281 0.0850 / 0.0422
16 0.0148 / 0.0093 301 0.0850 / 0.0376
17 0.0132 / 0.0093 321 0.0850 / 0.0376
18 0.0117 / 0.0093 341 0.0850 / 0.0298
19 0.0104 / 0.0093 361 0.0850 / 0.0265
20 0.0093 / 0.0093 381 0.0850 / 0.0236
21 n/a 401 0.0850 / 0.0210
22 n/a 421 0.0850 / 0.0187
23 n/a 441 0.0850 / 0.0166
24 n/a 461 0.0850 / 0.0148
25 n/a 481 0.0850 / 0.0132
26 n/a 501 0.0850 / 0.0117
27 n/a 521 0.0850 / 0.0104
28 n/a 541 0.0850 / 0.0093
104
5.5.2 Local and Global Deformation Analysis
The comparison of the model sizes between the conventional FEA method
and the NPH method are shown in Table 5.7. The FEA method required an
element number approximately a thousand times more than the NPH method; the
FEA employed 18,900 elements versus 16 elements used in the NPH model.
Table 5.7 Summary of Model Sizes for Case 10
FEA (ABAQUS)
Method
NPH
Method
Total Macro Elements 18,900 16
No. of Integration Points
per Element
4 (2 x 2) 16 (4 x 4)
Element Type Four (4) node Lagrange Nine (9) node Lagrange
Nodal DOF/ Total DOF 2/19176 6/486
The local and the global deformation results of the NPH method and the
FEA method are in Figure 5.15 and Figure 5.16. The NPH deformations
accurately followed the FEA results, and the global deformation U(x) matched the
fiber regions of the FEA curve. The deformation results of the entire structure are
shown in Figure 5.17. The highest deformation occurred at the center of the
geometry due to the lower density of horizontal fibers there.
105
Descending and Symmetric Matrix Case
at X1 = 1.995
0.0
0.5
1.0
1.6
2.1
0.000 0.002 0.004 0.006 0.008
Deformation
Distance
FEM Fine Mesh
NPH Global Def. U(x)
NPH Local Def.
Figure 5.15 Symmetric and Descending Structure: Deformation Displacement
at X
1
= 1.995
Detail View of the Descending and Symmetic Matrix Case
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.002 0.003 0.004 0.005 0.006 0.007 0.008
Deformation
Distance
FEA Fine Mesh
NPH Global Def. U(x)
NPH Local Def.
Figure 5.16 Detailed View of Symmetric and Descending Structure (Red Box
Region)
106
Figure 5.17 Deformation Results of the FEA Method (ABAQUS)
5.5.3 Local Stress Analysis
For this case study, the Von-Mises stress calculation algorithm was added
in the NPH procedure and the NPH stress results were compared with the results
of the commercially available FEA program such as ABAQUS. Figure 5.18
presents the FEA stress contours based on computed nodal point stresses. The
largest stresses occur in the horizontal fibers, not in the vertical fibers when the
force distribution is parallel to the horizontal fibers. The maximum stress was
located at the middle of the structure near the free edge and the Von-Mises stresses
value was 5.811. A detailed view of this high stress region is pictured in Figure
5.19.
107
Figure 5.18 Von-Mises Stress Results from FEA (ABAQUS)
Figure 5.19 Detailed View of the High Stress Region (Red Box
Region): Max. Stress = 5.811
108
NPH was used to compute the Von-Mises stress values in the same area
which was identified by the FEA as the maximum stress region. NPH computed
the stress values at the integration points rather than at the element nodal points.
Thus, the FEA stress values have been converted to the integration points. The
NPH algorithm used 16 integration points while the FEA used four integration
points per element. Figure 5.20 and Figure 5.21 show the stress results. The
summary of the Von-Mises stress results is shown in Table 5.8. It should be noted
that the FEA stress value at the nodal points were higher than at the integration
points.
Figure 5.20 Von-Mises Stress Results by NPH method (Red Box Region in
Figure 5.19): Max. Stress Value = 7.473
109
Figure 5.21 Von-Mises Stress Results by Commercially Available FEA
(ABAQUS): Max. Stress Value = 5.777
Table 5.8 Von-Mises Stress in Descending and Symmetric Case
FEA
Max. Von-Mises Stress
at Nodal Point
FEA
Max. Von-Mises Stress
at Integration Point
NPH
Max. Von-Mises Stress
at Integration Point
5.811 5.777 7.473
110
5.6 Case 11 – Descending Horizontal and Vertical Fiber Strips
5.6.1 Geometry Modeling
In this case, a high density nonperiodic structure was created to evaluate
the NPH method. The structure has vertical and horizontal fiber strips which are
linearly reduced in both X
1
and X
2
directions simultaneously. The model is
pictured in Figure 5.22. A total of sixteen macro elements were used to mesh the
entire geometry in the NPH method and nine matrix data collection points, per
finite element, were used. The material properties and the boundary conditions
were as presented in section 5.2.1.
Figure 5.22 Simultaneously Reducing Matrix Values in X
1
& X
2
Directions
(a) Linearly Descending Structure (b) 16 Macro Elements for NPH
F
X
2
X
1
Decreasing
111
The detailed view of the microcell structure is pictured in Figure 5.23.
The structure contains twenty microcells in each direction and a total of 400
microcells were used. These microcells were constructed by using imaginary lines,
which separate the nonperiodic matrix values. The structure was constructed in
similar manner to that used in section 5.5.2. For example, Cell 21 has the matrix
values hx
2
= 0.0850 and hy
2
= 0.0757. The fiber values were kept at the constant
value of 0.030. However, because of the symmetric condition requirement in the
microcell solution, a half fiber value of 0.015 was used in X
2
direction. The
matrix sizes hx
2
and hy
2
in the microcells are presented in Table 5.9.
Figure 5.23 Detailed View of the Microcell Structure
hx 2 =
0.0850
hx 2 =
0.0850
hx 2 =
0.0757
hx 2 =
0.0757
hx 2 =
0.0673
hx 2 =
0.0673
hy 2 =
0.0850
hy 2 =
0.0850
hy 2 =
0.0757
hy 2 =
0.0757
hy 2 =
0.0673
hy 2 =
0.0673
Cell 1 Cell 2 Cell 3
Cell 21
Cell 41
Imaginary
Lines
X
1
X
2
Fibers
112
Table 5.9 Matrix Sizes in Descending Horizontal and Vertical Strips Case
Cell No.
(X
1
- dir.)
Matrix Size
hx
2
/ hy
2
Cell No.
(X
2
- dir.)
Matrix Size
hx
2
/ hy
2
1 0.0850 / 0.0850 1 0.0850 / 0.0850
2 0.0757 / 0.0850 21 0.0850 / 0.0757
3 0.0673 / 0.0850 41 0.0850 / 0.0673
4 0.0599 / 0.0850 61 0.0850 / 0.0599
5 0.0533 / 0.0850 81 0.0850 / 0.0533
6 0.0475 / 0.0850 101 0.0850 / 0.0475
7 0.0422 / 0.0850 121 0.0850 / 0.0422
8 0.0376 / 0.0850 141 0.0850 / 0.0376
9 0.0335 / 0.0850 161 0.0850 / 0.0335
10 0.0298 / 0.0850 181 0.0850 / 0.0298
11 0.0265 / 0.0850 201 0.0850 / 0.0265
12 0.0236 / 0.0850 221 0.0850 / 0.0236
13 0.0210 / 0.0850 241 0.0850 / 0.0210
14 0.0187 / 0.0850 261 0.0850 / 0.0187
15 0.0166 / 0.0850 281 0.0850 / 0.0166
16 0.0148 / 0.0850 301 0.0850 / 0.0148
17 0.0132 / 0.0850 321 0.0850 / 0.0132
18 0.0117 / 0.0850 341 0.0850 / 0.0117
19 0.0104 / 0.0850 361 0.0850 / 0.0104
20 0.0093 / 0.0850 381 0.0850 / 0.0093
113
5.6.2 Local and Global Deformation Analysis
The model size comparison between the conventional FEA method and the
NPH method are summarized in Table 5.10. The FEA employed a total of 18,900
elements for the structure and the NPH used a total of sixteen elements.
Table 5.10 Summary of Model Sizes for Case 11
FEA (ABAQUS)
Method
NPH
Method
Total Macro Elements 17,807 16
No. of Integration Points
per Element
4 (2 x 2) 16 (4 x 4)
Element Type Four (4) node Lagrange Nine (9) node Lagrange
Nodal DOF/ Total DOF 2/18224 6/486
The NPH and the FEA results for the local and the global deformation are
shown in Figure 5.24. The material was stiffer at the top than the bottom of the
geometry due to more stiff fibers present at the top. Thus, the gradual deformation
results occurred from the top edge to the bottom edge of the structure. The
deformation results of the NPH method were computed at X
1
= 1.995 because the
microcell solutions were computed at Gauss integration points which were
adjacent to the vertical fiber in the microcell. The displacement results computed
by the FEA method are shown in Figure 5.25.
114
Descending Horizontal and Vertical Strips Case
at X1 = 1.995
0 . 0
0 . 5
1 . 0
1 . 5
2 . 0
0 . 0 0 2 0 . 0 0 6 0 . 0 1 0 0 . 0 1 4 0 . 0 1 8 0 . 0 2 2
D e f o r m a t i o n
Distance
NPH Local Disp
NPH U(x) Value
FEM(ABAQUS)
Figure 5.24 Deformation Results of Descending Horizontal and Vertical Strips
Figure 5.25 FEA (ABAQUS) Deformation Result
115
5.6.3 Local Stress Analysis
The Von-Mises Stresses were computed by the FEA and the NPH methods.
The overall stress contour is displayed in Figure 5.26. It can be seen that the
highest stress occurred at the right and the bottom of the geometry: Maximum
Stress = 11.35 at a nodal point. The detailed view of this particular location was
shown in Figure 5.27. The highest stresses occurred at the intersection of the
horizontal and vertical fibers when the force distribution was in action at the free
edge and parallel to the horizontal fibers.
Figure 5.26 Von-Mises Stress Results from FEA (ABAQUS)
116
Figure 5.27 Detailed View of Figure 5.26 (Red Box Region):
Maximum Stress Value is 11.35
The Von-Mises stress was computed using the NPH method. The stress
values were computed at the integration points. Thus, as in section 5.5.2, FEA
stress values, which were computed at the nodal points, have been converted to
integration point values. In this model, the FEA method used a total of four
integration points and the NPH method used a total of sixteen integration points.
The comparison of the local stress values by the FEA and the NPH methods are
presented in Figure 5.28 and Figure 5.29. The summary of the results is shown in
Table 5.11.
Table 5.11 Von-Mises Stress in Horizontal and Vertical Strips
FEA (ABAQUS)
Max. Von-Mises Stress
at Nodal Point
FEA (ABAQUS)
Max. Von-Mises Stress
at Integration Point
NPH
Max. Von-Mises Stress
at Integration Point
11.35 12.11 13.86
117
Figure 5.28 Von-Mises Stress Results by the NPH Method at
Integration Points: Maximum Stress = 13.86
Figure 5.29 Von-Mises Stress Results by FEA (ABAQUS) at
Integration Points: Maximum Stress = 12.11
118
Chapter 6
Summary and Conclusions
6.1 Summary
Based on the new theory of coupling the micro-macro structure, a
nonperiodic homogenized (NPH) algorithm has been developed to solve complex
Functionally Graded Materials (FGMs) problems. In order to verify the
performance and the accuracy of the NPH algorithm, local deformation values,
global deformation values and the Von-Mises stresses were computed. The results
from the NPH method were compared with the exact solutions in 1-D cases and
compared with commercially available FEA software results in 2-D cases. A total
of eleven independent FGMs cases were investigated: 6 cases with 1-D problems
and 5 cases with 2-D problems. And they are as follows:
(A) 1-D cases for FGMs:
Case 1 – Comparison between the NPH and the Homogenization Solution
Case 2 – Descending Low Density Microcell Structures
Case 3 – Descending High Density Microcell Structure
Case 4 – Descending and Ascending Microcell Structure
Case 5 – Descending Microcell Structure with a Sudden Jump
119
Case 6 – Rapidly Varying Descending, Ascending and Descending
Microcell Structures.
(B) 2-D cases for FGMs:
Case 7 – Periodic Microstructure
Case 8 – Descending Horizontal Fiber Strips in One Direction
Case 9 – Descending and Ascending FGMs with Square Fibers
Case 10 – Descending and Symmetric Matrix Structure
Case 11 – Descending Horizontal and Vertical Fiber Strips
The results of the NPH in 1-D cases indicated that the displacement in Case
3, High Density of Microcell Structure, provided better estimation results than in
Case 2, Low Density of Microcell Structure. For example, the accuracy of the
global deformation with respect to the analytical solution was 98 percent and 95
percent, in Case 3 and in Case 2, respectively. For the Case 5, due to geometric
discontinuity, the NPH method did not accurately follow the exact solution – the
NPH results were 93 percent compared to the exact solution. However, when the
microcell elements were continuously increased and/or decreased in the structure
such as Cases 4 and 6, the accuracy of the displacement was presented 99 percent
and 98 percent respectively.
For the 2-D Cases, local and global deformation results for the NPH
method and the FEA method were compared. The complex FGMs structures are
required to have massive number of macro elements and degree-of-freedom (DOF).
120
However, the NPH method used extremely low numbers of macro elements and
the NPH results accurately followed the results of the FEA displacement values.
In Case 2, the macro elements used 17,956 macro elements for the FEA and 16
elements for the NPH.
The local Von-Mises stresses at the integration points were computed in
Cases 10 and 11 based on the global displacement results. For the Case 10,
Descending and Symmetric Matrix Structure, the highest stresses occurred at the
middle of the structure in the horizontal fiber. The stress contour plots at the
microcell level were created to compare the NPH results with the FEA method.
The stress contour plots indicated that the NPH stress plots accurately compare
with the FEA results. The stress values of the NPH method indicated
approximately 29 percent more stress than the FEA estimation.
In Case 11, Descending Horizontal and Vertical Fiber Strips, the highest
stress occurred at the right and the bottom of the structure. The stress contours
were generated using the same method as Case 10. The results indicated that the
maximum stress value was 14 percent higher than the FEA estimation.
121
6.2 Conclusions
A new theory for nonperiodic materials has been developed and verified
for basic test cases. In particular, the NPH theory was developed and verified for
various complex cases in 1-D and 2-D problems. The NPH local and global
deformation results correctly followed the FEA solutions and the Von-Mises stress
values were computed. Two major critical factors were discovered in the cases
studied. One is the ratio of scale value ( / and x y x = ). As was shown in the
cases studied, a small number (high density of microcell structure) provided
better estimation results than the large number (low density of the microcell
structure). In all cases studied, the coefficient value is not a constant value in
the nonperiodic geometry cases. Thus, the coefficient values were varied from
0.018 to 0.180 in the structure. And the changes of the cell size were
approximately 10 percent. The other critical factor is the effectiveness of the NPH
global displacement method. The accuracy of the local NPH displacement is
depended on the global NPH displacement, U(x). Thus, the estimation of the
global displacement values is a critical process.
Overall, the NPH program demonstrated that it is a very efficient tool for
estimating the local and global displacements as well as computing the microcell
the Von-Mises stress levels. The NPH method requires a significantly less model
size compare to the conventional FEA method and reduces the computational time.
Considering the characteristic of the FGMs materials, which usually have a high
density oscillation among the microcell structures and continuously changing
122
volume, the NPH method is an ideal tool for modeling the complex nonperiodic
FGMs structures. Therefore, the NPH method is an attractive method for design
and estimation of the material behavior under various loading conditions.
123
Bibliography
[1] Aboudi, J., Pindera, M. J. and Arnold, S.M., “Higher-order theory for functionally
graded materials,” Composits Part B 30 (1999) 777-832
[2] Bathe, K. J., “Finite Element Procedures,” Prentice-Hall (1996)
[3] Bendsoe, M. P. and Kikuchi, N., “Generating Optimal Topologies in Structural
Design using a Homogenization Method,” Comp. Methods Appl. Mech. Engreg. 71
(1998) 197-224
[4] Berezovski, Engelbrecht, A., J. and Maugin, G., “ Numerical Simulation of Two-
dimensional Wave Propagation in Functionally Graded Materials,” European J.
Mechanics A/Solids 22 (2003) 257-265
[5] Budynas, R. G., “Advanced Strength and Applied Stress Analysis,” McGaw Hill
(1997)
[6] Chen, C. N. and Wellford, L. C. Jr., “Muti-Level Finite Element solution
Algorithms Based on Multiplicative and Additive Correction Procedures,” Int. J.
Numerical Engrg., 28 (1989) 27-41
[7] Cho, J.R. and Ha, D.Y., “Averaging and finite element discretization approaches in
the numerical analysis of functionally graded materials,” Mater. Sci. Engrg. A 302
(2001) 187-196
[8] Cook, R., Malkus, D. and Plesha M. E., “Concepts and Application of Finite
Element Analysis,” 3
rd
ed., John Wiley & Sons (1989)
[9] Fish, J., Shek, K., Pandheeradi, M., M. Shephard, S., “Computational Plasticity for
Composite Structures based on Methematrical Homogenization: Theory and practice,”
Comput. Methods Appl. Mech. Engrg., 148 (1997) 53 – 73
[10] Ghosh, S. and Mukhopadhyay, S.N., “A Material Based Finite Elemtent Analysis
of Heterogeneous Media involving Dirichlet Tessellations,” Comput. Methods Appl.
Mech. Eng. 104 (1993) 211-247
[11] Guedes, J. M. and Kikuchi, N., “Preprocessing and Postprocessing for Materials
Based on The Homogenization Method with Adaptive Finite Element Methods,
Comput. Methods Appl. Mech. Engrg. 83 (1990) 143-198
124
[12] Hassani, B., “A Direct Method to Derive the Boundary Conditions of the
Homogenization Equation for Symmetric Cells,” Comm. Int. Numer. Meth. Engrg., 12
(1996)185-196
[13] Hassani, B. and Hinton, E, “Homogenization and Structural Topology
Optimization”, Springer (1999)
[14] Heroux, M. A. and Thomas, J. W., “A Comparison of FAC and PCG Methods for
solving Composite Grid Problems,” Comm. Appl. Numerical Methods, 8 (1992) 573-
583
[15] Hill, R., “Elastic properties of reinforced properties: some theoretical principles,”
J. Mech. Phys. Solids 11 (1963) 357-372
[16] Kaminski, M., “Homogenized Properties of Periodic N-Component Composites,”
Int. J. Engrg. Sci. 38 (2000) 405-427
[17] Kaminski M., “Boundary Element Method Homogenization of the Periodic
Linear Elastic Fiber Composites,” Engng Analysis with Boundary Elements, 23
(1999) 815-823
[18] Kim, J. G. and Kim, Y. Y., “A New Higher-Order Hybrid-Mixed Curved Beam
Element,” Int. J. Numer. Meth. Engng 43 (1998) 925-940
[19] Kim, J. H. and Paulino, G. H., “Isoparametric Graded Finite Elements for
Nonhomogeneous Isotropic and Orthotropic Materials,” ASME J. Appl. Mech, 69
(2002) 502-514
[20] Lee, Y. D. and Erdogan, F., “Residual/Thermal Stresses in FGM and laminated
Thermal Barrier Coatings,” Int. J. Fract. 69 (1995) 145-165
[21] Markworth, A. J., Ramesh, K. S., and Parks, W.P. Jr., “Modeling Studies Applied
to Functionally Graded Materials,” J. Mater. Sci. 30 (1995) 2183-2193
[22] Maso, G. D. and Antonio, G.F. D., “Composite Media and Homogenization
Theory,” Progress in Nonliner Differential Equations and Their Applications (1991)
[23] Matache, A. M. and Schwab, Ch., “Homogenization via p-FEM for problems
with Microstructure,” Applied Numerical Mathematics 33 (2000) 43-59
[24] Oden, J. T. and Reddy, J.N., “Introduction to the Mathematical Theory of Finite
Elements,” John Wiley & Sons (1976)
125
[25] Pindera, M.J., Aboudi, J. and Arnold, S. M., “Thermomechanical Analysis of
functionally graded thermal barrier coating with different microstructural scales,” J.
Am. Ceram. Soc, 81 [6] (1998) 1525-36
[26] Raghavan, P., Moorthy, S., Ghosh, S., Pagano, N. J., “Revisiting the Composit
Laminate Problem with an Adaptive Muti-level Computational Model,” Composite
Science and Technology, 61 (2001) 1017-1040
[27] Reiter T., Dvorak G., Tvergaard, V., “Micromechanical models for graded
composite materials,” J. Mech. Phys. Solids 45 [8] (1997) 1282-1302
[28] Sadeghi, K. M., “Muti-level Adaptive Finite Element Analysis”, Ph.D.
Dissertation, University of Southern California, Aug. 1997
[29] Santare, M. H. and Lambros, J., “Use of Graded Finite Elements to Model the
Behavior of Nonhomogeneous Materials,” ASME J. Appl. Mech., 67 (2000) 819-822
[30] Shkoller, S., “An Approximate Homogenization Scheme for Nonperiodic
Materials,” Comput. Math. Applic. 33 [4] (1997) 15-34
[31] Sills, L., Eliasi, R. and Berlin, Y., “Modeling of functionally graded materials in
dynamic analyses,” Composites: Part B 33 (2002) 7-15
[32] Smit, R. J., Brekelmans, W. A. M. and Meijer, H. E. H., “Prediction of the
Mechanical Behavior of Nonliner Heterogeneous Systems by Multi-Level Finite
Element Modeling,” Comput. Methods Appl. Mech. Engrg. 155 (1998) 181-192
[33] Tabiei, A. and Ivanov I., “Computational micro-mechanical model of flexible
woven fabric for finite element impact simulation,” Int. J. Numer. Meth. Engng., 53
(2002) 1259-1276
[34] Vemaganti, K. and Deshmukh, P., “An Adaptive global-local approach to
modeling functionally graded materials,” Comput. Methods Appl. Mech. Engrg. 195
(2006) 4230-4243
[35] Wieckowski, Z., “Dual Finite Element Methods in Homogenization for Elastic-
Plastic Fibrous Composit Material,” Int. J. Plasticity 16 (2000) 199-221
126
Appendix A
Nonperiodic Homogenization (NPH) Cell Solution
If[]
[]
[]
[]
18 2
2
1
2
1
1
1
18 4
) (
0
0
0
0
0
0
) (
x
N
x
J
J
DJBX 3 0
3
3
3
3
3
#
#
#
#
#
#
#
#
#
#
$
%

#
$
%

=
,
,
(A.1)
For the equation (3.19),
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

N
N
H
N
N
H
N
N
H
N
N
H
C
X
u
C
X
u
C
X
u
C
X
u
2
2
2
2
1
2
1
2
1
1
1
1
) (
) (
) (
) (
) (
) (
) (
) (
x
x
x
x
x
x
x
x
0
0
0
0
= [] {}
1 18 18 4
2
2
1
1
) (
) ( 0 0 0
0 ) ( 0 0
0 0 ) ( 0
0 0 0 ) (
x
N
i x
H
H
H
H
C DJBX
u
u
u
u
3
3
3
3
3
#
#
#
#
#
#
$
%

(A.2)
127
For the equation (3.20),
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

N
N
N
H
N
N
N
H
N
N
N
H
N
N
N
H
C
X
u
C
X
u
C
X
u
C
X
u
2
2
) ( 2
2
1
) ( 2
1
2
) ( 1
1
1
) ( 1
) (
) (
) (
) (
) (
) (
) (
) (
x
x
x
x
x
x
x
x
0
0
0
0
=[] {}
1 18
18 18
9
2
9
1
1
2
1
1
18 4
0
.
.
0
) (
x
N
i
x
H
H
H
H
x
C
u
u
u
u
DJBX
#
#
#
#
#
#
#
#
$
%

3
(A.3)
For the equation (3.21),
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

N N
H
N N
H
N N
H
N N
H
C
X
u
C
X
u
C
X
u
C
X
u
2
2
2
2
1
2
1
2
1
1
1
1
) (
) (
) (
) (
) (
) (
) (
) (
x
x
x
x
x
x
x
x
0
0
0
0
=[] [] { }
N
i x
N
x
H
H
H
H
x
C
u
u
u
u
DJBX
18 2
2 18
9
2
9
1
1
2
1
1
18 4
) (
0
0
. .
. .
. .
0
0
) ( 3 0 3
#
#
#
#
#
#
#
#
#
$
%

(A.4)
For the equation (3.22),
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

N
N
l
H
k kl
N
N
l
H
k kl
N
N
l
H
k kl
N
N
l
H
k kl
C
X X
u
C
X X
u
C
X X
u
C
X X
u
2
2
2
2
1
2
1
2
1
1
1
1
) ( ) (
) , (
) ( ) (
) , (
) ( ) (
) , (
) ( ) (
) , (
x x
y x
x x
y x
x x
y x
x x
y x
0

0

0

0

(A.5)
128
[] {}
1 18 18 4
2
2
1
1
) (
) , ( 0
) (
) , (
) (
) , (
0
) (
) , (
x
N
i
X
l
H
k kl
l
H
k kl
l
H
k kl
l
H
k kl
C DJB
X
u
X
u
X
u
X
u
#
#
#
#
#
#
#
#
#
#
$
%

=
x
y x
x
y x
x
y x
x
y x

(A.6)
Note that, corresponding microcell approximation at the gauss point @ ij,
ij
x x = (for
natural coordinate system,
ij
3 ), we have the following;
) (
) ( ) (
) ( ) , (
) ( ) , ( ) , (
y
x y
x y x
x y x y x
ij
ij ij
kl
i
ij N kl
i
ij N
x x
N kl
i
x x
N
N kl
i
N kl
i

8
8
8
=
=
=
=
= =
(A.7)
[]
#
#
#
#
#
$
%

#
$
%

=
#
#
#
#
$
%

2
1
2
1
, 2
, 2
, 1
, 1
12
2
22
2
11
2
12
1
22
1
11
1
1 2
2
1
) (
) (
) (
) (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
) (
) , (
) (
) , (
X
H
X
H
X
H
X
H
x
l
H
k kl
l
H
k kl
u
u
u
u
A
X
u
X
u
x
x
x
x
y x y x y x
y x y x y x
x
y x
x
y x

(A.8)
[] {}
1 18 12
2
22
2
11
2
12
1
22
1
11
1
) (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
x
H
i
u DJBX A
#
$
%

#
$
%

= 3

y x y x y x
y x y x y x
(A.9)
129
[][ ]{}
1 18
12
2
22
2
11
2
12
1
22
1
11
1
) (
) ( ) ( ) (
) ( ) ( ) (
x
H
i
u DJBX A
ij ij ij
ij ij ij
3

#
#
$
%

=
y y y
y y y
(A.10)
But,
3 18
9 12
2
9 22
2
9 11
2
9 12
1
9 22
1
9 11
1
1 12
2
1 22
2
1 11
2
1 12
1
1 22
1
1 11
1
18 2
9 1
9 1
12
2
22
2
11
2
12
1
22
1
11
1
) ( ) ( ) (
) ( ) ( ) (
. . .
. . .
) ( ) ( ) (
) ( ) ( ) (
) ( 0 . . ) ( 0
0 ) ( . . 0 ) (
) ( ) ( ) (
) ( ) ( ) (
x
ij ij ij
ij ij ij
ij ij ij
ij ij ij
x
ij ij
ij ij
ij ij ij
ij ij ij
#
#
#
#
#
#
#
#
$
%

#
$
%

=
#
#
$
%

x x x
x x x
x x x
x x x
y y
y y
y y y
y y y

0 0
0 0

(A.11)
[]
3 18
11
18 2
) (
x
N
i x
ij
#
$
%

= 3 0 (A.12)
[] [][ ]{} [] { }
1 18 18 4
4 4
1 18
3 18
11
18 2
) ) ( ) ( (
x
i
X
x
x
H
i
x
N
i x
C DJB u DJBX A
ij
#
$
%

#
$
%

= 3 3 0 (A.13)
For the equation (3.23),
130
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
C
X
u
X
C
X
u
X
C
X
u
X
C
X
u
X
2
2
2
2
1
2
1
2
1
1
1
1
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
x
x y x
x
x y x
x
x y x
x
x y x
0

0

0

0

[] {}
1 18 18 2
2 4
2
2
1
2
2
1
1
1
) (
) ( ) , (
0
) ( ) , (
0
0
) ( ) , (
0
) ( ) , (
x
i
x
N
x
l
H
k
kl
l
H
k
kl
l
H
k
kl
l
H
k
kl
C x
X
u
X
X
u
X
X
u
X
X
u
X
0

#
#
#
#
#
#
#
#
#
#
$
%

=
x y x
x y x
x y x
x y x
(A.14)
For example,
if kl = 11 then,
l
H
k
i
kl
X
u
X

) ( ) , (
1
x y x
(i = 1, 2 )
1
1
2
11
2 11
2
2 1
1
2
11
2
1
1
1
11
2 11
2
1 1
1
1
11
2
1
1
2
11
1 11
1
2 1
1
2
11
1
1
1
1
11
1 11
1
1 1
1
1
11
1
) (
]
) (
) , ( ) (
) , (
[
) ( ) , (
) (
]
) (
) , ( ) (
) , (
[
) ( ) , (
) (
]
) (
) , ( ) (
) , (
[
) ( ) , (
) (
]
) (
) , ( ) (
) , (
[
) ( ) , (
X
u
X X X
u
X
X
u
X X X
u
X
X
u
X X X
u
X
X
u
X X X
u
X
H N
N N
N H
H N
N N
N H
H N
N N
N H
H N
N N
N H

+

=

+

=

+

=

+

=

x x
y x x
y x x y x
x x
y x x
y x x y x
x x
y x x
y x x y x
x x
y x x
y x x y x

8
8

8
8

8
8

8
8
(A.15)
131
#
#
#
#
#
#
#
#
#
#
$
%

+

+

+

+

=
1
1 11
2
2
2
2 2
2
2
1
1 11
2
1
2
2 1
2
2
1
1 11
1
2
2
2 2
2
2
1
1 11
1
1
2
2 1
2
2
) (
) ( ]
) (
) (
) , ( ) (
) (
) , (
[
) (
) ( ]
) (
) (
) , ( ) (
) (
) , (
[
) (
) ( ]
) (
) (
) , ( ) (
) (
) , (
[
) (
) ( ]
) (
) (
) , ( ) (
) (
) , (
[
X
u
X
hy
hy X
hx
hx
X
u
X
hy
hy X
hx
hx
X
u
X
hy
hy X
hx
hx
X
u
X
hy
hy X
hx
hx
H
N
N N
H
N
N N
H
N
N N
H
N
N N
x
x
x
x
y x x
x
y x
x
x
x
x
y x x
x
y x
x
x
x
x
y x x
x
y x
x
x
x
x
y x x
x
y x

8 8

8 8

8 8

8 8
1
1
2
2
2
11
2
2
2
2
11
2
1
2
2
11
2
1
2
2
11
2
2
2
2
11
1
2
2
2
11
1
1
2
2
11
1
1
2
2
11
1
) (
) (
) (
) ( ) (
) (
) (
) (
) (
) ( ) (
) (
) (
) (
) (
) ( ) (
) (
) (
) (
) (
) ( ) (
) (
) (
) , (
X
u
X
hy
hy X
hx
hx
X
hy
hy X
hx
hx
X
hy
hy X
hx
hx
X
hy
hy X
hx
hx
H
N N
N N
N N
N N
N

#
#
#
#
#
#
#
#
#
#
$
%

+

+

+

+

+
x
x
x
x x
x
x
x
x
x x
x
x
x
x
x x
x
x
x
x
x x
x
x
y x

8 (A.16)
Have 4x2 matrix form to 4x1 from the equation (A.14), then,
1 4
2
2
1
2
2
1
1
1
) ( ) , (
) ( ) , (
) ( ) , (
) ( ) , (
X
l
H
k
kl
l
H
k
kl
l
H
k
kl
l
H
k
kl
X
u
X
X
u
X
X
u
X
X
u
X
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

x y x
x y x
x y x
x y x

132
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

+

+

+

+

+

+

+

+

+

+

+

+

=
1
2
2
1
2
12
2
2
2
2
22
2
1
1
2
11
2
1
2
2
1
1
12
2
2
2
1
22
2
1
1
1
11
2
1
2
2
1
2
12
1
2
2
2
22
1
1
1
2
11
1
1
2
2
1
1
12
1
2
2
1
22
1
1
1
1
11
1
) ( ) ( ) , ( ) ( ) , ( ) ( ) , (
) ( ) ( ) , ( ) ( ) , ( ) ( ) , (
) ( ) ( ) , ( ) ( ) , ( ) ( ) , (
) ( ) ( ) , ( ) ( ) , ( ) ( ) , (
X
u
X
u
X X
u
X X
u
X
X
u
X
u
X X
u
X X
u
X
X
u
X
u
X X
u
X X
u
X
X
u
X
u
X X
u
X X
u
X
H H H H
H H H H
H H H H
H H H H
x x y x x y x x y x
x x y x x y x x y x
x x y x x y x x y x
x x y x x y x x y x

(A.17)
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

#
#
#
#
#
#
#
#
#
#
$
%

=
2
2
1
2
2
1
1
1
3 4
2
12
2
2
22
2
2
11
2
1
12
2
1
22
2
1
11
2
2
12
1
2
22
1
2
11
1
1
12
1
1
22
1
1
11
1
) (
) (
) (
) (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
X
u
X
u
X
u
X
u
A
X X X
X X X
X X X
X X X
H
H
H
H
X
x
x
x
x
y x y x y x
y x y x y x
y x y x y x
y x y x y x

(A.18)
since ) ( ) ( ) ( ) ( ) , ( ) , ( y x y x y x y x
x x x x x x
ij kl
i
ij N kl
i
ij N N kl
i
N kl
i
ij ij ij
8 8 = = =
= = =
[] {}
H
i
X
ij ij ij
ij ij ij
ij ij ij
ij ij ij
u DJBX A
X X X
X X X
X X X
X X X
)] ( [
) ( ) ( ) (
) ( ) ( ) (
) ( ) ( ) (
) ( ) ( ) (
3 4
2
12
2
2
22
2
2
11
2
1
12
2
1
22
2
1
11
2
2
12
1
2
22
1
2
11
1
1
12
1
1
22
1
1
11
1
3

#
#
#
#
#
#
#
#
#
#
$
%

=
y y y
y y y
y y y
y y y
(A.19)
133
Also, note that,
i
N kl
ij N N kl
i
i
ij N
i
ij kl
X X X

+

=

) (
) ( ) (
) ( ) (
1 1
x
y x
y y
8
8
#
$
%

+

+
#
$
%

+

=
i
ij
ij
N kl
i
ij
ij
N kl
ij N
N kl
i
i
ij
ij
ij N
i
ij
ij
ij N
X
hy
hy X
hx
hx
X
hy
hy X
hx
hx
) (
) (
) ( ) (
) (
) (
) (
) (
) (
) (
) ( ) (
) (
) (
2
2
1 2
2
1
2
2
2
2
x
x
x x
x
x
y
x
x
x
y x
x
y

8

8 8
(A.20)
Thus, the equation (A.19) can be expressed as following,
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
$
%

+

+

+

+

+

+

+

+

+

+

+

+

=
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
12
2
2
2
2
21
2
2 22
2
2
2
2
2
2
2 11
2
2
2
2
2
2
2
12
2
1
2
2
1
2
2 22
2
1
2
2
1
2
2 11
2
1
2
2
1
2
2
12
1
2
2
2
21
2
2 22
1
2
2
2
2
2
2 11
1
2
2
2
2
2
2
12
1
1
2
2
1
2
2 22
1
1
2
2
1
2
2 11
1
1
2
2
1
2
2
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
x
x
x
y
x
x
y
x
x
x
y
x
x
y
x
x
x
y
x
x
y
x
x
x
y
x
x
y
x
x
x
y
x
x
y
x
x
x
y
x
x
y
x
x
x
y
x
x
y
x
x
x
y
x
x
y
x
x
x
y
x
x
y
x
x
x
y
x
x
y
x
x
x
y
x
x
y
x
x
x
y
x
x
y

8
8

8
8

8
8

8
8

8
8

8
8

8
8

8
8

8
8

8
8

8
8

8
8
(A.21a)
134
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
$
%

+

+

+

+

+

+

+

+

+

+

+

+

+
2
2
2
12
2
2
2
2
12
2
2
2
2
22
2
2
2
2
22
2
2
2
2
11
2
2
2
2
11
2
1
2
2
12
2
1
2
2
12
2
1
2
2
22
2
1
2
2
22
2
1
2
2
11
2
1
2
2
11
2
2
2
2
12
1
2
2
2
12
1
2
2
2
22
1
2
2
2
22
1
2
2
2
11
1
2
2
2
11
1
1
2
2
12
1
1
2
2
12
1
1
2
2
22
1
1
2
2
22
1
1
2
2
11
1
1
2
2
11
1
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
) (
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
ij
ij
ij N
ij
ij
ij N
ij N
x
x
x
x
x
x
y
x
x
x
x
x
x
y
x
x
x
x
x
x
y
x
x
x
x
x
x
y
x
x
x
x
x
x
y
x
x
x
x
x
x
y
x
x
x
x
x
x
y
x
x
x
x
x
x
y
x
x
x
x
x
x
y
x
x
x
x
x
x
y
x
x
x
x
x
x
y
x
x
x
x
x
x
y

8

8

8

8

8

8

8

8

8

8

8

8
(A.21b)

And further more,
18 4
2
2
2
1
1
2
2
1
2
2
2
2
2
2
2
1
1
2
2
2
1
2
2
1
. . 0
) (
) (
) (
0
. . 0
) (
) (
) (
0
. .
) (
) (
) (
0
) (
) (
) (
. .
) (
) (
) (
0
) (
) (
) (
X
ij
ij
ij
ij
ij
ij
ij
ij
ij ij
ij
ij
ij
ij
ij ij
ij
ij
X
hx
hx
X
hx
hx
X
hx
hx X
hx
hx
X
hx
hx X
hx
hx
#
#
#
#
#
#
#
#
#
#
$
%

=
x
x
y
x
x
y
x
x
y x
x
y
x
x
y x
x
y
8
8
8 8
8 8
135
X
3 18
9 12
2
9 22
2
9 11
2
9 12
1
9 22
1
9 11
1
1 12
2
1 22
2
1 11
2
1 12
1
1 22
1
1 11
1
) ( ) ( ) (
) ( ) ( ) (
. . .
. . .
) ( ) ( ) (
) ( ) ( ) (
x
ij ij ij
ij ij ij
ij ij ij
ij ij ij
#
#
#
#
#
#
#
#
$
%

x x x
x x x
x x x
x x x

(A.22a)
18 4
2
2
2
1
1
2
2
1
2
2
2
2
2
2
2
1
1
2
2
2
1
2
2
1
. . 0
) (
) (
) (
0
. . 0
) (
) (
) (
0
. .
) (
) (
) (
0
) (
) (
) (
. .
) (
) (
) (
0
) (
) (
) (
X
ij
ij
ij
ij
ij
ij
ij
ij
ij ij
ij
ij
ij
ij
ij ij
ij
ij
X
hy
hy
X
hy
hy
X
hy
x hy
X
hy
hy
X
hy
hy X
hy
hy
#
#
#
#
#
#
#
#
#
#
$
%

+
x
x
y
x
x
y
x y x
x
y
x
x
y x
x
y
8
8
8 8
8 8
X
3 18
9 12
2
9 22
2
9 11
2
9 12
1
9 22
1
9 11
1
1 12
2
1 22
2
1 11
2
1 12
1
1 22
1
1 11
1
) ( ) ( ) (
) ( ) ( ) (
. . .
. . .
) ( ) ( ) (
) ( ) ( ) (
x
ij ij ij
ij ij ij
ij ij ij
ij ij ij
#
#
#
#
#
#
#
#
$
%

x x x
x x x
x x x
x x x

(A.22b)
+
36 4
1
1
1
2 1
. . 0 ) ( 0 0 0
. . 0 0 ) ( 0 0
. . 0 0 0 ) ( 0
. . ) ( 0 0 0 ) (
X
ij
ij
ij
ij ij
#
#
#
#
#
$
%

y
y
y
y y
8
8
8
8 8
136
X
3 36
2
2
2
1 12
2
2
2
2
1 22
2
2
2
2
1 11
2
1
2
2
1 12
2
1
2
2
1 22
2
1
2
2
1 11
2
2
2
2
1 12
1
2
2
2
1 22
1
2
2
2
1 11
1
1
2
2
1 12
1
1
2
2
1 22
1
1
2
2
1 11
1
. . .
. . .
. . .
) (
) (
) ( ) (
) (
) ( ) (
) (
) (
) (
) (
) ( ) (
) (
) ( ) (
) (
) (
) (
) (
) ( ) (
) (
) ( ) (
) (
) (
) (
) (
) ( ) (
) (
) ( ) (
) (
) (
X
ij
ij
ij ij
ij
ij ij
ij
ij
ij
ij
ij ij
ij
ij ij
ij
ij
ij
ij
ij ij
ij
ij ij
ij
ij
ij
ij
ij ij
ij
ij ij
ij
ij
X
hx
hx X
hx
hx X
hx
hx
X
hx
hx X
hx
hx X
hx
hx
X
hx
hx X
hx
hx X
hx
hx
X
hx
hx X
hx
hx X
hx
hx
#
#
#
#
#
#
#
#
#
#
#
#
#
#
$
%

x
x
x x
x
x x
x
x
x
x
x x
x
x x
x
x
x
x
x x
x
x x
x
x
x
x
x x
x
x x
x
x

(A.22c)
+
36 4
1
1
1
2 1
. . 0 ) ( 0 0 0
. . 0 0 ) ( 0 0
. . 0 0 0 ) ( 0
. . ) ( 0 0 0 ) (
X
ij
ij
ij
ij ij
#
#
#
#
#
$
%

y
y
y
y y
8
8
8
8 8
X
3 36
2
2
2
1 12
2
2
2
2
1 22
2
2
2
2
1 11
2
1
2
2
1 12
2
1
2
2
1 22
2
1
2
2
1 11
2
2
2
2
1 12
1
2
2
2
1 22
1
2
2
2
1 11
1
1
2
2
1 12
1
1
2
2
1 22
1
1
2
2
1 11
1
. . .
. . .
. . .
) (
) (
) ( ) (
) (
) ( ) (
) (
) (
) (
) (
) ( ) (
) (
) ( ) (
) (
) (
) (
) (
) ( ) (
) (
) ( ) (
) (
) (
) (
) (
) ( ) (
) (
) ( ) (
) (
) (
X
ij
ij
ij ij
ij
ij ij
ij
ij
ij
ij
ij ij
ij
ij ij
ij
ij
ij
ij
ij ij
ij
ij ij
ij
ij
ij
ij
ij ij
ij
ij ij
ij
ij
X
hy
hy X
hy
hy X
hy
hy
X
hy
hy X
hy
hy X
hx
hy
X
hy
hy X
hy
hy X
hy
hy
X
hy
hy X
hy
hy X
hy
hy
#
#
#
#
#
#
#
#
#
#
#
#
#
#
$
%

x
x
x x
x
x x
x
x
x
x
x x
x
x x
x
x
x
x
x x
x
x x
x
x
x
x
x x
x
x x
x
x

(A.22d)
137
) ( ) (
),
2
(
) ( ) (
),
2
(
) , ( ) , (
2 2
2
1
2 2
2
1
2 1
ij ij
ij ij
x hy x hy
hy
Y
x hx x hx
hx
Y
Y Y
2 2
2
3 3
3
2 3
, =

=
, =

=
= = y
(A.23)
)
) (
)( 1 )( 2 (
) (
) (
)
) (
)( 1 ( ) 1 2 (
2
1
) (
) (
)
) (
)( 1 ( ) 2 (
2
1
) (
) (
)
) (
)( 1 ( ) 1 2 (
2
1
) (
) (
)
) (
)( 1 ( ) 2 (
2
1
) (
) (
)
) (
)( 1 ( ) 1 2 (
4
1
) (
) (
)
) (
)( 1 ( ) 1 2 (
4
1
) (
) (
)
) (
)( 1 ( ) 1 2 (
4
1
) (
) (
)
) (
)( 1 ( ) 1 2 (
4
1
) (
) (
) (
) (
) (
) (
2
2
2
9
2
2
2
8
2 2
7
2
2
2
6
2 2
5
2 2
4
2 2
3
2 2
2
2 2
1
2
1
ij ij
ij
ij ij
ij
ij ij
ij
ij ij
ij
ij ij
ij
ij ij
ij
ij ij
ij
ij ij
ij
ij ij ij
ij
ij
ij
hx hx
hx hx
hx hx
hx hx
hx hx
hx hx
hx hx
hx hx
hx hx hx
x x
y
x x
y
x x
y
x x
y
x x
y
x x
y
x x
y
x x
y
x x
y
x
y
x
y
3
2 3
8
3
2 2 3
8
3
2 2 3
8
3
2 2 3
8
3
2 2 3
8
3
2 2 3
8
3
2 2 3
8
3
2 2 3
8
3
2 2 3
3
3
8 8
,
, , =

,
, , =

,
+ , =

,
, + =

,
, , =

,
+ , =

,
+ + =

,
, + =

,
, , =

=

(A.24)
138
)
) (
)( 2 )( 1 (
) (
) (
)
) (
)( 2 )( 1 (
2
1
) (
) (
)
) (
)( 1 2 )( 1 (
2
1
) (
) (
)
) (
)( 2 )( 1 (
2
1
) (
) (
)
) (
)( 1 2 )( 1 (
2
1
) (
) (
)
) (
)( 1 2 )( 1 (
4
1
) (
) (
)
) (
)( 1 2 )( 1 (
4
1
) (
) (
)
) (
)( 1 2 )( 1 (
4
1
) (
) (
)
) (
)( 1 2 )( 1 (
4
1
) (
) (
) (
) (
) (
) (
2
2
2
9
2 2
8
2
2
2
7
2 2
6
2
2
2
5
2 2
4
2 2
3
2 2
2
2 2
1
2
1
ij ij
ij
ij ij
ij
ij ij
ij
ij ij
ij
ij ij
ij
ij ij
ij
ij ij
ij
ij ij
ij
ij ij ij
ij
ij
ij
hy hx
hy hx
hy hx
hy hx
hy hx
hy hx
hy hx
hy hx
hy hy hy
x x
y
x x
y
x x
y
x x
y
x x
y
x x
y
x x
y
x x
y
x x
y
x
y
x
y
2
2 3
8
2
2 3 3
8
2
2 3
8
2
2 3 3
8
2
2 3
8
2
2 3 3
8
2
2 3 3
8
2
2 3 3
8
2
2 3 3
2
2
8 8
,
, , =

,
, , =

,
+ , =

,
, + =

,
, , =

,
+ , =

,
+ + =

,
, + =

,
, , =

=

(A.25)
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

+

+

+

+

+

+

+

+

+

+

+

+

=
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

....
) ( ) ( ) (
....
) ( ) ( ) (
....
) ( ) ( ) (
....
) ( ) ( ) (
) (
) (
) (
) (
3
2
2
3
2
2
2
2
1
2
2
1
3
2
1
3
2
2
1
2
1
2
1
1
3
2
2
3
2
2
2
2
1
2
2
1
3
2
1
3
2
2
1
2
1
2
1
1
2
2
1
2
2
2
1
2
hy
X
hy
X
hy
X
hy
X
hy
X
hy
X
hx
X
hx
X
hx
X
hx
X
hx
X
hx
X
X
hy
X
hy
X
hx
X
hx
ij ij ij
ij ij ij
ij ij ij
ij ij ij
ij
ij
ij
ij
x x x
x x x
x x x
x x x
x
x
x
x
8 8 8
8 8 8
8 8 8
8 8 8
(A.26)
139
[] [] {}
N
i X X
ij N
hx DJBX
hy
hx
hy
hx
hy
hx
X
X
X
X
18 4
9
2
9
2
2
2
2
2
1
2
1
2
18 2
2
1
2
1
) (
.
.
. ) (
0
0
0
0
3 0 =
&
&
&
&
&
&
'
&
&
&
&
&
&
(
)
&
&
&
&
&
&

&
&
&
&
&
&

#
#
#
#
#
#
#
#
#
$
%

= y (A.27)
Also,
2
) (
11
1 ) (
11
1
2
11
1 2 2 2
) ( ) (
) (
) (
hx hx
hx
N
hx hx
N
N
9
,
=

9 + x x
x x
x
x

(A.28)
2
) (
11
1 ) (
11
1
2
11
1 2 2 2
) ( ) (
) (
) (
hy hy
hy
N
hy hy
N
N
9
,
=

9 + x x
x x
x
x

(A.29)
[] {}
H
i
X
ij ij ij
ij ij ij
ij ij ij
ij ij ij
u DJBX A
X X X
X X X
X X X
X X X
)] ( [
) ( ) ( ) (
) ( ) ( ) (
) ( ) ( ) (
) ( ) ( ) (
3 4
2
12
2
2
22
2
2
11
2
1
12
2
1
22
2
1
11
2
2
12
1
2
22
1
2
11
1
1
12
1
1
22
1
1
11
1
3

#
#
#
#
#
#
#
#
#
#
$
%

y y y
y y y
y y y
y y y
140
[] []
[] {}
H
i
X
ij
ij
ij
kl
X
ij N
X
ij
ij
ij
kl
X
ij N
X
ij
N kl
i
X
ij
ij
ij
X
ij
N kl
i
X
ij
ij
ij
u DJBX A
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
)] ( [ x
) (
) (
) (
) , (
) (
) (
) (
) , (
) (
) (
) (
) , , (
) (
) (
) (
) , , (
3 36
1
2
2
N
1
36 4
3 36
1
2
2
N
1
36 4
3 18
18 4
1
2
2
1
3 18
18 4
1
2
2
1
3
3
3
3
2 3 8
3
3
3
2 3 8
3
3
3
2 3 8
3
3
3
2 3 8

#
#
#
$
%

+
#
#
#
$
%

+
#
$
%

#
#
#
$
%

+
#
$
%

#
#
#
$
%

=
x x
(A.30)
Therefore, the equation A.14 becomes:
Equation (A.14) =
[]
[]
[] {}
{}
1 18
18 2
3 36
1
2
2
N
1
36 4
3 36
1
2
2
N
1
36 4
3 18
18 4
1
2
2
1
3 18
18 4
1
2
2
1
) ( x
)] ( [
) (
) (
) (
) , (
) (
) (
) (
) , (
) (
) (
) (
) (
) (
) (
x
N
i
x
N
H
i
X
ij
ij
ij
kl
X
ij N
X
ij
ij
ij
kl
X
ij N
X
ij
N kl
i
X
ij
ij
ij
X
ij
N kl
i
X
ij
ij
ij
C
u DJBX A
X
hy
hy
X
hx
hx
X
hy
hy
X
hx
hx
#
$
%

&
&
&
&
&
&
&
&
'
&
&
&
&
&
&
&
&
(
)
&
&
&
&
&
&
&
&

&
&
&
&
&
&
&
&

#
#
$
%

+
#
#
$
%

+
#
$
%

#
#
$
%

+
#
$
%

#
#
$
%

3 0
3
3
3
3
2 3 8
3
3
3
2 3 8
3
3
3
8
3
3
3
8
(A.31)
141
For the equation (3.24),
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

N N
l
H
k kl
N N
l
H
k kl
N N
l
H
k kl
N N
l
H
k kl
C
X X
u
C
X X
u
C
X X
u
C
X X
u
2
2
2
2
2
1
2
2
1
2
2
1
1
1
2
1
) (
) (
) , (
) (
) (
) , (
) (
) (
) , (
) (
) (
) , (
x
x
y x
x
x
y x
x
x
y x
x
x
y x
0
0
0
0
{}
1 18
18 2
1 1
2 1
2 4
2
2
2
1
2
2
2
2
1
1
2
1
. . . ) ( 0 ) ( 0
. . . 0 ) ( 0 ) (
x
) (
) , ( 0
) (
) , ( 0
0
) (
) , (
0
) (
) , (
X
N
i
X
X
l
H
k kl
l
H
k kl
l
H
k kl
l
H
k kl
C
X X
u
X X
u
X X
u
X X
u
#
$
%

#
#
#
#
#
#
#
#
#
#
$
%

=
x x
x x
x
y x
x
y x
x
y x
x
y x
0 0
0 0

(A.32)
Assemble the matrix form,
142
&
&
'
&
&
(
)
&
&

&
&

+

+

+

+

+

+

=
&
&
'
&
&
(
)
&
&

&
&

1
2
2
1
1
12
2
2
2
1
22
2
1
1
1
11
2
1
2
2
1
1
12
1
2
2
1
22
1
1
1
1
11
1
1
2
2
1
2
1
) ( ) (
) , (
) (
) , (
) (
) , (
) ( ) (
) , (
) (
) , (
) (
) , (
) (
) , (
) (
) , (
X
u
X
u
X X
u
X X
u
X
X
u
X
u
X X
u
X X
u
X
X X
u
X X
u
H H H H
H H H H
l
H
k kl
l
H
k kl
x x
y x
x
y x
x
y x
x x
y x
x
y x
x
y x
x
y x
x
y x

[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

#
$
%

=
2
2
1
1
2
1
2
1
1
1
1
1
12
2
22
2
11
2
12
1
22
1
11
1
) (
) (
) (
) (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
X
u
X
X
u
X
X
u
X
X
u
X
A
H
H
H
H
x
x
x
x
y x y x y x
y x y x y x

[] []
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

#
$
%

=
2
2
1
1
2
1
2
1
1
1
1
1
3 18
11
18 2
) (
) (
) (
) (
) (
X
u
X
X
u
X
X
u
X
X
u
X
A
H
H
H
H
x
N
i x
ij
x
x
x
x
3 0 (A.33)
143
Also,
&
&
'
&
&
(
)
&
&

&
&

+

+

+

+

+

+

=
&
&
'
&
&
(
)
&
&

&
&

1
2
2
1
2
12
2
2
2
2
22
2
1
1
2
11
2
1
2
2
1
2
12
1
2
2
2
22
1
1
1
2
11
1
2
2
2
2
2
1
) ( ) (
) , (
) (
) , (
) (
) , (
) ( ) (
) , (
) (
) , (
) (
) , (
) (
) , (
) (
) , (
X
u
X
u
X X
u
X X
u
X
X
u
X
u
X X
u
X X
u
X
X X
u
X X
u
H H H H
H H H H
l
H
k kl
l
H
k kl
x x
y x
x
y x
x
y x
x x
y x
x
y x
x
y x
x
y x
x
y x

(A.34)
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

#
$
%

=
2
2
2
1
2
2
2
1
2
1
1
2
12
2
22
2
11
2
12
1
22
1
11
1
) (
) (
) (
) (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
X
u
X
X
u
X
X
u
X
X
u
X
A
H
H
H
H
x
x
x
x
y x y x y x
y x y x y x

[] []
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

#
$
%

=
2
2
2
1
2
2
2
1
2
1
1
2
3 18
11
18 2
) (
) (
) (
) (
) (
X
u
X
X
u
X
X
u
X
X
u
X
A
H
H
H
H
x
N
i x
ij
x
x
x
x
3 0 (A.35)
144
[] {}
1 18 18 2
2
1
2
1
1
2
2
1
1
2
1
2
1
1
1
1
1
) (
0
0
0
0
) (
) (
) (
) (
X
H
i X
H
H
H
H
u
X
X
X
X
X
X
u
X
X
u
X
X
u
X
X
u
X
x
x
x
x
x
0
#
#
#
#
#
#
#
#
#
$
%

=
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

[]
[]
[] {}
1 18
18 2
2
1
2
1
1
1
1
) (
0
0
0
0
0
0
X
H
i
x
N
u
J
J
X

#
#
#
#
#
#
#
#
#
#
$
%

#
$
%

=
,
,
3 0
3
3
3
3
(A.36)
[] {}
1 18 18 2
2
1
2
1
2
2
2
2
1
2
2
2
1
2
1
1
2
) (
0
0
0
0
) (
) (
) (
) (
X
H
i X
H
H
H
H
u
X
X
X
X
X
X
u
X
X
u
X
X
u
X
X
u
X
x
x
x
x
x
0
#
#
#
#
#
#
#
#
#
$
%

=
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

145
[]
[]
[] {}
1 18
18 2
2
1
2
1
1
1
2
) (
0
0
0
0
0
0
X
H
i
x
N
u
J
J
X

#
#
#
#
#
#
#
#
#
#
$
%

#
$
%

=
,
,
3 0
3
3
3
3
(A.37)
In order to compute the
1
X

and
2
X

, we will get from the following equation:
#
#
#
#
#
#
#
#
#
#
#
$
%

#
#
#
#
$
%

=
#
#
#
#
#
#
#
#
#
$
%

2
1
2
1
22 21
12 11
22 21
12 11
2
1
2
1
0
0
0
0
0 0
0 0
0 0
0 0
0
0
0
0
3
3
3
3
X
X
X
X
#
#
#
#
#
#
#
#
#
#
$
%

+

+

+

+

=
2
22
1
21
2
12
1
11
2
22
1
21
2
12
1
11
0
0
0
0
3 3
3 3
3 3
3 3
(A.38)
146
&
&

&
&

+

=

+

=

?
2
22
1
21
2
2
12
1
11
1
3 3
3 3
X
X
(A.39)
where
] [
1
], [
1
], [
1
], [
1
11 22 21 21 12 12 22 11
J
J
J
J
J
J
J
J
= , = , = = and ] det[J J = (A.40)
[]
[]
[] {}
1 18
18 2
2
1
2
1
1
1
1
) (
0
0
0
0
0
0
X
H
i
x
N
u
J
J
X

#
#
#
#
#
#
#
#
#
#
$
%

#
$
%

,
,
3 0
3
3
3
3
[]
[]
[] {}
1 18
18 2
2
1
2
1
1
1
2
12
1
11
) (
0
0
0
0
0
0
X
H
i
x
N
u
J
J

#
#
#
#
#
#
#
#
#
#
$
%

#
$
%

+

=
,
,
3 0
3
3
3
3
3 3
(A.41)
147
[]
[]
{}
1 18
18 2
2
2
2
12
2 1
2
11
2 1
2
12
2
1
2
11
2
2
2
12
2 1
2
11
2 1
2
12
2
1
2
11
1
1
) (
0
0
0
0
0
0
X
H
i
x
N
u
J
J
#
$
%

#
#
#
#
#
#
#
#
#
#
#
$
%

+

+

+

+

#
$
%

=
,
,
3 0
3 3 3
3 3 3
3 3 3
3 3 3
(A.42)
[ ]
[]
#
$
%

=
,
,
1
1
0
0
J
J
X
[] [] {}
1 18
18 2
2
2
2
12
2 1
2
12
2
2
2
12
2 1
2
12
18 2
2 1
2
11
2
1
2
11
2 1
2
11
2
1
2
11
) (
0
0
0
0
) (
0
0
0
0
X
H
i
x
N
x
N
u

#
#
#
#
#
#
#
#
#
#
#
$
%

+
#
#
#
#
#
#
#
#
#
#
#
$
%

3 0
3
3 3
3
3 3
3 0
3 3
3
3 3
3
(A.43)
[] ()
[]
[]
{}
1 18
18 2
2 4
2
12
18 2
2 4
2
11
4 4
1
) (
) (
X
H
i
x
N
X
i j
x
N
X
j i
X
u J

+

=
,
3 0
3 3
3 0
3 3
(A.44)
148
[]
[]
{}
1 18
18 2
2
1
2
1
1
1
1
) (
0
0
0
0
0
0
X
H
i
x
N
u
J
J
X

#
$
%

#
#
#
#
#
#
#
#
#
#
#
$
%

#
$
%

?
,
,
3 0
3
3
3
3
[] ()
[ ]
[]
{}
1 18
18 2
2 4
2
12
18 2
2 4
2
11
4 4
1
) (
) (
X
H
i
x
N
X
i j
x
N
X
j i
X
u J

+

=
,
3 0
3 3
3 0
3 3
(A.45)
[]
[]
{}
1 18
18 2
2
1
2
1
1
1
2
) (
0
0
0
0
0
0
X
H
i
x
N
u
J
J
X

#
$
%

#
#
#
#
#
#
#
#
#
#
#
$
%

#
$
%

?
,
,
3 0
3
3
3
3
[] ()
[ ]
[]
{}
1 18
18 2
2 4
2
22
18 2
2 4
2
21
4 4
1
) (
) (
X
H
i
x
N
X
i j
x
N
X
j i
X
u J

+

=
,
3 0
3 3
3 0
3 3
(A.46)
149
[] [][] () {}
[] [][] () {}
{}
1 18 18 2
2 4
1 18
18 2
2 4
2
22
18 2
2 4
2
21
4 4
1
3 18
11
18 2
1 18
18 2
2 4
2
12
18 2
2 4
2
11
4 4
1
3 18
11
18 2
) ( x
) (
) (
) (
) (
) (
) (
) 17 (
X
N
i
X
N
X
X
H
i
x
N
X
i j
x
N
X
j i
X
x
N
i x
X
H
i
x
N
X
i j
x
N
X
j i
X
x
N
i x
C
u J A
u J A
eq
ij
ij
3 0
3 0
3 3
3 0
3 3
3 0
3 0
3 3
3 0
3 3
3 0
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
#
$
%

#
$
%

+
#
$
%

#
$
%

#
$
%

+
#
$
%

#
$
%

=
,
,
(A.47)
For the equation (3.25),
150
{} [ ]
1 18 18 4
2
2
2
1
1
2
1
1
) (
) (
) (
) (
X
N
i X
N
N
N
N
N
N
N
N
u DJBX
u
X
u
X
u
X
u
X
=
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

x
x
x
x
0
0
0
0
(A.48)
For the equation (3.26),
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

2
2
2
2
2
1
1
2
1
2
2
2
1
2
1
1
1
1
1
1
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
X
Y
C
X
u
Y
X
Y
C
X
u
Y
X
Y
C
X
u
Y
X
Y
C
X
u
Y
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
x
x y x
x
x y x
x
x y x
x
x y x
0

0

0

0

(A.49)
[]{}
1 18 18 2
2
2
2
2
1
1
1
2
2
2
2
1
1
1
1
1
) (
) ( ) , (
0
) ( ) , (
0
0
) ( ) , (
0
) ( ) , (
X
i
X
l
H
k
kl
l
H
k
kl
l
H
k
kl
l
H
k
kl
C
X
Y
X
u
Y
X
Y
X
u
Y
X
Y
X
u
Y
X
Y
X
u
Y
x
x y x
x y x
x y x
x y x
0

#
#
#
#
#
#
#
#
#
#
$
%

= (A.50)
151
From the equation A.50, rearrange the matrix from 4x2 to 4x1 format. Then,
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

2
2
2
2
1
1
1
2
2
2
2
1
1
1
1
1
) ( ) , (
) ( ) , (
) ( ) , (
) ( ) , (
X
Y
X
u
Y
X
Y
X
u
Y
X
Y
X
u
Y
X
Y
X
u
Y
l
H
k
kl
l
H
k
kl
l
H
k
kl
l
H
k
kl
x y x
x y x
x y x
x y x

#
#
#
#
#
#
#
#
#
#
$
%

=
l
H
k
kl
l
H
k
kl
l
H
k
kl
l
H
k
kl
X
u
Y
X
u
Y
X
u
Y
X
u
Y
) ( ) , (
0 0 0
0
) ( ) , (
0 0
0 0
) ( ) , (
0
0 0 0
) ( ) , (
2
2
1
2
2
1
1
1
x y x
x y x
x y x
x y x

&
&
&
&
'
&
&
&
&
(
)
&
&
&
&

&
&
&
&

2
2
1
1
2
2
1
1
x
X
Y
X
Y
X
Y
X
Y
(A.51)
Reformat the equation A.51,
[]
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

#
#
#
#
#
#
#
#
#
#
$
%

=
#
#
#
#
#
#
#
#
#
#
$
%

2
2
1
2
2
1
1
1
2
12
2
2
22
2
2
11
2
1
12
2
1
22
2
1
11
2
2
12
1
2
22
1
2
11
1
1
12
1
1
22
1
1
11
1
2
2
1
2
2
1
1
1
) (
) (
) (
) (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
) ( ) , (
) ( ) , (
) ( ) , (
) ( ) , (
X
u
X
u
X
u
X
u
A
Y Y Y
Y Y Y
Y Y Y
Y Y Y
X
u
Y
X
u
Y
X
u
Y
X
u
Y
H
H
H
H
l
H
k
kl
l
H
k
kl
l
H
k
kl
l
H
k
kl
x
x
x
x
y x y x y x
y x y x y x
y x y x y x
y x y x y x
x y x
x y x
x y x
x y x

(A.52)
152
[][ ] {}
1 18 18 4
2
12
2
2
22
2
2
11
2
1
12
2
1
22
2
1
11
2
2
12
1
2
22
1
2
11
1
1
12
1
1
22
1
1
11
1
) , ( ) , ( ) , (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
) , ( ) , ( ) , (
X
H
i X
u DJBX A
Y Y Y
Y Y Y
Y Y Y
Y Y Y
#
#
#
#
#
#
#
#
#
#
$
%

=
y x y x y x
y x y x y x
y x y x y x
y x y x y x

(A.53)
[][ ] {}
1 18 18 4
3 18
9 12
2
9 22
2
9 11
2
9 12
1
9 22
1
9 11
1
1 12
2
1 22
2
1 11
2
1 12
1
1 22
1
1 11
1
18 4
2
2
2
1
1
2
1
1
2
2
2
1
1
2
1
1
x
) ( ) ( ) (
) ( ) ( ) (
. . .
. . .
) ( ) ( ) (
) ( ) ( ) (
...
) (
0
) (
0
...
) (
0
) (
0
... 0
) (
0
) (
... 0
) (
0
) (
X
H
i X
x
ij ij ij
ij ij ij
ij ij ij
ij ij ij
X
u DJBX A
Y Y
Y Y
Y Y
Y Y
#
#
#
#
#
#
#
#
$
%

#
#
#
#
#
#
#
#
#
#
$
%

=
x x x
x x x
x x x
x x x
y y
y y
y y
y y

8 8
8 8
8 8
8 8
(A.54)
[] [ ] [][ ] {}
1 18 18 4 3 18 1 18 2
2
1
2
1
) ( ) (
0
0
0
0
X
H
i X X
ij N kl
X
N
u DJBX A
Y
Y
Y
Y

#
#
#
#
#
#
#
#
#
$
%

= x y 0 (A.55)
[ ] [ ] ( ) [ ][ ] { } ) (
1 18 18 4 3 18 1 18 4 X
H
i X X
ij N kl
X
u DJBX A DJBY x =
(A.56)
153
&
&

&
&

+
= =

+
= =

=
=
2
1 2
2
2
2
2
1
1 2
1
1
1
1
2
2
X
hy hy
X
Y
X
Y
X
hx hx
X
Y
X
Y
ij
ij
x
x
x
x
(A.57)
&
&
&
&
&
'
&
&
&
&
&
(
)
&
&
&
&
&

&
&
&
&
&

2
2
2
2
2
1
1
2
1
2
2
2
1
2
1
1
1
1
1
1
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
) (
) ( ) , (
X
Y
C
X
u
Y
X
Y
C
X
u
Y
X
Y
C
X
u
Y
X
Y
C
X
u
Y
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
N N
l
H
k
kl
x
x y x
x
x y x
x
x y x
x
x y x
0

0

0

0

[][] ()[][ ] {} {} { }
1 18 18 2
1 4
1 18 18 4 3 18 1 18 4
) ( ) (
X
N
i
X
N
X
i
i
X
H
i X X
ij N kl
X
C
X
Y
u DJBX A DJBY 3 0
#
#
$
%

'
(
)

= x
(A.58)

Abstract (if available)

Abstract Functionally Graded Materials (FGMs) have a gradual material variation from one material character to another throughout the structure. Applications of these types of materials have significant advantages in civil and mechanical systems including thermal systems. Analyzing the FGMs at the microstructure level with the conventional Finite Element Method (FEM) takes enormous pre-processing and computational time due to the complex material characteristics at the microstructure level. Essentially, the model contains too many degrees of freedom to be solved economically.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Homogenization procedures for the constitutive material modeling and analysis of aperiodic micro-structures

PDF

Studies into vibration-signature-based methods for system identification, damage detection and health monitoring of civil infrastructures

PDF

Modeling of multiscale continuum-atomistic systems using homogenization theory

PDF

Analytical and experimental studies of modeling and monitoring uncertain nonlinear systems

PDF

Characterization of environmental variability in identified dynamic properties of a soil-foundation-structure system

PDF

Novel multi-stage and CTLS-based model updating methods and real-time neural network-based semiactive model predictive control algorithms

PDF

Towards an optimal design of feedback strategies for structures controlled with smart dampers

PDF

Analytical and experimental methods for regional earthquake spectra and structural health monitoring applications

PDF

Analytical and experimental studies in modeling and monitoring of uncertain nonlinear systems using data-driven reduced‐order models

PDF

Constitutive model of concrete confined by advanced fiber composite materials and applications in seismic retrofitting

PDF

Finite element model to understand the role of fingerprints in generating vibrations

PDF

Analytical and experimental studies in system identification and modeling for structural control and health monitoring

PDF

Modeling, analysis and experimental validation of flexible rotor systems with water-lubricated rubber bearings

PDF

Development of composite oriented strand board and structures

PDF

Efficient inverse analysis with dynamic and stochastic reductions for large-scale models of multi-component systems

PDF

Design optimization under uncertainty for rotor blades of horizontal axis wind turbines

PDF

Studies into computational intelligence approaches for the identification of complex nonlinear systems

PDF

Large-scale molecular dynamics simulations of nano-structured materials

PDF

Numerical study of shock-wave/turbulent boundary layer interactions on flexible and rigid panels with wall-modeled large-eddy simulations

PDF

A simplified building energy simulation tool: material and environmental properties effects on HVAC performance

Asset Metadata

Creator Rhee, Richard Sangwon (author)

Core Title Multi-scale modeling of functionally graded materials (FGMs) using finite element methods

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Civil Engineering (Structural Mechanics)

Publication Date 08/21/2009

Defense Date 08/07/2007

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag coupling micro-macro structural models,finite element method,functionally graded materials,nonperiodic microstructure,OAI-PMH Harvest

Language English

Advisor Wellford, Carter (committee chair), Anderson, James C. (committee member), Lee, Vincent W. (committee member), Masri, Sami F. (committee member), Shiflett, Geoffrey R. (committee member)

Creator Email richardrhe@hotmail.com

Permanent Link (DOI) https://doi.org/10.25549/usctheses-m795

Unique identifier UC1188327

Identifier etd-Rhee-20070821 (filename),usctheses-m40 (legacy collection record id),usctheses-c127-542891 (legacy record id),usctheses-m795 (legacy record id)

Legacy Identifier etd-Rhee-20070821.pdf

Dmrecord 542891

Document Type Dissertation

Rights Rhee, Richard Sangwon

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Repository Name Libraries, University of Southern California

Repository Location Los Angeles, California

Repository Email cisadmin@lib.usc.edu