# PiezoTensorNet: Crystallography informed multi-scale hierarchical machine learning model for rapid piezoelectric performance finetuning

Sachin Poudel, Rubi Thapa, Rabin Basnet, Anna Timofiejczuk, Anil Kunwar

► To cite this version:

Sachin Poudel, Rubi Thapa, Rabin Basnet, Anna Timofiejczuk, Anil Kunwar. PiezoTensorNet: Crystallography informed multi-scale hierarchical machine learning model for rapid piezoelectric performance finetuning. Applied Energy, 2024, 361, pp.122901. 10.1016/j.apenergy.2024.122901. hal-04508390

HAL Id: hal-04508390

https://hal.science/hal-04508390v1

Submitted on 18 Mar 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
Graphical Abstract

# PiezoTensorNet: Crystallography informed multi-scale hierarchical machine learning model for rapid piezoelectric performance finetuning

Sachin Poudel,Rubi Thapa,Rabin Basnet,Anna Timofiejczuk,Anil Kunwar

<figure><img src="image_2.png" /></figure>
## Highlights

# PiezoTensorNet: Crystallography informed multi-scale hierarchical machine learning model for rapid piezoelectric performance finetuning

Sachin Poudel,Rubi Thapa,Rabin Basnet,Anna Timofiejczuk,Anil Kunwar

* Crystallography informed machine learning models.

* Composition-based feature vector for allowing the entry of dopant elements of different amount.

* PiezoTensorNet able to design materials with quasi-effective piezoelectric coefficients.

* Crystal-rotated B₀.₃Er₀.₅Al₀.₂N alloy has larger piezoelectric response than AlN ceramic.

* **PiezoTensorNet** web app released online.
# PiezoTensorNet: Crystallography informed multi-scale hierarchical machine learning model for rapid piezoelectric performance finetuning

Sachin Poudel<sup>a,*</sup>, Rubi Thapa<sup>a</sup>, Rabin Basnet<sup>b</sup>, Anna Timofiejczuk<sup>a</sup> and Anil Kunwar<sup>a</sup>

## ARTICLE INFO

**Keywords:**

Piezoelectric tensors

Crystal rotation

Energy conversion

Quasi-effective coefficients

Electromechanical finetuning

## ABSTRACT

Piezoelectric devices offer numerous opportunities for sustainable harvesting of wasted mechanical energy, leading to a significant interest in data-driven research on these materials. This study presents the design of PiezoTensorNet, a comprehensive framework that encompasses a hierarchical classification neural network for crystal point group determination and modular ensembles of regression-based multi-dimensional models for predicting piezoelectric tensors. The machine learning models capable of forecasting tensors for dopant element alloying and various crystallographic transformations is integrated along with finite element analysis for electromechanical performance evaluation. The efficacy of integrated toolkit is demonstrated through the computational design and discovery of a lead-free microelectromechanical system based on AlN. The introduction of Boron and Erbium dopants in AlN enhances its piezoelectric performance, particularly when the crystal undergoes rotations along a preferred axis. Specifically, under a vertical loading of $5 \times 10^{-5}$ N/m$^2$ applied to a cantilever beam, the preferentially oriented $B_{0.3}Er_{0.5}Al_{0.2}N$ material generates a power 9.96 times larger than that of AlN ceramics.

## 1. Introduction

Energy harvesting has gained significant attention as a crucial technology for applications including but not limited to sensors and actuators [1], Internet of Things (IoT) [2], and e-health [3, 4]. Main types of micro and nano energy harvesters include piezoelectric, triboelectric, electromagnetic, electrostatic, and acoustic [5, 6]. Piezoelectric materials, including single crystals, ceramics, composites, and films, are indispensable components in electromechanical devices [7]. These functional materials possess the capability for efficient interconversion of mechanical and electrical energy due to the linear coupling between mechanical stress and electric polarization [8]. They manifest the direct piezoelectric effect, where mechanical stress induces alterations in surface charge density, as well as the inverse piezoelectric effect, where electric fields elicit material deformation [8]. As polarization represents a vector and stress corresponds to a second-rank tensor, the coupling phenomenon encompasses three distinct directions, thereby giving rise to a third-rank piezoelectric tensor [9]. In a Euclidean crystal reference frame because of the symmetricity, the piezoelectric tensor is typically represented as a $3 \times 6$ matrix [10]. The concise overview of the piezoelectric coefficients, their relations with electric and mechanical terms, and corresponding expressions are provided in Supplementary Information S1.

The functional characteristics of these materials are directly linked with their non-centrosymmetric structures [11]. Perovskite-structured piezoelectric materials, including barium titanate, lead titanate, and lead zirconate titanate, show

significant promise for a wide range of practical and industrial applications [12]. The European Union has imposed stringent restrictions on the utilization of lead-based materials, thereby instigating scientific communities to seek lead-free piezoelectric alternatives [13] that not only demonstrate environmental compatibility but also exhibit exceptional characteristics for the advancement of novel devices [14]. Unfortunately, the identification of suitable Pb-free substitutes, capable of fulfilling the requisite criteria, continues to elude the current state of research [15].

Several studies have made endeavors to apply machine learning techniques in the pursuit of high piezoelectric coefficient lead-free materials, showcasing their potential in this field [16, 17, 18]. These works have successfully improved piezoelectric performance through the identification of novel materials. However, a limitation of these studies is that they solely rely only on specific components of the piezoelectric tensor. Piezoelectric ceramics, characterized by $\infty mm$ symmetry [19], typically possess limited independent piezoelectric coefficients, namely longitudinal mode $d_{33}$, transversal mode $d_{31}$, and thickness-shear mode $d_{15}$ [20]. It is important to note that only these three coefficients may not fully capture the polarization distribution, which can be accurately represented by the entire tensor. The tensor encompasses contributions from both ionic and electronic sources, providing a more comprehensive understanding of the material's behavior [21, 22]. Moreover, studying crystal rotation in such functional materials is challenging without prior knowledge of the complete tensor elements [23].

The presence of a non-zero piezoelectric tensor in materials leads to a high electro-mechanical coupling polarization, which results in an increased figure of merit [24]. By introducing quasi-symmetry breaking in metamaterial design, it becomes possible to induce unconventional piezoelectric

<figure><img src="image_3.png" /></figure>

<footnote><sup>a</sup>Faculty of Mechanical Engineering, Silesian University of Technology, Konarskiego 18A, 44-100 Gliwice, Poland</footnote>

<footnote><sup>b</sup>School of Engineering, The Australian National University, 2601, Canberra, Australia</footnote>

<footnote>*Corresponding author</footnote>

<footnote>sachin.poudel@polsl.pl (S. Poudel)</footnote>

<footnote>ORCID(s): 0000-0003-1427-5550 (S. Poudel); 0000-0002-8973-331X (R. Thapa); 0000-0002-0406-6918 (R. Basnet); 0000-0003-2941-4955 (A. Timofiejczuk); 0000-0003-4295-5772 (A. Kunwar)</footnote>
coefficients in perovskite ceramics. These unconventional coefficients, referred to as "quasi-effective piezoelectric coefficients," when combined with conventional coefficients, form a 3 × 6 matrix with technically all elements being non-zero. This expanded coefficient matrix opens up a broader range of possibilities for developing devices based on piezoelectric ceramics [25, 26]. Hence, there is a need to develop methods that can predict piezoelectric tensors, enabling a more thorough investigation of polarization and facilitating the study of crystal rotation effects and their advancement in piezoelectric performance.

This study integrates alloying of base piezoelectric materials with crystallographic-supported data-driven computational modeling, incorporating the physics of crystal rotations based on tensor mathematics. This comprehensive approach offers a rapid screening tool for evaluating the piezoelectric performance of materials by establishing a link between the microscopic crystal structure and resulting piezoelectric polarization. The effectiveness of this approach is demonstrated using AlN-based materials in various scenarios, including pure AlN, single and paired dopant alloying, and crystal rotations along three Cartesian coordinate axes. Energy generation in these scenarios is investigated, and the machine learning-supported finite element method is employed for energy harvesting modeling, aiming to design materials that exhibit high voltage generation under similar loading conditions. The findings showcase the extensive capabilities of this integrated approach in piezoelectric material design and engineering applications.

# 2. Methodology

## 2.1. Data Acquisition, processing and Analytics

The piezoelectric properties of training materials consisting of 3301 were extracted from the Materials Project (MP) [27, 28]. A dataset of 1354 piezoelectric materials was used for machine learning model formulation. That harvested dataset consists of piezoelectric properties for various crystal structures. Selected materials comprising 367 cubic (23, 42m), 152 tetragonal 42m, 305 orthorhombic 222, and 195 hexagonal 6mm + tetragonal 4mm, and 335 orthorhombic mm2 piezoelectric materials were undertaken for further data analysis and machine learning model training purpose. Numerical featurization of physical properties played a pivotal role in the data-driven study. Therefore, the assessment of numerical descriptors for structural, physical, and thermo-mechanical characterization carries crucial in data processing and analytics.

Featurization involved transforming materials-related quantities into physically relevant numerical descriptors, capturing complex relationships between inputs and targeted piezoelectric tensors. Initial features were derived from chemical composition formulas using Materials Agnostic Platform for Informatics and Exploration (Magpie) [29]. While the Magpie featurization package includes 132 features encompassing elemental properties and stoichiometric attributes, it does not include thermodynamics, elastic, and

the structural properties described in [30]. Thus, Magpie featurization is complemented by 13 additional descriptors to account for missing descriptors.

Principal Component Analysis (PCA) was employed for dimensionality reduction, compressing large feature sets to streamline analysis and visualization. By decomposing the data matrix into the direction and magnitude components, PCA identified multiple data directions, eliminating less relevant features. The retained information's importance was quantified using explained variance (EV). Incorporating PCA-based feature selection enhanced analysis efficiency and interpretability, focusing on impactful data aspects [31]. Detailed information regarding each component of the harvested data, including its visualization, numerical feature calculations, and feature engineering, can be found in Supplementary Information S3.

## 2.2. Crystal point group clustering based on piezoelectric tensor symmetry

The lack of inversion symmetry in a crystal structure can give birth to a variety of unique characteristics including piezoelectricity [32]. Neumann's symmetry principle elucidates the profound interplay between the intrinsic symmetry of a material's structure and the resulting constraints imposed on tensors (a reference to the number of independent components) characterizing its physical properties. This fundamental principle asserts that the properties of a material cannot possess a symmetry lower than that exhibited by its crystal structure. The choice of crystallographic unit cell axes aligns with the principal symmetry elements, and basis coordinate systems for tensors are established accordingly. This alignment simplifies tensor descriptions by reducing the number of independent components. Neumann's principle highlights the interplay between crystal symmetry and material properties [33].

Crystal symmetry classes categorize the symmetry elements and operations of crystals, enabling analysis of their properties including impact on piezoelectric behavior [34]. Among 32 unique crystal classes, 20 point groups exhibit piezoelectricity [19]. The higher symmetry crystals such as cubic structures display simplified piezoelectric tensors due to their symmetricity. This higher symmetry results in a reduction in the number of independent non-zero piezoelectric coefficients [35]. Conversely, crystals with lower symmetry, as found in trigonal, tetragonal, or orthorhombic crystal systems, exhibit more intricate piezoelectric tensor forms. The number of independent piezoelectric coefficients is determined by the crystal's symmetry [36]. Supplementary Information S2 provides further elaboration on the relationship between crystallography and piezoelectricity.

While symmetry considerations can identify the potential for piezoelectricity, their ability to fully capture the magnitude and nature of the piezoelectric response is limited [37]. Incorporating all 20 point groups into a machine learning model is not feasible first because of the classification of higher nodal targets and higher dimensions of
<caption>**Table 1**</caption>

**Crystallographic Impact on Piezoelectric Tensors:** This table explores the relationship between crystallography and piezoelectric tensors through symmetry constraints and graphical representations. The category class serves as a pseudo-grouping mechanism in the hierarchical classification machine learning model, enhancing the efficiency of the parent classification process. Similar tensor structures cluster point groups into PG1, PG2, PG3, PG4, and PG5. The symbol $\#e_{ij}$ denotes the number of independent piezoelectric coefficients. Notably, the cubic (23, 42m) and hexagonal 6mm + tetragonal 4mm point groups share the same tensor structure, resulting in seven point groups represented by five sibling labels.

<table><thead><tr><th>Category Class</th><th>Crystal and Point groups</th><th>Crystal Structure representation</th><th>Piezoelectric tensor</th><th>#e<sub>ij</sub></th><th>Tensor visualization<sup>a</sup></th></tr></thead><tbody><tr><td rowspan="3">Class A</td><td>PG1<br>Cubic 23<br>Cubic 43m</td><td></td><td><table><tr><td>0</td><td>0</td><td>0</td><td>e<sub>14</sub></td><td>0</td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>e<sub>14</sub></td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>e<sub>14</sub></td></tr></table></td><td>1</td><td></td></tr><tr><td>PG2<br>Tetragonal 42m</td><td></td><td><table><tr><td>0</td><td>0</td><td>0</td><td>e<sub>14</sub></td><td>0</td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>e<sub>25</sub></td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>e<sub>36</sub></td></tr></table></td><td>3</td><td></td></tr><tr><td>PG3<br>Orthorhombic 222</td><td></td><td><table><tr><td>0</td><td>0</td><td>0</td><td>e<sub>14</sub></td><td>0</td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>e<sub>14</sub></td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>e<sub>36</sub></td></tr></table></td><td>2</td><td></td></tr><tr><td rowspan="2">Class B</td><td>PG4<br>Hexagonal 6mm<br>Tetragonal 4mm</td><td></td><td><table><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>e<sub>15</sub></td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>e<sub>15</sub></td><td>0</td><td>0</td></tr><tr><td>e<sub>31</sub></td><td>e<sub>31</sub></td><td>e<sub>33</sub></td><td>0</td><td>0</td><td>0</td></tr></table></td><td>3</td><td></td></tr><tr><td>PG5<br>Orthorhombic mm2</td><td></td><td><table><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>e<sub>15</sub></td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>e<sub>24</sub></td><td>0</td><td>0</td></tr><tr><td>e<sub>31</sub></td><td>e<sub>32</sub></td><td>e<sub>33</sub></td><td>0</td><td>0</td><td>0</td></tr></table></td><td>5</td><td></td></tr></tbody></table>

