Computational Modeling of Gold Nanoparticle Interacting with Molecules of Pharmaceutical Interest in Water

We derived a theory of biomolecule binding to the surface of Aun clusters and of the Au plane based on the hard soft acid base (HSAB) principle and the free electron metallic surface model. With the use of quantum mechanical calculations, the chemical potential (μ) and the chemical hardness (η) of the biomolecules are estimated. The effect of the gold is introduced via the empirical value of the gold chemical potential (−5.77 eV) as well as by using the expression (modified here) for the chemical hardness (η). The effect of an aqueous environment is introduced by means of the ligand molecular geometry influenced by the PCM field. This theory allows for a fast and low-cost estimation of binding biomolecules to the AuNPs surface. The predicted binding of thiolated genistein and abiraterone to the gold surface is about 20 kcal/mol. The model of the exchange reaction between these biomolecules and citrates on the Au surface corresponds well with the experimental observations for thiolated abiraterone. Moreover, using a model of the place exchange of linear mercaptohydrocarbons on 12-mercaptododecane acid methyl ester bound to the Au surface, the present results reflect the known relation between exchange energy and the size of the reagents.


Introduction
Gold nanoparticles are readily conjugated to biological molecules, e.g., protein, oligonucleotides, and antibodies, primarily due to the affinity of sulfhydryl (-SH) groups to the gold surface [1].Gold-biomolecule conjugates have been widely incorporated into diagnostic applications, where their bright red color is used in home and point-of-care tests, such as lateral flow assays [2,3].Gold nanoparticle coated with biomolecules can specifically target cancer cells [4].Recently, a considerably large interest in gold particles coated with thiol compounds has been noted.Such systems are being investigated, inter alia, for drug delivery systems in living organisms [3,5].Compounds chemically adsorbed on the surface of the gold nanoparticle core can be thought of as three-dimensional monolayer (3D SAM) systems.The great interest in nanostructures is the result of their chemical, electrical, mechanical and optical properties, which are often different and sometimes more favorable than larger structures [2].Gold nanostructures are gold nanoparticles (a diameter of 1 to 100 nm [6]) as well as nanocrystals ranging in size from one to several hundred nanometers [6,7].
Gold nanoparticles become an intermediate structure between the microscopic size of atoms or molecules and bulk materials.In such systems, one can distinguish a nanoparticle core and a stabilizing compound layer.Modifying nanoparticles with various compounds prevent their aggregation.In the case of gold nanoparticles stabilized with alkanethiols, the stabilizing layer is the alkanethiol monolayer.A specific property of gold nanoparticles modified with an organic monolayer is their ability to accumulate a charge; therefore, they are a material useful in the catalysis of electrode processes.Gold nanoparticles are already used in catalysis, optoelectronics, chemical technology, biosensors, and the pharmaceutical industry as systems for drug delivery in living organisms.In particular, it is worth noting the ease of modification of the gold nanoparticle monolayer, which may allow one to obtain materials with the desired properties.
The S-Au bond energy of gold nanoparticles is estimated to be about 44 ± 3 kcal/mol [8] or 28 kcal/mol (or 117 kJ/mol) [9].It was found that depending on different coverage thiol densities, the Au-S binding can vary, e.g., be stronger at higher densities and weaker at lower densities [10].In the monolayer on a flat surface (2D SAM), each thiol occupies the same area on the cluster.The packing of a thiolated ligand in the planar Au surface is almost twice the coverage of 27 nm AuNPs [11].This can be related to the curvature of the nanoparticle core surface, where gold atoms can interact with more than one sulfur atom of the thiol [12].
There are usually more defects in a nanoparticle monolayer than in a well-organized monolayer adsorbed on a flat (they are, e.g., steps, edges, vertices) [13].The surface of alkanethiol monolayers on 2D and 3D SAM surfaces are better understood compared to thiol-containing drug monolayers [14].
In our recent work, we have shown that the interaction between gold-adsorbed molecules, as well as the molecule binding to the gold surface, can be predicted theoretically using both 2D SAM and 3D SAM surfaces [15,16].The present work is a continuation of this innovative approach to research.The study focused on the modeling of new drug molecules' interaction with gold surfaces.Drugs containing thiol groups immobilized on 2D and 3D SAM surfaces were selected for the present study.It focused on comparing the properties of the Au n -drug complex, in particular the drug-gold surface interaction.
Molecules 2023, 28, x FOR PEER REVIEW 2 of 11 of gold nanoparticles modified with an organic monolayer is their ability to accumulate a charge; therefore, they are a material useful in the catalysis of electrode processes.Gold nanoparticles are already used in catalysis, optoelectronics, chemical technology, biosensors, and the pharmaceutical industry as systems for drug delivery in living organisms.In particular, it is worth noting the ease of modification of the gold nanoparticle monolayer, which may allow one to obtain materials with the desired properties.
The S-Au bond energy of gold nanoparticles is estimated to be about 44 ± 3 kcal/mol [8] or 28 kcal/mol (or 117 kJ/mol) [9].It was found that depending on different coverage thiol densities, the Au-S binding can vary, e.g., be stronger at higher densities and weaker at lower densities [10].In the monolayer on a flat surface (2D SAM), each thiol occupies the same area on the cluster.The packing of a thiolated ligand in the planar Au surface is almost twice the coverage of 27 nm AuNPs [11].This can be related to the curvature of the nanoparticle core surface, where gold atoms can interact with more than one sulfur atom of the thiol [12].
There are usually more defects in a nanoparticle monolayer than in a well-organized monolayer adsorbed on a flat (they are, e.g., steps, edges, vertices) [13].The surface of alkanethiol monolayers on 2D and 3D SAM surfaces are better understood compared to thiol-containing drug monolayers [14].
In our recent work, we have shown that the interaction between gold-adsorbed molecules, as well as the molecule binding to the gold surface, can be predicted theoretically using both 2D SAM and 3D SAM surfaces [15,16].The present work is a continuation of this innovative approach to research.The study focused on the modeling of new drug molecules' interaction with gold surfaces.Drugs containing thiol groups immobilized on 2D and 3D SAM surfaces were selected for the present study.It focused on comparing the properties of the Aun-drug complex, in particular the drug-gold surface interaction.
This paper contains the first theoretical and computational study of the recently synthesized new analog of drugs: "thiogenistein"-thiol derivative of genistein (TGE) [17] denotes the HS-CH2-COO-CH2-CH2-(-O)-genistein, and "thioabiraterone" denotes the HS-CH2-COO-abiraterone [18].The structures are presented in Figure 1.The present model becomes an extension of the former Gazquez et al. model [19,20], includes the water environment affecting the biomolecule-Au surface binding and includes elementary electrostatic concepts to approximate the chemical hardness (η) of the AuNPs.In this model, quantum mechanical calculations are only required for the The present model becomes an extension of the former Gazquez et al. model [19,20], includes the water environment affecting the biomolecule-Au surface binding and includes elementary electrostatic concepts to approximate the chemical hardness (η) of the AuNPs.In this model, quantum mechanical calculations are only required for the prediction of the chemical potential (µ) and the chemical hardness (η) of the biomolecules.The effect of the gold enters into the theoretical equations through the empirical value of the gold chemical potential (−5.77 eV), as well as by using the expression (modified here) for the chemical hardness (η).The effect of the aqueous environment is introduced by means of the ligand molecular geometry influenced by the PCM field [21].With the present modified Gazquez et al. model [19,20], it was possible to determine the upper bounds for the binding energy in the Au prediction of the chemical potential (µ) and the chemical hardness (η) of the biomolecules.The effect of the gold enters into the theoretical equations through the empirical value of the gold chemical potential (−5.77 eV), as well as by using the expression (modified here) for the chemical hardness (η).The effect of the aqueous environment is introduced by means of the ligand molecular geometry influenced by the PCM field [21].With the present modified Gazquez et al. model [19,20], it was possible to determine the upper bounds for the binding energy in the Au∷biomolecule complexes.

