Discovery of Novel Coumarin-Schiff Base Hybrids as Potential Acetylcholinesterase Inhibitors: Design, Synthesis, Enzyme Inhibition, and Computational Studies

To discover anti-acetylcholinesterase agents for the treatment of Alzheimer’s disease (AD), a series of novel Schiff base-coumarin hybrids was rationally designed, synthesized successfully, and structurally characterized using Fourier transform infrared (FTIR), Nuclear magnetic resonance (NMR), and High-Resolution Mass Spectrometry (HRMS) analyses. These hybrids were evaluated for their potential inhibitory effect on acetylcholinesterase (AChE). All of them exhibited excellent inhibitory activity against AChE. The IC50 values ranged from 87.84 to 515.59 μg/mL; hybrids 13c and 13d with IC50 values of 0.232 ± 0.011 and 0.190 ± 0.004 µM, respectively, showed the most potent activity as acetylcholinesterase inhibitors (AChEIs). The reference drug, Galantamine, yielded an IC50 of 1.142 ± 0.027 µM. Reactivity descriptors, including chemical potential (μ), chemical hardness (η), electrophilicity (ω), condensed Fukui function, and dual descriptors are calculated at wB97XD/6-311++ G (d,p) to identify reactivity changes of the designed compounds. An in-depth investigation of the natural charge pattern of the studied compounds led to a deep understanding of the important interaction centers between these compounds and the biological receptors of AChE. The molecular electrostatic surface potential (MESP) of the most active site in these derivatives was determined using high-quality information and visualization. Molecular docking analysis was performed to predict binding sites and binding energies. The structure-activity-property relationship studies indicated that the proposed compounds exhibit good oral bioavailability properties. To explore the stability and dynamic behavior of the ligand-receptor complexes, molecular dynamics simulations (MDS) were performed for 100 ns on the two best docked derivatives, 13c and 13d, with the AChE (4EY7) receptor. A popular method for determining the free binding energies (MM/GBSA) is performed using snapshots taken from the systems’ trajectories at 100 ns. These results revealed that the complex system of compound 13d acquired a relatively more stable conformation and exhibited better descriptors than the complex system of compound 13c and the Galantamine drug, suggesting its potential as an effective inhibiting drug. The binding free energy analysis revealed that the 13d-4EY7 complex exhibited greater stability with AChE receptors compared to other complexes.


Introduction
Alzheimer's disease (AD) is the main cause of dementia and one of the most prevalent neurodegenerative diseases that impact people around the world. It is characterized by progressive memory loss, severe behavioral abnormalities, and cognitive impairments [1]. By 2050, the number of AD patients will have tripled from the present level of more than 50 million worldwide [2]. The cholinergic hypothesis is one of the most significant hypotheses, despite the fact that the etiology of AD has not yet been fully explained. In the brains of AD patients, it was found that the level of the neuromodulator acetylcholine (ACh) was abnormal. Studies have shown that acetylcholine levels rise when acetylcholinesterase (AChE) is inhibited, thereby enhancing memory and cognitive function in patients [3]. Therefore, a major treatment strategy for AD involves inhibiting the acetylcholinesterase (AChE) enzyme, which catalyzes the hydrolysis of ACh neurotransmitters. Currently, drugs such as donepezil, rivastigmine, and galantamine, which belong to the class of acetylcholinesterase inhibitors (AChEIs), are approved and administered as therapeutic options for the management of AD [4]. However, while these medications can provide temporary improvements in memory and cognitive function, they are unable to prevent or reverse the progression of the disease. Given the limitations of current treatments and the urgent need to address both symptomatic relief and disease progression, researchers are continually seeking and developing potential agents that can effectively treat AD and potentially modify its course.Heterocyclic-based compounds have been reported as active inhibitors against the acetylcholinesterase enzyme (AChE). Medicinal chemists have developed and used a variety of coumarin cores to design novel treatments with a wide spectrum of pharmacological properties [4]. Mohammadi-Khanaposhtani et al. have investigated the inhibition of AChE with a series of coumarin-3-carboxamide-N-morpholine derivatives. The obtained IC 50 values ranged in micromolar levels (6.21-26.4 µM). Compounds 1 and 2 were found to be the most active among the series [5]. Kara et al. reported the synthesis and AChE inhibition activity of a series of coumarin-based compounds. Among the series, derivative 3, with an IC 50 value of 0.04 µM, displayed higher activity [6]. On the other hand, Schiff bases represent an important group of organic compounds with various biological activities [7][8][9][10][11]. The study of novel biologically significant Schiff bases has been drawing the attention of chemists and pharmacists [12][13][14][15]. Wang et al. have designed and synthesized a series of Schiff bases containing coumarin cores. The synthesized compounds were assessed for their anti-Alzheimer activity. The obtained results showed that derivative 4 showed higher inhibitory activities with IC 50 values of 0.673 and 0.711 µM for human monoamine oxidasehMAO-A and hMAO-B, respectively [16]. Compound 5 showed promising inhibitory activity against AChE with an IC 50 value of 4.12 µM [17]. In 2019, 16 hybrids of acrine-isatin Schiff bases were designed and synthesized by Riazimontazer and co-workers. The designed compounds were tested as potential anti-Alzheimer agents. Compound 6, with an IC 50 value of 0.42 nM, was the most potent (IC 50 = 38.72 nM). Interestingly, chloro-substituent compounds displayed stronger inhibitory activity than those with nitro and methoxy substituents [18]. Figure 1 shows that some coumarin and Schiff base-containing compounds were reported as AChE inhibitors. Our research group has reported chalcone-based coumarin scaffolds and psoralen derivatives as potential acetylcholinesterase inhibitors [19,20] in our ongoing search for a lead compound in drug discovery, as illustrated in Figure 1. This concept led us to design and synthesize a novel series of coumarin-Schiff base hybrids. This may contribute to discovering more effective inhibitors against AChE. Molecular modeling analyses of the synthesized compounds were performed to clarify their interaction modes and stabilities with the amino acid residues of AChE and explore their pharmacological effects. Finally, comprehensive reports were provided on the density functional theory (DFT) calculations, in silico assessments of local and global reactivity, and drug-likeness predictions using ADMET for the designed hybrids.