piezoelectric component prediction. To address this constraint, the application of crystal symmetry-guided clustering becomes imperative. Initially, this involves excluding monoclinic and triclinic crystals because of a higher number of independent coefficients from the scope of the regression-based ML model. Subsequently, the positions and structure of non-zero elements are systematically grouped to form distinct category classes. The organization of these classes is informed by the collected piezoelectric data, resulting in two main categories. These criteria capture the piezoelectric behaviour of 7 crystal point groups.

The crystal symmetry-based grouping is included in Table 1. In the first category, four crystal point groups together form three distinct sub-categories. Conversely, the

second category comprises three point groups organized into two sub-categories. This classification of crystal points operates on multiple layers, forming sub-classes within an overarching category class. Consequently, a data-driven approach necessitates a multi-level hierarchical classification in a specific sequence that is explained in a subsequent section.

## 2.3. Prognosticating Piezoelectric Tensors: Introducing PiezoTensorNet.

The integration of neural networks and crystallography discloses substantial potential for accelerated exploration of optimum materials design [38, 39]. In this study, this interdisciplinary fusion of pattern understanding ability of

<footnote><sup>a</sup>The tensor visualization was performed using the following tensor component values: For Class A: $e_{14} = 0.5$. For Class B Hexagonal 6mm and Tetragonal 4mm: $e_{15} = 0.0743$, $e_{31} = 0.0588$, $e_{33} = -0.0983$. For Class B Orthorhombic mm2: $e_{15} = 0.015$, $e_{15} = 0.0081$, $e_{31} = 0.014$, $e_{31} = 0.011$, $e_{33} = -0.02$. The unit of these coefficients is C/m<sup>2</sup>.</footnote>
<figure><img src="image_4.png" /><figcaption>**Figure 1:** A hierarchical multi-scale ML framework for piezoelectric material design: a) The featured piezoelectric database undergoes preprocessing, which involves PCA for feature compression. This step achieves a 96% explained variance by selecting 28-32 features. b) The PiezotensorNet operates at two scales: i) The hierarchical multi-level classifier, HierCrystalNet, categorizes materials in the initial layer into pseudo-parent classifications: class A or class B. Subsequent children's classifications assign crystal point groups, resulting in PG1, PG2, or PG3 for class A, and PG4 or PG5 for class B. ii) The regression-based ModularEnsembleNet ultimately provides the independent non-zero piezoelectric coefficient ($\#e_{ij}$) which is reordered into 3 × 6 tensor. The outputs from HierCrystalNet activate distinct modules within an ensemble artificial neural network, ensuring accurate prediction of the tensor.</figcaption></figure>

ML with crystallographic insights is considered which latter will be assisted by crystal rotation and Finite Element Analysis (FEM) in the design of functional material. The schematic depiction of this process is elucidated sequentially in Figure 1. The phenomenon of crystal point group classification unfolds within a two-tiered structure, activating a module-based ensemble neural network exclusively tailored for independent piezoelectric coefficient predictions. The final output of the ML model is further restructured to form a piezoelectric tensor in the 3×6 order. A detailed explanation of the crystal point group classification and $\#e_{ij}$ prediction process is explained in a subsequent sub-section within this section.

### 2.3.1. Hierarchical Multi-Level Classification of Crystal Groups in Piezoelectric Materials

Multi-label classification poses significant benefits as it involves the assignment of parent categories to closely related labels [40]. In this investigation, a two-layer classification was employed, focusing on the consistent patterns exhibited by mathematically homogeneous piezoelectric tensors. It is noteworthy that the child group is a true representation of the actual classification label which is the primary objective of this classification model. On the other hand, the parent groups are pseudo-clustering of the closely related point groups in assisting multi-level classification. The five sibling classes of seven crystal point groups are

characterized by their proximity to one another because of the same mathematical piezoelectric tensor as represented by piezoelectric tensor in Table 1. The initial step of the process involves the preliminary classification of piezoelectric materials, assigning them to either class A or class B. This pseudo-parent assignment serves as the first layer of classification, and the model achieves an accuracy of 68%. Subsequently, the predictions from the first layer are utilized in the second layer for specific point group classification. In the second stage, when the first layer results in class A, the focus of the classification process is on PG1, PG2, and PG3. Sub-classification within class A is carried out with an accuracy of 71.43%. Conversely, if the first layer results in class B, the attention of the second stage shifts to the classification of PG4 and PG5, achieving an accuracy of 82.22%, as illustrated in Fig. 1 b(i). This entails a finer level of categorization based on the pseudo-parent classification and then sibling point group recognition. The core of this hierarchical methodology lies in the precise identification of the material's crystal point group. Upon completion of the entire classification model, the ultimate outcomes are the crystal point groups denoted as PG1, PG2, PG3, PG4, and PG5. After the identification of the material's crystal point group, it serves as a trigger for activating the prediction process of independent piezoelectric coefficients. This prediction is accomplished through a multi-target regression-based ModularEnsembleNet.
### 2.3.2. Modular Ensemble Models to Foretell the Piezoelectric Coefficients

After the two successive classifications, this work introduces a well-engineered and intricate ensemble neural network model, comprising five meticulously designed modules. Each of them comprises an ensemble of two distinct ANN architectures. The collaborative synergy between these two ANN models involves arithmetic averaging their predicted values. This ensemble approach contributes to the generalization of the multi-dimensional regression model, resulting in a compact modular block. The multi-dimensional regression-based modular ensemble ANN seamlessly integrates with hierarchical classification, significantly improving the speed and efficiency of piezoelectric performance design.

The fusion of multi-layered classification HierCrystal-Net (Figure 1 bi) and multi-dimensional regression ModularEnsembleNet (Figure 1 b ii) is achieved from the unique activation mechanism of the ModularEnsembleNet to activate the specific modules among five different module set. This activation mechanism is designed to selectively choose a specific module from a set of five different modules. The crystal point groups are identified through the HierCrystalNet framework (Figure 1 b(i)). The activation of one among the five different modules is contingent on the outcome of the classification model, determining whether the result corresponds to PG1, PG2, PG3, PG4, or PG5. Given that a piezoelectric material possesses a unique point group, only one module is activated. This unified framework plays a crucial role in the model's architecture, where hierarchical classification guides the activation of ANN modules, as illustrated by ModularEnsembleNet in Figure 1 b(ii). The root mean squared error (RMSE) of each modular ensemble network is: 0.65, 0.85, 0.53, 1.05, 1.02 for PG1, PG2, PG3, PG4, and PG5 respectively. The outcome from the ML model is the prediction of piezoelectric coefficients in a 1D array having the same order as the $\#e_{ij}$ of the classified point group. Subsequently, the predicted results are reorganized into the format of a piezoelectric tensor with a 3×6 order.

By integrating the crystallography dominant HierCrystalNet into piezo-tensor focused ModularEnsembleNet framework, a powerful and versatile tool called PiezoTensorNet [41] is forged. This integrative approach enhances the model's capacity to first classify crystal point groups and then predict the piezoelectric coefficients, enhancing the field of piezoelectric material design and facilitating advancements in diverse scientific applications. For more detailed information on the architecture of the machine learning models, please refer to Supplementary Information S4, Tables A3, A4 and A5.

### 2.4. Crystal transformation of Piezoelectric tensor

Piezoelectric tensors, as 3rd rank tensors follow specific transformation laws. The orientation of a molecule or crystallographic lattice in relation to the laboratory coordinate system is determined by three Euler angles: $\phi$, $\theta$, and $\psi$ as visualized in Fig. 2. Two angles, $\theta$ and $\phi$, specify a direction

in space, while $\psi$ fixes the orientation of the other two sets in orthogonal axes. This can be visualized as three consecutive rotations. The direction cosine matrix is commonly used to derive the orientation of a general rotation.

$$ 
\begin{aligned}
a_1 &= \begin{bmatrix} \cos \phi & \sin \phi & 0 \\ -\sin \phi & \cos \phi & 0 \\ 0 & 0 & 1 \end{bmatrix} \\ 
a_2 &= \begin{bmatrix} 1 & 0 & 0 \\ 0 & \cos \theta & \sin \theta \\ 0 & -\sin \theta & \cos \theta \end{bmatrix} \\ 
a_3 &= \begin{bmatrix} \cos \psi & \sin \psi & 0 \\ -\sin \psi & \cos \psi & 0 \\ 0 & 0 & 1 \end{bmatrix}
\end{aligned} 
\quad (1) $$

The relationship between the $e_{ijk}$ tensor in the laboratory coordinate system obtained after eulerian rotation ($\phi$, $\theta$, $\psi$) of the $e^o_{lmn}$ tensor in the crystal coordinate system is given by the transformation for 3rd rank tensor as follow: [42]

$$ e'_{ijk} = (a_1)_{il} (a_2)_{jm} (a_3)_{kn} e^o_{lmn} \quad (2) $$

where, $[i, j, k, l, m, n] = 3$ ($a_1)_{il}$ ($a_2)_{jm}$ ($a_3)_{kn}$ are the direction cosines representing both the coordinate systems. In the representation of 3rd rank 3×3×3 in 3×6 piezoelectric tensors, the new tensor after transformation is given by [43]:

$$ e'_{ij} (\phi, \theta, \psi) = A_{il} e^o_{lm} N_{mj} \quad (3) $$

$[i, l] = 3$; $[j, m] = 6$. The Cartesian tensor notation and Einstein summation convention were used in the analytical solution of tensor mathematics and performing crystal rotations. In these engineering notations, repeated indices imply summation over those indices. The convention also allows for partial differentiation with respect to the coordinate associated with an index [44, 45]. These notations provide a concise and precise way to express and manipulate the piezoelectric tensors obtained from the multi-scale machine learning models. $A_{il}$ and $N_{mj}$ are 3×3 and 6×6 rotation matrix of direction cosines. The analytical calculation of crystal rotation guided by tensor mathematics can be found in Supplementary Information S5.

### 3. Discovering materials with enhanced piezoelectric performance

This section showcases the potential of PiezoTensor-Net in the rapid fine-tuning and discovery of piezoelectric materials. Aluminum Nitride (AlN) is a commonly used material in microelectromechanical systems (MEMS) due to its advantages of high-temperature stability, better thermal conductivity, and enhanced mechanical strength. However, its low piezoelectric coefficient has prompted studies on performance enhancement through alloying with transition metals (TM) [46, 47]. Previous investigations in this area
<figure><img src="image_5.png" /><figcaption>Figure 2: Crystallographic rotations for a crystal with three sets of Euler angles are described as follows: a) Counterclockwise rotation (φ) about the Z-axis, with a range of 0° to 360°. b) Counterclockwise rotation (θ) about the X'-axis, with a range of 0° to 180°. c) Counterclockwise rotation (ψ) about the Z''-axis, with a range of 0° to 360°. The transformation of the tensor by these rotations is given by a₁, a₂, and a₃, respectively, where each of them represents a 3×3 matrix of direction cosines. It is noteworthy that some of the zero value position of the piezoelectric coefficient transforms into non-zero coefficients.</figcaption></figure>

have focused on a limited number of TM. For instance, Uehara et al. [46] and Hirata et al. [47] explored the preparation and enhancement of AlN with specific TMs of a small subset in their analysis. Zha et al. [48] expanded the scope by investigating 28 sets of TM in AlN, alloying a single TM at a time. Despite the comprehensive approach, these studies lacked a data-driven methodology, resulting in limitations such as a focus on fewer alloying metals and primarily examining single alloying scenarios.

PiezoTensorNet is utilized in this study to investigate the different scenarios starting from pure AlN which is hexagolan 6mm crystal. The doping effect of foreign elements, both single alloying and binary configurations as visualized in Figure 3 (a) is done. Singular dopants [Ta, Zr, Y, Hf, Nb, K] result in higher e<sub>mod</sub> values compared to pure AlN, with Hf and Nb showing significantly larger 3.35 and 3.41 C/m² values at a doping composition of 0.525 and 0.725 respectively as represented in Figure 3 (b). This study focuses on observing the collective impact of the longitudinal, transverse, and shear modes within the piezoelectric coefficient. A statistical approach was implemented to deduce the overall influence of these coefficients by employing a reference coefficient modulus, denoted as e<sub>mod</sub> and explained by Eq. 4. The e<sub>mod</sub> solely serve as a criterion for selecting materials with superior longitudinal, transverse, and shear modes of piezoelectric components.

$$e_{\text{mod}_i} = \sqrt{\sum_{j=1}^{6} e_{ij}^2} \qquad (4)$$
$$e_{\text{mod}} = \max(e_{\text{mod}_i}) \quad [i = 1 \text{ to } 3]$$

Binary dopants exhibit superior piezoelectric modulus values, and among the seven binary dopant pairs studied [Zr, Hg], [Sc, Y], [Ti, Hf], [Zr, Hf], [Ta, Nb], [B, Er], and [Nb, Sc] which is illustrated in Figure 3 (c). The [B, Er]

pair achieves the optimal e<sub>mod</sub> values of 4.36 C/m² (e₃₃ 4.36 C/m² and e₃₁ 0.09 C/m²) which is followed by 3.74 C/m² (e₃₃ 3.54 C/m² and e₃₁ 0.86 C/m²) for [Nb, Sc]. The ternary surface representation in Figure 3 (d) demonstrates the distribution of modulus values across different compositions of B, Er, and AlN, with a maximum value of 4.36 C/m² attained at B₀.₃Er₀.₅Al₀.₂N. The doped piezoelectric material is also analyzed for crystal transformation using Euler angular transformations (φ, θ, ψ) as discussed in Section 2.4.

The impact of angular variation is clearly encapsulated by the plots of Figure 3 e) and f). It is noteworthy that the variation of angle φ shows no variation in tensor structure whereas there was a significant variation of tensor component diversification because of θ and ψ. The angular variations of ψ and θ depict that distribution of maximum values of all 18 components of the piezoelectric tensor falls in 91.52° and 268.47° value of ψ in the Figure 3 e) which is again validated from e*<sub>mod</sub> value lying for same values of ψ of Figure 3 f). Also, the maximum modulus is achieved for θ values of 180° and 0°. This analysis yields valuable insights into the comprehensive piezoelectric enhancement of amalgamated materials under rotational transformations.

