Ab-Initio Predictions of the Energy Harvesting Performance of L-Arginine and L-Valine Single Crystals

Biological piezoelectric materials are beginning to gain attention for their huge potential as eco-friendly energy harvesting materials. In particular, simple amino acid and peptide crystal assemblies are demonstrating large voltage outputs under applied force, and high sensitivity when detecting vibrations. Here we utilise Density Functional Theory (DFT) calculations to quantitatively predict the energy harvesting properties of two understudied proteinogenic amino acid crystals: L-Arginine and L-Valine. The work highlights the ability of quantum mechanical calculations to screen crystals as high-performance energy harvesters, and demonstrates the capability of small biological crystals as eco-friendly piezoelectric materials. L-Arginine is predicted to have a maximum piezoelectric voltage constant of gij = 274 mV m/N, with a Young’s Modulus of E = 17.1 GPa. L-Valine has a maximum predicted piezoelectric voltage constant of gij = 62 mV m/N, with a calculated Young’s Modulus of E = 19.8 GPa.

The more natural application avenue for biomolecular crystal piezoelectrics is in energy harvesting. For years this field has successfully utilised many techniques to modulate the dielectric constants of inorganic piezoelectric materials and increase their energy harvesting performance (Trolier-McKinstry et al., 2011;Roscow et al., 2019;Wang et al., 2021). The lower the dielectric constant of a material the higher its piezoelectric voltage constant g ij , and the higher the electric field generated per unit force. Biomolecular crystals have an average dielectric constant of 3, two orders of magnitude lower than ceramics (Zhang et al., 2007;Guerin et al., 2019). This means that even "weak" organic piezoelectrics can outperform ceramics and polymers in both output per unit force and sensitivity. We have recently demonstrated this using a flexible glycine-based energy harvester that can be used as a structural health monitor to detect leaks in water infrastructure networks (Okosun et al., 2021). Amino acids have also been utilised to detect ultra-low mechanical pressures (Bishara et al., 2020).
In this work, we add to the database of biomolecular crystals that can be used as low-cost, eco-friendly, high-performance energy harvesters, by predicting the electromechanical properties of the amino acids L-Arginine and L-Valine. We report the elastic, piezoelectric, and dielectric properties of these crystals, and rationalise their anisotropic response via supramolecular packing analysis. This methodology can be used to screen biomolecular crystals such as these for desirable electromechanical properties, which can then be experimentally validated using techniques such as impedance spectroscopy, scanning probe microscopy, nanoindentation, or the Berlincourt method.
L-Arginine is one of the most underutilised amino acids, with the molecule being used as a ligand in 1D-copper (II) coordination polymers (Alikhani et al., 2020). Arginine has been crystallized as L-arginine 4-nitrophenolate 4-nitrophenol dehydrate (LAPP) (Wang et al., 2011), L-arginine hydrofluoride (Pal and Kar, 2002), and L-arginine hydrochloride monohydrate (Kalaiselvi et al., 2008) for non-linear optical applications, as has L-Valine (Moitra and Kar, 2010). These biomolecular-crystal assemblies can be grown at room temperature with no byproducts, and do not require an external electric field to induce piezoelectricity, unlike PZT and other piezoceramics. Both the raw material and fabrication cost of arginine and valine are a tiny fraction of that of current commercial piezoelectrics which rely on heavy processing of heavy metals and their oxides, which carries a large financial and environmental burden (Ibn-Mohammed et al., 2018).

