A unified materials approach to mitigating optical nonlinearities in optical fiber. II. A. Material additivity models and basic glass properties

Funding information US Department of Defense High Energy Laser Joint Technology Office (HEL JTO), Grant/Award Number: N00014-171-2546; J. E. Sirrine Foundation Abstract The purpose of this paper, Part IIA in the trilogy (Int J Appl Glass Sci. 2018;9:263-277; Int J Appl Glass Sci. 2018 (in press); Int J Appl Glass Sci. 2018 (in press)), is to describe the continuum models employed to deduce the physical, acoustic, and optical properties of optical fibers that exhibit intrinsically low optical nonlinearities. The continuum models described herein are based on the additivity approaches of Winklemann and Schott (W-S). Initially developed over 120 years ago, W-S additivity works well for predicting the basic properties of bulk silicate glasses. While high-silica-content glasses are still the gold-standard for telecommunication and high energy laser fibers, the models have been systematically expanded to include deduction of the physical, thermophysical, and acoustic constants and coefficients that bear on parasitic nonlinearities. The stateof-the-art in W-S-based continuum materials models is reviewed here with specific examples provided based on canonical material systems suggested from the findings of Part I and treated in detail in Part III.


| INTRODUCTION
Without today's understanding of atomistic and quantum physics and chemistry and the benefits of modern high performance computing, our scientific ancestors developed continuum models for predicting the behavior of a range of materials. Amongst the original approaches was that of Adolph Winkelmann, who introduced an additivity approach for calculating the specific heat of (oxide) glasses of arbitrary composition based on the physical properties of the individual components. 1 That said, one does today have the benefit of high performance computing and numerous atomistic and quantum chemical models and codes with which to simulate the performance of glasses. [2][3][4] As a matter of fact, the simulation of glass properties and processing has become so effective, that glasses have been brought to commercial markets entirely designed and optimized based on modeling. 5 Accordingly, the Readers are then within their rights to wonder why continuum-level additivity models as described herein are employed. The simple answer is simplicity . . . and speed and versatility. In the work described in Companion Paper III, 6 there are entire binary, ternary, and quaternary glass systems whose properties need to be estimated in order to predict their physical, acoustic, and optical behavior as relates to the nonlinearities described in Paper I. 7 A macroscale continuum approach offers a rapid yet sufficiently accurate method for screening large computational spaces.
The mere fact that Hooke's "Ut tensio sic vis" 8 so well describes the linear elasticity of solids-ensembles of~10 23 atoms/molecules-is a great testament to the power of continuum approaches.
More specifically, the materials modeling described herein is not intended to provide specific insights into glass structure or chemistry. Its principal purposes are to (i) enable sufficiently accurate extrapolation of data measured on a few fibers to a broader range of compositions in that glass family; (ii) help quickly identify precursor core phases that will yield fibers of the desired properties; (iii) provide a simplified tool to present the enabling materials science to fiber designers (who are, at best, amateur glass scientists); and (iv) after fiber fabrication, provide guidance for further compositional optimization.

