X-ray Diffraction Analysis and Williamson-Hall Method in USDM Model for Estimating More Accurate Values of Stress-Strain of Unit Cell and Super Cells (2 × 2 × 2) of Hydroxyapatite, Confirmed by Ultrasonic Pulse-Echo Test

Taking into account X-ray diffraction, one of the well-known methods for calculating the stress-strain of crystals is Williamson-Hall (W–H). The W-H method has three models, namely (1) Uniform deformation model (UDM); (2) Uniform stress deformation model (USDM); and (3) Uniform deformation energy density model (UDEDM). The USDM and UDEDM models are directly related to the modulus of elasticity (E). Young’s modulus is a key parameter in engineering design and materials development. Young’s modulus is considered in USDM and UDEDM models, but in all previous studies, researchers used the average values of Young’s modulus or they calculated Young’s modulus only for a sharp peak of an XRD pattern or they extracted Young’s modulus from the literature. Therefore, these values are not representative of all peaks derived from X-ray diffraction; as a result, these values are not estimated with high accuracy. Nevertheless, in the current study, the W-H method is used considering the all diffracted planes of the unit cell and super cells (2 × 2 × 2) of Hydroxyapatite (HA), and a new method with the high accuracy of the W-H method in the USDM model is presented to calculate stress (σ) and strain (ε). The accounting for the planar density of atoms is the novelty of this work. Furthermore, the ultrasonic pulse-echo test is performed for the validation of the novelty assumptions.


Introduction
Young's modulus (E) plays an important role in the characterization of mechanical properties, and accurate knowledge of the engineering values of elastic modulus is vital for understanding materials' stiffness [1]. Recently, ceramic materials have been favored in industrial applications, because they exhibit good mechanical properties, such as high elastic moduli and high hardness [2]. One of the well-known bio ceramics is hydroxyapatite, which has bioactive properties very close to natural bone mineral and it has special biological significance [3]. The chemical formula of hydroxyapatite is Ca 10 (PO 4 ) 6 (OH) 2 and it differs very little from bone tissue [4,5]. Understanding the mechanical properties of hydroxyapatite during the crystallization and growth stages of the synthesis processes is

Methods
The synthesis route of hydroxyapatite powder is explained completely in part 2 of the Support Information (preparation of hydroxyapatite powder).

Structural Analyses of Synthesized Hydroxyapatite
The XRD pattern of synthesized hydroxyapatite powder is shown in Figure 1. The XRD pattern exhibits several diffraction peaks, which can be indexed as the hexagonal hydroxyapatite. The XRD pattern was evaluated based on X'pert and the pattern was in agreement with the standard XRD peaks of hydroxyapatite (ICDD 9-432). Similar observations were reported in References [11,12]. In addition, crystallographic parameters and details of each diffracted plane of hydroxyapatite were evaluated by X'pert and the values are tabulated in Tables 1 and S1.  According to Table S1, the values of the distance between planes are calculated by Equation (1). In this equation, h, k and l are indices of each plane, and a, c and d are lattice parameters and distance of planes, respectively [13].
Hydroxyapatite has a hexagonal system with a P63/m space group and has little deviation from stoichiometry [14]. Figure 2 shows a sketch of a unit cell of hexagonal hydroxyapatite and a cif file of synthesized hydroxyapatite. There are two different situations of calcium ions and, in total, 18 ions are closely packed to create the hexagonal structure. At each hexagonal corner, a calcium ion is restricted by 12 calcium ions shared with 3 hexagons. Void spaces between two hexagons are filled with three phosphate tetrahedral per unit cell. Ions in hydroxyapatite can be interchangeably replaced with biologically useful ions due to the inherent versatility of this crystal structure and can also be referred to as doping. In addition, the substitution of calcium, phosphate and/or hydroxyl ions is possible [15]. Notably, the specific feature of hydroxyapatite is related to the OH ions forming inner channels along the c axis. This property plays an important role in its mechanical and physical properties [16]. In addition, the Edax analysis of synthesized hydroxyapatite is presented in Figures S2. According to the EDX analysis, the value of the Ca/P ratio for hydroxyapatite obtained from bovine bone was recorded to be 1.60. In addition, thermal decomposition of hydroxyapatite into tricalcium phosphate and tetra cal-  According to Table S1, the values of the distance between planes are calculated by Equation (1). In this equation, h, k and l are indices of each plane, and a, c and d are lattice parameters and distance of planes, respectively [13].
Hydroxyapatite has a hexagonal system with a P6 3 /m space group and has little deviation from stoichiometry [14]. Figure 2 shows a sketch of a unit cell of hexagonal hydroxyapatite and a cif file of synthesized hydroxyapatite. There are two different situations of calcium ions and, in total, 18 ions are closely packed to create the hexagonal structure. At each hexagonal corner, a calcium ion is restricted by 12 calcium ions shared with 3 hexagons. Void spaces between two hexagons are filled with three phosphate tetrahedral per unit cell. Ions in hydroxyapatite can be interchangeably replaced with biologically useful ions due to the inherent versatility of this crystal structure and can also be referred to as doping. In addition, the substitution of calcium, phosphate and/or hydroxyl ions is possible [15]. Notably, the specific feature of hydroxyapatite is related to the OH − ions forming inner channels along the c axis. This property plays an important role in its mechanical and physical properties [16]. In addition, the Edax analysis of synthesized hydroxyapatite is presented in Figure S2. According to the EDX analysis, the value of the Ca/P ratio for hydroxyapatite obtained from bovine bone was recorded to be 1.60. In addition, thermal decomposition of hydroxyapatite into tricalcium phosphate and