RESULTS AND DISCUSSION
L-Arginine crystallises in the monoclinic space group P2 1 , allowing for eight non-zero piezoelectric tensor components.
The predicted piezoelectric charge constants (Table 1), denoted e ij , range from a minimum of e 14 0.0041 C/m 2 to a maximum of e 22 0.1031 C/m 2 . Dividing the charge tensor by the calculated elastic stiffness tensor ( Table 2) gives the piezoelectric strain constants d ij (Table 3), which for Arginine range from d 14 0.28 pC/N to d 16 6.52 pC/N. The elastic stiffness tensor itself shows high elastic anisotropy, with increased stability along the crystallographic a axis (c 11 29 GPa, c 44 14 GPa). The bulk Young's Modulus is predicted to be 17.1 GPa.
As mentioned, the key energy harvesting figures of merit (FoMs) are the piezoelectric voltage constants g ij ( Table 4). The predicted bulk dielectric constant of L-Arginine is 2.54, resulting in predicted voltage constants ranging from 12 mV m/N to 274 mV m/N. L-Arginine is predicted to have five voltage constants that exceed both the 40 mV/m N value for  KNN ceramics (Zhang et al., 2007), and the 64 mV m/N of organic polymer PVDF (Bernard et al., 2017). L-Valine is predicted to have a low piezoelectric response compared to other biological crystals. It too crystallises in the monoclinic P2 1 space group, with eight piezoelectric constants in each tensor. The piezoelectric charge constants of valine range from 0.0024 to 0.0151 C/m 2 . It is unique amongst biomolecular crystals in that the predicted electronic and ionic contributions to the piezoelectric tensor are both high but of opposite polarity, resulting in these reduced values. This also explains why valine has a low piezoelectric response, but high non-linear optical performance (Moitra and Kar, 2010), as we have recently demonstrated that the second harmonic generation capabilities of small biomolecular crystals correlates to the electronic contribution to the piezoelectric tensor only (Gleeson et al., 2020).
L-Valine is predicted to have lower elastic anisotropy than L-Arginine, but does demonstrate increased longitudinal stiffness along the crystallographic c axis (c 33 46.5 GPa). Its piezoelectric strain constants range from d 23 0.05 pC/N to d 34 1.8 pC/N, of similar value to quartz crystals. Still, we see that the predicted dielectric constant of 3.31 results in four piezoelectric voltage constants that exceed that of PZT (25 mV m/N), though they are all shear, i.e. voltage can only be produced under application of a shearing force. The predicted piezoelectric voltage constants of L-Valine range from 2 to 62 mV m/N.
Looking at the supramolecular packing of the L-Arginine unit cell (Figure 1), we can see that the high stiffness along the crystallographic a axis is due to both continuous directional hydrogen bonding and the molecular backbone, which runs perpendicular to the axis increasing mechanical stability. By symmetry and the presence of a 22 piezoelectric constant, the net dipole in the crystal will be predominantly in line with the b axis, which also contains a continuous directional hydrogen bond network. The zwitterionic molecules create this 2D hydrogen bond network in the ab plane, enhancing the piezoelectric polarisation both in equilibrium and under an applied longitudinal (along b only) or shear force.
L-Arginine has many features of a moderately piezoelectric molecular crystal, with medium to high density, small functional groups, low molecular weight, and a monoclinic angle of 97.4°. Comparing the response to other materials we see that the maximum piezoelectric strain constant of both crystals is low compared to many inorganic single crystals and ceramics ( Figure 2).  However Figure 3 shows that both crystals dramatically outperform ceramics in their energy harvesting responses, with the maximum piezoelectric voltage constants of L-Valine matching those of organic polymers PVDF and PLLA.
Examining the crystal structure of L-Valine (Figure 4), it can be seen that the low piezoelectric response is due to inhibited hydrogen bonding along the b axis (the only axis that can allow a longitudinal piezoelectric response by symmetry). The molecules pack as distinct sheets, with two layers of molecules forming hydrogen bonds along the a and c axes between NH 2 + and COOtermini. However the molecular sheets are separated by the valine methyl side chains, which cannot form hydrogen bonds and thus disrupt the formation of a long-range hydrogen bonding network that would increase the net polarisation in the unit cell and thus the piezoelectric response and energy harvesting performance.
Similarly to the L-Arginine single crystal, maximum mechanical stability is observed in the direction of the continuous hydrogen bond network that does form along the c axis, as well as perpendicular to the molecular backbone. The L-Valine unit cell also has a lower monoclinic angle than L-Arginine of 96°, which results in increased shear stiffness (c 55 c 66 5.3 GPa) and reduces the amount of ionic displacement that can occur under an applied force. This is also reflected in the larger predicted Young's Modulus of L-Valine of 19.8 GPa. Figure 5 shows that despite the low response of valine compared to the whole spectrum of piezoelectrics, the responses of both valine and arginine are average for biological materials overall, and far exceed the piezoelectric response of materials such as bone, tendon, calcite, wool, wood, and so on.

CONCLUSION
In this work, it is shown that both L-Arginine and L-Valine single crystals are predicted to have piezoelectric voltage constants that exceed those of piezoelectric ceramics, with L-Arginine having voltage constants that also far exceed piezoelectric polymers. This highlights that both crystals would be suitable as energy harvesting materials, with substantially reduced cost and environmental impact than currently used candidates. The higher response of L-Arginine is attributed to its continuous hydrogen bond network along the axis of polarisation and larger monoclinic angle, while the lower response of L-Valine is caused by inhibited long-range hydrogen bonding and higher shear stiffness. The predicted  elastic constants indicate that these amino acids could be used as rigid single crystal components, or as flexible thin films for a variety of energy harvesting applications.