| MATERIAL PROPERTY ADDITIVITY
The additivity models utilized here, and in Companion Paper II B, 9 treat the constituent species comprising the glass as (i) separable and fully independent of each other and (ii) individually well-mixed with no heterogeneities, such as phase separation or crystallization; see Figure 1 for a representative schematic. Put another way, the models predict ensemble averages. A case in point is that of alumina in an aluminosilicate glass where properties as a general function of alumina content, [Al 2 O 3 ], are deduced even though the aluminum in Al 2 O 3 is known to exist in multiple coordination environments. 10,11 Further, reactive species, such as Al 2 O 3 and P 2 O 5 , which can form AlPO 4 , an index-lowering compound relative to SiO 2 , 12 are treated as the single final reactant.
As such, and as noted above, this averaging removes direct physical insight into the glass structure and chemistry, but reflects the power of these additivity models as effective tools for the straight-forward property deduction over large compositional spaces. Therefore, they are useful for the rapid identification of compositions and compositional ranges for designer optical fibers with specific performance attributes, such as the intrinsically low nonlinearities considered in this trilogy.
Beginning with a simple binary system where m is the fractional volume* of component A and, therefore, (mÀ1) is the fractional volume of component B, one can compute the time-of-flight for light to propagate through this "unmixed" glass segment. If said piece of glass is defined as being of unit length and its constituents separated into independent volumes (ie, lengths in a one-dimensional approximation) that are proportional to their concentrations, then the time, t, taken for the light to traverse the glass is: where n A,B are the refractive indices of constituents A and B, respectively, and c is the speed of light in vacuum. Given the unit length of this representative binary glass, the time given in Equation 1 is equivalent to an inverse velocity of light. As such, the numerator in Equation 1 is the linear refractive index of the aggregate glass, n binary = mn A + (1Àm)n B . As will be discussed further later, to first order, mixed effects such as thermo-optic, stress-optic, and optical Kerr effect-related nonlinearities can be parametrically included. The (dimensionless) fractional volume, m, of component A can be determined from the known composition along with the mass density, q, and molar mass, M, of the 2 components (A and B). For a binary system, this is given as, 13 where [C A ] is the concentration of the component A (in mole fraction or percent). The ratio M/q is the molar volume, which can depend on the glass thermal history. While the refractive index is the logical place to begin when considering optical materials, the additivity approach may be generalized to other material properties. In its most basic form, the additivity models follow the governing equation, Here, x is the additivity parameter of component i (in this case fractional volume), g i is the physical property of component i, and G is the aggregate property value of the glass.
Provided below are examples for the additivity of the optical, acoustic, and thermal properties required to model the performance of fibers designed to materially mitigate the optical nonlinearities described in Companion Paper I. 7 Additional details can also be found in Ref.14.