The piezoelectric tensors obtained from the sequence of the above analysis will be used for the electro-mechanical characterization in the microelectromechanical self-energy sufficient system design.

$$e_o = \begin{bmatrix} 0 & 0 & 0 & 0 & -0.2893 & 0 \\ 0 & 0 & 0 & -0.2893 & 0 & 0 \\ -0.5801 & -0.5801 & 1.461 & 0 & 0 & 0 \end{bmatrix}$$

$$e^* = \begin{bmatrix} 0 & 0 & 0 & 0 & 0.2283 & 0 \\ 0 & 0 & 0 & 0.2283 & 0 & 0 \\ 0.0943 & 0.0943 & 4.363 & 0 & 0 & 0 \end{bmatrix}$$
<figure><img src="image_6.png" /><figcaption>Figure 3: Rapid materials design and discovery for achieving optimum piezoelectric performance by strategic incorporation of alloy doping and crystal rotation. a) The introduction of foreign elements, D1, D2, or a combination thereof, into the AlN base material initiates the alloying process. b) The influence of single dopants (Ta, Zr, Y, Hf, Nb, K) on the reference modulus of elasticity, e<sub>mod</sub>, is visually elucidated. Remarkably, all dopants yield higher e<sub>mod</sub> values compared to the base AlN, which itself exhibits a value of 1.67 C/m<sup>2</sup>. c) Comparative plots of e<sub>mod</sub>, e<sub>33</sub>, and e<sub>31</sub> for various binary dopant combinations unveil the superiority of the [B, Er] dopant pair, followed by [Nb, Sc], in terms of achieving elevated e<sub>mod</sub> values. d) The distribution of e<sub>mod</sub> values across diverse compositions of B, Er, and the AlN ceramic material is visualized through a ternary surface heatmap. Notably, a maximum value of 4.37 C/m<sup>2</sup> is attained at a composition of 0.3 B and 0.5 Er. Crystal rotation analysis. e) is conducted to investigate the impact of angular variations, specifically the ψ and θ angles, on the 18 components of the rotated piezoelectric coefficient tensor. The individual linear heatmap depicts the density distribution of each component as a function of varying ψ angles, with dark curves representing the density distribution of maximum values observed. Intriguingly, the density peaks at angles of 91.52° and 268.47°. f) Moreover, the interplay between e'<sub>mod</sub> and the angular variations of ψ and θ is explored, yielding valuable insights into the piezoelectric properties under rotational transformations.</figcaption></figure>

$$e_2^* = \text{Rotated piezoelectric tensor for alloyed AlN}, \\ \theta = 180^\circ, \psi = 268.47^\circ, 91.52^\circ; 0^\circ \le \phi \le 360^\circ$$

$$e_1^* = \begin{bmatrix} 0 & 0.08 & 0 & -0.0808 & 0.459 & -0.4554 \\ -0.1882 & -4.563 & -4.3232 & 8.7947 & -0.0807 & 0.146 \\ 0.1899 & 4.3775 & 4.3636 & -8.6464 & 0 & -0.0667 \end{bmatrix}$$

$$e_2^* = \begin{bmatrix} 1.1354 & 0.947 & 0.9856 & -1.825 & -2.005 & 1.9339 \\ 0.9564 & 1.098 & 0.9799 & -1.9783 & -1.8336 & 1.93 \\ -0.9946 & -0.9778 & -1.0156 & 1.9799 & 2.0 & -1.838 \end{bmatrix}$$

Where,

$e_o$ = Piezoelectric tensor for pure AlN

$e^*$ = Piezoelectric tensor for alloyed AlN,
$B_{0.3}Er_{0.5}Al_{0.2}N$

$e_1^* = \text{Rotated piezoelectric tensor for alloyed AlN}, \\ \theta = 0^\circ, \psi = 268.47^\circ, 91.52^\circ; 0^\circ \le \phi \le 360^\circ$

The tensor structure of pure AlN and doped material has a similar structure with five non-zero elements and has elevated piezoelectric effects because of a higher coefficient of 4.363 C/m<sup>2</sup> in $e^*$ in comparison to 1.461 C/m<sup>2</sup> of $e_o$. However, there is significant variation in tensor structure because of crystal rotation. The surface plot in Figure 4b) illustrates the tensor components $e_o$, $e^*$, $e_1^*$, and $e_2^*$ that are utilized in the FEM analysis.

### 3.1. Finite element analysis for micro-sized beams of ML-assigned piezoelectric tensors

With the in-silico toolkit to predict the piezoelectric tensors from new composition information, it can be very useful for beneficial applications if the properties can be
tested for energy generation using finite element models [49, 50]. A microelectromechanical model of a cantilever beam can be a substantially understandable and comprehensible formulation for assessing the property-performance characteristics of ML-designed novel piezoelectric materials. The governing equations for polar materials can be expressed with linear constitutive equations, which are close to those presented in Ref. [51]. For steady-state conditions the stress tensor ($T_{ij}$) of a piezoelectric beam of Fig. 4(a) undergoing elastic deformation is expressed by the equation 5,

$$-\nabla \cdot T_{ij} = \vec{f}_b \qquad (5)$$

The deformation of the piezoelectric solid material is quantified in terms of displacement vector ($\delta_i$), and it will be known that the stress tensor is a composite function of $\delta_i$. The term $\vec{f}_b$ in the right-hand side of the equation is the volumetric body force vector.

Variation of electric field potential ($\lambda$) in the micro-sized piezoelectric beam is mathematically defined by the Poisson equation 6:

$$-\nabla \cdot \varepsilon \nabla \lambda = \rho \qquad (6)$$

In Eq. 6, $\varepsilon$ is the permittivity of the piezoelectric material constituting the loaded cantilever beam, whereas $\rho$ is the volumetric charge density in the medium. Electric field intensity vector $E_k$ for the computational domain is computed as the gradient of the scalar field $\lambda$.

The electro-mechanical characterization within the simulation domain involves the interaction between mechanical stimuli (Eq. 5) and electrical response (Eq. 6). The piezoelectric tensor plays the crucial role of representing the tightly coupled physics of the electromechanical system such as the piezoelectric beams. As shown in Eq. 7, the linear piezoelectric constitutive equations are in the Stress-Electric Displacement form with strain components and electric field components as independent variables.

$$T_{ij} = c_{ijkl}S_{kl} - e_{kij}E_k \qquad (7)$$

$$D_i = e_{ijk}S_{kl} + \varepsilon_{ik}E_k \qquad (8)$$

In the equation, the strain tensor annotated as $S_{kl}$, is a function of displacement variables. $c_{ijkl}$ is of fourth-order stiffness tensor, $e_{kij}$ and $e_{ijk}$ are third order piezoelectric tensors, and the dielectric tensor $\varepsilon_{ik}$ is of second order. Applying materials symmetry conditions, the piezoelectric tensors can be reduced into the form of $e_{ki}$ and $e_{ij}$, and PiezoTensorNet is currently suited to work with the tensors of this format. The 3D surface representation of $e_{ij}$ tensors for four types of considered piezoelectric materials are graphically shown in Fig. 4(b). The tensor for undoped AlN i.e. $e_o$ as visualized in b(i) has its maximum peak value of 1.461 C/m² and maximum depth value of -0.5801 C/m². The corresponding tensors for the different doped variants are presented graphically in b(ii)- (iv). Without incorporating the rotation effect, the doping of AlN with B and Er in the composition B₀.₃Er₀.₅Al₀.₂N produces the 3D profile of e* tensor of b(ii). The maximum peak (maximum positive

valued coefficient = 4.363 C/m²) for this material is much higher than that of undoped AlN. However, each of the coefficients of the e* tensor is zero or larger, and so it loses the merit associated with the depth of the surface. For doped alloy (B₀.₃Er₀.₅Al₀.₂N), and rotated with $\psi = 91.52°$ or 268.47° and $\theta = 0°$, both the peak and depths are of larger magnitudes. Referring to b)(iii), the maximum peak of $e_1^{*'}$ tensor is of 8.7947 C/m² whereas the maximum depth is of -8.6464 C/m². When the rotational angle $\psi$ is maintained at 91.52° or 268.47° and $\theta$ is changed to 180°, the profile of $e_1^{*'}$ for B₀.₃Er₀.₅Al₀.₂N as shown in b) iv) has relatively smaller peak values (largest positive value of 2.0 C/m²) and shallower depth (most negative valued coefficient -2.005 C/m²). For the AlN ceramic materials and doped-materials, the relative permittivity of 9.1 is utilized. As illustrated earlier, the piezoelectric tensor for the 3D cantilever beam of AlN ceramics or doped alloys is represented by 18 coefficients in a 3 × 6 matrix. Applying symmetry conditions, the stiffness matrix for isotropic and orthotropic 3D beam materials can be represented with a 6 × 6 matrix with a total of 36 elements.

For isotropic assumption, the following stiffness tensor is used to represent the case of AlN piezoelectric material:

$$C_{ijkl}^{AlN} = 10^{11} \times \begin{bmatrix} 3.54 & 1.12 & 1.12 & 0 & 0 & 0 \\ 1.12 & 3.54 & 1.12 & 0 & 0 & 0 \\ 1.12 & 1.12 & 3.54 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1.21 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1.21 & 0 \\ 0 & 0 & 0 & 0 & 0 & 1.21 \end{bmatrix} N/m^2$$

Composition-based feature vector (CBFV) methodology, illustrated in Poudel et al. [31] has been utilized to compute the effective stiffness tensor for doped alloy. The stiffness tensor for B₀.₃Er₀.₅Al₀.₂N alloy (annotated as d-AlN) for elastic isotropy condition is:

$$C_{ijkl}^{d-AlN} = 10^{11} \times \begin{bmatrix} 2.65 & 0.84 & 0.84 & 0 & 0 & 0 \\ 0.84 & 2.65 & 0.84 & 0 & 0 & 0 \\ 0.84 & 0.84 & 2.65 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0.91 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0.91 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0.91 \end{bmatrix} N/m^2$$

The conditions for materials isotropy are utilized in this section to highlight the sole effect of alteration in piezoelectric coefficients in doped alloys. The values of Young's modulus ($Y_E$) and Poisson's ratio $\nu$ for AlN are 300 GPa and 0.24 respectively. Whereas, for B₀.₃Er₀.₅Al₀.₂N the corresponding stiffness tensor indicates $Y_E$ and $\nu$ as 225 GPa and 0.24 respectively. The FEM simulations corresponding to orthotropic elasticity are presented in Supplementary Information S6. The cantilever beam in Fig. 4(a) has a dimension of $l = 200$ µm, $b = 50$ µm and $h = 5$ µm. At the end face of the beam (AlN ceramics or its doped versions) attached to the fixed wall, fixed boundary condition (zero displacement with $\delta_x = \delta_y = \delta_z = 0$) is applied. At the face at the opposite end of this piezoelectric beam, a surface load ($F_z$) is applied in the downward direction. This stress load is applied with constant values ranging from $1 \times 10^{-5}$ to ($F_z$)$_{max} = 5 \times 10^{-5}$ N/m². The top and bottom faces
<figure><img src="image_7.png" /><figcaption>**Figure 4:** Finite element analysis for micro-sized beams. a) The cantilever beam (200 × 50 × 5) µm is subjected to varying loading conditions ranging from 1× 10⁻⁵ N/m² to a maximum (Fz)max = 5× 10⁻⁵ N/m². The dimensions are not according to the scale. The top and bottom surfaces are electrically grounded for boundary conditions. b) Deformation, voltage generation, and power results obtained from the electro-mechanical coupled FEA for five different loading conditions of four different tensor properties. c) Results of deformation distribution along the z-axis (δz) for a load (Fz)max in the pure AlN and alloyed beam. d) Electric potential (λ) for (Fz)max showing the maximum value of 7.8 V for the corresponding maximum deformation of 0.97 µm.</figcaption></figure>

of the beam are grounded electrically. With the materials properties and boundary conditions, the set of Eqs. 5 and 6 is solved using finite element method in Elmer Multiphysics software [52, 53].

As illustrated in Fig. 4 b), the potential distribution and deformation variation for B₀.₃Er₀.₅Al₀.₂N alloy (with the orientation described as ψ = 91.52° or 268.47° and θ = 0°) are presented. With the matrix for e₁*′ tensor being asymmetric, the voltage distribution pattern is characterized by spatial anisotropy. For a stress load of (Fz)max, maximum displacement with a magnitude of 0.97 µm is attained and potential field with a maximum value of 7.8 V is generated in this material type. It is further revealed from the curves of Fig. 4(d) that the value of maximum λ for this piezoelectric material (represented by the curve for Doped and rotated, e₁*′) is around 3.39 times that of the maximum value obtained in the AlN beam (material denoted by the curve named Pure AlN, e₀) at the stress loading of (Fz)max. The maximum value of potential in pure AlN for the above loading is 2.3 V. The maximum voltage for the other two variants of doped alloys with piezoelectric tensors of e* and e₂*′ are 3.9 and 4.7 V respectively for the same loading. Thus, the doped material in general offers better electromechanical merits

