# Integral Solution of the Interface Profile of Grain Boundary Grooving by Surface Diffusion

Authors: Caifang Wang; Hui Xing; Fengping Yao; Changle Chen; Jianyuan Wang; Fei Fan
Journal: Advances in Materials Science and Engineering
 (2012)

---

## Abstract

In polycrystalline thin films, formation of grain boundary grooving is investigated widely due to its effects on their properties. Infinite series solution of the interface of grain-boundary grooving by surface diffusion has been obtained. In this paper, with Fourier transform, a novel integral solution of the profile of grain boundary grooving is derived. Numerical examples are taken to test the performance of the integral solution. According to the numerical examples, we find that the novel solution is in accordance with the infinite series solution.


---

## Body

## 1. Introduction

In polycrystalline thin films, grain boundary grooving through the thickness of the film is a common failure mode that strongly affects their properties. Formation and development of grain boundary grooving have received a great deal of attention over decades [1–6]. The grooving forms and develops at the point of intersection when the grain boundary ends at a free boundary in order to reduce the total free energy. In absence of the external potential field, the grooving profile migrates via surface diffusion. Figure 1 shows the AFM image of grain boundary grooving of BiFeO3/LaNiO3 thin film. Mass transport by surface diffusion is driven by the curvature Laplacian. For convex surface, mass flows from high curvature parts to low parts, while in reverse for concave surfaces. The theory of the capillary-driven shape evolution was first provided by Mullins [1], and this kind of problem hence is so-called Mullins problem. In his study, he considered that the grooving is symmetric between the two grains in the small slope approximation. Assume that the surface energy is isotropic and the surface flux can be expressed as
(1)Js=-DsγΩnskT∇sκ,
where ∇sκ is the surface curvature gradient, Ds is the coefficient of surface diffusion, s is the arc length of the interface profile, γ is the isotropic surface energy, Ω is the atomic volume, ns is the number of mobile atoms in the surface layer per unit, and kT is with its usual meaning. Therefore the surface normal speed can be written as
(2)Vn=DsγΩ2nskT∇s2κ.
Defining y=y(x,t) as the interface profile and using κ=-y′′(1+y')-2/3, we could get
(3)yt=-B∂x[(1+Y′2)1/2∂x(1+y′)-2/3],
where B=DsγΩ2ns/kT. In small slope approximation, it is reduced to
(4)yt=-By′′′′.Figure 1
AFM image of grain boundary grooving of BiFeO3/LaNiO3 thin film.Wong et al. [7] have shown that the film profiles at different times are found to be self-similar. The wedge contact angle will change at initial time from wedge angle α to the equilibrium β. Therefore, the partial differential equation can be converted into an ordinary differential equation before solving it by the integral method. In this paper, we mainly seek the integral solution of (4) with proper initial and boundary conditions. The rest of the paper is organized as follows. In Section 2, we derive the main integral solution of the interface profile with Fourier transform. Then we change the integral solution into normalized solution with some variable substitution. In Section 3, we take numerical experiments to test the novel integral solution. In Section 4, we give the conclusion of this paper.


## 2. Integral Solution

In the following, we consider the symmetric case of the profiley(x,t). Since the symmetry is about y axis, we only consider the profile y(x,t) for x>0. As in [1], we assume that the initial surface is flat. Also, we assume that at the groove root, the flux is zero and there is a fixed dihedral angle. Hence, the groove profile y(x,t) satisfies the following initial boundary problem:
(5)yt+By′′′′=0,x>0,y(x,0)=0,y′(0,t)=m,y′′′(0,t)=0,
where m is the slope related to the dihedral angle enclosed by the surface.


### 2.1. Main Solution

According to the theories of partial differential equations and Fourier transform, we can get the integral solution of the initial boundary problem (5) as follows:
(6)y(x,t)=mx-12π×∫0∞‍∫-∞∞‍e-Bξ4t×(cos⁡ξ(x-y)+cos⁡ξ(x+y))·mydξdy.In the following, we use three steps to derive the integral solution (6).Letz(x,t)=mx-y(x,t). Then we only need to seek the solution z(x,t) of the following problems:
(7)zt+Bz′′′′=0,x>0,(8)z(x,0)=mx,(9)z′(0,t)=0,(10)z′′′(0,t)=0.Step 1.
We first consider the solution of(11)wt+Bw′′′′=0,-∞<x<∞,w(x,0)=g(x).
By taking Fourier transform with respect to x, we can obtain the solution of (11)
(12)w(x,t)=12π∫-∞∞‍g(y)dy∫-∞∞‍e-Bξ4te-iξ(x-y)dξ=12π∫-∞∞‍g(y)dy∫-∞∞‍e-Bξ4tcos⁡ξ(x-y)dξ.
Since a Schwartz function's Fourier transform is still a Schwartz function, we deduce that
(13)∫-∞∞‍e-Bξ4te-iξ(x-y)dξ
is a Schwarz function for the variable y, which makes the double integral sense in (12) when g(x) is a polynomial function. In fact, we can refer to the paper [8] for further discussions of elementary theory of (11).Step 2.
It is easy to verify that(14)z(x,t)=12π×∫0∞‍∫-∞∞‍e-Bξ4t×(cos⁡ξ(x-y)+cos⁡ξ(x+y))·mydξdy
satisfies (7)–(10). In fact, it is trivial to prove that z(x,t) satisfies (7), (9), and (10). Moreover, we can verify that z(x,t) also satisfies (8) by choosing
(15)g(x)={mx,x>0,-mx,x<0.Step 3.
Finally, we can obtain the solution of (5)
(16)y(x,t)=mx-z(x,t)=mx-12π×∫0∞‍∫-∞∞‍e-Bξ4t×(cos⁡ξ(x-y)+cos⁡ξ(x+y))·mydξdy.


### 2.2. Normalized Solution

Letx-=x(Bt)-1/4, y-=y(Bt)-1/4 and τ=ξ(Bt)1/4. With the variable substitution, the function z(x,t) can also be expressed as
(17)z(x-,t)=m(Bt)1/42π×∫0∞‍∫-∞∞‍e-τ4×(cos⁡τ(x--y-)+cos⁡τ(x-+y-))·y-dτdy-.Let(18)h(u)=∫-∞∞‍e-τ4cos⁡τudτ.
It is the Fourier transform of a Schwartz function e-τ4. According to the properties of Fourier transform, h(u) is still Schwartz function and z(x-,t) is integrable. With simple calculation, we arrive at
(19)z(x-,t)=m(Bt)1/42π×[2x-∫0x-‍h(u)du-2∫0x-‍h(u)udu+2∫0∞‍h(u)udu].
Finally, we get the normalized solution of the initial boundary problem (5)
(20)y(x-,t)=m(Bt)1/4x--z(x-,t)=m(Bt)1/4π[πx--x-∫0x-‍h(u)du+∫0x-‍h(u)udu-∫0∞‍h(u)udu].


## 2.1. Main Solution

According to the theories of partial differential equations and Fourier transform, we can get the integral solution of the initial boundary problem (5) as follows:
(6)y(x,t)=mx-12π×∫0∞‍∫-∞∞‍e-Bξ4t×(cos⁡ξ(x-y)+cos⁡ξ(x+y))·mydξdy.In the following, we use three steps to derive the integral solution (6).Letz(x,t)=mx-y(x,t). Then we only need to seek the solution z(x,t) of the following problems:
(7)zt+Bz′′′′=0,x>0,(8)z(x,0)=mx,(9)z′(0,t)=0,(10)z′′′(0,t)=0.Step 1.
We first consider the solution of(11)wt+Bw′′′′=0,-∞<x<∞,w(x,0)=g(x).
By taking Fourier transform with respect to x, we can obtain the solution of (11)
(12)w(x,t)=12π∫-∞∞‍g(y)dy∫-∞∞‍e-Bξ4te-iξ(x-y)dξ=12π∫-∞∞‍g(y)dy∫-∞∞‍e-Bξ4tcos⁡ξ(x-y)dξ.
Since a Schwartz function's Fourier transform is still a Schwartz function, we deduce that
(13)∫-∞∞‍e-Bξ4te-iξ(x-y)dξ
is a Schwarz function for the variable y, which makes the double integral sense in (12) when g(x) is a polynomial function. In fact, we can refer to the paper [8] for further discussions of elementary theory of (11).Step 2.
It is easy to verify that(14)z(x,t)=12π×∫0∞‍∫-∞∞‍e-Bξ4t×(cos⁡ξ(x-y)+cos⁡ξ(x+y))·mydξdy
satisfies (7)–(10). In fact, it is trivial to prove that z(x,t) satisfies (7), (9), and (10). Moreover, we can verify that z(x,t) also satisfies (8) by choosing
(15)g(x)={mx,x>0,-mx,x<0.Step 3.
Finally, we can obtain the solution of (5)
(16)y(x,t)=mx-z(x,t)=mx-12π×∫0∞‍∫-∞∞‍e-Bξ4t×(cos⁡ξ(x-y)+cos⁡ξ(x+y))·mydξdy.


## 2.2. Normalized Solution

Letx-=x(Bt)-1/4, y-=y(Bt)-1/4 and τ=ξ(Bt)1/4. With the variable substitution, the function z(x,t) can also be expressed as
(17)z(x-,t)=m(Bt)1/42π×∫0∞‍∫-∞∞‍e-τ4×(cos⁡τ(x--y-)+cos⁡τ(x-+y-))·y-dτdy-.Let(18)h(u)=∫-∞∞‍e-τ4cos⁡τudτ.
It is the Fourier transform of a Schwartz function e-τ4. According to the properties of Fourier transform, h(u) is still Schwartz function and z(x-,t) is integrable. With simple calculation, we arrive at
(19)z(x-,t)=m(Bt)1/42π×[2x-∫0x-‍h(u)du-2∫0x-‍h(u)udu+2∫0∞‍h(u)udu].
Finally, we get the normalized solution of the initial boundary problem (5)
(20)y(x-,t)=m(Bt)1/4x--z(x-,t)=m(Bt)1/4π[πx--x-∫0x-‍h(u)du+∫0x-‍h(u)udu-∫0∞‍h(u)udu].


## 3. Numerical Test

In this section, we take numerical experiments to show the performance of our integral solution of the interface profile of grain boundary grooving. All the numerical tests are implemented on a Dell desktop 4311s (Intel Core2 Duo CPU E7400 2.8 GHz with 2.00 GB RAM). The software Matlab is used to calculate all the numerical solutions.


### 3.1. Numerical Calculation ofh(u)

We first calculate Schwartz functionh(u). To facilitate the calculation, we expand the function cos⁡τu in (18) in a power series and integrate term by term. Then we get
(21)h(u)=∫-∞∞‍e-τ4cos⁡τudτ=12∑n=0+∞‍(-1)nΓ((2n+1)/4)(2n)!u2n,
where Γ(u) is the gamma function. Since h(u) is a series form, we use finite terms
(22)h2N(u)=12∑n=0N‍(-1)nΓ((2n+1)/4)(2n)u2n
to approximate h(u) in the numerical tests. The truncated remainder for h(u) is
(23)R2(N+1)(u)=|h(2N+2)(ξ)(2N+2)!u2N+2|ξ∈(0,u)≤Γ((2N+3)/4)(2N+2)!u2N+2.
We exam the dependence of truncated error and number of terms for different u in h(u) with decreasing truncated error. Table 1 illustrates the details and Table 2 lists the values of h2N(u) and uh2N(u) when truncated error e=1e-7. It is shown that the number of terms used exploded for large u as truncated error decreased. However, we do not have to calculate the value of h(u) for those u far away from zero since h(u) is a Schwarz function.Table 1
Number of terms to truncated error for differentu.
u
0.1
0.5
1
2
5
10
20
N
(
e
=
1
e
-
7
)
2
3
5
7
13
25
54
(
e
=
1
e
-
8
)
2
4
5
8
14
26
55
(
e
=
1
e
-
10
)
3
5
6
9
16
28
57
(
e
=
1
e
-
12
)
3
5
7
10
18
31
60Table 2
Values ofh2N(u) and uh2N(u) (e1=10-7).
u
0.1
0.5
1
2
5
10
20
N
2
3
5
7
13
25
54
h
2
N
(
u
)
1.80974330
1.73738665
1.52470975
0.85205558
-
0.17151472
-
0.00269424
-
0.00000168
u
h
2
N
(
u
)
0.18097433
0.86869332
1.52470975
1.70411115
-
0.85757362
-
0.02694244
-
0.00003369


### 3.2. Integral Solution

In the following, we provide the numerical result of our integral solution.h(u) used in the integral solution is approximated by h120(u). The part of infinite integration ∫0∞‍h(u)udu is truncated to finite integral over the interval from zero to 20. To calculate numerical integration, trapezium integral formula is applied. Since the parameter (Bt)1/4 only impacts the scale of the solution y(x-,t), here we select (Bt)1/4=1. We plot the normalized profile y/m in Figure 2. To test the performance of our integral solution, we also compare the integral solution with Mullins' series solution [1] at |x-|≤5. Here Mullins solution of the interface profile is displayed as
(24)y(x-,t)=m(Bt)1/4∑n=0+∞‍anx-n,
where an is calculated by the recursion relation
(25)an+4=ann-14(n+1)(n+2)(n+3)(n+4),(26)a0=-121/2Γ(5/4),a1=1,a2=-123/2Γ(3/4),a3=0.
Figure 2 displays the comparison results. Solid line is Mullins' solution and dash line is our integral solution. Monitoring the evolution of the depth and the width of the groove profile, we find that our solution agrees very well with Mullins' solution according to the figure.Figure 2
Normalized profile shape due to volume diffusion.Here we have not compared the groove profile at largex- because both solutions can easily become unstable due to the limited machine precision and the truncation. For Mullins' solution, we plot the profile between |x-|≤15 with different truncation of the infinite series in Figure 3. The truncation of the infinite series is chosen as n=25,26,27,28,41,42,43,44, respectively. From Figure 3, we find that for n=25,26,41, and 42, the curves of Mullins's solution grow to positive infinity and for n=27,28,43, and 44, the curves grow to negative infinity. This phenomenon appears in a four-term cycle, and it is in agreement with the coefficients of the infinite series in (25) and (26). However, for our integral solution, with the properties of the Fourier transform, we have
(27)∫0+∞‍h(u)du=12∫-∞+∞‍h(u)du=π(12π∫-∞+∞‍cos⁡(u·0)du∫-∞+∞‍e-τ4cos⁡τudτ)=πe-04=π.
Hence the groove profile y(x-,t)→0 at a fixed time when x-→+∞. Since h(u) is a Schwarz function, we can easily conclude that y(x-,t) is stable with the properties above, and we do not have to examine the profile at large x-.Figure 3
Normalized profile shape due to volume diffusion.


### 3.3. Properties of Integral Solution

At last, we conduct some properties of our integral solution (20). We compute
(28)dydx-=m(Bt)1/4π[π-∫0x-‍h(u)du]
and find that the maximum of the surface is at the first occurrence of x- with
(29)∫0x-‍h(u)du=π.
Combined with the data in Figure 2, the maximum of the surface is bout y(2.3)=0.194. The depth of the groove measured in y direction from the maximum of the surface to the grain boundary is given by
(30)d=0.9738m(Bt)1/4,
whereas the separation of the two maxima is
(31)s=4.6(Bt)1/4.
From these expressions, we have the time independent ratio
(32)sd=4.7236m.The spatial oscillation of the groove profile at a fixed time is shown in Figure4. The solution shows a spatial oscillation with a rapid exponential damping in the amplitude. Beyond this, we can also get the interval between the groove root and the beginning of the tail from the figure.Figure 4
The spatial oscillation of the groove profile.The curvature of the groove profiley(x,t) is defined as K=|y′′|/(1+y'2)3/2. However, in small slope approximation, we assume that |y'|≪1. Hence, the curvature of the normalized profile here is K=|y′′|=h(x-). The result is plotted in Figure 5. All these properties are in agreement with Mullins' solution.Figure 5
The curvature of the groove profile.


## 3.1. Numerical Calculation ofh(u)

We first calculate Schwartz functionh(u). To facilitate the calculation, we expand the function cos⁡τu in (18) in a power series and integrate term by term. Then we get
(21)h(u)=∫-∞∞‍e-τ4cos⁡τudτ=12∑n=0+∞‍(-1)nΓ((2n+1)/4)(2n)!u2n,
where Γ(u) is the gamma function. Since h(u) is a series form, we use finite terms
(22)h2N(u)=12∑n=0N‍(-1)nΓ((2n+1)/4)(2n)u2n
to approximate h(u) in the numerical tests. The truncated remainder for h(u) is
(23)R2(N+1)(u)=|h(2N+2)(ξ)(2N+2)!u2N+2|ξ∈(0,u)≤Γ((2N+3)/4)(2N+2)!u2N+2.
We exam the dependence of truncated error and number of terms for different u in h(u) with decreasing truncated error. Table 1 illustrates the details and Table 2 lists the values of h2N(u) and uh2N(u) when truncated error e=1e-7. It is shown that the number of terms used exploded for large u as truncated error decreased. However, we do not have to calculate the value of h(u) for those u far away from zero since h(u) is a Schwarz function.Table 1
Number of terms to truncated error for differentu.
u
0.1
0.5
1
2
5
10
20
N
(
e
=
1
e
-
7
)
2
3
5
7
13
25
54
(
e
=
1
e
-
8
)
2
4
5
8
14
26
55
(
e
=
1
e
-
10
)
3
5
6
9
16
28
57
(
e
=
1
e
-
12
)
3
5
7
10
18
31
60Table 2
Values ofh2N(u) and uh2N(u) (e1=10-7).
u
0.1
0.5
1
2
5
10
20
N
2
3
5
7
13
25
54
h
2
N
(
u
)
1.80974330
1.73738665
1.52470975
0.85205558
-
0.17151472
-
0.00269424
-
0.00000168
u
h
2
N
(
u
)
0.18097433
0.86869332
1.52470975
1.70411115
-
0.85757362
-
0.02694244
-
0.00003369


## 3.2. Integral Solution

In the following, we provide the numerical result of our integral solution.h(u) used in the integral solution is approximated by h120(u). The part of infinite integration ∫0∞‍h(u)udu is truncated to finite integral over the interval from zero to 20. To calculate numerical integration, trapezium integral formula is applied. Since the parameter (Bt)1/4 only impacts the scale of the solution y(x-,t), here we select (Bt)1/4=1. We plot the normalized profile y/m in Figure 2. To test the performance of our integral solution, we also compare the integral solution with Mullins' series solution [1] at |x-|≤5. Here Mullins solution of the interface profile is displayed as
(24)y(x-,t)=m(Bt)1/4∑n=0+∞‍anx-n,
where an is calculated by the recursion relation
(25)an+4=ann-14(n+1)(n+2)(n+3)(n+4),(26)a0=-121/2Γ(5/4),a1=1,a2=-123/2Γ(3/4),a3=0.
Figure 2 displays the comparison results. Solid line is Mullins' solution and dash line is our integral solution. Monitoring the evolution of the depth and the width of the groove profile, we find that our solution agrees very well with Mullins' solution according to the figure.Figure 2
Normalized profile shape due to volume diffusion.Here we have not compared the groove profile at largex- because both solutions can easily become unstable due to the limited machine precision and the truncation. For Mullins' solution, we plot the profile between |x-|≤15 with different truncation of the infinite series in Figure 3. The truncation of the infinite series is chosen as n=25,26,27,28,41,42,43,44, respectively. From Figure 3, we find that for n=25,26,41, and 42, the curves of Mullins's solution grow to positive infinity and for n=27,28,43, and 44, the curves grow to negative infinity. This phenomenon appears in a four-term cycle, and it is in agreement with the coefficients of the infinite series in (25) and (26). However, for our integral solution, with the properties of the Fourier transform, we have
(27)∫0+∞‍h(u)du=12∫-∞+∞‍h(u)du=π(12π∫-∞+∞‍cos⁡(u·0)du∫-∞+∞‍e-τ4cos⁡τudτ)=πe-04=π.
Hence the groove profile y(x-,t)→0 at a fixed time when x-→+∞. Since h(u) is a Schwarz function, we can easily conclude that y(x-,t) is stable with the properties above, and we do not have to examine the profile at large x-.Figure 3
Normalized profile shape due to volume diffusion.


## 3.3. Properties of Integral Solution

At last, we conduct some properties of our integral solution (20). We compute
(28)dydx-=m(Bt)1/4π[π-∫0x-‍h(u)du]
and find that the maximum of the surface is at the first occurrence of x- with
(29)∫0x-‍h(u)du=π.
Combined with the data in Figure 2, the maximum of the surface is bout y(2.3)=0.194. The depth of the groove measured in y direction from the maximum of the surface to the grain boundary is given by
(30)d=0.9738m(Bt)1/4,
whereas the separation of the two maxima is
(31)s=4.6(Bt)1/4.
From these expressions, we have the time independent ratio
(32)sd=4.7236m.The spatial oscillation of the groove profile at a fixed time is shown in Figure4. The solution shows a spatial oscillation with a rapid exponential damping in the amplitude. Beyond this, we can also get the interval between the groove root and the beginning of the tail from the figure.Figure 4
The spatial oscillation of the groove profile.The curvature of the groove profiley(x,t) is defined as K=|y′′|/(1+y'2)3/2. However, in small slope approximation, we assume that |y'|≪1. Hence, the curvature of the normalized profile here is K=|y′′|=h(x-). The result is plotted in Figure 5. All these properties are in agreement with Mullins' solution.Figure 5
The curvature of the groove profile.


## 4. Conclusion

In this work, we mainly seek a novel solution of the interface profile of grain-boundary grooving by surface diffusion with small slope approximation. With Fourier transform theory and properties of partial differential equation, we derive an integral solution (6) for the Mullins problem. Taking variable substitution, we get the normalized solution (20). According to the numerical tests, we find that our solution is identical with Mullins' solution. The rates of growth of the depth and width remain proportional to (Bt)1/4. Besides, we can get the properties of the groove profile y(x-,t) with theoretical analysis, and these properties are no longer affected by the machine precision and the truncation of infinite series.



---

*Source: 520967-2012-11-14.xml*