Planar Density of Unit Cell and Super Cells (2 × 2 × 2) of Hydroxyapat
According to the resulting list of planes by X-Ray diffraction (Fig planar density values of each diffracted plane in unit cell and supe hydroxyapatite are calculated in Figures 3, S6 and S7. To calculate the ues, diffracted planes were selected from a low angle to a high angle in mentioning that the matrix of all super cell lattices was considered to ing to the center of atoms, the planar density is calculated by the area plane divided by the total area of that plane [19].

Planar Density of Unit Cell and Super Cells (2 × 2 × 2) of Hydroxyapatite
According to the resulting list of planes by X-ray diffraction ( Figure 1, Table S1), the planar density values of each diffracted plane in unit cell and super cells (2 × 2 × 2) of hydroxyapatite are calculated in Figure 3, Figures S6 and S7. To calculate the planar density values, diffracted planes were selected from a low angle to a high angle in tandem. It is worth mentioning that the matrix of all super cell lattices was considered to be 2 × 2 × 2. According to the center of atoms, the planar density is calculated by the area of the atoms in the plane divided by the total area of that plane [19].

Elastic Stiffness Constant and Elastic Compliance of Hydroxyapatite
Hooke's law is shown in Equation (2); the stress corresponds to the strain for small displacements. It is the basic form, that this symmetry can be converted to the six items of stress (σ) and strain (ε) [20].
The elastic stiffness determines the response of crystal to an externally applied stress or strain and provides information about bonding characteristics and mechanical and structural stability [21]. The hydroxyapatite system has five elastic constants (Equation (4)). Therefore, the values of five independent C ij , can be named C 11 , C 12 , C 13 , C 33 , C 44 .

Planar Density of Unit Cell and Super Cells (2 × 2 × 2) of Hydroxyapatite
According to the resulting list of planes by X-Ray diffraction ( Figure 1, Table S1), the planar density values of each diffracted plane in unit cell and super cells (2 × 2 × 2) of hydroxyapatite are calculated in Figures 3, S6 and S7. To calculate the planar density values, diffracted planes were selected from a low angle to a high angle in tandem. It is worth mentioning that the matrix of all super cell lattices was considered to be 2 × 2 × 2. According to the center of atoms, the planar density is calculated by the area of the atoms in

Elastic Stiffness Constant and Elastic Compliance of Hydroxyapatite
Hooke's law is shown in Equation (2); the stress corresponds to the strain for small displacements. It is the basic form, that this symmetry can be converted to the six items of stress (σ) and strain (ε) [20]. The crystallographic nature of the hexagonal structure is shown in Figure S3. Furthermore, for conventional hexagonal systems, such as hydroxyapatite, the relationship between C ij and S ij is introduced in Equations (5)-(9) [22,23].
S 33 = C 11 + C 12  (9) According to the Equations (5)- (9), to obtain S values C values are needed. We have used two approaches to obtain C values: First, theoretical calculations were performed via the CASTEP model of materials studio software version 6.0 in the GGA level of theory with a PBE basis set ( Figure S4). The second approach is based on ultrasonic measurements (Table S2). The complete set of five elastic stiffness constant values (C 11 , C 12 , C 13 , C 33 and C 44 ) of the samples was found from ultrasonic measurements of the phase velocity anisotropy. The stiffness constant values were recorded by utilizing Equations (10)- (15). In these equations, ρ and V are density of sample and velocity, respectively [24][25][26][27][28].
For measuring the velocities, the standard ultrasonic pulse-echo ASTM E797/E797-M-15 was accomplished according to Reference [26]. In this study, the immersion procedure via water between probe and sample was utilized and the effect of pressure derived from a hand placed into the probe was decreased. In addition, the longitudinal frequency probe was 5.4 MHz, and to decrease the extension and depreciation of waves, especially the more energetic waves, a lens with a specific curvature was utilized according to the standard of Reference [28]. To measure the value of C66, a high amplitude of curve was carried out; therefore, changing the rotation angles was useful. To create the transverse waves, a probe of 2.3 MHz and a high viscos interface material (honey) were used and finally the pressure on the sample was adjusted by obtaining the better and smoother curve [29]. Taking into account each point of the samples, a three-dimensional axis, such as X 1 , X 2 and X 3, can be performed. In this case, according to Figure S5, X 1 is the radial coordinate, X 2 is the circumferential coordinate and X 3 is the axial coordinate. V i/j is the denoted velocity of an ultrasound wave that can be propagated in the X i direction with particle displacements in the Xj direction simultaneously. V i/j , with the same i and j, is longitudinal and with i = j being transverse waves. For measuring quasi-longitudinal or quasi-transverse velocity (Vij/ij), specimens should be cut (bezel) on the edges toward surfaces of perpendicular X directions ( Figure S5). The obtained values of velocities are shown in Table S2. On the other hand, C 11 is in good agreement with the longitudinal distortion and longitudinal compression/tension; thus, C 11 can be introduced as the hardness. Moreover, the transverse distortion depends on the C 12 , and C 12 comes from the transverse expansion related to the Poisson's ratio. Additionally, C 44 corresponds to the shear modulus, and C 44 is in the settlement with C 11 and C 12 [22,30]. Accordingly, the shear modulus is proportional to the Burgers vector and the Young's modulus; in addition, dislocation density is in agreement with the Young's modulus [31,32].
In the ultrasonic method, longitudinal and transverse waves were utilized for measuring the Young's modulus value [33,34]. According to this method (Equation (16)), based on the velocity of ultrasound waves and density of sample, the Young's modulus value was determined.
In Equation (16), ρ, c l and c t are density, velocity of longitudinal and transverse ultrasound waves tandemly. Furthermore, according to Equation (17), the velocity of longitudinal and transverse waves can be registered by determining the length of specimen and the differences between two echoes (t = t 2 − t 1 ) in the signals [35].
Here, L is the length of the sample and t is the difference between two echoes, and the density of the sample can be detected by measuring the mass and volume of the sample [36]. Additionally, with the substitution of Equations (16) and (17), the main equation for the calculation of Young's modulus is Equation (18).
t s and t l are differences between two echoes in longitudinal and transverse waves separately [37]. The results of the theoretical calculation, and the experimental measurements of elastic stiffness constant values of hydroxyapatite (from the literature and this study), are shown in Table 2. In addition, taking into account Equations (10)-(15), the elastic compliance values of this study were calculated and are recorded in Table 3. As a result, the theoretical values were in good agreement with the theoretical values extracted by References [38,39]. The replicated values (five times) for measuring the time of transverse and longitudinal waves of the samples are presented in Table 4. X-ray diffraction has provided data on diffracted planes and the location of atoms in each plane. In the previous section, the planar density of each diffracted plane was calculated and it could play an important role in the mechanical properties of each plane. The Young's modulus of each plane (E hkl ) of a hydroxyapatite lattice can be calculated as Equation (19) [5,43]. In this equation, h, k and l are the plane indices, a and c are the lattice parameters, C and S are the elastic stiffness constant and elastic compliance, respectively. The values of Young's modulus of 32 diffracted planes of hydroxyapatite in unit cells and super cells (2 × 2 × 2),E (h 1 k 1 l 1 ) , E (h 2 k 2 l 2 ) , E (h 3 k 3 l 3 ) . . . . . . . . . . . . ., E (h 32 k 32 l 32 ), related to the literature and the present study, are reported in Table S3.  (19) By using the least squares method between the Young's modulus and the planar density of diffracted planes (based on our proposed method), the Young's modulus value of hydroxyapatite was determined with high precision. Consequently, the calculation of C and S parameters for crystalline materials is essential for the application of this method. To show the feasibility and accuracy of our proposed method for determining the Young's modulus, the values of Young's modulus of each plane versus the planar density of the unit cell and supercells (2 × 2 × 2) are depicted in Figure 4. The Young's modulus values (intercept) of the unit cell and super cells (2 × 2 × 2) of hydroxyapatite, obtained from the least squares method, are tabulated in Table 5.  According to the uncertainty measurement (the measurements were replicated five times (Table 4)) and Equation (18), the Young's modulus value gained 113.08 ± 0.14 GPa  According to the uncertainty measurement (the measurements were replicated five times (Table 4)) and Equation (18), the Young's modulus value gained 113.08 ± 0.14 GPa by ultrasonic measurement. This value is in good agreement with the reported values of our study in Table 5. In this study, the difference between theory and experiment values for both unit cell and super cells (2 × 2 × 2) are identical. This difference for unit cell and super cells (2 × 2 × 2) are 10.47 GPa and 10.57 GPa, respectively. This means that the theoretical calculation is valid and, by reducing it by about~10 GPa, the experimental values can be obtained.

Positive and Negative Slope Values of Unit Cell and Super Cells (2 × 2 × 2)
The difference of Young's modulus values in the unit cell and super cells (2 × 2 × 2) is attributed to the locations of atoms. The calculated slope is a negative value in super cells (2 × 2 × 2). In the first aspect, it is clear that the slope is dependent on the planar density of diffracted planes, so the fitting based on the planar density of super cells with a matrix of eight unit cells could submit a better result of intercept [35]. It is because eight cells besides each other are completed and have more symmetry than two cells ( Figure S8) [35]. In the second aspect, the reason for the positive slope in the unit cell and the negative slope in the super cells (2 × 2 × 2) is related to the defects (imperfections) including point defects (vacancies, substitutional and interstitials), line defects (screw and edge dislocation), surface defects (grain boundaries) and volume defects (lack of order of atoms due to amorphous region in a very tiny area). The effect of these imperfections is more effective in super cells (2 × 2 × 2), while the unit cell is more ideal and less affected by these imperfections. This means that when the density of atoms in planes is increased, lower forces for dislocation motion are required and the strength will be decreased and, consequently, in super cells (2 × 2 × 2) the slope is negative. In the case of super cells, when the number of atoms is increased, by the increase of planar density, the effect of dislocation motion is increased; the strength and Young's modulus will be decreased so the slope is negative. The intercept of the fitting line is a value of Young's modulus, which can show the Young's modulus of the plane with zero planar density as a plane without any specific atom. Therefore, a discussion of defects and imperfections cannot be considered for such a plane without an atom at the origin and so, in this case, the Young's modulus value of the unit cell is less than that of the super cells (2 × 2 × 2) and consequently the plane without an atom is more realistic in a smaller area of a unit cell than in a wider area of the super cells. This means that the intercept in a unit cell is closer to the real values of the Young's modulus. As mentioned above, the planar density depends on the array and position of atoms into the plane and the situation of the planes. For illustration, according to Figure 4b, the planar density of (210) > (510) > (102) > (023) > (043), because sequencing of the planes is not related to the diffraction (XRD), but it depends on the planar density values. With this method, it can be possible to consider two or more planes with similar planar density values; for example, in super cells (2 × 2 × 2), planar density values of (040) and (331) are equal to 0.250. The control of the displacement and deformation process of the atoms in the planes are associated with the dislocation networks. The work or energy (W) for the movement of the atoms in each plane corresponds to Equation (20) [44,45]. W = Gb 2 l (20) In this equation, G, b and l are shear modulus, Burgers vector and dislocation length, respectively. In addition, by merging Equations (20) and (21), Equation (22) can be registered. In Equations (21) and (22), E and ν are Young's modulus and Poisson's ratio tandemly.
There are several studies on the properties of hydroxyapatite, especially with regard to biocompatibility and bioactivity, and these properties are dependent on the identification and recognition of hydroxyapatite structures [46,47] For example, when atoms are widely spaced, such as corner atoms (Ca) in (111) super cells (2 × 2 × 2), atoms require higher values of applied force for approaching, so knowledge on the planes of hydroxyapatite structure can be helpful for doping metals, polymers and other ceramics, for aims such as the controlled release of protein and also the fabrication of bioactive monolithic fragments in the biomaterials industry [48,49]. In addition, the hydroxyl ion is in the center of each unit cell. The hydroxyl group in the center of the hydroxyapatite lattice is surrounded by three calcium ions per hexagon, forming a ring (six calcium ions). A chord is formed by the structure, as these rings are responsible for many properties of hydroxyapatite, especially the biocompatibility and the position of the hydroxyl group in the planes have considerable importance for doping mechanisms [15,50].

Williamson-Hall Method in USDM Model
The W-H method has been used for determining different elastic properties. The best procedure is to mathematically reduce the errors and obtain the values of the Young's modulus by all the diffracted peaks, using the least squares method. The W-H method is a simplified integral expanse and, taking into account the peak width, strain-induced broadening is specified [51]. Taking into account the W-H method in the USDM model, Young's modulus values are examined. It is clear that in Equation (23) and Figure 5, the term of 4sin θ E hkl is along the X-axis and the term of β hkl .cos θ is along the Y-axis individually. (23) of each unit cell. The hydroxyl group in the center of the hydroxyapatite lattice is rounded by three calcium ions per hexagon, forming a ring (six calcium ions). A cho formed by the structure, as these rings are responsible for many properties of hydrox atite, especially the biocompatibility and the position of the hydroxyl group in the pl have considerable importance for doping mechanisms [15,50].

Williamson-Hall Method in USDM Model
The W-H method has been used for determining different elastic properties. The procedure is to mathematically reduce the errors and obtain the values of the You modulus by all the diffracted peaks, using the least squares method. The W-H metho a simplified integral expanse and, taking into account the peak width, strain-indu broadening is specified [51]. Taking into account the W-H method in the US model, Young's modulus values are examined. It is clear that in Equation (23) and Fi 5, the term of is along the X-axis and the term of β . cosθ is along the Y-axis i vidually.
β .cosϴ = ( ) + 4σ ϴ Here, β is the broadening peak from (hkl) plane and, in this study, the instrum tal broadening is taken as 0.02 degree for each diffracted peak and is subtracted from values before multiplying to to convert the degree to radian. In this model, the co tion is performed to calculate the strain and the average Young's modulus. The avera value has been calculated in the research and studies on using the W-H method, but subject to errors, because if the average values of the Young's modulus are consid (through Equation (19)), the final value would be far from the standard value of You modulus in each peak extracted by X-Ray diffraction. In addition, in some studies Young's modulus value is considered a value that is listed in the literature but is not a ciated with the prepared materials. As an illustration, based on the study in Refer [10], which refers to the use of the W-H method to calculate the crystal size and Sche analysis of hydroxyapatite, the average value of Young's modulus has a larger devia than the actual value of the Young's modulus of hydroxyapatite. The line broadenin the diffracted peak is extracted with Crystal Diffract 6.7.2.300 software. According to ure 5, the slope values are associated with the stress (σ). The values obtained were p tive, and the positive values of intrinsic strain and stress can prove the tensile stress strain, and if the values were negative, they are associated with compressive stress strain. Additionally, the resulting values of stress (σ) and strain (ε) by utilizing the W method in the USDM model are shown in Table 6. The value of σ is in good agreem with the values obtained from Reference [52].  Here, β hkl is the broadening peak from (hkl) plane and, in this study, the instrumental broadening is taken as 0.02 degree for each diffracted peak and is subtracted from β hkl values before multiplying to π 180 to convert the degree to radian. In this model, the condition is performed to calculate the strain and the average Young's modulus. The average E value has been calculated in the research and studies on using the W-H method, but it is subject to errors, because if the average values of the Young's modulus are considered (through Equation (19)), the final value would be far from the standard value of Young's modulus in each peak extracted by X-ray diffraction. In addition, in some studies the Young's modulus value is considered a value that is listed in the literature but is not associated with the prepared materials. As an illustration, based on the study in Reference [10], which refers to the use of the W-H method to calculate the crystal size and Scherrer analysis of hydroxyapatite, the average value of Young's modulus has a larger deviation than the actual value of the Young's modulus of hydroxyapatite. The line broadening of the diffracted peak is extracted with Crystal Diffract 6.7.2.300 software. According to Figure 5, the slope values are associated with the stress (σ). The values obtained were positive, and the positive values of intrinsic strain and stress can prove the tensile stress and strain, and if the values were negative, they are associated with compressive stress and strain. Additionally, the resulting values of stress (σ) and strain (ε) by utilizing the W-H method in the USDM model are shown in Table 6. The value of σ is in good agreement with the values obtained from Reference [52].

Conclusions
In this study, crystal of hydroxyapatite was successfully synthesized with the thermal treatment process. Moreover, the position of atoms and extracted planar density values of each diffracted plane (32 planes) of hydroxyapatite were determined and calculated for the first time. In addition, a new method based on the relationship between the measurement of elastic modulus and atomic planar density of crystalline hydroxyapatite, consisting of a unit cell and super cells (2 × 2 × 2), was fully presented. According to the method of this study, the slope of modules of elasticity against planar density in super cells (2 × 2 × 2) was negative; the reason was related to imperfections, and this means that when the density of atoms in the planes is increased, lower forces are required for the dislocation motion. As a result, a plane without an atom is more realistic in a smaller area of the unit cell than a larger area of the super cells (2 × 2 × 2), and this means that the intercept of the unit cell is closer to the real values of the Young's modulus in the hydroxyapatite lattice. Moreover, a comparison of theoretical and experimental data of the Young's modulus of hydroxyapatite showed that there is a small difference between the values for both unit cell and super cells (2 × 2 × 2), namely 10.47 GPa and 10.57 GPa, respectively; this means that the theoretical calculation is valid and, by decreasing by about 10 GPa, the experimental value can be obtained. Furthermore, the Young's modulus values of hydroxyapatite in the unit cell and super cells (2 × 2 × 2) were achieved at 108.15 and 121.17 GPa tandemly. Finally, one of the applications of the presented method was carried out in this study and the Williamson-Hall method in the USDM model can be used to minimize the errors in the least squares method and to obtain the correct elastic modulus of hydroxyapatite, which is much more accurate than the average value.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/10 .3390/ma14112949/s1, Figure S1: Images of production route of hydroxyapatite obtained from cow bones, Table S1: Crystallographic parameters of XRD pattern related to the natural hydroxyapatite obtained from cow bones, Figure S2: EDX spectrum and stoichiometric composition of hydroxyapatite obtained from cow bones, Figure S3: (a) Raw hexagonal structure, (b) fundamental axis, (c) principal crystallographic directions and (d) between the orthogonal axes and crystallographic directions, Figure S4: Transfer of hydroxyapatite unite cell in P63/m to P63 space group, Figure S5: Schematic of a design of cutted samples, Table S2: The values of Longitudinal and Transvers velocity of samples, Figure