compared to the AlN beam. Moreover, the doped candidate with rotation orientation resulting in e₁*′ is far superior to any of the other three candidates. The results of voltage and deformation distribution for a range of constant loads, varying from 1×10⁻⁵ N/m² to 5×10⁻⁵ N/m², as presented in Figure 4 (b) offer more insights than just the finetuning of electric potential for a given mechanical loading. The deformation curve for pure AlN material is of visibly shorter length than that of the other three variants all composed of the same doped material (B₀.₃Er₀.₅Al₀.₂N). Due to the compositional differences between AlN and B₀.₃Er₀.₅Al₀.₂N, their stiffness tensors are distinct. This variation in elastic response between the two materials leads to different deformations in pure AlN and its doped variants. For the vibrating beam of mass density (ρ) of 7500 Kg/m³ and moment of inertia (I = bh³/12), the natural frequency at first mode is given as $f = \frac{1.875^2}{2\pi} \sqrt{\frac{Y_E I}{pbhl^4}}$ [54]. From Fig. 4 c), Under (Fz)max, the curves demonstrate that a cantilever beam composed of AlN ceramic exhibits a maximum vertical displacement of only 0.73 µm at a first mode frequency of 629.11 kHz. In contrast, the B₀.₃Er₀.₅Al₀.₂N beams already deform by 0.97 µm at 544.82 kHz.
The average power (P) for the MEMS computed in ac-
cordance to the relation P = (2πf)CV² [55]. Where, f, C, and
V are the vibrational frequency of the beam, capacitance of
the connected capacitor, and generated voltage respectively.
For this system, considering the capacitance of 0.0492 pF
and for (F₂)max, the power is obtained as 2082.71 μW for
the rotated B₀.₃Er₀.₅Al₀.₂N (e₁') beam, whereas it is 209.1
μW for pure AlN (e₀). This represents a 9.96-fold increase
compared to pure AlN. In the case of the other system, the
maximum power is 520.82 μW for e* and 756.19 μW for e₂'.
In energy harvesting applications, the devices are designed
based on loading conditions so the inclusion of finite element
analysis will add to the flexibility of choosing the materials
that provide the optimum range of deflections for a given
applied force.

4. Conclusion

The integration of computational engineering, crystal-
lography, and electro-mechanical characterization in the de-
sign of piezoelectric materials is the focus of this study.
The introduced framework, PiezoTensorNet, is scientifically
robust and combines the interdisciplinary science of crystal
structure, tensor mathematics, and machine learning tech-
niques. It combines a hierarchical neural network-based
crystal point group classification with a modular ensemble
neural network for piezoelectric coefficient prediction. The
fusion of these interdisciplinary approaches results in an
efficient model that combines a multi-layer classification
model and a multi-dimensional regression-based modular
network. The crystal point group classification influences
the activation of specific blocks within the modular network.
This framework also incorporates the study of piezoelectric
behavior variations in the presence of single or binary for-
eign elements as dopants. Furthermore, it includes crystal
transformation and the enhancement of the piezoelectric
coefficient using Euler angular rotation by the combinations
of (φ, θ, ψ). The application of tensor mathematics supports
the physics of crystal rotation. This framework effectively
utilizes and enhances the piezoelectric response of the ma-
terial through electro-mechanical coupling. The integration
of computational engineering, crystallography, and electro-
mechanical characterization enables rapid design improve-
ments in piezoelectric materials.

In this research, the multidisciplinary integrated PiezoTen-
sorNet was applied to a MEMS based on Aluminum Nitride
(AlN). The study focused on the impact of singular and
binary foreign element dopants on the piezoelectric perfor-
mance of the system. It was determined that the composition
[B, Er] with a ratio of B₀.₃Er₀.₅Al₀.₂N exhibited optimal
piezoelectric properties. Using crystallography-guided ML
and FEM, the performance of a cantilever beam was an-
alyzed. The five different loading conditions which cause
varying deformation were studied. The results showed that
the power increased from 209.1 µW 520.67 µW in the doped
system when subjected to a case of maximum deformation

of 0.97 µm. Furthermore, crystallography rotation was per-
formed on the doped alloy, and it was found that rotating the
crystal through 91.52º along the vertical axis (Z-axis) from
either direction resulted in the highest power generation,
reaching a value of 2082.71 µW. This study demonstrates
the scientific application of PiezoTensorNet in the field of
piezoelectric materials, specifically in the context of MEMS
based on AlN. By exploring the effects of dopants and
crystallography rotation, significant improvements in the
piezoelectric performance of the system were achieved.

Declaration of Competing Interests

The authors declare that they have no known competing
financial interests or personal relationships that could have
appeared to influence the work reported in this paper.

CRediT author statement

**Sachin Poudel:** Writing–original draft, Visualization, Validation, Software, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. **Rubi Thapa:** Writing–review & editing, Visualization, Validation, Investigation, Formal analysis. **Rabin Basnet:** Writing–review & editing, Validation, Methodology. **Anna Timofiejczuk:** Writing–review & editing, Validation, Supervision, Methodology. **Anil Kunwar:** Writing–review & editing, Validation, Supervision, Software, Methodology, Investigation, Funding acquisition, Data curation, Conceptualization.

Acknowledgments

This work was supported by the National Science Centre,
Poland (UMO-2021/42/E/ST5/00339)

**Appendix A. Supplementary Information**

The Supplementary Material to this article is a document
consisting of 6 different sections. The page numbers within
the supplementary information document are indicated us-
ing Roman numerals.

The following are the six sections of the supplementary
material, each of which will be described in detail later in
the document.

S1. *Piezoelectric Coefficients: Mathematical Expressions:* contains **Table A1**, expressing piezoelectric coefficients in terms of electric field, electric displacement, mechanical stress, and strain.

S2. *Crystal Point Group and Piezoelectricity:* includes **Table A2**, with noncentrosymmetric point groups and piezoelectric tensor in 20 crystal classes.

S3. *Data Visualization and Feature Engineering:* comprises **Table A3**, Magpie Descriptors for Materials Analysis with added complementing features, alongside **Fig. A1**, Visualization and feature engineering of collected Piezoelectric datasets.
**S4. Machine Learning Descriptions:** consists of **Tables A4, A5, and A6** detailing hyperparameter selection and machine learning descriptions, along with **Fig. A2** presenting learning curves for training and validation accuracy.

**S5. Crystal Rotation Supported by Tensor Mathematics:**

**S6. FEM for Orthotropic Stiffness:** features **Fig. A3**, illustrating FEM results of generated voltage and corresponding deformation for the orthotropic stiffness matrix.

## Data Availability

The data used in this work is available in https://github.com/Sachinscnpdl/PiezoTensorNet-Data. The graphical user interface of PiezoTensorNet can be accessed via: https://piezotensornet.streamlit.app/

## References

[1] Y. Liu, J. Yang, S. Deng, Y. Zhang, Y. Zhang, S. Sun, L. Wang, X. Jiang, C. Huo, H. Liu, et al., Flexible polarization configuration in high-entropy piezoelectrics with high performance, Acta Materialia 236 (2022) 118115.

[2] G. Sui, X. Shan, C. Hou, H. Tian, J. Hu, T. Xie, An underwater piezoelectric energy harvester based on magnetic coupling adaptable to low-speed water flow, Mechanical Systems and Signal Processing 184 (2023) 109729.

[3] S. Chen, P. Zhu, L. Mao, W. Wu, H. Lin, D. Xu, X. Lu, J. Shi, Piezocatalytic medicine: An emerging frontier using piezoelectric materials for biomedical applications, Advanced Materials (2023) 2208256.

[4] C. Wang, Y. Hu, Y. Liu, Y. Shan, X. Qu, J. Xue, T. He, S. Cheng, H. Zhou, W. Liu, et al., Tissue-adhesive piezoelectric soft sensor for in vivo blood pressure monitoring during surgical operation, Advanced Functional Materials (2023) 2303696.

[5] H. Tian, X. Shan, X. Li, J. Wang, Enhanced airfoil-based flutter piezoelectric energy harvester via coupling magnetic force, Applied Energy 340 (2023) 120979.

[6] W. Sun, G. Ji, J. Chen, D. Sui, J. Zhou, J. Huber, Enhancing the acoustic-to-electrical conversion efficiency of nanofibrous membrane-based triboelectric nanogenerators by nanocomposite composition, Nano Energy 108 (2023) 108248.

[7] J. Deng, W. Chen, Y. Wang, S. Zhang, Y. Liu, Modeling and experimental evaluations of a four-legged stepper rotary precision piezoelectric stage, Mechanical Systems and Signal Processing 132 (2019) 153–167.

[8] J. M. Muganda, B. Jansen, E. Homburg, Y. van de Burgt, J. den Toonder, Influence function measurement technique using the direct and indirect piezoelectric effect for surface shape control in adaptive systems, IEEE Transactions on Automation Science and Engineering 19 (2021) 994–1002.

[9] J. Erhart, P. Púlpán, M. Pustka, Piezoelectric ceramic resonators, Springer, 2017.

[10] Y. Benveniste, G. Dvorak, Uniform fields and universal relations in piezoelectric composites, Journal of the Mechanics and Physics of Solids 40 (1992) 1295–1312.

[11] H. Khan, N. Mahmood, A. Zavabeti, A. Elbourne, M. A. Rahman, B. Y. Zhang, V. Krishnamurthi, P. Atkin, M. B. Ghasemian, J. Yang, et al., Liquid metal-based synthesis of high performance monolayer sns piezoelectric nanogenerators, Nature communications 11 (2020) 3449.

[12] M. C. Sekhar, E. Veena, N. S. Kumar, K. C. B. Naidu, A. Mallikarjuna, D. B. Basha, A review on piezoelectric materials and their applications, Crystal Research and Technology 58 (2023) 2200130.

[13] B. G. Baraskar, Y. D. Kolekar, B. R. Thombare, A. R. James, R. C. Kambale, C. Ramana, Enhanced piezoelectric, ferroelectric, and electrostrictive properties of lead-free (1-x) bczt-(x) bct electroceramics with energy harvesting capability, Small (2023) 2300549.

[14] J. Li, Z. Shen, X. Chen, S. Yang, W. Zhou, M. Wang, L. Wang, Q. Kou, Y. Liu, Q. Li, et al., Grain-orientation-engineered multilayer ceramic capacitors for energy storage applications, Nature materials 19 (2020) 999–1005.

[15] M. Pomeroy, Encyclopedia of materials: technical ceramics and glasses, Elsevier, 2021.

[16] J. He, C. Yu, Y. Hou, X. Su, J. Li, C. Liu, D. Xue, J. Cao, Y. Su, L. Qiao, et al., Accelerated discovery of high-performance piezocatalyst in batio3-based ceramics via machine learning, Nano Energy 97 (2022) 107218.

[17] R. Yuan, D. Xue, Y. Xu, D. Xue, J. Li, Machine learning combined with feature engineering to search for batio3 based ceramics with large piezoelectric constant, Journal of Alloys and Compounds 908 (2022) 164468.

[18] S.-H. V. Oh, W. Hwang, K. Kim, J.-H. Lee, A. Soon, Using feature-assisted machine learning algorithms to boost polarity in lead-free multicomponent niobate alloys for high-performance ferroelectrics, Advanced Science 9 (2022) 2104569.

[19] A. L. Kholkin, N. A. Pertsev, A. V. Goltsev, Piezoelectricity and crystal symmetry, Piezoelectric and acoustic materials for transducer applications (2008) 17–38.

[20] M. Abbasipour, R. Khajavi, A. H. Akbarzadeh, A comprehensive review on piezoelectric polymeric and ceramic nanogenerators, Advanced Engineering Materials 24 (2022) 2101312.

[21] K. Choudhary, K. F. Garrity, V. Sharma, A. J. Biacchi, A. R. Hight Walker, F. Tavazza, High-throughput density functional perturbation theory and machine learning predictions of infrared, piezoelectric, and dielectric responses, NPJ Computational Materials 6 (2020) 64.

[22] G. Rao, H. Fang, T. Zhou, C. Zhao, N. Shang, J. Huang, Y. Liu, X. Du, P. Li, X. Jian, et al., Robust piezoelectricity with spontaneous polarization in monolayer tellurene and multilayer tellurium film at room temperature for reliable memory, Advanced Materials 34 (2022) 2204697.

[23] A. Bartasyte, G. Clementi, Q. Micard, I. Labbaveettil, A. S. L. Moreira, S. Boujnah, M. Ouhabaz, A. Verma, A. Ichangi, G. Mallandrino, et al., Material strategies to enhance the performance of piezoelectric energy harvesters based on lead-free materials, Journal of Micromechanics and Microengineering 33 (2023) 053001.

[24] Z. J. Wegert, A. P. Roberts, V. J. Challis, Novel properties of multi-poled piezoelectric network structures, Smart Materials and Structures 30 (2021) 105016.

[25] J. Yang, Z. Li, X. Xin, X. Gao, X. Yuan, Z. Wang, Z. Yu, X. Wang, J. Zhou, S. Dong, Designing electromechanical metamaterial with full nonzero piezoelectric coefficients, Science Advances 5 (2019) eaax1782.

[26] J. Wu, Z. Hu, X. Gao, M. Cheng, X. Zhao, W. Su, L. Wang, M. Guan, Y. Du, R. Mao, et al., Unconventional piezoelectric coefficients in perovskite piezoelectric ceramics, Journal of Materiomics 7 (2021) 254–263.

[27] A. Jain, S. P. Ong, G. Hautier, W. Chen, W. D. Richards, S. Dacek, S. Cholia, D. Gunter, D. Skinner, G. Ceder, K. A. Persson, Commentary: The Materials Project: A materials genome approach to accelerating materials innovation, APL Materials 1 (2013). 011002.

[28] M. De Jong, W. Chen, H. Geerlings, M. Asta, K. A. Persson, A database to enable discovery and design of piezoelectric materials, Scientific data 2 (2015) 1–13.

[29] L. Ward, A. Agrawal, A. Choudhary, C. Wolverton, A general-purpose machine learning framework for predicting properties of inorganic materials, npj Computational Materials 2 (2016) 1–7.

[30] G. Deffrennes, K. Terayama, T. Abe, R.amura, A machine learning-based classification approach for phase diagram prediction, Materials & Design 215 (2022) 110497.
[31] S. Poudel, U. Subedi, M. O. Hamid, K. Gyanwali, N. Moelans, A. Timofiejczuk, A. Kunwar, Alloymanufacturingnet for discovery and design of hardness-elongationsynergy in multi-principal element alloys, Engineering Applications of Artificial Intelligence (2024) 107902.

[32] Y.-R. Chang, R. Nanae, S. Kitamura, T. Nishimura, H. Wang, Y. Xi-ang, K. Shinokita, K. Matsuda, T. Taniguchi, K. Watanabe, et al., Shift current photovoltaics based on a noncentrosymmetric phase in in-plane ferroelectric ns, Advanced Materials (2023) 2301172.

[33] R. E. Newnham, Properties of materials: anisotropy, symmetry, structure, Oxford University Press on Demand, 2005.

[34] A. Kelly, K. M. Knowles, Crystallography and crystal defects, John Wiley & Sons, 2020.

[35] M. Smith, S. Kar-Narayan, Piezoelectric polymers: Theory, challenges and opportunities, International Materials Reviews 67 (2022) 65-88.

[36] M.-M. Yang, Z.-D. Luo, Z. Mi, J. Zhao, S. P. E. M. Alexe, Piezoelectric and pyroelectric effects induced by interface polar symmetry, Nature 584 (2020) 377-381.

[37] S. Trolier-McKinstry, Crystal chemistry of piezoelectric materials, in: Piezoelectric and Acoustic Materials for Transducer Applications, Springer, 2008, pp. 39-56.

[38] K. Veasna, Z. Feng, Q. Zhang, M. Knezevic, Machine learning-based multi-objective optimization for efficient identification of crystal plasticity model parameters, Computer Methods in Applied Mechanics and Engineering 403 (2023) 115740.

[39] G. Trezza, L. Bergamasco, M. Fasano, E. Chiavazzo, Minimal crystallographic descriptors of sorption properties in hypothetical mofs and role in sequential learning optimization, npj Computational Materials 8 (2022) 123.

[40] H. Huang, H. Liu, Feature selection for hierarchical classification via joint semantic and structural information of labels, Knowledge-Based Systems 195 (2020) 105655.

[41] S. Poudel, R. Thapa, A. Kunwar, PiezoTensorNet: A machine learning based tooklit for exploratory design and discovery of piezoelectric materials, https://piezotensornet.streamlit.app/, 2023. Accessed: 03.07.2023.

[42] S. V. Kalinin, B. J. Rodriguez, S. Jesse, J. Shin, A. P. Baddorf, P. Gupta, H. Jain, D. B. Williams, A. Gruverman, Vector piezoresponse force microscopy, Microscopy and Microanalysis 12 (2006) 206-220.

[43] Y. Jin, X. Lu, J. Zhang, Y. Kan, H. Bo, F. Huang, T. Xu, Y. Du, S. Xiao, J. Zhu, Studying the polarization switching in polycrystalline bifeo3 films by 2d piezoresponse force microscopy, Scientific reports 5 (2015) 12237.

[44] J. Coope, R. Snider, F. McCourt, Irreducible cartesian tensors, The Journal of Chemical Physics 43 (1965) 2269-2275.

[45] A. P. Harrison, D. Joseph, Numeric tensor framework: Exploiting and extending einstein notation, Journal of computational science 16 (2016) 128-139.

[46] M. Uehara, Y. Amano, S. A. Anggraini, K. Hirata, H. Yamada, M. Akiyama, Preparation of ybaln piezoelectric thin film by sputtering and influence of yb concentration on properties and crystal structure, Ceramics International 47 (2021) 16029-16036.

[47] K. Hirata, H. Yamada, M. Uehara, S. A. Anggraini, M. Akiyama, Enhancement of piezoelectric property in mgtmal (tm= cr, mo, w): First-principles study, Journal of Physics and Chemistry of Solids 152 (2021) 109913.

[48] X.-H. Zha, X. Ma, J.-T. Luo, C. Fu, Enhanced piezoelectric response of aln via alloying of transitional metals, and influence of type and distribution of transition metals, Nano Energy 111 (2023) 108390.

[49] E. Herth, E. Algre, J.-Y. Rauch, J.-C. Gerbedoen, N. Defrance, P. Delobelle, Modeling and characterization of piezoelectric beams based on an aluminum nitride thin-film layer, physica status solidi (a) 213 (2016) 114-121.

[50] M. Hara, L. Minh, H. Oguchi, H. Kuwano, Power estimation for piezoelectric energy harvesters in flexure mode with large displacement amplitude, in: Journal of Physics: Conference Series, volume 476, IOP Publishing, p. 012024.

[51] C. Zhang, X. Wang, W. Chen, J. Yang, An analysis of the extension of a zno piezoelectric semiconductor nanofiber under an axial force, Smart Materials and Structures 26 (2017) 025030.

[52] T. Sikanen, T. Zwinger, S. Tuomikoski, S. Franssila, R. Lehtiniemi, C.-M. Fager, T. Kotiaho, A. Pursula, Temperature modeling and measurement of an electrokinetic separation chip, Microfluidics and Nanofluidics 5 (2008) 479-491.

[53] T. Gustafsson, P. Råback, J. Videman, Mortaring for linear elasticity using mixed and stabilized finite elements, Computer Methods in Applied Mechanics and Engineering 404 (2023) 115796.

[54] L. Barboni, G. Gillich, C. Chioncel, C. Hamat, I. Mituletu, A method to precise determine the young's modulus from dynamic measurements, in: IOP Conference Series: Materials Science and Engineering, volume 416, IOP Publishing, p. 012063.

[55] G.-W. Kim, J. Kim, J.-H. Kim, Flexible piezoelectric vibration energy harvester using a trunk-shaped beam structure inspired by an electric fish fin, International Journal of Precision Engineering and Manufacturing 15 (2014) 1967-1971.
# Appendix A. Supplementary Information:

## PiezoTensorNet: Crystallography informed multi-scale hierarchical machine learning model for rapid piezoelectric performance finetuning

### S1. Piezoelectric Coefficients: Mathematical Expressions

The piezoelectric response in a material is directly influenced by the piezoelectric tensors. In total, there are four piezoelectric coefficients, $d_{ij}$, $e_{ij}$, $g_{ij}$, and $h_{ij}$, which is represented in Table A1 along with their corresponding representations. The electric field strength, $E$, relates to the electric flux density, $D$, while the linearized strain, $S$, is connected to the stress, $T$. In this study, the focus is on analyzing the piezoelectric stress tensor, $e_{ij}$, to understand and study the piezoelectric response of the material. By examining this tensor, insights can be gained into the material's behavior and its ability to convert mechanical stress into electrical signals or vice versa.

<table><caption>**Table A1:** This table presents the four piezoelectric coefficients in terms of the electric field, electric displacement, mechanical stress, and strain. The expressions for direct and indirect piezoelectric effects are included, along with their respective units. The symbols have meaning throughout this article as electric field strength $E$, electric flux density $D$, the linearized strain $S$, stress $T$</caption><thead><tr><th>Symbol</th><th>Expression</th><th>Units</th></tr></thead><tbody><tr><td rowspan="2">d<sub>ij</sub></td><td>(∂D<sub>i</sub>/∂T<sub>j</sub>) E</td><td>[C/N]</td></tr><tr><td>(∂S<sub>j</sub>/∂E<sub>i</sub>) T</td><td>[mV/N]</td></tr><tr><td rowspan="2">e<sub>ij</sub></td><td>(∂D<sub>i</sub>/∂S<sub>j</sub>) E</td><td>[C/m<sup>2</sup>]</td></tr><tr><td>-(∂T<sub>j</sub>/∂E<sub>i</sub>) S</td><td>[N/V·m]</td></tr><tr><td rowspan="2">g<sub>ij</sub></td><td>-(∂E<sub>i</sub>/∂T<sub>j</sub>) D</td><td>[V/m·N]</td></tr><tr><td>(∂S<sub>j</sub>/∂D<sub>i</sub>) T</td><td>[m<sup>2</sup>/C]</td></tr><tr><td rowspan="2">h<sub>ij</sub></td><td>-(∂E<sub>i</sub>/∂S<sub>j</sub>) D</td><td>[V/m]</td></tr><tr><td>(∂T<sub>j</sub>/∂D<sub>i</sub>) S</td><td>[N/C]</td></tr></tbody></table>

### S2. Crystal point group and piezoelectricity

The relationship outlined by Neumann's law connects a material's physical characteristics to its crystallographic point group symmetry. Among the 32 groups, only 21 lack a central inversion, facilitating piezoelectric behavior. However, centrosymmetric structures, exemplified by point group 432, hinder odd-rank tensor properties, thereby eliminating piezoelectricity (Trolier-McKinstry (2008)). Conversely, the 20 noncentrosymmetric groups hold promise for demonstrating piezoelectric properties. These groups encompass crystal systems with Laue class, noncentrosymmetric point groups housing independent piezoelectric tensor elements across 20 crystal classes endowed with piezoelectricity. From a mathematical perspective, the crystal class is classified based on its correlation with the minimum number of independent components required to fully characterize the material's piezoelectric response. The numeric structure of the piezoelectric tensor relies on the symmetry displayed by the crystal structure and its corresponding point group symmetry (Prince (2004)). The crystal system with laue class and point group and the number of independent piezoelectric coefficients is as tabulated in Table A2. It is to be noted that $\#e_{ij}$ represents the number of independent coefficients from the tensor, for example in cubic crystal there are three non-zero coefficients and all are equal to $e_{14}$, so, cubic system will have 1 independent coefficient instead of 3.
**Table A2**

Crystal systems with Laue class, noncentrosymmetric point group with independent piezoelectric tensor elements in 20 crystal classes with piezoelectric property. #e<sub>ij</sub> indicates the number of independent coefficients in the fourth-order piezoelectric tensor.

<table><thead><tr><th>Crystal Structure</th><th>Laue class</th><th>Noncentrosymmetric Point group</th><th># e<sub>ij</sub></th></tr></thead><tbody><tr><td rowspan="2">Cubic</td><td>2/m<sup>3</sup></td><td>23</td><td rowspan="2">1</td></tr><tr><td>4/m<sup>3</sup>2/m</td><td><sup>4</sup>3m</td></tr><tr><td rowspan="4">Hexagonal</td><td rowspan="2">6/m</td><td>6</td><td>4</td></tr><tr><td><sup>6</sup>6</td><td>2</td></tr><tr><td rowspan="2">6/m2/m2/m.</td><td>622</td><td>1</td></tr><tr><td>6mm</td><td>3</td></tr><tr><td rowspan="2">Trigonal</td><td><sup>3</sup>3</td><td><sup>6</sup>2m</td><td>1</td></tr><tr><td><sup>3</sup>2/m</td><td>3</td><td>6</td></tr><tr><td rowspan="3">Tetragonal</td><td rowspan="2">4/m</td><td>32</td><td>2</td></tr><tr><td>3m</td><td>4</td></tr><tr><td>4/m</td><td>4</td><td>4</td></tr><tr><td rowspan="4">Orthorhombic</td><td rowspan="2">4/m2/m2/m</td><td>422</td><td>1</td></tr><tr><td>4mm</td><td>3</td></tr><tr><td rowspan="2">2/m2/m2/m</td><td><sup>4</sup>2m</td><td>2</td></tr><tr><td>222</td><td>3</td></tr><tr><td>Monoclinic</td><td>2/m</td><td>mm2</td><td>5</td></tr><tr><td rowspan="2">Triclinic</td><td rowspan="2">1</td><td>2</td><td>8</td></tr><tr><td>m</td><td>10</td></tr><tr><td></td><td></td><td>1</td><td>18</td></tr></tbody></table>

**S3. Data visualization and Feature Engineering**

The graph in Figure A1a provides a visual comparison of the occurrence frequencies for 58 unique elements in the collected piezoelectric datasets. Each element is represented by a bar, and the length of each bar corresponds to the number of times that element appears in the datasets. From the graph, it can be observed that 21 elements occur less than 25 times, 20 occur between 25 to 75 times, and 17 elements occur more than 75 times. Notably, oxygen exhibits the highest occurrence frequency among all the elements, indicating its significant presence in the collected piezoelectric materials. Similarly, Figure A1b showcases the constituent piezoelectric alloys found in the datasets, with each alloy represented by different numbers of elements. The data collected encompasses a wide range of elemental compositions, varying from alloys consisting of a single element to those composed of up to six different elements. Specifically, there are 3 alloys consisting of a single element, 190 alloys with 2 elements, 749 alloys with 3 elements, 339 alloys with 4 elements, 72 alloys with 5 elements, and 3 alloys with 6 elements. This representation highlights the diversity in the compositions of piezoelectric materials and underscores the potential complexity involved in their design and synthesis.

Magpie descriptors are a collection of composition-based materials attributes that encompass stoichiometric attributes, elemental properties, electronic structure attributes, and ionic compound attributes. These descriptors have found extensive usage in constructing machine learning models for predicting properties of composition-based materials (Ward et al. (2016)). In our study, we have utilized the statistics of elemental properties as the foundational dataset for predicting the space group and crystal system. The Magpie element property statistics consist of 22 features, including Atomic Number, Mendeleev Number, Atomic Weight, Melting Temperature, Periodic Table Row and Column, Covalent Radius, Electronegativity, the number of Valence electrons in each orbital (s, p, d, f, total), the
<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /><figcaption>**Figure A1:** Visualization of Constituent Elements in the Collected Piezoelectric Datasets a) Graphical comparison depicting the number of occurrences of 58 unique elements in the collected piezoelectric datasets. Oxygen exhibits the highest frequency among the elements. b) Representation of the constituent piezoelectric alloys with varying numbers of elements. The collected data encompasses elemental compositions ranging from one to six elements. c) The feature compression with PCA and corresponding principal component for each model training.</figcaption></figure>