COMPUTATIONAL DETAILS
Classical piezoelectricity manifests itself in both a direct and a converse effect. The direct piezoelectric effect is a linear coupling between mechanical stress and electrical polarization. The converse piezoelectric effect is a linear coupling between mechanical strain and an applied electric field. Both behaviours are described by the same set of piezoelectric constants, most commonly the piezoelectric strain constant d ij , which is measured in pC/N for the direct piezoelectric effect, and pm/V for the converse piezoelectric effect. Mathematically, this piezoelectric response can be described by a third rank tensor in the form of a 3 × 6 matrix:  (Argaman and Makov, 2000) calculations on single crystals using the VASP (Hafner, 2007) code. For a full overview on DFT readers are directed to more thorough overviews (Fiolhais et al., 2003;Koch and Holthausen, 2015). A number of multiscale modelling techniques can be used to evaluate and simulate electromechanical properties (Sass et al., 2004), and within DFT other suitable softwares include CP2K (Kühne et al., 2020), ABINIT (Gonze et al., 2016), and CASTEP (Clark et al., 2005). Alternative methods to DFT for prediction and calculation of material properties without requiring higherorder parameters include deterministic artificial intelligence, which has recently been proposed for material property utilization in control (Sands, 2020). Whiplash compensation, as applied to flexible space robotics has also been used for stiffness estimation (Sands, 2019), stemming from a similar problem formulation utilizing the Hamiltonian-akin to the many-electron time-independent Schrödinger equation used for piezoelectric constant estimation in DFT.
Electronic structures were calculated using the PBE functional (Perdew et al., 1992) with Grimme-D3 dispersion corrections (Grimme et al., 2010) and projector augmented wave (PAW) pseudopotentials (Kresse and Joubert, 1999). Calculations were carries out using Gaussian smearing, and a plane wave cut-off of 600 eV. Piezoelectric charge constants, e ij, were calculated using density functional perturbation theory (DFPT) (Wu et al., 2005). For this a 2 × 2 × 2 Γ-centred k-point grid was also used, The theoretical methodology herewith consists of four steps: -Geometry optimisation -Calculation of piezoelectric charge tensor -Calculation of elastic constants -Calculation of static dielectric tensor

Geometry Optimisation
The first step is to take the experimental crystal structures of interest, and allow their lattice parameters and atomic positions to relax in order to obtain ideal ground state structures. Experimental structures were downloaded from the Cambridge Crystallographic Data Centre (CCDC). All crystal structures were optimised using the conjugate gradient algorithm (Štich et al., 1989). For a periodic system, integrals in real space over the infinitely extended system are replaced by integrals over the finite first Brillouin zone in reciprocal space, in accordance with Bloch's theorem. These integrals are performed at a finite number of points in the Brillouin zone, called the k-point mesh, or grid. A 4 × 4 × 4 Γ-centred k-point grid was used for geometry optimisations; with a plane wave energy cutoff of 600 eV (all plane waves with a kinetic energy less than this are included in the basis set). Γ-centred k-point grids are generally recommended for non-centrosymmetric unit cells. These values were obtained after energy convergence tests with respect to N (where N number of subdivisions of the reciprocal lattice vector), and plane wave cut off energy. Minimum converged are chosen to balance calculation accuracy with computational expense.

Elastic Constants
The elastic stiffness constants c kj , are required to calculate the piezoelectric strain constants d ik . The elastic constants are also important experimentally. The elastic compliance can easily be derived from the stiffness and measured using impedance spectroscopy or nanoindentation, as can the Young's Modulus, which is important for device applications. Using the piezoelectric charge coefficients, e ij , which are calculated directly by VASP, and the elastic stiffness constants, c kj , we can calculate the more useful piezoelectric strain coefficient, d ik , using the relationship The elastic constants are calculated in the form of the stiffness tensor, C, presented as a 6 × 6 matrix. VASP outputs this tensor in Voigt notation in kB, so post analysis is required to produce C in matrix notation in GPa.
For this work, it is only necessary to extract the six primary diagonal matrix components, shown in bold in Eq. (2). The number of non-zero elements in both the elastic and piezoelectric matrices will vary according to the symmetry of the crystal being studied. Young's Moduli were derived from the stiffness and its inverse compliance matrix components. Values are presented as Voigt-Reuss-Hill averages (Chung and Buessem, 1967;Zuo et al., 1992) as calculated by the ELATE software tool (Gaillac et al., 2016). Crystal structures were visualised using VESTA (Momma and Izumi, 2011). Currently, VASP only supports the calculation of the ionic contribution to the elastic stiffness constants using DFPT, and so cannot be used to derive the full elastic stiffness tensor. For more information on these constants, their matrix representations and their units, readers are directed to Nye (Nye, 1985).

Dielectric Tensor
DFPT calculations can also predict the static dielectric tensor of the crystal being studied. Using the dielectric tensor we can extract the final piezoelectric constant: the voltage constant, g ik . This is an important figure of merit (FoM) for energy harvesting applications, and in motion and pressure sensing. To obtain these g ik values we divide the corresponding piezoelectric strain constant, d ik, by the relevant dielectric constant ε ii , as shown in Eq. (3). These constants are measured in V m/N.
The accuracy of the computational methods used in the simulation of permittivity, elastic stiffness, and piezoelectricity, has been extensively benchmarked and validated in previous publications. Initially the methodology was benchmarked with respect to three well-known inorganic piezoelectric materials; namely aluminium nitride (AlN), zinc oxide (ZnO) and α-quartz (SiO 2 ). This was then extended to the proteinogenic amino acids (Guerin et al., 2018a), biominerals (Guerin et al., 2018b), co-crystals (Ji et al., 2020) and peptides (Bera et al., 2021), with deviations from experiment ranging from 1-20%, which is highly accurate for identifying highperformance materials. The upper limit is observed in highly flexible materials with individual stiffness constants of less than 5 GPa.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.