| A note on the use of crystal versus glass properties
In some cases, the properties of individual glass formers or of multicomponent glasses are known directly and the effective properties of a given component can be deduced by interpolation or extrapolation. However, in some cases, only property values for the crystalline precursor phase is *Fractional volume is employed here since it is the refractive index that is being used as an example, which is a volumetric material property. Where additivity in mole fraction is more appropriate, such as in the section on Debye Temperature and Specific Heat, this will be specified.
known. This section treats the estimation of the properties of an amorphous analog to those of the known crystalline state. These "glassy" † values can then be combined using the additivity approaches described herein.
The method relies on homogenizing (ie, "relaxing" in technical vernacular) the single crystal constants, for each substance, in a self-consistent manner to obtain equivalent isotropic values. With the elastic constants as an example, the starting point is homogenization (amorphization) of the elasticities of each crystal constituent using a simple extension of the Voigt-Reuss-Hill (VRH) procedure. [15][16][17][18] The extension, denoted herein as VRHx, yields a single set of relaxed quantities for each material, rather than a range of values. An instructive and typical comparison of scenarios (A) and (B) may be made by using the admirable experimental results of Nassau et al, 19 on silica, along with the VRHx procedure applied to the a-alumina values of Goto et al. 20 Addition of alumina to silica is shown to result in linear variations of density and acoustic velocities (shear and pressure waves), at least up to 10 mole percent alumina. 19 From these relations, using the 10 mole percent values, one finds, for the mixture, density <q> = 2.284 [g/cm 3 ], and elastic stiffness values <c 11 > = 90. 16 and Subjecting the a-alumina values of Goto 20 to the VRHx algorithm gives: <c 11x > = 471.23, <c 44x > = 163.12 [GPa]. Alumina density is taken as q = 3.982 [g/cm 3 ]; this value is that of the crystal. When incorporated into a glass, this value will change in a manner dependent on interatomic forces and molecular sizes; determination of densities of component species in any given situation is a matter of some complexity, 21 and is one reason for the loss of accuracy when using the VRHx procedure in conjunction with the density of the starting crystal.
When the VRHx alumina values, along with the crystal density, are used in the W-S additivity formula at 10 mole % alumina in silica, one finds: <q (mix) > = 2.368 (3.7%); <c 11 (mix) > = 87.75 (2.7%) and <c 44 (mix) > = 34.73 (2.3%). Use of the word "glassy" here is not meant to imply the existence of a glass transition, as is one of the conventional requirements for a material to be considered a glass (the other being that the material is a non-crystalline solid). Its use is simply a short-hand way to denote an amorphous analog to the crystalline composition.
The values in parentheses are the errors compared to Ref. 19 at 10 mole percent alumina in silica. It is reasonable to assume that the accuracy could be improved by use of an improved alumina density value.
Typical scenario (C) examples would be VRHx-reduced a-alumina mixed with one of the VRHx-reduced polymorphs of silica. As indicated above, the resultant mixture predictions are less accurate, compared to the aluminosilicate glass results in Ref. 19. For b-cristobalite, for example, one obtains results within 10%-15% for mixtures of a few mole % alumina in silica.
Voigt 15 averaged elastic stiffnesses over all directions in space whereas Reuss 16 did the same for elastic compliances, which are the matrix inverse of stiffness values. These produce values for non-physical, hypothetical, isotropic substances. Hill 17 took arithmetical averages of these upper (Voigt) and lower (Reuss) bounds, finding that the averages for Young's modulus and shear modulus then agreed better with experiment than either the Voigt or Reuss bound. The Hill averages still do not produce completely consistent values when the matrix inverses are compared. The VRHx procedure, which uses repeated averaging of the averages and inverses, does just this.
The VRHx approach to estimating the elastic properties of a glass from the corresponding crystal begins with the Voigt stiffness matrix: c V = c V (0) and the Reuss compli- ) and s RV (n+1) = ½(s RV (n) + s VR (n) ) rapidly converge to a consistent set of "best" isotropic matrices, <c x > and <s x > = <c x > À1 . From these coefficients and the mass density, the shear (S) and pressure (P) wave speeds for component i are found from, Using both the S and P wave-speeds so determined, one finds the isotropic elastic constants of the mixture. Key here is the use of slownesses (V Pi À1 and V Si À1 ) in determining average elasticities. Slowness is proportional to transit-time, and has a direct physical meaning, whereas the meanings of velocity averaging, or other criteria, are less obvious.
Mass density of the mixture (q (mix) ) is determined from Equation 3, using the component densities. Wave speeds (V P (mix) & V S (mix) ) are found from Equation 3 using . Introducing the Lam e constants, k and l, the elasticities of the multicomponent glass are found using the following relationships: Isotropic stiffnesses: Isotropic compliances: Figure 3 provides a comparison of the density and acoustic velocities for aluminosilicate glasses. The computed values employ the VRHx amorphization approach applied to crystalline b-cristobalite (for silica) and sapphire (for amorphous alumina) end-members coupled with W-S additivity for the binary compositions.
The experimental values are based on measured data from 19,22,23 and trend-lines provided by Nassau et al. 19 ‡ Elasticity data for b-cristobalite employed for the determination of the acoustic velocities is from Ref. 24,25 whereas that for sapphire is from Refs. 18,20. The use of b-cristobalite as the crystalline phase from which the properties of silica glass are computed is because, on those occasions when SiO 2 crystallizes from the melt, b-cristobalite is the phase that generally forms. 26 Further, it possesses the density of silica glass at ambient conditions (2.200 g/cm 3 ) and is crystallographically cubic (m 3m). 18,27 Over the range from 0 to about 50 mole percent, where aluminosilicate glasses have been formed, the modeling and experimental results match more than adequately for the purposes of this trilogy. Perhaps most importantly is that the slope, representing the change in the property, (density or acoustic velocity), with alumina concentration, is consistent with those reported in the literature.

| Refractive index (n)
The refractive index is the central optical property of a material and enters into nearly all considerations of its linear and nonlinear performance.
As noted in Ref. 28 a power of the additivity model can be seen in its remarkable resemblance to the often-used Sellmeier model for a mixed system. The Sellmeier model, which can be directly derived from the physics of the macroscopic refractive index, is given for a unary material by, 29,30 Here, n is the refractive index at (free space) wavelength k, N is the number of oscillators used to fit the refractive index data (typically 3 as in Ref. 29 or 4 as in Ref. 30; this is often referred to as the "term" of the Sellmeier fit; ie, a 3-term Sellmeier), A i and k i , are the oscillator strength and its resonant wavelength, respectively. The binary form of Equation 6 can be found in Ref. 29 where additivity in the mole fraction is usually employed.
As an illustrative example of the simplicity and power of the additivity model, consider the comparison provided in Figure 4. Plotted there is the refractive index of the SiO 2 -GeO 2 binary system as calculated using the additivity model and the Sellmeier model (data taken from Ref. 29). The data used for the additivity model results are the mass densities (2200 kg/m 3 for SiO 2 and 3650 kg/m 3 for GeO 2 ), molar masses (60.08 g/mole for SiO 2 and 104.59 g/mole for GeO 2 ), and end-member refractive index values (1.4442 for SiO 2 and 1.5873 for GeO 2 taken for a wavelength of 1534 nm). The index values between the 2 models differ only in the 5th decimal place.