number of unfilled electrons in each orbital (s, p, d, f, total), Ground State Volume, Ground State Band Gap Energy, Ground State Magnetic Moment, and the Space Group Number of elements. The primary features of the Magpie feature set are obtained by calculating the mean, average deviation, range, mode, minimum, and maximum of the aforementioned elemental properties, with weights determined by the fraction of each element in the composition (Li et al. (2021)). This transformation of raw materials data enables compatibility with machine learning algorithms. Additionally, a new set of descriptors that contribute to the improvement of the piezoelectric system were introduced to complement the magpie features.
<caption>**Table A3**</caption>

Elemental Properties and Predictors in Magpie Descriptors for Materials Analysis and added feature to complement the magpie features

<table><thead><tr><td>Element Property statistics of Magpie</td><td>Added Predictors in this work</td></tr></thead><tbody><tr><td>Atomic Number</td><td>Average Atom Number</td></tr><tr><td>Mendeleev Number</td><td>Average Thermal Conductivity</td></tr><tr><td>Atomic Weight</td><td>Average Bulk Modulus</td></tr><tr><td>Melting Temperature</td><td>Average Bulk Modulus Asymmetry</td></tr><tr><td>Periodic Table Row and Column</td><td>Average Shear Modulus</td></tr><tr><td>Covalent Radius</td><td>Average Shear Modulus Asymmetry</td></tr><tr><td>Electronegativity</td><td>Average Atomic Radius Difference</td></tr><tr><td>Number of Valence e in each Orbital (s, p, d, f, total)</td><td>Average Electronegativity</td></tr><tr><td>Number of unfilled e in each orbital (s, p, d, f, total)</td><td>Average Melting Temperature</td></tr><tr><td>Ground State Volume</td><td>Average Melting Temperature Asymmetry</td></tr><tr><td>Ground State Band Gap Energy</td><td>Valence Electron Concentration</td></tr><tr><td>Ground State Magnetic Moment</td><td>Mixing of Entropy</td></tr><tr><td>Space Group Number</td><td>Lambda Parameter in Entropy</td></tr></tbody></table>