Chemistry
The Pechmann condensation was performed between 3-methoxyphenol 7 and ethyl 4-chloroacetoacetate 8 in the presence of concentrated sulfuric acid. This reaction resulted in the formation of 4-(chloromethyl)-7-methoxy-2H-chromen-2-one 9 as a white powder with a yield of 94.0% (Scheme 1). Schiff bases (12a-j) were synthesized via an acidcatalyzed condensation reaction between commercially available hydroxybenzaldehydes (11a-d) and various aniline substituents (10a-d) in ethanol for 24 h. These compounds were synthesized in moderate to good yields (56.7-89.6%). Our research group has reported chalcone-based coumarin scaffolds and psoralen derivatives as potential acetylcholinesterase inhibitors [19,20] in our ongoing search for a lead compound in drug discovery, as illustrated in Figure 1. This concept led us to design and synthesize a novel series of coumarin-Schiff base hybrids. This may contribute to discovering more effective inhibitors against AChE. Molecular modeling analyses of the synthesized compounds were performed to clarify their interaction modes and stabilities with the amino acid residues of AChE and explore their pharmacological effects. Finally, comprehensive reports were provided on the density functional theory (DFT) calculations, in silico assessments of local and global reactivity, and drug-likeness predictions using ADMET for the designed hybrids.
The formation of Schiff bases (12a-j) is illustrated in Scheme 1. The presence of the expected CH=N group in these imines is strong evidence for the formation of target compounds and could be confirmed by IR and 1 H NMR. The IR spectra of compounds (12a-j) showed characteristic C=N bands from (1603-1627 cm −1 ). Moreover, their 1 H NMR spectra exhibited singlet signals in the range (δ 8. 44-8.54). Thus, they are in agreement with previous studies [21][22][23]. Finally, the newly synthesized 4-(chloromethyl)-7-methoxy-2Hchromen-2-one 9 and Schiff bases (12a-j) were combined to form coumarin-imine hybrids (13a-j) following the procedure outlined in Scheme 2. The formation of Schiff bases (12a-j) is illustrated in Scheme 1. The presence of the expected CH=N group in these imines is strong evidence for the formation of target compounds and could be confirmed by IR and 1 H NMR. The IR spectra of compounds (12a-j) showed characteristic C=N bands from (1603-1627 cm −1 ). Moreover, their 1 H NMR spectra exhibited singlet signals in the range (δ 8. 44-8.54). Thus, they are in agreement with previous studies [21][22][23]. Finally, the newly synthesized 4-(chloromethyl)-7methoxy-2H-chromen-2-one 9 and Schiff bases (12a-j) were combined to form coumarinimine hybrids (13a-j) following the procedure outlined in Scheme 2.

Anti-Acetylcholinesterase of Coumarin-Imine Hybrids (13a-j)
The in vitro activity of the coumarin-imine series (13a-j) was mentioned in Table 1. The study found that all hybrids (except for compound 13h) had excellent inhibition activity. It was observed that the compounds with chloro substituent, which reduced electron density on the aromatic ring [24], exhibited better biological activity in comparison to those without a substituent or with a methoxy group at C-4 of the ring-B ( Figure 2). As for chloro-hybrids, 13c and 13d showed activity at 0.232 and 0.190 μM, respectively. These two compounds were more potent than galantamine by ~5-fold. However, the IC50 value calculated for coumarin 13h (1.175 μM) was lower than the positive control. The in vitro activity of the coumarin-imine series (13a-j) was mentioned in Table 1. The study found that all hybrids (except for compound 13h) had excellent inhibition activity. It was observed that the compounds with chloro substituent, which reduced electron density on the aromatic ring [24], exhibited better biological activity in comparison to those without a substituent or with a methoxy group at C-4 of the ring-B ( Figure 2). As for chloro-hybrids, 13c and 13d showed activity at 0.232 and 0.190 µM, respectively. These two compounds were more potent than galantamine by~5-fold. However, the IC 50 value calculated for coumarin 13h (1.175 µM) was lower than the positive control.  Furthermore, the introduction of a methoxy group on the ring-B, as in compounds 13b, 13f, 13g, and 13h, decreased anti-AChE activity, whereas the inhibitory activity against the enzyme was increased by the incorporation of a chloro group, as in compounds 13c, 13d, 13i and 13j. That is worth mentioning; the connection of the imine Furthermore, the introduction of a methoxy group on the ring-B, as in compounds 13b, 13f, 13g, and 13h, decreased anti-AChE activity, whereas the inhibitory activity against the enzyme was increased by the incorporation of a chloro group, as in compounds 13c, 13d, 13i and 13j. That is worth mentioning; the connection of the imine moiety to the coumarin core at 4-O was preferred. The structure-activity relationship (SAR) study of hybrids (13a-j) was summarized as follows ( Figure 2): (1) The anti-AChE activity of the target compounds, which substituent group on imine ring-B was the electron-withdrawing group was better than that of the electron-donating group (-Cl > H > OCH 3 ). (2) The bioassay outcome pointed out that the presence of chloro on C-4 of ring-B together with methoxy substituents on phenyl ring-A at the third position of the imine scaffold enhanced the AChE inhibition.
(3) The coumarin core attached to the imine moiety via the O atom at C-4 was more active than those attached at C-3.