| Thermo-optic coefficient (TOC)
The refractive index as a function of temperature in a linear regime may be expressed simply as n(T) = n 0 + dn/dT (TÀT 0 ) where the subscript "0" refers to some starting temperature (T 0 ) such as room temperature. This expression may be substituted into Equation 1 to arrive at the following simple expression for the additivity of the TOC in a binary material, Generally, for a multicomponent glass Equation 7 may be used with g i ?dn i /dT.

| Acoustic velocity (V A )
The acoustic velocity can be estimated using a similar time-of-flight consideration for an acoustic phonon as the refractive index was calculated using a photon. In this case, both the aggregate acoustic velocity, V a , and the total acoustic attenuation, a, through the unit fiber can be deduced. 13,31 The determination of the acoustic attenuation is of interest to this work since multiplying it by (V a /p) yields the Brillouin spectral width, Dm B , which will be used later to compute the Brillouin gain coefficient, BGC. For an n-component glass, these quantities are given by, The factor F i (m) in Equation 9 typically is given by where a i are the phonon damping coefficients (m -1 ) of the oxide constituent(s), which typically are proportional to the square of the acoustic frequency. 32 For completeness, however, it is worth noting that this relationship generally holds for glasses, such as silica, where the maximum viscoelastic damping occurs far below the glass transition temperature, T g . 33 For glasses such as those doped with B 2 O 3 , where the dynamic viscoelastic damping processes peak at temperatures above T g , the acoustic attenuation is not proportional to m 2 . 34 Further, since there is dispersion in the acoustic attenuation (a = a(m)), a i should be specified at some reference frequency, m ref . For ease, a convenient reference frequency is 11 GHz, which corresponds to the Brillouin scattering frequency of standard single mode telecommunications fiber, when measured at a wavelength of 1534 nm. 13 In summary, since m and therefore a (at a fixed optical wavelength) change with composition, F i (m) takes this into consideration by scaling the a i reference values accordingly.

| Photoelastic coefficient (p)
As shown and discussed in Companion Paper I, 7 the Pockels photoelastic coefficient, p 12 , factors into both Brillouin scattering and into the contribution of density fluctuations to Rayleigh scattering. In order to deduce the p 12 coefficients, consider the axial (z-axis) elongation of a fiber due to an applied strain, e z . The transverse electric field components of the propagating (assume x-polarized) optical mode will experience a concomitantly modified refractive index, n x , given by, where m is Poisson's ratio and the subscript "0" corresponds to the nominal (unstrained) value. 35 This equation is also valid for strain applied in the transverse (y) direction also perpendicular to the optical polarization. If a transverse strain is applied to the fiber in the direction of polarization (e x ), Equation 10 becomes, 36 n x ðeÞ ¼ n 0;x À 1 2 e x n 3 0;x ½p 11 À 2mp 12 : With these in hand, and assuming that the photoelastic constants can be parametrically included into the additivity model through the refractive index, then deduction of the photoelastic coefficients are simplified using the transverse (ie, strain transverse to polarization direction, e y,z ) strainoptic coefficient (eOC), which is defined as eOC = p 12 Àm (p 11 + p 12 ) and the parallel (ie, strain parallel to the polarization direction, e x ) stress-optic coefficient (rOC), which is defined here as rOC = p 11 À2mp 12 . These definitions are selected to disambiguate these values from the usual stress optic coefficient, C, which is related to p 44 = (p 11 Àp 12 )/2 for an isotropic material, by C = À(n 3 /2c 44 ) p 44 .
Accordingly, then, for a mixed binary glass of Components A and B, the two effective strain-optic coefficients are, where n 0 is the zero-strain refractive index of the binary glass and, as above, m is the fractional volume of Component 1. Given the transverse and parallel strain-optic components, the p 12 photoelastic constant for the aggregate glass is: 36 Figure 5 provides the calculated p 12 photoelastic coefficient as a function of additive into silica for a number of crystal-derived all-glass optical fibers, 22,[37][38][39][40][41] along with corresponding measured data. For the cases of BaO, Al 2 O 3 , and YAG, the additive model was fit to fiber-based measurements and values of p 12 for the non-silica additive were subsequently determined. The calculation in Figure 5 provides an extrapolation to points beyond compositions available for characterization, which should be taken with caution as any structural changes to the glass outside the available compositional range will cause the fittings to break down. In the available compositional ranges, and in all cases, the fittings are found to be very good. For the alkali metal oxides, the model was fit to bulk data found in Ref. 42, which is also reproduced in Figure 5. In the case of spinel-derived fiber, the acoustically anti-guiding behavior precluded the determination of p 12 directly, so the value for MgO found in Ref. 43 (p 12 = À0.07) was instead used in the modeling. Finally, in the case of SrO, its p 12 was deduced from a ternary glass (strontium aluminosilicate) and that value was utilized in the modeling of the binary system.