Principal Component Analysis (PCA) was employed for feature selection in this study. PCA revealed the explained variance and the number of features required to capture 96% of the explained variance. In the classification process, 1354 datasets were used for crystal point groups. Hierarchical classification determined that 29 principal components captured 96% of the data variance. For modular tensor prediction, different crystal structures were considered, with varying numbers of principal components selected. For the cubic structure (42m), 33 principal components accounted for 96% of the variance. For the tetragonal structure (42m), 24 principal components were identified. The orthorhombic structure (222) utilized 28 principal components, while the combined hexagonal (6mm) and tetragonal (4mm) structures employed 30 principal components. The orthorhombic structure (mm2) also utilized 30 principal components. Figure A1c illustrates the results of the feature selection and the corresponding variances.

## S4. Machine Learning Descriptions

In our study, seven point groups are initially considered, and they are segmented into five distinct labels. Recognizing the challenge posed by having multiple models and the associated complexity, we implemented a clustering mechanism to group these labels into two category classes, namely Class A and B. This clustering helps create a more manageable structure, improving the overall effectiveness of the fault classification. The classification process begins by determining whether the piezoelectric material falls into Class A or B. Subsequently, based on this outcome, the second layer is activated to further classify the material. The two-layer classification system focuses on consistent patterns exhibited by mathematically homogeneous piezoelectric tensors. The first layer involves three distinct children categories, while the second layer includes two children categories. An initial categorization of parent groups is conducted before classifying the children. The hierarchical approach captures sibling crystal point group relationships and structural inter-class dependencies, with child groups representing the actual crystal point groups and parent groups serving as pseudo-clusters. The response emphasizes that this multi-layer classification approach effectively identifies five distinct sibling crystal point groups.

The dataset is initially divided using the standard train-test split with a ratio of 90:10, allocating 90% of the data for training and 10% for testing. Subsequently, within the training data, a further division is implemented, extracting
<table><caption>**Table A4**

Hyper-parameters used to randomly tune for the best combination.</caption><thead><tr><th>SN</th><th></th><th>Hyper parameters</th></tr></thead><tbody><tr><td>1</td><td>Epochs</td><td>[100, 150, 200]</td></tr><tr><td>2</td><td>batch size</td><td>[30, 40, 50, 60]</td></tr><tr><td>3</td><td>layers</td><td>[4,6]</td></tr><tr><td>4</td><td>neuron size</td><td>[8, 12, 16, 24, 32,48,64]</td></tr><tr><td>5</td><td>optimizer</td><td>[RMSprop, Adam, Nadam]</td></tr><tr><td>6</td><td>activation</td><td>[ReLU, selu, relu, PReLU, LeakyReLU]</td></tr><tr><td>7</td><td>loss function</td><td>[mse, msle]</td></tr><tr><td>8</td><td>drops</td><td>[0.05, 0.075, 0.08, 0.1]</td></tr><tr><td>9</td><td>learning rate</td><td>[0.0001, 0.0005, 0.001, 0.005, 0.01, 0.015, 0.031]</td></tr></tbody></table>

<figure><img src="image_11.png" /><figcaption>**Figure A2:** Learning curve for classification machine learning models</figcaption></figure>

a portion for validation purposes in a 90:10 split. Throughout the training process, a sequential model architecture is employed to train all of the models.

## S4.1. Hyperparameter Tuning Description

In this section, we present the results of hyperparameter tuning for the neural network model. The hyperparameters considered include batch size, optimizer, dropout rate, learning rate, weight initialization, and weight constraint. The randomized search was conducted with 10 iterations and 10-fold cross-validation. The hyperparameter grid explored during the randomized search is as in Table A4:

## S4.2. Neural Network Model Architecture

The NN architecture for the HierCrystalNet is presented in Table A5. The classification model consists of numerous layers on which the dense layer and dropout layer follow sequentially. The model is checked for accuracy metrics. The loss function for two classifications like that of pseudo-parent is binary crossentropy and for actual sub-category classification, it is categorical crossentropy. The learning curve of training versus validation accuracy is presented in Fig. A2.
<table><caption>**Table A5**

Hyperparameters of Classification Models: Multi-level hierarchical machine learning framework</caption><thead><tr><th>Classification Model</th><th>Hyperparameters</th></tr></thead><tbody><tr><td>Pseudo-Parent:<br>Class A<br>or<br>Class B</td><td>Architecture: Sequential<br>Layers: Dense(32), Dropout, Dense(24), Dropout, Dense(16),<br>Dropout, Dense(12), Dropout, Dense(8), Dropout, Dense(1)<br>Dropout Rate: 0.05<br>Learning Rate: 0.005<br>Activation: selu, relu, LeakyReLU, PReLU, relu, softmax<br>Optimizer: adam<br>Loss Function: binary_-crossentropy<br>Epochs: 200<br>Batch Size: 60</td></tr><tr><td>Class A:<br>PG 1<br>or<br>PG2<br>or<br>PG3</td><td>Architecture: Sequential<br>Layers: Dense(48), Dropout(0.05), Dense(16),<br>Dropout(0.05), Dropout(0.1), Dense(8), Dense(3)<br>Activation Functions: relu, relu, LeakyReLU, softmax<br>Optimizer: Adam<br>Learning Rate: 0.015<br>Loss Function: categorical_-crossentropy<br>Epochs: 200<br>Batch Size: 50</td></tr><tr><td>Class B:<br>PG4<br>or<br>PG5</td><td>Architecture: Sequential<br>Layers: Dense(48), Dropout(0.08), Dense(16),<br>Dropout(0.08), Dense(8), Dense(2),<br>Activation Functions: relu, relu, LeakyReLU, softmax<br>Optimizer: Adam<br>Learning Rate: 0.031<br>Loss Function: categorical_-crossentropy<br>Epochs: 200<br>Batch Size: 30</td></tr></tbody></table>

## S5. Crystal rotation supported by tensor mathematics

Neumann's principle asserts that the symmetry of a material's structure limits the number of independent components in tensors describing its physical properties. These properties cannot possess lower symmetry than the crystal structure itself. In crystallography, the choice of standard unit cells and their axes is based on the principal symmetry elements of the crystal system. The alignment of XYZ coordinates used for tensors is established in accordance with the standard crystallographic unit cell axes (Newnham (2004)). Three angles ($\phi$, $\theta$, $\psi$) are required to specify the mutual orientation of two sets of orthogonal axes.

Piezoelectric tensors, being 3rd rank tensors, follow specific transformation laws. The orientation of a molecule or crystallographic lattice in relation to the laboratory coordinate system is determined by three Euler angles: $\phi$, $\theta$, and $\psi$. Two angles, $\theta$ and $\phi$, specify a direction in space, while $\psi$ fixes the orientation of the other two axes. This can be visualized as three consecutive rotations. The direction cosine matrix is commonly used to derive the orientation of a general rotation (Newnham (2004)).

1. First rotation: a counterclockwise rotation of ($\phi$) about Z:

$$a_1 = \begin{bmatrix} \cos \phi & \sin \phi & 0 \\ -\sin \phi & \cos \phi & 0 \\ 0 & 0 & 1 \end{bmatrix}$$
<table><caption>**Table A6**

Neural Network Architectures for five Different modular based ensemble models</caption><thead><tr><th>ModularEnsembleNet</th><th>Model I Architecture</th><th>Model II Architecture</th></tr></thead><tbody><tr><td>PG1</td><td>Dense(128, activation='relu')<br>Dense(64, activation='relu')<br>Dense(1)</td><td>Dense(32, activation='relu')<br>Dense(16, activation='relu')<br>Dense(1)</td></tr><tr><td>PG2</td><td>Dense(128, activation='relu')<br>Dense(64, activation='relu')<br>Dense(3)</td><td>Dense(32, activation='relu')<br>Dense(16, activation='relu')<br>Dense(3)</td></tr><tr><td>PG3</td><td>Dense(128, activation='relu')<br>Dense(64, activation='relu')<br>Dense(2)</td><td>Dense(32, activation='relu')<br>Dense(16, activation='relu')<br>Dense(2)</td></tr><tr><td>PG4</td><td>Dense(128, activation='relu')<br>Dense(64, activation='relu')<br>Dense(3)</td><td>Dense(32, activation='relu')<br>Dense(16, activation='relu')<br>Dense(3)</td></tr><tr><td>PG5</td><td>Dense(28, activation='relu')<br>Dense(18, activation='LeakyReLU')<br>Dropout(0.1)<br>Dense(12, activation='selu')<br>Dense(8, activation='PReLU')<br>Dense(6, activation='relu')<br>Dropout(0.1)<br>Dense(5)</td><td>Dense(30, activation='relu')<br>Dense(16, activation='relu')<br>Dropout(0.1)<br>Dense(12, activation='selu')<br>Dropout(0.1)<br>Dense(5)</td></tr></tbody></table>

2. First rotation: a counterclockwise rotation of $(\theta)$ about X':

$$a_2 = \begin{bmatrix} 1 & 0 & 0 \\ 0 & \cos \theta & \sin \theta \\ 0 & -\sin \theta & \cos \theta \end{bmatrix}$$

3. First rotation: a counterclockwise rotation of $(\psi)$ about Z":

$$a_3 = \begin{bmatrix} \cos \psi & \sin \psi & 0 \\ -\sin \psi & \cos \psi & 0 \\ 0 & 0 & 1 \end{bmatrix}$$

$$A_{ij} = (a_{ij}) = a_3 \cdot a_2 \cdot a_1$$

$$A_{ij} = (a_{ij}) = \begin{bmatrix} (\cos \phi \cos \psi - \cos \theta \sin \phi \sin \psi) & (\cos \phi \sin \psi - \cos \theta \cos \phi \sin \psi) & \sin \theta \sin \psi \\ (-\cos \theta \cos \psi \sin \phi - \cos \phi \sin \psi) & (\cos \theta \cos \phi \cos \psi - \sin \phi \sin \psi) & \cos \theta \sin \psi \\ \sin \theta \sin \phi & -\cos \phi \sin \theta & \cos \theta \end{bmatrix}$$

$$N_{ij} = \begin{bmatrix} a_{11}^2 & a_{21}^2 & a_{31}^2 & 2a_{21}a_{31} & 2a_{31}a_{11} & 2a_{11}a_{21} \\ a_{12}^2 & a_{22}^2 & a_{32}^2 & 2a_{22}a_{32} & 2a_{32}a_{12} & 2a_{12}a_{22} \\ a_{13}^2 & a_{23}^2 & a_{33}^2 & 2a_{23}a_{33} & 2a_{33}a_{13} & 2a_{13}a_{23} \\ a_{12}a_{13} & a_{22}a_{23} & a_{32}a_{33} & a_{22}a_{33} + a_{32}a_{23} & a_{12}a_{33} + a_{32}a_{13} & a_{22}a_{13} + a_{12}a_{23} \\ a_{13}a_{11} & a_{23}a_{21} & a_{33}a_{31} & a_{21}a_{33} + a_{31}a_{23} & a_{31}a_{13} + a_{11}a_{33} & a_{11}a_{23} + a_{21}a_{13} \\ a_{11}a_{12} & a_{21}a_{22} & a_{31}a_{32} & a_{21}a_{32} + a_{31}a_{22} & a_{31}a_{12} + a_{11}a_{32} & a_{11}a_{22} + a_{21}a_{12} \end{bmatrix}$$
The relationship between the $e_{ijk}$ tensor in the laboratory coordinate system and the $e_{ijk}^0$ tensor in the crystal coordinate system (Kalinin et al. (2006)) is

$$e'_{ijk} = A_{il} A_{jm} A_{kn} e_{lmn}^0 \qquad (1)$$

For the second rank 3×6 representations, the rotation relations is given by:

$$e'_{ij} (\phi, \theta, \psi) = A_{ik} e_{kl}^o N_{lj} \qquad (2)$$

## S5.1. Analytical tensor rotation illustration for cubic crystals

The analytical tensor rotation for the cubic system is represented here and all other remaining system follows a similar approach. For a cubic system where the non-zero element in the laboratory coordinate system is $e_{14}$ there is only one set of $(k, l)$ where $k = 1, l = 4$.

$$e'_{11} = A_{11} e_{14}^o N_{41} \qquad e'_{12} = A_{11} e_{14}^o N_{41} \qquad (3)$$

The Einstein summation convention simplifies the notation further by omitting the summation symbol $(\Sigma)$. In this convention, when an index appears as both a subscript and a superscript, it indicates summation.

In the expression $e'_{ij} = \sum_{k=1}^{3} \sum_{l=1}^{6} A_{ik} \cdot e_{kl}^o \cdot N_{lj}$, the indices of A and N should be based on the indices of $e'_{ij}$, $e_{kl}^o$, and N.

For matrix A:

- The first index of A (i) corresponds to the first index of $e'_{ij}$ (i).

- The second index of A (j) corresponds to the first index of $e_{kl}^o$ (k).

For matrix N:

- The first index of N (l) corresponds to the second index of $e_{kl}^o$ (l).

- The second index of N (j) corresponds to the second index of $e'_{ij}$ (j).