Results and Discussion
The present theoretical approach was designed to estimate the binding of molecules of pharmaceutical interest with the gold nanoparticles by means of fast and low-cost theoretical calculations.A goal to attain is a qualitative rather than quantitative assessment of molecule binding.It may help to answer questions from laboratory experimenters interested in having a prediction about whether there can or cannot be binding.An essential quantity of the present theory is the interaction energy between involved species.Here, the interaction energy modeling is based on the Gazquez et al. (1994) approach derived from the Density Functional Theory and the µ-chemical potential, η-hardness, (1/η)-softness concepts and fAi-condensed Fukui's functions of the i-th atom in the Amolecule.The present theory incorporates the free electron metallic surface model and is also modified to include the effect of the water medium.
The µ and η parameters are obtained from the B3LYP calculations with the moderatesize basis sets and uses the three-point finite difference approximation to the first-and second derivative of the total molecular energy.The water effect is introduced implicitly via the PCM model in the course of the molecular geometry optimization.

The Interation Energy
The term ∆Eν is the energy that corresponds to the charge transfer process between A and B species arising from the chemical potential equalization principle at a constant external potential ν [20].
The term ∆Eµ is the energy that corresponds to a reshuffling of the charge distribution, and it is basically a manifestation of the maximum hardness principle [20] and takes place at a constant chemical potential ν [19].
Following Gazquez et al. [19,20], the expressions for the interaction energy in Equation ( 1) can be written as: where Ne [22] can be interpreted as the square of the effective number of valence electrons participating in the interaction between A and B, while fi is the condensed Fukui's function of the i-th of A and is equal to the derivative of the charge qi of the i-th atom with respect to the charge of the molecule q.

