- Email: [email protected]
A numerical technique for the analysis of coupled surface and grain-boundary diffusion
Acta metall, mater. Vol. 43, No. 4, pp. 1395-1406, 1995
Pergamon
0956-7151(94)00365-3
Copyright ~ 1995ElsevierScienceLtd Printed in Great Britain.All rights reserved 0956-7151/95 $9.50+ 0.00
A NUMERICAL TECHNIQUE FOR THE ANALYSIS OF COUPLED SURFACE A N D G R A I N - B O U N D A R Y DIFFUSION J. PANt and A. C. F. COCKS:~ Engineering Department, Cambridge University, Trumpington Street, Cambridge CB2 1PZ, England (Received 20 July 1994)
Abstract--The continuity conditions across a junction between a grain-boundary and a free surface are discussed for the coupled surface and grain-boundary diffusion problem. An appropriate treatment of the continuity is proposed, which forms the basis of a numerical procedure for an arbitrary network of grains containing pores of arbitrary shape and size. The analysis technique, when combined with a suitable time integration algorithm, is able to simulate physical processes such as powder sintering, diffusional void growth and creep crack propagation etc. which are dominated by the coupled diffusion mechanism under many practical circumstances. Currently only two-dimensional problems are discussed. The basic principles proposed in this paper, however, are also valid for three-dimensional situations. Simple numerical examples of powder sintering are given while the full potential of the numerical technique for the above physical problems will be exploited in a forthcoming paper.
1. INTRODUCTION At elevated temperatures, many physical processes which occur in polycrystalline materials require a coupling between surface and grain-boundary diffusion. Take the sintering of a ceramic powder compact for example: driven by the reduction of the associated total free surface energy, matter flows from gi'ainboundaries into pores, resulting in shrinkage of the pores and therefore the densification of the powder compact. Under many practical circumstances [1], matter is transported to the junctions between grainboundaries and pores by grain-boundary diffusion and is taken away from the pore tips and deposited onto the pore surfaces by surface diffusion. Another example is diffusional void growth in which the above process is reversed, i.e. matter flows from the pores into the grain-boundaries, driven by an external tensile stress state [2]. The solid state diffusion of matter is driven by the gradient of chemical potential of the diffusing species. This gradient of chemical potential is induced by the gradient of the free surface curvature in the case of surface diffusion and by the gradient of stress acting normal to the grain-boundaries in the case of grainboundary diffusion [3]. When a grain-boundary meets a free surface, the two processes are coupled. Because of the complex geometric problems involved, analytical solutions of the governing equations for this coupled diffusion problem, especially for the surface tCurrent address: Department of Mechanical Engineering, University of Surrey, Guildford, Surrey GU2 5XH, England. :~To whom all correspondence should be addressed.
diffusion part, are difficult to obtain. A given physical problem can often be simplified by first examining two extreme situations: (a) the fast surface diffusion limit; and (b) the fast grain-boundary diffusion limit. At the extreme of fast surface diffusion, it is assumed that surface diffusion is much faster than boundary diffusion and matter transported to the pore tip is instantly taken away and distributed along the pore surface (or an abundant source of matter is available at the pore tip to be taken away by grainboundary diffusion). As a result, an equilibrium pore shape is maintained and the pores simply serve as boundary conditions for the grain-boundary diffusion problem in the form of capillarity stresses. For geometrically symmetric problems, analytical solutions for the grain-boundary part of the problem can be easily obtained. This approach has been adopted by many researchers to model void growth [2], creep crack propagation [4] and powder sintering [ll. At the extreme of fast grain-boundary diffusion, more attention is paid to the evolution of pore profile while the grain-boundaries only act as sources or sinks of matter and as a boundary condition for the chemical potential which links the pore tip curvature with the remote stress. Chuang and Rice [5] obtained analytical solutions at this extreme under a steady state assumption for a one pore and one boundary system. Chuang et al. [6] also investigated the coupled problem and developed a class of self-similar solutions for the one pore one boundary system. For more general situations, analytical solutions are often too difficult to obtain and numerical methods are often adopted. Pharr and Nix [7] studied cavity
1395
1396
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
growth controlled by surface diffusion using a finite difference scheme. Similar schemes were used by Binh and Uzan [8], and Takahashi and Nishiguchi [9] etc. for problems of solid drop formation at the tip of a solid wedge and void growth controlled by surface diffusion. A simplified treatment of surface diffusion during sintering has been carried out by Ashby [10] and Swinkels and Ashby [11]. Recently Svoboda and Riedel [12] developed similarity solutions at both extremes for small scale diffusion. They also suggested a variational solution for the full coupled situation under the small scale diffusion assumption. Bouvard and McMeeking [13] and Svoboda and Riedel [12] used a scheme of combining the analytical solution for grain-boundary diffusion with the finite difference procedure for surface diffusion for a uniform array of particles. This scheme has provided useful results for the investigated sintering problems but it can only be applied to situations where analytical solutions for the grain-boundary part of the problem are available and where the pore shape is highly symmetric. For more general situations, two issues have to be dealt with in order to solve the coupled diffusion problem: (a) solution of the grain-boundary part of the problem for any arbitrary boundary network and (b) the continuity conditions at the junctions between pore surfaces and grain-boundaries. The authors have developed a finite element type of numerical technique which can solve grain-boundary diffusion problems for any arbitrary boundary network with pores of equilibrium shape [14, 15]. A natural extension of this technique is to combine it with the finite difference procedure for the surface diffusion part of the problem so that the coupled problem can be solved. However, the continuity conditions have not been fully examined in the literature. In this paper, we first discuss the continuity conditions, clarify some areas of confusion that often arise in the literature and then propose a numerical scheme which can solve the coupled diffusion problems for any arbitrary grain-boundary network with pores of arbitrary shape. Simple examples of powder sintering are given while the full potential of the numerical technique will be exploited in a forthcoming paper.
¢
¢
E~
E~
Fig. 1. Polycrystalline material with pores. According to Fick's law, the diffusive flux j, defined as volume of matter flowing across unit area perpendicular to the flux direction per unit time, is linearly dependent on the gradient of the chemical potential # of the diffusing species Dt$ t3/a
j
k T ~s
where D is diffusivity, 6 is the thickness of the layer through which the material diffuses, k is Boltzman's constant, T is the absolute temperature and s is a local coordinate along the diffusion path.
2.1. Grain-boundary diffusion Along a grain-boundary, the gradient of the atomic chemical potential is induced by the gradient of stress o acting normal to the boundary [3] # = - f~o
A schematic example of the coupled diffusion problem is shown in Fig. 1. In the current paper we confine our attention to two-dimensional problems. Grain-boundaries are assumed to be straight and grains are assumed to be rigid. Grains either side of a grain-boundary are free to slide with respect to each other. Their relative motion in the direction normal to the grain-boundary, however, has to be accommodated by grain-boundary diffusion.
(2.2)
where f~ is the atomic volume, therefore
Db6bf] ~0" k T ~3s
A- - -
(2.3)
in which D b is the grain-boundary diffusivity and 6b is the grain-boundary thickness. As matter is deposited onto (or removed from) a grain-boundary, a normal velocity Vb between grains either side of the boundary results. Matter conservation requires that
eJb
~ s + Vb= O. 2. THE GOVERNING EQUATIONS
(2.1)
(2.4)
At each junction where a number, n, of grain-boundaries meet, the flux Jb must satisfy the condition that ~ Jb,i " ti = 0
(2.5)
i=l
where ti is a unit vector in the direction of the ith boundary and pointing away from the junction.
2.2. Surface diffusion Along a free surface, the gradient of the atomic chemical potential is induced by the gradient of the
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION free surface curvature [3]
1397
where ?s is the surface free energy and x is the principal curvature of the surface. It is positive for a concave free surface. So we have
Furthermore, the surface tension ),~ and grainboundary tension ~'gbmust be in equilibrium as shown in Fig. 2. This equilibrium requirement leads to a discontinuity in the tangent of free surface at the junction which is often referred to as the dihedral angle 20
D~6~7~ ~ j~ - - kT ~s
cos 0 = - - . 2?.,
# = -- ~?s x
(2.6)
in which D~ is the free surface diffusivity. As matter is deposited onto (or removed from) a free surface a normal velocity v~ of the free surface results. Matter conservation requires that
~j~
0~ + t,~ = 0
~gb
(2.7)
Finally, the chemical potential must be continuous (otherwise there would be an unbounded flux j there). Because of equations (2.2) and (2.6), the continuity of the chemical potential is always expressed as O'tip = ?s/('tip
(2.8)
with v~ taken as positive when matter is added to the surface. The surface velocity v~ due to matter redistribution is imaginary, i.e. it is not a velocity of any material element. It only makes sense when a point-to-point correspondence is defined between the surfaces before and after matter redistribution. Equation (2.8) implies that this point to point correspondence is defined by the normal of the surface before matter redistribution, i.e. v~ represents the 'thickening' rate of the surface as matter is deposited onto the surface.
2.3. Global equilibrium The stresses, a, normal to the grain-boundaries, the free surface tension, ?s, and the applied remote stress must be in equilibrium.
(Jb)tip = -- (J + )tip + (]" s )tip
(2.11)
in the literature. However, (2.10) implies that xtip does not exist. We examine this feature of the continuity requirements further in Section 3.
2.5. Nondimensionalization Let d be a characteristic diffusion distance (the grain-size for example) and O-0 a suitable reference stress (the applied remote stress for example). We can then define a reference rate
Eb--
k~
d3
(2.12)
which has the dimensions of strain-rate, and write the governing equations in terms of dimensionless quantities
= xd, 6 = O-/o-o,f=J/(Ebd2), f = V / ( d b d ) , t ' = ~b t
2.4. Continuity When a grain-boundary meets a free surface, several continuity conditions must be satisfied. First matter must be conserved at the pore tip
(2.10)
(2.13)
with any lengths in the physical description scaled by d. Equations (2.3) and (2.7) then become = ~--~
(2.14)
js = C~z_
(2.15)
Ds3s ?s Cs - - = Of~ Db3b O-od
(2.16)
Ds6s 7s 0 = 0 - ~ and ~s = ~ .
(2.17)
(2.9)
where + and - indicate the two free surfaces joining the grain-boundary and the sign convention for these local fluxes are defined in Fig. 2.
and OS
where
i.e.
Fig. 2. A junction between a grain-boundary and a free surface. AM 43/4~.5
Other equations keep their original form in terms of the corresponding nondimensionalized variables. Three groups of non-dimensionalized material constants emerge from the above analysis: O, which represents the relative importance of the two diffusion processes; ?s/?gb, which determines the dihedral angle at the pore tips and 7s/o-0d, which represents the relative significance of the applied stress and the surface tension.
1398
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION 3. THE CONTINUITY CONDITIONS
For the purpose of discussion we consider the simple problem of Fig. 3(a) to examine the coupling between grain-boundary and surface diffusion. Equilibrium shaped pores with spherical caps of radius, R, are periodically distributed along an infinitively long grain-boundary, which is subjected to a remote normal stress, E~. Because of symmetry, only the part of the problem shown in Fig. 3(b) needs to be considered. At the pore tip, continuity of potential, equation (2.11), requires that O'ti p = ~)s/R and a potential gradient is set up to drive material away from the tip at a rate 1 fs ~ = 2~s(~ + sin 0 ) - - T= - 2 l ' ( 1 - ~ )
(3.1)
where the reference stress a0 has been identified with the applied stress E® and the characteristic diffusion distance d with the half ligament length a. Our assumption of an equilibrium shape for the pore, however, means that~ = 0 at the pore tip and we are unable to satisfy equation (2.9). This inconsistency indicates that within the framework of steady state diffusion theory employed here we cannot completely define the pore profile. A pore tip profile (or, at least, a chemical potential profile, which need not be the same thing) must therefore develop which allows both equations (2.11) and (2.9) to be satisfied. The chemical potential at the pore tip can therefore only be determined by considering the coupled diffusion process. The chemical potential has to be continuous everywhere otherwise there would be an infinite flux. As a result, the curvature and therefore the local profile of the pore in the region very close to the pore tip are also determined by the coupled diffusion process and can not be arbitrarily prescribed. We refer to this local profile close to the pore tip as the steady state tip profile. Any departure from this profile will be adjusted to the steady state instantaneously and it is this adjustment that kicks off the diffusion process along the pore surface in the above example. E
0
(a)
t
t
0
tz
t
t
In analytical solutions [5, 6], this problem is avoided by directly seeking the steady state profile. The solution then is not designed to satisfy any initial pore profile condition. In numerical procedures, however, each time step has to start from a pore profile determined from integration of all the previous timesteps. This profile is generally not the exact steady state one because of integration errors and node re-arrangement as the pore profile evolves. The numerical solution also has to start from an initial profile, a circular cap for example, which is not necessarily the steady state profile. Finding the correct chemical potential at the pore tip at each time step is then the critical step for the numerical solution of the coupled diffusion problem. This will be further discussed in the next section. When analysing a symmetric coupled diffusion problem using the finite difference method Svoboda and Riedel [12] determined the 'pore tip curvature' by drawing a circle through the pore tip and its neighbouring finite difference node, ensuring that the angle between the two symmetric curves satisfies the dihedral angle relationship of equation (2.10) at the pore tip. The resulting pore tip curvature provides a surface flux that violates matter conservation. This problem was solved by shifting the pore tip coordinate iteratively until matter conservation was satisfied. This procedure found the correct pore tip chemical potential but it is difficult to apply the technique to other situations. If we consider an asymmetric pore profile as shown in Fig. 2, then we have two such circles obtained from either side of a pore tip giving curvatures of opposite sign. The procedures described in the following section were designed to allow this more general situation to be evaluated. 4. A COUPLED FINITE ELEMENT AND FINITE DIFFERENCE METHOD
4.1. Finite element formulation for grain-boundary diffusion For the grain-boundary diffusion part of the problem, Cocks [16] extended an earlier result of Needleman and Rice [17] to establish the following variational principle: among all the kinematically admissible fields .~ [i.e. those satisfying matter conservation requirements (2.4) and (2.5)], the true field [i.e. that also satisfying (2.3) and global equilibrium] makes the functional 1"10 a minimum, where ]'I 1 (" -2 -T 0 = ~ | Jb dA q"~.~tipm'Jb
~JAb
-- ftiP ~' " ~" dS + ~-~fs ~Tbsin 0 JST
(b)
r
[ .L.
a
Fig. 3. Periodically distributed pores along an infinitely long grain-boundary.
(4.1)
tip
and the summation is performed over all the pore tips, m is a unit vector along a grain-boundary towards the pore tip, ST is the part of the specimen boundary where tractions T are prescribed, ~" is the velocity of ST, )7 is the nondimensionalized surface tension and 0 is the dihedral angle.
1399
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION Cocks [14] introduced the constraint of (2.5) into 1"10and obtained a modified functional
~ l'I=l~10+ ~
(2~m.]b ~
junction\
gb
a
(4.2)
/
where 2 represents a set of Lagrange multipliers. When constructing approximate solutions using (4.2), does not have to automatically satisfy (2.5). Surface diffusion and grain-boundary diffusion are only coupled through the pore tip. For the coupled diffusion problem, the variational principle (4.2) still represents the governing equations for the grainboundary part of the problem except that fft~p is unknown and the variation of)~ must satisfy the extra constraint of equation (2.9). This constraint provides the extra equation required to determine ¢Ttip' A 'grain-boundary element' proposed by the authors [15] is shown in Fig. 4. For a straight grain-boundary, its opening rate Vb can be directly related to the velocities of its two associated grains ~/1, ~7~, and 05~ and t~2, 172,052 ~Tb= [/~ ][00]
(4.3)
[/~]=[--nl --n2 ( £ + f l ) n l n2 - - ( g + ~ ) ]
(4.4)
where
and
j
~ram-boundary "
B -
v....
~
j~
Fig. 5. Discretized free surface. where [(~] = [n I s n2,~ --(0.5s 2 -k-~l)
--nlg--nS(O.5gz +~)l] [/-'~] ~---[/gl91 (~1 /~2 ~TZO)2J~,0]T
(4.8) (4.9)
and J~,0 is the matter flux across the origin of the local co-ordinate system. Substituting the matrix expressions ofj~ and Vbinto (4.2), we have the discrete form of f'I l=l = 1 I~R~ T + ~
,~ ,v
-~ ~tip~"C ["T "-~~s]~) T -- 1~'~ T (4.10) [00] = [ffl vl °51 u2 f2 °52]T.
(4.5)
In (4.4) g is the local co-ordinate along the grainboundary, with its origin located at the mid-point of the boundary, n i, i = 1, 2, is the normal to the grain-boundary, the length ~ is the negative value of the local co-ordinate where a line drawn parallel to n through the centre of the first grain intersects the grain-boundary facet, Fig. 4, and similarly T2 is the negative value of the local co-ordinate where a line drawn parallel to n through the centre of the second grain intersects the grain boundary facet. From equation (2.4), we have Jb =
fo
-- [B][O0] dg +J~,0
(grain 2)
$
16
61']=0
(4.11)
(4.6)
which gives the following set of linear simultaneous equations AX = R (4.12)
(4.7)
where A is assembled from 1~, ~ ' and (~", R is assembled from Fsl] and 1~, and finally X is the global matrix of unknowns consisting of 1~, .~ and :~tip. Since ~ p does not join the variation, (4.11) has not provided enough equations for all the unknowns. The extra equations are provided by the surface diffusion part of the problem.
which can be written as = [(~ 1[01
where 10, I~' and Q" are the global assemblies of [U ], 1~6 and l~g are the global assemblies of [L~0], and/~, ~tip and ~ are the global assemblies of 2, atip, and ~" respectively. 1~, (~', (~" and n are the global stiffness and constraint matrices, with their detailed expressions given in [15]. The true grain centre velocities •, q, 05, the true matter fluxes across the origins of the grain-boundary elementsj% and the Lagrange multipliers ~'introduced at the grain-boundary junctions are those that make the functional l"l stationary, i.e.
4.2. Finite differenceformulation for surface diffusion
(grain 1) Fig. 4. Grain-boundary element.
For surface diffusion, we adopt a finite difference scheme. The pore surface is discretized into knot points as shown in Fig. 5. We shall refer to these points defining the pore profile as P-points. The curvature x and the centre of curvature of the surface at a P-point are determined by constructing a circle
1400
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
through the P-point and its two neighbouring Ppoints. For a given pore profile, expressed in terms of the P-point co-ordinates, the curvature is known at all P-points except at the pore tips where curvature is not defined. A point bisecting two neighbouring P-points is referred to as a B-point. Matter flux f along the surface across any B-point can be determined according to (2.15) (4.13)
Js,ib = Cs (/~ib,+P -- Kib,-P) ASib
It proves convenient to define an imaginary pore tip curvature through (2.6) K'tip --
(4.16)
which has nothing to do with the pore tip curvature prescribed by the initial geometry, because #tip is unknown. The matter fluxes across the two B-points next to the pore tip are given by -
(-, /('tip -- •+2
Js'+1=vs
except for those points next to the pore tips, where gib.+V and x~b,-e are curvatures of the two P-points neighbouring the ith B-point, and Asib represents the distances between these two neighbouring P-points. The normal velocity of the pore surface at any P-point can then be determined according to (2.8)
/~tip
~')Ts
-
AS+l---- ' Js'-l=
C K'-2 -- -Ktip
s
AS_l
in which g±2 are known. The velocities of the two P-points next to the pore tip are given by
V'S,+2~----2 J~,+t-A,+2 AS+I "[- As+2
~,,-2 = - 2 AL,-2 ~ _ ~-L,-1 T~_ 2
17s,ip=
Js,ip,+O -- Js,ip,-B Asip
in which ~,±2 are known. Combining equations (2.11), (4.7), (4.17), (4.18) with (4.15), we obtain [C]fip[U] + Q1 t~tip + Q2 = 0
With reference to the pore tip shown in Fig. 5, we refer to ~,,ip as the matter flux into the pore tip from the grain boundary, f,+~ and f , - i the matter fluxes across the two B-points next to the pore tip, and ~7~,~p,i v~+2 and ~,-2 the velocities of the pore tip and the two P-points next to the pore tip. The dihedral angle is enforced between the line linking the pore tip and its neighbouring P-point and the grain-boundary. We further assume that matter deposits uniformly on the two elements adjacent to the pore tip. As a result there is no rotation of these elements or the boundary that meets the tip, although this boundary might be required to migrate. We examine this assumption further in Section 4.3.2. Matter conservation requires that the amount of material deposited on to the surface between the two B-points is equal to the amount flowing into this region and we have As+ l
As_ 1
Jb,tip + f.+ 1 -- f,_ 1 = vs,+2 -----~ + ~Ts,_2 ---~
(4.19)
where [(~ ]tip is the matrix of equation (4.8) at the pore tip, [U] is a matrix which contains the degrees of freedom associated with the grain-boundary [equation (4.9)] and QI and Q2 are known coefficients calculated from g_+2,J~,+2, AS±I, As±2 and Cs
QI = C ,
+ A-~_I+ AS+l + As+2
1
-~ As l+As_2 Q2 = - - C s Y K+2 --+
LAs+l
4.3. Coupling between surface and boundary diffusion
_
(4.18)
(4.14)
except for the pore tips and those P-points next to the pore tips. Here f,ip,+B and f,ip.-B are fluxes across the two B-points neighbouring the ith P-point, and Asip is the mean of the distances between the ith P-point and its two neighbouring P-points. The direction of ~ is taken as that along the direction of the normal of a P-point defined by the P-point and its centre of curvature. From the resulting values of q~ and their directions, the pore profile can be updated by a suitable time integration scheme. However, es at the pore tips and their immediate neighbouring P-points cannot be determined without considering the boundary diffusion part of the problem.
-
(4.17)
l K+2
t~-2
A--~_i+ As+l+As+2
re_2
As_~-AS_~] As,
-
AS+~
-~ I~S_1 -t-'-* IxS-2 Js-Z--As+l + "
~
] _l
By appending equation (4.19) to equation (4.12), we obtain the appropriate number of equations to determine all the degrees of freedom and the values of the chemical potential at the pore tips. 4.3.1. Application of the continuity condition. As an example, we consider the problem of Fig. 2 discussed in Section 3. Consider, as before, the situation where the initial pore profile is a circular cap. If it is uniformly divided into elements of the same length, As, we find, from equations (3.1) to (3.4) and (4.19) /('tip --
+
~
-- sin 0
(4.20)
1+-O (4.15)
where As± l are the two distances between the pore tip and its two neighbouring P-points.
where the half length of the grain-boundary between the two neighbouring pores a has again been taken as the characteristic diffusion distance d and E~ has been interpreted as the reference stress tr0. /~ is the
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION normalized radius of the circular cap a n d / the half distance between the two neighbouring pore centres. The flux across the B-point next to the pore tip is then given by, equation (4.17)
('
J s = C ~ ss ~ - Xtip
) As 1+-O
(4.21)
It proves instructive to examine these expressions in certain limits. We would expect the discrete description here to reduce to a continuum description in the limit As ---, 0. In this limit equations (4.20) and (4.21) become l
f~ip = ~
(4.22)
and
js (l+sin0) r -
"L
) 423,
i.e. the tip curvature is the same as that for the rest of the pore, and the fluxf towards the pore tip (along the two surfaces) is equal to the flux~ away from the tip along the grain-boundary [compare equations (4.23) and (3.1)]. It is also interesting to note that the surface flux at this initial stage is independent of the magnitude of the surface diffusivity. We can gain further insight by allowing As to increase slightly. Then
~,i~=~ + g [ ~ - s i n 0 - ~
,] +O
,
l
.(4.24)
This equation demonstrates that flow of material away from the pore tip by boundary diffusion generates a small perturbation in the pore tip profile which drives the diffusive flux of material to the tip along the pore surface. As the void shrinks, this perturbation propagates along the surface, allowing for internal redistribution of material. If As~® is small then equations (4.24) and (4.21) reduce to (4.22) and (4.23). If ® is large, i.e. if the surface diffusivity is much greater than the boundary diffusivity, then the continuum response can be well approximated using a coarse mesh. In this limit only small perturbations in the pore profile will be set up to drive the diffusive flux of material around the surface throughout the diffusion process. For moderate values of O, a suitable initial node spacing can be obtained by requiring that As/® is small. We now turn to the limit of ® ~ 0 for finite As. Equation (4.20) then becomes
This is precisely the tip boundary condition employed by Pharr and Nix [7] for their computation in the limit of fast boundary diffusion, and simply reflects the fact that in this limit the potential gradient along the boundary is small, with the normal stress a n approximately equal to the net section stress along the entire length of boundary. In this limit equation (4.21) becomes ~=~s[)Ts(l+sin0)-/"]
Vs,tip" n+ = vs,+2 • n+
(4.27)
vs.,ip • n_ = vs,_2 • n_
(4.28)
where n+ and n_ are the normals of the two straight lines linking the pore tip and its two neighbouring P-points, see Fig. 6. The dihedral angle comes into play through the dot product between vs.tip and n_+. It should be emphasized that not only the magnitude, but also the direction of the pore tip velocity, is determined from the coupled diffusion process. If a pore profile is not symmetric about a grain-boundary, then the above equations give a pore tip velocity which is not along the grain-boundary direction. Since, in the derivation of the above equations, we have not allowed the pore tip profile to rotate, the grain-boundary must migrate to retain its orientation with respect to the pore tip. This assumption is valid provided the characteristic rate of grain-boundary migration is much faster than the characteristic rates of the two diffusion processes considered here. If this is not true then the finite element formulation for
V,.e
/3*,÷~
i.e. #tip = - - ( E ~ ! - y~sin 0 ) fL
(4.25)
(4.26)
i.e. the surface flux now depends on the surface diffusivity and the flux at the first B-point depends on the length of the first element. This reflects the fact that there is significant rearrangement of material next to the pore tip and the amount of material flowing to the tip decreases with increasing distance from the tip. 4.3.2. Pore tip velocity. As soon as grip is known, the velocities of the two P-points next to a pore tip can be determined through equations (4.17) and (4.18). Because the dihedral angle is enforced between the pore tip and these two P-points, the pore tip velocity is related to the velocities of these two P-points in the following way
/g~ip= - - sin 0 7~
1401
Fig. 6. Pore tip velocity.
1402
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
grain-boundary diffusion has to be extended to allow for curved grain-boundaries and the coupling conditions (4.15), (4.27) and (4.28) have to be further generalized. This is beyond the scope of the current paper, although there is no difficulty in principle in performing such an extension and generalization within the framework of the current numerical technique. 5. UPDATING OF THE PORE PROFILE
Equations (4.12) and (4.28) can be solved by any standard routine for linear simultaneous equations. The solution provides velocities of each individual grain and each P-point of the pore surface. The grain-boundary network and the pore profile can then be updated according to these velocities. The updating procedures for the grain-boundary network has been discussed in [15], therefore only the procedures for updating the pore profile are discussed here. Euler's method is used to perform the time integration of the pore profile
xi(t -q- At ) = xi(t ) + VsnxA i yi(t + A t ) = yi(t ) + VsnyAt i i
(5.1)
where x i, yi are P-point co-ordinates, vis the P-point velocities, n i, n~ the P-point normals and At the time step. The time step At is determined from the maximum surface velocity vmax experienced by the Ppoints at the current time step, such that the displacement increment due to At and v max is kept constant for all time step increments. If we define the mean distance between all the neighbouring P-points of the initial pore profile as dmean, then this constant displacement increment can be set to a prescribed percentage, EI, of d. . . . 13max
At = Etdmea---~'
(5.2)
During the numerical simulation, the pore profile experiences large scale evolution. The finite difference mesh has to be modified accordingly. A minimum critical distance between any two neighbouring Ppoints is set in terms of a prescribed percentage, Emin, of dm.... a P-point will be deleted if the distance between the P-point and its neighbouring P-point is shorter than e~ndme~.. Similarly, a maximum critical
Fig. 7. The sintering of a hexagonal array of circular cylinders.
Fig. 8. Detailed initial profile of a pore tip. distance, £maxdmean, is set so that a B-point will be converted into a P-point if the distance between its two neighbouring P-points becomes longer than this critical distance. Like all numerical procedures, the density of the finite difference mesh has to be determined by experience. The three adjustable parameters, Et, Emi, and Em~, have also to be determined by trial and error. This will be further discussed with regard to the numerical examples presented in the next section.
6. NUMERICAL EXAMPLES In this section we use the procedure described in the previous two sections to analyse four simple problems that relate to the sintering of fine grained materials. First, we consider the pressureless sintering of a close packed array of circular cylinders. This problem has been examined by Svoboda and Riedel [18], and we are able to compare our results directly with theirs. We then examine the situation where the array of cylinders is subjected to a non-hydrostatic stress state, and follow the irregular evolution of the pore shape. Next, we consider the co-sintering of two cylindrical particles of different radius. Finally, we examine the free sintering of a row of cylinders of different sizes and follow the process of dynamic grain-growth. Further examples of the technique are presented in a forthcoming paper.
6.1. Sintering of a close packed array of cylinders As our first example, we use the numerical technique established in this paper to simulate the sintering of a hexagonal array of circular cylinders. The representative element used in the simulation is shown in Fig. 7. For the loading cases considered in this paper, a smaller representative element can be identified. The representative element of Fig. 7 is used deliberately to check the symmetry of the numerical results. All the lengths in the example are scaled by the initial radius of the cylinder. Each free surface between the two contacts are uniformly divided into 60 nodes. The pore tips are given an initial profile as shown in Fig. 8. First, we consider the case of pressureless sintering. Figure 9 shows the evolution of the pore shape and Fig. 10 shows the variation of the relative density and 'radius' of the pore tip as a function of the normalized time. The corresponding stages between Figs 9 and 10 are marked by (a), (b), (c) . . . and so on. The nondimensionalized material parameters used in the simulation are Ds6s/Db6b= 1, )'gb/(27,)=0.5 and ys/(dtro) = 0.01. After several test runs of the program, suitable parameters controlling the timestep
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
(al (b)
{c)
(a)
(e)
(£)
Fig. 9. Evolution of pore shape during pressureless sintering. length and node elimination were determined: et---0.02 and Emin = 0.05. The computer simulated pore shape evolution follows the pattern observed experimentally [19]. The imaginary pore tip radius calculated from the pore tip chemical potential initially increases as a result of material depositing onto the pore tip; then, after the pore has reached its equilibrium shape, it decreases as a result of pore shrinkage. There is a certain degree of noise in the pore tip chemical potential due to its high sensitivity to the detailed local profile. This noise was smeared out when creating Fig. 10. Matter conservation requires that the total area (excluding the pores) should remain constant. The actual change of this area 1.0. 0.98
Relative density ....... "Poretip nadius" J • • Results of Svoboda & R i e d e l /
during the simulation can therefore be taken as a measure of computing accuracy. Over the range of density variation shown by Fig. 10, about 0 . 2 o of the total area was lost. Also shown in Fig. 10 are the results obtained by Svoboda and Riedel [18]. Although different methods were used to deal with the pore-tip continuity in these two analyses, the two simulations agree almost exactly with each other. The method used by Svoboda and Riedel is, however, unable to deal with more general situations like non-hydrostatic loading because it assumes the development of a symmetric pore profile. The method developed in this paper was used to simulate sintering under a remote stress state Y,x x = - 0 2 1 and : ~ y y = - 0 . 1 2 for the situation D s r s / D b ~ b = 5 with 7gb/(2?s) = 0.5 and ?J(da0) = 0.01 as before. Figure 11 shows the simulated densification curve with the evolution of pore shape. This macroscopic stress field produces a mean compressive stress across the angled boundaries of Fig. 7 and a mean tensile stress across the vertical boundary of the repeating cell. The effect of the compressive stress across the angled boundaries is to accelerate the flow of material into the pore from these boundaries, while flow of material into the pore from the vertical boundaries is retarded by the microscopic tensile stress. As a result, in the early stage of densification, the pore develops an anisotropic shape as it shrinks more perpendicular to the direction of maximum compressive stress. In the later stages, material can readily flow around the small pore, allowing it to adopt an equiaxed equilibrium profile. Overall, the superposition of this net hydrostatic stress state accelerates the rate of densification (compare Figs 10 and 11), as one would expect. Reducing the magnitude of Y,xx to zero, increases the magnitude of the tensile stress across the vertical boundary to a level which now causes material to flow away from the pore and to plate out on this boundary. The grains either side of the boundary move apart and the contact patch shrinks in size.
098
f' .0.8
1403
13
0.97 ~ .~ 0.96
"~
J 0.96
t"
.0.6
"'"..
~2
"~
0.94
_0.4 o.94 _
~ 0.93
• ..........
?. 0,92_ b: £
...... ""..
02 ~ "''.
0.92 0.91
I
I
I
t
I
r
5
10
15
20
25
30
x 10 ~ N o r m a l i z e d time Normalized time
Fig. 10. Variations of relative density and 'pore-tip' during pressureless sintering.
Fig. 11. Variation of relative density and evolution of pore shape under a non-hydrostatic stress state of :~xx= -0.01 and ~yy= -0.12.
1404
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
Eventually the grains separate as they are wedged apart and a grain switching event occurs. The implication of these results for developing the constitutive models of sintering is that non-hydrostatic stress states can trigger anisotropic behaviour of the material in the early stage of densification, but in the later stages surface diffusion is sufficiently rapid to redistribute material around the surface of the pores, allowing them to adopt equilibrium profiles and the material response to become isotropic.
6.2. Sintering of cylinders of two different sizes Most commercial powders consist of particles with a range of sizes. This difference in particle size can promote grain-growth as a compact densities during sintering. In the later stages of densification (stage II), when the pores are isolated, the grain-growth kinetics is determined by the pinning effect and rate of migration of the pores. Constitutive models for this process have been described by Du and Cocks [20], based on the grain-growth studies of Hillert [21], Shewmom [22], Brook [23] and Ashby [1]. These models predict a much faster rate of grain-growth during pressureless sintering than HIPing; which is consistent with experimental observations [24]. In the early stages of densification (stage I), however, when
(a). [ = 0
(c). ~ = 6.80 • 10 -3
(e), f = 2.30* I0-I
(b). [ = 4.31 * 10 -4
(d). [ = 1 . 3 9 , 1 0 - 2
(f). ~ =
5.63
Fig. 12. Sintering of two cylinders of different radii.
J 0'12"1"~~ . . . . .
L
0.1
~
of ~
cyfed~
-3.26
/]
"~-,
0.06
-3.24
-3.22 ~.,
-3.2 , %',
I~ ~ ~
-3.~
5
10
t5
20
25
5O
55
4O
45
5O
55
6O
65
Normolized time
Fig. 13. Variation of cross-sectional areas of the two cylinders vs non-dimensionalized time. distinct necks exist between the particles a correlation between grain-growth and densification is generally observed [25] which is not captured in existing models [20]. Kellett and Lange [25, 26] have obtained expressions for the free energy of simple arrays of irregular particles and demonstrated how material redistribution, driven by changes in this free energy, can result in both densification and grain-growth, but full kinetic equations have not been developed for this process. As a preliminary study, to aid the development of appropriate kinetic equations for the coupling between grain-growth and densification in stage I, we examine two simple representative situations: first we simulate the sintering of two cylinders with a radius ratio of 5, Fig. 12(a); then we consider the situation where these two cylinders represent a repeating cell within a row of cylinders which is allowed to sinter. An initial neck size of 0.1 times the radius of the small particle was used in all the simulations. Other material constants and parameters used w e r e : Ds~s/Db~ b = 1, ygt,/(2ys)= 0.5, 7J(dtro)= 0.0l, Et = 0 . 0 2 , £min ~--" 0 . 0 5 and Emax= 1.2. Figure 12 shows the evolution of geometry for the two sintering particles. There is approximately an order of magnitude difference in time between the different stages of Fig. 12. Initially, the diffusion distance is short, as material flows from the neck and deposits in the crevice between the two particles, and the profile evolves at a fast rate. Material transported away from the boundary between the two particles is deposited non-uniformly on the surfaces of the particles, with most material depositing on the larger particle (it is possible to view this process in an alternative way: one can assume that the material removed from the boundary region deposits uniformly onto the two particles, and superimposed on this there is a net flow of material from the smaller particle to the larger particle, driven by the difference in the curvature between the two particles). This non-uniformity of deposition results in a pore-tip velocity in the direction of the smaller particle. As shown in Fig. 13, there is a gradual increase in the volume of the larger particle and a corresponding
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION decrease in the volume of the smaller particle. Two factors contribute to this: the net flow of material between the two particles; and grain-boundary migration, which results from the pore-tip motion. The grain-boundary gradually moves towards the top of the small particles, Fig. 12(c), and eventually meets the free surface, Fig. 12(d), leaving a bump on the large cylindrical particle (numerically, it is not possible to completely follow this process; in the simulation the grain-boundary was deleted when the area of the smaller particle became less than 0.2% of the total area of the two particles). Material contained within this bump subsequently redistributes around the surface of the particle, Fig. 12(d-f). As it does so the characteristic diffusion distance gradually increases, and the process slows down as the particle asymptotes towards a circular equilibrium profile. Equivalent results for the row of cylinders are shown in Figs 14 and 15. Now, the material from the smaller particle is transferred to the larger particles either side of it. As it does so, the two boundaries between the particles move towards each other, Fig. 15(b,c), Until they meet. At this point we could assume that the boundaries annihilate each other, but the larger particles will generally have different crystallographic orientations and it is more appropriate to assume that the boundaries combine to form a new grain-boundary. As before we cannot completely follow this process, and in the simulations the bound-
(a) t = 0
(b). ( = 8 , 6 7 * 10 -~
0.12" ~\
.....
0.1.
~
1405 f3.26
A r ~ of =rrd
3.24
3.2
~'
338 o
316
5
10
15
20
25
50
.10" Normolized
time
Fig. 15. Variation of cross-sectional areas of the large and small cylinders as a function of non-dimensionalized time before the small particle disappears. aries were replaced by a single boundary, located midway between the original boundaries, when the area of the smaller particle became less than 0.2% of the total area of the two particles, as before. Since a small particle is now sandwiched between two larger particles material is transferred from both ends of the particle and the time taken for the two boundaries to meet is roughly half the time required for the boundary to disappear from the small particle of Fig. 12 (compare Figs 13 and 15). The problem is symmetric after the smaller particle disappears and there is no further migration of the grain boundary. Material continues to flow from this boundary and deposit adjacent to the pore tip, causing the neck between the particles to grow, Fig. 14(d-f). The distance between the centres of the two large particles throughout the sintering process is plotted as a function of the dimensionless time in Fig. 16. Initially the rate of approach of the two centres is rapid as the small particle disappears. This rate gradually decreases as the characteristic diffusion distance increases and the configuration approaches equilibrium, i.e. as the driving force for material redistribution decreases. For the current arrangement equilibrium is reached when the normalized centreto-centre distance equals 0.9.
2.5( ¢ ) i = 3,30 * 10 - 3
(d). { = 1 9 4 * 10 -~
o= ~ 2.25o 22-
~
2.15o
~
2.1_
2.05"
, 02 (e). t = 5.23 * 10 -L
(f). t = 1,38
Fig. 14. Sintering of a series of cylinders with two different radii,
i 0.4
i 0.6
i o.8
Norm~ized
i 1.0 time
= 1.2
i 1.4
1.6
Fig. 16. Variation of the distance between the centres of the two large cylinders either side of a small cylinder.
1406
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
In each of the simulations examined above the total change in volume was less than 0.06%.
4. 5.
7. CONCLUSIONS
6.
In general, the chemical potential at a junction between a grain-boundary and a free surface has to be determined from the coupled diffusion process. As a result, the local profile of the free surface is also determined by the coupled diffusion process; and we refer to this profile as the steady state local profile. Any departure from the steady state profile is unstable. When developing numerical techniques for the coupled process it is necessary to develop procedures that allow for non-steady state local profiles, which can arise due to the various approximations and time integration schemes employed in the analysis. The correct chemical potential at the junction can be found by allowing the pore tip potential to be an extra u n k n o w n in the discretized equations. Based on this concept, a numerical technique has been developed which is able to analyse the coupled diffusion problem in an arbitrary grain-boundary network with pores of arbitrary shape. Initial numerical examples have demonstrated that this concept provides an efficient and effective way of analysing practical diffusion problems.
7.
Acknowledgements--Financial support from the SERC and
the European Union are gratefully acknowledged.
8. 9. 10. I 1. 12. 13. 14. 15. 16. 17. 18. 19. 20. 21. 22.
REFERENCES
1. M. F. Ashby, HIP 6.0 Background reading, Engineering Department, Cambridge University (1991). 2. A. C. F. Cocks and A. A. Searle, Acta metall, mater. 38, 2493 (1990), 3. C. Herring, The Physics of Powder Metallurgy (edited
23. 24. 25. 26.
by W. E. Kingston). McGraw-Hill, New York (1951). D. S. Wilkinson and V. Vitek, Acta metall. 30, 1723 (1982). T. J. Chuang and J. R. Rice, Acta metall. 21, 1625 (1973). T. J. Chuang, K. I. Kagawa, J. R. Rice and L. B. Sills, Acta metall. 27, 265 (1979). G. M. Pharr and W. D. Nix, Acta metall. 27, 1615 (1979). V. T. Birth and R. Uzan, Surf. Sci. 179, 540 (1987). Y. Takahashi, F. Ueno and K. Nishiguchi, Acta metall. 36, 3007 (1988). M. F. Ashby, Acta metall. 22, 275 (1974). F. B. Swinkels and M. F. Ashby, Acta metall. 29, 259 (1981). J. Svoboda and H. Riedel, New solutions describing the formation of interparticle necks in solid-state sintering. Submitted for publication (1994). B. Bouvard and R. M. McMeeking, The deformation of interparticle necks by diffusion controlled creep. Submitted for publication (1994). A. C. F. Cocks, Applied Solid Mechanics 3 (edited by I. M. Allison and C. Ruiz). Elsevier, North Holland (1989). J. Pan and A. C. F. Cocks, Comput. Mater. Sci. 1, 95 (1993). A. C. F. Cocks, Acta metall, mater. 42, 2191 (1994). A. Needleman and J. R. Rice, Acta metall. 28, 1315 (1980). J. Svoboda and H. Riedel, A two-dimensional sintering model for coupled grain-boundary and surface diffusion. Submitted for publication (1994). N. Charles and G. Weaver, Materials Principles and Practice, 185pp. The Open University, ButterworthHeinemann (1990). Z-Z. Du and A. C. F. Cocks, Acta metall, mater. 40, 1969 (1992). M. Hillert, Acta metall. 13, 227 (1965). P. G. Shewmon, Trans. metall. Soc. A.I.M.E. 230, 1134 (1964). R. J. Brook, J. Am. Ceram. Soc. 52, 56 (1969). Miyamoto and T. Miyashita, J. Am. Ceram. Soe. 73, 74 (1990). B. J. Kellett and F. F. Lange, J. Am. Ceram. Soc. 72, 725 (1989). F. F. Lange and B. J. Kellett, J. Am. Ceram. Soc. 72, 735 (1989).
Pergamon
0956-7151(94)00365-3
Copyright ~ 1995ElsevierScienceLtd Printed in Great Britain.All rights reserved 0956-7151/95 $9.50+ 0.00
A NUMERICAL TECHNIQUE FOR THE ANALYSIS OF COUPLED SURFACE A N D G R A I N - B O U N D A R Y DIFFUSION J. PANt and A. C. F. COCKS:~ Engineering Department, Cambridge University, Trumpington Street, Cambridge CB2 1PZ, England (Received 20 July 1994)
Abstract--The continuity conditions across a junction between a grain-boundary and a free surface are discussed for the coupled surface and grain-boundary diffusion problem. An appropriate treatment of the continuity is proposed, which forms the basis of a numerical procedure for an arbitrary network of grains containing pores of arbitrary shape and size. The analysis technique, when combined with a suitable time integration algorithm, is able to simulate physical processes such as powder sintering, diffusional void growth and creep crack propagation etc. which are dominated by the coupled diffusion mechanism under many practical circumstances. Currently only two-dimensional problems are discussed. The basic principles proposed in this paper, however, are also valid for three-dimensional situations. Simple numerical examples of powder sintering are given while the full potential of the numerical technique for the above physical problems will be exploited in a forthcoming paper.
1. INTRODUCTION At elevated temperatures, many physical processes which occur in polycrystalline materials require a coupling between surface and grain-boundary diffusion. Take the sintering of a ceramic powder compact for example: driven by the reduction of the associated total free surface energy, matter flows from gi'ainboundaries into pores, resulting in shrinkage of the pores and therefore the densification of the powder compact. Under many practical circumstances [1], matter is transported to the junctions between grainboundaries and pores by grain-boundary diffusion and is taken away from the pore tips and deposited onto the pore surfaces by surface diffusion. Another example is diffusional void growth in which the above process is reversed, i.e. matter flows from the pores into the grain-boundaries, driven by an external tensile stress state [2]. The solid state diffusion of matter is driven by the gradient of chemical potential of the diffusing species. This gradient of chemical potential is induced by the gradient of the free surface curvature in the case of surface diffusion and by the gradient of stress acting normal to the grain-boundaries in the case of grainboundary diffusion [3]. When a grain-boundary meets a free surface, the two processes are coupled. Because of the complex geometric problems involved, analytical solutions of the governing equations for this coupled diffusion problem, especially for the surface tCurrent address: Department of Mechanical Engineering, University of Surrey, Guildford, Surrey GU2 5XH, England. :~To whom all correspondence should be addressed.
diffusion part, are difficult to obtain. A given physical problem can often be simplified by first examining two extreme situations: (a) the fast surface diffusion limit; and (b) the fast grain-boundary diffusion limit. At the extreme of fast surface diffusion, it is assumed that surface diffusion is much faster than boundary diffusion and matter transported to the pore tip is instantly taken away and distributed along the pore surface (or an abundant source of matter is available at the pore tip to be taken away by grainboundary diffusion). As a result, an equilibrium pore shape is maintained and the pores simply serve as boundary conditions for the grain-boundary diffusion problem in the form of capillarity stresses. For geometrically symmetric problems, analytical solutions for the grain-boundary part of the problem can be easily obtained. This approach has been adopted by many researchers to model void growth [2], creep crack propagation [4] and powder sintering [ll. At the extreme of fast grain-boundary diffusion, more attention is paid to the evolution of pore profile while the grain-boundaries only act as sources or sinks of matter and as a boundary condition for the chemical potential which links the pore tip curvature with the remote stress. Chuang and Rice [5] obtained analytical solutions at this extreme under a steady state assumption for a one pore and one boundary system. Chuang et al. [6] also investigated the coupled problem and developed a class of self-similar solutions for the one pore one boundary system. For more general situations, analytical solutions are often too difficult to obtain and numerical methods are often adopted. Pharr and Nix [7] studied cavity
1395
1396
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
growth controlled by surface diffusion using a finite difference scheme. Similar schemes were used by Binh and Uzan [8], and Takahashi and Nishiguchi [9] etc. for problems of solid drop formation at the tip of a solid wedge and void growth controlled by surface diffusion. A simplified treatment of surface diffusion during sintering has been carried out by Ashby [10] and Swinkels and Ashby [11]. Recently Svoboda and Riedel [12] developed similarity solutions at both extremes for small scale diffusion. They also suggested a variational solution for the full coupled situation under the small scale diffusion assumption. Bouvard and McMeeking [13] and Svoboda and Riedel [12] used a scheme of combining the analytical solution for grain-boundary diffusion with the finite difference procedure for surface diffusion for a uniform array of particles. This scheme has provided useful results for the investigated sintering problems but it can only be applied to situations where analytical solutions for the grain-boundary part of the problem are available and where the pore shape is highly symmetric. For more general situations, two issues have to be dealt with in order to solve the coupled diffusion problem: (a) solution of the grain-boundary part of the problem for any arbitrary boundary network and (b) the continuity conditions at the junctions between pore surfaces and grain-boundaries. The authors have developed a finite element type of numerical technique which can solve grain-boundary diffusion problems for any arbitrary boundary network with pores of equilibrium shape [14, 15]. A natural extension of this technique is to combine it with the finite difference procedure for the surface diffusion part of the problem so that the coupled problem can be solved. However, the continuity conditions have not been fully examined in the literature. In this paper, we first discuss the continuity conditions, clarify some areas of confusion that often arise in the literature and then propose a numerical scheme which can solve the coupled diffusion problems for any arbitrary grain-boundary network with pores of arbitrary shape. Simple examples of powder sintering are given while the full potential of the numerical technique will be exploited in a forthcoming paper.
¢
¢
E~
E~
Fig. 1. Polycrystalline material with pores. According to Fick's law, the diffusive flux j, defined as volume of matter flowing across unit area perpendicular to the flux direction per unit time, is linearly dependent on the gradient of the chemical potential # of the diffusing species Dt$ t3/a
j
k T ~s
where D is diffusivity, 6 is the thickness of the layer through which the material diffuses, k is Boltzman's constant, T is the absolute temperature and s is a local coordinate along the diffusion path.
2.1. Grain-boundary diffusion Along a grain-boundary, the gradient of the atomic chemical potential is induced by the gradient of stress o acting normal to the boundary [3] # = - f~o
A schematic example of the coupled diffusion problem is shown in Fig. 1. In the current paper we confine our attention to two-dimensional problems. Grain-boundaries are assumed to be straight and grains are assumed to be rigid. Grains either side of a grain-boundary are free to slide with respect to each other. Their relative motion in the direction normal to the grain-boundary, however, has to be accommodated by grain-boundary diffusion.
(2.2)
where f~ is the atomic volume, therefore
Db6bf] ~0" k T ~3s
A- - -
(2.3)
in which D b is the grain-boundary diffusivity and 6b is the grain-boundary thickness. As matter is deposited onto (or removed from) a grain-boundary, a normal velocity Vb between grains either side of the boundary results. Matter conservation requires that
eJb
~ s + Vb= O. 2. THE GOVERNING EQUATIONS
(2.1)
(2.4)
At each junction where a number, n, of grain-boundaries meet, the flux Jb must satisfy the condition that ~ Jb,i " ti = 0
(2.5)
i=l
where ti is a unit vector in the direction of the ith boundary and pointing away from the junction.
2.2. Surface diffusion Along a free surface, the gradient of the atomic chemical potential is induced by the gradient of the
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION free surface curvature [3]
1397
where ?s is the surface free energy and x is the principal curvature of the surface. It is positive for a concave free surface. So we have
Furthermore, the surface tension ),~ and grainboundary tension ~'gbmust be in equilibrium as shown in Fig. 2. This equilibrium requirement leads to a discontinuity in the tangent of free surface at the junction which is often referred to as the dihedral angle 20
D~6~7~ ~ j~ - - kT ~s
cos 0 = - - . 2?.,
# = -- ~?s x
(2.6)
in which D~ is the free surface diffusivity. As matter is deposited onto (or removed from) a free surface a normal velocity v~ of the free surface results. Matter conservation requires that
~j~
0~ + t,~ = 0
~gb
(2.7)
Finally, the chemical potential must be continuous (otherwise there would be an unbounded flux j there). Because of equations (2.2) and (2.6), the continuity of the chemical potential is always expressed as O'tip = ?s/('tip
(2.8)
with v~ taken as positive when matter is added to the surface. The surface velocity v~ due to matter redistribution is imaginary, i.e. it is not a velocity of any material element. It only makes sense when a point-to-point correspondence is defined between the surfaces before and after matter redistribution. Equation (2.8) implies that this point to point correspondence is defined by the normal of the surface before matter redistribution, i.e. v~ represents the 'thickening' rate of the surface as matter is deposited onto the surface.
2.3. Global equilibrium The stresses, a, normal to the grain-boundaries, the free surface tension, ?s, and the applied remote stress must be in equilibrium.
(Jb)tip = -- (J + )tip + (]" s )tip
(2.11)
in the literature. However, (2.10) implies that xtip does not exist. We examine this feature of the continuity requirements further in Section 3.
2.5. Nondimensionalization Let d be a characteristic diffusion distance (the grain-size for example) and O-0 a suitable reference stress (the applied remote stress for example). We can then define a reference rate
Eb--
k~
d3
(2.12)
which has the dimensions of strain-rate, and write the governing equations in terms of dimensionless quantities
= xd, 6 = O-/o-o,f=J/(Ebd2), f = V / ( d b d ) , t ' = ~b t
2.4. Continuity When a grain-boundary meets a free surface, several continuity conditions must be satisfied. First matter must be conserved at the pore tip
(2.10)
(2.13)
with any lengths in the physical description scaled by d. Equations (2.3) and (2.7) then become = ~--~
(2.14)
js = C~z_
(2.15)
Ds3s ?s Cs - - = Of~ Db3b O-od
(2.16)
Ds6s 7s 0 = 0 - ~ and ~s = ~ .
(2.17)
(2.9)
where + and - indicate the two free surfaces joining the grain-boundary and the sign convention for these local fluxes are defined in Fig. 2.
and OS
where
i.e.
Fig. 2. A junction between a grain-boundary and a free surface. AM 43/4~.5
Other equations keep their original form in terms of the corresponding nondimensionalized variables. Three groups of non-dimensionalized material constants emerge from the above analysis: O, which represents the relative importance of the two diffusion processes; ?s/?gb, which determines the dihedral angle at the pore tips and 7s/o-0d, which represents the relative significance of the applied stress and the surface tension.
1398
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION 3. THE CONTINUITY CONDITIONS
For the purpose of discussion we consider the simple problem of Fig. 3(a) to examine the coupling between grain-boundary and surface diffusion. Equilibrium shaped pores with spherical caps of radius, R, are periodically distributed along an infinitively long grain-boundary, which is subjected to a remote normal stress, E~. Because of symmetry, only the part of the problem shown in Fig. 3(b) needs to be considered. At the pore tip, continuity of potential, equation (2.11), requires that O'ti p = ~)s/R and a potential gradient is set up to drive material away from the tip at a rate 1 fs ~ = 2~s(~ + sin 0 ) - - T= - 2 l ' ( 1 - ~ )
(3.1)
where the reference stress a0 has been identified with the applied stress E® and the characteristic diffusion distance d with the half ligament length a. Our assumption of an equilibrium shape for the pore, however, means that~ = 0 at the pore tip and we are unable to satisfy equation (2.9). This inconsistency indicates that within the framework of steady state diffusion theory employed here we cannot completely define the pore profile. A pore tip profile (or, at least, a chemical potential profile, which need not be the same thing) must therefore develop which allows both equations (2.11) and (2.9) to be satisfied. The chemical potential at the pore tip can therefore only be determined by considering the coupled diffusion process. The chemical potential has to be continuous everywhere otherwise there would be an infinite flux. As a result, the curvature and therefore the local profile of the pore in the region very close to the pore tip are also determined by the coupled diffusion process and can not be arbitrarily prescribed. We refer to this local profile close to the pore tip as the steady state tip profile. Any departure from this profile will be adjusted to the steady state instantaneously and it is this adjustment that kicks off the diffusion process along the pore surface in the above example. E
0
(a)
t
t
0
tz
t
t
In analytical solutions [5, 6], this problem is avoided by directly seeking the steady state profile. The solution then is not designed to satisfy any initial pore profile condition. In numerical procedures, however, each time step has to start from a pore profile determined from integration of all the previous timesteps. This profile is generally not the exact steady state one because of integration errors and node re-arrangement as the pore profile evolves. The numerical solution also has to start from an initial profile, a circular cap for example, which is not necessarily the steady state profile. Finding the correct chemical potential at the pore tip at each time step is then the critical step for the numerical solution of the coupled diffusion problem. This will be further discussed in the next section. When analysing a symmetric coupled diffusion problem using the finite difference method Svoboda and Riedel [12] determined the 'pore tip curvature' by drawing a circle through the pore tip and its neighbouring finite difference node, ensuring that the angle between the two symmetric curves satisfies the dihedral angle relationship of equation (2.10) at the pore tip. The resulting pore tip curvature provides a surface flux that violates matter conservation. This problem was solved by shifting the pore tip coordinate iteratively until matter conservation was satisfied. This procedure found the correct pore tip chemical potential but it is difficult to apply the technique to other situations. If we consider an asymmetric pore profile as shown in Fig. 2, then we have two such circles obtained from either side of a pore tip giving curvatures of opposite sign. The procedures described in the following section were designed to allow this more general situation to be evaluated. 4. A COUPLED FINITE ELEMENT AND FINITE DIFFERENCE METHOD
4.1. Finite element formulation for grain-boundary diffusion For the grain-boundary diffusion part of the problem, Cocks [16] extended an earlier result of Needleman and Rice [17] to establish the following variational principle: among all the kinematically admissible fields .~ [i.e. those satisfying matter conservation requirements (2.4) and (2.5)], the true field [i.e. that also satisfying (2.3) and global equilibrium] makes the functional 1"10 a minimum, where ]'I 1 (" -2 -T 0 = ~ | Jb dA q"~.~tipm'Jb
~JAb
-- ftiP ~' " ~" dS + ~-~fs ~Tbsin 0 JST
(b)
r
[ .L.
a
Fig. 3. Periodically distributed pores along an infinitely long grain-boundary.
(4.1)
tip
and the summation is performed over all the pore tips, m is a unit vector along a grain-boundary towards the pore tip, ST is the part of the specimen boundary where tractions T are prescribed, ~" is the velocity of ST, )7 is the nondimensionalized surface tension and 0 is the dihedral angle.
1399
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION Cocks [14] introduced the constraint of (2.5) into 1"10and obtained a modified functional
~ l'I=l~10+ ~
(2~m.]b ~
junction\
gb
a
(4.2)
/
where 2 represents a set of Lagrange multipliers. When constructing approximate solutions using (4.2), does not have to automatically satisfy (2.5). Surface diffusion and grain-boundary diffusion are only coupled through the pore tip. For the coupled diffusion problem, the variational principle (4.2) still represents the governing equations for the grainboundary part of the problem except that fft~p is unknown and the variation of)~ must satisfy the extra constraint of equation (2.9). This constraint provides the extra equation required to determine ¢Ttip' A 'grain-boundary element' proposed by the authors [15] is shown in Fig. 4. For a straight grain-boundary, its opening rate Vb can be directly related to the velocities of its two associated grains ~/1, ~7~, and 05~ and t~2, 172,052 ~Tb= [/~ ][00]
(4.3)
[/~]=[--nl --n2 ( £ + f l ) n l n2 - - ( g + ~ ) ]
(4.4)
where
and
j
~ram-boundary "
B -
v....
~
j~
Fig. 5. Discretized free surface. where [(~] = [n I s n2,~ --(0.5s 2 -k-~l)
--nlg--nS(O.5gz +~)l] [/-'~] ~---[/gl91 (~1 /~2 ~TZO)2J~,0]T
(4.8) (4.9)
and J~,0 is the matter flux across the origin of the local co-ordinate system. Substituting the matrix expressions ofj~ and Vbinto (4.2), we have the discrete form of f'I l=l = 1 I~R~ T + ~
,~ ,v
-~ ~tip~"C ["T "-~~s]~) T -- 1~'~ T (4.10) [00] = [ffl vl °51 u2 f2 °52]T.
(4.5)
In (4.4) g is the local co-ordinate along the grainboundary, with its origin located at the mid-point of the boundary, n i, i = 1, 2, is the normal to the grain-boundary, the length ~ is the negative value of the local co-ordinate where a line drawn parallel to n through the centre of the first grain intersects the grain-boundary facet, Fig. 4, and similarly T2 is the negative value of the local co-ordinate where a line drawn parallel to n through the centre of the second grain intersects the grain boundary facet. From equation (2.4), we have Jb =
fo
-- [B][O0] dg +J~,0
(grain 2)
$
16
61']=0
(4.11)
(4.6)
which gives the following set of linear simultaneous equations AX = R (4.12)
(4.7)
where A is assembled from 1~, ~ ' and (~", R is assembled from Fsl] and 1~, and finally X is the global matrix of unknowns consisting of 1~, .~ and :~tip. Since ~ p does not join the variation, (4.11) has not provided enough equations for all the unknowns. The extra equations are provided by the surface diffusion part of the problem.
which can be written as = [(~ 1[01
where 10, I~' and Q" are the global assemblies of [U ], 1~6 and l~g are the global assemblies of [L~0], and/~, ~tip and ~ are the global assemblies of 2, atip, and ~" respectively. 1~, (~', (~" and n are the global stiffness and constraint matrices, with their detailed expressions given in [15]. The true grain centre velocities •, q, 05, the true matter fluxes across the origins of the grain-boundary elementsj% and the Lagrange multipliers ~'introduced at the grain-boundary junctions are those that make the functional l"l stationary, i.e.
4.2. Finite differenceformulation for surface diffusion
(grain 1) Fig. 4. Grain-boundary element.
For surface diffusion, we adopt a finite difference scheme. The pore surface is discretized into knot points as shown in Fig. 5. We shall refer to these points defining the pore profile as P-points. The curvature x and the centre of curvature of the surface at a P-point are determined by constructing a circle
1400
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
through the P-point and its two neighbouring Ppoints. For a given pore profile, expressed in terms of the P-point co-ordinates, the curvature is known at all P-points except at the pore tips where curvature is not defined. A point bisecting two neighbouring P-points is referred to as a B-point. Matter flux f along the surface across any B-point can be determined according to (2.15) (4.13)
Js,ib = Cs (/~ib,+P -- Kib,-P) ASib
It proves convenient to define an imaginary pore tip curvature through (2.6) K'tip --
(4.16)
which has nothing to do with the pore tip curvature prescribed by the initial geometry, because #tip is unknown. The matter fluxes across the two B-points next to the pore tip are given by -
(-, /('tip -- •+2
Js'+1=vs
except for those points next to the pore tips, where gib.+V and x~b,-e are curvatures of the two P-points neighbouring the ith B-point, and Asib represents the distances between these two neighbouring P-points. The normal velocity of the pore surface at any P-point can then be determined according to (2.8)
/~tip
~')Ts
-
AS+l---- ' Js'-l=
C K'-2 -- -Ktip
s
AS_l
in which g±2 are known. The velocities of the two P-points next to the pore tip are given by
V'S,+2~----2 J~,+t-A,+2 AS+I "[- As+2
~,,-2 = - 2 AL,-2 ~ _ ~-L,-1 T~_ 2
17s,ip=
Js,ip,+O -- Js,ip,-B Asip
in which ~,±2 are known. Combining equations (2.11), (4.7), (4.17), (4.18) with (4.15), we obtain [C]fip[U] + Q1 t~tip + Q2 = 0
With reference to the pore tip shown in Fig. 5, we refer to ~,,ip as the matter flux into the pore tip from the grain boundary, f,+~ and f , - i the matter fluxes across the two B-points next to the pore tip, and ~7~,~p,i v~+2 and ~,-2 the velocities of the pore tip and the two P-points next to the pore tip. The dihedral angle is enforced between the line linking the pore tip and its neighbouring P-point and the grain-boundary. We further assume that matter deposits uniformly on the two elements adjacent to the pore tip. As a result there is no rotation of these elements or the boundary that meets the tip, although this boundary might be required to migrate. We examine this assumption further in Section 4.3.2. Matter conservation requires that the amount of material deposited on to the surface between the two B-points is equal to the amount flowing into this region and we have As+ l
As_ 1
Jb,tip + f.+ 1 -- f,_ 1 = vs,+2 -----~ + ~Ts,_2 ---~
(4.19)
where [(~ ]tip is the matrix of equation (4.8) at the pore tip, [U] is a matrix which contains the degrees of freedom associated with the grain-boundary [equation (4.9)] and QI and Q2 are known coefficients calculated from g_+2,J~,+2, AS±I, As±2 and Cs
QI = C ,
+ A-~_I+ AS+l + As+2
1
-~ As l+As_2 Q2 = - - C s Y K+2 --+
LAs+l
4.3. Coupling between surface and boundary diffusion
_
(4.18)
(4.14)
except for the pore tips and those P-points next to the pore tips. Here f,ip,+B and f,ip.-B are fluxes across the two B-points neighbouring the ith P-point, and Asip is the mean of the distances between the ith P-point and its two neighbouring P-points. The direction of ~ is taken as that along the direction of the normal of a P-point defined by the P-point and its centre of curvature. From the resulting values of q~ and their directions, the pore profile can be updated by a suitable time integration scheme. However, es at the pore tips and their immediate neighbouring P-points cannot be determined without considering the boundary diffusion part of the problem.
-
(4.17)
l K+2
t~-2
A--~_i+ As+l+As+2
re_2
As_~-AS_~] As,
-
AS+~
-~ I~S_1 -t-'-* IxS-2 Js-Z--As+l + "
~
] _l
By appending equation (4.19) to equation (4.12), we obtain the appropriate number of equations to determine all the degrees of freedom and the values of the chemical potential at the pore tips. 4.3.1. Application of the continuity condition. As an example, we consider the problem of Fig. 2 discussed in Section 3. Consider, as before, the situation where the initial pore profile is a circular cap. If it is uniformly divided into elements of the same length, As, we find, from equations (3.1) to (3.4) and (4.19) /('tip --
+
~
-- sin 0
(4.20)
1+-O (4.15)
where As± l are the two distances between the pore tip and its two neighbouring P-points.
where the half length of the grain-boundary between the two neighbouring pores a has again been taken as the characteristic diffusion distance d and E~ has been interpreted as the reference stress tr0. /~ is the
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION normalized radius of the circular cap a n d / the half distance between the two neighbouring pore centres. The flux across the B-point next to the pore tip is then given by, equation (4.17)
('
J s = C ~ ss ~ - Xtip
) As 1+-O
(4.21)
It proves instructive to examine these expressions in certain limits. We would expect the discrete description here to reduce to a continuum description in the limit As ---, 0. In this limit equations (4.20) and (4.21) become l
f~ip = ~
(4.22)
and
js (l+sin0) r -
"L
) 423,
i.e. the tip curvature is the same as that for the rest of the pore, and the fluxf towards the pore tip (along the two surfaces) is equal to the flux~ away from the tip along the grain-boundary [compare equations (4.23) and (3.1)]. It is also interesting to note that the surface flux at this initial stage is independent of the magnitude of the surface diffusivity. We can gain further insight by allowing As to increase slightly. Then
~,i~=~ + g [ ~ - s i n 0 - ~
,] +O
,
l
.(4.24)
This equation demonstrates that flow of material away from the pore tip by boundary diffusion generates a small perturbation in the pore tip profile which drives the diffusive flux of material to the tip along the pore surface. As the void shrinks, this perturbation propagates along the surface, allowing for internal redistribution of material. If As~® is small then equations (4.24) and (4.21) reduce to (4.22) and (4.23). If ® is large, i.e. if the surface diffusivity is much greater than the boundary diffusivity, then the continuum response can be well approximated using a coarse mesh. In this limit only small perturbations in the pore profile will be set up to drive the diffusive flux of material around the surface throughout the diffusion process. For moderate values of O, a suitable initial node spacing can be obtained by requiring that As/® is small. We now turn to the limit of ® ~ 0 for finite As. Equation (4.20) then becomes
This is precisely the tip boundary condition employed by Pharr and Nix [7] for their computation in the limit of fast boundary diffusion, and simply reflects the fact that in this limit the potential gradient along the boundary is small, with the normal stress a n approximately equal to the net section stress along the entire length of boundary. In this limit equation (4.21) becomes ~=~s[)Ts(l+sin0)-/"]
Vs,tip" n+ = vs,+2 • n+
(4.27)
vs.,ip • n_ = vs,_2 • n_
(4.28)
where n+ and n_ are the normals of the two straight lines linking the pore tip and its two neighbouring P-points, see Fig. 6. The dihedral angle comes into play through the dot product between vs.tip and n_+. It should be emphasized that not only the magnitude, but also the direction of the pore tip velocity, is determined from the coupled diffusion process. If a pore profile is not symmetric about a grain-boundary, then the above equations give a pore tip velocity which is not along the grain-boundary direction. Since, in the derivation of the above equations, we have not allowed the pore tip profile to rotate, the grain-boundary must migrate to retain its orientation with respect to the pore tip. This assumption is valid provided the characteristic rate of grain-boundary migration is much faster than the characteristic rates of the two diffusion processes considered here. If this is not true then the finite element formulation for
V,.e
/3*,÷~
i.e. #tip = - - ( E ~ ! - y~sin 0 ) fL
(4.25)
(4.26)
i.e. the surface flux now depends on the surface diffusivity and the flux at the first B-point depends on the length of the first element. This reflects the fact that there is significant rearrangement of material next to the pore tip and the amount of material flowing to the tip decreases with increasing distance from the tip. 4.3.2. Pore tip velocity. As soon as grip is known, the velocities of the two P-points next to a pore tip can be determined through equations (4.17) and (4.18). Because the dihedral angle is enforced between the pore tip and these two P-points, the pore tip velocity is related to the velocities of these two P-points in the following way
/g~ip= - - sin 0 7~
1401
Fig. 6. Pore tip velocity.
1402
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
grain-boundary diffusion has to be extended to allow for curved grain-boundaries and the coupling conditions (4.15), (4.27) and (4.28) have to be further generalized. This is beyond the scope of the current paper, although there is no difficulty in principle in performing such an extension and generalization within the framework of the current numerical technique. 5. UPDATING OF THE PORE PROFILE
Equations (4.12) and (4.28) can be solved by any standard routine for linear simultaneous equations. The solution provides velocities of each individual grain and each P-point of the pore surface. The grain-boundary network and the pore profile can then be updated according to these velocities. The updating procedures for the grain-boundary network has been discussed in [15], therefore only the procedures for updating the pore profile are discussed here. Euler's method is used to perform the time integration of the pore profile
xi(t -q- At ) = xi(t ) + VsnxA i yi(t + A t ) = yi(t ) + VsnyAt i i
(5.1)
where x i, yi are P-point co-ordinates, vis the P-point velocities, n i, n~ the P-point normals and At the time step. The time step At is determined from the maximum surface velocity vmax experienced by the Ppoints at the current time step, such that the displacement increment due to At and v max is kept constant for all time step increments. If we define the mean distance between all the neighbouring P-points of the initial pore profile as dmean, then this constant displacement increment can be set to a prescribed percentage, EI, of d. . . . 13max
At = Etdmea---~'
(5.2)
During the numerical simulation, the pore profile experiences large scale evolution. The finite difference mesh has to be modified accordingly. A minimum critical distance between any two neighbouring Ppoints is set in terms of a prescribed percentage, Emin, of dm.... a P-point will be deleted if the distance between the P-point and its neighbouring P-point is shorter than e~ndme~.. Similarly, a maximum critical
Fig. 7. The sintering of a hexagonal array of circular cylinders.
Fig. 8. Detailed initial profile of a pore tip. distance, £maxdmean, is set so that a B-point will be converted into a P-point if the distance between its two neighbouring P-points becomes longer than this critical distance. Like all numerical procedures, the density of the finite difference mesh has to be determined by experience. The three adjustable parameters, Et, Emi, and Em~, have also to be determined by trial and error. This will be further discussed with regard to the numerical examples presented in the next section.
6. NUMERICAL EXAMPLES In this section we use the procedure described in the previous two sections to analyse four simple problems that relate to the sintering of fine grained materials. First, we consider the pressureless sintering of a close packed array of circular cylinders. This problem has been examined by Svoboda and Riedel [18], and we are able to compare our results directly with theirs. We then examine the situation where the array of cylinders is subjected to a non-hydrostatic stress state, and follow the irregular evolution of the pore shape. Next, we consider the co-sintering of two cylindrical particles of different radius. Finally, we examine the free sintering of a row of cylinders of different sizes and follow the process of dynamic grain-growth. Further examples of the technique are presented in a forthcoming paper.
6.1. Sintering of a close packed array of cylinders As our first example, we use the numerical technique established in this paper to simulate the sintering of a hexagonal array of circular cylinders. The representative element used in the simulation is shown in Fig. 7. For the loading cases considered in this paper, a smaller representative element can be identified. The representative element of Fig. 7 is used deliberately to check the symmetry of the numerical results. All the lengths in the example are scaled by the initial radius of the cylinder. Each free surface between the two contacts are uniformly divided into 60 nodes. The pore tips are given an initial profile as shown in Fig. 8. First, we consider the case of pressureless sintering. Figure 9 shows the evolution of the pore shape and Fig. 10 shows the variation of the relative density and 'radius' of the pore tip as a function of the normalized time. The corresponding stages between Figs 9 and 10 are marked by (a), (b), (c) . . . and so on. The nondimensionalized material parameters used in the simulation are Ds6s/Db6b= 1, )'gb/(27,)=0.5 and ys/(dtro) = 0.01. After several test runs of the program, suitable parameters controlling the timestep
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
(al (b)
{c)
(a)
(e)
(£)
Fig. 9. Evolution of pore shape during pressureless sintering. length and node elimination were determined: et---0.02 and Emin = 0.05. The computer simulated pore shape evolution follows the pattern observed experimentally [19]. The imaginary pore tip radius calculated from the pore tip chemical potential initially increases as a result of material depositing onto the pore tip; then, after the pore has reached its equilibrium shape, it decreases as a result of pore shrinkage. There is a certain degree of noise in the pore tip chemical potential due to its high sensitivity to the detailed local profile. This noise was smeared out when creating Fig. 10. Matter conservation requires that the total area (excluding the pores) should remain constant. The actual change of this area 1.0. 0.98
Relative density ....... "Poretip nadius" J • • Results of Svoboda & R i e d e l /
during the simulation can therefore be taken as a measure of computing accuracy. Over the range of density variation shown by Fig. 10, about 0 . 2 o of the total area was lost. Also shown in Fig. 10 are the results obtained by Svoboda and Riedel [18]. Although different methods were used to deal with the pore-tip continuity in these two analyses, the two simulations agree almost exactly with each other. The method used by Svoboda and Riedel is, however, unable to deal with more general situations like non-hydrostatic loading because it assumes the development of a symmetric pore profile. The method developed in this paper was used to simulate sintering under a remote stress state Y,x x = - 0 2 1 and : ~ y y = - 0 . 1 2 for the situation D s r s / D b ~ b = 5 with 7gb/(2?s) = 0.5 and ?J(da0) = 0.01 as before. Figure 11 shows the simulated densification curve with the evolution of pore shape. This macroscopic stress field produces a mean compressive stress across the angled boundaries of Fig. 7 and a mean tensile stress across the vertical boundary of the repeating cell. The effect of the compressive stress across the angled boundaries is to accelerate the flow of material into the pore from these boundaries, while flow of material into the pore from the vertical boundaries is retarded by the microscopic tensile stress. As a result, in the early stage of densification, the pore develops an anisotropic shape as it shrinks more perpendicular to the direction of maximum compressive stress. In the later stages, material can readily flow around the small pore, allowing it to adopt an equiaxed equilibrium profile. Overall, the superposition of this net hydrostatic stress state accelerates the rate of densification (compare Figs 10 and 11), as one would expect. Reducing the magnitude of Y,xx to zero, increases the magnitude of the tensile stress across the vertical boundary to a level which now causes material to flow away from the pore and to plate out on this boundary. The grains either side of the boundary move apart and the contact patch shrinks in size.
098
f' .0.8
1403
13
0.97 ~ .~ 0.96
"~
J 0.96
t"
.0.6
"'"..
~2
"~
0.94
_0.4 o.94 _
~ 0.93
• ..........
?. 0,92_ b: £
...... ""..
02 ~ "''.
0.92 0.91
I
I
I
t
I
r
5
10
15
20
25
30
x 10 ~ N o r m a l i z e d time Normalized time
Fig. 10. Variations of relative density and 'pore-tip' during pressureless sintering.
Fig. 11. Variation of relative density and evolution of pore shape under a non-hydrostatic stress state of :~xx= -0.01 and ~yy= -0.12.
1404
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
Eventually the grains separate as they are wedged apart and a grain switching event occurs. The implication of these results for developing the constitutive models of sintering is that non-hydrostatic stress states can trigger anisotropic behaviour of the material in the early stage of densification, but in the later stages surface diffusion is sufficiently rapid to redistribute material around the surface of the pores, allowing them to adopt equilibrium profiles and the material response to become isotropic.
6.2. Sintering of cylinders of two different sizes Most commercial powders consist of particles with a range of sizes. This difference in particle size can promote grain-growth as a compact densities during sintering. In the later stages of densification (stage II), when the pores are isolated, the grain-growth kinetics is determined by the pinning effect and rate of migration of the pores. Constitutive models for this process have been described by Du and Cocks [20], based on the grain-growth studies of Hillert [21], Shewmom [22], Brook [23] and Ashby [1]. These models predict a much faster rate of grain-growth during pressureless sintering than HIPing; which is consistent with experimental observations [24]. In the early stages of densification (stage I), however, when
(a). [ = 0
(c). ~ = 6.80 • 10 -3
(e), f = 2.30* I0-I
(b). [ = 4.31 * 10 -4
(d). [ = 1 . 3 9 , 1 0 - 2
(f). ~ =
5.63
Fig. 12. Sintering of two cylinders of different radii.
J 0'12"1"~~ . . . . .
L
0.1
~
of ~
cyfed~
-3.26
/]
"~-,
0.06
-3.24
-3.22 ~.,
-3.2 , %',
I~ ~ ~
-3.~
5
10
t5
20
25
5O
55
4O
45
5O
55
6O
65
Normolized time
Fig. 13. Variation of cross-sectional areas of the two cylinders vs non-dimensionalized time. distinct necks exist between the particles a correlation between grain-growth and densification is generally observed [25] which is not captured in existing models [20]. Kellett and Lange [25, 26] have obtained expressions for the free energy of simple arrays of irregular particles and demonstrated how material redistribution, driven by changes in this free energy, can result in both densification and grain-growth, but full kinetic equations have not been developed for this process. As a preliminary study, to aid the development of appropriate kinetic equations for the coupling between grain-growth and densification in stage I, we examine two simple representative situations: first we simulate the sintering of two cylinders with a radius ratio of 5, Fig. 12(a); then we consider the situation where these two cylinders represent a repeating cell within a row of cylinders which is allowed to sinter. An initial neck size of 0.1 times the radius of the small particle was used in all the simulations. Other material constants and parameters used w e r e : Ds~s/Db~ b = 1, ygt,/(2ys)= 0.5, 7J(dtro)= 0.0l, Et = 0 . 0 2 , £min ~--" 0 . 0 5 and Emax= 1.2. Figure 12 shows the evolution of geometry for the two sintering particles. There is approximately an order of magnitude difference in time between the different stages of Fig. 12. Initially, the diffusion distance is short, as material flows from the neck and deposits in the crevice between the two particles, and the profile evolves at a fast rate. Material transported away from the boundary between the two particles is deposited non-uniformly on the surfaces of the particles, with most material depositing on the larger particle (it is possible to view this process in an alternative way: one can assume that the material removed from the boundary region deposits uniformly onto the two particles, and superimposed on this there is a net flow of material from the smaller particle to the larger particle, driven by the difference in the curvature between the two particles). This non-uniformity of deposition results in a pore-tip velocity in the direction of the smaller particle. As shown in Fig. 13, there is a gradual increase in the volume of the larger particle and a corresponding
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION decrease in the volume of the smaller particle. Two factors contribute to this: the net flow of material between the two particles; and grain-boundary migration, which results from the pore-tip motion. The grain-boundary gradually moves towards the top of the small particles, Fig. 12(c), and eventually meets the free surface, Fig. 12(d), leaving a bump on the large cylindrical particle (numerically, it is not possible to completely follow this process; in the simulation the grain-boundary was deleted when the area of the smaller particle became less than 0.2% of the total area of the two particles). Material contained within this bump subsequently redistributes around the surface of the particle, Fig. 12(d-f). As it does so the characteristic diffusion distance gradually increases, and the process slows down as the particle asymptotes towards a circular equilibrium profile. Equivalent results for the row of cylinders are shown in Figs 14 and 15. Now, the material from the smaller particle is transferred to the larger particles either side of it. As it does so, the two boundaries between the particles move towards each other, Fig. 15(b,c), Until they meet. At this point we could assume that the boundaries annihilate each other, but the larger particles will generally have different crystallographic orientations and it is more appropriate to assume that the boundaries combine to form a new grain-boundary. As before we cannot completely follow this process, and in the simulations the bound-
(a) t = 0
(b). ( = 8 , 6 7 * 10 -~
0.12" ~\
.....
0.1.
~
1405 f3.26
A r ~ of =rrd
3.24
3.2
~'
338 o
316
5
10
15
20
25
50
.10" Normolized
time
Fig. 15. Variation of cross-sectional areas of the large and small cylinders as a function of non-dimensionalized time before the small particle disappears. aries were replaced by a single boundary, located midway between the original boundaries, when the area of the smaller particle became less than 0.2% of the total area of the two particles, as before. Since a small particle is now sandwiched between two larger particles material is transferred from both ends of the particle and the time taken for the two boundaries to meet is roughly half the time required for the boundary to disappear from the small particle of Fig. 12 (compare Figs 13 and 15). The problem is symmetric after the smaller particle disappears and there is no further migration of the grain boundary. Material continues to flow from this boundary and deposit adjacent to the pore tip, causing the neck between the particles to grow, Fig. 14(d-f). The distance between the centres of the two large particles throughout the sintering process is plotted as a function of the dimensionless time in Fig. 16. Initially the rate of approach of the two centres is rapid as the small particle disappears. This rate gradually decreases as the characteristic diffusion distance increases and the configuration approaches equilibrium, i.e. as the driving force for material redistribution decreases. For the current arrangement equilibrium is reached when the normalized centreto-centre distance equals 0.9.
2.5( ¢ ) i = 3,30 * 10 - 3
(d). { = 1 9 4 * 10 -~
o= ~ 2.25o 22-
~
2.15o
~
2.1_
2.05"
, 02 (e). t = 5.23 * 10 -L
(f). t = 1,38
Fig. 14. Sintering of a series of cylinders with two different radii,
i 0.4
i 0.6
i o.8
Norm~ized
i 1.0 time
= 1.2
i 1.4
1.6
Fig. 16. Variation of the distance between the centres of the two large cylinders either side of a small cylinder.
1406
PAN and COCKS: COUPLED SURFACE AND GRAIN BOUNDARY DIFFUSION
In each of the simulations examined above the total change in volume was less than 0.06%.
4. 5.
7. CONCLUSIONS
6.
In general, the chemical potential at a junction between a grain-boundary and a free surface has to be determined from the coupled diffusion process. As a result, the local profile of the free surface is also determined by the coupled diffusion process; and we refer to this profile as the steady state local profile. Any departure from the steady state profile is unstable. When developing numerical techniques for the coupled process it is necessary to develop procedures that allow for non-steady state local profiles, which can arise due to the various approximations and time integration schemes employed in the analysis. The correct chemical potential at the junction can be found by allowing the pore tip potential to be an extra u n k n o w n in the discretized equations. Based on this concept, a numerical technique has been developed which is able to analyse the coupled diffusion problem in an arbitrary grain-boundary network with pores of arbitrary shape. Initial numerical examples have demonstrated that this concept provides an efficient and effective way of analysing practical diffusion problems.
7.
Acknowledgements--Financial support from the SERC and
the European Union are gratefully acknowledged.
8. 9. 10. I 1. 12. 13. 14. 15. 16. 17. 18. 19. 20. 21. 22.
REFERENCES
1. M. F. Ashby, HIP 6.0 Background reading, Engineering Department, Cambridge University (1991). 2. A. C. F. Cocks and A. A. Searle, Acta metall, mater. 38, 2493 (1990), 3. C. Herring, The Physics of Powder Metallurgy (edited
23. 24. 25. 26.
by W. E. Kingston). McGraw-Hill, New York (1951). D. S. Wilkinson and V. Vitek, Acta metall. 30, 1723 (1982). T. J. Chuang and J. R. Rice, Acta metall. 21, 1625 (1973). T. J. Chuang, K. I. Kagawa, J. R. Rice and L. B. Sills, Acta metall. 27, 265 (1979). G. M. Pharr and W. D. Nix, Acta metall. 27, 1615 (1979). V. T. Birth and R. Uzan, Surf. Sci. 179, 540 (1987). Y. Takahashi, F. Ueno and K. Nishiguchi, Acta metall. 36, 3007 (1988). M. F. Ashby, Acta metall. 22, 275 (1974). F. B. Swinkels and M. F. Ashby, Acta metall. 29, 259 (1981). J. Svoboda and H. Riedel, New solutions describing the formation of interparticle necks in solid-state sintering. Submitted for publication (1994). B. Bouvard and R. M. McMeeking, The deformation of interparticle necks by diffusion controlled creep. Submitted for publication (1994). A. C. F. Cocks, Applied Solid Mechanics 3 (edited by I. M. Allison and C. Ruiz). Elsevier, North Holland (1989). J. Pan and A. C. F. Cocks, Comput. Mater. Sci. 1, 95 (1993). A. C. F. Cocks, Acta metall, mater. 42, 2191 (1994). A. Needleman and J. R. Rice, Acta metall. 28, 1315 (1980). J. Svoboda and H. Riedel, A two-dimensional sintering model for coupled grain-boundary and surface diffusion. Submitted for publication (1994). N. Charles and G. Weaver, Materials Principles and Practice, 185pp. The Open University, ButterworthHeinemann (1990). Z-Z. Du and A. C. F. Cocks, Acta metall, mater. 40, 1969 (1992). M. Hillert, Acta metall. 13, 227 (1965). P. G. Shewmon, Trans. metall. Soc. A.I.M.E. 230, 1134 (1964). R. J. Brook, J. Am. Ceram. Soc. 52, 56 (1969). Miyamoto and T. Miyashita, J. Am. Ceram. Soe. 73, 74 (1990). B. J. Kellett and F. F. Lange, J. Am. Ceram. Soc. 72, 725 (1989). F. F. Lange and B. J. Kellett, J. Am. Ceram. Soc. 72, 735 (1989).