The 3 × 6 piezoelectric tensor after crystal rotation is represented for each of the 18 components as below:

$$e'_{11} = A_{11} \cdot e_{11}^o \cdot N_{11} + A_{12} \cdot e_{12}^o \cdot N_{21} + A_{13} \cdot e_{13}^o \cdot N_{31} + A_{11} \cdot e_{14}^o \cdot N_{41} + A_{12} \cdot e_{15}^o \cdot N_{51} + A_{13} \cdot e_{16}^o \cdot N_{61} \\ + A_{11} \cdot e_{21}^o \cdot N_{11} + A_{12} \cdot e_{22}^o \cdot N_{21} + A_{13} \cdot e_{23}^o \cdot N_{31} + A_{11} \cdot e_{24}^o \cdot N_{41} + A_{12} \cdot e_{25}^o \cdot N_{51} + A_{13} \cdot e_{26}^o \cdot N_{61} \\ + A_{11} \cdot e_{31}^o \cdot N_{11} + A_{12} \cdot e_{32}^o \cdot N_{21} + A_{13} \cdot e_{33}^o \cdot N_{31} + A_{11} \cdot e_{34}^o \cdot N_{41} + A_{12} \cdot e_{35}^o \cdot N_{51} + A_{13} \cdot e_{36}^o \cdot N_{61}$$

$$e'_{12} = A_{11} \cdot e_{11}^o \cdot N_{21} + A_{11} \cdot e_{12}^o \cdot N_{21} + A_{11} \cdot e_{13}^o \cdot N_{21} + A_{11} \cdot e_{14}^o \cdot N_{21} + A_{11} \cdot e_{15}^o \cdot N_{21} + A_{11} \cdot e_{16}^o \cdot N_{21} \\ + A_{12} \cdot e_{21}^o \cdot N_{21} + A_{12} \cdot e_{22}^o \cdot N_{21} + A_{12} \cdot e_{23}^o \cdot N_{21} + A_{12} \cdot e_{24}^o \cdot N_{21} + A_{12} \cdot e_{25}^o \cdot N_{21} + A_{12} \cdot e_{26}^o \cdot N_{21} \\ + A_{13} \cdot e_{31}^o \cdot N_{21} + A_{13} \cdot e_{32}^o \cdot N_{21} + A_{13} \cdot e_{33}^o \cdot N_{21} + A_{13} \cdot e_{34}^o \cdot N_{21} + A_{13} \cdot e_{35}^o \cdot N_{21} + A_{13} \cdot e_{36}^o \cdot N_{21}$$

$$e'_{13} = A_{11} \cdot e_{11}^o \cdot N_{31} + A_{11} \cdot e_{12}^o \cdot N_{31} + A_{11} \cdot e_{13}^o \cdot N_{31} + A_{11} \cdot e_{14}^o \cdot N_{31} + A_{11} \cdot e_{15}^o \cdot N_{31} + A_{11} \cdot e_{16}^o \cdot N_{31} \\ + A_{12} \cdot e_{21}^o \cdot N_{31} + A_{12} \cdot e_{22}^o \cdot N_{31} + A_{12} \cdot e_{23}^o \cdot N_{31} + A_{12} \cdot e_{24}^o \cdot N_{31} + A_{12} \cdot e_{25}^o \cdot N_{31} + A_{12} \cdot e_{26}^o \cdot N_{31} \\ + A_{13} \cdot e_{31}^o \cdot N_{31} + A_{13} \cdot e_{32}^o \cdot N_{31} + A_{13} \cdot e_{33}^o \cdot N_{31} + A_{13} \cdot e_{34}^o \cdot N_{31} + A_{13} \cdot e_{35}^o \cdot N_{31} + A_{13} \cdot e_{36}^o \cdot N_{31}$$

$$e'_{14} = A_{11} \cdot e_{11}^o \cdot N_{41} + A_{11} \cdot e_{12}^o \cdot N_{41} + A_{11} \cdot e_{13}^o \cdot N_{41} + A_{11} \cdot e_{14}^o \cdot N_{41} + A_{11} \cdot e_{15}^o \cdot N_{41} + A_{11} \cdot e_{16}^o \cdot N_{41} \\ + A_{12} \cdot e_{21}^o \cdot N_{41} + A_{12} \cdot e_{22}^o \cdot N_{41} + A_{12} \cdot e_{23}^o \cdot N_{41} + A_{12} \cdot e_{24}^o \cdot N_{41} + A_{12} \cdot e_{25}^o \cdot N_{41} + A_{12} \cdot e_{26}^o \cdot N_{41} \\ + A_{13} \cdot e_{31}^o \cdot N_{41} + A_{13} \cdot e_{32}^o \cdot N_{41} + A_{13} \cdot e_{33}^o \cdot N_{41} + A_{13} \cdot e_{34}^o \cdot N_{41} + A_{13} \cdot e_{35}^o \cdot N_{41} + A_{13} \cdot e_{36}^o \cdot N_{41}$$

$$e'_{15} = A_{11} \cdot e_{11}^o \cdot N_{51} + A_{11} \cdot e_{12}^o \cdot N_{51} + A_{11} \cdot e_{13}^o \cdot N_{51} + A_{11} \cdot e_{14}^o \cdot N_{51} + A_{11} \cdot e_{15}^o \cdot N_{51} + A_{11} \cdot e_{16}^o \cdot N_{51} \\ + A_{12} \cdot e_{21}^o \cdot N_{51} + A_{12} \cdot e_{22}^o \cdot N_{51} + A_{12} \cdot e_{23}^o \cdot N_{51} + A_{12} \cdot e_{24}^o \cdot N_{51} + A_{12} \cdot e_{25}^o \cdot N_{51} + A_{12} \cdot e_{26}^o \cdot N_{51} \\ + A_{13} \cdot e_{31}^o \cdot N_{51} + A_{13} \cdot e_{32}^o \cdot N_{51} + A_{13} \cdot e_{33}^o \cdot N_{51} + A_{13} \cdot e_{34}^o \cdot N_{51} + A_{13} \cdot e_{35}^o \cdot N_{51} + A_{13} \cdot e_{36}^o \cdot N_{51}$$

$$e'_{16} = A_{11} \cdot e_{11}^o \cdot N_{61} + A_{11} \cdot e_{12}^o \cdot N_{61} + A_{11} \cdot e_{13}^o \cdot N_{61} + A_{11} \cdot e_{14}^o \cdot N_{61} + A_{11} \cdot e_{15}^o \cdot N_{61} + A_{11} \cdot e_{16}^o \cdot N_{61} \\ + A_{12} \cdot e_{21}^o \cdot N_{61} + A_{12} \cdot e_{22}^o \cdot N_{61} + A_{12} \cdot e_{23}^o \cdot N_{61} + A_{12} \cdot e_{24}^o \cdot N_{61} + A_{12} \cdot e_{25}^o \cdot N_{61} + A_{12} \cdot e_{26}^o \cdot N_{61} \\ + A_{13} \cdot e_{31}^o \cdot N_{61} + A_{13} \cdot e_{32}^o \cdot N_{61} + A_{13} \cdot e_{33}^o \cdot N_{61} + A_{13} \cdot e_{34}^o \cdot N_{61} + A_{13} \cdot e_{35}^o \cdot N_{61} + A_{13} \cdot e_{36}^o \cdot N_{61}$$
$$e'_{21} = A_{21} \cdot e_{11}^{0} \cdot N_{11} + A_{21} \cdot e_{12}^{0} \cdot N_{21} + A_{21} \cdot e_{13}^{0} \cdot N_{31} + A_{21} \cdot e_{14}^{0} \cdot N_{41} + A_{21} \cdot e_{15}^{0} \cdot N_{51} + A_{21} \cdot e_{16}^{0} \cdot N_{61} \\ + A_{22} \cdot e_{21}^{0} \cdot N_{11} + A_{22} \cdot e_{22}^{0} \cdot N_{21} + A_{22} \cdot e_{23}^{0} \cdot N_{31} + A_{22} \cdot e_{24}^{0} \cdot N_{41} + A_{22} \cdot e_{25}^{0} \cdot N_{51} + A_{22} \cdot e_{26}^{0} \cdot N_{61} \\ + A_{23} \cdot e_{31}^{0} \cdot N_{11} + A_{23} \cdot e_{32}^{0} \cdot N_{21} + A_{23} \cdot e_{33}^{0} \cdot N_{31} + A_{23} \cdot e_{34}^{0} \cdot N_{41} + A_{23} \cdot e_{35}^{0} \cdot N_{51} + A_{23} \cdot e_{36}^{0} \cdot N_{61}$$

$$e'_{22} = A_{21} \cdot e_{11}^{0} \cdot N_{12} + A_{21} \cdot e_{12}^{0} \cdot N_{22} + A_{21} \cdot e_{13}^{0} \cdot N_{32} + A_{21} \cdot e_{14}^{0} \cdot N_{42} + A_{21} \cdot e_{15}^{0} \cdot N_{52} + A_{21} \cdot e_{16}^{0} \cdot N_{62} \\ + A_{22} \cdot e_{21}^{0} \cdot N_{12} + A_{22} \cdot e_{22}^{0} \cdot N_{22} + A_{22} \cdot e_{23}^{0} \cdot N_{32} + A_{22} \cdot e_{24}^{0} \cdot N_{42} + A_{22} \cdot e_{25}^{0} \cdot N_{52} + A_{22} \cdot e_{26}^{0} \cdot N_{62} \\ + A_{23} \cdot e_{31}^{0} \cdot N_{12} + A_{23} \cdot e_{32}^{0} \cdot N_{22} + A_{23} \cdot e_{33}^{0} \cdot N_{32} + A_{23} \cdot e_{34}^{0} \cdot N_{42} + A_{23} \cdot e_{35}^{0} \cdot N_{52} + A_{23} \cdot e_{36}^{0} \cdot N_{62}$$

$$e'_{23} = A_{21} \cdot e_{11}^{0} \cdot N_{13} + A_{21} \cdot e_{12}^{0} \cdot N_{23} + A_{21} \cdot e_{13}^{0} \cdot N_{33} + A_{21} \cdot e_{14}^{0} \cdot N_{43} + A_{21} \cdot e_{15}^{0} \cdot N_{53} + A_{21} \cdot e_{16}^{0} \cdot N_{63} \\ + A_{22} \cdot e_{21}^{0} \cdot N_{13} + A_{22} \cdot e_{22}^{0} \cdot N_{23} + A_{22} \cdot e_{23}^{0} \cdot N_{33} + A_{22} \cdot e_{24}^{0} \cdot N_{43} + A_{22} \cdot e_{25}^{0} \cdot N_{53} + A_{22} \cdot e_{26}^{0} \cdot N_{63} \\ + A_{23} \cdot e_{31}^{0} \cdot N_{13} + A_{23} \cdot e_{32}^{0} \cdot N_{23} + A_{23} \cdot e_{33}^{0} \cdot N_{33} + A_{23} \cdot e_{34}^{0} \cdot N_{43} + A_{23} \cdot e_{35}^{0} \cdot N_{53} + A_{23} \cdot e_{36}^{0} \cdot N_{63}$$

$$e'_{24} = A_{21} \cdot e_{11}^{0} \cdot N_{14} + A_{21} \cdot e_{12}^{0} \cdot N_{24} + A_{21} \cdot e_{13}^{0} \cdot N_{34} + A_{21} \cdot e_{14}^{0} \cdot N_{44} + A_{21} \cdot e_{15}^{0} \cdot N_{54} + A_{21} \cdot e_{16}^{0} \cdot N_{64} \\ + A_{22} \cdot e_{21}^{0} \cdot N_{14} + A_{22} \cdot e_{22}^{0} \cdot N_{24} + A_{22} \cdot e_{23}^{0} \cdot N_{34} + A_{22} \cdot e_{24}^{0} \cdot N_{44} + A_{22} \cdot e_{25}^{0} \cdot N_{54} + A_{22} \cdot e_{26}^{0} \cdot N_{64} \\ + A_{23} \cdot e_{31}^{0} \cdot N_{14} + A_{23} \cdot e_{32}^{0} \cdot N_{24} + A_{23} \cdot e_{33}^{0} \cdot N_{34} + A_{23} \cdot e_{34}^{0} \cdot N_{44} + A_{23} \cdot e_{35}^{0} \cdot N_{54} + A_{23} \cdot e_{36}^{0} \cdot N_{64}$$

$$e'_{25} = A_{21} \cdot e_{11}^{0} \cdot N_{15} + A_{21} \cdot e_{12}^{0} \cdot N_{25} + A_{21} \cdot e_{13}^{0} \cdot N_{35} + A_{21} \cdot e_{14}^{0} \cdot N_{45} + A_{21} \cdot e_{15}^{0} \cdot N_{55} + A_{21} \cdot e_{16}^{0} \cdot N_{65} \\ + A_{22} \cdot e_{21}^{0} \cdot N_{15} + A_{22} \cdot e_{22}^{0} \cdot N_{25} + A_{22} \cdot e_{23}^{0} \cdot N_{35} + A_{22} \cdot e_{24}^{0} \cdot N_{45} + A_{22} \cdot e_{25}^{0} \cdot N_{55} + A_{22} \cdot e_{26}^{0} \cdot N_{65} \\ + A_{23} \cdot e_{31}^{0} \cdot N_{15} + A_{23} \cdot e_{32}^{0} \cdot N_{25} + A_{23} \cdot e_{33}^{0} \cdot N_{35} + A_{23} \cdot e_{34}^{0} \cdot N_{45} + A_{23} \cdot e_{35}^{0} \cdot N_{55} + A_{23} \cdot e_{36}^{0} \cdot N_{65}$$

