Document Sample

Numerical Solution to the Drift-Diﬀusion Equation of Point Defects under Spherical Symmetry Mohammed Suleiman Hussein Suleiman (suleiman@aims.ac.za) African Institute for Mathematical Sciences (AIMS) Supervised by: Prof. D.T. Britton and Prof. M. Harting University of Cape Town, South Africa 22 May 2008 Submitted in partial fulﬁllment of a postgraduate diploma at AIMS Abstract The diﬀusion of point defects is considered as a basis for reactions in solids which change important physical properties of materials. In this work, after reviewing some essential concepts on the subject, a ﬁrst approximation to the problem of clustering of point defects in cubic crystallines has been made by considering certain spherical symmetries. Then, a numerical solution to the obtained drift-diﬀussion equation is acquired using a ﬁnite diﬀerence technique. The obtained results are discussed within the conceptual framework introduced. Declaration I, the undersigned, hereby declare that the work contained in this essay is my original work, and that any work done by others or by myself previously has been acknowledged and referenced accordingly. Mohammed Suleiman Hussein Suleiman, 22 May 2008 i Contents Abstract 1 Introduction 2 Crystalline Solids 2.1 2.2 2.3 2.4 2.5 2.6 Space Lattice . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Unit and Primitive Cell . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Crystal Systems and Bravais Lattices . . . . . . . . . . . . . . . . . . . . . . . . . . . . Interstices in Crystal Structures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Crystal Binding . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Crystal Defects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . i 1 2 2 2 3 3 4 5 7 7 7 7 9 10 3 Diﬀusion and Diﬀusion Equations 3.1 3.2 3.3 3.4 Fick’s First Law . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Fick’s Second Law . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Forced Diﬀusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . The Fokker-Planck (FP) Equation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 Diﬀusion in Crystalline Solids 4.1 4.2 4.3 4.4 Diﬀusivity of Crystalline Solids . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 D Tensor of Cubic Crystals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 D Tensor of Non-Cubic Crystals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 Mechanisms of Diﬀusion in Crystalline Solids . . . . . . . . . . . . . . . . . . . . . . . 12 14 5 Diﬀusion of Point Defects in Cubic Crystals Under Spherical Symmetry 5.1 5.2 5.3 5.4 Applicability of the FP Equation to the Problem of Point Defects Transport in CubicCrystals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 Spherical Symmetries . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 The FP Equation in Spherical Polar Coordinates . . . . . . . . . . . . . . . . . . . . . . 15 Numerical Solution to the FP Equation . . . . . . . . . . . . . . . . . . . . . . . . . . 16 20 6 Results and Discussion 6.1 Numerical Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 ii 6.1.1 6.1.2 6.1.3 6.2 Results for True Diﬀusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 Results for Diﬀusion with an Outward Driving Force . . . . . . . . . . . . . . . . 21 Results for Diﬀusion with an Inward Driving Force . . . . . . . . . . . . . . . . . 21 Results Justiﬁcation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 23 25 7 Conclusion References iii 1. Introduction The subject of diﬀusion in solids has been an active area of research for a long time. Physical properties of materials are strongly related to their internal structure. A real understanding of changes of many physical properties must be based on knowledge of diﬀusion, since most changes in the structure of materials occur by diﬀusion [She63]. The purpose of this work is to introduce some basic notions of diﬀusion in crystalline solids, then to try to apply these essentials to model and solve a problem of early stage clustering of point defects in a crystalline solid bulk. The treatment is very general, that is, no real system is considered. A basic assumption and some simpliﬁcations are made to solve the problem. Thus, no attempt is made at completeness or rigor. Rather it is a ﬁrst step or approximation. Since crystalline solids are to be considered, in our problem, as the mediums of diﬀusion, their internal geometry is discussed in chapter 2. The nature of the diﬀusing particles, point defects, are also included. In chapter 3, the general phenomenon of diﬀusion of particles through a medium is introduced brieﬂy and in a very general sense. That is, in a way that the general derived diﬀerential equation can be equally applied to any physical diﬀusion phenomenon. Diﬀusivity D and mobility η are two basic quantities which contain the physics of any diﬀusion system [BH02]. Their mathematical nature is given in chapter 4. Then this will be followed by a brief discussion about how a particle diﬀuses in a crystal lattice. Chapter 5 starts with introducing the basic assumption of this work: the diﬀusion of point defects can be treated by ignoring the crystalline structure of the medium. Then, to simplify the problem, three spherical symmetries are assumed. According to this assumption and to the simpliﬁcations made, it is shown that the problem becomes a matter of solving the obtained partial diﬀerential equation. Then, because of the geometry of the problem, this obtained equation is rewritten in spherical polar coordinates, and a numerical method of solution is investigated. The obtained results are presented and justiﬁed in chapter 6, followed by a conclusion and recommendations for future work in the last chapter. 1 2. Crystalline Solids According to their internal geometry, solids can be classiﬁed into two main types: crystalline and noncrystalline or amorphous [All01]. Each type is characterized by the size of an ordered region within the material. In crystals, atoms or molecules have a very regular geometric periodicity over extended regions, i.e. over a large number of atoms or molecules. In contrast, amorphous materials have order only over a few atomic or molecular dimensions [Tab00][Nea03]. Structural and chemical defects exist in both structures, but for clarity the discussion below is focused on defects in crystals. 2.1 Space Lattice A crystal is constructed by repetition, in the three dimensions, of identical atoms or groups of atoms which are called bases. Then we deﬁne the lattice, or space lattice, which is a mathematical abstraction, as the set of mathematical points to which bases are attached [Kit05]. 2.2 Unit and Primitive Cell A small parallelepiped of the crystal that can be used to to produce the entire crystal is known as a unit cell. It is not a unique entity. The smallest unit cell is called the primitive cell. In many cases, it is more convenient to use a unit cell that is not primitive[Nea03].A generalized primitive unit cell is shown in ﬁgure 2.1. Figure 2.1: A generalized primitive unit cell. Modiﬁed from: [Nea03] The relationship between the unit primitive cell and the lattice can be given by the three primitive 2 Section 2.3. Crystal Systems and Bravais Lattices Page 3 c γ α β b a Figure 2.2: A generalized unit cell. translation vectors: a1 , a2 and a3 , such that every two lattice points r and r in the crystal can be related to each other by ′ r = r + n1 a1 + n2 a2 + n3 a3 , (2.1) where n1 , n2 and n3 are integers[Nea03][Kit05]. So, the set of points r given by 2.1 for all n1 , n2 and n3 deﬁnes the lattice, and the volume of the unit primitive cell is given by [Kit05]: v = a1 · a2 × a3 . (2.2) ′ ′ 2.3 Crystal Systems and Bravais Lattices Crystal lattices can be carried into themselves by the lattice translation T = n1 a1 + n2 a2 + n3 a3 and by various other symmetry operations. A typical symmetry operation is that of the rotation about an axis that passes through a lattice point. Lattices can be found such that one-, two-, three-, four- and six-fold rotation axes carry the lattice into itself, crresponding to rotation by 2π, 2π/2, 2π/3, 2π/4 and 2π/6 radians and by integral multiples of these rotations. For example, the a cubic lattice has three four-fold rotation axes, four three-fold rotation axes and six two-fold rotation axes [Kit05]. It was found that - by applying equation 2.1 and the other possible symmetry operations - there are only seven ways to draw a unit cell. These possible cells are called crystal systems, and when points are distributed on these systems it produces fourteen possibilities called Bravais lattices. Each cell is characterized by its axial relations, that is the lenghts of its three axes a, b, c and the angles α, β, γ between them as shown in ﬁgure 2.2 [All01]. The Bravais unit cells are given in table 2.1 [Kit05]. The three cubic crystals are the simple cubic (sc) lattice, the body-centered cubic (bcc) and the facecentered cubic (fcc) lattice [Kit05]. Their unit cells are shown in ﬁgure 2.3 and they will be considered in chapter 4. 2.4 Interstices in Crystal Structures By knowing the crystal structure of a material and its lattice dimensions, several charateristics of the crystal can be determined such as the volume density of atoms. Two factors are always used: (i) the packing density is deﬁned as number of atoms per unit cell [Nea03] (ii) the packing fraction is the Section 2.5. Crystal Binding Table 2.1: Crystal Systems Restrictions on Axes’s Lengths Restrictions on Angles a=b=c α=β=γ a=b=c α = β = 90◦ = γ a=b=c α = β = γ = 90◦ a=b=c α = β = γ = 90◦ a=b=c α = β = γ = 90◦ a=b=c α = β = γ = 90◦ a=b=c α = β = 90◦ γ = 120◦ Page 4 System Triclinic Monoclinic Orthorhombic Tetragonal Cubic Trigonal Hexagonal Number of Lattices 1 2 4 2 3 1 1 Figure 2.3: The three cubic cells: (a) simple ubic, (b) body-centered cubic, (c) face-centered cubic. Source: [Nea03] maximum proportion of the available volume -of the unit cell- that can be ﬁlled with hard spheres [Kit05]. The unoccupied spaces in the lattice are called interstices. Interstices have certain geometrical shapes and certain number in a unit cell [All01]. For example, the interstitial sites of an fcc lattice also form an fcc lattice [She63]p.43. The existence of the interstices in the crystal play important roles such as hosting the interstitial atoms in the interstitial solid solutions [All01] and in the diﬀusion of point defects. 2.5 Crystal Binding The cohesive energy of a crystal is deﬁned as the energy required to seperate the crystal components into neutral free atoms at rest, at inﬁnite seperation, with the same electronic conﬁguration [Kit05].The types of bonds in crystals are: 1. Ionic bonds: Oppositely charged ions experience a coloumb attraction and form a bond refered to as ionic bond. If the ions were to get too close, a repulsive force would become dominant, so an equilibrium distance results between these two ions [Nea03]. The term lattice energy is used in the discussion of ionic crystals and is deﬁned as the energy that must be added to the crystal to seperate its component ions into free ions at rest at inﬁnite seperation [Kit05]. Section 2.6. Crystal Defects Page 5 2. Covalent bonds: Here atoms tend to achieve closed-valence energy shells by sharing their outer electrons [Nea03][FLS75]. 3. Metallic bonds: Positive metallic ions might be thought of as being sorrounded by a sea of negative electrons, the solid being held together by the electrostatic forces [Nea03][FLS75]. 4. Van der waals bonds: In some molecules the eﬀective centre of the positive charge of the molecule is not the same as the eﬀective centre of the negative charge. This nonsymmetry in the charge distribution results in a small electric dipole that can interact with a dipole of other molecules [Nea03]. In many materials the bonds between atoms in a crystal are mixed. That is, bonds are not purely ionic, purely covalent or purely metallic. Instead, it is a mixture of two or more bond types [All01]. The diﬀerent strength of some crystals along diﬀerent directions depends on the kind of interatomic binding in the diﬀerent directions [FLS75]. 2.6 Crystal Defects The term defect can be interpreted to mean any equilibrium deviation from the perfect homogeneity and the geometric periodicity of the crystal structure. According to this deﬁnition, ﬂuctuations are not considered. Crystal defects can be classiﬁed, according to their geometrical shape and according to the dimensions that they take, into four main classes [All01]: 1. Point Defects: The word point is used to denote the localized nature of the defects. Theoretically, point defects are zero dimensional, but in reality involve single atoms or single-atomic locations. The commom point defects in crystals are: (i) foreign atoms, or impurities, may be located at normal lattice sites or between normal sites and they are called substitutional impurities and interstitial impurities, respectively [Nea03], (ii) vacancies, that is vacant lattice sites and (iii) self-interstitials, which are extra atoms of the same type of the crystal but not in regular lattice positions. Although point defects are generally introduced by non-equilbrium processes, such as ion implantation [HNR+ 07] or irradiation [PHI91], there is always thermal equilibrium concentrations of point defects present. If the thermal energy of an atom in a crystal exceeds the energy required to remove an atom from its lattice site then that atom can leave its site. As a result, a number of point defects can be formed in the structure: (a) Isolated vacancy: If the atom released and migrates to the surface, or to a grain boundary, it will leave an isolated vacancy inside the crystal. This defect is also known as a Schottky defect, and is illustrated in ﬁgure 2.4(a). The common Schottky defect occurs in ionic crystals, where two vacancies are formed in two diﬀerent lattice sites because of the imigration of a positive ion “cation” and a negative ion “anion”. This type of defect is also called the ion-pair defect, in which the crystal is roughly kept electrostaticaly neutral on a local scale. (b) Isolated interstitial: An atom can diﬀuse in from the surface of the crystal and transfer to an interstitial position, that is a position not normally occupied by an atom, forming an isolated interstitial (ﬁgure 2.4(b)). Section 2.6. Crystal Defects Page 6 (a) (b) (c) Figure 2.4: Two dimensional illustration of the formation of: (a) an isolated vacancy, (b) an isolated interstitial, and (c) a vacancy-interstitial pair. (c) Frenkel defect: If the atom released and migrates only a short distance from the vanacy that it created, the two form an interstitial-vacancy pair called a Frenkel defect, as illustrated in ﬁgure 2.4(c). In higly-packed structures the energy required to introduce an interstitial atom of the element is so large and that make their occurence extremely unlikely. In contrast to that, interstitials are more easily accommodated in structures with low packing fracture [Wal92][All01][Kit05]. 2. Line Defects (Dislocations): The point defects involve single atoms or single-atoms locations. In contrast, dislocation is a deviation from the perfect lattice that involve an entire row of atoms or atoms sites. There are many types of dislocations such as edge and screw dislocation. Dislocations aﬀect most of the physical properties of the crystal [Nea03] [All01]. 3. Planar Defects: The previously described crystals are known as single-crystal materials, that is, they have a high degree of order and regular geometric periodicity throughout the entire volume of the material. In other materials the regular geometric periodicity may extend only over many atomic or molecular dimensions forming single-crystal regions which usually vary in size and in orientation with respect to one another. These single-crystal regions are called grains and they are seperated from one another by surfaces called grain boundaries [Nea03][All01]. 4. Bulk Defects: When vacancies cluster together they form “voids”, and when foreing atoms cluster together they are called “precipitates”. These clusters are called bulk defects. That is, bulk defects are clusters of point defects. These defects aﬀect the physical properties of the perfect crystals, and they can interact with each other forming complex defects [Kit54][Nea03]. 3. Diﬀusion and Diﬀusion Equations Although later the treatment will be restricted speciﬁcally to the diﬀusion of certain speicies in certain mediums, in this chapter, the general phenomenon of diﬀusion of particles through a medium will be brieﬂy introduced. That is, the following discusion can be equally applied to gas diﬀusion and charge transport in semiconductors. Hence, the term particle will be used in a very general sense, which includes even crystal vacancies. 3.1 Fick’s First Law True diﬀusion, or simply diﬀusion, of a certain species in a medium is the outcome of a random motion which, by virtue of a concentration gradient of that species results in a net ﬂux of that species away from the highest concentration. The net diﬀusive ﬂux JD is related to the gradient of the concentration C(x, t) by a famous phenomenological relation called Fick’s ﬁrst law: JD = −D∇C(x, t) , (3.1) where C(x, t) is the concentration of particles at a position x and at time t, (∇C) is the concentration gradient, and D is called the diﬀusion coeﬃcient or diﬀusivity, which is deﬁned as the net ﬂux of particles per unit concentration gradient [She63] [BH02]. If there are other types of gradient in the system, other terms are added to Fick’s ﬂux equation [She63] [BH02]. This will be discussed later. 3.2 Fick’s Second Law Equation 3.1 can only be used in the steady state of ﬂow, that is, when the net ﬂux J is independent of time: J = −D∇C(x, t) = const. (3.2) In the time-dependent case, the rate of change of particles concentration at any point x in the distribution can be related to the net ﬂux J by ∂ C(x, t) = −∇ · J . (3.3) ∂t Equation 3.3 is known as Fick’s second law and can be derived from equation 3.1. It simply states that a net ﬂow out of a volume results in a decreased density inside the volume. This equation appears in a wide variety of physical problems, ranging from a probability density in quantum mechanics to neutrons leakage in a nuclear reactor [Arf85] [She63]. 3.3 Forced Diﬀusion The ﬂow of particles can in fact be due to two causes: one is the mentioned concentration gradient, and the other is the action of a driving force. That is, if we consider a single particle moving in a potential 7 Section 3.3. Forced Diﬀusion ﬁeld V (x); the potential gradient exerts a force F on the particle given by the equation F = −∇V (x) . Page 8 (3.4) It is found empirically that under the inﬂuence of this force, particles move with a certain average terminal elocity vd . This fact is mathematically expressed by the equation vd = ηF , (3.5) where η, is called mobility, and can be deﬁned as the ratio between the average terminal velocity of the particles and the applied force. It is a measure of the response of the particles of the medium to an applied force [She63] [PHI91]. For a concentration distribution of particles C(x, t) under the inﬂuence of a driving force F(x), vd gives rise to a ﬂux: JF = C(x, t)vd = C(x, t)ηF. (3.6) Accordingly, one can deﬁne the forced diﬀusion as a directed motion that is superimposed onto the random motion of the particles. Forced diﬀusion is also refered to as drift [BH02]. Three kinds of driving force that can drive particles, including point defects, in a solid medium are brieﬂy mentioned below. Each of them can be exspressed as a potential gradient [Wip76][All01] [PHI91]: 1. Electric ﬁeld: The motion of charged particles under an electric ﬁeld E is called electroimigration or electrotransport. The electric force Fe is given by Fe = q ∗ E , (3.7) where q ∗ is the eﬀective electric charge of the particle, and E can be expressed in terms of the electrical potential Φ as E = −∇Φ(x) . (3.8) 2. Temperature gradient: A temperature gradient will result in a force given by Fth = − Q∗ ∇T (x) , T (3.9) where Q∗ is called migration enthalpy or heat of transport, and the motion of particles under temperature gradient is called thermoimigration or thermotransport. 3. Stress gradient: Stress gradient will result in a force given by Fσ = −∇U (x) , (3.10) where U is the elastic interaction energy in the stress ﬁeld σ(x). A uniform stress cannot generate a ﬂux (Curie-Neumann principle); thus the driving force must be the stress gradient and its eﬀect must be considered whenever the elastic interaction energy of the particle with the stress ﬁeld is large enough. A typical case is Gorski eﬀect: imigration of atoms in an elastically deformed sample. Section 3.4. The Fokker-Planck (FP) Equation Page 9 3.4 The Fokker-Planck (FP) Equation The total net ﬂux J in equation 3.3 is due to both contributions of JD and JF deﬁned by equations 3.1 and 3.6 respectively. That is J = JD + JF = −D∇C(x, t) + ηC(x, t)F. Substiuting for J from equation 3.11 into equation 3.3 yields ∂ C(x, t) = ∇ · (D∇C(x, t)) − ∇ · (ηC(x, t)F). ∂t (3.12) (3.11) Equation 3.3, or equation 3.12, is a continuity equation, that is, it is valid only for species which obey a conservation law [PHI91]. If one deals with entities which are not conserved, then an additional term S(x, t) which is equal to the rate of production or destruction of these entities per unit volume should be added to the right hand side of equation 3.3, or equation 3.12 [She63] [BH02]. Thus ∂ C(x, t) = −∇ · J + S(x, t) , ∂t or (3.13) ∂ C(x, t) = ∇ · (D∇C(x, t)) − ∇ · (ηC(x, t)F) + S(x, t) . (3.14) ∂t Equation 3.14 is the general form of the drift-diﬀusion equation [BH02], and it is known in the literature as a Fokker-Planck (FP) equation1 [GT88]. 1 A Fokker-Planck (FP) equation is a partial diﬀerential equation for the time evolution of the probability density of a particle (or ﬁeld). That is, for a diﬀusing particle, one writes a diﬀerential equation for C(x, t), where C(x, t) is to be deﬁned as the probability that the diﬀusing particle is to be found in the interval x + dx at time t [LMN07]. Equation 3.14 is one of the simplest Fokker-Planck equations [Ris96]. 4. Diﬀusion in Crystalline Solids While in chapter 3 the diﬃsivity D was deﬁned and its units were mentioned, but nothing more about the mathematical nature of D or its dependence on the structure of the meduim was mentioned. This will be given in the ﬁrst sections of this chapter in order to gustify the applicability of the diﬀusion equations introduced in chapter 3 to the problem of the diﬀusion in crystalline solids. Then this will be followed by a brief discussion about how a particle diﬀuses in a lattice. 4.1 Diﬀusivity of Crystalline Solids In equation 3.1, the two vectors J and ∇C are, in the most general case, not parallel, and then D cannot be considered as a constant. Instead, for the general case, it is assumed that a given component of the ﬂux is inﬂuenced by each of the components of the gradient. Thus, in the cartesian coordinates, this situation is described by the equations: Jx = −D11 Jy Jz ∂C ∂C ∂C − D12 − D13 ∂x ∂y ∂z ∂C ∂C ∂C − D22 − D23 = −D21 ∂x ∂y ∂z ∂C ∂C ∂C − D32 − D33 = −D31 dx dy ∂z , , . (4.1) (4.2) (4.3) where the three scalar ﬂuxes Jx , Jy and Jz are parallel to the three cartesian coordinates, and the set of nine numbers designated Dij is called a second-order tensor and is deﬁned by the above equations. The same argument can be made for the mobility η in equation 3.6 [She63][FLS75]. The above system of equations can be rewritten in the matrix form as D11 D12 D13 ∂C/∂x Jx − D21 D22 D23 ∂C/∂y = Jy D31 D32 D33 ∂C/∂z Jz , (4.4) 4.2 D Tensor of Cubic Crystals Let us consider a cube of a material whose lattice is a cubic one. Taking a single crystal and setting up its axes parallel to the cartesian coordinates, and a concentration gradient is established such that ∂C/∂x = α = 0 and ∂C/∂y = ∂C/∂z = 0. Equation 4.1 then gives the component of the ﬂux parallel to the gradient as ∂C . (4.5) Jx = −D11 ∂x If this gradient is removed and replaced by an equal gradient along the y axis, that is, ∂C/∂y = α = 0 and ∂C/∂x = ∂C/∂z = 0, then the component of the ﬂux parallel to the gradient will be Jy = −D22 ∂C ∂y . (4.6) 10 Section 4.2. D Tensor of Cubic Crystals The same argument can be made for Jz and Jz = −D33 ∂C ∂z . Page 11 (4.7) In a cubic lattice crystal x, y and z axes are indistinguishable. And since the gradient were of equal magnitudes in the three cases, the cubic symmetry requires that Jx = Jy = Jz in the three experiments. This implies that D11 α = D22 α = D33 α, that is D11 = D22 = D33 . (4.8) Again, let us establish a gradient such that dC/dy = α = 0 but ∂C/∂x = ∂C/∂z = 0. By measuring the component of the ﬂux in the x axis, one can determine D12 , since equation 4.1 gives Jx = −D12 ∂C = −D12 α . ∂y (4.9) Now, let us imagine that the source and sink which established this gradient are removed and the crystal rotated 180◦ about its x axis, relative to the source and sink. Since the lattice x direction has four-fold rotational symmetry, the original and ﬁnal positions of the lattice will be indistinguishable. If the source and the sink are again applied, the gradient is again established along the y axes. However, since the y axis is ﬁxed in the crystal, the 180◦ has interchanged the +y and the −y axes, so the gradient will be negative of its previous value. This means that now ∂C/∂y = −α, and ∂C/∂x = ∂C/∂z = 0. This gives ∂C = D12 α . (4.10) Jx = −D12 ∂y But by symmetry, Jx from equation 4.9 must equal Jx from equation 4.11, that is D12 α = −D12 α . (4.11) or more simply This requirement of symmetry can be satisﬁed only if D12 = 0. If the gradient had been along the z axis, the same rotation could have been used to show that D13 = 0. Since the y and the z axes also have two-fold rotation symmetry, it follows that all oﬀ-diagonal terms are zero. Accordingly, for cubic crystals D11 0 0 D = 0 D11 0 , (4.12) 0 0 D11 D = D(x). (4.13) The same argument can be made for the mobility η in equation 3.6, that is η = η(x). (4.14) for cubic crystals [She63] [FLS75]. Equations 4.13 and 4.14 simply show that cubic crystals are isotropic. Section 4.3. D Tensor of Non-Cubic Crystals Page 12 4.3 D Tensor of Non-Cubic Crystals and for orthohrombic crystals It is shown in reference [She63] that equation 4.12 is no longer true for non-cubic crystals. For example, for tetragonal and hexagonal lattices D is given by D11 0 0 D = 0 D11 0 , (4.15) 0 0 D33 D11 0 0 D = 0 D22 0 0 0 D33 . (4.16) That is, non-cubic crystals are inisotropic. 4.4 Mechanisms of Diﬀusion in Crystalline Solids In both cases, true and forced diﬀusion, for a particle to diﬀuse in a crystal, it must overcome the potential energy barrier presented by its neighbours. If the barrier height is E and the particle characteristic vibrational frequency is ν, then the propability ρ that during a unit time this particle will have enough thermal energy to pass over the barrier is ρ = νe−E/kT , (4.17) that is, each unit time the particle makes ν jumps, with a probability e−E/kT of overcoming the barrier in each try. ρ is called the jump frequency [Kit05]. Most mechanisms of diﬀusion require the presence of point defects discussed in section 2.6. There are three possible elementary mechanism of diﬀusion of particles in a solid meduim [Kit05][Wal92]: 1. Ring Mechanism: For an atom in a lattice site, it is easy to conceive a mechanism of atomic motion by direct interchange by a rotation about a midpoint (ﬁgure 4.1a), or by rotation of a ring of three or more atoms 1 [PHI91][She63]. 2. Interstitial Mechanism: An atom is jumping from interstitial site (discussed in section 2.4) to interstitial site is a point defect, and is said to diﬀuse by interstitial mechanism. This is illustrated, in two dimensional lattice, in ﬁgure 4.1b. An appreciable local dilation of the lattice must occur before the jump can occur. It is this dilation or distortion which constitutes the barrier to an interstitial atom changing sites. A variety of interstitial mechanisms is possible. In the direct interstitial mechanism, the atom passes without permanently displacing any of the matrix atoms and also there is no dependence between the successive jumps. In contrast, in the indirect interstitial mechanism (also known as interstitialcy mechanism) a relatively large interstitial atom pushes one of its nearest-neighbour atoms into an interstitial position and occupied the lattice site previously occupied by the displaced atom [PHI91] [She63]. The ring mechanism is not known to operate in any metal or alloy, but it has been suggested as a mechanism which would explain some apparent anomalies in D for bcc metals [She63]. 1 Section 4.4. Mechanisms of Diﬀusion in Crystalline Solids Page 13 3. Vacancy Mechanism: If a lattice site is not occupied, a nearest neighbour atom can jump onto this site, and the vacancy will appear on the site that the atom has just vacated, the atom is said to have diﬀused by a vacancy mechanism [PHI91][She63]. This mechanism is illustrated in ﬁgure 4.1c in two dimensional lattice. (a) ( b) ( c ) Figure 4.1: Three basic mechanisms of diﬀusion : (a) Interchange by rotation about a common midpoint. (b) Migration through interstitial sites. (c) Exchange of position with a vacancy. More complex mechanisms are possible, and several elementary mechanisms can be operative simultaneously. And if diﬀerent types of particles participate in the diﬀusion process (e.g. host atoms, vacancies and interstitials), then the succesive jumps of the particles are dependent [PHI91]. 5. Diﬀusion of Point Defects in Cubic Crystals Under Spherical Symmetry In this chapter, the problem of diﬀussion of point defects in cubic crystals will be treated. After introducing some basic assumptions to simplify the problem, the Fokker-Planck equation, equation 3.14, will be rewritten, for such homogenous isotropic solid mediums, in spherical coordinates. Then a numerical method of solution, given an initial distribution C(x, 0) and deﬁned driving ﬁelds F(x), will be investigated. 5.1 Applicability of the FP Equation to the Problem of Point Defects Transport in Cubic-Crystals Fick’s law (equation 3.1), continuity equation (equation 3.3) and the consequent derived equations are valid for continuum mediums and continous C(x, t) [PHI91]. Fick treated the diﬀusion problem by assuming that the medium is continuum, neglecting the discontinuity of the real material structure [PHI91] [She63] and by ignoring the diﬀusion mechanisms [She63]. We are going to make the same assumption since all subsequent developments in the theory of solids have in no way aﬀected the validity of Fick’s approach [She63]. Accordingly, under this assumption, at no point in the following discussion the structural nature of the diﬀusion medium or the diﬀusion mechanism will enter the problem. Moreover, for cubic crystals, equations 4.13 and 4.14 become D = const. , η = const. That is, neither D nor η has spatial dependence, and the medium is said to be homogeneous. (5.1) (5.2) 5.2 Spherical Symmetries To simplify the problem more, we consider a single cubic crystal with a spherical boundary, or a polycrystalline material which has cubic structure, then imagining one of the grains has a spherical boundary. Also a Gaussian initial distribution is assumed, that is C0 = C(r, θ, φ, 0) = N e− for some constant N 1 , and with boundary condition: C(R, θ, φ, 0) = N e− C(R, θ, φ, t) = 0, (R−r0 )2 2σ 2 (r−r0 )2 2σ 2 , (5.3) , t=0 t > 0, (5.4) (5.5) where R is the radius of the bulk or the grain, 0 < r0 < R is the expected value with 0 < r0 < R, and σ is the standard deviation. The prefactor N can be chosen as a normalization factor in order for C(r, θ, φ, 0) to satisfy the deﬁnition of the FP equation (equation 3.14). That is C(r, θ, φ, t) must be redeﬁned as a probability density. 1 14 Section 5.3. The FP Equation in Spherical Polar Coordinates Page 15 Also, the case when there is only true diﬀusion will be considered ﬁrst, then a driving force with spherical symmetry, that is 1 F ∼ ± 2 er (5.6) r will be introduced. 5.3 The FP Equation in Spherical Polar Coordinates Equation 3.14 yields ∂ C = D∇2 C + ∇D · ∇C − C∇η · F − η∇C · F − Cη∇ · F + S ∂t . (5.7) And since we assumed that the medium which we are dealing with is homogenous and isotropic, then both ∇D and ∇η will vanish, since D and η are constants, and equation 5.7 becomes ∂ C = D∇2 C − η(∇C) · F − Cη∇ · F + S ∂t . (5.8) Also, since the bulk or the grain boundary, the driving force and the initial distribution all have spherical symmetry, it is better then to write equation 5.8 in spherical polar coordinates. In this coordinates we have [Arf85]: ∇C(r, θ, φ, t) = er ∇2 C(r, θ, φ, t) = + ∇ · F(r, θ, φ) = 1 ∂ 1 ∂ ∂ C(r, θ, φ, t) + eθ C(r, θ, φ, t) + eφ C(r, θ, φ, t) , ∂r r ∂θ rsinθ ∂φ ∂ ∂ ∂ 1 ∂ (sinθ C(r, θ, φ, t)) (5.9) sinθ (r 2 C(r, θ, φ, t)) + 2 sinθ r ∂r ∂r ∂θ ∂θ 1 ∂2 C(r, θ, φ, t) , sinθ ∂φ2 ∂ ∂ ∂ 1 . sinθ (r 2 Fr (r, θ, φ)) + r (sinθFθ (r, θ, φ)) + r Fφ (r, θ, φ) r 2 sinθ ∂r ∂θ ∂φ Substituting for ∇C(r, θ, φ, t), ∇2 C(r, θ, φ, t) and ∇ · F(r, θ, φ) from 5.9 into equation 5.8 yields ∂ C(r, θ, φ, t) = D ∂t ∂2C ∂2C 2 ∂C 1 ∂2C cotθ ∂C 1 + + 2 2 + 2 + 2 2 ∂r 2 r ∂r r ∂θ r ∂θ r sin θ ∂φ2 1 ∂C 1 ∂C ∂C + eθ + eφ · (er Fr + eθ Fθ + eφ Fφ ) − η er ∂r r ∂θ rsinθ ∂φ ∂Fr 2 1 ∂Fθ cotθ 1 ∂Fφ − Cη + Fr + + Fθ + ∂r r r ∂θ r rsinθ ∂φ + S(r, θ, φ, t) . (5.10) This is the FP equation in the spherical coordinates, which must be solved for a given initial distribution and a well realistic deﬁned driving force F to model the clustering of point defects. Section 5.4. Numerical Solution to the FP Equation Page 16 5.4 Numerical Solution to the FP Equation Except for a few special cases, there is no straightforward analytical solutions to such drift-diﬀusion equation [Bri94] [Ris96]. Generally, however, it is diﬃcult to obtain analytical solutions especially if no seperation of variables is possible, and its numerical solution is a non-trivial problem in numerical analysis [GT88]. The clearest numerical method is to follow the dynamical processes step by step using a ﬁnite diﬀerencing technique [Bri94]. In this technique, ﬁrst we represent C(r, θ, φ, t) by its values at discrete set of points: ri = i△r, θj = j△θ, φk = k△φ, tn = n△t, i = 0, 1, 2, ..., I j = 0, 1, 2, ..., J k = 0, 1, 2, ..., K n = 0, 1, 2, ..., N (5.11) where △r, △θ and △φ are the chosen spatial spacings, and △t is the temporal one[Pea02]. The Explicit (FTCS) Method: The simplest ﬁnite diﬀerence method is to calculate C(r, θ, φ, t) explicitly. That is, at each time step C(r, θ, φ, t) can be calculated explicitly from its values in the previous time step [Pea02]. In this scheme, the numerical form of equation 5.10 can be written as Section 5.4. Numerical Solution to the FP Equation Page 17 C(r, θ, φ, t + △t) − C(r, θ, φ, t) △t C(r + △r, θ, φ, t) − 2C(r, θ, φ, t) + C(r − △r, θ, φ, t) (△r)2 2 C(r + △r, θ, φ, t) − C(r − △r, θ, φ, t) + r 2△r 1 C(r, θ + △θ, φ, t) − 2C(r, θ, φ, t) + C(r, θ − △θ, φ, t) + r2 (△θ)2 cotθ C(r, θ + △θ, φ, t) − C(r, θ − △θ, φ, t) + r2 2△θ C(r, θ, φ + △φ, t) − 2C(r, θ, φ, t) + C(r, θ, φ − △φ, t) 1 + 2 sin2 θ r (△φ)2 C(r + △r, θ, φ, t) − C(r − △r, θ, φ, t) −η Fr (r, θ, φ) 2△r Fθ (r, θ, φ) C(r, θ + △θ, φ, t) − C(r, θ − △θ, φ, t) + r 2△θ Fφ (r, θ, φ) C(r, θ, φ + △φ, t) − C(r, θ, φ − △φ, t) + rsinθ 2△φ Fr (r + △r, θ, φ) − Fr (r − △r, θ, φ) 2 −η C(r, θ, φ, t) + Fr (r, θ, φ) 2△r r Fθ (r, θ + △θ, φ) − Fθ (r, θ − △θ, φ) cotθ 1 + Fθ (r, θ, φ) + r 2△θ r 1 Fφ (r, θ, φ + △φ) − Fφ (r, θ, φ − △φ) + rsinθ 2△φ + S(r, θ, φ, t). (5.12) = D n If we write Cijk for C(ri , θj , φk , tn ), Fr,ijk for Fr (ri , θj , φk ) ... and so on, equation 5.12, in this notations, becomes Section 5.4. Numerical Solution to the FP Equation Page 18 n+1 n Cijk − Cijk △t = + + −η + D n n n C(i+1)jk − 2Cijk + C(i−1)jk (△r)2 n n 2 C(i+1)jk − C(i−1)jk + ri 2△r n n n n n cotθj Ci(j+1)k − Ci(j−1)k 1 Ci(j+1)k − 2Cijk + Ci(j−1)k + 2 2 (△θ)2 2△θ ri ri n n + Cn Cij(k+1) − 2Cijk 1 ij(k−1) 2 ri sin2 θ (△φ)2 2△r + n n Fθ,ijk Ci(j+1)k − Ci(j−1)k ri 2△θ Fr,ijk n n C(i+1)jk − C(i−1)jk n n Fφ,ijk Cij(k+1) − Cij(k−1) ri sinθj 2△φ n −η Cijk + + Fr,(i+1)jk − Fr,(i−1)jk 2Fr,ijk + 2△r ri cotθj Fθ,ijk 1 Fθ,i(j+1)k − Fθ,i(j−1)k + ri 2△θ ri 1 Fφ,ij(k+1) − Fφ,ij(k−1) n + Sijk ri sinθj 2△φ (5.13) This ﬁnite diﬀerence approximation to equation 5.10 is known as FTCS representation (Forward Time Centered Space). It is an easy algorithm to derive, takes little storage in the computer and executes quickly [Pea02]. But this method is said to be unstable if round oﬀ errors or other errors grow rapidly as the computations proceed. The necessity of using very small time steps is the principal fault of this method [ZC97]. Following this scheme, a code has been written in Python programming language to solve equation 5.13 n+1 for Cijk at each time step. The concentration at each spatial point is calculated as follows: 1. Using 5.11 the space and time were conﬁgurated with N = 100, I = 60, J = 50 and K = 50. 2. An N × I × J × K matrix was created to store the computed values of C 2 . 3. The radius R of the bulk and the total time T of the observation were both chosen to be unity for some arbitrary units. 4. r = 0 was avoided, since it represents a singularity in equation 5.10. 5. The values D = 1 × 10−2 and η = 1 × 10−4 were chosen. 6. The initial distribution 5.3 was stored with an expected value at r0 = 0.4R and standard deviation σ = 0.1R. n 7. The Dirichlet boundary condition: CIjk = 0 for n > 0 was introduced. n 8. Neumann boundary conditions were applied for the angular dependence. That implies: Ci0k = n = C n , and C n n n , Cn = Cn Ci1k iJk ij1 ijK = CijK−1 . i(J−1)k , Cij0 9. Radial concentration distributions at some selected time steps were plotted. 2 It could be thought of as N matrices each of order I × J × K. Section 5.4. Numerical Solution to the FP Equation Page 19 The Implicit (BTCS) Method: In this scheme the the spatial derivatives on the right-hand side are evaluated at time step n + 1. That is n+1 n Cijk − Cijk △t = + + −η + D n+1 n+1 n+1 C(i+1)jk − 2Cijk + C(i−1)jk (△r)2 n+1 n+1 2 C(i+1)jk − C(i−1)jk + ri 2△r n+1 n+1 n+1 n+1 n+1 cotθj Ci(j+1)k − Ci(j−1)k 1 Ci(j+1)k − 2Cijk + Ci(j−1)k + 2 2 (△θ)2 2△θ ri ri n+1 n+1 n+1 Cij(k+1) − 2Cijk + Cij(k−1) 1 2 (△φ)2 ri sin2 θ Fr,ijk n+1 n+1 C(i+1)jk − C(i−1)jk 2△r + n+1 n+1 Fθ,ijk Ci(j+1)k − Ci(j−1)k ri 2△θ n+1 n+1 Fφ,ijk Cij(k+1) − Cij(k−1) ri sinθj 2△φ n+1 −η Cijk + + Fr,(i+1)jk − Fr,(i−1)jk 2Fr,ijk + 2△r ri Fθ,i(j+1)k − Fθ,i(j−1)k cotθj Fθ,ijk 1 + ri 2△θ ri 1 Fφ,ij(k+1) − Fφ,ij(k−1) n+1 + Sijk ri sinθj 2△φ (5.14) Such a scheme is called fully implicit or backward time. To solve equation 5.14 one has to solve a set n+1 of simultaneous linear equations at each time step for Cijk [Pea02]. This method does not suﬀer from stability problems [ZC97], but, compared to the explicit scheme, it takes more storage in the computer and more time to execute. 6. Results and Discussion In this chapter, the obtained results are presented, followed by a brief discussion. 6.1 Numerical Results Graphic results for the three cases: F = 0, F = + r12 er and F = − r12 er are shown below. For each obtained ﬁgure the computer runs for about 40 minutes. 6.1.1 Results for True Diﬀusion Figure 6.1 represents some results obtained in the absence of a driving force. The four diﬀerent curves show the radial distribution of point defects at diﬀerent times: (i) at t = 0 the curve represents the given initial Gaussian distribution, (ii) at a later time t = 0.10T the peak’s decrease shows that the diﬀusion starts to take place, (iii) as the time runs, at t = 0.60T and at t = 0.99T , the decrease of the peak and the broadening of the distribution becomes clearer. Figure 6.1: Radial distribution of point defects under true diﬀusion at four diﬀerent times. 20 Section 6.2. Results Justiﬁcation Page 21 6.1.2 Results for Diﬀusion with an Outward Driving Force Figure 6.2 represents some results obtained in the presence of an outward driving force. In addition to the same behaviour shown in ﬁgure 6.1 one can notice that: (i) the values of C at r = 0 are greater than those in 6.1 at the same time steps, (ii) there is a formation of a new peak at r = 0. Figure 6.2: Radial distribution of point defects driven by an outward radial force at four diﬀerent times. 6.1.3 Results for Diﬀusion with an Inward Driving Force Figure 6.3 represents some results obtained in the presence of an inward driving force. In addition to the same behaviour shown in ﬁgure 6.1 one can notice that: (i) the values of C at r = 0 are lower than those in ﬁgure 6.1 at the same time steps, (ii) there is a formation of a new peak at r = 0. 6.2 Results Justiﬁcation In all three cases: (i) The decreasing of the concentration at r = r0 implies that some point defects have transported to somewhere else. (ii) The movement of the peak toward the centre r = 0 and the increase of C values at r = 0 show the stage of accumulation of point defects at the centre of our spherical bulk. (iii) The increase in C values in the region r0 << r < R shows that there is also a diﬀusion in the outward direction but still less than the inward one. (iv) The zero values of C at r = R (at t > 0) satisﬁes the boundary condition 5.5 and guarantee that the boundary works as a sink for the Section 6.2. Results Justiﬁcation Page 22 Figure 6.3: Radial distribution of point defects driven by an inward radial force at four diﬀerent times. diﬀusing point defects. At r = 0 : (i) There is no concentration gradient at any time, and in the very early stage r = 0 is a stable point. Once the peak at r = 0 is formulated, then it is no longer stable, that is, each particle that reaches r = 0 will be reﬂected (in our macro scale) by the negative gradient. At this stage, the higher positive gradient in the left side of the old peak will gain and drives more particles toward the centre than the reﬂected ones. This gives rise to the accumulation at that point. (ii) As the time runs, we expect this peak to keep increasing, while the other two moving peaks1 are decreasing. This will continue untill these two moving peaks diminish and the peak at r = 0 reaches its maximum2 . Then we expect the true diﬀusion to be totally in the outward direction, enhancing the outward force and retarding the inward one. In the presence of the driving forces: (i) One can notice that the true diﬀusion behaviour is always dominant. That is, there is no big diﬀerence between the velocity of the corresponding peaks in the presence and in the absence of a driving force. This indeed depends on the relative values of D and η, and also on |F|. (ii) We expect the true diﬀusion toward the centre to be enhanced by the inward force and competed by the outward one. Thus the two behaviours shown in ﬁgures 6.2 and 6.3 are the exact opposite of what is expected. It should be because of a minor algebraic error in the code. Unfortunatley, time lacked to check this before the deadline. 1 It is meant by the two moving peaks: the moving one shown in the given ﬁgures and the other one on the other half which is not shown here and which should be moving in the opposite direction to the ﬁrst one. 2 It was not possible to achieve this moment because of the breakdown of the stability of our numerical explicit scheme. That is, the ﬁrst phase of the dynamics is only studied. 7. Conclusion In this work: 1. Some essentials on the subject of diﬀusion in solids have been reviewed. 2. The problem of clustering of point defects has been simpliﬁed and modeled succesfully. 3. The algorithm and the achieved code work and give stable numerical solutions within the range of the chosen values for the diﬀerent variables. 4. In terms of the obtained results, within the conceptual framework introduced in the earlier chapters, and based on the simpliﬁcations made in chapter 5: a great partial success has been achieved. 5. Simpliﬁcations introduced in chapter 5 allowed to solve the problem and help to gain a better understanding of the problem of clustering of point defects in solids. Thus, the obtained results are to be considered as ﬁrst approximation. Future Work In future work: 1. Asymmetric problems should be considered. That is, asymmetric forces, arbitrary initial distributions and arbitrary boundaries. 2. Realistic forces, diﬀusivity and mobility values (and forms) should be inserted into the numerical form of the drif-diﬀusion equation. 3. Non-isotropic (non-cubic) and inhomogenous crystals are to be studied. 4. The stability of the followed explicit scheme must by studied to determine the stability condition for the solution to the FP equation. 5. It is better to follow a stable implicit numerical scheme, which provides more ﬂexibility in applying the realistic values and realistic forms. 23 Acknowledgements “WHO DOES NOT THANK PEOPLE WILL NEVER THANK ALLAH” First of all, I thank my god for the uncountable granted blessings and for the help in completion of this work. Without Him, I would not have ﬁnished neither AIMS diploma nor this essay. I am indeed indebted to the guidance I received from my supervisors, Prof. David T. Britton and Prof. Margit Harting, who very patiently provided me with adequate materials and encouraged me to go a head with this work in spite of many set backs. I am grateful to the PhD students of the Solid State and Materials Physics Research Group at UCT: Abdulraﬁu Raji, Shadrack Nsengiyumva, and Lloyd Corker, for their assistance and encouragement. I am at loss when I think of expressing my sincere gratitude to the super tutor Matteo Smerlak (a very special man to me) who helped me in this work and during the whole year with his assistance in various forms at various times. To him I owe many things. I also appreciate being helped, in the numerical part of this work, by Dr. Aziz Ouhinou. His help was priceless. My thanks also to everyone gave me a hand with this essay, specially my mates Sara, Felix, Bewketu and Esra. I express my profound gratitude to all AIMS Council, donars and staﬀ; in particular Prof. Neil Turok and Prof. Fritz Hahne for the opportunity that they have given me to study at AIMS. Many special thanks to our father Igsaan Kamalie, Ms. Lynne Teixeira, sister Aida Mpofu, Mr. Jan Groenewald, Ms. Gudrun Schirge, Mr. Emmanuel Kongolo and all Mothers in the kitchen for the wonderful generosity and kindness they showed during our stay at AIMS. I will be failing if I do not mention the angels of AIMS: the lecturers who volunteered themselves for teaching us. Their knowledge, their help, their kindness, their smiles and their patience have left lasting fervent impressions on me. I extend my warm thanks to all our tutors. Special thanks to Dr. Laure Gouba who was always there whenever I needed her. My many thanks are also to our English teacher, Ms. Frances Aron, for her unlimited support. I would like to express my sincere appreciation to my colleagues at AIMS. You guys are all really wonderful. My thanks to all of you for the wonderful atmosphere we shared. May our friendship continue to grow. I wish you all every success in the future. May I ask each and every one of you for forgiveness. And may my god guide all of us to His heaven. I heartily thank my family, to whom I owe everything, for their love, their prayers, their support, their constant encouragement, and for their patience. Finally, I would like to dedicate this essay to my lovely mother and to my two lovely daughters. 24 References [All01] [Arf85] [BH02] [Bri94] [FLS75] [GT88] M. Talb Allah, The science of engineering materials, in Arabic, Omdurman Islamic University, 2001. George Arfken, Mathematical methods for physicsists, third ed., Academic Press, Inc., 1985. D.T. Britton and M. Harting, The inﬂuence of strain on point defect dynamics, Advanced Engineering Materials 4 (2002), 629–635. D. T. Britton, Time-dependent diﬀusion models of the annihilation of low energy positrons implanted in solids, The Royal Society 445 (1994), 57–67. Richard P. Feynman, Robert B. Lighton, and Matthew Sands, The feynman lectures on physics, Addison-Wesley Publishing Company, 1975. Harvey Gould and Jan Tobochnik, An introduction to computer simulation methods, application to physical systems, Addison-Wesley Publishing Company, 1988. [HNR+ 07] M. Hartin, S. Nsengiyumva, A. T. Raji, G. Dollinger, P. Sperr, S. R. Naido, T. E. Derry snd C. M. Comrie, and D. T. Britton, Near surface stress determination in kr-implanted polycrystalline titanium by x-ray sin2 φ − method, Surface and Coatings Technology 201 (2007), 8237–8241. [Kit54] [Kit05] [LMN07] [Nea03] [Pea02] [PHI91] [Ris96] [She63] [Tab00] [Wal92] [Wip76] [ZC97] Charles Kittel, Introduction to solid state physics, John Wiley and Sons, Inc., 1954. C. Kittel, Introduction to solid state physics, eigth ed., John Wiley and Sons, Inc., 2005. T. B. Liverpool and K. K. Muller-Nedebock, Random walks, polymer and biopolymers: statics and dynamics, second ed., AIMS course study guide, 2007. Donald A. Neamen, Semiconductors physics and devices: Basic principles, third ed., McGrawHill Higher Education, 2003. William H. Press and et al, Numerical recipes in c, the art of scientiﬁc computing, 2 ed., Cambridge University Press, 2002. ´ Jean PHILIBERT, Atom movements diﬀusion and mass transport in solids, Les Editions de Physique, 1991. H. Risken, The fokker-planck equation, methods of solution and applications, 2 ed., Springer, 1996. Paul G. Shewmon, Diﬀusion in solids, McGraw-Hill Book Company, 1963. D. Tabor, Gases, liquids and solids and other states of matter, third ed., Cambridge University Press, 2000. Alan J. Walton, Three phases of matter, 2 ed., Oxford University Press, 1992. H. Wipf, The gorsky eﬀect, electrotransport and the thermotransport of hydrogen in metals, Journal of Less-Common Metals 49 (1976), 291–307. Dennis G. Zill and Michael R. Cullen, Diﬀerential equations with boundary-value problems, International Thomson Publishing Inc., 1997. 25

DOCUMENT INFO

Shared By:

Categories:

Tags:
Numerical, Solution, Drift-Diffusion, Equation, Point

Stats:

views: | 337 |

posted: | 12/19/2009 |

language: | English |

pages: | 29 |

Description:
Numerical Solution to the Drift-Diffusion Equation of Point

OTHER DOCS BY monkey6

How are you planning on using Docstoc?
BUSINESS
PERSONAL

By registering with docstoc.com you agree to our
privacy policy and
terms of service, and to receive content and offer notifications.

Docstoc is the premier online destination to start and grow small businesses. It hosts the best quality and widest selection of professional documents (over 20 million) and resources including expert videos, articles and productivity tools to make every small business better.

Search or Browse for any specific document or resource you need for your business. Or explore our curated resources for Starting a Business, Growing a Business or for Professional Development.

Feel free to Contact Us with any questions you might have.