Interface problems in nonlocal diffusion and sharp transitions between local and nonlocal domains

Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

Contents lists available at ScienceDirect

Comput. Methods Appl. Mech. Engrg. journal homepage: www.elsevier.com/locate/cma

Interface problems in nonlocal diffusion and sharp transitions between local and nonlocal domains Pablo Seleson a,⇑, Max Gunzburger b, Michael L. Parks c a b c

Institute for Computational Engineering and Sciences, 201 East 24th St, Stop C0200, Austin, TX 78712-1229, United States Department of Scientific Computing, Florida State University, 400 Dirac Science Library, Tallahassee, FL 32306-4120, United States Sandia National Laboratories, Computing Research Center, P.O. Box 5800, MS 1320, Albuquerque, NM 87185-1320, United States

a r t i c l e

i n f o

Article history: Received 3 April 2012 Received in revised form 14 March 2013 Accepted 24 May 2013 Available online 7 June 2013 Keywords: Nonlocal diffusion Interface problems Interface conditions Multiscale modeling

a b s t r a c t We investigate interface problems in nonlocal diffusion and demonstrate how to reformulate and generalize the classical treatment of interface problems in the presence of nonlocal interactions. Through formal derivations, we show that nonlocal diffusion interface problems converge to their classical local counterparts, in the limit of vanishing nonlocality. A central focus of this paper is local/nonlocal interface problems, or interface problems with sharp transitions between local and nonlocal domains. Such problems can be cast as instances of a nonlocal interface problem, with a finite horizon in certain regions and a vanishing horizon in other regions. We derive a local/nonlocal interface problem and utilize conservation principles to obtain local/nonlocal interface conditions. Comparisons between nonlocal, local, and local/nonlocal interface problems are presented, analytically and numerically, with a focus on multiscale aspects of nonlocal models induced by their inherent length scales. Ó 2013 Elsevier B.V. All rights reserved.

1. Introduction Nonlocal models have been utilized in many fields, including image processing [1,2], nonlocal Dirichlet forms [3,4], kinetic equations [5,6], phase transition [7–9], nonlocal solid mechanics [10–17], nonlocal heat conduction [18], and nonlocal diffusion [19–22]. In this paper, we consider the development, analysis, and implementation of interface problems in nonlocal diffusion. We refer to the boundary separating domains with different model or material parameters as the material interface. We focus specifically on material interfaces defined by sharp discontinuities in model or material parameters. Local diffusion models, which utilize partial differential equations (PDEs), are undefined on such material interfaces. As such, local diffusion models require separate interface conditions to be defined on material interfaces. Special techniques are used to derive these interface conditions [23]. For a local model, the material interface is also the interaction interface. Nonlocal diffusion models, based upon integral equations (IEs), involve direct interactions between points separated by finite distances. In nonlocal models, the interaction interface is a volumetric region surrounding the material interface, meaning that these interfaces are not coincident. As such, nonlocal interface problems may be regarded as a generalization of local interface problems. We note that special interface conditions are not needed for nonlocal models. Their governing equations are valid everywhere in the domain, even at material discontinuities, since they do not assume spatial differentiability of material or model parameters. We show that, under certain assumptions, nonlocal diffusion models converge to classical local diffusion models in the limit of vanishing nonlocality or local limit. Nonlocal models introduce length scales and thus are suitable for multiscale modeling. Of particular interest are problems that feature significantly different length scales in different subdomains. We formulate interface problems for systems with sharp transitions between local and nonlocal regions. A system of this type can be cast as an instance of a nonlocal interface ⇑ Corresponding author. E-mail addresses: [email protected] (P. Seleson), [email protected] (M. Gunzburger), [email protected] (M.L. Parks). 0045-7825/$ - see front matter Ó 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.cma.2013.05.018

186

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

problem, with a finite horizon on one side of the material interface and a vanishing horizon on the other side. Subdomain governing equations are derived for such systems and conservation principles are used to obtain local/nonlocal interface conditions along material interfaces. The analysis and the numerical simulations of local/nonlocal systems for static and diffusion problems provide valuable insights toward concurrent multiscale modeling. Methods of this type in solid mechanics appear in the atomistic-to-continuum (AtC) coupling literature, where an interface connects discrete (possibly nonlocal) and continuum (local) subdomains [24]. Multiscale approaches in solid mechanics using the nonlocal continuum peridynamics theory appear, e.g., in [25–32]. This paper is organized as follows. In Section 2, we introduce nonlocal diffusion models and the generalization of domain boundaries. In Section 3, we describe nonlocal interface problems, emphasizing differences between material and interaction interfaces. We also introduce a nonlocal flux, showing that governing equations in nonlocal diffusion are instances of conservation principles. The convergence of nonlocal interface problems to local interface problems is presented in Section 4, resulting in an explicit relation between local and nonlocal diffusion tensors. We also show in Appendix A the convergence of the nonlocal flux to a local flux. In Section 5, we present local/nonlocal interface problems for diffusion, deriving governing equations and an interface condition at the material interface. In Section 6, we present numerical simulations for one-dimensional static and diffusion problems, comparing nonlocal diffusion models with local and local/nonlocal diffusion models. Conclusions are given in Section 7. 2. A nonlocal diffusion model This section describes a nonlocal diffusion model. Models of this type can be found in [19,20,22]. Related models for nonlocal heat conduction appear, e.g, in [18]. Those models are special instances of the model presented in this section. In this paper, we use lower-case letters for scalars, e.g., u; q; c, lower-case boldface letters for vectors, e.g., x, and upper-case blackboard letters for second-order tensors, e.g., K; D [33]. Double indices are used to denote components of tensors, e.g., ðKÞ‘k . The nonlocal diffusion problem is to find the diffusing material density u such that

8 @u nl ðaÞ > < @t ðx; tÞ ¼ LX ðuÞðx; tÞ þ qðx; tÞ; in X  ð0; TÞ; uðx; tÞ ¼ 0; in B X  ð0; TÞ; ðbÞ > : in X; ðcÞ uðx; 0Þ ¼ u0 ðxÞ;

ð1Þ

where

Lnl X ðuÞðx; tÞ :¼

1 2

Z

fT ðuÞ½x; thx0  xi  T ðuÞ½x0 ; thx  x0 igdV x0

ð2Þ

X

is an integral (nonlocal) operator representing material transport due to some long-range interaction mechanism defined by T . The operator T ðuÞ½x; thi : Rd ! R (d is the dimension) is a general mapping, defined at the point x at time t, depending on the state variable u.1 Data for the initial condition u0 and source q are provided. Nonlocal interactions involve action-at-a-distance. Given a domain X  Rd , nonlocal boundary conditions are imposed on a boundary layer denoted BX  Rd defined as [34]

BX :¼ supp ðT Þ n X:

ð3Þ