$$e'_{26} = A_{21} \cdot e_{11}^{0} \cdot N_{16} + A_{21} \cdot e_{12}^{0} \cdot N_{26} + A_{21} \cdot e_{13}^{0} \cdot N_{36} + A_{21} \cdot e_{14}^{0} \cdot N_{46} + A_{21} \cdot e_{15}^{0} \cdot N_{56} + A_{21} \cdot e_{16}^{0} \cdot N_{66} \\ + A_{22} \cdot e_{21}^{0} \cdot N_{16} + A_{22} \cdot e_{22}^{0} \cdot N_{26} + A_{22} \cdot e_{23}^{0} \cdot N_{36} + A_{22} \cdot e_{24}^{0} \cdot N_{46} + A_{22} \cdot e_{25}^{0} \cdot N_{56} + A_{22} \cdot e_{26}^{0} \cdot N_{66} \\ + A_{23} \cdot e_{31}^{0} \cdot N_{16} + A_{23} \cdot e_{32}^{0} \cdot N_{26} + A_{23} \cdot e_{33}^{0} \cdot N_{36} + A_{23} \cdot e_{34}^{0} \cdot N_{46} + A_{23} \cdot e_{35}^{0} \cdot N_{56} + A_{23} \cdot e_{36}^{0} \cdot N_{66}$$

$$e'_{31} = A_{31} \cdot e_{11}^{0} \cdot N_{11} + A_{31} \cdot e_{12}^{0} \cdot N_{21} + A_{31} \cdot e_{13}^{0} \cdot N_{31} + A_{31} \cdot e_{14}^{0} \cdot N_{41} + A_{31} \cdot e_{15}^{0} \cdot N_{51} + A_{31} \cdot e_{16}^{0} \cdot N_{61} \\ + A_{32} \cdot e_{21}^{0} \cdot N_{11} + A_{32} \cdot e_{22}^{0} \cdot N_{21} + A_{32} \cdot e_{23}^{0} \cdot N_{31} + A_{32} \cdot e_{24}^{0} \cdot N_{41} + A_{32} \cdot e_{25}^{0} \cdot N_{51} + A_{32} \cdot e_{26}^{0} \cdot N_{61} \\ + A_{33} \cdot e_{31}^{0} \cdot N_{11} + A_{33} \cdot e_{32}^{0} \cdot N_{21} + A_{33} \cdot e_{33}^{0} \cdot N_{31} + A_{33} \cdot e_{34}^{0} \cdot N_{41} + A_{33} \cdot e_{35}^{0} \cdot N_{51} + A_{33} \cdot e_{36}^{0} \cdot N_{61}$$

$$e'_{32} = A_{31} \cdot e_{11}^{0} \cdot N_{12} + A_{31} \cdot e_{12}^{0} \cdot N_{22} + A_{31} \cdot e_{13}^{0} \cdot N_{32} + A_{31} \cdot e_{14}^{0} \cdot N_{42} + A_{31} \cdot e_{15}^{0} \cdot N_{52} + A_{31} \cdot e_{16}^{0} \cdot N_{62} \\ + A_{32} \cdot e_{21}^{0} \cdot N_{12} + A_{32} \cdot e_{22}^{0} \cdot N_{22} + A_{32} \cdot e_{23}^{0} \cdot N_{32} + A_{32} \cdot e_{24}^{0} \cdot N_{42} + A_{32} \cdot e_{25}^{0} \cdot N_{52} + A_{32} \cdot e_{26}^{0} \cdot N_{62} \\ + A_{33} \cdot e_{31}^{0} \cdot N_{12} + A_{33} \cdot e_{32}^{0} \cdot N_{22} + A_{33} \cdot e_{33}^{0} \cdot N_{32} + A_{33} \cdot e_{34}^{0} \cdot N_{42} + A_{33} \cdot e_{35}^{0} \cdot N_{52} + A_{33} \cdot e_{36}^{0} \cdot N_{62}$$

$$e'_{33} = A_{31} \cdot e_{11}^{0} \cdot N_{13} + A_{31} \cdot e_{12}^{0} \cdot N_{23} + A_{31} \cdot e_{13}^{0} \cdot N_{33} + A_{31} \cdot e_{14}^{0} \cdot N_{43} + A_{31} \cdot e_{15}^{0} \cdot N_{53} + A_{31} \cdot e_{16}^{0} \cdot N_{63} \\ + A_{32} \cdot e_{21}^{0} \cdot N_{13} + A_{32} \cdot e_{22}^{0} \cdot N_{23} + A_{32} \cdot e_{23}^{0} \cdot N_{33} + A_{32} \cdot e_{24}^{0} \cdot N_{43} + A_{32} \cdot e_{25}^{0} \cdot N_{53} + A_{32} \cdot e_{26}^{0} \cdot N_{63} \\ + A_{33} \cdot e_{31}^{0} \cdot N_{13} + A_{33} \cdot e_{32}^{0} \cdot N_{23} + A_{33} \cdot e_{33}^{0} \cdot N_{33} + A_{33} \cdot e_{34}^{0} \cdot N_{43} + A_{33} \cdot e_{35}^{0} \cdot N_{53} + A_{33} \cdot e_{36}^{0} \cdot N_{63}$$

$$e'_{34} = A_{31} \cdot e_{11}^{0} \cdot N_{14} + A_{31} \cdot e_{12}^{0} \cdot N_{24} + A_{31} \cdot e_{13}^{0} \cdot N_{34} + A_{31} \cdot e_{14}^{0} \cdot N_{44} + A_{31} \cdot e_{15}^{0} \cdot N_{54} + A_{31} \cdot e_{16}^{0} \cdot N_{64} \\ + A_{32} \cdot e_{21}^{0} \cdot N_{14} + A_{32} \cdot e_{22}^{0} \cdot N_{24} + A_{32} \cdot e_{23}^{0} \cdot N_{34} + A_{32} \cdot e_{24}^{0} \cdot N_{44} + A_{32} \cdot e_{25}^{0} \cdot N_{54} + A_{32} \cdot e_{26}^{0} \cdot N_{64} \\ + A_{33} \cdot e_{31}^{0} \cdot N_{14} + A_{33} \cdot e_{32}^{0} \cdot N_{24} + A_{33} \cdot e_{33}^{0} \cdot N_{34} + A_{33} \cdot e_{34}^{0} \cdot N_{44} + A_{33} \cdot e_{35}^{0} \cdot N_{54} + A_{33} \cdot e_{36}^{0} \cdot N_{64}$$

$$e'_{35} = A_{31} \cdot e_{11}^{0} \cdot N_{15} + A_{31} \cdot e_{12}^{0} \cdot N_{25} + A_{31} \cdot e_{13}^{0} \cdot N_{35} + A_{31} \cdot e_{14}^{0} \cdot N_{45} + A_{31} \cdot e_{15}^{0} \cdot N_{55} + A_{31} \cdot e_{16}^{0} \cdot N_{65} \\ + A_{32} \cdot e_{21}^{0} \cdot N_{15} + A_{32} \cdot e_{22}^{0} \cdot N_{25} + A_{32} \cdot e_{23}^{0} \cdot N_{35} + A_{32} \cdot e_{24}^{0} \cdot N_{45} + A_{32} \cdot e_{25}^{0} \cdot N_{55} + A_{32} \cdot e_{26}^{0} \cdot N_{65} \\ + A_{33} \cdot e_{31}^{0} \cdot N_{15} + A_{33} \cdot e_{32}^{0} \cdot N_{25} + A_{33} \cdot e_{33}^{0} \cdot N_{35} + A_{33} \cdot e_{34}^{0} \cdot N_{45} + A_{33} \cdot e_{35}^{0} \cdot N_{55} + A_{33} \cdot e_{36}^{0} \cdot N_{65}$$

$$e'_{36} = A_{31} \cdot e_{11}^{0} \cdot N_{16} + A_{31} \cdot e_{12}^{0} \cdot N_{26} + A_{31} \cdot e_{13}^{0} \cdot N_{36} + A_{31} \cdot e_{14}^{0} \cdot N_{46} + A_{31} \cdot e_{15}^{0} \cdot N_{56} + A_{31} \cdot e_{16}^{0} \cdot N_{66} \\ + A_{32} \cdot e_{21}^{0} \cdot N_{16} + A_{32} \cdot e_{22}^{0} \cdot N_{26} + A_{32} \cdot e_{23}^{0} \cdot N_{36} + A_{32} \cdot e_{24}^{0} \cdot N_{46} + A_{32} \cdot e_{25}^{0} \cdot N_{56} + A_{32} \cdot e_{26}^{0} \cdot N_{66} \\ + A_{33} \cdot e_{31}^{0} \cdot N_{16} + A_{33} \cdot e_{32}^{0} \cdot N_{26} + A_{33} \cdot e_{33}^{0} \cdot N_{36} + A_{33} \cdot e_{34}^{0} \cdot N_{46} + A_{33} \cdot e_{35}^{0} \cdot N_{56} + A_{33} \cdot e_{36}^{0} \cdot N_{66}$$
## S6. FEM for orthotropic stiffness

Orthotropic materials exhibit unique and independent mechanical properties along three mutually perpendicular directions. This study employs a stiffness matrix for the characterization of orthotropic materials, mirroring the considerations and setup outlined in Section 3.1 of the main article. Specifically, a load of $(F_z)_{max} = 5 \times 10^{-5} \text{ N/m}^2$ is applied solely to investigate the piezoelectric performance under orthotropic assumptions. Young's modulus along the applied load i.e., along the z-axis $(C_{33})$ from the orthotropic stiffness matrix was used in frequency calculations. The stiffness $C_{ijkl}^{d-AIN}$ is used for all three cases of alloyed and rotated crystals.

Figure A3 a) illustrates the voltage generation results in volts for four different tensors: pure AlN, doped $B_{0.3}Er_{0.5}Al_{0.2}N$, doped AlN with Type-1 rotation, and doped AlN with Type-2 rotation. For pure AlN, the deformation and first mode frequency are 0.79 µm and 70.71 kHz, respectively. Conversely, for $B_{0.3}Er_{0.5}Al_{0.2}N$, these values are 0.88 µm and 56 kHz, as depicted in Fig A3 b). The results encompassing power and voltage for all four setups are presented in Fig A3 c).

### Tensor Values Used

Doped: $B_{0.3}Er_{0.5}Al_{0.2}N \rightarrow e^*$

$$ \begin{bmatrix} 0 & 0 & 0 & 0 & 0.2283 & 0 \\ 0 & 0 & 0 & 0.2283 & 0 & 0 \\ 0.0943 & 0.0943 & 4.363 & 0 & 0 & 0 \end{bmatrix} $$

Rotated- Type1 : ($\theta = 0, \psi = 268.47, 91.52; 0 \le \phi \le 360$) $\rightarrow e_1^{*'}$

$$ \begin{bmatrix} 0 & 0.08 & 0 & -0.0808 & 0.459 & -0.4554 \\ -0.1882 & -4.563 & -4.3232 & 8.7947 & -0.0807 & 0.146 \\ 0.1899 & 4.3775 & 4.3636 & -8.6464 & 0 & -0.0667 \end{bmatrix} $$

Rotated-Type 2 : ($\theta = 180, \psi = 268.47, 91.52; \phi = 0 \le \phi \le 360$) $\rightarrow e_2^{*'}$

$$ \begin{bmatrix} 1.1354 & 0.947 & 0.9856 & -1.825 & -2.005 & 1.9339 \\ 0.9564 & 1.098 & 0.9799 & -1.9783 & -1.8336 & 1.93 \\ -0.9946 & -0.9778 & -1.0156 & 1.9799 & 2.0 & -1.838 \end{bmatrix} $$

### Orthotropic Stiffness matrix

Pure AlN: $E_1 = 275 \text{ GPa}, E_2 = 300 \text{ GPa}, E_3 = 325 \text{ GPa}; v_{12} = 0.26, v_{23} = 0.24, v_{31} = 0.22$

$$ C_{ijkl}^{AlN} = \begin{bmatrix} 3.24 & 1.17 & 1.15 & 0 & 0 & 0 \\ 1.08 & 3.57 & 1.21 & 0 & 0 & 0 \\ 0.972 & 1.12 & 3.79 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1.31 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1.13 & 0 \\ 0 & 0 & 0 & 0 & 0 & 1.19 \end{bmatrix} \times 10^{11} \text{N/m}^2 $$

$B_{0.3}Er_{0.5}Al_{0.2}N: E_1 = 250 \text{ GPa}, E_2 = 225 \text{ GPa}, E_3 = 275 \text{ GPa}; v_{12} = 0.22, v_{23} = 0.24, v_{31} = 0.26$

$$ C_{ijkl}^{d-AlN} = \begin{bmatrix} 2.95 & 0.795 & 0.783 & 0 & 0 & 0 \\ 0.884 & 2.63 & 0.744 & 0 & 0 & 0 \\ 0.979 & 0.837 & 2.38 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0.806 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0.992 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0.922 \end{bmatrix} \times 10^{11} \text{N/m}^2 $$
<figure><img src="image_12.png" /><figcaption>**Figure A3:** FEM results of generated voltage and corresponding deformation for the orthotopic stiffness matrix. The four different tensors considered are pure AlN, doped $B_{0.3}Er_{0.5}Al_{0.2}N$, doped AlN with Type-1 rotation, and doped AlN with Type-2 rotation.</figcaption></figure>

## References

Kalinin, S. V., Rodriguez, B. J., Jesse, S., Shin, J., Baddorf, A. P., Gupta, P., Jain, H., Williams, D. B., & Gruverman, A. (2006). Vector piezoresponse force microscopy. *Microscopy and Microanalysis*, **12**, 206–220.

Li, Y., Dong, R., Yang, W., & Hu, J. (2021). Composition based crystal materials symmetry prediction using machine learning with enhanced descriptors. *Computational Materials Science*, **198**, 110686.

Newnham, R. E. (2004). *Properties of Materials: Anisotropy, Symmetry, Structure*. Oxford University Press.

Prince, E. (Ed.) (2004). *International Tables for Crystallography, Volume C: Mathematical, Physical and Chemical Tables, 3rd Edition*. Wiley.

Trolier-McKinstry, S. (2008). Crystal chemistry of piezoelectric materials. In A. Safari, & E. K. Akdoğan (Eds.), *Piezoelectric and Acoustic Materials for Transducer Applications* (pp. 39–56). Springer.

Ward, L., Agrawal, A., Choudhary, A., & Wolverton, C. (2016). A general-purpose machine learning framework for predicting properties of inorganic materials. *npj Computational Materials*, **2**, 16028.