The Hardness of the Aun Nanoparticles in the Metallic Sphere Approximation
The chemical hardness η is related to the capacitance C through the equation [23,24]: biomolecule complexes.

Results and Discussion
The present theoretical approach was designed to estimate the binding of molecules of pharmaceutical interest with the gold nanoparticles by means of fast and low-cost theoretical calculations.A goal to attain is a qualitative rather than quantitative assessment of molecule binding.It may help to answer questions from laboratory experimenters interested in having a prediction about whether there can or cannot be binding.An essential quantity of the present theory is the interaction energy between involved species.Here, the interaction energy modeling is based on the Gazquez et al. (1994) approach derived from the Density Functional Theory and the µ-chemical potential, η-hardness, (1/η)-softness concepts and f Ai -condensed Fukui's functions of the i-th atom in the A-molecule.The present theory incorporates the free electron metallic surface model and is also modified to include the effect of the water medium.
The µ and η parameters are obtained from the B3LYP calculations with the moderatesize basis sets and uses the three-point finite difference approximation to the first-and second derivative of the total molecular energy.The water effect is introduced implicitly via the PCM model in the course of the molecular geometry optimization.

The Interation Energy ∆E
The term ∆E ν is the energy that corresponds to the charge transfer process between A and B species arising from the chemical potential equalization principle at a constant external potential ν [20].
The term ∆E µ is the energy that corresponds to a reshuffling of the charge distribution, and it is basically a manifestation of the maximum hardness principle [20] and takes place at a constant chemical potential ν [19].
Following Gazquez et al. [19,20], the expressions for the interaction energy in Equation ( 1) can be written as: where N e [22] can be interpreted as the square of the effective number of valence electrons participating in the interaction between A and B, while f i is the condensed Fukui's function of the i-th of A and is equal to the derivative of the charge q i of the i-th atom with respect to the charge of the molecule q.

The Hardness of the Au n Nanoparticles in the Metallic Sphere Approximation
The chemical hardness η is related to the capacitance C through the equation [23,24]: The factor of 1/2 in the original definition of the global hardness has been omitted here for convenience, and e is the electronic charge [25].
The capacitance C of a spherical capacitor with radius r is given by [25]: From Equations ( 4), it follows that the chemical hardness (η) is: where α is a geometrical correcting factor that takes into account the fact that even if, macroscopically, an AuNP surface is a sphere within a good approximation, the gold atoms are not geometrical points on a surface (e.g., the Au atomic radius is about 0.146 nm).