This is in contrast to classical local problems involving PDEs, for which boundary conditions are imposed on a boundary surface @ X  Rd1 . An illustration of a boundary layer BX for a nonlocal problem with a single interaction range appears in Fig. 1(a). In this case, the boundary layer is of uniform width. A nonlocal problem with different interaction ranges has a boundary layer of variable width, as illustrated in Fig. 1(b). In (1b), homogeneous Dirichlet-like nonlocal boundary conditions are imposed. The nonlocal closure of X is defined by X :¼ X [ BX [34]. For simplicity, we choose the linear constitutive relation

T ðuÞ½x; thx0  xi ¼ cðx0 ; xÞðuðx0 ; tÞ  uðx; tÞÞ; d

d

ð4Þ

þ

0

with cð; Þ : R  R ! R [ f0g a two-point scalar-valued function [35]. The mapping cðx ; xÞ is defined at the point x for its interaction with the point x0 ; this mapping is not necessarily symmetric. Nonlocal models usually assume a finite interaction range by introducing a length scale 0 < e < 1, referred to as the horizon of the nonlocal model. In nonlocal diffusion, we assume direct mass transport only occurs between a given point x 2 X and its neighborhood

Hðx; eÞ :¼ fx0 2 Rd : kx0  xk < eg:

ð5Þ

We thus have finite support in (4), i.e., cðx ; xÞ ¼ 0 for kx  xk P e. By using (4) in (2) we obtain 0

Lnl ðuÞðx; tÞ ¼ X

Z

0

csym ðx0 ; xÞðuðx0 ; tÞ  uðx; tÞÞdV x0 ;

X

where the symmetric part of cðx0 ; xÞ is

1

Operators analogous to T in peridynamic solid mechanics are called force states, and act as generalizations of stress tensors [16].

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

187

Fig. 1. Two-material system illustration. Material 1 is dark gray and material 2 is light gray. The boundary layer BX is uncolored. In (a), we illustrate the different interaction types in a nonlocal system. These may be intra-material interactions (involving pairs of points of the same material) or inter-material interactions (involving pairs of points of dissimilar materials). We assume a single-horizon system, so that the boundary layer is of uniform width. In (b), we show the subdomains X1 ; X2 , the interaction interface C  Rd , the material interface C0  Rd1 , and the boundary layer BX with non-uniform width due to the different horizons in the system. We assume that material 1 and material 2 are characterized by constant horizons e1 and e2 , respectively, with e1 > e2 . We define  :¼ maxðe1 ; e2 Þ.

csym ðx0 ; xÞ :¼

cðx0 ; xÞ þ cðx; x0 Þ : 2

ð6Þ

We now introduce a model capable of representing heterogeneity and anisotropy. Anisotropic behavior is especially important in orthotropic materials. Our model can be extended to nonlocal solid mechanics to describe fiber-reinforced composite laminates such as those presented in [36–38]. Let the two-point scalar-valued function cðx0 ; xÞ be (cf. [33]) a cðx0 ; xÞ ¼ veðxÞ ðx0  xÞ cd; eðxÞ

ðx0  xÞ  Kðx; x0  xÞðx0  xÞ kx0  xk2þa

;

ð7Þ

a with cd; eðxÞ depending on the dimension d, the parameter a, and possibly on x through the horizon eðxÞ, but independent of 0 n :¼ x  x. The characteristic function, which ensures a finite support for cðx0 ; xÞ, is defined as

ve ðnÞ :¼



1 knk < e; 0

otherwise:

ð8Þ

Furthermore, the nonlocal diffusion tensor K is assumed to satisfy

lim e!0

Z

Hð0;eÞ

cd;e a 2knk2þa

ðn  nÞKðx; 0Þðn  nÞdV n < 1:

ð9Þ

We assume a < d þ 2, which ensures the existence of the integral in (9) (cf. Section 4.2) [39]. Remark 2.1. Given two points x; x0 2 X, we may have eðxÞ – eðx0 Þ and Kðx; x0  xÞ – Kðx0 ; x  x0 Þ, so that cðx0 ; xÞ – cðx; x0 Þ. Regardless, the integrand in (2) is antisymmetric.

3. Nonlocal interface problems We define interface problems to be problems involving non-overlapping regions, with a nonzero intersection between their closures, described by different model or material parameters. These regions are said to be connected by interfaces. For such problems concerning nonlocal diffusion, different regions may be characterized by different horizons and diffusion tensors. In the case of two-material nonlocal systems, two types of interactions are present: intra-material interactions and inter-material interactions. The former consist of interactions involving pairs of points belonging to the same material, either material 1 or material 2. The latter consist of interactions involving pairs of points of dissimilar materials (see Fig. 1(a)). Due to the presence of long-range interactions, we can identify two types of interfaces in nonlocal systems. The first one, referred to as the material interface, is a boundary located at the discontinuities of the model or material parameters (e.g., C0  Rd1 in Fig. 1(b)). The second one, referred to as the interaction interface, is a volumetric layer containing all inter-material interactions (e.g., C  Rd in Fig. 1(b)). Assume a two-material system, where material 1 and material 2 are characterized by constant horizons e1 and e2 , respectively, and define

 :¼ maxðe1 ; e2 Þ: We define the subdomains Xi ; i ¼ 1; 2, and the interaction interface C (see Fig. 1(b)) as

ð10Þ

188

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

n

o

Xi :¼ x 2 X : e and K continuous in Hðx; Þ \ X; Kðx; Þ ¼ Ki ðx; Þ; eðxÞ ¼ ei ; n

i ¼ 1; 2;

o

C :¼ x 2 X : K or e discontinuous in Hðx; Þ \ X :

ð11Þ ð12Þ

The subdomains X1 and X2 are decoupled, i.e., no direct interactions exist between these two domains. Their interaction is carried through the interaction interface C . The definitions (11) and (12) can be generalized to the case of multi-material systems with variable horizons. Interface problems and their analysis are based on the concept of a flux and a conservation principle. In diffusion, for instance, we use the conservation of mass, which states that the rate of change of mass within any control volume is equal to the production of mass, by possible sources, in that volume minus the transport of mass out of it. The transport of mass is determined by the flux. The choice of conservation principle depends on the given physical problem. For example, in heat conduction, one uses the principle of conservation of energy. A summary of conservation principles and their corresponding fluxes for different physical problems is provided in [23, p. 44]. The nonlocal flux is identified as

RXA ;XB :¼

Z

Z

XA

XB

1 fT ðuÞ½x; thx0  xi  T ðuÞ½x0 ; thx  x0 igdV x0 dV x ; 2

ð13Þ

which represents the net nonlocal mass transport from XB into XA . As opposed to a local flux, a nonzero nonlocal flux between XA and XB may exist even when their closures have an empty intersection [21,40]. In Appendix A, we show that the nonlocal flux (13) converges to a local flux in the limit of vanishing nonlocality. The nonlocal flux possesses the following properties (cf. [21,35,40]):

RXA ;XB ¼ RXB ;XA ;

ð14Þ

RXA ;XA ¼ 0;

ð15Þ

RXA [XC ;XB ¼ RXA ;XB þ RXC ;XB ;

ð16Þ

RXA ;XB [XC ¼ RXA ;XB þ RXA ;XC :

ð17Þ

Properties (15) and (17) give the important result

RXA ;XB nXA ¼ RXA ;XB nXA þ RXA ;XA ¼ RXA ;XB :

ð18Þ

Let x  X denote a control volume. By the principle of conservation of mass,

d dt

Z

uðx; tÞdV x ¼

x

Z

qðx; tÞdV x þ

x

Z Z x

Xnx

1 fT ðuÞ½x; thx0  xi  T ðuÞ½x0 ; thx  x0 igdV x0 dV x : 2

Using property (18), we can write

 Z  Z @u 1 ðx; tÞ  qðx; tÞ  fT ðuÞ½x; thx0  xi  T ðuÞ½x0 ; thx  x0 igdV x0 dV x ¼ 0: x @t X 2 Because x is arbitrary, we obtain the nonlocal diffusion equation (1a) at each point in X. This demonstrates that (1a) is an instance of a conservation principle. Because (1a) does not assume that model or material parameters are spatially differentiable it holds everywhere in X, i.e., on or off material discontinuities. This is in contrast to classical local models where governing PDEs are valid for points in smooth regions but do not hold at material interfaces. Interface conditions are thus required in local models, in addition to subdomain equations. The nonlocal diffusion interface problem is given by (1) and does not require any additional interface conditions. 4. Convergence of the nonlocal interface problem to the classical local interface problem In this section, using formal derivations, we demonstrate that the nonlocal interface problem (1) converges, under suitable assumptions, to the classical local interface problem in the limit of vanishing nonlocality or local limit, i.e.,  ! 0 (cf. (10)). 4.1. Convergence of the nonlocal diffusion equation to the local diffusion equation Assume the linear constitutive relation (4) with (7). Then, the governing equation (1a) in Xi ; i ¼ 1; 2; t > 0, can be written as

@u ðx; tÞ ¼ @t

Z Hð0;ei Þ

cd;ei a

n knk2þa

  Ki ðx; nÞ þ Ki ðx þ n; nÞ n ðuðx þ n; tÞ  uðx; tÞÞdV n þ qðx; tÞ: 2

ð19Þ

For the purpose of connecting local and nonlocal diffusion models, we assume u and Ki are smooth enough so that the following Taylor expansions

189

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

  1 uðx þ n; tÞ ¼ uðx; tÞ þ ðn  rÞuðx; tÞ þ ðn  rÞðn  rÞuðx; tÞ þ O knk3 ; 2   Ki ðx; nÞ ¼ Ki ðx; 0Þ þ ðn  rn ÞKi ðx; 0Þ þ O knk2 ;

ð21Þ

  Ki ðx þ n; nÞ ¼ Ki ðx; 0Þ þ ðn  rÞKi ðx; 0Þ  ðn  rn ÞKi ðx; 0Þ þ O knk2 ;

ð22Þ

ð20Þ

exist for x 2 Xi and knk < ei , where r ¼ rx . By using (21) and (22), we obtain

  Ki ðx; nÞ þ Ki ðx þ n; nÞ 1 ¼ Ki ðx; 0Þ þ ðn  rÞKi ðx; 0Þ þ O knk2 : 2 2 Using, in addition, (20) we can express the nonlocal diffusion equation (19) as  Z      cd;ei a @u 1 1 2 3 ðx;tÞ ¼ n  Ki ðx; 0Þ þ ðn  rÞKi ðx;0Þ þ O knk n ðn  rÞuðx;tÞ þ ðn  rÞðn  rÞuðx; tÞ þ O knk dV n þ qðx;tÞ: 2þ a @t 2 2 Hð0;ei Þ knk Collecting all the contributions to the leading terms and noting that the integration range is spherically symmetric, we get

)

@2u @ðKi Þ‘k @u nn n‘ ðKi Þ‘k ðx; 0Þnk nm ðx; tÞ þ nn n‘ ðx; 0Þnk nm ðx; tÞ dV n þ O e2i þ qðx; tÞ 2þa @x @x @x @x n m n m Hð0;ei Þ 2knk

¼ Di ðxÞ : ðr  rÞuðx; tÞ þ ðr  Di ðxÞÞ  ruðx; tÞ þ O e2i þ qðx; tÞ

¼ r  ½Di ðxÞruðx; tÞ þ O e2i þ qðx; tÞ; Z

@u ðx; tÞ ¼ @t

(

cd;ei a

where

Di ðxÞ :¼

Z

cd;ei a

Hð0;ei Þ

ðn  nÞKi ðx; 0Þðn  nÞdV n

2knk2þa

ð23Þ

a 2þda and cd; is a scaling factor chosen so that Di ðxÞ does not depend on ei / 1=ei cover the classical local subdomain diffusion equation

ei (cf. Section 4.2). In the limit ei ! 0, we re-

@u ðx; tÞ ¼ r  ½Di ðxÞruðx; tÞ þ qðx; tÞ; @t

ð24Þ

where Di is the diffusion coefficient tensor for material i. Remark 4.1. The connections established here between local and nonlocal models use their strong forms. Similar results can be obtained for their weak forms with weaker assumptions about the smoothness of u and Ki (cf. [33,34]). 4.2. Explicit relation between Di and Ki We now write (23) in component form:

ðDi Þnm ðxÞ ¼

1 d;a c ðKi Þ‘k ðx; 0Þ 2 ei

Z

nn n‘ nk nm

Hð0;ei Þ

knk2þa

dV n :

ð25Þ

In the following, we calculate the integrals in (25) for different dimensions. In 3D (a < 5) and 2D (a < 4),

Z Hð0;ei Þ

nn n‘ nk nm 2þa

knk

dV n ¼ ½dn‘ dkm ð1  dnm Þ þ dnk d‘m ð1  dk‘ Þ þ dnm d‘k ð1  dn‘ Þ

Z Hð0;ei Þ

n21 n22 2þa

knk

dV n þ dnm d‘k dn‘

Z Hð0;ei Þ

n41 knk2þa

dV n

a d;a d;a ¼ ðdn‘ dkm ð1  dnm Þ þ dnk d‘m ð1  dk‘ Þ þ dnm d‘k ð1  dn‘ ÞÞId; ei þ 3dnm d‘k dn‘ Iei ¼ ðdn‘ dkm þ dnk d‘m þ dnm d‘k ÞIei a 2;a 5a 4a with dnm the Kronecker’s delta, d ¼ 2; 3 the dimension, I3; ei :¼ 4pei =15ð5  aÞ and Iei :¼ pei =4ð4  aÞ. In 1D (a < 3),

Z

Hð0;ei Þ

nn n‘ nk nm knk

2þa

dV n ¼

Z

ei

ei

n4 jnj

2þa

dn ¼

Z

ei ei

jnj2a dn ¼ 2

Z

ei 0

a n2a dn ¼ 3I1; ei ;

a 3a with I1;  aÞ. ei :¼ 2ei =3ð3 d;a 1 a Take cd; e ¼ ðI e Þ , which can be written as

cd;e a ¼

1 ðd þ 2Þðd þ 2  aÞ ; e2a jHð0; eÞj

ð26Þ

with jHð0; eÞj the size of the neighborhood, i.e., its volume, area, or length, for d ¼ 3; 2; or 1, respectively. Then, (9) is satisfied and we obtain

190

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

Di ðxÞ ¼

Ki ðx; 0Þ þ KTi ðx; 0Þ 1 þ trðKi Þðx; 0Þ I; 2 2

ð27Þ

with I a unit tensor. In 1D, we can obtain the alternative simpler relation Di ðxÞ ¼ K i ðx; 0Þ, by choosing instead 1

c1;e a ¼ 2ðI1;e a Þ =3. Relation (27) can be inverted by observing that trðKi Þðx; 0Þ ¼

2 trðDi ÞðxÞ 2þd

with d the dimension. If Ki ¼ KTi , we obtain

Ki ðx; 0Þ ¼ Di ðxÞ 

1 trðDi ÞðxÞ I; 2þd

ð28Þ

which is similar, but not identical to the deviatoric part of Di . 4.3. The local limit of the nonlocal interface problem In Section 4.1, we demonstrated that the nonlocal diffusion equation (1a) converges to the classical local diffusion equation (24) under the assumptions (20)–(22). These assumptions do not hold, in general, for n ¼ x0  x with x and x0 points of dissimilar materials, due to discontinuities in the model horizon or in the diffusion tensor along C0 . Eq. (24) is thus only valid in the subdomains X1 and X2 , but not in the interaction interface C . However, in the limit as  ! 0, the interaction interface C reduces to the material interface C0 , so that (24) becomes valid everywhere in X n C0 . To ensure conservation of mass in the system we thus need an interface condition along C0 . A classical approach can be used to derive the interface condition (29d) below (cf. [23]). The nonlocal interface problem, in the local limit, is to find the diffusing material density u such that

8 @u > ðx; tÞ ¼ Ll ðuÞðx; tÞ þ qðx; tÞ; > @t > > < uðx; tÞ ¼ 0; > uðx; 0Þ ¼ u0 ðxÞ; > > > : D ðxÞruðÞ ðx; tÞ  n ¼  D ðxÞruðþÞ ðx; tÞ  n ; 1 1 2 2

in X1 [ X2  ð0; TÞ; ðaÞ on @ X  ð0; TÞ;

ðbÞ

in X;

ðcÞ

on C0  ð0; TÞ;

ðdÞ

ð29Þ

with n1 and n2 outward unit normals to X1 and X2 , respectively,

Ll ðuÞðx; tÞ :¼ r  ½DðxÞruðx; tÞ; and

 DðxÞ :¼

D1 ðxÞ x 2 X1 ; D2 ðxÞ x 2 X2 :

ð30Þ

ð31Þ

In (29d), we use the standard notation of ruðÞ ðx; tÞ and ruðþÞ ðx; tÞ to represent the evaluation of the gradient of u from the left and right sides of the material interface, respectively. 5. Local/nonlocal interface problems We study two-material systems, where one material is characterized by nonlocal interactions and the other by local interactions. Interactions between these two materials occur through a local/nonlocal interaction interface. We observe that such a system is an instance of a nonlocal interface problem where the interactions are nonlocal with a finite horizon  > 0 in one subdomain and local in the other subdomain, which we represent as ‘‘nonlocal’’ interactions with a vanishing horizon. We approach the formal analysis of local/nonlocal interface problems from this perspective. In this section, we formulate interface problems for these types of systems and derive a local/nonlocal interface condition. The geometry of local/nonlocal interface problems deserves special attention. Following Section 3, we expect the subdomain X1 (containing points in the nonlocal region) to be described by a nonlocal diffusion model and the subdomain X2 (containing points in the local region) to be described by a classical local diffusion model. In contrast, governing equations for points in C which contains inter-material interactions, will combine local and nonlocal contributions. The width of the interaction interface is determined by the horizon of the nonlocal model, and is equal to 2. In addition, part of the boundary of the domain is a volumetric layer of width  and part of it is a classical boundary surface. An illustration is given in Fig. 2. 5.1. Subdomain equations in local/nonlocal interface problems for diffusion We derive governing equations for subdomains in a local/nonlocal interface problem (cf. Fig. 2). We begin with the nonlocal framework presented in Section 3, assuming a nonlocal interface problem with two finite horizons: e1 and e2 (cf. Fig. 1(b)), i.e.,

(

eðxÞ ¼

e1 x 2 Xf1g ; e2 x 2 Xf2g :

ð32Þ

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

191

Fig. 2. Two-material interface problem with a sharp transition across a material interface C0 between local and nonlocal subdomains. Material 1 (darkgray) is described by nonlocal interactions with horizon  > 0. Material 2 (light-gray) is described by local interactions. The boundary layer is uncolored. The boundary of the entire system is the union of the boundary layer BXf1g n C (cf. (34)) and the boundary surface @ X2 n ‘2 . The lines ‘i are the boundary surfaces between C and Xi ; i ¼ 1; 2, respectively. Further details for the subdomains are given in Fig. 1(b).

The material subdomains are defined as

n

o

Xfig :¼ x 2 X : Kðx; Þ ¼ Ki ðx; Þ; eðxÞ ¼ ei ; fig

X

:¼ fx 2 X : Kðx; Þ ¼ Ki ðx; Þ; eðxÞ ¼ ei g;

i ¼ 1; 2;

ð33Þ

i ¼ 1; 2:

ð34Þ

Note that Xfig (Xfig . We also define the interaction interface material subdomains

Cfig  :¼ fx 2 C n C0 : Kðx; Þ ¼ Ki ðx; Þ; eðxÞ ¼ ei g;

i ¼ 1; 2:

To derive subdomain equations in a local/nonlocal interface problem for diffusion, we proceed as follows. First, we write down the governing equation for the subdomain of interest in a reference nonlocal diffusion model, assuming 0 < e2 < e1 ¼  (cf. (10)) (without loss of generality, we assume X2 represents, in the limit, a local domain). Second, we take the limit of e2 ! 0 while keeping e1 finite. In a reference nonlocal model, the governing equation at x 2 X; t > 0, can be written as

@u ðx; tÞ ¼ @t

Z

1

X

2knk

n 2þa





a a 0 0 veðxÞ ðnÞcd;eðxÞ Kðx; nÞ þ veðx0 Þ ðnÞcd; eðx0 Þ Kðx ; nÞ nðuðx ; tÞ  uðx; tÞÞdV x0 þ qðx; tÞ;

ð35Þ

where n :¼ x0  x. Eq. (35) can be written as

@u ðx; tÞ ¼ @t

where

Z

  d;a d;a 0 n  v ðnÞ c Kðx; nÞ þ v ðnÞ c K ðx ; nÞ nðuðx0 ; tÞ  uðx; tÞÞdV x0 1 e e ðxÞ e e ðxÞ 2þa 1 1 Xf1g 2knk Z   1 a d;a 0 0 þ n  veðxÞ ðnÞcd; eðxÞ Kðx; nÞ þ ve2 ðnÞce2 K2 ðx ; nÞ nðuðx ; tÞ  uðx; tÞÞdV x0 þ qðx; tÞ; 2þa f2g 2knk X 1

( Kðx; nÞ ¼

K1 ðx; nÞ x 2 Xf1g ; K2 ðx; nÞ x 2 Xf2g :

ð36Þ

ð37Þ

We use (36) in different subdomains, to obtain governing equations in the e2 ! 0 limit. Subdomain X1 . Let x 2 X1 . Then kx0  xk > e1 ; e2 for all x0 2 Xf2g . The second term on the right-hand side of (36) thus vanishes. In addition, X1  Xf1g . We can then express the governing equation at x as

@u ðx; tÞ ¼ @t

Z

ve1 ðnÞcd;e1a

X1

2knk2þa

n  ðK1 ðx; nÞ þ K1 ðx0 ; nÞÞnðuðx0 ; tÞ  uðx; tÞÞdV x0 þ qðx; tÞ

¼ Lnl ðuÞðx; tÞ þ qðx; tÞ; X1

ð38Þ

with Lnl X ðuÞðx; tÞ given in (2). This equation does not depend on e2 . It thus remains unaltered in the limit as e2 ! 0. Subdomain X2 . Let x 2 X2 . Then kx0  xk > e1 ; e2 for all x0 2 Xf1g . The first term on the right-hand side of (36) vanishes and, in addition, Hðx; e2 Þ  Xf2g . The governing equation at x can then be expressed as

@u ðx; tÞ ¼ @t

Z Hðx;e2 Þ

cd;e2a 2knk2þa

n  ðK2 ðx; nÞ þ K2 ðx0 ; nÞÞnðuðx0 ; tÞ  uðx; tÞÞdV x0 þ qðx; tÞ:

We assume u and K2 are smooth enough to satisfy (20)–(22), for x with knk < e2 . Then, following Section 4.1, we obtain in the limit as e2 ! 0 the governing equation

192

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

@u ðx; tÞ ¼ Ll ðuÞðx; tÞ þ qðx; tÞ @t

ð39Þ

with Ll ðuÞðx; tÞ given in (30), and DðxÞ ¼ D2 ðxÞ given in (23). f1g f2g Subdomain Cf1g ¼ ;. Then,  . Let x 2 C with Hðx; e2 Þ \ X Thus, the governing equation at x is

Z

@u ðx; tÞ ¼ @t

ve1 ðnÞcd;e1a

n  ðK1 ðx; nÞ þ K1 ðx0 ; nÞÞnðuðx0 ; tÞ  uðx; tÞÞdV x0

2knk2þa

Xf1g

Z

ve2 ðnÞ ¼ 0 in the second term on the right-hand side of (36).

ve1 ðnÞcd;e1a

n  K1 ðx; nÞnðuðx0 ; tÞ  uðx; tÞÞdV x0 þ qðx; tÞ 2knk2þa Z 1 ¼ Lnlf1g ðuÞðx; tÞ þ T ðuÞ½x; thx0  xidV x0 þ qðx; tÞ; X 2 Xf2g þ

Xf2g

ð40Þ

f1g with Lnl X ðuÞðx; tÞ given in (2). In the limit as e2 ! 0, this equation holds for all x 2 C . f2g f2g f2g Subdomain C . Let x 2 C with Hðx; e2 Þ  X . Then, veðxÞ ðnÞ ¼ 0 in the first term on the right-hand side of (36). Thus, the governing equation at x is

Z

@u ðx; tÞ ¼ @t

ve1 ðnÞcd;e1a 2knk2þa

Xf1g

þ

Z

n  K1 ðx0 ; nÞnðuðx0 ; tÞ  uðx; tÞÞdV x0

cd;e2a Hðx;e2 Þ

2knk2þa

n  ðK2 ðx; nÞ þ K2 ðx0 ; nÞÞnðuðx0 ; tÞ  uðx; tÞÞdV x0 þ qðx; tÞ:

We assume u and K2 are smooth enough to satisfy (20)–(22), for x with knk < e2 . Then, following Section 4.1, we obtain in the limit as e2 ! 0 the governing equation

@u 1 ðx; tÞ ¼  @t 2

Z

T ðuÞ½x0 ; thx  x0 idV x0 þ Ll ðuÞðx; tÞ þ qðx; tÞ:

ð41Þ

Xf1g

This equation holds for all x 2 Cf2g  in the limit as

e2 ! 0.

5.2. Interface condition for local/nonlocal interface problems for diffusion In Section 5.1, we obtained governing equations for points in X n C0 . We now derive an interface condition on C0 , which ensures conservation of mass in X. For this purpose, we integrate Eqs. (38)–(41) over their respective domains and then add the resulting equations. Using, in addition, the divergence theorem, we obtain Z Z Z Z Z @u 1 ½D2 ðxÞruðx;tÞ  ndSx ðx;tÞdV x ¼ fT ðuÞ½x;thx0  xi  T ðuÞ½x0 ;thx  x0 igdV x0 dV x þ ½D2 ðxÞruðþÞ ðx;tÞ  n2 dSx þ X @t Xf1g 2 Xf1g ‘2 @ X2 n‘2 Z Z Z Z Z  1 1 þ T ðuÞ½x;thx0  xidV x0 dV x  T ðuÞ½x0 ;thx  x0 idV x0 dV x þ D2 ðxÞruðþÞ ðx;tÞ  n2 dSx 2 Cf1g 2 Cf2g Xf2g Xf1g C0   Z Z Z  þ D2 ðxÞruðÞ ðx;tÞ  n1 dSx þ ½D2 ðxÞruðx;tÞ  ndSx þ qðx;tÞdV x : ð42Þ f2g

@ C nðC0 [‘2 Þ

‘2

X

We now proceed as follows: 1. Use property (18) to rewrite the first term on the right-hand side of (42). 2. Use a local interface condition (cf. (29d)) over ‘2 to cancel the second and seventh terms on the right-hand side of (42). 3. Cancel the contributions of nonlocal fluxes between material 1 and material 2 within X, from the fourth and fifth terms on the right-hand side of (42). We first write2

1 2

Z

Z

Cf1g 

T ðuÞ½x; thx0  xidV x0 dV x ¼ Xf2g

1 2 þ



1 2

Z

Z Cf2g 

Z

Z

Cf1g 

1 2

Cf2g 

Z

1 2



1 2

Xf1g

Z

Cf1g 

T ðuÞ½x0 ; thx  x0 idV x0 dV x ¼ 

T ðuÞ½x; thx0  xidV x0 dV x

Xf2g nXf2g

Z Z

Z Cf2g 

Cf2g 

Z

Cf1g 

T ðuÞ½x; thx0  xidV x0 dV x ;

T ðuÞ½x0 ; thx  x0 idV x0 dV x

Xf1g nXf1g

T ðuÞ½x0 ; thx  x0 idV x0 dV x

and observe that the first terms on the right-hand sides of the above equations cancel each other upon addition. 2

fjg fjg Given a point x 2 Cfig  , we have Hðx; Þ \ ðX n C0 Þ ¼ Hðx; Þ \ C , i; j ¼ 1; 2; i – j.

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

193

4. Combine the third and the eighth terms on the right-hand side of (42). We now can write (42) as

d dt

Z

uðx; tÞdV x ¼

X

Z

Z Z 1 ½D2 ðxÞruðx; tÞ  ndSx fT ðuÞ½x; thx0  xi  T ðuÞ½x0 ; thx  x0 igdV x0 dV x þ Xf1g 2 Xf1g nXf1g @ Xf2g nC0 Z Z Z Z 1 1 þ T ðuÞ½x; thx0  xidV x0 dV x  T ðuÞ½x0 ; thx  x0 idV x0 dV x f2g 2 Cf1g 2 f2g nXf2g f1g nXf1g X C X  Z  Z  þ D2 ðxÞruðþÞ ðx; tÞ  n2 dSx þ qðx; tÞdV x : ð43Þ C0

X

The left-hand side of (43) represents the rate of change of mass in X. Noticing that the first four terms on the right-hand side of (43) include all net outward mass transport (either local or nonlocal) from X and that the last term contains the mass sources, we conclude that mass conservation in X imposes the constraint

Z



D2 ðxÞruðþÞ ðx; tÞ  n2 dSx ¼ 0:

ð44Þ

C0

Repeating the above exercise for arbitrary subdomains, we get an equivalent condition to (44) with an integral over an arbitrary subregion c0  C0 . Thus, the local/nonlocal interface condition is

 D2 ðxÞruðþÞ ðx; tÞ  n2 ¼ 0;

x 2 C0 :

ð45Þ

Remark 5.1. Condition (45) states that there is no normal local flux across C0 , and that mass transport through C0 occurs only by nonlocal diffusion mechanisms. 5.3. A local/nonlocal interface problem for diffusion Following the derivations in Section 5.1 and Section 5.2, the local/nonlocal interface problem for diffusion is to find the diffusing material density u such that

8 @u > ðx; tÞ ¼ Lnl ðuÞðx; tÞ þ qðx; tÞ; > @t > X1 > > > l > @u > ðx; tÞ ¼ L ðuÞðx; tÞ þ qðx; tÞ; > @t > > > > @u ðx; tÞ ¼ Lnl ðuÞðx; tÞ þ 1 R > T ðuÞ½x; thx0  xidV x0 þ qðx; tÞ; < @t 2 Xf2g Xf1g R l @u ðx; tÞ ¼ L ðuÞðx; tÞ  12 f1g T ðuÞ½x0 ; thx  x0 idV x0 þ qðx; tÞ; > > X > @t > > > uðx; tÞ ¼ 0; > > > > > uðx; 0Þ ¼ u0 ðxÞ; > > > : D2 ðxÞruðþÞ ðx; tÞ  n2 ¼ 0;

with

in X1  ð0; TÞ;

ðaÞ

in X2  ð0; TÞ;

ðbÞ

in Cf1g   ð0; TÞ; ðcÞ in Cf2g   ð0; TÞ; ðdÞ in BX  ð0; TÞ; ðeÞ in X;

ðfÞ

on C0  ð0; TÞ;

ðgÞ

ð46Þ

    BX ¼ BXf1g n C [ @ Xf2g n C0 :

6. Numerical results This section presents numerical studies of one-dimensional static and diffusion interface problems. We would like to investigate and understand computationally the principal features characterizing both nonlocal and local/nonlocal interface problems, and how they compare to classical local interface problems. The computational studies here also demonstrate analytical results presented in previous sections. 6.1. Static problems The nonlocal static problem is to find the material density u such that

(



R X

0

csym ðx0 ; xÞðuðx0 Þ  uðxÞÞdx ¼ qðxÞ; in X;

uðxÞ ¼ gðxÞ;

ðaÞ

in BX; ðbÞ

ð47Þ

with g a function providing Dirichlet-like nonlocal boundary conditions, csym ðx0 ; xÞ given in (6),

cðx0 ; xÞ ¼ veðxÞ ðx0  xÞ

3a 3a

ðeðxÞÞ

1 DðxÞ; jx0  xja

a < 3;

ð48Þ

194

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

Fig. 3. One-dimensional geometry for the numerical simulations. The domain of interest X ¼ ð0; LÞ and its boundary layer BX ¼ BXð1Þ [ BXð2Þ ¼ ðe1 ; 0 [ ½L; L þ e2 Þ. The point x0 represents a material interface.

with ve ðnÞ given in (8), and the source q is given. Eq. (48) is a one-dimensional specialization of (7). In (47), X ¼ ð0; LÞ and BX ¼ BXð1Þ [ BXð2Þ ¼ ðe1 ; 0 [ ½L; L þ e2 Þ (see Fig. 3) with e1 and e2 two possibly different horizons. The corresponding classical local problem is a one-dimensional static version of (29a) with Dirichlet boundary conditions uð0Þ ¼ gð0Þ and uðLÞ ¼ gðLÞ. We assume a two-material system, i.e.,

DðxÞ ¼



D1 D2

x < x0 ; x P x0 ;

eðxÞ ¼



e1 x < x0 ; e2 x P x0 ;

ð49Þ

with x0 representing a material interface, as illustrated in Fig. 3. Our computations are based on finite element discretizations of weak forms, using linear shape functions. We implement a continuous Galerkin method for the discretized weak form of the classical local diffusion problem. In the case of nonlocal diffusion, numerical simulations suggest the presence of a discontinuity in the state variable u at material interfaces. We thus use a discontinuous Galerkin method at the material interface, defining two separate nodes at that point, and a continuous Galerkin method elsewhere. Continuous and discontinuous Galerkin methods in nonlocal models are used in [41,42]. For related weak formulations and well-posedness of nonlocal problems see, e.g., [21,33,34,43–46]. 6.1.1. Patch-test consistency in static nonlocal interface problems In nonlocal interface problems, different horizons represent different length scales, thus producing multiscale nonlocal models. The study of patch-test consistency (cf. Definition 6.1) is an important component in the design of novel concurrent multiscale methods. In solid mechanics, examples appear in, e.g., [24,32,47–49]. We seek to determine the conditions under which the nonlocal problem (47) is patch-test consistent. Definition 6.1. A problem is patch-test consistent if its static solution is given by a linear function, provided zero sources and linear Dirichlet-like nonlocal boundary conditions. Let q  0. We can write (47a) as

Z

X

Z

ð3  aÞDðxÞ 0 ðuðx0 Þ  uðxÞÞdx ðeðxÞÞ3a jx0  xja Z 1 ð3  aÞDðx0 Þ 0 þ ðuðx0 Þ  uðxÞÞdx ¼ 0: 2 X\Hðx;eðx0 ÞÞ ðeðx0 ÞÞ3a jx0  xja

cðx0 ; xÞ þ cðx; x0 Þ 1 0 ðuðx0 Þ  uðxÞÞdx ¼ 2 2

X\Hðx;eðxÞÞ

ð50Þ

We assume uðxÞ ¼ mx þ n with m and n constants, and check whether it satisfies the nonlocal equilibrium equation (50). Note that uðx0 Þ  uðxÞ ¼ uðx þ nÞ  uðxÞ ¼ mn. In the following derivations, we assume e2 6 e1 and a < 2. Let x 2 X1 ¼ ð0; x0  e1 Þ. In this case, DðxÞ ¼ Dðx0 Þ ¼ D1 and eðxÞ ¼ eðx0 Þ ¼ e1 . Then,

Z

X

cðx0 ; xÞ þ cðx; x0 Þ 0 ðuðx0 Þ  uðxÞÞdx ¼ 2

Let x 2 ðx0  e2 ; x0 Þ  Cf1g  . Then,

Z

Z

e1

 e1

Z Z ð3  aÞD1 1 1 x0 x ð3  aÞD1 1 1 e2 ð3  aÞD2 1 mndn a mndn þ a mndn þ 3a 3 a a 2  e1 2 x0 x e3 e1 e1 jnj jnj jnja  e1 2 ( " "  2a #  2a #) 1 3  a D1 x0  x D2 x0  x  : ¼ m 1 1 2 2  a e1 e1 e2 e2

cðx0 ; xÞ þ cðx; x0 Þ 1 0 ðuðx0 Þ  uðxÞÞdx ¼ 2 2 X

Z

ð3  aÞD1 1 m ndn ¼ 0: a e3 jnja 1

e1

The analysis for points in other subdomains of X, follows a similar derivation. We conclude: – Far from the material interface the solution has a linear profile. – Close to the material interface, a linear profile is obtained if uðxÞ ¼ constant (m ¼ 0) or if we have a single-material system (D1 ¼ D2 and e1 ¼ e2 ). Table 1 Static simulation parameters.

X

x0

a

u0

uL

Xfine

hfine

Xcoarse

hcoarse

(0, 100)

50

0

0.0

0.5

[48,52]

0.005

½e1 ; 30 [ ½70; 100 þ e2 

0:5

195

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

0.26

0.26

Num. Nonlocal Soln. Exact Local Soln.

0.258

0.26

Num. Nonlocal Soln. Exact Local Soln.

0.256

0.254

0.254

0.254

0.252

0.252

0.252

0.25

u(x)

0.256

0.25 0.248

0.248

0.246

0.246

0.246

0.244

0.244

0.244

0.242

0.242 48.5

49

49.5

50

50.5

51

51.5

0.24 48

52

0.242 48.5

49

49.5

x

50

50.5

51

51.5

0.24 48

52

0.012

Num. Nonlocal Soln. Exact Local Soln.

0.008

0.008

0.006

du / dx (x)

0.008

du / dx (x)

0.01

0.006

0.004

0.004

0.002

0.002

0.002

50

50.5

51

51.5

52

0 48

50

50.5

51

51.5

52

48.5

49

49.5

x

50

50.5

51

51.5

52

50.5

51

51.5

52

Num. Nonlocal Soln. Exact Local Soln.

0.006

0.004

49.5

49.5

0.012

Num. Nonlocal Soln. Exact Local Soln. 0.01

49

49

x

0.01

48.5

48.5

x

0.012

0 48

Num. Nonlocal Soln. Exact Local Soln.

0.25

0.248

0.24 48

du / dx (x)

0.258

0.256

u(x)

u(x)

0.258

0 48

48.5

49

49.5

x

50

x

Fig. 4. Results for static interface problems with a single diffusion coefficient. The top row shows the solution in Xfine and the bottom row shows the first derivative of the solution in the same region. The results include: (a, d) a nonlocal diffusion model with a single horizon e, (b, e) a nonlocal diffusion model with two different horizons: e1 and e2 , and (c, f) a local/nonlocal diffusion model with e as the horizon of the nonlocal model. Results are compared with the exact solution in classical local diffusion (black-dashed line).

6.1.2. Numerical results for static problems We solve numerically problem (47) with q  0 and gðxÞ given by the classical local solution, i.e.,

gðxÞ ¼

u  u  L 0 x þ u0 ; L

in BX;

for D1 ¼ D2 . Table 1 contains the simulation parameters. To resolve important features around the material interface, we implement a non-uniform finite element mesh as follows. Around the material interface, in Xfine ,3 we use a fine mesh with mesh size hfine . Far from the material interface, in Xcoarse , we use a coarse mesh with mesh size hcoarse (cf. Table 1). Transitional meshes with a monotonically increasing mesh size are used to connect the fine-resolution region Xfine to the coarse-resolution region Xcoarse . In Fig. 4, we present results for systems with a single diffusion coefficient, i.e., D1 ¼ D2 . The top row shows the numerical solution in Xfine and the bottom row shows the first derivative of the solution in the same region. The results include a nonlocal diffusion model with a single horizon e in (a,d), a nonlocal diffusion model with two different horizons e1 – e2 in (b,e), and a local/nonlocal diffusion model with e as the horizon of the nonlocal model in (c,f). In each case, a comparison to the exact solution of classical local diffusion is shown. The main results from Fig. 4 are summarized below. 1. Patch-test consistency. The results in (a,d) demonstrate that nonlocal single-material diffusion models with a single horizon are patch-test consistent. In this case, nonlocal and local models have the same solution. In contrast, problems with two different horizons, as in (b,e) and (c,f), are not patch-test consistent. This is in agreement with Section 6.1.1. 2. Multiple discontinuities in solution gradients. In (b,e), results suggest that two distinct horizons in a nonlocal model create five points of discontinuities in the first derivative of the solution, at x0 ; x0 e1 , and x0 e2 . This implies that indeed three interfaces may be present within the interaction interface C , in contrast to Fig. 1(b). 3. Discontinuity in solution at the material interface. A discontinuity in the solution at x0 appears in (c). This can, in principle, be observed in (b) as well. In contrast to classical local models, discontinuous solutions are thus possible in nonlocal models. 4. Local/nonlocal interface condition. In (c,f), the right-sided first derivative of the solution at x0 is zero. This is consistent with the local/nonlocal interface condition (46g). 5. Local interface condition. In (f), no discontinuity in the first derivative of the solution is observed at x0 þ e, as opposed to result 2 above. This is in agreement with a local interface condition at that point (cf. (29d)). 3

We choose Xfine ¼ ½x0  2 ; x0 þ 2  with x0 the material interface and

 ¼ 1 the maximum horizon value in many simulations.

196

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

Fig. 5 shows results from systems with two diffusion coefficients, i.e., D1 – D2 . The first and third rows show the solution in Xfine , and the second and fourth rows show the first derivative of the solution in the same region. Results are compared with the exact solution from classical local diffusion. In (a, d) and (b, e) we show results for nonlocal diffusion models with a single horizon e. In (c, f) and (g, j), we show results for a nonlocal diffusion model with two different horizons e1 – e2 , alternating the choice of region with largest horizon. In (h,k) and (i,l), we show results for local/nonlocal diffusion models with e as the horizon of the nonlocal diffusion model. The local diffusion model is on the right half of the domain in (h, k) and on the

0.34

0.34

Num. Nonlocal Soln. Exact Local Soln.

0.338

0.34

Num. Nonlocal Soln. Exact Local Soln.

0.336

0.334

0.334

0.334

0.332

0.332

0.332

0.33

u(x)

0.336

0.33 0.328

0.328

0.326

0.326

0.326

0.324

0.324

0.324

0.322

0.322 48.5

49

49.5

50

50.5

51

51.5

0.32 48

52

0.322 48.5

49

49.5

x

51

51.5

0.32 48

52

0.008

0.008

0.008

0.006

du / dx (x)

0.01

0.006

0.004

0.004

0.002

0.002

0.002

49.5

50

50.5

51

51.5

0 48

52

48.5

49

49.5

x

50

50.5

51

51.5

0 48

52

0.34

Num. Nonlocal Soln. Exact Local Soln.

0.338

0.338

0.334

0.334

0.332

0.332

0.332

u(x)

0.334

u(x)

0.336

0.33 0.328

0.328

0.326

0.326

0.326

0.324

0.324

0.324

0.322

0.322 49

49.5

50

50.5

51

51.5

0.32 48

52

49

49.5

50

50.5

51

51.5

0.32 48

52

0.012

0.008

du / dx (x)

0.008

du / dx (x)

0.008

0.006

0.006

0.004

0.004

0.002

0.002

0.002

50

x

50.5

51

51.5

52

0 48

50.5

51

51.5

52

50.5

51

51.5

52

50.5

51

51.5

52

Num. Nonlocal Soln. Exact Local Soln.

48.5

49

49.5

50

0.006

0.004

49.5

50

Num. Nonlocal Soln. Exact Local Soln. 0.01

49

49.5

0.012

Num. Nonlocal Soln. Exact Local Soln. 0.01

48.5

49

x

0.01

0 48

48.5

x

Num. Nonlocal Soln. Exact Local Soln.

52

0.322 48.5

x

0.012

51.5

0.33

0.328

48.5

51

0.34

Num. Nonlocal Soln. Exact Local Soln.

0.336

0.33

50.5

x

0.336

0.32 48

50

Num. Nonlocal Soln. Exact Local Soln.

x

0.34 0.338

49.5

0.006

0.004

49

49

0.012

Num. Nonlocal Soln. Exact Local Soln. 0.01

48.5

48.5

x

0.01

du / dx (x)

du / dx (x)

50.5

0.012

Num. Nonlocal Soln. Exact Local Soln.

u(x)

50

x

0.012

0 48

Num. Nonlocal Soln. Exact Local Soln.

0.33

0.328

0.32 48

du / dx (x)

0.338

0.336

u(x)

u(x)

0.338

48.5

49

49.5

50

x

50.5

51

51.5

52

0 48

48.5

49

49.5

50

x

Fig. 5. Results for static interface problems with two diffusion coefficients. The first and third rows show the solution in Xfine . The second and fourth rows show the first derivative of the solution in Xfine . In (a, d) and (b, e), results are for nonlocal diffusion models with a single horizon e. In (c, f) and (g, j) we present results for nonlocal diffusion models with two horizons: e1 and e2 . Subfigures (h, k) and (i, l) show results for local/nonlocal diffusion models with e as the horizon of the nonlocal model. Results are compared with the exact solution in classical local diffusion (black-dashed line).

197

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

0.334

0.334

Num. Nonlocal Soln. Exact Local Soln.

0.334

Num. Nonlocal Soln. Exact Local Soln.

0.3338

0.3336

0.3334

0.3334

0.3334

0.3332

0.3332

0.3332

u(x)

0.3336

0.333

0.333

0.333

0.3328

0.3328

0.3328

0.3326

0.3326

0.3326

0.3324

0.3324

0.3324

0.3322

0.3322

0.332 49.8

49.85

49.9

49.95

50

50.05

50.1

50.15

50.2

0.3322

0.332 49.8

49.85

49.9

49.95

50

x

50.1

50.15

50.2

0.332 49.8

8

8

6

4

du / dx (x)

6

du / dx (x)

6

5

5

4

3 50

50.05

50.1

50.15

50.2

49.8

50.05

50.1

50.15

50.2

x 10

5

4

3 49.95

50

Num. Nonlocal Soln. Exact Local Soln. 7

49.9

49.95

Num. Nonlocal Soln. Exact Local Soln. 7

49.85

49.9

−3

x 10

7

49.8

49.85

x

−3

x 10

Num. Nonlocal Soln. Exact Local Soln.

du / dx (x)

50.05

x

−3

8

Num. Nonlocal Soln. Exact Local Soln.

0.3338

0.3336

u(x)

u(x)

0.3338

3 49.85

49.9

49.95

x

50

50.05

50.1

50.15

50.2

x

49.8

49.85

49.9

49.95

50

50.05

50.1

50.15

50.2

x

Fig. 6. Results for static interface problems with two diffusion coefficients, using a uniform horizon e in the nonlocal diffusion models. The top row shows the solution in ½49:8; 50:2 and the bottom row shows the first derivative of the solution in the same region. Results are compared with the exact solution in classical local diffusion (black-dashed line).

left half of the domain in (i, l). We obtain the same results 2–5 as for Fig. 4. In addition, we observe that even when the horizon is uniform throughout the domain, a discontinuity in the diffusion coefficient at x0 generates discontinuities in the first derivative of the solution at x0 e (cf. (d, e)). Fig. 6 shows results from nonlocal diffusion models with a single horizon e and two diffusion coefficients D1 – D2 . A larger horizon is used in the left column (a, d) and a very small horizon is used in the right column (c, f). This allows one to see the effects of the value of the horizon upon the solution. The top row shows the solution in ½49:8; 50:2 and the bottom row shows the first derivative of the solution in the same region. The results demonstrate that the solution of the nonlocal diffusion model converges to the solution of the classical local diffusion model in the limit as e ! 0. Furthermore, we observe that the first derivative of the classical local solution at the material interface satisfies ðÞ

du dx

ðþÞ

ðx0 Þ ¼ 2

du dx

ðx0 Þ;

in agreement with (29d), assuming n1 ¼ n2 ¼ 1. 6.2. Diffusion problems The nonlocal diffusion problem is to find the material density u such that

R 8 @u 0 ðx; tÞ ¼ csym ðx0 ; xÞðuðx0 ; tÞ  uðx; tÞÞdx þ qðx; tÞ; > @t > X > R > 0 < c ðx0 ; xÞðuðx0 ; tÞ  uðx; tÞÞdx ¼ sðx; tÞ; BX sym > uðx; 0Þ ¼ u0 ðxÞ; > > > R :1 R uðx; tÞdx ¼ jX1 j X u0 ðxÞdx; jXj X

in X  ð0; TÞ; ðaÞ in X  ð0; TÞ; ðbÞ in X;

ðcÞ

in ð0; TÞ;

ðdÞ

ð51Þ

with sðx; tÞ a function providing Neumann-like nonlocal boundary conditions, u0 ðxÞ an initial condition, csym ðx0 ; xÞ given in (6) with cðx0 ; xÞ in (48), and the source q is given. The additional condition (51d) ensures well-posedness of the problem. In (51), X ¼ ð0; LÞ. The corresponding classical local diffusion problem is a one-dimensional version of (29) with Neumann boundary conditions (cf. Appendix A) and the additional condition (51d). If s  0 and q  0, we can express (51a) as [21]

d dt

Z X

u2 ðx; tÞdx ¼ 

Z Z X

2

0

csym ðx0 ; xÞðuðx0 ; tÞ  uðx; tÞÞ dx dx < 0; X

198

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

Table 2 Diffusion simulation parameters.

X

x0

a

Dt

N elem

T

(0, 1)

0.5

0

0.001

1000

0.25

showing that the nonlocal model indeed represents diffusion. Due to (51d), we can write

d dt

Z

 Þ2 dx ¼ ðuðx; tÞ  u

X

d dt

Z

 :¼ u2 ðx; tÞdx < 0 with u

X

1 jXj

Z

u0 ðxÞdx

ð52Þ

X

 , the mean of the field, is time-independent. We see that the variance of the field gets smaller in time. where u As in the static simulations, computations are based on finite element discretizations of weak forms using linear shape functions. We use a continuous Galerkin method in classical local diffusion. In nonlocal diffusion, we use a discontinuous Galerkin method at the material interface and a continuous Galerkin method elsewhere. We discretize the domain with a uniform mesh with N elem elements. We assume no external sources, i.e., q  0, and choose s  0 in (51b). Accordingly, we impose zero Neumann boundary conditions in the classical local diffusion problem, i.e., @u ð0; tÞ ¼ @u ðL; tÞ ¼ 0 for t 2 ð0; TÞ. @x @x The system is evolved in time using a backward Euler scheme in v-form (cf. [50, pg. 460]) with time step Dt and final time T. Table 2 describes the simulation parameters. The initial condition is chosen as

u0 ðxÞ ¼



0

x < x0 ;

1 x P x0 :

1.2

1

1

1

0.8

0.8

0.8

0.6

0.6

0.6

0.4

0.4

0.2

0.2

0

0

−0.2 0

0.2

0.4

0.6

0.8

−0.2 0

1

1.2

u(x)

1.2

u(x)

u(x)

In Fig. 7, we present diffusion simulations with a single diffusion coefficient, i.e., D1 ¼ D2 . We compare (a) a classical local diffusion model, (b) a nonlocal diffusion model with a single horizon e, (c, d) nonlocal diffusion models with two horizons: e1 and e2 , and (e,f) local/nonlocal diffusion models with e as the horizon of the nonlocal model. In (e) the local model is on the right half of the domain and in (f) on the left half. Colored lines represent material density profiles at different times [20]. A time of 0:02 separates the lines, with the initial profile (in green) at t ¼ 0, and the final profile (in black) at t ¼ 0:24.

0.2 0

0.2

0.4

x

0.6

0.8

−0.2 0

1

1

1

0.8

0.8

0.8

0.6

0.6

0.6

0.4

0.4

0.2

0.2

0

0

0.4

0.6

x

0.8

1

−0.2 0

0.6

0.8

1

0.6

0.8

1

1.2

u(x)

1.2

1

0.2

0.4

x

1.2

−0.2 0

0.2

x

u(x)

u(x)

0.4

0.4 0.2 0

0.2

0.4

0.6

x

0.8

1

−0.2 0

0.2

0.4

x

Fig. 7. Field values for one-dimensional diffusion simulations with single diffusion coefficients. We compare (a) the classical local diffusion model, (b) the nonlocal diffusion model with a single horizon e, (c, d) nonlocal diffusion models with two horizons: e1 and e2 , and (e,f) local/nonlocal diffusion models with e as the horizon of the nonlocal model, which is on the left half in (e) and on the right half in (f). Each colored line represents a different time [20]. Dotted vertical lines indicate the material interface x0 , the points x0 e in (b, e, f), and the points x0 e1 and x0 e2 in (c, d). Simulation choices are described in Table 2.

199

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

The main results of Fig. 7 are summarized below. 1. Discontinuity preservation. In nonlocal diffusion the initial discontinuity is preserved in time (see (b)). This is in contrast to classical local diffusion (a), where initial discontinuities are immediately smoothed out. The size of the jump in the nonlocal solution decreases in time, however. 2. Discontinuities in solution gradients. As in static simulations, nonlocal diffusion models with two horizons have discontinuities in the first derivative of the solution at x0 e1 and x0 e2 , in addition to the discontinuity at x0 . Such discontinuities are clearly observed on the side having the largest horizon, i.e., on the left side in (c) and on the right side in (d). 3. Local/nonlocal interface condition. In (e,f), the one-sided first derivative of the solution at x0 , on the local side (right side in (e) and left side in (f)), is zero. This is in agreement with the local/nonlocal interface condition (46g) and holds for all times. 4. Local interface condition. In (e,f), no discontinuities in the first derivative of the solution seem to appear at a point e away from x0 , in the local side (right side in (e) and left side in (f)). This is in agreement with a classical local interface condition at that point (cf. (29d)). In Fig. 8, we reproduce the simulations of Fig. 7 with two different diffusion coefficients D1 – D2 . Very similar results are obtained. The effect of disparate diffusion coefficients is however observed in the line slopes around the material interface, as seen by comparing Figs. 7(a) and 8(a). The convergence of the nonlocal diffusion model to the classical local diffusion model is demonstrated in Fig. 9 for two diffusion coefficients D1 – D2 . Analogous convergence studies for single-material nonlocal diffusion models appear in [20]. We compare a classical local diffusion model (a) to nonlocal diffusion models with a single horizon e (b, c, d, e, f). As in Figs. 8 and 7, colored lines represent material density profiles at different times [20]. The main results of Fig. 9 are summarized below.

1.2

1

1

1

0.8

0.8

0.8

0.6

0.6

0.6

0.4

0.4

0.2

0.2

0

0

−0.2 0

0.2

0.4

0.6

0.8

−0.2 0

1

1.2

u(x)

1.2

u(x)

u(x)

1. Convergence in the local limit. For small horizons, e.g., in (b), the solution of the nonlocal diffusion model converges to the solution of the classical local diffusion model (a). 2. Horizon-dependent diffusion rate. Following (52), we expect the variance of the solution to decrease in time until the system reaches a steady-state solution, in this case u  0:5. The diffusion rate in nonlocal diffusion seems to depend on the horizon. Models with larger horizons take larger times to reach the steady-state profile (compare (b) through (f)). 3. Horizon-dependent smoothing rate. As in Figs. 7 and 8, a smoothing effect is observed, even though nonlocal diffusion models tend to preserve discontinuities in time. The smoothing rate depends upon the horizon. The larger the horizon, the smaller the smoothing rate (compare (b) through (f)).

0.2 0

0.2

0.4

x

0.6

0.8

−0.2 0

1

1

1

0.8

0.8

0.8

0.6

0.6

0.6

0.4

0.4

0.2

0.2

0

0

0.4

0.6

x

0.8

1

−0.2 0

0.6

0.8

1

0.6

0.8

1

1.2

u(x)

1.2

1

0.2

0.4

x

1.2

−0.2 0

0.2

x

u(x)

u(x)

0.4

0.4 0.2 0

0.2

0.4

0.6

x

0.8

1

−0.2 0

0.2

0.4

x

Fig. 8. Field values for one-dimensional diffusion simulations with two diffusion coefficients. These results extend the results of Fig. 7. Refer to the caption of Fig. 7 for a detailed explanation.

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

1.2

1

1

1

0.8

0.8

0.8

0.6

0.6

0.6

0.4

0.4

0.2

0.2

0

0

−0.2 0

0.2

0.4

0.6

0.8

−0.2 0

1

1.2

u(x)

1.2

u(x)

u(x)

200

0.2 0

0.2

0.4

x

0.6

0.8

−0.2 0

1

1

1

0.8

0.8

0.8

0.6

0.6

0.6

0.4

0.4

0.2

0.2

0

0

0.4

0.6

x

0.8

1

−0.2 0

0.6

0.8

1

0.6

0.8

1

1.2

u(x)

1.2

1

0.2

0.4

x

1.2

−0.2 0

0.2

x

u(x)

u(x)

0.4

0.4 0.2 0

0.2

0.4

0.6

x

0.8

1

−0.2 0

0.2

0.4

x

Fig. 9. Field values for one-dimensional diffusion simulations with two different diffusion coefficients. We compare the classical local diffusion model (a) to nonlocal diffusion models with a single horizon e (b, c, d, e, f). Each colored line represents a different time [20]. Dotted vertical lines indicate the material interface x0 and the points x0 e in (b, c, d). Simulation choices are described in Table 2.

7. Conclusions In this paper we explored interface problems in nonlocal diffusion and demonstrated how to reformulate the classical treatment of interface problems in the presence of long-range interactions. We focused on connections to classical local interface problems for diffusion and on the analysis of systems with sharp transitions between local and nonlocal domains. In nonlocal interface problems we must differentiate between a material interface and an interaction interface. The former is the boundary separating domains with different model or material parameters, and the latter is a volumetric region surrounding the material interface that contains all long-range inter-material interactions. These two interfaces coincide in classical local interface problems. Domain boundaries are also treated differently in local and nonlocal problems. In classical local problems, a surface boundary is used, where usually Dirichlet or Neumann boundary conditions are imposed. In nonlocal problems, however, a boundary layer is used to impose nonlocal boundary conditions. In classical local diffusion interface problems, subdomain governing equations are not valid at material interfaces, and special interface conditions must be introduced to satisfy conservation principles. In contrast, nonlocal diffusion equations are valid everywhere in the domain, regardless of the presence of material discontinuities. This makes nonlocal models naturally applicable to interface problems without any additional interface conditions. We showed that nonlocal diffusion interface problems converge to classical local diffusion interface problems, in the limit of vanishing nonlocality or local limit. This analysis was performed using strong forms of governing equations. We obtained an explicit relation between local and nonlocal diffusion tensors. In Appendix A, we also demonstrate that the nonlocal flux converges to its local counterpart in the local limit. A principal result in this paper concerns the formulation of interface problems for diffusion with sharp transitions between local and nonlocal domains. We cast such systems as nonlocal interface problems with two horizons, one finite and the other zero. By assuming a model with two finite horizons and then letting one of them go to zero, we derived governing equations for each of the relevant subdomains in the system. We found that an interface condition is required at local/ nonlocal material interfaces, and we used a conservation principle to derive it. That interface condition suggests that in this type of local/nonlocal interface problem, mass transport across material interfaces only occurs through nonlocal mechanisms. We demonstrated computationally, through one-dimensional static and diffusion examples, the analytical results presented in this paper. In particular, we showed numerically the convergence of nonlocal diffusion models to local diffusion models, and the interface conditions for both local and local/nonlocal interface problems. Numerical simulations also suggested that, in contrast to classical local models, discontinuous solutions may appear in nonlocal interface

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

201

problems. Furthermore, we demonstrated that discontinuities in diffusion coefficients or horizons in nonlocal models induce discontinuities in solution gradients. Diffusion simulations also showed that nonlocal models can preserve discontinuities, as opposed to classical local ones, with smoothing and diffusion rates depending on the value of the horizon. Important insights toward multiscale modeling are obtained in this paper. Diffusion simulations demonstrated that we can represent systems with part of the domain governed by a nonlocal model and the other part by a local model. However, analytical and numerical results for static problems showed that a nonlocal system with two distinct horizons does not in general pass a patch test. This suggests a mismatch of length scales typical of local/nonlocal coupling methods in multiscale material modeling, such as AtC coupling methods in solid mechanics. Development of multiscale descriptions of single materials is of major interest in the concurrent multiscale modeling community. This is a topic for future work. Acknowledgements This research was supported by DOE Grant DE-SC0004970 at Florida State University, by DOE Grant DE-FG02-05ER25701 at the University of Texas, and by the Laboratory Directed Research and Development program at Sandia National Laboratories. Sandia is a multiprogram laboratory operated by Sandia Corporation, a Lockheed Martin Company, for the United States Department of Energy under contract DE-AC0494-AL85000. We acknowledge helpful discussions with David Littlewood, Stewart Silling, Jakob Ostien, and Qiang Du. Pablo Seleson acknowledges support from the ICES Postdoctoral Fellowship Program and useful discussions with Serge Prudhomme and Leszek Demkowicz. Appendix A. Convergence of the nonlocal flux to the classical local flux We derive the convergence of the nonlocal flux to its classical local counterpart, in the local limit. Following (13) and (4), we express the nonlocal flux from XB into XA as

RXA ;XB ¼

Z XA

Z

csym ðx0 ; xÞðuðx0 ; tÞ  uðx; tÞÞdV x0 dV x :

ðA:1Þ

XB

Let XA and XB be two non-overlapping connected domains such that XA \ XB ¼ C0 . For the purpose of connecting local and nonlocal fluxes, we assume that, for x 2 XA and x0 2 XB \ Hðx; eÞ, the following holds:

  uðx0 ; tÞ ¼ uðx; tÞ þ ððx0  xÞ  rÞuðx; tÞ þ O kx0  xk2 ; 0

ðA:2Þ

0

Kðx; x  xÞ ¼ Kðx; 0Þ þ Oðkx  xkÞ;

ðA:3Þ

Kðx0 ; x  x0 Þ ¼ Kðx; 0Þ þ Oðkx0  xkÞ;

ðA:4Þ

where, for simplicity, we assume a single-material model with horizon e. Under these assumptions and following similar algebraic manipulations as in Section 4.1, we obtain (cf. (7))

RXA ;XB ¼

Z XA

cd;e a ðKÞ‘k ðx; 0Þ

@u ðx; tÞ @xj

Z Hðx;eÞ\XB







x0‘  x‘ x0k  xk x0j  xj kx0  xk2þa

dV x0 dV x þ OðeÞ:

ðA:5Þ

To simplify the derivations, we take XA as a rectangular prism (in 3D), a rectangle (in 2D), or a line segment (in 1D), although results can be derived for general geometries. We introduce the coordinate system illustrated in Fig. 10, where we identify x1 with x; x2 with y, and x3 with z. We denote by z0 the z-coordinate of the surface C0 connecting XA and XB .

Fig. 10. Coordinate system for the convergence analysis of the nonlocal flux to the classical local flux. The z-axis is oriented along the x3 -direction. The angle h is the inclination and / the azimuthal angle.

202

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

We now can write (A.5) as

RXA ;XB

Z

@u ¼ ce ðKÞ‘k ðxc ; 0Þ ðxc ; tÞ @xj Sðzc Þ d;a

"Z

0





x‘  x‘ x0k  xk x0j  xj

Z

z0

z0  e

kx0  xk2þa

Hðx;eÞ\XB

# dV x0 dz dSx þ OðeÞ;

ðA:6Þ

where xc :¼ ðx; y; zc Þ; zc 2 ½z0  e; z0 , and SðzÞ is a surface in the x–y plane corresponding to the coordinate z. In (A.6), we used the first mean value theorem for integration with respect to the z coordinate, assuming the following: 1. The functions Kðx; 0Þ and ruðx; tÞ are continuous in XA . 2. The integrals inside the square brackets in (A.6) are bounded, i.e.,

Z

z0

z0 e



Z





x0‘  x‘ x0k  xk x0j  xj kx0  xk2þa

Hðx;eÞ\XB

dV x0 dz < 1

ðA:7Þ

and the integrand in (A.7) does not change sign. The lack of sign change in Property 2 above holds due to the symmetry of the integration domain around the z-axis; in particular, integrals involving odd powers of ðx0j  xj Þ; j ¼ 1; 2, vanish. In fact, we only need to compute the following two integrals:

Z

z0

z0 e

Z

n33

Hðx;eÞ\XB

knk

2þa

Z

dV x0 dz;

Z

z0

z0 e

n21 n3

Hðx;eÞ\XB

knk2þa

dV x0 dz;

with n :¼ x0  x, because, due to symmetry, we have

Z

z0

z0 e

Z

Z

n21 n3

Hðx;eÞ\XB

knk

dV x0 dz ¼ 2þa

Z

z0 z0 e

n22 n3

Hðx;eÞ\XB

knk2þa

dV x0 dz:

ðA:8Þ

The calculations below assume a < d þ 2 with d the dimension (cf. Section 4.2). Derivations in 3D. The first integral is calculated directly as follows:

Z

z0

z0 e

Z

Z

n33

Hðx;eÞ\XB

knk

dV x0 dz ¼ 2þa

Z

z0 z0 e

Z 2p Z

e

z0 z

0

cos1 ð

0

z0 z r

Þ ðr cos hÞ3 2p e5a r 2 sin hdhd/drdz ¼ : 2þ a r 5 ð5  aÞ

To calculate the second integral, we observe that

Z

z0

z0 e

Z

n21 n3

Hðx;eÞ\XB

1 dV x0 dz ¼ 2 knk2þa

(Z

Z

z0

z0 e

Hðx;eÞ\XB

n3 dV x0 dz  knka

Z

z0 z0 e

Z

n33

Hðx;eÞ\XB

knk2þa

) dV x0 dz ;

where (A.8) is used. Using similar calculations as for the first integral, we compute the first term inside the curly brackets. We then obtain

Z

z0

z0 e

Z



1 2p e5a 2p e5a 2p e5a : ¼  2 3 ð5  aÞ 5 ð5  aÞ 15 ð5  aÞ

n21 n3

Hðx;eÞ\XB

knk

dV x0 dz ¼ 2þa

Derivations in 2D. We consider a similar coordinate system as in Fig. 10, but using instead polar coordinates (no azimuthal angle /). We take the ‘‘z-direction’’ as the x-direction (thus, replacing z0 by x0 ). The following results are obtained for a ¼ 0; 1; 2; 3.4 We start by computing the second integral:

Z

x0

x0 e

Z Hðx;eÞ\XB

n1 n22

dV x0 dx ¼ knk2þa

Z

Z

x0 x0 e

Z

e

x0 x

x0 x r

cos1 ð

Þ r cos hðr sin hÞ2 pe4a rdhdrdx ¼ : 2þa x0 x r 8ð4  aÞ 1  cos ð r Þ

To obtain the first integral, we first compute

Z

x0

x0 e

Z Hðx;eÞ\XB

n1 dV x0 dx ¼ knka

Z

Z

x0

Z

e

x0  e

x0 x

Z

Z

cos1 ð

x0 x r

Þ r cos h pe4a rdhdrdx ¼ : a x x r 2ð4  aÞ  cos1 ð 0r Þ

Using this result, we get

Z

x0

x0 e

Z Hðx;eÞ\XB

n31

dV x0 dx ¼ knk2þa

x0 x0 e

Hðx;eÞ\XB

n1 dV x0 dx  knka

Z

x0

x0 e

Z Hðx;eÞ\XB

n1 n22 2þa

knk

dV x0 dx ¼

3pe4a : 8ð4  aÞ

4 The calculation appears to be more complicated for the two-dimensional integrals than for the three-dimensional integrals. Therefore, we present results for specific values of a, i.e., a ¼ 0; 1; 2; 3, rather than for a general a.

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

203

Derivations in 1D. Here, we only need to compute the following integral:

Z

Z

x0

x0 e

xþe

x0

n e3a 0 : a dx dx ¼ 3a jnj

In the following, we perform calculations in 3D. Results in 2D and 1D are obtained in an analogous way. We can now write (A.6) as a RXA ;XB ¼ c3; e

(

Z

K33 Sðzc Þ

@u @x3

Z

z0 z0 e

Z

n33 Hðx;eÞ\XB

knk

2þa

!  !)  Z z0 Z @u @u @u n21 n3 dV x0 dz þ ðK11 þ K22 Þ þ ðK13 þ K31 Þ þ ðK23 þ K32 Þ dV x0 dz dSx 2þ a @x3 @x1 @x2 z0 e Hðx;eÞ\XB knk

þ OðeÞ  Z 



Z 5a 1 @u 1 @u 1 @u @u @u @u a 4p e dSx þ OðeÞ ¼ dSx þ OðeÞ ð2K33 þ trðKÞÞ þ ðK13 þ K31 Þ þ ðK23 þ K32 Þ D33 þ D31 þ D32 ¼ c3; e 15 ð5  aÞ Sðzc Þ 2 @x3 2 @x1 2 @x2 @x3 @x1 @x2 Sðzc Þ Z Z ðDðxc Þruðxc ; tÞÞ  ndSx þ OðeÞ! ðDðxÞruðx; tÞÞ  ndSx ; ¼ e!0 C0

Sðzc Þ

@u @xi

@u @xi

where Kij :¼ ðKÞij ðxc ; 0Þ, :¼ ðxc ; tÞ, Dij :¼ ðDÞij ðxc Þ, and n ¼ ð0; 0; 1Þ is an outward unit normal vector to XB . We also a used relations (27) for D and (26) for c3; e . Recall that the local flux vector is defined as

rðx; tÞ :¼ DðxÞruðx; tÞ:

ðA:9Þ

The minus sign in (A.9) indicates that particles move towards regions of lower concentration. We then obtain

lim RXA ;XB ¼ e!0

Z

rðx; tÞ  ndSx ;

ðA:10Þ

C0

which demonstrates that the total nonlocal mass transport from XB into XA across C0 reduces, in the local limit, to the classical local flux across C0 from XB into XA . References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32]

A. Buades, B. Coll, J.M. Morel, Image denoising methods. A new nonlocal principle, SIAM Rev. 52 (1) (2010) 113–147. G. Gilboa, S. Osher, Nonlocal operators with applications to image processing, Multiscale Model. Simul. 7 (2008) 1005–1028. D. Appelbaum, Lévy Processes and Stochastic Calculus, Cambridge studies in advance mathematics, vol. 93, Cambridge University Press, 2004. M.T. Barlow, R.F. Bass, Z.-Q. Chen, M. Kassmann, Non-local Dirichlet forms and symmetric jump processes, Trans. Am. Math. Soc. 361 (4) (2009) 1963– 1999. C. Bardos, R. Santos, R. Sentis, Diffusion approximation and computation of the critical size, Trans. Am. Math. Soc. 284 (2) (1984) 617–649. J.-G. Liu, L. Mieussens, Analysis of an asymptotic preserving scheme for the linear kinetic equations in the diffusion limit, SIAM J. Numer. Anal. 48 (4) (2010) 1474–1491. P.W. Bates, A. Chmaj, A discrete convolution model for phase transitions, Arch. Ration. Mech. Anal. 150 (1999) 281–305. P.W. Bates, A. Chmaj, An integrodifferential model for phase transitions: stationary solutions in higher space dimensions, J. Stat. Phys. 95 (1999) 1119– 1139. D.B. Duncan, M. Grinfeld, I. Stoleriu, Coarsening in an integro–differential model of phase transitions, Eur. J. Appl. Math. 11 (2000) 561–572. A.C. Eringen, Nonlocal Continuum Field Theories, Springer, New York, 2002. E. Kröner, Elasticity theory of materials with long range cohesive forces, Int. J. Solids Struct. 3 (5) (1967) 731–742. I.A. Kunin, Elastic media with microstructure I: one-dimensional models, Springer Series in Solid State Sciences, vol. 26, Springer-Verlag, Berlin, 1982. I.A. Kunin, Elastic media with microstructure II: three-dimensional models, Springer Series in Solid State Sciences, vol. 44, Springer-Verlag, Berlin, 1983. D. Rogula, Introduction to nonlocal theory of material media, in: D. Rogula (Ed.), Nonlocal Theory of Material Media, Springer-Verlag, Berlin, 1982, pp. 125–222. S.A. Silling, Reformulation of elasticity theory for discontinuities and long-range forces, J. Mech. Phys. Solids 48 (2000) 175–209. S.A. Silling, M. Epton, O. Weckner, J. Xu, E. Askari, Peridynamic states and constitutive modeling, J. Elast. 88 (2007) 151–184. Z.P. Bazˇant, M. Jirásek, Nonlocal integral formulations of plasticity and damage: survey of progress, J. Engrg. Mech. 128 (11) (2002) 1119–1149. F. Bobaru, M. Duangpanya, The peridynamic formulation for transient heat conduction, Int. J. Heat Mass Transfer 53 (2010) 4047–4059. F. Andreu, J.M. Mazón, J.D. Rossi, J. Toledo, A nonlocal p-Laplacian evolution equation with nonhomogeneous Dirichlet boundary conditions, SIAM J. Math. Anal. 40 (5) (2009) 1815–1851. N. Burch, R. Lehoucq, Classical, nonlocal, and fractional diffusion equations on bounded domains, Int. J. Multiscale Comput. Engrg. 9 (2011) 661–674. Q. Du, M. Gunzburger, R.B. Lehoucq, K. Zhou, Analysis and approximation of nonlocal diffusion problems with volume constraints, SIAM Rev. 54 (4) (2012) 667–696. P. Fife, Some nonclassical trends in parabolic and parabolic-like evolutions, in: M. Kirkilionis, S. Krömker, R. Rannacher, F. Tomi (Eds.), Trends in nonlinear analysis, Springer, Berlin, 2003, pp. 153–191. Ch. 3. J.T. Oden, E.B. Becker, G.F. Carey, Finite Elements: An Introduction, Prentice Hall, Volume I, 1981. W.A. Curtin, R.E. Miller, Atomistic/continuum coupling in computational materials science, Model. Simul. Mater. Sci. Engrg. 11 (2003) R33–R68. F. Bobaru, Y.D. Ha, Adaptive refinement and multiscale modeling in 2D peridynamics, Int. J. Multiscale Comput. Engrg. 9 (6) (2011) 635–660. R.W. Macek, S.A. Silling, Peridynamics via finite element analysis, Finite Elements in Anal. Des. 43 (2007) 1169–1178. W. Liu, J.-W. Hong, A coupling approach of discretized peridynamics with finite element method, Comput. Methods Appl. Mech. Engrg. 245–246 (2012) 163–175. B. Kilic, E. Madenci, Coupling of peridynamic theory and the finite element method, J. Mech. Mater. Struct. 5 (5) (2010) 707–733. G. Lubineau, Y. Azdoud, F. Han, C. Rey, A. Askari, A morphing strategy to couple non-local to local continuum mechanics, J. Mech. Phys. Solids 60 (2012) 1088–1102. F. Han, G. Lubineau, Coupling of nonlocal and local continuum models by the Arlequin approach, Int. J. Numer. Methods Engrg. 89 (2012) 671–685. P. Seleson, Peridynamic multiscale models for the mechanics of materials: constitutive relations, upscaling from atomistic systems, and interface problems, Ph.D. thesis, Florida State University, 2010. P. Seleson, S. Beneddine, S. Prudhomme, A force-based coupling scheme for peridynamics and classical elasticity, Comput. Mater. Sci. 66 (2013) 34–49.

204

P. Seleson et al. / Comput. Methods Appl. Mech. Engrg. 266 (2013) 185–204

[33] M. Gunzburger, R.B. Lehoucq, A nonlocal vector calculus with application to nonlocal boundary value problems, Multiscale Model. Simul. 8 (5) (2010) 1581–1598. [34] B. Aksoylu, M.L. Parks, Variational theory and domain decomposition for nonlocal problems, Appl. Math. Comput. 217 (2011) 6498–6515. [35] Q. Du, M. Gunzburger, R.B. Lehoucq, K. Zhou, A nonlocal vector calculus, nonlocal volume-constrained problems, and nonlocal balance laws, Math. Models Methods Appl. Sci. 23 (2013) 493–540. [36] E. Oterkus, E. Madenci, Peridynamics for failure prediction in composites, in: 53th AIAA/ASME/ASCE/AHS/ASC Structures, structural dynamics, and materials conference, Honolulu, HI, no. 2012–1692, 2012. [37] E. Oterkus, E. Madenci, O. Weckner, S. Silling, P. Bogert, Combined finite element and peridynamic analyses for predicting failure in a stiffened composite curved panel with a central slot, Composite Struct. 94 (2012) 839–850. [38] J. Xu, A. Askari, O. Weckner, S. Silling, Peridynamic analysis of impact damage in composite laminates, J. Aerosp. Engrg., SPECIAL ISSUE: Impact Mechanics of Composite Materials for Aerospace Application, 21 (2008) 187–194. [39] Q. Du, K. Zhou, Mathematical analysis for the peridynamic nonlocal continuum theory, ESAIM: Math. Model. Numer. Anal. 45 (2011) 217–234. [40] Q. Du, J.R. Kamm, R.B. Lehoucq, M.L. Parks, A new approach for a nonlocal, nonlinear conservation law, SIAM J. Appl. Math. 72 (1) (2012) 464–487. [41] X. Chen, M. Gunzburger, Continuous and discontinuous finite element methods for a peridynamics model of mechanics, Comput. Methods Appl. Mech. Engrg. 200 (2011) 1237–1250. [42] Q. Du, L. Ju, L. Tian, K. Zhou, A posteriori error analysis of finite element method for linear nonlocal diffusion and peridynamic models, Math. Comput. 82 (2013) 1889–1922. [43] B. Hinds, P. Radu, Dirichlet’s principle and wellposedness of solutions for a nonlocal pLaplacian system, Appl. Math. Comput. 219 (2012) 1411–1419. [44] B. Aksoylu, T. Mengesha, Results on nonlocal boundary value problems, Numer. Funct. Anal. Optim. 31 (2010) 1301–1317. [45] F. Andreu, J.M. Mazón, J.D. Rossi, J. Toledo, Nonlocal diffusion problems, Mathematical Surveys and Monographs, vol. 165, American Mathematical Society, 2010. [46] E. Emmrich, O. Weckner, On the well-posedness of the linear peridynamic model and its convergence towards the Navier equation of linear elasticity, Commun. Math. Sci. 5 (4) (2007) 851–864. [47] S. Badia, M. Parks, P. Bochev, M. Gunzburger, R. Lehoucq, On atomistic-to-continuum coupling by blending, Multiscale Model. Simul. 7 (1) (2008) 381– 406. [48] B.V. Koten, M. Luskin, Analysis of energy-based blended quasi-continuum approximations, SIAM J. Numer. Anal. 49 (2011) 2182–2209. [49] P. Seleson, M. Gunzburger, Bridging methods for atomistic-to-continuum coupling and their implementation, Commun. Comput. Phys. 7 (4) (2010) 831–876. [50] T.J.R. Hughes, The finite element method: linear static and dynamic finite element analysis, Prentice-Hall, Englewood Cliffs, N.J., 1987.