| Coefficient of thermal expansion (CTE)
In general, thermal effects in fiber-based lasers are known to induce a variety of performance-modifying phenomena. 44 While the coefficient of thermal expansion (CTE) does not directly influence the optical nonlinearities that are the focus of this trilogy, it does influence the effective thermo-optic coefficient, as first suggested by Prod'homme, 45 which is believed to induce transverse mode instabilities (TMI). Additionally, differential thermal expansion between the core and clad (or cladding and coating) can lead to polarization effects, 46 increased transmission losses, 47 and phase shifts, the latter causing jitter in transmission systems and noise in interferometric sensors. 48 It is worth noting that, in addition to differential thermal expansion, thermal stresses also can be generated from differences in core/clad viscosity. 49 Core/clad CTE mismatches can also lead to some novel and useful phenomena not present in bulk glasses, such as temperature independent (athermal) and strain independent (atensic) Brillouin frequency shifts. 40,50,51 Many empirical models have been developed regarding the additivity of the coefficient of thermal expansion in glasses. [52][53][54][55][56] Recently, the authors (MC, PD, JB) have used the additivity of density, which is governed by Equation 3, to develop an additive model for CTE in binary silicate glasses, 57 based on a derivation of the CTE for isotropic solids such as glass. The expression of CTE takes the following form, Here, CTE i and q i are the coefficient of thermal expansion and the density for each individual constituent of the glass, respectively, whereas y i is defined as the volume fraction per structural unit for each glass compound.
As shown in Figure 6, Equation 15 fits the data extremely well, not just when considering the addition of ubiquitous glass formers in silica ( Figure 6A), but also for some others dopants such as the alkaline earth family (Figure 6B). [58][59][60][61][62][63] For completeness, it is noted that a proper additivity relationship needs to be rooted in some physical basis. For example, linear fits to the data in Figure 6 would likely lead to acceptably good fits over a certain range, but would have no physical meaning. The time-of-flight approach and volume additivity for the refractive index (Figure 4), and this CTE additivity ( Figure 6) based on the additivity of density both have their origins in the nature of the systems which they are modeling.

