**7. Bridging- Scale Modeling**

Where *u***˜** is displacement field. Using this equation, we could find displacement fields. By the use of displacement field we are able to find the forces exerted on the discrete dislocation as

*c*

*i* ¶*<sup>E</sup>* = - ¶

At this point, an iterative procedure involving the discrete dislocation positions, FE positions

**•** simulating the passage of defect from the atomistic to continuum in three dimensions

**•** annular dislocations that reside in both atomistic and continuum regions at the same time

Bridging Domain is a dynamic multiple scale method that couples MD with continuum. In this method the boundary that overlays the MD region and surrounding continuum region

Within overlaying bridging domain, the Hamiltonian is defined as a linear combination of the

 a

To make compatibility between the molecular and the continuum regions, we involve

( ) 0 *N ud i i* æ ö = -= ç ÷ è ø *I J JJ I* <sup>å</sup>*<sup>J</sup>*

Where *gI* are the Lagrange Multipliers, *uiJ* are the FE nodal displacements, and *diJ* are the MD

*M ext int* && =-- *<sup>L</sup>*

<sup>M</sup> = - H+ (19)

*g X* (20)

*II I I I uF F F* (21)

(1 ) *<sup>C</sup> H H* a

where *α* acts to scale the contribution of each domain to the total Hamiltonian.

*<sup>d</sup>* (18)

*i*

*P*

and atomic positions is solved until all degrees of freedom are at equilibrium.

However, we have some challenge in this method as follows:

**•** extending this approach to dynamic problems

has varying size, termed the bridging domain.

Lagrange multipliers in the overlap region as

The coupled equations of motion will be as follows:

molecular and continuum Hamiltonians:

**6.5. Bridging Domain**

16 Graphene - New Trends and Developments

displacements.

The bridging scale method that we want to introduce is an incorporation of MD and continu‐ um. The total displacement as our solution is decomposed to fine and coarse- scale as follows:

$$
\mu\left(\mathbf{x}\right) = \overline{\mu}\left(\mathbf{x}\right) + \mu'\left(\mathbf{x}\right) \tag{23}
$$

where *u*¯ is a coarse- scale solution that will be state in terms of finite element shape functions and *u* ' is a fine- scale solution. These two solutions are orthogonal, that is, one projection onto another is zero.

If *x*=*x*α (we use Greek indices for atoms and Roman indices for coarse- scale nodes) was the initial position for typical atoms of interesting body we can rewrite (23) as:

$$\mathfrak{u}\left(\boldsymbol{\chi}\_{a}\right) = \widetilde{\mathfrak{u}}\left(\boldsymbol{\chi}\_{a}\right) + \mathfrak{u}\left(\boldsymbol{\chi}\_{a}\right) \tag{24}$$

u¯(Xα) is an average behavior of body since it's interpolation between atoms so we anticipate it is continues function. We can write:

$$\overline{\mathfrak{u}}\left(\mathbf{X}\_{a}\right) = \sum\_{I} \mathbf{N}\_{I}^{a} \mathbf{d}\_{I} \tag{25}$$

where, NI *<sup>α</sup>* =NI (X*α*) is a shape function of node I at initial atomic position *X*α and *d*I is the displacement degree of freedom at node I.