The Model of the Interaction Energy in Water
One obtains an expression that is a good approximation of the (∆E) Ai : (∆E It follows that the first term in Equation ( 7) represents the interaction energy of the biomolecule with a gold plane and the second term represents a correction for the molecule A interacting through the i-th atom with the Au cluster of the radius r: where (∆E ∞ ) A is the interaction energy with a gold plane and the molecule A, and µ G is the chemical potential of gold (−5.77eV).The total interaction energy ∆E A between the gold nanoparticles Au n and the molecule A can be obtained by replacing the condensed Fukui's function set f Ai with a value of 1 in Equations ( 7) and ( 8); it follows that:

The Model of Reaction Energy in Water
We consider the following exchange chemical reaction between a molecule ligand1 bounded with the gold nanocluster Au n and a molecule ligand2: The reaction energy ∆E can be derived, in this simple model, as the difference between the interaction energy (Au) n -ligand2 and (Au) n -ligand1: Equation ( 10) can be rewritten as: where: ) It follows that we can make the following approximations: More details on the present theory, in particular the derivation of all the expressions, can be found in the Supplementary Information.

Biomolecules' Hardness and Chemical Potential in Water Medium
Since biochemical reactions occur mainly in water solution, taking into account the changes in reactivity indices, going from the gas to solvent phase permits the elucidation of the actual reactivity of the compounds within DFT for neutral and charged systems using the polarizable continuum model (PCM).The Taylor expansion of the total electronic energy E(∆N) as a function of the number of extra electrons ∆N [26,27] is: where, E 0 is the energy of the neutral species, µ is the chemical potential, and η is the chemical hardness.It follows that: Summing and subtracting Equations ( 14) and ( 15), it follows that: The calculated values, in water medium, of the chemical potential µ and of the chemical hardness η of different ligands are reported in Table 1.

Ligand Exchange Energies
Based on the data included in Table 1, one can predict the energy output ∆E of the hypothetical ligand exchange reaction (R1) according to Equations ( 11) and (12).The reaction energy can be simply calculated using the following Equation: the reaction (R1) is energetically favorable if, in Table 2, the row of the ligand2 is above the row of the ligand1.The energy values in Table 2 are ordered by increasing energies values.
Based on this result, one can predict that a hypothetical ligand exchange will occur when the citrates are doubly negatively charged.Such a result, coming from two independent theories, corresponds with the experimentally observed ligand exchange in the course of the formation of the AuNPs dressed with "thioabiraterone" [16].[19,20].Energies in kcal/mol, distances in Angstroms.

