**4. Conclusions**

A detailed analysis of interactions of proteins and ligands represents an optimized method for the rational design of new compounds. Binding energy helps to identify a lead compound. MD, free energy perturbation, meta-dynamics, and other methods are consistently used for studying drug–target binding. Assessment of these methods may help to achieve the most definite target-optimized affinity with improved drug efficacy [43]. Integrating MD simulation with binding-free-energy calculation by means of g\_mmpbsa to determine interaction free energy between a ligand and protein is an efficient method for distinguishing between active and inactive molecules [44]. The comparative analysis of PF-07321332, α-ketoamide, lopinavir, and ritonavir via MD simulation provided a detailed insight into the interactions of these compounds with 3CLpro. We believe that our findings revealed the binding mechanism of PF-07321332 and α-ketoamide and further explained the inability of lopinavir and ritonavir to cause clinical improvement in severe COVID-19 patients by due to their low affinity towards 3CLpro. We believe the dynamic interaction and binding energy of PF-07321332 and α-ketoamide will be helpful in designing potent inhibitors for 3CLpro.

#### **5. Materials and Methods**

#### *5.1. System Preparation for Molecular Docking and MD Simulations*

The crystal structure of SARS-CoV-2 3CLpro (Protein Data Bank [PDB] ID: 6M03) and structures of 3CLpro in complex with 3CLpro (PDB ID 6Y2K) [45,46] were downloaded from the PDB [47]. The 2-D structure of PF-07321335 was drawn in Chemdraw and the threedimensional (3D) coordinates of lopinavir (PubChem CID: 92727) and ritonavir (PubChem CID: 392622) were retrieved from the PubChem database [48] for molecular docking. All hydrogens and missing atoms were added to the proteins, while crystal water was removed. The MOE software from Chemical Computing Group Canada (Montreal, QC, Canada) [49] was used for energy minimization of 3CLpro and both drugs to remove any steric clashes and to improve the bond lengths and bond angles. The default parameters were used for energy minimization, while the gradient was set to 0.01 rms kcal/mol.

#### *5.2. Molecular Docking of 3CLpro to Lopinavir and Ritonavir*

PF-07321335, α-ketoamide, lopinavir, and ritonavir were docked into the binding pocket of 3CLpro using the MOE software [49]. The key residues involved in ligand binding reported in the literature were selected for docking: His41, Asn142, Cys145, His164, Met165, Glu166, Gln189, and Thr190. The triangle matcher placement method was selected with the London dG scoring function, and the induced-fit method was applied for refinement of the binding poses of ligands. A total of 100 conformations were generated for each ligand, and the best pose with the lowest binding energy was selected. Complexes of 3CLpro with PF-07321335, α-ketoamide, lopinavir, and ritonavir were analyzed and saved for MD simulations.

#### *5.3. Building the Ligand Topology*

Prior to MD simulations, the ligand topology was generated on the CHARMM general force field (CGenFF) server (University of Maryland, Baltimore, MD, USA) [50]. Ligands, along with their 3D coordinates, were extracted from their respective complexes and saved in mol2 file format, and all hydrogens were added. The mol2 file was uploaded to the CGenFF server, which generated a CHARMM-compatible stream file comprising ligand information, such as atom types, charges, and bond parameters. A Python script was executed to convert the stream files into Gromacs-compatible [51,52] files. These files were used for the MD simulations.

#### *5.4. MD Simulations of Complexes*

The dataset for MD simulations contained apo-3CLpro and four 3CLpro complexes. All five MD simulations were carried out using the Gromacs (University of Groningen, Groningen, Netherlands) software for 100 ns each, and the CHARMM36 forcefield [53] was applied to the systems. Periodic boundary conditions were used, and the protein was placed inside a 10 Å cubic box. For solvation of the system, the TIP3P water model [54] was used, and an appropriate amount of counterions was added to neutralize the system. The solvated electroneutral system was subjected to energy minimization by the steepest-descent algorithm, followed by temperature and pressure equilibration steps. The production run was carried out for 100 ns in each system, and data analysis was performed by means of Gromacs built-in tools, MOE, and PyMOL (Schrödinger, LLC. DeLano Scientific, San Francisco, CA, USA) [55,56].

## *5.5. Binding-Free-Energy Calculation*

The molecular mechanics Poisson–Boltzmann surface area (MMPBSA) method was used to compute the binding free energy of all complexes. A thousand frames were extracted from each 100 ns simulation trajectory by means of the *gmx trjconv* utility of Gromacs. The extracted frames were used for free-binding-energy calculation with the help of the g\_mmpbsa tool [44]. The total binding free energy (ΔGbinding) was calculated as

$$
\Delta \mathbf{G}\_{\text{binding}} = \Delta \mathbf{G}\_{\text{complex}} - (\Delta \mathbf{G}\_{\text{protein}} + \Delta \mathbf{G}\_{\text{ligand}}),
$$

where Gcomplex is the average free energy of the complex, and Gprotein and Gligand are the average energy values of the protein and ligand.

**Author Contributions:** Conceptualization, B.A. and M.B.; methodology, B.A., Q.u.A., and M.B.; formal analysis, B.A., Q.u.A., and M.B.; data curation, B.A., Q.u.A., and M.B.; writing—original draft preparation, B.A., Q.u.A., and M.B.; writing—review and editing, B.A. and M.B.; supervision, S.C.; project administration, M.S.K. and S.C.; funding acquisition, M.S.K. and S.C. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work was supported by the National Research Foundation of Korea (NRF- 2020R1F1A 1071517, 2019M3D1A1078940, and 2019R1A6A1A11051471).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** The structure models and simulation trajectories are available upon request (sangdunchoi@ajou.ac.kr). All remaining data are contained within the manuscript.

**Conflicts of Interest:** The authors declare no conflict of interest.