As we said, we want to use MD method in this approach, to do this we have to utilize *q*α. *q*α is a displacement solution we derive from MD simulation (*q*<sup>α</sup> has the same role as *u*(*X*α) in other words, we are simulating *u*(*X*α) by the use of *q*α. Undoubtedly, *q*<sup>α</sup> will have projection onto u¯ (coarse- scale solution) and *u* (fine- scale solution), and with respect to their orthogonality, each projection will not represent the other one, so we can easily find the fine- scale solution, by subtracting the projection of *q*α onto coarse- scale solution, from *q*α, that is,

$$
\hat{u}' = q - \text{NW} \tag{26}
$$

or with indices representation:

$$
\mu^{\dagger} \left( \mathbf{X}\_{a} \right) = \mathbf{q}\_{a} - \sum\_{\mathbf{l}} \mathbf{N}\_{\mathbf{l}}^{\alpha} \mathbf{W}\_{\mathbf{l}} \tag{27}
$$

∑ I NI *<sup>α</sup>*WI is a coarse- scale projection. In other words, to find coarse- scale projection, we have projected *q*α onto shape functions NI. However, there is one burning question: how can we find WI? The answer is laid in least squares weighted residual method. First, we define *J* function as follows:

$$\mathbf{J} = \sum\_{a} m\_{a} \left( \mathbf{u}\_{a} - \sum\_{I} \mathbf{N}\_{I}^{a} \mathbf{w}\_{I} \right)^{2} \tag{28}$$

where *m*α is the atomic mass. The choice of *J* (called the mass-weighted square of the fine scale) allows the coarse and fine- scale kinetic energies to be decoupled. Minimizing of *J* will result in WI. If we gather all masses in a matrix as follows:

$$\mathbf{M}\_A = \begin{bmatrix} m\_1 \mathbf{I}\_{3 \times 3} \\ & m\_2 \mathbf{I}\_{3 \times 3} \\ & & \ddots \\ & & & \ddots \\ & & & & \ddots \end{bmatrix} \tag{29}$$

*J* can be writhed as:

$$J = \left(q - N\boldsymbol{\omega}\right)^{\mathsf{T}} \boldsymbol{\mathcal{M}}\_A \left(q - N\boldsymbol{\omega}\right) \tag{30}$$

Solving for *w* to minimize this quantity gives

$$\mathbf{w} = \mathbf{M}^{-1} \mathbf{N}^{\top} \mathbf{M}\_A \mathbf{q} \tag{31}$$

where the coarse- scale mass matrix *M* is

$$\mathbf{M} = \mathbf{N}^{\mathrm{T}} \mathbf{M}\_{\mathrm{A}} \mathbf{N} \tag{32}$$

Substituting *W* in Equation 26 with Equation 31, we will have

A Review on Modeling, Synthesis, and Properties of Graphene http://dx.doi.org/10.5772/61564 19

$$
\dot{u}' = q - Pq \tag{33}
$$

where the matrix *P* has been defined as

$$P = N M^{-1} N^{\dagger} M\_A \tag{34}$$

So we can write the total displacement *u*α (with respect to Equation 23-25 and 33) as

$$
\mu = Nd + q - Pd\tag{35}
$$

The last term in the above expression is called the "bridging scale"; it is that part of the molecular dynamics calculation that must be subtracted from the total in order to achieve a complete separation of scales. Our sake about the "scales" is fine and coarse scales, and our sake about "complete separation" is the orthogonality of these scales.

Indeed, *P* was the operator for the projection onto the coarse scale; with respect to Equation 35, we could represent another more perfect operator for this projection as

$$
\underline{\mathbf{Q}} = \overline{I} - \mathbf{P} \tag{36}
$$

or

'

( ) '

a

I I

a

is a coarse- scale projection. In other words, to find coarse- scale projection, we have

a

2

*J uw* (28)

*M* (29)

*J q Nw M q Nw* (30)


*<sup>A</sup>* = *<sup>T</sup> M NMN* (32)

I

projected *q*α onto shape functions NI. However, there is one burning question: how can we find WI? The answer is laid in least squares weighted residual method. First, we define *J* function

where *m*α is the atomic mass. The choice of *J* (called the mass-weighted square of the fine scale) allows the coarse and fine- scale kinetic energies to be decoupled. Minimizing of *J* will result

> 2 33

*I*

´

O O

*m*

ë û

=- - ( ) *<sup>A</sup>* ( ) *<sup>T</sup>*

1

*A*

é ù ê ú

*u* X q NW

*m N <sup>a</sup>*

æ ö = - ç ÷ è ø å å *I I α I*

a

1 33

*I*

´

*m*

 

*A*

Substituting *W* in Equation 26 with Equation 31, we will have

=

in WI. If we gather all masses in a matrix as follows:

Solving for *w* to minimize this quantity gives

where the coarse- scale mass matrix *M* is

 a

or with indices representation:

18 Graphene - New Trends and Developments

∑ I NI *<sup>α</sup>*WI

as follows:

*J* can be writhed as:

*u q NW* = - (26)

= -å (27)

$$\underline{\mathbf{Q}} = \mathbf{I} - \mathbf{N} \mathbf{M}^{-1} \mathbf{N}^{\mathsf{T}} \mathbf{M}\_{\mathcal{A}} \tag{37}$$

Thus, we can rewrite 35 as

$$
\mu = Nd + \underline{Q}q\tag{38}
$$

#### **7.1. Equations of motion**

In this step, we will derive the coupled MD and FE equations of motions. First, we adopt the Lagrangian equation definition we submitted in MD discussion for multiscale condition as "multiscale Lagrangian."

*7.1.1. Multiscale Lagrangian*

$$L\left(\mathbf{u}, \dot{\mathbf{u}}\right) = K\left(\dot{\mathbf{u}}\right) - U\left(\mathbf{u}\right) - f\_{\text{ext}}^{\text{T}}\mathbf{u} \tag{39}$$

$$K\left(\dot{\boldsymbol{\mu}}\right) = \frac{1}{2}\dot{\boldsymbol{\mu}}^T \boldsymbol{M}\_A \dot{\boldsymbol{\mu}} \tag{40}$$

$$
\dot{\mu} = N\dot{d} + \mathbb{Q}\dot{q} \tag{41}
$$

Substituting *u*˙ in Equation 40 with Equation 41 we will have

$$K\left(\dot{\mu}\right) = \frac{1}{2}\dot{\mu}^T M\_{\dot{\mathcal{A}}} \dot{\mu} = \frac{1}{2}\dot{d}^T N^T M\_{\dot{\mathcal{A}}} N\dot{d} + \dot{d}^T N^T M\_{\dot{\mathcal{A}}} \mathcal{Q}\dot{q} + \frac{1}{2}\dot{q}^T \mathcal{Q}^T M\_{\dot{\mathcal{A}}} \mathcal{Q}\dot{q} = \frac{1}{2}\dot{d}^T M \dot{d} + \frac{1}{2}\dot{q}^T M \dot{q} \tag{42}$$

The second term in the right- hand side of abvoe equation is zero because

$$N^{\mathcal{T}} \mathcal{M}\_A \mathcal{Q} = N^{\mathcal{T}} \mathcal{M}\_A \left( I - N \mathcal{M}^{-1} N^{\mathcal{T}} \mathcal{M}\_A \right) = N^{\mathcal{T}} \mathcal{M}\_A - M \mathcal{M}^{-1} N^{\mathcal{T}} \mathcal{M}\_A = 0$$

The fine- scale mass matrix *M* in equation 42 defined as

$$M \equiv \underline{\mathbf{Q}}^T M\_A \underline{\mathbf{Q}} = M\_A \underline{\mathbf{Q}} = \underline{\mathbf{Q}}^T M\_A \tag{43}$$

This will be proven by the complete expression for Q (Equation 36) and coarse- scale mass matrix M (Equation 30).

Thus, we can write Lagrangian 36 as

$$L\left(d, \dot{d}; q, \dot{q}\right) = \frac{1}{2}\dot{d}^T M \dot{d} + \frac{1}{2}\dot{q}^T M \dot{q} - \mathcal{U}\left(Nd + \mathcal{Q}q\right) + f\_{\text{ext}}^T Nd + f\_{\text{ext}}^T \mathcal{Q}q \tag{44}$$

This is the desired equation for us because we have eliminated cross terms such as *d* and *q* in the kinetic energy. The bridging scale is responsible for this elimination. Using the bridging scale makes us capable of decomposing total kinetic energy into the sum of the kinetic energy in the coarse scale plus that in the fine scale.

#### *7.1.2. Multiscale equations of motion*

We can derive equations of motion by the use of Lagrangian equation:

$$
\frac{d}{dt} \left( \frac{\partial \mathcal{L}}{\partial \dot{\mathbf{d}}} \right) - \frac{\partial \mathcal{L}}{\partial \mathbf{d}} = 0 \tag{45}
$$

$$\frac{d}{dt}\left(\frac{\partial \mathbf{L}}{\partial \dot{\mathbf{q}}}\right) - \frac{\partial \mathbf{L}}{\partial \mathbf{q}} = \mathbf{0} \tag{46}$$

Substituting the Lagrangian Equation 44 into Equations 45 and 46 gives us

= +& *u Nd Qq* & & (41)

*AAA A* & && *T T T T* & & & && *T T* & & *T T<sup>M</sup>* (42)

*M* º == *T T QMQ MQ QM AA A* (43)

*d d* (45)

(46)

Substituting *u*˙ in Equation 40 with Equation 41 we will have

The fine- scale mass matrix *M* in equation 42 defined as

,;, 2 2

in the coarse scale plus that in the fine scale.

*7.1.2. Multiscale equations of motion*

matrix M (Equation 30).

20 Graphene - New Trends and Developments

Thus, we can write Lagrangian 36 as

( ) 1 1 1 11 2 2 2 22 *K u u M u d N M Nd d N M Qq q Q M Qq d Md q q* & && = = + + =+ *<sup>T</sup>*

( ) 1 1 0 - - = - =- = *T T TT T N M Q N M I NM N M N M MM N M A A AA A*

This will be proven by the complete expression for Q (Equation 36) and coarse- scale mass

This is the desired equation for us because we have eliminated cross terms such as *d* and *q* in the kinetic energy. The bridging scale is responsible for this elimination. Using the bridging scale makes us capable of decomposing total kinetic energy into the sum of the kinetic energy

<sup>0</sup> *<sup>d</sup>*

<sup>0</sup> *<sup>d</sup>*

æ ö ¶ ¶ ç ÷ - = è ø ¶ ¶ & *L L q q*

æ ö ¶ ¶ ç ÷ - = è ø ¶ ¶ & *L L*

*M U ext ext* = + - ++ + & && & && *T T T T L d d q q d Md q q Nd Qq f Nd f Qq* (44)

The second term in the right- hand side of abvoe equation is zero because

( ) ( ) 1 1

We can derive equations of motion by the use of Lagrangian equation:

*dt*

*dt*

$$\mathcal{M}\ddot{\mathbf{d}} = -\frac{\partial \mathcal{U}(\mathbf{d}, \mathbf{q})}{\partial \mathbf{d}}\tag{47}$$

$$M\ddot{q} = -\frac{\partial \mathcal{U}(d, q)}{\partial q} \tag{48}$$

For simplicity, we have ignored the external force. Using the chain rule, we can rewrite Equations 47 and 48 as

$$\mathbf{M}\ddot{\mathbf{d}} = \frac{\partial \mathbf{U}}{\partial \mathbf{u}} \frac{\partial \mathbf{u}}{\partial \mathbf{d}}\tag{49}$$

$$
\hbar M \,\dot{\mathfrak{q}} = -\frac{\partial \mathcal{U} \,\mathcal{U}}{\partial \mathfrak{u}} \frac{\partial \mathfrak{u}}{\partial \mathfrak{q}}\tag{50}
$$

As we see, the derivation of potential energy *U* couples scales in both equations. We know that the variation of energy to displacement is the same type as force. These means we can conclude that the derivation of potential energy is the total forces on atoms:

$$f = -\frac{\partial \mathcal{U}(\mathfrak{u})}{\partial \mathfrak{u}}\tag{51}$$

Using Equations 51 and 38, we can rewrite Equations 49 and 50 as

$$\mathbf{M}\ddot{\mathbf{d}} = -\left(\frac{\partial \mathbf{u}}{\partial \mathbf{d}}\right)^{T}\frac{\partial \mathcal{U}}{\partial \mathbf{u}} = \mathbf{N}^{T}\mathbf{f} \tag{52}$$

$$\mathcal{M}\ddot{\vec{q}} = -\left(\frac{\partial \mathfrak{u}}{\partial \mathfrak{q}}\right)^{\mathsf{T}} \frac{\partial \mathcal{U}}{\partial \mathfrak{u}} = \mathbf{Q}^{\mathsf{T}} \mathcal{F} \tag{53}$$

where *M* =*Q*<sup>T</sup>*M*A. If the external forces in Equation 44 are kept in formulation, we will have

$$\mathbf{M}\ddot{\mathbf{d}} = \mathbf{N}^{\mathrm{T}} \left( \mathbf{f} + \mathbf{f}\_{\mathrm{ext}} \right) \tag{54}$$

$$M\ddot{\boldsymbol{\eta}} = \mathbf{Q}^{\mathrm{T}} \left( \boldsymbol{f} + \boldsymbol{f}\_{\mathrm{ext}} \right) \tag{55}$$

Equation 55 can be rewritten as

$$\mathbf{Q}^{\sf T} \mathbf{M}\_A \ddot{q} = \mathbf{Q}^{\sf T} \left( f + f\_{\rm ext} \right) \tag{56}$$

*Q* is a singular since multiplying a field by *Q* gives the fine- scale part of the field, and many different fields can have same fine- scale part. Therefore the solution of and equation like *Qu* =*u* ' for *u* is nonunique, and *Q* must be singular. It is concluded that Equation 56 does not have a unique solution, and we are free to choose any *q* that has the proper fine- scale part and thus satisfies Equation 56. One *q* that comply these conditions is that which exactly satisfies the following equation:

$$\mathbf{M}\_A \ddot{\mathbf{q}} = \mathbf{f} + \mathbf{f}\_{\text{ext}} \tag{57}$$

This is nothing but a molecular dynamics simulation for the atomic displacements, with the atomic forces given by Equation 51. So we are able to solve this equation using MD simulation.

Anyway, now we have the equation of motions as

$$
\mathbf{M}\_A \ddot{\mathbf{q}} = \mathbf{f} + \mathbf{f}\_{\text{ext}} \tag{58}
$$

$$\mathbf{M}\ddot{\mathbf{d}} = \mathbf{N}^{\mathrm{T}} \left( \mathbf{f} + \mathbf{f}\_{\mathrm{ext}} \right) \tag{59}$$

There are some relevant comments:


$$Nd = Pq\tag{60}$$

It means that the coarse scale is a projection of *q*, and instead of solving the FE equation, we can use *q*. However, as we will discuss, since we will eliminate the fine scale from the coarse scale, *q* will be limited to fine scale (not the entire domain), and thus it is not possible to calculate the coarse- scale solution everywhere via direct the projection of the MD displacements.


### **7.2. Removing fine- scale degree of freedom in coarse- scale region**

*Mq*&& = + ( *ext* ) *<sup>T</sup> Qf f* (55)

*<sup>A</sup>q*&& = + ( *ext* ) *T T QM Q f f* (56)

*<sup>A</sup> ext M ff q*&& = + (57)

*<sup>A</sup> ext M ff q*&& = + (58)

= + ( *ext* ) && *<sup>T</sup> Md N f f* (59)

*Q* is a singular since multiplying a field by *Q* gives the fine- scale part of the field, and many different fields can have same fine- scale part. Therefore the solution of and equation like *Qu* =*u* ' for *u* is nonunique, and *Q* must be singular. It is concluded that Equation 56 does not have a unique solution, and we are free to choose any *q* that has the proper fine- scale part and thus satisfies Equation 56. One *q* that comply these conditions is that which exactly satisfies the

This is nothing but a molecular dynamics simulation for the atomic displacements, with the atomic forces given by Equation 51. So we are able to solve this equation using MD simulation.

**•** Equation 58 is a fine- scale equation of motion and a standard MD solver can be used to

**•** Equation 59 is a coarse- scale equation of motion. This is an FE momentum equation, and to solve it, we should use finite elements method. It should be noted that for the coarse scale,

**•** As we said before, *q* (MD solution) simulate *u* (real total displacement). Although we proved Equation 58 for *q*, it is obvious that *u* as real total displacement have to satisfy Equation 58 in each atomic position (*X*α), this is found from Newton's second law. *u* is not necessarily equal to *q*; they just have the same behavior like solutions of a differential equations. Indeed, initial conditions determine their equality. Thus, if these two quantities have the same initial

conditions, then (which they should) then they are equal forever.

Equation 55 can be rewritten as

22 Graphene - New Trends and Developments

following equation:

Anyway, now we have the equation of motions as

There are some relevant comments:

obtain the MD displacements *q*.

**•** Equality of *q* and *u* and 38 gives

the summation sign turns into integral sign.

In the previous section, we saw the fine- scale equation of motion (Equation 58) and coarsescale equation of motion (Equation 59) in areas among this two- scale overlay together, and this is in contrast with the orthogonality of these two scale, so we should eliminate the effect of fine scale from coarse scale. However why do we eliminate fine scale from coarse scale and not vice versa? The answer is because eliminating one scale from the other will enter another terms in the eliminated scale's equation of motion. Since we want to avoid wave reflection back into the MD domain and avoid heat retention within this region, we should add some terms to the fine- scale equation of motion, and it is the best opportuni‐ ty for us to achieve this gold by eliminating the fine- scale effect from the coarse- scale domain. Thus, we limit our studies to the MD region and its equation of motion (Equa‐ tion 58) with assumption *f*ext=0:

$$\mathcal{M}\_A \ddot{q} = f\left(q\right) \tag{61}$$

To find the area that shall be removed, we should consider the behavior of these two scales. In fine scale, we have a harmonic or nonlinear potential in atomic interactions. However, at coarse scale, we have a harmonic or linear potential energy in atomic interaction. Thus, we can conclude that there is a transition area among these two scales that belong in the MD region. This is the area that shall be eliminated, and instead of this elimination, we should add some boundary conditions for the MD region. Thus, we divide the MD region in two sections: linear and nonlinear. This division is called the "linearized MD equation of motion." Two find MD boundary condition, since we are transient from fine to coarse scale in the transition area, we use FE method and then Lattice mechanics. This means we are treating with transition area, like a coarse scale, and maybe we could say we are dividing the MD area into fine and coarse sections. We first study a chain of atoms and then consider the general state.

Consider a one-dimensional chain of atoms according to the following figure:

From the lattice mechanics, we have

$$f\_u = -\frac{\partial \mathcal{U}(u)}{\partial u\_u}|\_{u=0} \tag{62}$$

$$\left.K\_{u \to u^{\cdot}} = \frac{\partial f\_{u}}{\partial u\_{n^{\cdot}}}\right|\_{u \to 0} = -\frac{\partial^{2} \mathcal{U}\{\boldsymbol{u}\}}{\partial \boldsymbol{u}\_{n} \partial \boldsymbol{u}\_{n^{\cdot}}}\big|\_{u \to 0} \tag{63}$$

where *f <sup>n</sup>* denotes the total force on atom (unit cell) *n* and *Kn*−*n*' is the springy factor. Using the Lagrange equation results in the equation of motion as follows:

$$\min \ddot{u}\_n\left(t\right) - \sum\_{n^+=n-\nu}^{n+\nu} K\_{n-n^+} u\_n\left(t\right) = f\_n^{ext}\left(t\right) \tag{64}$$

This equation is validated only for harmonic lattice, and *f <sup>n</sup> ext*(*t*) is the total external force acting on the unit cell *n*.

By dividing the MD region into fine and coarse scales, we have

$$
\overline{q} = \overline{\mu} + \mu^\* \tag{65}
$$

where *u*¯ denotes the coarse scale and *u* ′ denotes the fine scale of the MD region.

For the left side of Equation 61, we can write

$$M\_A \ddot{q} = M\_A \ddot{\overline{\mu}} + (M\_A \dddot{\mu}') = f(\overline{\mu}) + \left(Ku' + f'^{ext}\right) \tag{66}$$

where

$$\mathbf{K} = \frac{\partial \mathbf{f}}{\partial \mathbf{u}} \Big|\_{\mathbf{u}' = 0} \tag{67}$$

Moreover, *f ext* is the total external force acting on the harmonic section of the MD by an unharmonic section of the MD (note that *<sup>f</sup> ext* is different from *<sup>f</sup> ext* we neglected in Equation 48; *<sup>f</sup> ext* takes place within the MD section while *<sup>f</sup> ext* lies outside of this section).

A comparison between Equations 66 and 61 gives

From the lattice mechanics, we have

24 Graphene - New Trends and Developments

( ) <sup>0</sup> | *n u n U u*

'0 0 ' ' | | *<sup>n</sup> n n u u n n n f U u*

where *f <sup>n</sup>* denotes the total force on atom (unit cell) *n* and *Kn*−*n*' is the springy factor. Using the

( ) ' ' () ()

*n n nn n*


u

*mu t K u t f t* u

*u uu* -= =

*u* <sup>=</sup>

( ) <sup>2</sup>

*ext*

&& - = å (64)

*quu* = + ' (65)

denotes the fine scale of the MD region.

*AA A <sup>q</sup>* = + = ++ *<sup>M</sup> <sup>u</sup> <sup>M</sup> uM <sup>f</sup> <sup>u</sup> Ku <sup>f</sup>* && && && (66)

( ) ( )' ( ) ' *ext*

<sup>0</sup> | <sup>=</sup> ¶ <sup>=</sup> ¶ *' <sup>u</sup> <sup>f</sup> <sup>K</sup> u*

*ext*(*t*) is the total external force acting

(67)

¶ = - ¶ (62)

¶ ¶ = = - ¶ ¶¶ (63)

*f*

*K*

Lagrange equation results in the equation of motion as follows:

This equation is validated only for harmonic lattice, and *f <sup>n</sup>*

where *u*¯ denotes the coarse scale and *u* ′

For the left side of Equation 61, we can write

By dividing the MD region into fine and coarse scales, we have

on the unit cell *n*.

where

'

*n n*

= -

*n*

+

$$f\left(q\right) = f\left(\overline{u}\right) + \left(Ku^\prime + f^{ext}\right) \tag{68}$$

In Equation 66, we have grouped terms based on fine scale and coarse scale. Because of the orthogonality of coarse and fine scales and that the timescale for the coarse scale is much larger than that of atomic vibrations present in the fine scale, it is reasonable to equate the corre‐ sponding parts from each side of equality sign, that is,

$$M\_A \ddot{\overline{\mu}} = f(\overline{\mu}) \tag{69}$$

$$(M\_A \ddot{\mu}') = \left(K\mu' + f^{ext}\right) \tag{70}$$

Equation 70 is the one we were looking for; in other words, by removing the fine- scale degree Of freedom in the coarse- scale region, Equation 61 turns into Equation 70. Thus we use Equation 70 instead of Equation 61. However we still have not terminated; we should calculate *f ext* as an issue problem. To do this, we use lattice mechanics.

In lattice mechanics, Equation 70 turns into Equation 71:

$$\dot{\boldsymbol{u}}\_{n}^{\cdot}\left(\boldsymbol{t}\right) = \sum\_{\boldsymbol{n}^{\cdot} = \boldsymbol{n} - 1}^{\boldsymbol{n}} \boldsymbol{m}\_{A}^{-1} \mathbb{K}\_{\boldsymbol{n} = \boldsymbol{n}} \dot{\boldsymbol{u}}\_{n^{\cdot}}\left(\boldsymbol{t}\right) + \boldsymbol{m}\_{A}^{-1} \boldsymbol{f}\_{n}^{ext}\left(\boldsymbol{t}\right) \tag{71}$$

Equation 70 is the one we were looking for; in other words, by removing the fine- scale degree Of freedom in the coarse- scale region, Equation 61 turns into Equation 70. Thus we use Equation 70 instead of Equation 61. However we still have not terminated; we should calculate *f ext* as an issue problem. To do this, we use lattice mechanics.

In lattice mechanics, Equation 70 turns into Equation 71:

$$
\vec{\mu}\_n^{\cdot} \left( t \right) = \sum\_{n^-=n^-1}^n m\_A^{-1} \mathcal{K}\_{n-n^-} \vec{\mu\_{n^-}} \left( t \right) + m\_A^{-1} f\_n^{ext} \left( t \right) \tag{71}
$$

where *u*n is the displacement of unit cell *n*.

Again, consider the one- dimensional chain of atoms like the following figure:

*n*=0 is the position that elimination has took place, and the result of that is *f <sup>n</sup> ext*. We can write the following:

$$f\_{n}^{\text{ext}}\left(t\right) = \mathcal{S}\_{u,0} f\_{0}^{\text{ext}}\left(t\right) \tag{72}$$

Applying discrete Fourier transform (DFT) to Equation 71, we have

$$
\hat{\tilde{\mu}}^{\prime}(p,t) = m\_A^{-1} \hat{\mathcal{K}}(p) \hat{\mu}^{\prime}(p,t) + m\_A^{-1} f\_0^{\text{ext}}(t) \tag{73}
$$

Using Laplace transform (LT) for Equation 73 and solving for resulting displacements *U* ^ gives

$$
\widehat{\mathbf{U}}\left(p,s\right) = \hat{\mathbf{G}}\left(p,s\right) \left(m\_A^{-1} \mathbf{F}\_0^{ext}\left(s\right) + s\hat{\mathbf{u}}\left(0,p\right) + \hat{\mathbf{u}}\left(0,p\right)\right) \tag{74}
$$

where

$$
\hat{\mathbf{G}}\left(\boldsymbol{p},\mathbf{s}\right) = \left(\mathbf{s}^2 - \hat{A}(\boldsymbol{p})\right)^{-1} \tag{75}
$$

and

$$
\hat{A}(p) = m\_A^{-1} \hat{\mathcal{K}}(p) \tag{76}
$$

For *f* <sup>0</sup> *ext* we can write:

$$\left(f\_0^{ext}\left(t\right) = \mathcal{K}\_{-1}\mu'\_1(t)\right) \tag{77}$$

Taking LT from Equation 77 gives impedance force:

A Review on Modeling, Synthesis, and Properties of Graphene http://dx.doi.org/10.5772/61564 27

$$\mathbf{f}\_0^{\text{imp}}\left(\mathbf{s}\right) = \mathbf{K}\_{-1}\mathbf{U}\_{-1}^{\prime}\left(\mathbf{s}\right) \tag{78}$$

So to find *f* <sup>0</sup> *ext*, we should calculate *U*<sup>1</sup> ´ (*s*), This requires to omit *F*<sup>0</sup> *ext*(*s*) from Equation 74, so we first take the inverse discrete Fourier transform of Equation 61 in terms of *p*:

$$\mathbf{U}\_{\boldsymbol{n}}\left(\mathbf{s}\right) = \mathbf{G}\_{\boldsymbol{n}}\left(\mathbf{s}\right)m\_{A}^{-1}F\_{0}^{ext}\left(\mathbf{s}\right) + R\_{\boldsymbol{n}}^{d}\left(\mathbf{s}\right) \tag{79}$$

where

*ext*. We can write

^ gives

(72)

*A A u pt m K pu pt m f t* - - && = + (73)

*<sup>A</sup> p s p s m s su p u p* - *UG F* = ++ & (74)

*G* = - (75)


<sup>0</sup> ( ) 1 1'() *ext f t Ku t* - = (77)

Again, consider the one- dimensional chain of atoms like the following figure:

*n*=0 is the position that elimination has took place, and the result of that is *f <sup>n</sup>*

Applying discrete Fourier transform (DFT) to Equation 71, we have

() () ,0 0 *ext ext n n ft ft* = d

( ) ( ) ( ) () 1 1 <sup>0</sup> <sup>ˆ</sup> <sup>ˆ</sup> <sup>ˆ</sup> ' , ' , *ext*

Using Laplace transform (LT) for Equation 73 and solving for resulting displacements *U*

µ( ) ( ) () ( ) ( ) ' 1 <sup>0</sup> <sup>ˆ</sup> , , ˆ' 0, ' <sup>ˆ</sup> (0, ) *ext*

( ) ( ) <sup>1</sup> <sup>2</sup> <sup>ˆ</sup> , () <sup>ˆ</sup> *ps s Ap* -

() () <sup>1</sup> ˆ ˆ *Ap m Kp <sup>A</sup>*

the following:

26 Graphene - New Trends and Developments

where

and

For *f* <sup>0</sup>

*ext* we can write:

Taking LT from Equation 77 gives impedance force:

( ) ' ' ' ' 1 1 2 2 ' ' ' ' / 2 / 2 ( ) (0) ( ) (0) *N N d n nn n n n n n N n N s s G su G su* - - - - =- =- *<sup>R</sup>* = + å å & (80)

and *N* is the total number of unit cells in the domain. Now we write Equation 66 for *n*=0,1:

$$\mathbf{U}\_0^\dagger(\mathbf{s}) = \mathbf{G}\_0(\mathbf{s}) m\_A^{-1} \mathbf{F}\_0^{\text{ext}}(\mathbf{s}) + \mathbf{R}\_0^d(\mathbf{s}) \tag{81}$$

$$\mathbf{M}\_1^\cdot \left(\mathbf{s}\right) = \mathbf{G}\_1(\mathbf{s}) m\_A^{-1} \mathbf{F}\_0^{ext} \left(\mathbf{s}\right) + \mathbf{R}\_1^d(\mathbf{s}) \tag{82}$$

Solving these two equations simultaneously can omit *F*<sup>0</sup> *ext*(*s*) and gives

$$\mathbf{U}\_1^\dagger(\mathbf{s}) = \mathbf{G}\_1(\mathbf{s})\mathbf{G}\_0^{-1}(\mathbf{s}) \left( \mathbf{U}\_0^\dagger(\mathbf{s}) - \mathbf{R}\_0^d(\mathbf{s}) \right) + \mathbf{R}\_1^d(\mathbf{s}) \tag{83}$$

We can rewrite Equation 83 as

$$\mathbf{U}\_1^\cdot(\mathbf{s}) = \mathbf{Q}(\mathbf{s}) \left( \mathbf{U}\_0^\cdot(\mathbf{s}) - \mathbf{R}\_0^d(\mathbf{s}) \right) + \mathbf{R}\_1^d(\mathbf{s}) \tag{84}$$

where

$$\mathbf{Q}\begin{pmatrix}\mathbf{s}\end{pmatrix} = \mathbf{G}\_1\begin{pmatrix}\mathbf{s}\end{pmatrix}\mathbf{G}\_0^{-1}\begin{pmatrix}\mathbf{s}\end{pmatrix} \tag{85}$$

Now taking the inverse Laplace transform (ILT) of Equations 78 and 84 and convoluting the displacement, we obtain

$$f\_0^{imp}\left(t\right) = \underset{0}{\dot{\boldsymbol{\Theta}}} \theta(t-\tau) \left(\dot{\boldsymbol{\mu}}\_0\left(\tau\right) - \boldsymbol{\mathcal{R}}\_0^d(\tau)\right) d\tau + \boldsymbol{\mathcal{R}}\_0^f(t) \tag{86}$$

We know

$$
\overline{u}\_0'(t) = q\_0(t) - \overline{u}\_0(t) \tag{87}
$$

And then

$$f\_0^{imp}\left(t\right) = \int\_0^t \Theta(t-\tau) \left(q\_0\left(\tau\right) - \overline{u}\_0\left(\tau\right) - \mathbf{R}\_0^d(\tau)\right) d\tau + \mathbf{R}\_0^f(t) \tag{88}$$

Where the time history kernel *θ(t)* and random force *R*<sup>0</sup> *<sup>f</sup>* (*t*) are defined to be

$$\theta(t) = \mathcal{L}^{-1}\left(\mathcal{K}\_{-1}\mathcal{Q}(s)\right) \tag{89}$$

$$\mathbf{R}\_0^\prime \left( t \right) = \mathbf{K}\_{-1} \mathbf{R}\_1^\prime \left( t \right) \tag{90}$$

Other terms in Equation 71 is equal with *f* (*u* ′ ), and we can rewrite Equation 71 as

$$M\_A \ddot{\mu}' = f\left(\mu'\right) + f\_0^{imp}\left(t\right) \tag{91}$$

We can turn *u* ′ to *q* as

$$\mathcal{M}\_{\mathcal{A}}\ddot{q} = f\left(q\right) + f\_0^{imp}\left(t\right) \tag{92}$$

Thus, the final coupled forms of the equations of motion are

$$\mathbf{M}\ddot{\mathbf{d}} = \mathbf{N}^T f(\mathbf{u})\tag{93}$$

$$\mathbf{M}\_{\mathcal{A}}\ddot{\mathbf{q}} = \mathbf{f}\left(\mathbf{q}\right) + f\_0^{imp}\left(\mathbf{t}\right) \tag{94}$$

A Review on Modeling, Synthesis, and Properties of Graphene http://dx.doi.org/10.5772/61564 29

$$f\_0^{imp}\left(t\right) = \bigwedge\_{0}^{t} \theta(t-\tau) \left(q\_0\left(\tau\right) - \overline{u}\_0\left(\tau\right) - \mathbf{R}\_0^d(\tau)\right)d\tau + \mathbf{R}\_0^f(t) \tag{95}$$

where the impedance force *f* <sup>0</sup> *imp* acts only on the boundary atom 0.

#### *7.2.1. Elimination of fine- scale degrees of freedom: 3D generalization*

( ) ( ( ) ) ' 0 00 0

() () '

<sup>0</sup> ( ) ( 00 0 0 ( ) ( ) )

 t

*imp <sup>d</sup> <sup>f</sup> ft t q u d t* =- - - +

() ( ) <sup>1</sup> 1

( ) () <sup>0</sup>

( ) () <sup>0</sup>

*tL s*( ) -

( ) ( ) ()

 tt

ò *R R* (88)

*<sup>f</sup>* (*t*) are defined to be


0 11 ( ) ( ) *<sup>f</sup> <sup>d</sup> t t* - *R KR* = (90)

), and we can rewrite Equation 71 as

( ) () <sup>0</sup> '' *imp MA*&& = + *f u f tu* (91)

*imp q*&& = + *f t M fq <sup>A</sup>* (92)

( ) *<sup>T</sup>* = *Md N f u* && (93)

*imp q*&& = + *f t M fq <sup>A</sup>* (94)

 t

*imp <sup>d</sup> <sup>f</sup> ft t u d t* =- - +

 t

( ) ( ) ()

 tt

ò *R R* (86)

<sup>000</sup> *u t q t ut* = - ( ) (87)

0

0

Where the time history kernel *θ(t)* and random force *R*<sup>0</sup>

Other terms in Equation 71 is equal with *f* (*u* ′

Thus, the final coupled forms of the equations of motion are

to *q* as

q

*t*

q t

We know

28 Graphene - New Trends and Developments

And then

We can turn *u* ′

*t*

q t

> In 3D state, each unit cell can be labeled with three indices, *l, m*, and *n* indicating the position along axes in the direction of the three primitive vectors of the crystal structure. The equation of motion can be rewritten as

$$\dot{\boldsymbol{m}}\_{l,m,n}^{\cdot}\left(\boldsymbol{t}\right) = \sum\_{l-l-1:n}^{l} \sum\_{-m-\mu n-\mu n-\nu}^{m\*\mu} \sum\_{-n-\nu}^{n\*\nu} \mathbf{M}\_{\boldsymbol{A}}^{-1} \mathbf{K}\_{\boldsymbol{l}-l,m-m;\,\boldsymbol{n}-n} \dot{\boldsymbol{\mu}}\_{l,m;\,\boldsymbol{n}}^{\cdot}\left(\boldsymbol{t}\right) + \boldsymbol{M}\_{\boldsymbol{A}}^{-1} \boldsymbol{f}\_{l,m;\boldsymbol{n}}^{\text{ext}}\left(\boldsymbol{t}\right) \tag{96}$$

where *μ* and *υ* represent the range of the forces in the *m* and *n* coordinate directions. We consider the one dimensional boundary condition in our study. Our boundary is located in *l*=0, like the following figure:

Like the one-dimensional chain, we can write

$$f\_{l,m,n}^{ext}\left(t\right) = \delta\_{l,0} f\_{0,m,n}^{ext}\left(t\right) \tag{97}$$

By taking LT and DFT of Equation 96 and solving for *U* ^ ', we obtain

$$
\hat{\mathbf{U}}^{\prime}(p,q,r,\mathbf{s}) = \hat{\mathbf{G}}(p,q,r,\mathbf{s}) \left( \mathbf{M}\_A^{-1} \hat{\mathbf{F}}\_0^{\text{ext}} \left( p,q\_\prime r,\mathbf{s} \right) + \mathbf{s} \hat{\mathbf{u}}^{\prime} \left( p,q\_\prime r,0 \right) + \dot{\mathbf{u}}^{\prime}(p,q,r,0) \right) \tag{98}
$$

where

$$
\ddot{A}\left(p,q,r\right) = \mathbf{M}\_A^{-1}\mathbf{K}\_{p,q,r} \tag{99}
$$

and

$$\hat{\mathbf{G}}\left(p,q,r,s\right) = \left(\mathbf{s}^2 I - \hat{\mathbf{A}}(p,q,r)\right)^{-1} \tag{100}$$

For *<sup>f</sup>* 0,*m*,*<sup>n</sup> ext* , we can write

$$f\_{0,m,n}^{ext}\left(t\right) = \sum\_{m'-m-\mu\stackrel{\circ}{u}-n\stackrel{\circ}{u}-n\stackrel{\circ}{v}}^{m'-m+\mu\stackrel{\circ}{u}-n\stackrel{\circ}{v}} K\_{-1,m-m',n-n'}\mu\stackrel{\circ}{1}\_{1,m',n'}(t)\tag{101}$$

Taking LT from Equation 101 gives impedance force:

$$\mathbf{M}\_{w,u}^{imp}\left(\mathbf{s}\right) = \sum\_{i^{-}=m-\mu^{i}\mu^{i}\mathbf{i}=n+\nu}^{m^{i}-m^{i}\mu^{i}}\mathbf{K}\_{-1,m-m^{i},n-n^{i}}\mathbf{U}\_{1,m^{i},n^{i}}^{\cdot}\left(\mathbf{s}\right)\tag{102}$$

Thus, to find *<sup>f</sup>* 0,*m*,*<sup>n</sup> ext* we should calculate *<sup>U</sup>* ' 1,*m* ' ,*n* ' (*s*). This requires to omit *F* ^ 0 *ext* from Equation 74, so we first take inverse discrete Fourier transform of (98) in terms of p:

$$
\tilde{\mathbf{U}}\_{l}^{\prime}(q\_{\prime}r\_{\prime}s) = \mathbf{M}\_{A}^{-1}\tilde{\mathbf{G}}\_{l}(q\_{\prime}r\_{\prime}s)\hat{\mathbf{F}}\_{0}^{\mathrm{ext}}\left(q\_{\prime}r\_{\prime}s\right) + \tilde{\mathbf{R}}\_{l}^{d}(q\_{\prime}r\_{\prime}s) \tag{103}
$$

where

$$
\tilde{\mathbf{R}}\_{l}^{d}\left(q\_{l},r,s\right) = \sum\_{l'=-L/2}^{\frac{L}{2}-1} \tilde{\mathbf{G}}\_{l-l'}\left(q\_{l'}r,s\right)\left(s\mathbf{U}\_{l'}\left(q\_{l'}r,0\right) + \dot{\mathbf{U}}\_{l'}\left(q\_{l'}r,0\right)\right) \tag{104}
$$

Here *L* is the total number of unit cells in the x-direction. Letting *l*=0,1 and eliminating *F* ^ 0 *ext* we have:

$$
\tilde{\mathbf{U}}'\_1(q,r,s) = \tilde{\mathbf{Q}}(q,r,s) \Big( \tilde{\mathbf{U}}\_0(q,r,s) - \tilde{\mathbf{R}}\_0^d(q,r,s) \Big) + \tilde{\mathbf{R}}\_1^d(q,r,s) \tag{105}
$$

where

where

and

where

have:

For *<sup>f</sup>* 0,*m*,*<sup>n</sup> ext* , we can write

30 Graphene - New Trends and Developments

( ) <sup>1</sup> , , <sup>ˆ</sup> , , *A pqr M KA pqr*

( ) ( ) <sup>1</sup> <sup>2</sup> <sup>ˆ</sup> ,,, (,,) <sup>ˆ</sup> *pqrs s pqr* -

0, , 1, , ' 1, , '( )

, 1, , ' 1, , '( )

1,*m* ' ,*n* '

0 ˆ ,, ,, ,, (,,) *ext <sup>d</sup>*

( ) ' ( ( ) )

*qrs qrs s qr qr*

Here *L* is the total number of unit cells in the x-direction. Letting *l*=0,1 and eliminating *F*

*l l l l l*

( ) ( )( ) ( ) ( ) '


, , ( , , ) , ,0 ( , ,0)


*m n mmnn mn*

*fs K U s*

u

u

*m n mmnn mn*

*ft K ut* m

( )

*ext*

Taking LT from Equation 101 gives impedance force:

*imp*

Thus, to find *<sup>f</sup>* 0,*m*,*<sup>n</sup> ext* we should calculate *<sup>U</sup>* '

*d*

( )

' 1

'

*l L*

=-

1 2


*L*

/ 2

''

''

' '

*mm nn*

74, so we first take inverse discrete Fourier transform of (98) in terms of p:

( ) ( )( )

m u

= - =-

*m m n n*

m

= + = +

' '

*mm nn*

m u

= - =-

*m m n n*

= + = +


*G IA* = - (100)

<sup>=</sup> å å (101)

<sup>=</sup> å å (102)

(*s*). This requires to omit *F*

^ 0

*ext* from Equation

^ 0 *ext* we

' '

' '

'

*<sup>l</sup> A l <sup>l</sup> qrs qrs qrs qrs* - *<sup>U</sup>* = + % % *MG F R*% (103)

' ' ' '

<sup>1</sup> <sup>0</sup> <sup>10</sup> ,, ,, ,, ,, (,,) *d d U'* % %% *qrs qrs qrs qrs qrs* = -+ *Q RR <sup>U</sup>* % % (105)

*R G* % = + å % *U U*& (104)


'

$$
\tilde{\mathbf{Q}}\left(q,r,s\right) = \tilde{\mathbf{G}}\_1(q,r,s)\tilde{\mathbf{G}}\_0^{-1}(q,r,s) \tag{106}
$$

Taking IFT of 105 and using convolution property of the DFT, we gain

$$\mathbf{U}'\_{1,m,n}(\mathbf{s}) = \sum\_{m'=-M/2}^{\frac{M}{2}-1} \sum\_{n'=-N/2}^{\frac{N}{2}-1} \mathbf{Q}\_{m-m',n-n'}(\mathbf{s}) \left( \mathbf{U}'\_{0,m',n'}(\mathbf{s}) - \mathbf{R}\_{0,m',n'}^{d}(\mathbf{s}) \right) + \mathbf{R}\_{1,m,d}^{d}(\mathbf{s}) \tag{107}$$

where upon plane *l* we can write:

$$\mathbf{R}\_{l,m,n}^{d}\left(\mathbf{s}\right) = \sum\_{l=-L/2m}^{\underline{L}-1} \sum\_{-M/2n}^{\underline{M}-1} \sum\_{-N/2}^{\underline{N}-1} \mathbf{G}\_{l-l,m-m,n-n}\left(\mathbf{s}\right) \left(\mathbf{s}\mathbf{U}\_{\underline{l},m,n}^{\underline{i}}\left(\mathbf{0}\right) + \dot{\mathbf{U}}\_{\underline{l},m,n}^{\underline{i}}\left(\mathbf{0}\right)\right) \tag{108}$$

Taking ILT of Equation 108 gives

$$\mathbf{R}\_{l,m,n}^{d}\left(t\right) = \sum\_{l=-L/2}^{\frac{L}{2}-1} \sum\_{m=-M/2}^{\frac{M}{2}-1} \sum\_{i=-N/2}^{\frac{N}{2}-1} \left(\dot{\mathbf{g}}\_{l-l,m-m,n-n}(t)\dot{\boldsymbol{u}}\_{l,m,n}^{\cdot}(0) + \mathbf{g}\_{l-l,m-m,n-n}(t)\dot{\boldsymbol{u}}\_{l,m,n}^{\cdot}(0)\right) \tag{109}$$

where

$$\mathbf{g}\_{l,w,n}\left(t\right) = L^{-1}\left(\mathbf{G}\_{l,w,n}(s)\right) \tag{110}$$

$$\dot{\mathbf{g}}\_{l,w,u}\left(t\right) = L^{-1}\left(\mathbf{sG}\_{l,w,u}(\mathbf{s})\right) \tag{111}$$

Substituting Equation 107 into Equation 102 and applying ILT, the impedance force boundary condition obtains as

$$\mathbf{f}\_{w,n}^{\text{imp}}\left(\mathbf{t}\right) = \sum\_{\substack{m=-M/2\\m=-M/2}}^{\frac{M}{2}-1} \sum\_{n=-N/2}^{N} \left\{ \boldsymbol{\theta}\_{m-\dot{m},n\to\dot{n}}(t-\tau) \left(\boldsymbol{\mu}\_{0,\dot{m},\dot{n}}^{\cdot}\left(\tau\right) - \boldsymbol{\mathcal{R}}\_{0,\dot{m},\dot{n}}^{d}\left(\tau\right)\right) d\tau + \boldsymbol{\mathcal{R}}\_{0,\dot{m},n}^{f}(t) \right\} \tag{112}$$

where the time history kernel *θ*(*t*) is defined to be

$$\Theta\_{m,u}\left(t\right) = L^{-1}\left(\Theta\_{m,u}\left(s\right)\right) \tag{113}$$

$$\Theta\_{m,n}\left(s\right) = \sum\_{n'=m-\mu\nmid \, -n-\nu}^{m\ast\mu} \sum\_{-1, \, m-m', \, n-n'}^{n\ast\nu} \mathcal{Q}\_{m',n'}\left(s\right) \tag{114}$$

and

$$\mathbf{R}\_{0,m,n}^{\prime}\left(t\right) = \sum\_{l}^{\frac{L}{2}-1} \sum\_{m=-M/2}^{M} \sum\_{i=-N/2}^{\frac{N}{2}-1} \mathbf{K}\_{-1,m-m^{\prime},n-m} \mathbf{R}\_{l,m,n}^{d}\left(t\right) \tag{115}$$

where *R <sup>d</sup>* is given by Equation 109.

For simplicity, we rewrite Equation 112 as

$$\mathbf{f}\_{m,u}^{imp}(t) = \sum\_{i=-n\_c}^{n\_c} \sum\_{i=-n\_c}^{n\_c} \left[ \theta\_{m-u',n-i}(t-\tau) \left( \dot{\boldsymbol{u}}\_{0,u',n'}^{\cdot}(\tau) - \mathbf{R}\_{0,u',i}^{d}(\tau) \right) d\tau + \mathbf{R}\_{0,m,u}^{f}(t) \tag{116} \right]$$

where *nc* refers to a maximum number of atomic neighbors, which will be used to compute the impedance force. Certainly, all of this *nc* atoms are located in the MD region with linear or harmonic behavior.

As we evolved in Equation 92, we can write the final coupled form of the equations of motion as:

$$\mathbf{M}\_A \ddot{\mathbf{q}} = \mathbf{f}\left(\mathbf{q}\right) + f\_{0,m,u}^{imp} \tag{117}$$

$$\mathbf{M}\ddot{\mathbf{d}} = \mathbf{N}^{\mathsf{T}}f(\mathsf{u})\tag{118}$$

where

$$\mathbf{f}\_{m,n}^{imp}(t) = \sum\_{i=-n\_i}^{n\_c} \sum\_{u=-n\_c}^{n\_c} \left[ \theta\_{m-i, u=n}(t-\tau) \left( \mathbf{q}\_{0,n',u} \left( \tau \right) - \overline{\mathbf{u}}\_{0,n',u}(\tau) - \mathbf{R}\_{0,m',u}^d(\tau) \right) \right] d\tau + \mathbf{R}\_{0,m,n}^f(t) \tag{119}$$

Knowing that

$$\stackrel{\cdot}{\boldsymbol{\mu}}\_{0,m^{\cdot},n^{\cdot}}\left(\boldsymbol{\tau}\right) = \boldsymbol{q}\_{0,m^{\cdot},n^{\cdot}}\left(\boldsymbol{\tau}\right) - \overline{\boldsymbol{u}}\_{0,m^{\cdot},n^{\cdot}}\left(\boldsymbol{\tau}\right) \tag{120}$$

Note:

( ) ' ' ''

*m n mmnn mn mm nn*

( ) ' '

*t t*

'

where *nc* refers to a maximum number of atomic neighbors, which will be used to compute the impedance force. Certainly, all of this *nc* atoms are located in the MD region with linear or

As we evolved in Equation 92, we can write the final coupled form of the equations of motion

( ) 0, , *imp*

 t

= --+ å å ò *<sup>f</sup> uR R* (116)

*t d t* - -

, 0, , , 0, , ' 0, ,

*imp d f m n mmnn m n m n m n*

 t -- -

*<sup>R</sup>* <sup>=</sup> ååå *K R* (115)

( ) ( ) ()

 tt

*<sup>A</sup> m n M fq q*&& = + *f* (117)

*d N*= ( ) && *<sup>T</sup> M fu* (118)

 tt

= - (120)

( ) ( ) ( ) ()

 t

*t d t* - -

tt

<sup>111</sup> 22 2 0, , 1, , , , /2 /2 /2


*LM N f d m n mmnn lmn*

u

, 1, , , ( )

*s s*


*Θ KQ* <sup>=</sup> å å (114)

( )

' '

'' '

*l Lm Mn N*

( ) ' ' ( ' ( ) ' ' )

( ) ' ' ( ' ' ( ) ' ' ' ' )

 t

'

t

'

, 0, , , 0, , 0, , 0, ,

() () ' ' ' '

0, , ' 0, , 0, , ( ) *m n m n m n uqu*

 t

<sup>=</sup> å å - -- + ò *<sup>f</sup> q uR R* (119)

*imp d f m n mmnn m n m n m n m n*

=- =- =-

and

where *R <sup>d</sup>* is given by Equation 109.

32 Graphene - New Trends and Developments

*t*

' '

*t*

Knowing that

*m nn n*

=- =-

*c c*

*n n t*

*c c*

0

q

harmonic behavior.

as:

where

For simplicity, we rewrite Equation 112 as

' '

*m nn n*

=- =-

*c c*

*n n t*

*c c*

0

q

m u

= - =-

*m n*

+ +

m

