Recent Advances in Fragment Molecular Orbital-Based Molecular Dynamics (FMO-MD) Simulations

Fragment molecular orbital (FMO)-based molecular dynamics simulation (MD), hereafter referred to as "FMO-MD," is an ab initio MD method (Komeiji et al., 2003) based on FMO, a highly parallelizable ab initio molecular orbital (MO) method (Kitaura et al., 1999). Like any ab initio MD method, FMO-MD can simulate molecular phenomena involving electronic structure changes such as polarization, electron transfer, and reaction. In addition, FMO's high parallelizability enables FMO-MD to handle large molecular systems. To date, FMOMD has been successfully applied to ion-solvent interaction and chemical reactions of organic molecules. In the near future, FMO-MD will be used to handle the dynamics of proteins and nucleic acids.


Introduction
Fragment molecular orbital (FMO)-based molecular dynamics simulation (MD), hereafter referred to as "FMO-MD," is an ab initio MD method (Komeiji et al., 2003) based on FMO, a highly parallelizable ab initio molecular orbital (MO) method (Kitaura et al., 1999). Like any ab initio MD method, FMO-MD can simulate molecular phenomena involving electronic structure changes such as polarization, electron transfer, and reaction. In addition, FMO's high parallelizability enables FMO-MD to handle large molecular systems. To date, FMO-MD has been successfully applied to ion-solvent interaction and chemical reactions of organic molecules. In the near future, FMO-MD will be used to handle the dynamics of proteins and nucleic acids.
In this chapter, various aspects of FMO-MD are reviewed, including methods, applications, and future prospects. We have previously published two reviews of the method (Komeiji et al., 2009b;chapter 6 of Fedorov & Kitaura, 2009), but this chapter includes the latest developments in FMO-MD and describes the most recent applications of this method. Fig. 1. Schematics of the FMO-MD method exemplified by an ion solvation with four water molecules. The atomic nuclei are represented by black circles (the large one for the ion, medium ones for Oxygens, and small ones for Hydrogens) and the electron cloud by a grey shadow. The electronic structure is calculated by FMO to give force (F) and energy (E), which are then used to update the 3D coordinates of nuclei (r) by MD, i.e., by solving the classical equation of motion.
performed by the PEACH/ABINIT-MP software system composed of the PEACH MD program (Komeiji et al., 1997) and the ABINIT-MP 1 (F)MO program (Nakano et al., 2000). We have revised the system several times (Komeiji et al., 2004(Komeiji et al., , 2009a, but here we describe the latest system, which has not yet been published. In the latest system, the PEACH program prepares the ABINIT-MP input file containing the list of fragments and 3D atomic coordinates, executes an intermediate shell script to run ABINIT-MP, receives the resultant FMO energy and force, and updates the coordinates by the velocity-Verlet integration algorithm. This procedure is repeated for a given number of time steps.
The above implementation of FMO-MD, referred to as the PEACH/ABINIT-MP system, has both advantages and disadvantages. The most important advantage is the convenience for the software developers; both FMO and MD programmers can modify their programs independently from each other. Also, if one wants to add a new function of MD, one can first write and debug the MD program against an inexpensive classical force field simulation and then transfer the function to FMO-MD, a costly ab initio MD. Nonetheless, the PEACH/ABINIT-MP system has several practical disadvantages as well, mostly related to the use of the systemcall command to connect the two programs. For example, frequent invoking of ABINIT-MP from PEACH sometimes causes a system error that leads to an abrupt end of simulations. Furthermore, use of the systemcall command is prohibited in many supercomputing facilities. To overcome these disadvantages, we are currently

FMO
FMO, the essential constituent of FMO-MD, is an approximate ab initio MO method (Kitaura et al., 1999). FMO scales to N 1-2 , is easy to parallelize, and retains chemical accuracy during these processes. A vast number of papers have been published on the FMO methodology, but here we review mainly those closely related to FMO-MD. To be more specific, those on the FMO energy gradient, Energy Minimization (EM, or geometry optimization), and MD are preferentially selected in the reference list. Thus, those readers interested in FMO itself are referred to Fedorov & Kitaura (2007b) for comprehensive reviews of FMO. Also, one can find an extensive review of fragment methods in Gordon et al. (2011), where FMO is re-evaluated in the context of its place in the history of the general fragment methods.
First, the molecular system of interest is divided into N f fragments. Second, the initial electron density, ρ I (r), is estimated with a lower-level MO method, e.g., extended Hückel, for all the fragments. Third, self-consistent field (SCF) energy, E I , is calculated for each fragment monomer while considering the electrostatic environment. The SCF calculation is repeated until all ρ I (r)'s are mutually converged. This procedure is called the self-consistent charge (SCC) loop. At the end of the SCC loop, monomer electron density ρ I (r) and energy E I are obtained. Finally, an SCF calculation is performed once for each fragment pair to obtain dimer electron density ρ IJ (r) and energy E IJ . Total electron density ρ(r) and energy E are calculated using the following formulae: (2 ) In calculation of the dimer terms, electrostatic interactions between distant pairs are approximiated by simple Coulombic interactions (dimer-ES approximation, Nakano et al., 2002). This approximation is mandatory to reduce the computation cost from O(N 4 ) to O(N 2 ).

www.intechopen.com
Molecular Dynamics -Theoretical Developments and Applications in Nanotechnology and Energy 6 The total energy of the molecular system, U, is obtained by adding the electrostatic interaction energy between nuclei to E, namely, where r ij denotes the distance between nuclei i and j and Z i and Z i their charges, respectively.
Force (F i ) acting on atomic nucleus i can be obtained by differenciation of eq. (3) by r i as follows: Analytical formulation of eq. (4) was originally derived for the HF level by Kitaura et al. (2001) and used in several EM calculations (for example, Fedorov et al., 2007a) and in the first FMO-MD simulation (Komeiji et al., 2003). Later on, the HF gradient was made fully analytic by Nagata et al. (2009Nagata et al. ( , 2010Nagata et al. ( , 2011a.

FMOn
The procedure described in the previous subsection is called FMO2, with "2" indicating that the energy is expanded up to 2-body terms of fragments. It is possible to improve the precision of FMO by adding 3-body, 4-body, ..., and n-body terms (FMOn) at the expense of the computation cost of O(1). FMO3 has been implemented in both GAMESS and ABINIT-MP. The improvement by FMO3 is especially apparent in FMO-MD, as exemplified by a simulation of proton transfer in water . Recently, FMO4 was implemented in ABINIT-MP (Nakano et al., 2012), which will presumably make it possible to regard even a metal ion as a fragment.

Second-order Moeller-Plesset perturbation (MP2)
The HF calculation neglects the electron correlation effect, which is necessary to incorporate the so-called dispersion term. The electron correlation can be calculated fairly easily by the second-order Moeller-Plesset perturbation (MP2). Though the MP2/FMO energy formula was published as early as 2004 (Fedorov et al., 2004b;Mochizuki et al., 2004ab), the energy gradient formula for MP2/FMO was first published in 2011 by Mochizuki et al. (2011) and then by . In Mochizuki's implementation of MP2 to ABINIT-MP, an integral-direct MP2 gradient program module with distributed parallelism was developed for both FMO2 and FMO3 levels, and a new option called "FMO(3)" was added, in which FMO3 is applied to HF but FMO2 is applied to MP2 to reduce computation time, based on the relatively short-range nature of the electron correlation compared to the range of the Coulomb or electrostatic interactions.
The MP2/FMO gradient was soon applied to FMO-MD of a droplet of water molecules . The water was simulated with the 6-31G* basis set with and without MP2, and the resultant trajectories were subjected to calculations of radial distribution functions (RDF). The RDF peak position of MP2/FMO-MD was closer to the experimental 7 value than that of HF/FMO-MD was. This result indicated the importance of the correlation energy incorporated by MP2 to describe a condensed phase.

Configuration Interaction Singles (CIS)
CIS is a useful tool to model low-lying excited states caused by transitions among near HOMO-LUMO levels in a semi-quantitative fashion (Foresman et al., 1992). A tendency of CIS to overestimate excitation energies is compensated for by CIS(D) in which the orbital relaxation energy for an excited state of interest as well as the differential correlation energy from the ground state correlated at the MP2 level (Head- Gordon et al., 1994). Both CIS and CIS(D) have been introduced to multilayer FMO (MFMO; Fedorov et al., 2005) in ABINIT-MP (Mochizuki et al., 2005a(Mochizuki et al., , 2007a. Very recently, Mochizuki implemented the parallelized FMO3-CIS gradient calculation, based on the efficient formulations with Fock-like contractions (Foresman et al., 1992). The dynamics of excited states is now traceable as long as the CIS approximation is qualitatively correct enough. The influence of hydration on the excited state induced proton-transfer (ESIPT) has been attracting considerable interest, and we have started related simulations for several pet systems such as toropolone.

Unrestricted Hartree-Fock (UHF)
UHF is the simplest method for handling open-shell molecular systems, as long as care for the associated spin contamination is taken. The UHF gradient was implemented by preparing and β-density matrices. Simulation of hydrated Cu(II) has been underway at the FMO3-UHF level, and the Jahn-Teller distortion of hexa-hydration has been reasonably reproduced (Kato et al., in preparation). The extension to a UMP2 gradient is planned as a future subject, where the computational cost may triple the MP2 gradient because of the three types of transformed integrals, (,), (,), and (,) (Aikens et al., 2003).

Model Core Potential (MCP)
Heavy metal ions play major roles in various biological systems and functional materials. Therefore, it is important to understand the fundamental chemical nature and dynamics of the metal ions under physiological or experimental conditions. Each heavy metal element has a large number of electrons to which relativistic effects must be taken into account, however. Hence, the heavy metal ions increase the computation cost of high-level electronic structure theories. A way to reduce the computation is the Model Core Potential (MCP; Sakai et al., 1987;Miyoshi et al., 2005;Osanai et al., 2008ab;Mori et al., 2009), where the proper nodal structures of valence shell orbitals can be maintained by the projection operator technique. In the MCP scheme, only valence electrons are considered, and core electrons are replaced with 1-electron relativistic pseudo-potentials to decrease computational costs. The MCP method has been combined with FMO and implemented in ABINIT-MP (Ishikawa et al., 2006), which has been used in the comparative MCP/FMO-MD simulations of hydrated cis-platin and trans-platin (see subsection 3.6). Very recently, the 4f-in-core type MCP set for trilvalent lanthanides has been developed and made available (Fujiwara et al., 2011).

Periodic Boundary Condition (PBC)
PBC was finally introduced to FMO-MD in the TINKER/ABINIT-MP system by Fujita et al. (2011). PBC is a standard protocol for both classical and ab intio MD simulations.
Nonetheless, partly due to the complexity of PBC in formulation but mostly due to its computation cost, FMO-MD simulations reported in the literature had been performed under a free boundary condition, usually with a cluster solvent model restrained by a harmonic spherical potential. This spherical boundary has the disadvantage of exposing the simulated molecular system to a vacuum condition and altering the electronic structure of the outer surface . Hence, PBC is expected to avoid the disadvantage and to extend FMO-MD to simulations of bulk solvent and crystals. For PBC simulations to be practical, efficient approximations in evaluating the ESP matrix elements will need to be developed. A technique of multipole expansion may be worth considering.

Miscellaneous
Analytic gradient formulae have been derived for several FMO methods and implemented in the GAMESS software, including those for the adaptive frozen orbital bond detachment scheme (AFO; ), polarizable continuum model method (PCM; Li et al., 2010), time-dependent density functional theory (TD-DFT; Chiba et al., 2009), MFMO with active, polarisable, and frozen sites , and effective fragment potential (EFP; Nagata et al., 2011c). Also, Ishikawa et al. (2010) implemented partial energy gradient (PEG) in their software PACIS. These gradients have been used for FMO-EM calculations of appropriate molecules. Among them, the EFP gradient has already been applied successfully to FMO-MD (Nagata et al., 2011c), and the others will be combined with FMO-MD in the near future.

MD
The MD portion of FMO-MD resembles the conventional classical MD method, but several algorithms have been introduced to facilitate FMO-MD.

Dynamic Fragmentation (DF)
DF refers to the redefinition of fragments depending on the molecular configuration during FMO-MD. For example, in an H + -transfer reaction (AH + + B → AHB + → A + BH + ), AH + and B can be separate fragments before the reaction but should be unified in the transition state AHB + , and A and BH + may be separated after the reaction. The DF algorithm handles this fragment rearrangement by observing the relative position and nuclear species of the constituent atoms at each time step of a simulation run.
The need for DF arose for the first time in an FMO-MD simulation of solvated H 2 CO (Mochizuki et al. 2007b; see subsection 3.1). During the equilibration stage of the simulation, an artifactual H + -transport frequently brought about an abrupt halt of the simulation. To avoid the halt by the H + -transport, T. Ishikawa developed a program to unite the donor and acceptor of H + by looking up the spatial formation of the water molecules. This program was executed at each time step of the simulation. This was the first implementation of the DF algorithm (see Komeiji et al., 2009a, for details The DF algorithm was generalized later to handle arbitrary molecular systems . The algorithm requires each atom's van der Waals radius and instantaneous 3D coordinate, atomic composition and net charge of possible fragment species, and certain threshold parameters.
Presently, PEACH has four fragmentation modes, as follows: Mode 0: Use the fragmentation data in the input file throughout the simulation. Mode 1: Merge covalently connected atoms, namely, those constituting a molecule, into a fragment.
Mode 2: Fragments produced by Mode 1 are unified into a larger fragment if they are forming an H-bond. Mode 3: Fragments produced by Mode 2 are unified if they are an ion and coordinating solvent molecules.
The modes are further explained as follows. Heavy atoms located significantly close to each other are united as a fragment, and each H atom is assigned to its closest heavy atom (Mode 1). Then, two fragments sharing an H atom are unified (Mode 2). Finally, an ion and surrounding molecules are united (Mode 3). See Figure 2 for typical examples of DF. Usually, Mode 1 is enough, but Mode 2 or 3 sometimes become necessary. The DF algorithm gracefully handles molecular systems consisting of small solute and solvent molecules, but not those containing large molecules such as proteins and DNA, which should be fragmented at covalent bonds. Currently, Mode 0 is the only choice of fragmentation for these large molecules, in which the initial fragmentation should be used throughout and no fragment rearrangement is allowed (Nakano et al., 2000;Komeiji et al., 2004). This limitation of the DF algorithm will be abolished soon by the introduction of a mixed algorithm of DF and a static fragmentation.

Blue moon ensemble
The blue moon ensemble method (Sprik & Ciccotti, 1998) is a way to calculate the free energy profile along a reaction coordinate (RC) while constraining RC to a specified value.
The method was implemented in FMO-MD (Komeiji, 2007) and was successfully applied to drawing a free energy profile of the Menschutkin reaction (Komeiji et al., 2009a).

Path Integral Molecular Dynamics (PIMD)
The nuclei were handled by the classical mechanics in most of the FMO-MD simulations performed to date (Fig. 1), but PIMD (Marx & Parrinello, 1996) has been introduced into FMO-MD to incorporate the nucleic quantum effect (Fujita et al., 2009). FMO-PIMD consumes tens of times more computational resource than the classical FMO-MD does but is necessary for a better description of, for example, a proton transfer reaction.

Miscellaneous
Miscellaneous MD methods implemented in the PEACH/ABINIT-MP system include the Nosé-Hoover (chain) thermostat, RATTLE bond constraint, RC constraint, spherical solvent boundary, and so on (Komeiji et al., 2009a). Another research group has implemented the Hamiltonian Algorithm (HA) to FMO-MD to enhance conformation sampling of, for example, polypeptides (Ishimoto et al., 2004(Ishimoto et al., , 2005Tamura et al., 2008).

Applications of FMO-MD
FMO-MD has been extensively applied to hydrated small molecules to simulate their solvation and chemical reactions. Some benchmark FMO-MD simulations were described briefly in the previous section. In this section, we review genuine applications of FMO-MD in detail.

Excitation energy of hydrated formaldehyde
FMO-MD and MFMO-CIS(D) were combined to evaluate the lowest n* excitation energy of hydrated formaldehyde (H 2 CO) molecules (Mochizuki et al., 2007b). The shift of excitation energy of a solute by the presence of a solvent, known as solvatochronism, has drawn attention of both experimentalists and theorists and has been studied by various computational methods, mostly by the quantum mechanics and molecular mechanics (QM/MM) method. Alternatively, Mochizuki et al. (2007b) tried a fully ab initio approach, in which FMO-MD sampled molecular configurations for excited calculations. In the configuration sampling, H 2 CO was solvated within a droplet of 128 water molecules ( Fig. 3 left), and the molecular system was simulated by FMO-MD at the FMO2-HF/6-31G level to generate a 2.62-ps trajectory at 300 K. From the last 2-ps portion of the trajectory, 400 conformations were chosen and were subjected to MFMO-CIS(D) calculations at the FMO2/HF/6-31G* level. In MFMO, the chromophore region contained H 2 CO and several water molecules and was the target of CIS(D) calculation. The calculated excitation energy was averaged over the 400 configurations (Fig. 3 right). A similar protocol was applied to an isolated H 2 CO molecule to calculate the excitation energy in a vacuum. The blue-shift by solvatochromism thus estimated was 0.14 eV, in agreement with preceding calculations.
The solvatochromism of H 2 CO is frequently challenged by various computational methods, but this study distinguishes itself from preceding studies in that all the calculations were fully quantum, without classical force field parameters.

Hydrolysis of a methyl diazonium ion
The hydrolysis of the methyl-diazonium ion (CH 3 N 2 + ) is an S N 2-type substitution reaction that proceeds as follows: Traditionally, this reaction is believed to occur in an enforced concerted mechanism in which a productive methyl cation after N 2 leaving is too reactive to have a finite lifetime, and consequently the attack by H 2 O and the bond cleavage occur simultaneously. This traditional view was challenged by Sato et al. (2008) using FMO-MD. The FMO-MD simulations exhibited diverse paths, showing that the chemical reaction does not always proceed through the lowest energy paths.
This reaction was simulated as follows. FMO-MD simulations were conducted at the FMO2/HF/6-31G level. CH 3 -N 2 + was optimized in the gas phase and then hydrated in a sphere of 156 water molecules. The water was optimized at 300 K for 0.5 ps with the RATTLE bond constraint. The temperature of the molecular system was raised to 1000 K, and the simulation was continued for 5 ps. From the 1000 K trajectory, 15 configurations were taken and subjected to a further run at 700 K without any constraint. Ten trajectories out of fifteen produced the final products (CH 3 -OH 2 + +N 2 ). The ten productive trajectories were classified into three groups: tight S N 2, loose S N 2, and intermediate. Trajectory A in Fig. 4 is of the tight S N 2 type, in which the attack by H 2 O and C-N bond cleavage, i.e. release of N 2 , occur concertedly. Trajectory B is of the loose S N 2 type, which shows a two-stage process in which C-N bond cleavage precedes the attack by H 2 O.
The difference between trajectories A and B was further analyzed by the configuration analysis for fragment interaction (CAFI; Mochizuki et al., 2005b), and the results are plotted in Fig. 5. Charge-transfer (CT) interaction between the two fragments increases rapidly when the C-N distance increases to 1.6 Å for trajectory A, but for trajectory B the CT increased only when R CN was 2.4 Å or longer. In trajectory B, the C-N bond cleavage and O-C bond formation events take place in a two-stage fashion. The CT interaction energy is larger for trajectory B than for A at R C-O = 2.6 Å, because at the same C-O distance the C-N bond is cleaved to a larger extent, and hence the CH 3 moiety has more positive charge for trajectory B than for trajectory A.
Most of the other productive trajectories exhibited intermediate characteristics between those of trajectories A and B. The diversity of the reaction path can be illustrated by the twodimensional R C-N -R O-C plot (Fig. 6). The existence of different paths indicates that the reaction does not always proceed through the lowest energy pathway with optimal solvation.
In summary, this series of simulations illustrated for the first time how the atoms in reacting molecules, from reactant to product, behave in solution at the molecular level. This was made possible by the advent of the full ab initio FMO-MD method.
www.intechopen.com  Fig. 6. R C-N -R O-C plot of the ten trajectories that resulted in product formation. Those trajectories that proceeded along the diagonal line are regarded as tight S N 2, in which attack by water and the exit of N 2 occurred simultaneously, while a trajectory that deviated from the diagonal line is regarded as loose S N 2, in which N 2 left before the attack by water. Reproduced from Sato et al. (2008) by permission. Sato et al. (2010) tackled the reaction mechanism of the amination of H 2 CO by FMO-MD simulations. In particular, they focused on whether the reaction proceeds via a zwitterion (ZW) intermediate (Fig. 7). The results indicated that the reaction proceeds through a stepwise mechanism with ZW as a stable intermediate. The FMO-MD simulations were designed as follows. RC was defined as R N-C -R N-H . With RC constrained, structural changes of the reactant (RT) molecules in MD simulations are confined to the line that has the slope=1 and intercept=RC in a More O'Ferrall-Jencks-type diagram (Fig. 8). This diagram allows the reader to distinguish between the stepwise process and the concerted one. for RC = 0.9, 1.0,..., 1.8 Å starting from configuration B. For each RC value the configuration was equilibrated for 0.3 ps and sampled for a further 0.3 ps.

Amination of formaldehyde
The diagram thus obtained clearly favored the stepwise mechanism over the concerted mechanism (Fig. 9). Nevertheless, there remained a possibility of the MD trajectory being trapped in a local minimum. To investigate the possibility, we conducted additional FMO-MD simulations starting from configuration C, the concerted TS-like one. These additional trajectories all diverted from the TS-like structure toward the trajectory of the stepwise path (see Sato et al., 2010, for details), thus confirming the validity of the stepwise mechanism. In summary, the constraint FMO-MD simulations indicated that the H 2 CO amination in water solvent occurs by the stepwise mechanism, not by the concerted one.

Hydration of Zn(II)
The divalent zinc ion, Zn(II), plays bio-chemically relevant roles, e.g., as the reaction center of superoxide dismutase. By using a droplet model of the Zn(II) ion with 64 water molecules, FMO2-and FMO3-MD simulations were performed at the HF/6-31G level, supposing that the electrostatic and coordination interactions are dominant in this system (Fujiwara et al., 2010b). The Zn-O peak positions at the first hydration shell were investigated, and a better accuracy of FMO3-MD than that of FMO2-MD was demonstrated, where the FMO3 value of 2.05 Å agreed well with the experimental value of 2.06±0.02 Å (Fig. 10). The coordination number of the first hydration shell was 6 consistently. Additionally, the charge fluctuations on the Zn atom were evaluated by the natural population analysis (NPA) as well as the conventional Mulliken population analysis (MPA).
The NPA results showed a consistent picture with the coordination bond with reasonable fluctuation (around a net charge of 1.8), while MPA yielded an artificially enhanced fluctuation with a larger extent of electron donation (net charge of 1.3-1.4). Discussion with NPA was found to be preferable for hydrated metal ions.

Hydration of Ln(III)
The lanthanide contraction and the gadolinium break have attracted considerable attention in the inorganic chemistry. As an application of the 4f-in-core MCP (Fujiwara et al., 2011), a series of FMO3-MD simulations on droplet model of Ln(III) plus 64 water molecules have been underway at the HF level (Fujiwara et al., in preparation). The RDF peal positions for La(III) (nona-hydration) and Lu(III) (octa-hydration) were estimated to be 2.59 Å and 2.31 Å, respectively, and they were comparable to the corresponding experimental values of 2.54 Å and 2.31 Å. Interestingly, the octa-and nona-hydration results for Gd(III) were evaluated as 2.46 Å and 2.53 Å, respectively. The former value is in closer agreement with the experimental value of 2.42 Å, suggesting that the octa-hydration is preferable.

Comparison on hydration dynamics of cis-and trans-platin
FMO-MD has also given important insight into the difference in the hydration dynamics of cis-and trans-platin (Mori et al., 2012). Since cis-platin (cis-[Pt II Cl 2 (NH 3 ) 2 ]) is recognized as an anticancer substance, quite a few studies have been devoted to the biochemical functions of its derivatives. Particularly interesting in the pharmaceutical research field of Pt-based anticancer drugs is the behaviour of its geometrical isomer, trans-platin, which only shows very low anticancer activity (Fig. 11). Trans-platin had not been considered to form DNA adducts that lead to anticancer activity. However, trans-type Pt-complexes that shows antitumor activities was found recently. Despite the extensive research on both cis-and trans-platin, the origin of their difference in biochemical activity still remains unclear. The final step of the antitumor treatment is the combination of cis-platin and DNA leading modifications of the DNA structure. Meanwhile, some earlier steps, such as solvation before reaching the final target, are also believed to play important roles in the efficacy of drugs.
Their hydration should be investigated to understand the difference in the medical application between cis-and trans-platins. FMO-MD simulations were performed for hydrated cis-and trans-platins. The simulation conditions were set as described below. Each platin complex was hydrated with a spherical droplet of water centred at the Pt atom with a diameter of 10.5 Å. This diameter was determined to include up to the second solvation shell, so that the physicochemical properties of the first shell should be reproduced. In the FMO-MD simulations, the electronic states of the hydrated platin complexes were described by FMO(3)-MP2. The basis sets were MCPdz for Pt, MCPdzp for Cl, and 6-31G(d) for the others, respectively. The MCP basis sets were applied for heavy elements (see subsection 2.2.6). The central platin and each of the water molecules were regarded as independent fragments. DF was applied to allow for the generation of proton-transferred species during the production MD runs. For each cis-and trans-platin system, a 1-ps equilibration and a subsequent 2-ps production MD run were performed using the Nose-Hoover Chains NVT ensemble at 300 K. NPA was also performed during the FMO-MD run to analyze the differences in charge fluctuations between cis-and trans-platin, illuminating the differences in the hydration environment around polarized Pt  + -Cl  -bonds, which should be cleaved by the nucleophilic attack of a solvent water molecule.
The time evolution of the natural charge on each ligand in cis-and trans-platin, and that of Pt-Cl bond lengths are shown in Fig. 12. Relatively larger charge fluctuations were observed on the Pt/Cl sites than on the NH 3 sites in both platins. This difference among the sites was attributable to the fact that NH 3 has no amplitude in the highest occupied molecular orbital. A close comparison of the left and right graphs in Fig. 12 revealed a correlation between fluctuation of the Pt/Cl sites and that of the Pt-Cl bond. By applying the Fourier transform technique to the charge fluctuation, we calculated the frequency of the fluctuation to be 334 cm -1 . This frequency can also be assigned to the Pt-Cl stretching mode coupled with intermolecular vibrations between the solute platin and solvent water molecules. The correlation observed in charge fluctuation on Pt and Cl sites means that there is a CT interaction between them. Since the frontier MO that participates in the CT process is a Pt-Cl antibonding orbital, the CT interaction coupled with the fluctuation of the solvent water should induce a Pt-Cl bonds fluctuation. Since trans-platin has inversion symmetry, the dipole moment of trans-platin is much smaller than that of cis-platin. This means that the number of water molecules which coordinates to the platin complex is larger for cis-platin than for trans-platin. Thus, the CT interaction coupled with the solvent motion is stronger in cis-platin than in trans-platin. As a result, the Pt-Cl bonds are easier to elongate for the cleavage in the hydrated cis-platin than in the hydrated trans-platin. Thus, by using FMO-MD simulations, we obtained new quantum chemical insight into the solvation of platin complexes.

Prospects and conclusion
As reviewed so far, FMO-MD has been applied to various chemical phenomena in the presence of explicit solvents and has given realistic molecular pictures of the phenomena. We are planning to extend the field of FMO-MD by introduction of new capabilities, as follows.
The so-called QM/MM scheme will enhance the target size of FMO-MD. QM/MM has attracted great interest in simulating condensed-phase systems as well as proteins. In this scheme, the chemically relevant region is subjected to QM calculations while the environmental effects are incorporated through a set of MM parameters. MFMO has a conceptual similarity to QM/MM, and hence we have a plan to implement a general QM/MM ability in conjunction with MFMO.
The improvement of accuracy in FMO gradient evaluations may be a future subject. Nagata's reformulation, including the supplemental response terms of monomers (Nagata et al., 2011a) as well as the BDA-related residual contributions , are of interest for implementation at the HF level.
Another important issue is the extraction of more information from FMO-MD trajectories. From a series of configurations, the time-dependent fluctuations in electronic densities can be derived, some of which are correlated with the creation and destruction of bonding interactions. For example, the Fourier transform-based analyses may shed light on the detailed dynamical picture of nucleophilic attack reactions.
In conclusion, FMO-MD is a highly-parallelizable ab initio MD method. FMO-MD has advanced rapidly by improvement of both the FMO and MD portions of the method and has been successfully applied to various chemical phenomena in solution. We are planning to extend the methodology and application of FMO-MD by incorporating several new features.

Acknowledgment
Thanks are due to Dr. Makoto Sato, Mr. Takayuki  Molecular Dynamics is a two-volume compendium of the ever-growing applications of molecular dynamics simulations to solve a wider range of scientific and engineering challenges. The contents illustrate the rapid progress on molecular dynamics simulations in many fields of science and technology, such as nanotechnology, energy research, and biology, due to the advances of new dynamics theories and the extraordinary power of today's computers. This first book begins with a general description of underlying theories of molecular dynamics simulations and provides extensive coverage of molecular dynamics simulations in nanotechnology and energy. Coverage of this book includes: Recent advances of molecular dynamics theory Formation and evolution of nanoparticles of up to 106 atoms Diffusion and dissociation of gas and liquid molecules on silicon, metal, or metal organic frameworks Conductivity of ionic species in solid oxides Ion solvation in liquid mixtures Nuclear structures