Model for Small Radius of the AuNPs
The total interaction energy ∆E: between the gold nanoparticles Au n , the ligand can be obtained by replacing the condensed Fukui's function f Ai with a value of 1 in Equation (2); it follows that: e e 2 8πε 0 α r r e 2 4πε 0 αrη A + 1 (19) In order to evaluate the interaction energy ∆E 0 between the gold nanoparticles Au n and the ligand for a small number n of Au n gold atoms, taking into account Equation (17), we take the limit for the radius r of the AuNP that tends to zero, and it follows that: The lower limit of ∆E in Equation ( 20) corresponds to r = 0 and becomes the right-hand side of Equation (19).The upper limit of ∆E corresponds to r = ∞ and becomes the first term in the nominator of Equation (19).The dependence of ∆E on r in Equation ( 19) is illustrated in Figure 2.
The calculations according to Equations ( 19) and ( 20) correspond to the case with solvent.The solvent, water in our case, implicitly enters into Equations (19) and (20) through the modification of the molecular geometry under the action of the PCM molecular field.The molecular geometry without water is somewhat different from the geometry under the PCM field; therefore, the calculated (with water) chemical potential µ A and chemical hardness η A , as well as Equations ( 19) and (20), lead to somewhat different results in comparison to the case when water was absent.
first term in the nominator of Equation (19).The dependence of ΔE on r in Equation ( 19) is illustrated in Figure 2. The calculations according to Equations ( 19) and ( 20) correspond to the case with solvent.The solvent, water in our case, implicitly enters into Equations ( 19) and ( 20) through the modification of the molecular geometry under the action of the PCM molecular field.The molecular geometry without water is somewhat different from the geometry under the PCM field; therefore, the calculated (with water) chemical potential µA and chemical hardness ηA, as well as Equations ( 19) and (20), lead to somewhat different results in comparison to the case when water was absent.
The theoretical (ΔE(∞))ligand binding between the ligand and an infinite gold plane is obtained by using Equation (9).
The energy output of the reaction (R1) is favorable if ΔE(∞) in the row of Table 2 corresponding to the ligand2 is more negative than ΔE(∞) in the row with the ligand1.Consequently, the complex Au-ligand2 is energetically favorable in comparison to the complex Au-ligand1; in other words, ligand1 can be replaced by ligand2.Let us consider three examples of the ΔE citrate exchange on the biomolecules using the binding energies from the Table 2  The theoretical (∆E(∞)) ligand binding between the ligand and an infinite gold plane is obtained by using Equation (9).
The energy output of the reaction (R1) is favorable if ∆E(∞) in the row of Table 2 corresponding to the ligand2 is more negative than ∆E(∞) in the row with the ligand1.Consequently, the complex Au-ligand2 is energetically favorable in comparison to the complex Au-ligand1; in other words, ligand1 can be replaced by ligand2.Let us consider three examples of the ∆E citrate exchange on the biomolecules using the binding energies from the Table 2 The ligand replacement should also occur, as was evidenced experimentally by Stolarczyk et al. [16].
The prediction of the citrate exchange of citrate on thiogenistein or thioethanegenistein on the AuNP does not yet have an experimental counterpart in our laboratory.
It is worth noting that the citrate replacement should not occur in cases where the citrate is triply ionized, i.e., when the citrate bears the (−3) charge.These results correspond to the experimental observation of a non-complete displacement of citrate capping the AuNPs by thiolated ligands [23,29].
One can also note that in some cases (e.g., for acetate clusters) presented in Table 2, the DFT results for large Au n : acetate clusters [9] are closer to ∆E(∞) than to ∆E(r), which suggests that finite Au clusters (for a small n) seem to be a poor representation of the biomolecule binding to the gold surface.
The experimental estimations of the Au-S(thiol) binding of about 45 kcal/mol [24], 28 kcal/mol [9] and 30 kcal/mol for Au-S-CH 3 [21] are not far from the estimated Au-S bond energy of about 30 kcal/mol from our models of Au:ethanethiol, Au:butanethiol and Au:octanethiol complexes (c.f.Table 2).
Another interesting observation from the data shown in Table 2 is the place exchange of the AuNP-thioalkyl conjugates, e.g., [30].Based on the ∆E(∞) values from Table 2, one can predict that the thioalkyls with a longer CH 2 chain can replace the thioalkyls with shorter CH 2 chains, e.g., octanethiols (∆E(∞) = −30.9kcal/mol) can replace ethanethiols (∆E(∞) = −29.4kcal/mol) in AuNP conjugates.A similar trend is predicted for the place exchange of linear mercaptohydrocarbons on 12-mercaptododecane acid methyl ester on the Au surface.These predictions follow the trend of experimental results [30,31].More details are provided in the Supplementary Information.

Computational Methods
We used the quantum mechanical approach using the hybrid Density Functional B3LYP method with the split valence basis set and d, p polarization functions 6-311G(d,p) (DFT:B3LYP/6-311G(d,p)) [32] for calculations of the ligands' µ and η values.
The geometry of the models of biomolecules have been DFT: B3LYP/6-311G(d,p)optimized.The harmonic frequencies were calculated using Gaussian 16 software [33], taking into account the solute-solvent interactions within the DFT, using the polarizable continuum model (PCM) [34] that is the default in Gaussian 16 software [33].The water medium has been chosen due to its importance in biochemistry.