| Specific heat and Debye temperature
The specific heat was shown in Companion Paper I 7 to factor into both stimulated thermal Rayleigh scattering, which facilitates transverse mode instabilities (TMI), and the density-related component of (classical) Rayleigh scattering. Accordingly, its estimation through simple additivity relationships is of practical value. Toward this end, it is convenient first to estimate the Debye temperature, Θ D , 64 which arises in the consideration of the specific heats of elastic solids and provides a quantitative estimate demarcating classical (Dulong-Petit 65 ) and quantum regimes. The Debye temperature can also be informative as to the thermal conductivity of a solid, 66 which also factors into TMI through STRS per the Dong formalism. 67 F I G U R E 5 Comparison of the calculated and measured/deduced Pockels p 12 photoelastic coefficient as a function of non-silica content for selected fibers deduced using the additivity approach described in the text. *Because of the acoustic anti-guiding nature of the spinelderived optical fiber, the p 12 photoelasticity for MgO, taken from Ref. For a single component, linearly-elastic solid, Θ D is given by where h is Planck's constant, k is Boltzmann's constant, q is the number of atoms in the molecule, q is mass density, N is Avogadro's number, and M is molecular mass. V is the average acoustic wave speed given by V = ffiffiffiffiffiffiffiffi ffi C=q 2 p , where C is either c D , (Θ = Θ D ), the true average elastic stiffness, or <c mix >, (Θ = Θ Dmix ), the extended Voigt-Reuss-Hill (VRHx) approximate average elastic stiffness (see Section 2.1 "A note on the use of crystal vs glass properties").
For glasses comprised of a mixture of components, Equation 16 is modified as follows. With f i (mol) being the mole fraction of the ith component, the expression for Debye temperature is: , whereas the previous additivity modeling generally employed volume fractions.
The elastic wave speed averaging procedure associated with Debye temperature calculations affords a simple comparison of the VRHx results with a more accurate computation. 68 The Debye procedure globally averages the 3 plane acoustic wave speeds (V i ) in a crystal over all propagation directions The waves are labeled i = a, b, and c, denoting, respectively, the quasi-longitudinal (pressure) wave, the fast and slow quasi-shear (transverse) waves. The elastic stiffness associated with each wave velocity is c i = qÁV i 2 . In terms of the spaceaveraged elastic stiffness, the requisite Debye expression is: The VRHx algorithm also produces a global space averaging; the resultant equivalent elastic stiffness corresponding to c D is: For VRHx averaging procedures to render accurate results it is necessary that c x % c D . As a practical computational matter, the procedure necessitates covering a sphere with a grid of finite patches, computing the velocities at a convenient point within each patch, then summing the result over all patches to approximate the true values that would be obtained in the limit as the patch sizes diminish to zero. Depending on the grid size and shape, this has, in the past, resulted in either a loss of accuracy, or excessive computational effort. The essential difficulty is the impossibility of uniformly tessellating a sphere. However, a precise and efficient algorithm based on Fibonacci sequences 69 is available to arrive quickly at accurate results for c D to compare with <c mix >. It is found that these generally agree within a few percent.
The molar specific heat capacity, c V , of the glass at a given temperature then can be determined from knowledge of the Debye temperature using the well-known relationship: x 4 e x ðe x À 1Þ 2 dx: The Debye temperature for silica is about 365K in the range from 70K to 370K. 70 Accordingly, one is neither in the c v / (T/Θ D ) 3 quantum nor C v /(3ÁN A Ák) % 1 Dulong-Petit regimes. Tabulated values of the Debye integral are, however, readily available.

| CONCLUSIONS
Provided herein is a compilation of simple additivity models for the deduction of the material properties of multicomponent silicate glasses, central to optical nonlinearities, from the properties of amorphous or crystalline end-member compounds. More specifically, additivity as relates to the refractive index, density, acoustic velocity, photoelasticity, stressand strain-optic coefficients, thermal expansion coefficients, thermo-optic coefficients, Debye temperature, and specific heat are discussed and examples provided. In general, over reasonable glass forming ranges in which the glass is presumed to remain compositionally uniform, homogeneous, and structurally similar, the models are found to be sufficiently accurate in order to guide the materials development of optical fibers that exhibit reduced optical nonlinearities. Use of these material properties to model the nonlinear coefficients will be provided in the second half of this paper.

ACKNOWLEDGMENTS
The Authors wish to thank Drs. Siva Mani, Matt Leigh, and Sarwat Chappell of the US Department of Defense High Energy Laser Joint Technology Office (HEL JTO) for posing the original question that led to this trilogy. The Authors gratefully acknowledge the financial support from the JTO through sustained support over many years including contracts: W911NF-05-1-0517, FA9550-07-1-0566, W911NF-12-1-0602, FA9451-15-D-0009/0001 and 0002, and N00014-17-1-2546. The J. E. Sirrine Foundation also is gratefully acknowledged for supporting the efforts of Authors (JB and MC).