|USFDC Home | USF Electronic Theses and Dissertations||| RSS|
This item is only available as the following downloads:
xml version 1.0 encoding UTF-8 standalone no
record xmlns http:www.loc.govMARC21slim xmlns:xsi http:www.w3.org2001XMLSchema-instance xsi:schemaLocation http:www.loc.govstandardsmarcxmlschemaMARC21slim.xsd
leader nam Ka
controlfield tag 001 001925260
007 cr mnu|||uuuuu
008 080130s2007 flu sbm 000 0 eng d
datafield ind1 8 ind2 024
subfield code a E14-SFE0002163
RKEM implementation for strain gradient theory in multiple dimensions
h [electronic resource] /
by Abhishek Kumar.
[Tampa, Fla.] :
b University of South Florida,
ABSTRACT: The Reproducing Kernel Element Method (RKEM) implementation of the Fleck-Hutchinson phenomenological strain gradient theory in 1D, 2D and 3D is implemented in this research. Fleck-Hutchinson theory fits within the framework of Touplin- Mindlin theories and deals with first order strain gradients and associated work conjugate higher-order stress. Theories at the intrinsic or material length scales find applications in size dependent phenomena. In elasticity, length scale enters the constitutive equation through the elastic strain energy function which depends on both strain as well as the gradient of rotation and stress. The displacement formulation of the Touplin Mindlin theory involve diffrential equations of the fourth order, in conventional finite element method C1 elements are required to solve such equations, however C1 elements are cumbersome in 2D and unknown in 3D. The high computational cost and large number of degrees of freedom soon place such a formulation beyond the realm of practicality. Recently, some mixed and hybrid formulations have developed which require only C0 continuity but none of these elements solve complicated geometry problems in 2D and there is no problem yet solved in 3D. The large number of degrees of freedom is still inevitable for these formulation. As has been demonstrated earlier RKEM has the potential to solve higher-order problems, due to its global smoothness and interpolation properties. This method has the promise to solve important problems formulated with higher order derivatives, such as the strain gradient theory.
Thesis (M.S.)--University of South Florida, 2007.
Includes bibliographical references.
Text (Electronic thesis) in PDF format.
System requirements: World Wide Web browser and PDF reader.
Mode of access: World Wide Web.
Title from PDF of title page.
Document formatted into pages; contains 40 pages.
Advisor: Daniel Simkins, Ph.D.
Higher order diffrential equation.
t USF Electronic Theses and Dissertations.
RKEM Implementation for Strain Gradient Theory in Multiple Dimensions by Abhishek K umar A thesis submitted in partial fulllment of the requirements for the de gree of Master of Science in Engineering Department of Ci vil and En vironmental Engineering Colle ge of Engineering Uni v ersity of South Florida Major Professor: Daniel Simkins,Ph.D. Andres T ejada Martinez, Ph.D. Muhammad Mustazur Rahman, Ph.D. Date of Appro v al: July 11, 2007 K e yw ords: couple stress, bimaterial, higher order dif frential equation, messless, isogeometric c r Cop yright 2007, Abhishek K umar
Dedication This w ork is dedicated to my mother Late Usha Sharma, whose presence i al w ays felt in whate v er i did
Ackno wledgments I w ould lik e to thank Prof. Simkins for gi ving me an opportunity to w ork under his guidance. Prof. Simkins has played a v ery important role in helping me with my thesis w ork, without his time endless encouragement and patience this w ork w ould not ha v e been possible. I thank Prof. T ezada and Prof Rahman for time the y took to re vie w this thesis and there v aluable suggestion to impro v e it. I am indebted to the support and encouragement I ha v e got from my f amily friends and colleagues.
T able of Contents List of Figures ii Abstract iii Chapter 1 Introduction 1 Chapter 2 Re vie w of Linear Elastic Strain Gradient Theory 4 Chapter 3 Re vie w of RKEM 7 3.1 Concept of RKEM 7 3.2 Global P artition Polynomial 10 3.2.1 The L4P3I1 Element 10 3.2.2 The T9P2I1 Element 12 3.3 Salient Features of RKEM 13 Chapter 4 Galerkin F ormulation for Strain Gradient Problems 16 Chapter 5 1D Examples 21 5.1 One-dimensional T oupin-Mindlin Strain Gradient Theory 21 5.2 Shear Layer Problem with the T oupin-Mindlin Theory 22 5.2.1 Model Problem 24 5.2.2 Con v er gence Study 24 Chapter 6 2D Examples 26 6.1 Numerical Examples in T w o Dimension 26 6.1.1 Boundary Layer Analysis 26 6.1.2 Analytical Solution 26 6.2 An Innite Plate W ith a Hole 33 Chapter 7 Conclusions 38 References 39 i
List of Figures Figure 1. The T9P2I1 element with v ariable associated at each nodes 11 Figure 2. The parent triangle for T9P2I1 element 11 Figure 3. The global shape function of T9P2I1 element (a) (00)1 (b) (10) 1 (c) (01) 1 14 Figure 4. Shear layer attached on the left side (x =0) with traction acting on the right side (x=L) 23 Figure 5. Exact solution vs RKEM solution for T oupin-Mindlin shear layer model problem 24 Figure 6. Error plot for T oupin-Mindlin shear layer model problem 25 Figure 7. Con v er gence rate for RKEM L2P3I2 element 25 Figure 8. Geometry of a bimaterial under uniform shear 27 Figure 9. Mesh for b undary layer analysis 29 Figure 10. Shear strain plot by strain gradient theory 30 Figure 11. Shear strain plot by con v entional theory 31 Figure 12. Comparison of strain gradient vs con v entional theory 32 Figure 13. Comparison of e xact solution vs RKEM solution 32 Figure 14. Con v er gence plot for RKEM solution) 33 Figure 15. Notation and geometry of an innite plate subjected to bi-axial remote tension 34 Figure 16. mesh for plate problem 35 Figure 17. V ariation of u r for the plate with a hole 36 Figure 18. V ariation of r r for the plate with a hole 37 Figure 19. V ariation of for the plate with a hole 37 ii
RKEM Implementation for Strain Gradient Theory in Multiple Dimensions Abhishek K umar ABSTRA CT The Reproducing K ernel Element Method (RKEM) implementation of the Fleck-Hutchinson phenomenological strain gradient theory in 1D, and higher dimension is implemented in this research. Fleck-Hutchinson theory ts within the frame w ork of T ouplin-Mindlin theories and deals with rst order strain gradients and associated w ork conjugate higher -order stress. Theories at the intrinsic or material length scales nd applications in size dependent phenomena. In elasticity length scale enters the constituti v e equation through the elastic strain ener gy function which depends on both strain as well as the gradient of rotation and stress. The displacement formulation of the T ouplin Mindlin theory in v olv es dif frential equations of the fourth order In con v entional FEM C 1 elements are required to solv e such equations. C 1 elements are cumbersome in 2D and unkno wn in 3D. The high computational cost and lar ge number of de grees of freedom soon place such formulation be yond the realm of practicality Recently some mix ed and hybrid formulations ha v e been de v eloped which require only C 0 continuity b ut none of these elements solv e complicated geometry problems in 2D and there is no problem yet solv ed in 3D. The lar ge number of de grees of freedom is still ine vitable for these formulations. As has been demonstrated earlier RKEM has the potential to solv e higher -order problems, the de gree of freedom consist of nodal displacement and their deri v ati v es. This method has the promise to solv e important problems formulated with higher order deri v ati v es, such as The strain gradient theory iii
Chapter 1 Introduction Classical (local) continuum constituti v e models possess no material/intrinsic length scale. The typical dimensions of length that appear in the corresponding boundary v alue problems are associated with the o v erall geometry of the domain under consideration. In spite of the f act that classical theories are quite suf cient for most applications, there is ample e xperimental e vidence which indicates that, in certain applications, there is signicant dependence on additional length/size parameters. Some of these instances, as selected from the literature, include the dependence of the initial o w stress upon particle size, the dependence of hardness on size of the indenter the ef fect of wire thickness on torsional response, the de v elopment and e v olution of damage in concrete, etc. All these in v estigations highlight the inadequac y of local continuum models in e xplaining the observ ed phenomena, thereby moti v ating the need to introduce non-local continuum models that ha v e length scales present in them. An e xtensi v e summary of this e xperimental e vidence is gi v en in a recent re vie w article by Fleck and Hutchinson . In gradient-type plasticity theories, length scale is introduced through the coef cient of spatial gradients of one or more internal v ariables. In elasticity length scale enters the constituti v e equation through the elastic strain ener gy function, which in this case depends not only on the strain tensor b ut also on gradients of the rotation and strain tensor Con v entional continuum mechanics theories assume that stress at a material point is a function of state' v ariables, such as stress at the same point. This local assumption has long been pro v ed to be adequate when the w a v elength of the deformation eld is much lar ger than the micro-structural length scale of the material. Ho we v er when the tw o length scales are comparable, the assumption is questionable as the material beha vior at a point is inuenced by deformation of neighboring points. Starting from the pioneering Cosserat couple stress theory , v arious non-local or strain gradient continuum theories ha v e been proposed. In the full Cosserat theory , an independent rotation quantity is dened in addition to the material displacement u ; couple stresses (bending moment per unit area) are introduced as the w ork conjugate to the micro-curv ature (that is, the spatial gradient of ). Later T oupin  and Mindlin  proposed a more general theory which includes not only micro-curv ature, b ut also gradients of normal strain. Both the Cosserat and T oupin-Mindlin theories were de v eloped for linear elastic materials. Afterw ards, non-local theories for plastic materials 1
w as de v eloped by among others, Aif antis  and Fleck and Hutchinson [8, 9]. Fleck-Hutchinson strain gradient plasticity theories f all within the T ouplin -Mindlin frame w ork. Interest in non-local continuum plasticity theories ha v e been rising recently due to an increasing number of observ ed size ef fects in plasticity phenomena. V ariants of this strain gradient plasticity theory ha v e also appeared [11, 14, 31]. These theories ha v e been widely applied to studying length scale-dependent deformation phenomena in metals. Polar and higher -order continuum theories ha v e been applied to layered materials, composites and granular media, in addition to polycrystalline metals. The solution of the initial and boundary v alue problems posed in terms of the higher -order theories is not straight forw ard: the go v erning dif ferential equation and boundary conditions are complicated  and analytical solutions are restricted to the simplest cases. Computational dif culties also arise. While boundary conditions are easier to treat in the v ariational setting, requirements of re gularity dictates that the displacement must be a C 1 function o v er the domain. The de grees of freedom include nodal displacement and displacement gradients. The situation is partially analogous to classical Bernoulli-Euler beam and Poisson-Kirchhof f plate theories in one and tw o spatial dimension respecti v ely A nite element formulation incorporating C 1 displacement elds are therefore a natural rst choice for strain gradient theories. F or e xample, the use of Specht' s triangular element  for the special case of couple stress theory w as e xamined . The element contains displacement deri v ati v es as e xtra nodal de grees of freedom and C 1 continuity is satised only in a weak a v eraged sense along each side of the element; therefore the element is not a strict C 1 element. Furthermore, the element f ails to deli v er an accurate pressure distrib ution for an incompressible, non linear solid. There is a rectangular C 1 element , b ut its shape and number of de grees of freedom are strong limitation for its implementation in tw o and higher dimensional problems. As is easily appreciated, the computational cost is high; the lar ge number of de grees of freedom soon place such formulation be yond realm of practicality The lack of rob ust C 1 -continuous elements then moti v ated the de v elopment of v arious C 0 -continuous elements for couple stress theory in recent decades( , ,   ,   etc). Finite element formulations for the Fleck-Hutch strain gradient plasticity theory ha v e been de v eloped with plate elements as a basis, b ut were generally found to perform poorly . Mix ed and hybrid formulations ha v e also been de v eloped in the same w ork and else where the w ork of  introduced some C 0 element types, where nodal de grees of freedom include nodal displacement and corresponding gradients, and the kinematic constraint between displacement and displacement gradient are enforced via the Lagrange multiplier method. Their lo west order triangular element requires 28 unkno wns per element, and their lo west order quadrilateral ele2
ment 38 unkno wns; Amanatidou and Ara v as  proposed mix ed C 0 -continuity nite element formulations, where e v ery element includes around 70 nodal de grees of freedom in 2-D problems. G.Engel et. al.  tries to solv e the same problem using a continuous/discontinuous nite element approximation. In this w ork the y ha v e presented rate of con v er gence and error estimates deri v ed for ener gy and L 2 norms. But this has been limited to 1D problems and to the best of the author' s kno wledge there is no e xtension of this formulation to tw o and higher dimension. From this information it is e vident that currently no nite element method nor e xtensions of nite element methods is a v ailable for strain gradient theory formulation in higher dimensions. During the last tw o decades, the technique of meshless interpolation of trial and test functions has been attracting great attention. Meshfree methods such as the hp-cloud method , the method of nite spheres , the particle partition of unity method , reproducing k ernel particle methods [ ? ] the element free Galerkin method , the nite point method , the dif fuse element method , the modied local Petro v-Galerkin method , smooth particle hydrodynamics , among others, of fer an attracti v e alternati v e for solution of man y classes of problem that are dif cult or e v en not feasible to solv e using the nite elements because these methods posess intrinsic non-local properties. Unlik e a typical nite element method, the nonlocal properties of meshfree approximations confer an arbitrary de gree of smoothness on solutions and ha v e been applied to v arious problems. The Reproducing K ernel Element Method(RKEM) rst presented by Liu et al.is the rst globally compatible, minimal de gree of freedom, arbitrary de gree of continuity basis function on a gi v en mesh. The salient feature of RKEM is the hybridization of nite element shape functions with the mesh free k ernel function. Detailed formulation for RKEM is gi v en in section 3. Use of RKEM to solv e plate and shell structures problem which contain higher order dif ferential deri v ati v e has been demonstrated In this w ork, the detail re vie w of strain gradient theory is gi v en in chapter 2. The application of RKEM to strain gradient theory is discussed in detail in chapter 3 and 4. T o study the accurac y of the present method, con v er gence test are carried out and se v eral problems in one and tw o dimensions are analyzed in chapter 5 and 6. From these tests, the RKEM method is found to gi v e quite accurate results. The remarkable accurac y in these numerical simulations sho w promising characteristics for solving general problems for materials whose constituti v e la w in v olv e strain-gradients. 3
Chapter 2 Re vie w of Linear Elastic Strain Gradient Theory T ouplin  and Mindlin  de v eloped a theory of linear elasticity whereby the strain ener gy density per unit v olume( w ) depends upon both the strain ij ( u i;j + u j;i ) = 2 and strain gradient ij k u k ;ij Here u is the displacement eld and comma represents partial dif ferentiation with respect to a Cartesian co-ordinate. In addition to the Cauchy stress ij this theory also tak es into account higher order stress ij k which is w ork conjugate to ij k Due to the symmetric property of strain tensor ij = j i The strain ener gy function w is assumed to be a con v e x function, with respect to its ar gument( ; ) for each point x of a solid of v olume V The total ener gy W stored in the solid is determined by the displacement eld u (x) within V W ( u ) Z V w ( ( u ) ; ( u ); x ) dx (1) with being deri v ed from u as discussed abo v e. The ener gy v ariation of the solid due to an arbitrary v ariation of the displacement u is : W = Z V ( ij ij + ij k ij k ) dx: (2) Mindlin  sho wed that for a general isotropic linear hyper -elastic material the solid strain ener gy per unit v olume( w ) can be e xpressed as w = 1 2 ii j j + ii j j + a 1 ij j ik k + a 2 iik k j j + a 3 iik j j k + a 4 ij k ij k + a 5 ij k j k i (3) in terms of the in v ariants of the second-order strain tensor and the third-order strain gradient tensor Here and are the standard Lame constants and the v e a n are additional constants of dimensions of stress times length squared. 4
From the constituti v e la w the stress ij and the higher order stress ij k for an elastic solid is deri v ed as : ij = @ w @ ij and ij k = @ w @ ij k F or Linear Elasticity ij and ij k can be written as ij = @ w @ ij = 2 ij + k k ij (4) ij k = @ w @ ij = a 1 ( ipp j k + j pp ik + a 2 ( k pp ij + 1 2 ppi j k + 1 2 ppj ik + a 3 (2 ppk ij ) + a 4 (2 ij k ) + a 5 ( j k i + ik j ) (5) The principal of virtual w ork can be e xpressed as  Z V [ ij ij + ij k ij k ] dV = Z v [ b k u k ] dV + Z s [ f k u k + r k D u k ] dS (6) for an arbitrary displacement increment u Here b k is the body force per unit v olume of the body V while f k and r k are the traction and double stress traction per unit area of the surf ace S, respecti v ely The y are in conjunction with the Cauchy stress ij and higher order stress ij k according to b k + ( ik j ik ;j ) ;i = 0 (7) f k = n i ( ik j ik ;i ) + i j ij k ( D p p ) D j ( i ij k ) (8) and r k = i j ij k (9) In the abo v e equation D j ( ) = ( j k j k ) @ ( ) @ x k is a surf ace gradient operator and D( ) = k @ ( : ) @ x k is surf ace normal-gradient operator n i is the ith componenent of the unit surf ace normal operator The strain gradient theory dictates that, in order to ha v e a unique solution, six boundary conditions must be independently 5
prescribed at an y point on the surf ace of the body i.e. f k or u k and r k or D u k The e xtra boundary conditions on r k or D u k are characteristic of the strain gradient theory and, as can be seen later in one particular e xample, the y imply strain continuity across the interf ace between tw o dissimilar materials. 6
Chapter 3 Re vie w of RKEM A ne w class of method, the Reproducing K ernel Element Method(RKEM)   ,  has been recently de v eloped. This method tak es salient features from the FEM and meshfree methods. The essence of RKEM is to e xtend the notion of the FEM element shape functions to the entire domain and use a meshfree k ernel to combine them in such a w ay as to ensure global compatibility partition of unity and global reproducing properties. 3.1 Concept of RKEM The RKEM method is illustrated briey in this section. As the name implies, RKEM is an element based method not lik e meshless methods. An element is a subset of the domain of the problem. The solution is represented as a linear combination of shape or basis functions on the domain. The shape functions are associated with nodal de grees of freedom (DOF) at each node, which are related to the primary unkno wns and its v arious deri v ati v es. F or each DOF associated global partition polynomials are constructed. Unlik e the construction of high order nite element shape functions, the RKEM shape function does not need an y e xtra de grees of freedom either on the side of the element or in the interior The requirement on the number of de grees of freedom is absolutely minimum. The RKEM is de v eloped to achie v e these follo wing objecti v es: 1. No numerically induced discontinuity between elements. Inter -element boundary continuities limit the smoothness of FEM shape function. F or e xample, to solv e a fourth order dif ferential equation, one needs C 1 elements in a standard conforming method. In practice it is v ery dif cult to get continuity between tw o such element in tw o or higher dimensional domains. RKEM is de v eloped to o v ercome this dif culty 2. No special treatment required for enforcing essential boundary conditions. F or most meshfree methods, the treatment of Dirichlet boundary conditions is problematic due to loss of the Kroneck er delta 7
property of meshfree shape functions. The treatment of Dirichlet boundary conditions in RKEM is straightforw ard. 3. High order smooth interpolation in arbitrary domains of multiple dimensions. RKEM shape functions are formed by combining global partition polynomials and the meshfree reproducing k ernel functions. Global partition polynomials (gpp) for RKEM ha v e the same properties as those of FEM, b ut the y are dened globally Detailed formulation of gpp' s are gi v en in the ne xt section in this chapter It should be noted that the gpp' s are al w ays constructed to be Hermite polynomials. K ernel functions are used to truncate the gpp' s outside a compact support such that the reproducing conditions are satised. The term reproducing condition refers to ability of shape function to reproduce a constant, as part of the requirements for con v er gence. An RKEM interpolation eld is dened as: I f ( x ) = N elem X e =1 Z n e K ( x y : x ) dy nno des X i = 1 e ; i ( x ) f ( x e ; i ) !# (10) where f( x ) is the interpolated function; K is the meshfree k ernel function; N elem is the number of elements in the domain; nnodes is the number of nodes for a particular element, e.g. nnodes =3 for a triangular element; n e is element domain, e;i ( x ) are the global partition polynomials; and I is the interpolation operator After nodal inte gration Eq. 10 can be written as I h f ( x ) = A e 2 V E 24 0@ nnodes X j =1 K ( x x e;j : x ) V e;j 1A nnodes X i =1 e;i ( x ) f ( x e;i ) 35 (11) where V e;j is the nodal inte gration weight, the Lobatto quadrature rule   is used to assign this weight. The symbol A e 2 V E is the assembly operator and denotes the summation o v er the set of elements, ( V E of the mesh. The reproducing k ernel function is a function with compact support and is chosen to ha v e the form K = 1 de;j w x x e ; j e;j b ( x ) (12) where is the support radius, d is the spatial dimension, x is the point at which the k ernel is e v aluatedd, w is a compactly-supported smooth windo w function, and the f actor b ( x ) is used for normalization. A smooth 8
windo w function C n (n) is chosen to serv e as the core of the k ernel. In the present w ork, w in one dimension is tak en to be a cubic spline, which gi v es C 2 continuity If f ( x ) = 1 in Eq. 11 then we ha v e 1 = A e 2 V E 24 0@ nnodes X j =1 1 de;j w x x e;j e;j b ( x ) V e ; j 1A nnodes X i =1 e;i ( x ) 35 (13) Since the global partition polynomials for a single element are a partitions of unity n e i.e: nnodes X i =1 e;i ( x ) = 1 we ha v e A e 2 V E 8<: 0@ nnodes X j =1 1 de;j w x x e;j e;j b ( x ) V e;j 1A 9=; = 1 8 x 2 n The e xpression for the normalizer b( x ) can be written using abo v e equation as : b ( x ) := 8<: A e 2 V E 24 0@ nnodes X j =1 1 de;j w x x e;j e;j V e;j 1A 35 9=; 1 (14) Using the connecti vity relation, the shape function ( I ( x ) ) for RKEM at node I can be written as I ( x ) = l Xk =1 0@ X j 2 V e k 1 de;j w x x e;j e;j V e;j 1A b ( x ) e k ;i k ( x ) (15) and the reproducing k ernel element interpolation in term of RKEM shape function can be written as I h f ( x ) = N nodes X I =1 I ( x ) f I (16) where N nodes is total number of nodes for domain. The RKEM shape functions are Generalized Hermite Functions, these are a set of functions, not necessarily polynomials which satises the interpolation condition P k ( i ) = y i ; k = 0 ; 1 ; ; n i 1 ; 1 = 0 ; 1 ; ; m (17) 9
where P is polynomial and y i is v alue of polynomial at i i.e. Hermite interpolation interpolates not only the primary v ariable, b ut also the rst n i 1 deri v ati v es. W e will denote the RKEM domain by n throughout this thesis. 3.2 Global P artition Polynomial The construction of global partition polynomials for RKEM is similar to that of FEM, e xcept that it is dened globally instead of locally in FEM. The systematic procedure to construct the global partition polynomial is gi v en in . In this section the global partition polynomial for the element used in this w ork, the so called L4P3I1 and T9P2I1 will be discussed. The nomenclature L4P3I1 stands for Linear element which has 4 de grees of freedom, globally reproduces 3rd order Polynomials and Interpolates 1st order deri v ati v es and T9P2I1 stands for T riangular element which has 9 de grees of freedom, globally reproduces 2nd order Polynomials and Interpolates 1st order deri v ati v es. 3.2.1 The L4P3I1 Element F or one dimensional elements, the global partition polynomials are particularly easy Depending on the order of interpolation, one uses the appropriate Hermite interpolant. The L4P3I1 element is constructed using the rst order Hermite polynomials. First order Hermite polynomials on an interv al of length x e with 2 [0,1] are: H 0 1 ( ) = 1 3 2 + 2 3 (18) H 0 2 ( ) = 3 2 2 3 ) (19) H 1 1 ( ) = x e ( 2 2 + 3 ) (20) H 1 2 ( ) = x e ( 3 2 ) (21) (22) The global partition polynomials for this element are the same as the standard FEM beam element shape functions. The shape function 00 and 10 for this element is sho wn in Fig. ?? and Fig. ?? 10
[ u 1 ; u 1 ;x ; u 1 ;y ] [ v 1 ; v 1 ;x ; v 1 ;y ] [ v 3 ; v 3 ;x ; v 3 ;y ] [ u 2 ; u 2 ;x ; u 2 ;y ] [ v 2 ; v 2 ;x ; v 2 ;y ] [ u 3 ; u 3 ;x ; u 3 ;y ] Figure 1. The T9P2I1 element with v ariable associated at each nodes (0,0) (1,0) (1,1) Figure 2. The parent triangle for T9P2I1 element 11
3.2.2 The T9P2I1 Element The triangular element used in this w ork is the T9P2I1. In Fig. 1 the solid circle in triangle represents the node, nodal de grees of freedom associated with each node is gi v en at each v erte x, u and v represent displacement in x and y direction respecti v ely As mentioned pre viously in this section, the smoothness or continuity of the global RKEM shape function is determined by the continuity of the reproducing k ernel, a third order spline is used as the windo w function for this element. Therefore, the element is a global C 2 compatible element. The element has 9 de grees of freedom, the global partition polynomials are constructed using the parametric approach. In the parametric approach the geometric(physical) triangle sho wn in Fig. 1 is mapped to a parent triangle sho wn in Fig. 2. The reason for choosing this parent triangle is e xplained in . The shape function for this parent triangle are: N 1 ( s; t ) = 1 s (23) N 1 ( s; t ) = s t (24) N 1 ( s; t ) = t (25) The global partition polynomials are deri v ed by considering the follo wing element interpolation eld in a three node triangular parent element (e), I f = 3 X i =1 (00) e;i ( x ) f e;i + (10) e;i ( x ) @ f @ x e;i + (01) e;i ( x ) @ f @ y e;i = T f (26) where is denoted as the local shape function array and v ector f is nodal data array i.e. T := 00 e; 1 ; 10 e; 1 ; 01 e; 1 ; 00 e; 2 ; 10 e; 2 ; 01 e; 2 ; 00 e; 3 ; 10 e; 3 ; 01 e; 3 (27) f T := f e; 1 ; f e; 1 ;x ; f e; 1 ;y ; f e; 2 ; f e; 2 ;x ; f e; 2 ;y ; f e; 3 ; f e; 3 ;x ; f e; 3 ;y (28) 12
Based on  the nal result for the global partition polynomials for this parent element are: ~ (00) 1 = 2 s 3 3 s 2 + 1 (29) ~ (10) 1 = s 3 2 s 2 + s (30) ~ (01) 1 = 1 2 ( st 2 + s 2 t t 2 3 st ) + t (31) ~ (00) 2 = 2 t 3 2 s 3 + 3 s 2 3 t 2 (32) ~ (10) 2 = 1 2 ( t 2 st 2 s 2 t + st ) + s 3 s 2 (33) ~ (01) 2 = t 3 1 2 ( st 2 + s 2 t + 3 t 2 3 st ) (34) ~ (00) 3 = 3 w 2 2 w 3 (35) ~ (10) 3 = 1 2 ( sw 2 + s 2 w w 2 w s ) (36) ~ (01) 3 = 1 2 ( tw 2 + t 2 w w 2 w t ) (37) (38) The global RKEM interpolation eld is constructed based on the Eq. 15. Quadratic Polynomials are reproduced using these shape functions. The shape and prole of the global RKEM shape function of this element is displayed in Fig. 3. It is important to note that each plot is scaled for visibility thus relati v e magnitudes cannot be determined from the plots. The shape of the global interpolation function depends upon the mesh and node, whether it is a boundary node or it is in the interior of the domain. The shape functions plotted in Fig. 3 are computed for the middle node of the mesh sho wn in Fig. 9. 3.3 Salient Features of RKEM The follo wing are some salient features of RKEM which dif ferentiates it from FEM and other meshless methods: 1. The shape functions are Generalized Hermite interpolants. As such, at each node the primary v ariable and v arious of its deri v ati v es are interpolated. This property is useful in interpolating deri v ati v es of primary v ariable, e.g. stresses in elasticity problems. 2. The shape functions possess the Higher -order Kroneck er property: D ( ) I x = x J = I J ; x I ; x J 2 n ; j j ; j j m; (39) 13
Figure 3. The global shape function of T9P2I1 element (a) (00)1 (b) (10) 1 (c) (01) 1 14
Inde x means which deri v ati v e of shape function, e.g. in one dimension 00 corresponds to = 0 and 10 correspond to = 1 is an inte ger the maximum v alue of which depends on highest order to which RKEM shape function is dif ferentiable. I and x j are node and coordinates respecti v ely at which RKEM shape function are e v aluated. This property is helpful for treatment of Dirichlet boundary conditions, the majority of meshless methods ha v e dif culty in imposing these boundary conditions. 3. The shape functions form a P artition of Unity: P I I ( x ) = 1 8 x 2 n; j j = 0 This helps to ensure consistenc y for RKEM k ernel function 4. The shape functions possess the global reproducing property: X I ( X I ( x )( D x ) j x = x I ) = x ; 8 x 2 n (40) 5. Each shape function has compact support with x ed size, re gardless of reproducing order The support size is dependent solely on the topology of the RKEM mesh.This property mak e the stif fness matrix sparse,lending to economic formation and solution. 6. Smoothness: Gi v en a k ernel function that is C n continuous, the resulting RKEM shape functions are also C n This property is useful for smooth geometry representation in computational geometry ( ). T o guarantee the abo v e properties, RKEM meshes must satisfy a quasi-uniformity condition as discussed in . In simple terms, this condition places some restriction on the aspect ratio of indi vidual elements and the gradation of element sizes within the mesh. 15
Chapter 4 Galerkin F ormulation for Strain Gradient Problems The plane stress/strain assumption is tak en in 2D e xamples in this research. In general, the displacement eld for tw o dimensional problems are u i = u i ( x 1 ; x 2 ) i = 1 ; 2 u 3 = 0 strain and strain gradients are ij = 1 2 ( u i;j + u j;i ) i = 1 ; 2 : j = 1 ; 2 (41) ij r = u r ;ij r = 1 ; 2 (42) and the corresponding con v entional and higher order stresses are ij = @ w @ ij = 2 ij + k k ij (43) ij k = a 1( ipp j k + j pp ik + a 2( k pp ij + 1 2 ppi j k + 1 2 ppj ik + a 3(2 ppk ij ) + a 4(2 ij k ) + a 5( j k i + ik j ) (44) The equilibrium equation for strain gradient theory is gi v en by Eq. 7, and the corresponding weak form for the problem is Z V [ ij ij + ij k ij k ] dV = Z v [ b k u k ] dV + Z s [ f k u k + r k D u k ] dS (45) F or con v enience, we assume the body force and body double force to be zero, due to this assumption the weak form becomes Z V [ ij ij + ij k ij k ] dV = Z s [ f k u k ] dS (46) 16
The unkno wn v ariable u(displacement in x direction) can be written in terms of RKEM shape functions as : u ( x ) = N nodes X I =1 00I u I + 00I u I + 00I u I (47) It is con v enient to recast this in matrix form as u ( x ) = Nu where N = [ N 1 N 1 ;x N 1 ;y N 2 N 2 ;x N 2 ;y N N nodes N N nodes ;x N N N odes ;y ] and u T = [ u 1 u 1 ;x u 1 ;y u 2 u 2 ;x u 2 ;y u N nodes u N nodes ;x u N N odes ;y ] Similarly u ( x ) = u T N T F or tw o dimensional problems strain and strain-gradient can be written in v ector form as = 266664 xx y y xy 377775 = 266664 u ;x v ;y ( v ;x + u ;y ) 377775 (48) = 266666666664 xxx xy x y y x xxy xy y 377777777775 = 2666666666666664 u ;xx u ;xy u y y v xx v xy v y y 3777777777777775 (49) where v is displacement in y direction. F or the T9P2I1 element, the strain-displacement matrix ( B i ) for node i is gi v en as B i = 266664 00 i;x 10 i;x 01 i;x 0 0 0 0 0 0 00 i;y 10 i;y 01 i;y 00 i;y 10 i;y 01 i;y 00 i;x 10 i;x 01 i;x 377775 (50) 17
such that ( x ) = B u and strain gradient-displacement matrix ( B SG ) i is B SG i = 2666666666666664 00 i;xx 10 i;xx 01 i;xx 0 0 0 00 i;xy 10 i;xy 01 i;xy 0 0 0 00 i;y y 10 i;y y 01 i;y y 0 0 0 0 0 0 00 i;xx 10 i;xx 01 i;xx 0 0 0 00 i;xy 10 i;xy 01 i;xy 0 0 0 00 i;y y 10 i;y y 01 i;y y 3777777777777775 (51) so that ( x ) = B S G u where 00 10 and 01 are the three RKEM shape functions for node i Stress and conjugate stress can be written in v ector form as = 266664 xx y y xy 377775 (52) = 266666666664 xxx xy x y y x xxy xy y 377777777775 (53) F or plane stress, the elastic modulus matrix ( D ) that relates stress ( ) with strain( ) is D = E 1 2 266664 1 0 1 0 0 0 1 2 377775 (54) where E and are Y oung' s Modulus and Poisson' s ratio, respecti v ely 18
F or 2D problem the matrix( D SG that relate with ) is gi v en as D SG = 266666666666666666666666666666666666664 2( a 1 + a 2 + 0 ( a 2 + 2 a 3 ) 0 ( a 2 + 2 a 3 ) 0 a 3 + a 4 + a 5 ) 0 ( a 1 + a 5 0 ( a 5 + 0 : 5 a 2 ) 0 ( a 1 + 0 : 5 a 2 ) +2 a 4 ) 2( a 2 + 2 a 3 ) 0 2( a 3 + a 4 ) 0 ( a 2 + 2 a 5 ) 0 0 ( a 2 + 2 a 5 ) 0 2( a 3 + a 4 ) 0 ( a 2 + 2 a 3 ) ( a 1 + 0 : 5 a 2 ) 0 ( a 5 + 0 : 5 a 2 ) 0 ( a 1 + 2 a 4 + 0 a 5 ) 0 (2 a 1 + 2 a 2 ) 0 ( a 2 + 2 a 3 ) 0 2( a 1 + a 2 + a 3 + a 4 + a 5 ) 377777777777777777777777777777777777775 (55) The matrix equi v alent of the weakform is gi v en as follo ws Ku = P (56) Where u is a v ector of nodal unkno wns, P represents a set of loads applied at the nodes, and K is the stif fness matrix. Both K and P in v olv e inte grals o v er the problem domain. The inte grands in v olv e the basis functions or v arious of their deri v ati v es and products of such functions.The stif fness matrix K can be written as sum of tw o stif fness matrix terms K1 and K2 K1 corresponds to the rst term on the left hand side of Eq. (46) while K2 correspond to the second term on the left hand side of the same equation. The stif fness matrix term K1 can be written as K1 = Z B B T ( x ; y ) D B ( x ; y ) dxdy : (57) The stif fness matrix term K2 can be written in this form K2 = Z B B SG T ( x ; y ) D SG B SG ( x ; y ) dxdy : (58) 19
The load v ector P i for node i because of traction load is gi v en as P = Z n N t dn : (59) where N is 6X2 matrix gi v en as N = 2666666666666664 00 0 10 0 01 0 0 00 0 10 0 01 3777777777777775 (60) and t = 264 t x t y 375 (61) where t x and t y traction in x and y direction.Point load can be added directly into Load V ector as in FEM. 20
Chapter 5 1D Examples 5.1 One-dimensional T oupin-Mindlin Strain Gradient Theory T oupin and Mindlin included higher -order stresses and strains in their theory of linear elasticity which serv es today as the foundation of more adv anced strain gradient plasticity formulations   .Let us introduce a one-dimensional problem follo wing their concepts. The solution for problem choosen in this section has been solv ed by C/DG method  Let n be an open, bounded domain and its boundary Let g q t and r denote the displacement, displacement gradient, couple stress and traction boundaries, respecti v ely The strong form of the problem can be written as ;x ;xx + f = 0 in n ; (62) = g on g (63) ;x n = q on q (64) = r on r (65) ( ;x n = t on t (66) W e ha v e g S t = g T t = ; q S r = q T r = ; and f,g,q,r and t are gi v en data. The constituti v e equations for the stress and higher order stress can be e xpressed as = ;x ; (67) = l 2 ;xx ; (68) 21
where is a material parameter and l is a length scale. W e can write eq(62) to eq(66) with eq(67) and eq(68) as ( ;x ) ;x ( l 2 ;xx ) ;xx + f = 0 in n ; (69) = g on g (70) ;x n = q on q (71) l 2 ;xx = r on r (72) ;x ( l 2 ;xx ) ;x n = t on t (73) From a mathematical point of vie w the one-dimensional T oupin-Mindlin strain gradient theory is a generalization of the Bernoulli-Euler beam theory in v olving in addition to the fourth-order deri v ati v e also a second order deri v ati v e. 5.2 Shear Layer Problem with the T oupin-Mindlin Theory T oupin-Mindlin shear layer problem is considered as a model problem to v alidate our method(RKEM) for a strain gradient theory .Model problem its e xact solution and con v er gence study is presented in this section. W e w ant to simulate a shear -deformable body which is x ed on its left and upon which a traction acts on its right as sho wn in g(4), in this problem the length from the attachment to where the traction acts is assumed to be L. The problem can be formulated as ( ;x ) ;x ( l 2 ;xx ) ;xx + f = 0 in ]0; L [ ; (74) (0) = 0 (75) ;x (0) = 0 (76) ;x ( L ) = 0 (77) ;x ( L ) l 2 ;xx ( L ) = t (78) 22
Figure 4. Shear layer attached on the left side (x =0) with traction acting on the right side (x=L) 23
5.2.1 Model Problem The e xact solution to this problem can be e xpressed as ( x ) = tl ( e L l + l ) 1 e L l + e L x l e x l + t x (79) The comparison of e xact solution with RKEM solution is sho wn in Fig. (5) and plot of error(Exact solution RKEM solution) is sho wn in Fig.(6) 5.2.2 Con v er gence Study The con v er gence rates for interpolation can be seen in Fig.(7), where we ha v e N el = 1 h Rate of con v er gence in L 2 Norm is 3.5. Our numerical observ ations conrm our analytical results and order of interpolation considered lead to good rate of con v er gence in the L 2 Norm. 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0 0.2 0.4 0.6 0.8 1 Fx RKEM (20 elements) Exact Solution Figure 5. Exact solution vs RKEM solution for T oupin-Mindlin shear layer model problem 24
-2.5e-07 -2e-07 -1.5e-07 -1e-07 -5e-08 0 5e-08 1e-07 1.5e-07 2e-07 2.5e-07 3e-07 0 0.2 0.4 0.6 0.8 1 errorx Figure 6. Error plot for T oupin-Mindlin shear layer model problem -8.5 -8 -7.5 -7 -6.5 -6 -5.5 -5 -4.5 -4 -3.5 -3 -2.5 -2 -1.5 log(error)log(h) Convergence plot Figure 7. Con v er gence rate for RKEM L2P3I2 element 25
Chapter 6 2D Examples 6.1 Numerical Examples in T w o Dimension Due to the comple xity and dif culity of a gradient theory the obtainable analytical solutions are restricted to some simple problems.In this section we will focus on tw o problems : 1. Boundary layer analysis; and 2. The stress eld analysis in an innite plate, with a hole, subject to a bi-axial tension p at innity under a plane stress assumption. 6.1.1 Boundary Layer Analysis Higher -order gradient theories predict the e xistence of a boundary layers adjacent to inhomogeneties such as interf aces. Consider for e xample, a bimaterial composed of tw o perfectly bonded half planes of elastic strain gradient solids, subjected to remote shear stress 1 12 as sho wn in Fig(8). 6.1.2 Analytical Solution An analytical solution is presented for a bimaterial, consisting of tw o perfectly bonded half planes of dissimilar linear elastic strain gradient solids. The bimaterial is subjected to a remote uniform shear stress 1 12 as sho wn in Fig(8). Here, we assume that material 1 lying belo w interf ace has a shear modulus 1 and an internal length scale l 1 while Material 2, lying abo v e interf ace has a shear modulus 2 and an internal length scale l 2 F or this bimaterial system, con v entional elasticity theory dictates that the shear stress is uniform and the shear strain jumps in magnitude at the interf ace from 12 = 1 12 2 1 in material 1 to 12 = 1 12 2 2 to material 2. By including strain gradient ef fects, a continuously distrib uted shaer strain can be obtained. In this problem, the only non-zero displacement, strain, stress and higher order stress are u 1 1 2 12 and 221 respecti v ely and the y are functions of the co-ordinate x 2 only From constituti v e equation and , it follo ws 26
Material # 2Material # 1 o 21 21 x 1 x 2 Figure 8. Geometry of a bimaterial under uniform shear 27
that 12 = 2 12 and 221 = 2 l 2 i 221 = 4 l 2 i @ 12 @ x 2 (80) in material i. Substitution of the abo v e relation into the equilibrium equation leads to @ 12 @ x 2 ^ l i 2 @ 3 12 @ x 2 3 = 0 (81) where ^ l i = p 2 l i The general solution to the abo v e ordinary dif frential equation is 12 = d 1 + d 2 e x 2 ^ l 1 + d 3 e x 2 ^ l 1 f or x 2 < 0 (82) and 12 = d 4 + d 5 e x 2 ^ l 2 + d 6 e x 2 ^ l 2 f or x 2 > 0 (83) Here d 1 to d 6 are 6 constants yet to be determined. The general solution is subjected to the follo wing boundary condtions 1. 12 1 12 2 1 and x 2 1 and 12 1 12 2 2 as x 2 1 and at the interf ace 2. continuity of traction: ( 21 221 ; 2 ) j x 2 0 = ( 21 221 ; 2 ) j x 2 0 + ; 3. continuity of double stress traction 221 ; 2 j x 2 0 = 221 ; 2 j x 2 0 + ; 4. continuity of strain 12 j x 2 0 = 12 j x 2 0 + The particular solution satisfying all these conditions is 12 = 1 12 2 1 ( 1 + 1 2 2 2 ^ l 2 1 ^ l 1 + 2 ^ l 2 e x 2 ^ l 1 ) f or x 2 < 0 (84) 12 = 1 12 2 2 ( 1 + 2 1 1 1 ^ l 1 1 ^ l 1 + 2 ^ l 2 e x 2 ^ l 2 ) f or x 2 > 0 (85) In a specic quantati v e e xample, we shall mak e the follo wing arbitrary choice of constituti v e parameters. The shear modulus of material 1 is tak en as to be twice that of material 2. F or each material, the constants a 3 and a 4 as dened in eq(5) are equal to 1 2 l 2 while a 1 a 2 and a 5 v anish. Here l is usually called the internal length scale of materials with strain gradient ef fects. 28
Mesh for a typical problem is sho wn in Figure(9),size of domain is(1 X 1) or ( 50 l X 50 l ), where l is length scale boundary condition is of pure shear ,beside this force boundary condition we set u=v=0 at left-bottom corner and v=0 at right-bottom corner to a v oid the rigid mo v ement. Plot of shear strain 12 for whole domain is gi v en in Fig(10),Fig(11) is shear strain plot by con v entional theory A clear continious band along the interf ace is clear from comparison of these tw o plots.T o mak e point more clear we ha v e sho wn plot of shear strain along x=0.5 by both strain gradient and con v entional theory in Fig(12) Fig(13) is comparison of analytical result and result obtained by RKEM, from this plot it is quite e vident that solution by RKEM matches closely with analytical solutionThe L2 error for this problem were computed and plotted in Fig(14).The slope of line as determined by re gression is 0.857 50 l 50 lt1= 0 t2= 21t1= 0 t2= 21t 1 = 21 t 2 = 0 t 1 = 21 t 2 = 0 Figure 9. Mesh for b undary layer analysis 29
Figure 10. Shear strain plot by strain gradient theory 30
Figure 11. Shear strain plot by con v entional theory 31
0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 -0.4 -0.2 0 0.2 0.4 exyx plot of strain( e xy ) at x=0.5 Strain gradient theory Convemtional theory Figure 12. Comparison of strain gradient vs con v entional theory 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 -0.4 -0.2 0 0.2 0.4 exyx RKEM Exact Solution Figure 13. Comparison of e xact solution vs RKEM solution 32
-2.2 -2.1 -2 -1.9 -1.8 -1.7 -1.6 -1.5 -1.2 -1.1 -1 -0.9 -0.8 -0.7 -0.6 -0.5 Log(error)Log(h) Figure 14. Con v er gence plot for RKEM solution) 6.2 An Innite Plate W ith a Hole In this section we consider the problem of an innite plate with a hole of radius a subjected to a biaxial tension p at innity under plane strain condition as sho wn in Fig(15). The constituti v e model presented in chapter 2 is used to describe the mechancical response of the elastic material. The problem is axially symmetric and the displacement eld is of the form u r = u(r), u = u z = 0, where(r ,z) are c ylindrical coordinates centered at center of the hole. The e xact solution of this problem has been de v eloped by Exadaktylos  and is of the form u ( r ) = p 2 G (1 2 ) r + a 2 r + l c h a r K 1 ( a l ) (1 2 ) K 1 ( r l ) i (86) where c = 1 2 2 K 0 ( a l ) + 1 2 4 l a + a l K 1 ( a l ) (87) and K n ( x ) are the well kno wn modied Bessel functions of the second kind. The problem is solv ed numerically using the RKEM with T9P2I1 element. One quarter of the plate is analyzed because of symmetry; the 33
1 22 = p 1 22 = p1 11= p 1 11= px 1 x 2 a Figure 15. Notation and geometry of an innite plate subjected to bi-axial remote tension 34
u1= 0 t1= pt 2 = p u 2 = 0 Figure 16. mesh for plate problem 35
mesh used for solving the problem is sho wn in Fig(16). The numerical calculations are carried out for = 0.3 and radius of hole(a) =3l, l is internal length scale of material. The size L of the domain analyzed is tak en to be L = 10a. Since L is lar ge compared to both a and l, the solution of this problem is e xpected to be close to that of an innite plate. Fig(17) is comparison of analytical result and result obtained by RKEM for plate problem, from this plot it is quite e vident that solution by RKEM matches closely with analytical solution gi v en in [ ? ]. In Fig(18) and Fig(19) v ariation of r r and has been sho wn respecti v ely these stress v alues are calculated along a radial line = 5 0 From these result we can mak e conclusion for stress concentration 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 urr RKEM Exact Solution Figure 17. V ariation of u r for the plate with a hole v alues, classical elasticity solution with l=0 predicts a stress concentration j r = a = 2p and a v alue of r r j r = a = 0, whereas the present elasticity solution with a =3l predicts j r = a = 1.94p and r r j r = a = 0.16p 36
0 0.2 0.4 0.6 0.8 1 1 2 3 4 5 6 7 8 9 10 srr/pr/a Figure 18. V ariation of r r for the plate with a hole 1 1.2 1.4 1.6 1.8 2 1 2 3 4 5 6 7 8 9 10 sq q/pr/a Figure 19. V ariation of for the plate with a hole 37
Chapter 7 Conclusions The Reproducing k ernel Element Method(RKEM) has been de v eloped for materials with the T oupinMindlin frame w ork of strain gradient type constituti v e theory The good accurac y in these numerical simulations sho ws promising characteristics of RKEM for general problems of material in elasticity where strain gradient ef fect may be important. The concepts and methods de v eloped in this thesis are easily generalizable to other situations and of fer the opportunity to deri v e ne w formulations for the ne w classes of problems and for problems abandoned in the past due to their comple xity 38
References  G. T ouzot B. Nayroles and P V illon. Generalizing the nite element method: Dif fuse approximation and dif fuse elements. Computational Mec hanics 10:307318, 1992.  JT Oden C.Duarte. Hp-cloud a meshless method to solv e boundary-v alue problems. Computer Methods in Applied Mec hanics and Engineering 139:237262, 1996.  E. Cosserat and F Cosserat. Theorie des corps deformables. Hermann et F ils, P aris 1909.  K. J. De, S.; Bathe. The method of nite spheres. Computational Mec hanics 25:329345, 2000.  F Perazzo1 E. Oate and J. Miquel. A nite point method for elasticity problems. Computer s and Structur es 79:21512163, 2001.  E.Amanatidou and N.Ara v as. Mix ed nite element formulation of strain gradient elasticity problems. Comput. Methods Appl. Mec h. Engr g 191:17231751, 2002.  E.C.Aif antis. On the microstructural origin of certain inelastic models. T r ans. ASME J Engng Mater T ec h. 106:326330, 1984.  N. A. Fleck and J. W Hutchinson. A phenomenological theory for strain gradient ef fects in plasticity J Mec h. Phys.Solids 41:18251857, 1993.  N. A. Fleck and J. W Hutchinson. Strain gradient plasticity Adv Appl. Mec h. 33:295361, 1997.  G.E.Exadaktylos. On the problem of the circular hole in an elastic material with microstructure,pri v ate communication.  Y .Huang H.Gao and W .D.Nix. Mechanism based strain gradient plasticity ii. theory J ournal of Mec hanics and Physics of Solids 47:12391263, 1999.  T .J.R. Hughes. The F inite Element Method Do v er Mineola, Ne w Y ork, 2000.  JM Melenk I.Bab uska. The partition of unity method. International J ournal for Numerical Methods in Engineering 40:727758, 1997.  H.Y uan J.Chen and F .H.W ittmann. Computational simulations of micro-indention tests using gradient plasticity CMES: Computer Modeling in Engineering and Science 3(6):743754, 2002.  J.Y .Shu and N.A.Fleck. Prediction of a size ef fect in micro indention. International J ournal of Solids and Structur es 35:13631383, 1998.  S. Li, H. Lu, W Han, W K. Liu, and D. C. Simkins, Jr Reproducing k ernel element method, P art II. Global conforming I m / C n hierarchy Computer Methods in Applied Mec hanics and Engineering 193:953987, 2004.  W K. Liu, W Han, H. Lu, S. Li, and J. Cao. Reproducing k ernel element method: P art I. Theoretical formulation. Computer Methods in Applied Mec hanics and Engineering 193:933951, 2004. 39
 L.T .T enek and E.C.Aif antis. A tw o-dimensional nite element implementation of a special form of gradient elasticity Computer Modelling in Engineering Science .(CMES) 3(6):731742, 2002.  H. Lu, S. Li, D. C. Simkins, Jr ., W K. Liu, and J. Cao. Reproducing k ernel element method P art III. Generalized enrichment and applications. Computer Methods in Applied Mec hanics and Engineering 193:9891011, 2004.  L.B. Luc y A numerical approach to the testing of the ssion hypothesis. The Astr onomical J ournal 82:10131024, 1977.  G.Engel K.Garikipati T .J.R. Hughes M.G. Larson L. Mazzei and R.L. T aylor Continuous/discontinuous nite element approximations of fourth-order elliptic problems in structural and continuum mechanics with applications to thin beams and plates, and strain gradient elasticity Computer Method in Engineering Science 191:36693750, 2002.  R.D. Mindlin. Micro-structure in linear elasticity Ar c h. Rational Mec h. Anal. 16:5178, 1964.  O.C.Zienkie wicz and R.L.T aylor The F inite Element, 4th Edition McGra w-Hill, Ne w Y ork, 1994.  D. C. Simkins, Jr ., S. Li, H. Lu, and W K. Liu. Reproducing k ernel element method P art IV. Globally compatible C n ( n 1) triangular hierarchy Computer Methods in Applied Mec hanics and Engineering 193:10131034, 2004.  Daniel C. Simkins, Jr Gener al Repr oducing K ernel Element Hier ar c hies PhD thesis, Uni v ersity of California, Berk ele y CA, May 2004.  T .Zhu S.N.Atluri. A ne w meshless local petro v-galerkin(mlpg) approach in computational mechanics. Computational Mec hanics 22:117127, 1998b  B. Specht. Modied shape functions for the three node plate bendinf element passing the patch test. International J ournal of Numerical Methods and Engineering 26:705715, 1988.  L.Gu T .Belytschk o, YY Lu. Element-free galerkien methods. International J ournal for Numerical Methods in Engineering 37:229256, 1994.  R. A. T oupin. Elastic materials with couple stresses. Ar c h. Rational Mec h. Anal. 11:385414, 1962.  J.Y .Shu W .E.King and N.A.Fleck. Finite elements for material with strain gradient ef fects. International J ournal for Numerical Methods in Engineering 44:373391, 1999.  W .D.Nix Y .Huang, H.Gao and J.W .Hutchinson. Mechanism based strain gradient plasticity ii. analysis. J ournal of Mec hanics and Physics of Solids 48:99128, 2000.  Z.C.Xia and J.W .Hutchinson. Crack tip eld in strain gradient plasticity J ournal of Mec hanics and Physics of Solids 44:16211648, 1996.  S.Chen Z.T ang and S.N.Atluri. Analysis of materials with strain-gradient ef fects:mlpg approach with nodal displacement only Computer Modelling in Engineering Science .(CMES) 4(1):177196, 2003. 40