Conclusions
Although approximate in nature, the results of Au-biomolecule binding energies in water medium, derived here on the hard soft acid base (HSAB) principle and the free electron metallic surface model, are comparable to the more advanced quantum mechanical calculations for small Au n ::biomolecule systems that are available.The present theory of biomolecule binding to the surface of the Au n clusters and of the Au plane is a reasonable tool for the fast and low-cost estimation of problems that are hard to investigate experimentally.The present theory predicts a binding of about 20 kcal/mol for thiolated genistein and thiolated abiraterone with the gold surface.The predicted exchange of citrate on thiolated abiraterone occurring on the Au surface corresponds to the experimental observations.The predicted exchanges of citrate on thiogenistein and thioethanegenisten still anticipate an experimental counterpart.The predicted place exchange of linear mercaptohydrocarbons on 12-mercaptododecane acid methyl ester on the Au surface corresponds reasonably well with the known relation linking a more efficient reaction with compounds possessing a larger number of carbon atoms.

( 1 )
∆E(r = ∞) (upper limit, also denoted as ∆E(∞)) corresponds to the binding of the AuNP(planar)-biomolecule complexes, Equation(9); ∆E(r = 2.88 Å) corresponds to the binding of AuNP-biomolecule complexes for r equal to the closest Au-Au distance in the Au(111) plane, Equation The effect of water was simulated with the PCM field modifying the molecular geometry in the course of its optimization as implemented in the Gaussian 16 suite of programs.(2)"thioethanegenistein" denotes the HS-CH 2 -CH 2 -COO-CH 2 -CH 2 -(-O)-genistein; "thioabiraterone" denotes the HS-CH 2 -COO-abiraterone; "thiogenistein" denotes the HS-CH 2 -COO-CH 2 -CH 2 -(-O)-genistein; citrates code: both external carboxyl groups ionized: (-H-), one external and the central carboxyl group ionized: (--H), one external carboxyl group ionized: (-HH), the central carboxyl group ionized: (H-H).(3)The DFT calculations for the Au 25 ::biomolecule complexes were performed with the DFT/6-31G(d,p) method supplemented with the lanl2dz effective potential for the Au atoms.The fixed planar geometry of Au25 cluster (5 × 5) corresponds to the fragment of the Au(111) plane with the closest Au-Au distance of r = 2.88 Å.The optimization of biomolecule geometry was carried out with the Berny's algorithm implemented in the Gaussian G16 software until the geometry plateau was reached.The interaction energy was estimated by subtracting the total energy of the monomers (Au 25 , biomolecule) from the energy of the dimer (Au 25 ::biomolecule).Two positions of biomolecule relative to the Au 25 plane were studied, i.e., thiol-moiety close to the edge Au atom and almost perpendicular to the Au 25 plane.Another DFT calculations of Al-Johani et al. [28] include acetate interacting with the 110 Au atoms forming three layers of 5 × 5, 6 × 6, and 7 × 7 atoms.a biomolecule bound to the Au 25 plane edge; b biomolecule bound perpendicular to the Au 25 plane; c B3LYP/6-31G(d,p) and lanl2dz for (thioabiraterone) 2 -Au 7 estimation from ref. [16]; d model of acetate binding to Au 110 [28].

Figure 2 .
Figure 2. The plot of ΔE(r) in Equation (19) for the model of AuNP interacting with thioabiraterone.ΔE in kcal/mol, r in Angstroms.

Figure 2 .
Figure 2. The plot of ∆E(r) in Equation (19) for the model of AuNP interacting with thioabiraterone.∆E in kcal/mol, r in Angstroms.

Table 2 .
The predicted AuNP-biomolecule binding energies in water, ∆E (r), following Gazquez et al.'s theory, modified here