Molecular Docking
The proposed molecular docking scores and binding mode with representation keys for the type of interaction docking in the active site of the acetylcholinesterase (4EY7) protein of newly synthesized hybrids (13a-j) studied in this work are presented in Table 2 and Figure S1. Upon careful inspection of these results, derivatives (13a-j) displayed binding affinity values ranging from −10.60 to −13.2 kcal/mol, which was the highest as compared to the positive compound, galantamine (−9.60 kcal/mol). As shown in Figure S1, the coumarin moiety, located in the catalytic anionic site (CAS) of the receptor, established hydrogen bonds (HBs) between carbonyl oxygen and amino acid residues. Compounds 13a, 13c, and 13d formed HBs with GLY121 and SER203, but in cases 13c and 13d (Figure 3), another HBs was observed with TYR337, while hybrid 13a formed HBs with amino acid residues of GLY122 and PHE295 with affinity energies of −11.9, −13.2, and −13.2 kcal/mol, respectively. Hybrids 13b and 13e formed another HBs with GLY122, but 13b formed one HBs with SER203 while 13e formed with TYR13, and their affinity energies values were −11.7 and −13.1 kcal/mol, respectively.  Compounds 13f, 13h, 13i, and 13j formed HBs with the TYR132, SER203, TYR132, and TYR337 amino acid residues, respectively, and the corresponding binding energies were −12.7, −10.6, −11.7, and −12.3 kcal/mol, respectively, with more than one HBs for compound 13j with TYR124 amino acid residues ( Table 2). In the case of hybrids 13c and 13d, the aryl ring of the coumarin core formed two π-π T-shaped interactions with residues of TRP86. Additionally, another aryl ring was established through two π-π stacking interactions with TRP286 amino acid residues. The CH 3 group of coumarin interacted with residues of TYR74, TYR124, TRP86, and TRP286 through π-alkyl interactions. Careful inspection of the binding site pattern and binding energy between the AChE enzyme and the studied hybrids (13a-j) indicated that both hybrids (13c and 13d) could be good candidates for the treatment of AChE, which is correlated with experimental results. We will be focusing on both (13c and 13d) in our next discussion of stability complexes by molecular dynamics simulation. According to the findings in Table 3, the size (volume) and molecular weight of proposed hybrids are often proportional to polarizability data, molecular refractivity, and surface area grid, such as hybrid 13f, which has the highest volume value (1263. 19

Structural Activity Relationships (SARs)
In this study, the physicochemical properties such as molar volume (V), hydration energy (HE), molar refractivity (MR), surface area grid (SAG), and polarizability (Pol) for hybrids (13a-j) were calculated (Table 3) and discussed using HyperChem (v8.0.7). The molecular polarizability (Pol) characteristics of a compound are determined based on how efficiently its electronic system controls itself in response to the presence of an external electric field of light. The importance of molecular polarizability is that it plays a crucial role in simulating a variety of compound characteristics and bioactivities [25]. Molecule volume, which controls things like blood-brain barrier permeability and intestinal absorption, is the main factor that influences molecular polarizability. Thus, molecular volume must be used in QSAR investigations to simulate molecular characteristics and bioactivities. A further SAR parameter is molar refractivity (MR), a steric characteristic that is dependent on the spatial arrangement of the phenyl ring in the compounds under evaluation. The spatial arrangement is significant because it is crucial to understanding how drug molecules interact with biological receptors. The London dispersive force, which is greatly involved in the interaction between drug molecules and receptors, is another factor that influences molar refractivity in addition to its dependence on molecular volume. Table 3. The physicochemical properties analysis and QSAR properties of the synthesized compounds (13a-j) against AChE.
The obtained results in Table 3 exhibit an increase in the values of hydrophobicity, causing a decline in hydration energy. The hydration energy determines the various molecular conformations' stability in aqueous solutions [26,27]. The change in the hydration energy value is affected by the increase or decrease in the number of hydrogen bonds (acceptors and donors). Table 3 illustrates the absolute values of hydration energy ordered as 13f < 13h < 13d < 13i < 13b < 13e < 13g < 13c < 13j < 13a with values of (−12. Lipophilicity is a major determinant of many ADME properties. Log P expresses the portioning of the drug molecules between the aqueous medium outside the cell membrane and the lipid nature of the cell membrane. This means that compounds with a lower Log P are more polar and have poorer lipid bilayer permeability, whereas hybrids with a higher Log P are more nonpolar and poorly soluble in water [28]. For that reason, all compounds except compounds 13c and 13j have good aqueous solubility. Furthermore, Log P values of compounds 13f = 13h < 13d = 13i < 13b = 13e = 13g < 13a are in the field of optimal values (0 < Log P < 5) [29]. It can be concluded that these hybrids have optimal biological activity and good oral bioavailability. Hybrids 13c and 13j need a drug delivery carrier to deposit them on the surface of a suitable nanomaterial with specific properties to enhance oral bioavailability.

Molecular Dynamics Simulation and System Stability
The conformational stability of the complex interaction is influenced by the molecular interactions and the solvent conditions around the receptor. As the initial structure, the best-docked pose of the most active compounds (13c and 13d) with the highest binding affinity was chosen. Moreover, MD was performed to explore the interaction modes and stability of these compounds [30,31]. For this reason, a long-range MD simulation of 100 ns was investigated to study structural stability and conformational stability as well as the dynamics of protein-ligand complexes. Herein, the root-mean-square deviation (RMSD) throughout the 100 ns simulations was used to determine the stability of the systems. An RMSD value lower than 3.0 Å was the most acceptable; it shows that the system is the most stable [20]. For all frames of the AChE protein, ligands (13c or 13d), and ligand-protein complex systems, as presented in Figure 4, the average RMSD values were 1.530, 1.691, and 2.364 Å in the 13d complex, whereas in the 13c complex, they were 1.631, 2.229 and 2.396 Å, respectively. The standard deviation of the average RMSD values were 0.199, 0.296, and 0.218 Å, in the 13d complex, whereas in the 13c complex, they were 0.208, 0.328, and 0.262 Å, respectively. These results revealed that, compared to the 13c-AChE complex, the other complex (13d-AChE) established a significantly more stable conformation. During MD simulation production, examining amino acid residue behavior and its interaction with the compound necessitates assessing receptor structural flexibility upon compound binding [20]. average RMSD values were 0.199, 0.296, and 0.218 Å, in the 13d complex, whereas in the 13c complex, they were 0.208, 0.328, and 0.262 Å, respectively. These results revealed that, compared to the 13c-AChE complex, the other complex (13d-AChE) established a significantly more stable conformation. During MD simulation production, examining amino acid residue behavior and its interaction with the compound necessitates assessing receptor structural flexibility upon compound binding [20]. The residue changes were evaluated through the root-mean-square fluctuation (RMSF) method to study the effect of ligand binding to the relevant targets throughout the simulations (100 ns). The calculated average RMSF value for the 13c-AChE complex was 0.85 Å, whereas its value for the 13d-AChE complex was 0.83 Å to protein systems. Figure 5 illustrates the overall amino acid residue fluctuations of both complex systems. These findings indicate that the inhibition of the 13d-AChE complex system is lower than that of the complex (13c-AChE) system, which will reflect well on the complex stability. Figure 6 depicts the number of hydrogen bond interactions between hybrids (13c and 13d) and the target protein (AChE) with an angle cut of 10 degrees and an rcut of 3.0 Å, plotted against simulation time (100 ns). The average number of hydrogen bonds per timeframe was calculated to be 1.09 for 13c-AChE and 1.005 for 13d-AChE. It could be observed that the interactions dramatically increased the number of hydrogen bonds per trajectory analysis from 1 to 5 HBs. Therefore, the obtained results show that the system 13d complex acquired a relatively more stable conformation than the other system, the 13c complex. The radius of gyration (Rg) provides evidence for both simulation stability and protein structural compactness. Rg values for the studied complexes were 23.077 and 23.198 Å, The residue changes were evaluated through the root-mean-square fluctuation (RMSF) method to study the effect of ligand binding to the relevant targets throughout the simulations (100 ns). The calculated average RMSF value for the 13c-AChE complex was 0.85 Å, whereas its value for the 13d-AChE complex was 0.83 Å to protein systems. Figure 5 illustrates the overall amino acid residue fluctuations of both complex systems. These findings indicate that the inhibition of the 13d-AChE complex system is lower than that of the complex (13c-AChE) system, which will reflect well on the complex stability. Figure 6 depicts the number of hydrogen bond interactions between hybrids (13c and 13d) and the target protein (AChE) with an angle cut of 10 degrees and an rcut of 3.0 Å, plotted against simulation time (100 ns). The average number of hydrogen bonds per timeframe was calculated to be 1.09 for 13c-AChE and 1.005 for 13d-AChE. It could be observed that the interactions dramatically increased the number of hydrogen bonds per trajectory analysis from 1 to 5 HBs. Therefore, the obtained results show that the system 13d complex acquired a relatively more stable conformation than the other system, the 13c complex. The radius of gyration (Rg) provides evidence for both simulation stability and protein structural compactness. R g values for the studied complexes were 23.077 and 23.198 Å, respectively, as shown in Figure 7. R g of 13c-AChE showed a more rigid structure than 13d-AChE.   The MD simulation in ligand-bound conditions was utilized to determine the solventaccessible surface area (SASA) of the receptor. As shown in Figure 8, when the ligand bound to the target receptor, the SASA values changed. The analysis shows that the 13d-protein complex is more stable upon ligand binding than the 13c-protein complex in regards to protein folding states and stability. To estimate the binding between the ligands 13c or 13d and AChE complexes, the contactFreq.tcl module on VMD and a cutoff of 4 Å were used to perform a contact frequency (CF) study, as shown in Figure 9a,b. In the simulation 13c complex case study, the following amino acid residues exhibited higher CF values: TYR72, THR75, LEU76, TYR124, TRP286, and TYR341. Whereas in the 13d complex, there were good contact surfaces with the protein pocket. Through this study, the following residues displayed higher CF values: TRP86, GLY121, TYR124, TRP286, LEU289, PRO290, GLN291, GLU292, ARG296, PHE297, TYR337, PHE338, and TYR341.    The MD simulation in ligand-bound conditions was utilized to determine the solvent-accessible surface area (SASA) of the receptor. As shown in Figure 8, when the ligand bound to the target receptor, the SASA values changed. The analysis shows that the 13d-protein complex is more stable upon ligand binding than the 13c-protein complex in regards to protein folding states and stability. To estimate the binding between the ligands 13c or 13d and AChE complexes, the contactFreq.tcl module on VMD and a cutoff of 4 Å were used to perform a contact frequency (CF) study, as shown in Figures 9a,b. In the simulation 13c complex case study, the following amino acid residues exhibited higher CF values: TYR72, THR75, LEU76, TYR124, TRP286, and TYR341. Whereas in the 13d of 4 Å were used to perform a contact frequency (CF) study, as shown in Figures 9a,b. In the simulation 13c complex case study, the following amino acid residues exhibited higher CF values: TYR72, THR75, LEU76, TYR124, TRP286, and TYR341. Whereas in the 13d complex, there were good contact surfaces with the protein pocket. Through this study, the following residues displayed higher CF values: TRP86, GLY121, TYR124, TRP286, LEU289, PRO290, GLN291, GLU292, ARG296, PHE297, TYR337, PHE338, and TYR341.

Binding Free Energy by MM/GBSA Methods
The molecular mechanics energy approach (MM/GBSA), which combines surface area continuum solvation and generalized Born, is a commonly used approach for calculating the free binding energies of small molecules to biological macromolecules. This approach might be more reliable compared to docking scores [32]. Therefore, to validate the docking scores that were anticipated by molecular docking studies for hybrids 13c and 13d against the AChE receptor, the binding free energy of the simulated ligand-protein was determined. The binding free energy of the simulated complex was computed to revalidate the inhibitor affinity predicted by docking simulation studies for the hybrids 13c and 13d with the AChE receptor. The MolAICal tool was used to take snapshots of the system trajectories being investigated in order to determine the binding free energy [33]. All the reported computed energy solvation components are presented in Table 4. Most negative values demonstrate favorable interactions. The binding free energy of both complexes 13c and 13d was calculated to be −23.645 and −36.042 kcal/mol, respectively. A close inspection of the individual energy contributions displays that the van der Waal energy of both complexes was found to be 13c-AChE with −44.913 kcal/mol, which had less binding affinity. Whereas 13d-AChE (−51.081 kcal/mol) exhibited strong binding affinity. The electrostatic energy for both complexes has considerable and moderate values. The binding free energy shows that the 13d-AChE complex was found to have more stability than the 13d-AChE complex.    The optimized geometrical parameters (dihedral angles, bond angles, and bond lengths), natural charges, natural population of the nucleus of proposed derivatives, the energetics of the ground state, molecular electrostatic potential maps, and reactivity descriptors of synthesized hybrids were calculated and analyzed utilizing the spectroscopic data and elemental analysis.
Ground State Geometric Parameters (S281-S299) The optimized geometry, dihedral angles, bond angles, bond lengths, vector of the dipole moment, and numbering system of coumarin hybrids (13a-j) are presented in Table  S1 and Figure S2. In the present work, the gas-phase wB97XD/6-311++G(d,p) has been compared to the available crystal data of 7-acetoxy-coumarin (ref. CCDC 1113091) [34] to evaluate the geometrical parameters. The mean absolute errors (MAE) that were calculated for the bond lengths and angles of the coumarin nucleus are given in Table S2. Careful inspection shows that MAEs range from 0.001 to 0.051 Å in bond length, from 0.014 to 0.562 degrees in bond angles, and from 0.004 to 1.354 degrees in dihedral angles in long-range corrected hybrid functionals (wB97XD), which give complete reducibility in predicting bond lengths and angles with experimental results with respect to experimental results and computational time and power uses. Therefore, the wB97XD/6-311++G (d,p) level of theory was selected for geometry optimizations and all calculations. The majority of the calculated bond lengths show underestimated values with percentages ranging from 0.25 to 1.1% in O1-O11 and an overestimation with values ranging from 1.1 to 5.1% in C8-C17. Generally, there is no major change. Inspection of the values of the dihedral angles compiled in Table S1 shows that almost all molecules are planar except the N36-phenyl moiety, which is out of plane and perpendicular in all the studied compounds (13a-j) with dihedral angle values of −41.8 to −47.6 degrees. The bond angle values calculated are shown in Table S1. Results vary from 109.7 to 125.0 degrees, which nicely compares to a regular SP 2 hybridization geometry.

Natural Charges and Natural Population
Natural charge analysis performed on the electronic structures of synthesized compounds (13a-j) clearly describes the distribution of electrons in various subshells of their atomic orbitals. The charge analysis carried out for all compounds using wB97XD/6-311++G(d,p) level of calculation is presented in Table 5. In Table 5, the most electronegative charges for hybrids (13a-j) are accumulated on O11, O24, O16, O1, N36, and C7, respectively, from −0.563 e to −0.313 e. According to an electrostatic point of view, these electronegative atoms tend to have active sites for electrons. However, the most electropositive atoms, such as C2, C10, C8, and C25, from +0.774e to +0.271e, respectively, tend to accept electrons at their active sites. Going from 13a to 13j is a minor change in natural charge with the same pattern of electrostatic mapping with the order. An in-depth study of the natural charge pattern of hybrids is extremely beneficial for gaining insight into the crucial interactions between title hybrids and biological receptors of AChE, which enhances the investigation of cytotoxicity activity.

Frontier Molecular Orbitals (FMOs) Analysis
As presented in Table 6, among all synthesized hybrids, 13e displayed higher stability and less reactivity with an energy gap value of 7.96 eVm, whereas hybrid 13b (∆E gap = 7.38 eV) showed the lowest stability and highest reactivity [35][36][37][38]. The energy gaps of the rest of the hybrids were ordered as follows: 13j < 13i < 13d < 13c < 13a < 13f < 13h < 13g. Due to the significance of the parameters such as I (potential ionization) and A (electron affinity), their calculations enable us to determine the global reactivity descriptors. The I and A parameters are related to the one-electron orbital energies of the HOMO and LUMO. The obtained results (Table 6) exhibited that hybrids 13j had the highest values of I (8.29 eV) and A (0.37 eV). The I value of the other hybrids are 13d > 13c > 13e > 13i > 13a > 13g > 13f > 13b > 13h. Whereas in the case of A values, the order is as follows: 13d > 13c > 13a > 13b > 13i > 13g > 13e > 13f > 13h. Then we can predict that derivatives 13j, 13d, and 13c are the best candidates for interaction with other biological AChE receptors. All of the hybrids have almost similar HOMO and LUMO isodensity dispersion as depicted in Figure S3, except for 13c and 13d (Figure 10), which have slightly different dispersion on the coumarin rings and the (E)-N-benzylideneaniline moiety. The direction of the electronic charge transfer motion is represented by the dipole moment vector with the order norm vector; synthesized compounds are ordered as 13j > 13i > 13f > 13e > 13a > 13d > 13c > 13h > 13b > 13g. Values are mean ± SD duplicate assays.

Global Reactivity Descriptors
The density functional theory (DFT) uses the chemical system's electron density to explain several basic ideas about how chemicals react [39]. In chemistry, understanding the nature of chemical interactions and predicting the chemical reactivity of molecules, atoms, or ions are the two most challenging problems. Herein, we studied the reactivity of novel synthesized derivatives (13a-j). Table 7 shows the values of the important reactivity descriptors that help us figure out how reactive and stable hybrids (13a-j) are. Among all compounds, coumarin 13j, with the highest value of η = 3.98 eV, is the chemically hardest compound, while coumarin 13b, with the lowest value (η = 3.69 eV), is chemically soft and more reactive. The chemical hardness of other hybrids is ordered as 13j > 13i > 13d > 13c > 13a > 13f > 13h > 13g. A general idea of charge transfer in any molecule's ground state can be obtained from the electronic chemical potential (V) value. In terms of chemical potential, hybrid 13h has the greatest value (−3.8 eV), whereas hybrid 13j has the lowest (−4.33 eV), and other coumarin hybrids are in the following order: 13f > 13b > 13g > 13e > 13i > 13a > 13c > 13d (Table 7).

Global Reactivity Descriptors
The density functional theory (DFT) uses the chemical system's electron density to explain several basic ideas about how chemicals react [39]. In chemistry, understanding the nature of chemical interactions and predicting the chemical reactivity of molecules, atoms, or ions are the two most challenging problems. Herein, we studied the reactivity of novel synthesized derivatives (13a-j). Table 7 shows the values of the important reactivity descriptors that help us figure out how reactive and stable hybrids (13a-j) are. Among all compounds, coumarin 13j, with the highest value of η = 3.98 eV, is the chemically hardest compound, while coumarin 13b, with the lowest value (η = 3.69 eV), is chemically soft and more reactive. The chemical hardness of other hybrids is ordered as 13j > 13i > 13d > 13c > 13a > 13f > 13h > 13g. A general idea of charge transfer in any molecule's ground state can be obtained from the electronic chemical potential (V) value. In terms of chemical potential, hybrid 13h has the greatest value (−3.8 eV), whereas hybrid 13j has the lowest (−4.33 eV), and other coumarin hybrids are in the following order: 13f > 13b > 13g > 13e > 13i > 13a > 13c > 13d (Table 7).   Values are mean ± SD triplicate assays.
A thermodynamic parameter that is represented by the electrophilicity index (ω) estimates the energy changes that occur when a chemical system reaches saturation with the addition of more electrons. This is very beneficial in determining a system's chemical reactivity. As shown in Table 7, coumarin hybrid 13h is nucleophilic in nature with the lowest electrophilicity index value equal to 1.9 eV, while hybrid 13j is strongly electrophilic in nature (ω = 2.37 eV). The lowest electrophilicity index (ω) order of other coumarins is the same as electron affinity (A). The electronegativity (X) describes the tendency of an atom in a covalent bond to draw electrons towards it. From the obtained electronegativity of the synthesized hybrids, coumarin 13j is the best electron acceptor with X = 4.33 eV, which had the highest electronegativity value as compared to the other hybrids. The electronegativity (X) order of other coumarins shows the same behavior in both the electron affinity and electrophilicity indexes. In terms of global softness, hybrid 13b displayed the highest reactivity and softness values (0.14 eV), whereas the other hybrids exhibited almost the same values of softness (0.126 to 0.133 eV) with order 13g > 13h > 13f > 13a > 13c > 13d > 13i > 13j > 13e.

Local Reactivity Descriptor
To study the site selectivity and chemical reactivity of a molecule, the principles of local reactivity descriptors have been frequently applied [40,41]. The Fukui function is a local descriptor that can be used to study molecular site selectivity [42]. It is the first derivative of the electronic density ρ(r) in relation to the electron numbers (N) of a system at a constant external potential ν(r) [43], as represented in the following equation: Based on the changes in electrical density throughout a reaction process, we can calculate Fukui functions to identify the active sites. As shown in the following equation, for the three different environments of chemicals, the Fukui functions f + (r), f -(r), and f 0 (r) are determined [44][45][46]: where q k (N), q k (N + 1) and q k (N − 1) are the atomic populations on the k th atom for the neutral molecule, anionic and cationic species, respectively. Tables S3 and S4 represent the descriptor values of coumarin hybrids (13a-j) computed at the wb97xd/6-311++G (d, p) level. In addition to knowing how an atomic site in a molecule could be electrophilic or nucleophilic, Labbe et al. [47] suggested an additional dual descriptor (∆ f (r)) that is provided by the following equation: The obtained results indicate that the most electrophilic reactivity is on the imine moiety, which is mostly found on the atoms O24, N36, C37, C40, and C42, while the nucleophilic active site on the coumarin moiety is located on O1, C2, C3, C4, C6, C8, O11, and O16 localized. From Tables S3 and S4, when considering the dual descriptor ∆ f (r) for the nucleophilic and electrophilic attacks, as well as the philicity indices, the same result could be obtained. The high electronegativity of atoms N and O led to an electron density redistribution, in addition to the effect of -OCH 3 insertion groups in R1 and R2 substitutions, that causes these characteristic changes. These findings are in agreement with the analysis of the natural population using calculated HOMO and LUMO. In 2004, Chattaraj et al. proposed the generalized philicity concept; with the help of corresponding condensed-to-atom Fukui function variations, they developed a local quantity known as philicity coupled with a site k in a molecule ( f α k ), as given in the following equation [48].
where α represents the local philic quantities, radical (α = 0), nucleophilic (α = +), and electrophilic (α = −) attacks. According to the aforementioned equation, the most electrophilic property has the highest value of ω α k . Moreover, different local softness was proposed by Lee et al. to define the reactivity of a molecule [49] as in the following Equation.
In the aforementioned equation, α is represented by local softness quantities for nucleophilic (α = +), electrophilic (α = −), and radical attacks (α = 0). In order to complete the picture, the software package Multiwfn (v. 3.7) determined the local electrophilicity and nucleophilicity index, condensed local softness, and relative electrophilicity/nucleophilicity for each atom in the compounds from a CDFT point of view [50]. A close inspection would show that all the compounds had the donating and the back-donating processes at the center of their active sites, in agreement with the Fukui functions and also with the frontier orbital, as shown in the results obtained, which are boldly given in Tables S5 and S6. These findings showed the studied compounds to have several active sites, making them able to interact with the surface of pocket proteins via donating electrons. Lastly, the aforementioned local descriptors show that the experimental data in this study are in agreement with the theoretical variation of the compound's efficiency.

Molecular Electrostatic Potential (MEP)
In several fields of chemistry, the electrostatic potential (ESP) on molecular surfaces has become one of the most effective tools for identifying, analyzing, and understanding trends [51,52]. It is related to electronic density, which is an excellent descriptor for describing the charge distributions on a molecule, identifying regions that are differently charged, and identifying the sites where hydrogen-bonding interactions, electrophilic properties, and nucleophilic properties are most likely to take place [53]. Electrostatic potentials (ESPs) are essential for predicting and understanding intermolecular interactions [52]. The crucial interactions between the synthesized hybrids and biological targets must be better understood with the help of an in-depth analysis of their ESPs. The ESP, as denoted by V(r) (in a.u.) at a given point r (x,y,z) in the molecule's vicinity, is a calculation of the electrostatic energy that a positive unit test charge would experience at that point. Negative and positive ESPs corresponded to attractive and repulsion interactions, respectively. The following equation defines the ESP as the interaction energy between a proton at r and the electrical charge produced by the electrons and nuclei.
where Z A is the charge, R A is the position of nucleus A, and ρ (r') is the electron density at position r'. The electrostatic potentials-mapped surfaces of the hybrids (13a-j) are presented in Figure S4. The overall van der Waals surface can be divided into several fragments by the quantitative molecular surface analysis module of the Multiwfn package, and this capability enables us to analyze the features of the ESP distribution. For hybrid 13a, the surface displays a high negative value of ESP at the O11, O16, O24, and N36 (−43.97, −16.9, −11.0, and −28.01 kcal/mol), respectively, with the positive charge being distributed among various active sites. As for the potential of various derivatives (13b-j) to redistribute electrons, the global minima of ESPs around O11, O16, O24, and N36 centers increase, reaching minimum values in case 13f with sequence (−46.9, −18.5, −22.0, and −27.9 kcal/mol) due to the electron-donating OCH 3 group but in case 13c and 13d (Figure 11), Cl derivatives show that the negative decrease due to electron drawing behavior as sequence (−43.3, −16.3, −9.8 and −24.5 kcal/mol) and (−45.9, −18.5, −6.5 and −23.5 kcal/mol), respectively. The global maxima of ESPs on the derivatives (13b-j) surfaces are located on the carbon with the proton of these derivatives, which vary from +36.2, +35.6, +37.4, +29.4, +28.3, +27.6, +34.1, +26.8, +28.7 and +35.9 kcal/mol for hybrids (13a-j), respectively. This indicates that electrostatic or hydrogen bonding will be the main interaction between +28.3, +27.6, +34.1, +26.8, +28.7 and +35.9 kcal/mol for hybrids (13a-j), respectively. This indicates that electrostatic or hydrogen bonding will be the main interaction between hybrids and their target receptors. The ability to generate hydrogen bond interaction and intramolecular charge transfer (ICT) is confirmed by a careful examination of these ESP values, indicating that they can act as therapeutics. These values show the same findings from the analyses of the NBO population and local reactivity descriptors reported in the previous section.

General Procedures
Reagents and solvents used in the research were commercially available. The reagents 3-Methoxyphenol, ethyl 4-chloroacetoacetate, 4-hydroxy-3methoxybenzaldehyde, 3-hydroxy-4-methoxybenzaldehyde, 4-chloroaniline, and 4methoxyaniline were purchased from Sigma-Aldrich. The other materials were purchased from Acros Organics, such as aniline, glacial acetic acid, anhydrous potassium carbonate, acetone, chloroform, dichloromethane, ethyl acetate, methanol, ethanol, and n-hexane. Analytical grade (AR) solvents such as ethanol absolute AR (99.9%) were purchased from Fisher Scientific. The purity of these chemicals was 90-99.9%, and they were used without further purification. Thin-layer chromatography (TLC) was used to detect compounds present in the products. The TLC plates used were the thin aluminum plates from Merck, pre-coated with silica gel F254 with 0.2 mm thickness. The spots were visualized under UV light at 254 nm or 365 nm. Melting points (uncorrected) were determined using a Barnstead Electrothermal 9100 melting point. The infrared IR spectra were recorded using a Perkin Elmer FTIR spectrometer. Samples were prepared as KBr discs. The 1 H NMR (400 MHz) and 13  3-Methoxyphenol (30 mmol) was added to a cooled aqueous solution of 70% H2SO4 (60 mL) and followed by adding ethyl 4-chloroacetoacetate (35 mmol). The solution was stirred in an ice-bath for a period of 9 h and continuously stirred for 96 h. The reaction mixture was poured into cold water. The precipitate was filtered off, washed several times with cold water, dried, and recrystallized from ethanol to afford 4-(chloromethyl)-7methoxy-2H-chromen-2-one (9) as a white powder (6.34 g, 94%); Rf = 0.51 (hexane: EtOAc

General Procedures
Reagents and solvents used in the research were commercially available. The reagents 3-Methoxyphenol, ethyl 4-chloroacetoacetate, 4-hydroxy-3-methoxybenzaldehyde, 3-hydroxy-4methoxybenzaldehyde, 4-chloroaniline, and 4-methoxyaniline were purchased from Sigma-Aldrich. The other materials were purchased from Acros Organics, such as aniline, glacial acetic acid, anhydrous potassium carbonate, acetone, chloroform, dichloromethane, ethyl acetate, methanol, ethanol, and n-hexane. Analytical grade (AR) solvents such as ethanol absolute AR (99.9%) were purchased from Fisher Scientific. The purity of these chemicals was 90-99.9%, and they were used without further purification. Thin-layer chromatography (TLC) was used to detect compounds present in the products. The TLC plates used were the thin aluminum plates from Merck, pre-coated with silica gel F 254 with 0.2 mm thickness. The spots were visualized under UV light at 254 nm or 365 nm. Melting points (uncorrected) were determined using a Barnstead Electrothermal 9100 melting point. The infrared IR spectra were recorded using a Perkin Elmer FTIR spectrometer. Samples were prepared as KBr discs. The 1 H NMR (400 MHz) and 13 (9) 3-Methoxyphenol (30 mmol) was added to a cooled aqueous solution of 70% H 2 SO 4 (60 mL) and followed by adding ethyl 4-chloroacetoacetate (35 mmol). The solution was stirred in an ice-bath for a period of 9 h and continuously stirred for 96 h. The reaction mixture was poured into cold water. The precipitate was filtered off, washed several times with cold water, dried, and recrystallized from ethanol to afford 4-(chloromethyl)-7-methoxy-2Hchromen-2-one (9) as a white powder (6.

General Procedure of Synthesis of Schiff bases (12a-j)
Glacial acetic acid (1 mL) was added dropwise into a solution of hydroxybenzaldehyde (10 mmol) and substituted aniline (10 mmol) in 95% EtOH (20 mL). The reaction mixture was refluxed for one hour at 80 • C and left to stir at room temperature overnight. After that, the mixture was poured into iced water (25 mL). The solvent was removed by vacuum filtration to form a solid. The solids were purified by crystallization from ethanol to afford the desired Schiff bases (12a-j).

Biological Evaluation of Hybrids against AChE
The procedure described by Ellman [64], slightly modified, was performed to assess the AChE inhibitory activity of synthesized compounds, according to reported procedures in the literature [20,[65][66][67]. Detailed protocol is provided in supplementary materials.

Molecular Docking Study
The 3D structures of all synthesized hybrids were previously geometry optimized and then saved in PDBQT file format utilizing Autodock Tools (v1.5.6rc3) [68]. This involves the following steps: developing a torsion tree by detecting and choosing a root and saving it as a pdpqt file for mapping prior to molecular docking simulation. The crystallographic structure of AChE (PDB ID: 4EY7) was downloaded from the protein data bank website (https://www.rcsb.org/structure/4EY7, accessed on 15 June 2023). The protein structure was cleaned from heteroatoms, and the Swiss-PdbViewer was used for molecular mechanics energy minimization [69]. The protein structure was prepared for molecular docking study according to our previous studies and saved in PDBQT format [28,70]. Finally, molecular docking simulation was performed, visualized, and analyzed as previously reported [71][72][73].

Molecular Dynamics Simulation (MDS)
MD simulations were investigated to calculate binding energies and study the relative stabilities of the interactions between the two top hybrids (13c and 13d), consensus docking scores, and target receptors. The simulation studies were carried out utilizing the NAMD (v2.13) suit [74] and the CHARMM36 [75] force field. Furthermore, the latest CHARMM/CGenFF force field was used to generate the parameters and topological files for the selected compounds [19,76]. The complex structure (ligand-protein) was placed in the middle of a box that solvated molecules of water with a TIP3P explicit solvation model. A 0.15 molar solution of (145 Na + and 135 Cl -) ions, to mimic the physiological salt concentrations, was added to provide electrostatic screening and charge neutralization which extended 20 Å from the protein. CHARMM and the periodic boundary conditions were set with dimensions of a rectangle-cubic system of 117.0, 117.0, and 117.0 Å in x, y, and z directions, respectively.
The MD studies include minimization, equilibration, and data analysis. There were no atoms restricted in the MD simulations. The isothermal-isobaric (NPT) ensemble and a 2 fs time integration step were selected in this study, and the pressure was set at 1 atm utilizing the Nose'-Hoover Langevin piston barostat [77]. The Langevin thermostat has been used to set the temperature at 300.0 K [78]. The force-field parameters were assumed in order to minimize and equilibrate the complexes in the system, which have a scaling of 1.0 Å. The preliminary energy of the complex was minimized through 2000 steps at 300 K. The temperature, kinetic energy, and/or pressure of the system were controlled by Langevin dynamics simulation through another 144,000 steps. The 500,000 minimization steps were used to equilibrate the solvated system, and 50,000,000 runs for 100 ns. The VMD package was used for the analysis of the output data [79]. A 20.0 Å for the distance cut-off was used for Lennard Jones interactions, and short-range non-bonded interactions with a pair list distance of 12 Å were smoothly truncated at 8.0 Å. The long-range electrostatic interactions were analyzed and visualized using the particle-mesh Ewald (PME) procedure [80], where simulated cells were placed in a grid box with 1.0 Å.

Binding Energy Calculations
The relative binding energy calculations were carried out by the one-average molecular mechanics generalized Born surface area (MM/GBSA) method [81]. The MolAICal [33] was employed to calculate the MM/GBSA for complex (ligands and receptors) based on molecular dynamical (MD) simulated results by NAMD, in which the ligand (L) binds to the target receptor (R) to generate the receptor-ligand (RL). We were only interested in relative binding energies based on MM/GBSA calculations, which were the Gibbs relative binding energy, given by:

Geometry DFT Optimization
During this study, packages of programs were used to run the molecular modeling calculations of all synthesized hybrids utilizing the Gaussian 09W software package [82]. The molecular structures of hybrids were geometrically optimized using density functional theory (DFT) with long-range corrections functional wB97XD (DFT/wB97XD) [83], which includes empirical dispersion with basis set 6-311++G (d,p) [84]. During the geometry optimization, no symmetry restrictions were used [85,86]. The choice of long-range corrections functional wB97XD with a large basis set was due to the accuracy, consistency, flexibility, and better performance of Grimme's D2 dispersion model, which includes empirical dispersion [83,87]. The same level of theory has been applied to calculate the vibrational frequency for each compound, and the molecular structures of each hybrid were found to correspond to real minima of the potential energy surface [88]. In order to identify the reactive site of the molecules, the DFT/wB97XD was employed to describe reactivity descriptors and molecular stability. A descriptor of local reactivity was computed using the Fukui function and the dual descriptor [40,42,47].
Furthermore, the quantum chemical descriptors from conceptual density functional theory (CDFT) were calculated by utilizing the Multiwfn (v3.7) package [41]. The electrostatic potential (ESP) of the molecules was rendered by the Visual Molecular Dynamics package (VMD 1.9 program) based on the data outputted by the Multiwfn program [50,79]. Natural bond orbital (NBO) analyses have been calculated utilizing NBO 3.1, which is provided in the Gaussian 09W program. The GaussView (v6.1) [89] and ChemCraft (v1.6) packages [90] were used to visualize the optimized structure and molecular orbitals. The QSAR features included in the HyperChem program (v8.0.7) [91] were used to determine the SAR properties of all synthesized compounds.

Conclusions
A series of Schiff base-coumarin hybrids (13a-j) were designed with the assistance of a molecular docking study. The affinity energy of hybrids (13a-j) on the active site of acetylcholinesterase (AChE) was in the range of −10.6 to −13.2 kcal/mol, while the affinity energy of positive controls, donepezil and galantamine, was recorded at −11.4 and −9.6 kcal/mol, respectively. The obtained results showed that most of the designed hybrids were more active than the positive controls. Hybrids were synthesized and bioevaluated against AChE. Results showed that most of them could potentially inhibit the target enzyme. Compounds 13c and 13d with IC 50 values of 0.232 ± 0.011 and 0.190 ± 0.004 yielded the lowest IC 50 values and were 5-fold stronger than the positive control. The detailed insights into the electronic structure properties, its link to drug-like properties, and the structure-activity relationships of the designed compounds were provided by applying the high level of computational approaches of DFT/wB97XD methods by using 6-311++G(d,p) as the basis set wave function. Moreover, the molecular docking studies revealed that the synthesized compounds exhibit molecular interactions via hydrogen bonds with GLY121, GLY122, TYR132, SER203, PHE295, and TYR337 amino acid residues of the target enzyme. Careful inspection of the pattern of binding sites and binding energy indicated 13c and 13d compounds could be good candidates for AChE inhibitors, which is correlated with experimental results. The drug likeness and QSAR descriptors show that compound 13f has the highest value, while compound 13a has the lowest values in all descriptors for other derivatives. All target compounds except 13c and 13j have good aqueous solubility. Furthermore, Log P values of compounds 13f = 13h < 13d = 13i < 13b = 13e = 13g < 13a are in the field of optimal values (0 < Log P < 5). Based on these, it can be stated that these compounds have good oral bioavailability and optimal biological activity. The stability of the systems was measured using RMSD during the 100 ns simulations. Obtained results revealed that the 13d complex system acquired a relatively more stable conformation than the 13c complex. The calculated average RMSF values for the 13c-AChE and 13d-AChE complexes to protein systems were 0.85 and 0.83 Å, respectively, which demonstrate that the system inhibition of the 13d-AChE protein complex system is lower than the other systems inhibition, which will reflect well on the complex stability. The binding free energy technique (MM/GBSA) of the simulated complex was computed, and it was found to be −23.645 and −36.042 kcal/mol for complexes 13c and 13d, respectively. Thus, it is suggested that the stability of the 13d-4EY7 complex over the other is better. The overall study indicates that 13d has comparable and even better descriptors than Galantamine and can be a potential AChE candidate.
Supplementary Materials: The following supporting information can be downloaded at: The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/ ph16070971/s1, Detailed description of biological evaluation protocol and IC 50 calculation, Inhibition diagrams for the compounds, FTIR spectra of the compounds, 1H and 13C spectra of the compounds, HRMS spectra of the compounds, Figure S1: The proposed binding mode of synthesized hybrids (13a-j) docked in the active site of acetylcholinesterase (4EY7) protein; (2D and 3D ligand-receptor interactions), Figure S2: The optimized geometry, numbering system, vector of dipole moment of the synthesized compounds, Figure S3: Frontier molecular orbitals of the synthesized compounds, Figure  S4: MEP surfaces of the synthesized compounds, Table S1: The selected bond length (Ao), bond angles and dihedral angles, (degree) of the compounds, Table S2: Mean absolute errors computed for selected bond lengths (Ao) and angles (degree) of synthesized compounds, Tables S3 and S4: Values of the Fukui functions and dual descriptor of compounds, Table S5: Values of the condensed local softnesses (Hartree*e) of compounds, Table S6: Values of the condensed local electrophilicity (ElP)/nucleophilicity (NuP) index (e*eV) of compounds.