Many marvelous biological functions in proteins and DNA and their profound dynamic mechanisms, such as switch between active and inactive states [58, 59] , cooperative effects [60] , allosteric transition [61, 62] , intercalation of drugs into DNA [63] , and assembly of microtubules [64] , can be revealed by studying their internal motions [65] . Likewise, to really in-depth understand the action mechanism of receptor-ligand binding, we should consider not only the static structures concerned but also the dynamical information obtained by simulating their internal motions or dynamic process. To realize this, the MD simulation is one of the feasible tools. Initial coordinates for the protein atoms were taken from the wild type (WT) and mutant forms of both L-LTA4, L-RAR complex structures. Mutations were introduced at E271Q, R563A, and K565A of the enzyme based on the previous experimental reports to investigate the single active site that catalyzes two different functions upon diverse substrates [14, 45] . The protonation states of all ionizable residues were set to their normal states at pH 7. Eight MD simulations were performed for systems including WT and mutant forms of L-LTA4 and L-RAR complexes ( Table 1 ). All MD simulations were performed with GROMOS96 forcefield using GROMACS 4.5.3 package running on a high performance linux cluster computer [66, 67] . During the MD simulations, all the protein atoms including divalent metal ion (Zn 2+ ) were surrounded by a cubic water box of SPC3 water molecules that extended 10 Å from the protein and periodic boundary conditions were applied in all directions. The systems were neutralized with Na + and Cl 2 counter ions replacing the water molecules and energy minimization was performed using steepest descent algorithm for 10,000 steps. A 100 ps position restrained MD simulations were performed for every system followed by 5 ns production MD simulations with a time step of 2 fs at constant pressure (1 atm), temperature (300 K). The electrostatic interactions were calculated by the PME algorithm and all bonds were constrained using LINCS algorithm. A twin range cutoff was used for long-range interactions including 0.9 nm for van der Waals and 1.4 nm for electrostatic interactions. The snapshots were collected at every 1 ps and stored for further analyses of MD simulations. The system stability and behavior of the catalytic structural components present in every system were analyzed using the tools available with GROMACS 4.0.5 and PyMol programs.