Molecular Dynamics Simulations of Nimodipine ... - ACS Publications

Subscriber access provided by Access provided by University of Liverpool Library

C: Surfaces, Interfaces, Porous Materials, and Catalysis

Molecular Dynamics Simulations of Nimodipine Confined in an Ordered Mesoporous Silica Matrix Aleksandra Pajzderska, Miguel Angel Gonzalez, and Jan Wasicki J. Phys. Chem. C, Just Accepted Manuscript • Publication Date (Web): 10 Apr 2018 Downloaded from http://pubs.acs.org on April 10, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Molecular Dynamics Simulations of Nimodipine Confined in an Ordered Mesoporous Silica Matrix A. Pajzderska1*, M. A. Gonzalez2, J. Wąsicki1,3

1

Faculty of Physics, A. Mickiewicz University, ul. Umultowska 85, 61-614 Poznan, Poland Institute Laue Langevin, B.P. 156x, 38042 Grenoble Cedex 9, France 3 The NanoBioMedical Centre in Poznań, Umultowska 85, Poznan, Poland 2

* corresponding author: [email protected]

Abstract The structural and dynamical properties of crystalline nimodipine and nimodipine confined in mesoporous SBA-15 have been studied by means of molecular dynamics simulations. As those have been motivated by nuclear magnetic resonance measurements the percentage of filling of the silica matrix was chosen to be comparable with the experimental value, i.e. only about 25% of the silica surface was covered by nimodipine. We find that nimodipine molecules connect to the silica surface by different types of hydrogen bonds, which are formed, broken and rebuilt continuously during the simulation. These interactions inhibit translational diffusion of nimodipine, but not the reorientational dynamics. Direct comparison of spin-lattice relaxation times T1 extracted from the simulations with experimental ones shows a very good agreement for both systems. MD simulations indicate that the NMR data for the crystal can be perfectly explained in terms of the rotation of the five methyl groups present in the nimodipine molecule, while when confined in the silica matrix the molecules of nimodipine present additional degrees of reorientational freedom and that they must be taken into account in order to reproduce the NMR signal.

ACS Paragon Plus Environment

1

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 26

1. Introduction Porous materials are of great interest because of their ability to interact with a variety of substances like isolated atoms or molecules, liquids, nanoparticles, etc. Due to their tunable surface area, porosity, pore volume, pore shape and framework compositions, mesoporous materials are currently under intense investigation. Some of the main areas where these materials are applied include catalysis, adsorption, or drug delivery.1-3 Among the wide variety of existing mesoporous systems, ordered mesoporous silica materials are of particular interest. The best known examples are MCM-41 and SBA-15, which have an ordered hexagonal structure with mesopores of sizes ranging from 2 to 10 nm.4 They can exhibit large surface areas (up to about 1000 m2/g), large pore volumes (close to 1 cm3/g), small dispersion of pore diameters and high thermal stability. For the above-mentioned reasons, porous silica materials are very good adsorbents and can be used as hosts for the adsorption of guest molecules. It is also well known that the properties of the guest molecules will differ from those in the bulk state and the changes in the dynamics of confined liquids (e.g. water, benzene, toluene, acetonitrile or liquid crystals) relative to the dynamics in the bulk have been extensively studied.512

In 2001 it was shown for the first time that MCM-41 can be applied as a drug delivery

system.13 Since then, many experimental studies have shown the possibilities of adsorption of different therapeutic substances.14-20 We have recently studied the dynamics of nimodipine (a calcium channel blocker21,22 from the group of 1,4-dihydropyridinium) confined in ordered silica mesoporous SBA-15 by means of

13

C and 1H solid state nuclear magnetic resonance (NMR) methods23. The spin-lattice

relaxation time T1 and second moment measurements of the NMR line suggested that nimodipine molecules have additional degrees of rotational freedom with respect to the crystal. Motivated by these findings, we have studied the dynamics of crystalline and confined nimodipine by means of molecular dynamics (MD) simulations and present here the results of such work.

2. Methods MD simulations were performed using DL_POLY_classic24 for crystalline nimodipine (the molecular geometry is schematically shown in fig. 1) and nimodipine confined in a silica matrix. Crystalline nimodipine was simulated using a supercell constructed on the basis of X-ray


2

Page 3 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


diffraction (XRD) data25 of dimensions 41.9 Å x 44.0 Å x 44.6 Å and containing 144 nimodipine molecules. In order to simulate the confined system, we carved a single cylindrical pore of diameter 62 Å corresponding to the average dimensions of the pore calculated from the adsorption−desorption method23 in an amorphous silica simulation box following the same procedure as in our previous paper.26 Briefly, using graphical tools built in the Materials Studio 2016 package27 we generated a block of amorphous silica and all atoms within a cylinder of diameter 62 Å along the z-axis were removed. Then, all silica atoms which did not have a complete tetrahedral environment were also removed. Finally, non-bridging oxygen atoms were saturated with hydrogen atoms. The surface obtained in such a way is irregular and rough and we can estimate that the concentration of hydroxyl groups is 7.5 per nm2. In the next step the Adsorption Locator module in Materials Studio 2016 was used to add 50 nimodipine molecules to the simulated pore, corresponding approximately to the percentage of filling obtained experimentally. Using the generated configuration as the starting point, MD simulations were then performed. NPT (for crystalline cluster) and NVT (for nimodipine in silica matrix) simulations were done at 3 different temperatures (250 K, 225 K, and 200 K), using Berendsen’s thermostat and barostat28 with a relaxation constant of 1 ps for the temperature and 5 ps for the pressure. Periodic boundary conditions were applied in all directions. The silica frame was kept rigid by freezing the positions of O and Si atoms, while rotations of the H atoms around Si−O bonds were allowed and the O−H distance was kept fixed at 0.95 Å using the SHAKE algorithm. The GAFF force field was used for nimodipine29, while for the silica matrix we used the Lennard-Jones potential parameters proposed by30 and geometric mixing rules were applied. The partial charges for Si was 2.4 e and for O was equal to -1.2 e in silica matrix. This silica forcefield has been already successfully used in comparable cases, such as to study water dynamics in silica nanopores31 and in MCM-4126, as well as the behaviour of confined glucose solutions32. A cutoff distance of 20 Å (crystalline cluster) and 30 Å (silica matrix) was applied for the van der Waals forces and the electrostatic interactions were treated using the Ewald summation method with the same cutoff in real space. In all cases a time step of 1 fs was used and the systems were equilibrated over 2 ns. Then the trajectory was saved every 5 ps for a total simulation time of 10 ns. The analysis of the trajectories was performed with nMoldyn 333 and self-written programs.


3


Page 4 of 26

3. Results 3.1. Distribution of nimodipine molecules in the silica matrix The first step consisted in determining the distribution of nimodipine molecules inside the silica pore. Figures 1 and 2 show two snapshots corresponding to the pure crystal and the confined system, respectively. In the latter the amount of molecules inside the pore is insufficient to form one complete layer and, as expected, long-range order is lost. The molecular distribution analysis shows that nimodipine molecules are not distributed homogeneously all around the pore surface. It should be noted that after the simulated adsorption step, nimodipine molecules are randomly distributed, but that during the minimization they tend to aggregate and we observed the formation of clusters containing between 5 and 7 molecules which are attached to the SBA-15 surface and therefore remain close to the pore wall. We can also estimate that only 25 % of the surface of the pore is covered.

Fig. 1. Left: Snapshot of the crystal of nimodipine obtained from the simulation performed at 250 K. Gray spheres represent C, red O, blue N and white H atoms. Right: Sketch of the molecular geometry of nimodipine along with the atoms and methyl groups notation adopted throughout his work.


4

Page 5 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


a)

b)

Fig. 2. Snapshot of a nimodipine/SBA configuration obtained from the simulation performed at 250 K. Yellow and red lines correspond to Si and O atoms in the silica matrix. Gray, red, blue, and white spheres represent C, O, N, and H in nimodipine. a) View along z axis, b) View along y axis (the silica matrix is presented in semi-transparent mode) The average distances over the trajectory between molecules in the crystal and in confined nimodipine were also analysed. We chose atom C4 as a reference (as it is nearly the centre of mass of the molecule) and we calculated the distance between C4 atoms belonging to two different molecules. In the crystal the distance between two nearest molecules is equal to 7.7 Å and oscillates about ±0.6 Å. In the silica matrix most of the neighboring molecules are separated from each other by 12 - 15 Å. This distance does not change with time, and only small oscillations are visible, but for some pairs of molecules we can observe occasional jumps of 2-3 Å. Only a few pairs of nimodipine molecules show C4-C4 distances that are smaller than in the crystal, approximately 6 Å. The lack of translational mobility of the nimodipine molecules is also visible in Fig. 3, which shows a projection of the trajectory of C4 atoms on the xy and xz planes. All the nimodipine molecules in the matrix are close to the silica surface. The strong interaction with it inhibits long-range diffusion and as shown in Fig. 3, only a restricted motion of about 2 Å is observed. A similar behavior has been found for a small amount of confined water in MCM26.


5


Page 6 of 26

Fig 3. Left: Projection in the xy plane of the positions of all C4 atoms (red) of nimodipine during the full trajectory obtained at 250 K in the confined system. The blue points show the individual trajectories of atom C4 for one representative nimodipine molecule. Right: Same as in the left panels, but projected on the xz plane.

To determine how nimodipine molecules are oriented relative to the silica matrix, the angle between the plane defined by the 1,4-dihydropyridinium ring (atoms N1-C2-C3-C4-C5C6) and the axis of the channel (z-axis) was calculated for all the molecules as a function of time. For 80% of the molecules (40 molecules out of 50) the determined angle is in the range of −30 to 20 deg. On this basis we can say that the 1,4-dihydropyridinium ring tends to align parallel or almost parallel to the silica surface. For half of the molecules, this angle does not change with time and only fast oscillations of approximately ± 10 deg are observed. For the other half, instantaneous (in a simulation time-scale) jumps of approximately 20 degrees are observed, indicating a reorientation of the entire molecule in a limited angular range. In figure 4 we show an example of such reorientational jumps. It is worth noting that similar reorientations are not possible in the crystalline box.


6

Page 7 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Fig. 4. Left and right panels show one selected nimodipine molecule at 2.5 ns and 4.5 ns respectively. The middle panel shows the time dependence of the angle between the plane of the 1,4-dihydropyridinium ring and the OZ axis (axis of the channel). Hydrogen bonds between nimodipine and the silica OH groups are also shown as blue dotted lines. In the left figure (2.5 ns) the molecule of nimodipine is bonded to the surface by O9a…H-O, O15a…H-O, O16a…H-O and O19a…H-O HBs; while in the right one (4 ns) we have and N1-H1a…O, O15a…H-O, O19a…H-O HBs. 3.2. Hydrogen bond analysis and interactions with the surface The next step of our study consisted in the analysis of the interactions between nimodipine molecules and the silica surface. To identify hydrogen bonds (HBs) between neighboring nimodipine molecules and/or between nimodipine and the silica surface we used the geometric criterion.34,35 According to this criterion, an HB exists when the distance d between the hydrogen and the acceptor heavy atoms is less than 2.5 Å, and the angle between the hydrogen atom, the donor heavy atom and the acceptor heavy atom is less than 30 degrees. In the nimodipine molecule, only the N1 atom can act as a HB donor and in the crystal structure, molecules are bound by N-H…O (N1-H1a_O19a) HBs. By analyzing how this HB behaves as a function of time during the simulation of the crystal, we can say that its length (distance from hydrogen to oxygen) oscillates by about ± 0.4 Å around its average value of 2.0 Å and that no HB breakings are observed in our trajectory. As mentioned above, nimodipine molecules confined in the SBA matrix form small clusters. Thus, HBs can be formed between neighbouring nimodipine molecules as well as between nimodipine and the SBA surface. However from the analysis of the whole set of saved configurations we found that only a single N-H ... O HB formed between two molecules and lasted only for 1.6 ns. Thus we conclude that the formation of HBs between neighboring nimodipine molecules in the confined system is negligible, at least at the low filling fractions


7


Page 8 of 26

explored here. This is consistent with the large distances that on average we found between the centers of gravity of neighboring molecules. On the other hand, nimodipine has seven oxygen atoms (noted as O9a, O9b, O15a, O15b, O16, O16b, O19a) that can act as acceptors to form an HB with an hydroxyl group located on the silica surface. In addition, the nitrogen atom in the 1,4-dihydropyridinium ring may also donate its hydrogen H1a to form an HB with SBA oxygen atoms. Figure 5 illustrates the hydrogen bonds behavior as a function of time for a selected molecule. In this case, during most of the trajectory the molecule forms three stable HBs (O16a…H-O, O15a…H-O and O16a…H-O), with some extra HBs forming occasionally but not lasting for long. As an example, figure 3 shows the selected molecule of nimodipine at 2.5 ns and 4 ns. Initially two HBs exist (O19a…H-O and N1H1a…O), while after reorientation N1-H1a…O is broken and two new HB appear: O16b…H-O and O15a…H-O. A more detailed analysis reveals that on average each nimodipine molecule forms 2.3 HBs with the surface, but there is a large variation with about 20% of them forming on average less than 1.5 HBs and another 20% forming more than 3 HBs. Preferred HBs are formed between O16a and O15a atoms and the hydroxyl groups while the less common HBs are O16b…H-O and O15b…H-O. It is interesting to mention that ab initio molecular dynamics simulations for ibuprofen confined in silica matrix36 have also revealed that hydrogen bonds between drug and silica are continuous forming, breaking and reforming during simulations.


8

Page 9 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Fig. 5. Distances between hydrogen atom and acceptor as a function of simulation time for 8 hydrogen bonds for a selected molecule at 250 K. Only bonds O16a…H-O, O15a…H-O and O19a…H-O remain intact during the full trajectory. Other HBs are formed, broken and rebuilt continuously during the simulation.

3.3 Torsional angles distribution


9


Page 10 of 26

The earlier analysis, especially the variation of the angle between the plane of the molecule and the axis of the channel and the breaking and formation of HBs suggest that confined molecules maintain some degree of rotational freedom and can perform either whole rotations around a given axis or center of gravity or partial reorientations of molecular fragments around the individual C-O and C-C bonds. This conclusion is consistent with the NMR results (measurements of the second moment of NMR line and T1 relaxation time) which showed that the mobility of nimodipine in the matrix was significantly more complex than that in the crystal.23 However, experimental results did not allow an unambiguous identification of these reorientations. Thus, in order to identify them, in the next part of this paper we analyze carefully the geometry and time scale of these motions. First we focus on the local methyl group dynamics. There are five different methyl groups in the nimodipine molecule and for each of them the dihedral angle along the trajectory was computed from the positions of methyl H, C13, C2 and C3 (I); methyl H, C14, C6 and C5 (II); methyl H, C21, C17 and O15b (III); methyl H, C22, C17 and O15b (IV), and methyl H, C20, O19a and C19 (V), respectively. Figure 6 shows their time dependence for the crystal and the confined system for a randomly selected molecule. For the five types of methyl groups we have 3 possible torsional angles, showing fast oscillations around the equilibrium value and occasional instantaneous jumps between the three angles. Jump rates are different for different types of CH3, the fastest being groups I and II and the slowest IV and V. Their different time scales are analyzed later. A more global information about the geometry of this reorientation is given by the orientational distribution function P(Θ) calculated as:

P(Θ ) =

1 1 N step N group

N step

∑ N (Θ, Θ + dΘ) i

(1)

i =1

where Ni(Θ,Θ+dΘ) is the number of considered groups found with an orientation in the range (Θ,Θ+dΘ) at step i. The results are also shown in Fig. 6. In all cases, 3 preferred positions with almost equivalent populations are clearly visible, both for the pure crystal and nimodipine confined in the matrix.


10

Page 11 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Fig. 6. Time-dependence of torsional angles calculated for five methyl groups, top – in the crystal (blue points), middle – confined molecules in the matrix (red points). Bottom orientational distribution function P(Θ). Torsional angles were also calculated for the atoms belonging to isopropyl and methoxyethyl ester chains. For the isopropyl ester chain, the angular distribution of the following three dihedrals was calculated: C2, C3, C15 and O15a; C3, C15, O15b and C17; and C15, O15b, C17 and C21 (Fig. 7a). As expected, the first two show a stable planar conformation (180 deg) both in the crystal and in confined nimodipine (see both insets in Fig. 7b), so the orientational distribution function P(Θ) shows a single peak at 180°. But the third dihedral exhibits a different behaviour in the crystal and the matrix. In the crystal a single peak at 73° is observed and none of the molecules composing the simulated crystal exhibits any conformational change during the full length of our simulations. On the other hand, in the simulation of confined nimodipine, we observe that some molecules perform occasional jumps between two positions, one similar to the angle observed in the crystal (~75°) and a second one at ~164°, giving naturally two distinct peaks in the P(Θ) distribution (Figure 7c).


11


Page 12 of 26

For the methoxyethyl ester chain, 4 different dihedral angles were studied, formed by the following atoms: C5, C16, O16b and C18; C16, O16b, C18 and C19; O16b, C18, C19 and O19a; and C18, C19, O19a and C20. The first angle distribution shows a single peak at about 160° and no torsional jumps are observed during the simulation neither in the crystal nor in the confined system. Clear differences between both systems can be observed for the C16-O16b-C18-C19 dihedral distribution, shown in Fig. 9. In the crystal the distribution is significantly blurred, spanning values in the range between 55 and 165°. For confined nimodipine, instead, two clearly different peaks are observed, one centered at ~72° and another at ~177°, and jumps between both positions are occasionally observed as shown in the inset of Fig. 7d. It is also interesting to compare the distribution for the O16b-C18-C19-O19a torsion. In the crystal it shows two peaks at −63 and 65°, while in the matrix there is an additional one at 179°. It should also be noted that in the crystal the conformation at −63° is clearly preferred, and that jumps between both positions during the trajectory are scarce. On the other hand, nimodipine molecules in the silica matrix exhibit an increased mobility. Not only the additional position appears, but also jumps between the three possible conformations occur more often than in the crystal (compare insets in Figure 7e). It was also checked that for molecules for which we observe such jumps the hydrogen bond O19a ... O-H with the silica surface is either absent or it is not formed for longer than 20% of the trajectory time. Finally the distributions corresponding to the torsional angle formed by C18, C19, O19a and C20 do not show significant differences between nimodipine in the crystal and in the matrix. In both cases, a planar configuration is preferred, giving rise to a main peak centered at 180°. Only in the matrix we can observe two additional peaks at −78 and 63°, due to single jumps of a few molecules, but they have a very low occupation, as those jumps are usually followed by a quick return to the more stable planar conformation (see inset of Fig. 7f).


12

Page 13 of 26

a)

b) 50 240

45 angle Θ [deg]

180

40 35

120 60 0 -60

-120 0

1

2

3

4

5

3

4

5

time [ns]

25 240

20

180 angle Θ [deg]

P(Θ) [a.u.]

30

15 10

120 60 0 -60 -120 0

5

1

2 time [ns]

0 -120

-60

0

60

120

180

240

angle Θ [deg]

c)

d) 25

25 240

180

180

angle Θ [deg]

20

60

20

0

angle Θ [deg]

240

120

120 60 0

-60

-60 -120

15

1

2

3

4

5

P(Θ) [a.u.]

180

10

120 60

0

1

2

3

4

4

5

2 3 time [ns]

4

5

120 60 0

-120 0

0 -60

3

-60

time [ns]

-120

2

180

10

5

5

1

time [ns]

0

-120

0

240

-60

5

-120

15

time [ns] 240

angle Θ [deg]

0

angle Θ [deg]

P (Θ) [a.u.]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


0

60

120

180

1

0

240

-120

-60

angle [deg]

e)

0

60

120

180

240

angle Θ [deg]

f)


13


180

120

50 240

45

120

60 0

60

-60

-120

-120 0

1

2

3

4

40

0

-60 5

0

1

time [ns]

2

3

4

35

5

180 angle Θ [deg]

240

180

P(Θ) [a.u.]

15

10

60 0

30

-120

25

240

0

1

2

3

4

5

4

5

time [ns]

180

20 15

5

120

-60

time [ns]

angle Θ [deg]

20

240

angle Θ [deg]

angle Θ [deg]

25

P(Θ) [a.u.]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 26

120 60 0 -60

10

-120 0

5 0

1

2

3

time [ns]

0 -120

-60

0

60

120

180

240

-120

-60

angle [deg]

0

60

120

180

240

angle Θ [deg]

Fig. 7: Angular distribution function P(Θ) calculated for different dihedral angles of the nimodipine molecule in the crystal (blue squares) and confined in a silica matrix (red circles). Calculated distributions correspond to the average over all the nimodipine molecules present in the system and the insets show the time dependence along the trajectory of the corresponding dihedral angle for a single molecule. (a) 2D sketch of the nimodipine molecule showing the different dihedral angles whose distribution is presented in graphs b-f.(b) C2-C3-C15-O15a distribution. A similar behaviour is observed for C5-C16-O16b-C18 and C3-C15-O15b-C17 torsions. (c) C15-O15b-C17-C21 distribution. (d) C16-O16b-C18-C19 distribution.(e) O16bC18-C19-O19a distribution.(f) C18-C19-O19a-C20 distribution. 3.4 Time scale of the motions The time scale of the different intramolecular reorientational motions presented in the previous section was analyzed through the following angular correlation function (ACF):

ACFl (t ) = Pl [r (t0 )r (t )]

(1)

where ACFl(t) is the correlation function of order l, Pl(x) is the Legendre’s polynomial of order l and r (t ) is the vector connecting two atoms. Since the second order function (l=2) is proportional to the Fourier transform of the spectral density J(ω) measured in a NMR experiment we calculated ACF2(t) for the following motions: (i) reorientation of methyl groups, (ii) reorientation of the entire molecule and (iii) reorientation of fragments of isopropyl and methoxyethyl ester chains around specific axis. The ACF2(t) can then be fitted using a single exponential function (ACF2(t) ~ exp(−t/τ)), where τ is the correlation time, so the characteristic reorientational correlation times were extracted from the analysis of the autocorrelation functions. First, we analyzed the results of the reorientation of the methyl groups in the crystal at 250K (Fig. 8). As it can be seen, two methyl groups are characterized by short correlation times


14

Page 15 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(groups I and II), one by intermediate (group III) and two by long times (IV and V) (the corresponding correlation times at 250 K are 5.5 ps, 7.0 ps, 19 ps, 170 ps and 177 ps, respectively). A similar analysis of ACF2(t) was performed as a function of temperature and correlation times are collected in Table 1. As expected, with decreasing temperature, the methyl reorientation slows down. The temperature dependence of the average correlation time allowed us (on the basis of the Arrhenius relation) to obtain the activation energy. The reorientation of the methyl groups I and II is characterized by a low activation barrier (~5.0 kJ/mol), the barrier for the methyl group III shows a middle value (~11.5 kJ/mol), and the slowest groups IV and V present the highest activation barrier (~15.0 kJ/mol). Experimental values for the nimodipine crystal obtained on the basis of a NMR experiment supported by DFT calculations are 3.64 kJ/mol (groups I and II), 7.11 kJ/mol (group III) and 12.50 kJ/mol (IV and V) kJ/mol23,37. The agreement is reasonably good, thus validating the method of calculation, as well as the partial charges and force field employed to describe the dynamic properties of nimodipine. Table 1. The correlation times [ps] obtained from fitted ACF2(t) using a single exponential function for reorientation of methyl groups in the crystal group I group II group III group IV group V 250 K

5.5 ± 0.2

7.0 ± 0.2

19.0 ± 0.6

170 ± 5

177 ± 5

225 K

7.1 ± 0.2

10.2 ± 0.3

37.2 ± 1.1

433 ± 12

503 ± 13

200 K

10.2 ± 0.3

17.0 ± 0.5

79.2 ± 2.4

1355 ± 39

1420 ± 39

The same kind of analysis was done for nimodipine confined in the silica matrix and the correlation times are collected in Table 2, while Figure 8 shows ACF2(t) for the five methyl groups at 250 K. Comparing the correlation times, it can be seen that the one of group III is comparable with the crystal, while the other four methyl groups exhibit slightly slower dynamics than in the crystal. This is possibly due to the interaction of nimodipine molecules with the silica surface. The activation energies obtained are ~11.5 kJ/mol (groups I and II), ~10 kJ/mol (group III), and ~15 kJ/mol (groups IV and V).

Table 2. The correlation times [ps] obtained from fitted ACF2(t) using a single exponential function for reorientation of methyl groups and for rotations of fragments of chains around axis O15b-C17 and C18-C19 of nimodipine molecules in the silica matrix group I group II group III group IV group V O15b-C17 C18-C19


15


Page 16 of 26

250 K 38.0 ± 1.7

47.0 ± 1.4

19.2 ± 0.6

304 ± 9

363 ± 11

650 ± 16

450 ± 14

225 K 64.2 ± 1.9

70.0 ± 2.1

48.7 ± 1.4

715 ± 21

505 ± 15

1624 ± 49

592 ± 18

200 K 130.1 ± 3.9

96.2 ± 3.2

83.4 ± 3.3

2101 ± 45

1002 ± 31

3158 ± 69

878 ± 26

As shown above, the main difference in the mobility of confined molecules in comparison with the crystal is the possibility of having reorientations of the whole molecule in a limited angular range, as well as rotations of fragments of chains around specific axis. Therefore, the ACF2 was first calculated for the reorientation of the whole molecule around different axis and/or around the center of mass. The correlation time of this motion was found to be about two orders of magnitude longer than the correlation time of the slowest of the methyl groups. Next, ACF2 was calculated for the reorientation of the vectors connecting two atoms perpendicular to the axis O15b-C17 and C18-C19. These two axes were chosen because the torsional angles analysis showed that reorientation around theses axes are possible for confined molecules. The correlation time for first motion is equal to 650 ps, while for the second reorientation the correlation time is equal to 450 ps (Fig. 12). Analogous calculations as a function of temperature allowed to estimate the activation energy of this motion, which is equal to 6.5 kJ/mol. It should be noted that τc extracted from the ACF2 at 250 K for the same reorientation in the crystal is much longer. Decreasing the temperature causes (as opposed to confined nimodipine) these reorientations to slow down very much. Therefore, there is no possibility to determine the activation energy of these motions in the crystal.


16

Page 17 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Fig. 8. ACF2 corresponding to the reorientation of the five methyl groups and of a fragment of the nimodipine molecule around the O15b-C17 and C18-C19 axis in the crystal (blue lines) and in the confined system (red lines) at 250 K.

3.5. Comparison with experimental data


17


Page 18 of 26

NMR relaxation time measurements of T1 spin-lattice as a function of temperature have showed that T1 for confined nimodipine is much shorter than for nimodipine bulk and, also importantly, has a different slope of lnT1(1000/T). The activation energy derived from this slope for crystalline nimodipine is 9.8 kJ/mol, while for the confined system one gets 3.7 kJ/mol. Similar temperature dependences of the relaxation-time T1 were detected for other confined systems: nonpolar liquids in porous silica glasses,38,39 benzene confined in cylindrical mesoporous materials MCM-41 and SBA-15,5 liquid crystals confined in controlled pore glasses or even polymers in nanoscopic pores.40,41 It is now interesting to compare the experimental values with estimations of T1 obtained from the calculations. Spin-lattice relaxation times T1 were calculated as linear combinations of spectral densities J(w), according to equation42:

1 1 n = C ∑ = C ∑ i [J i (ω ) + J i (4ω )] T1 i T1i i N

Where:

J i (ω ) =

τ ci 1 + ω02τ ci2

, τci is the correlation time of the considered motion, ω0 stands for the

resonant frequency, C is a proportional constant, ni is the number of rotating hydrogens, and N is the total number of hydrogens in the molecule. The reciprocal of total spin-lattice relaxation time T1 (the relaxation rate) is a sum of reciprocal of spin-lattice relaxation times T1i of individual motions. Using the above equation, ω0 = 2π*59.8 MHz (spectrometer frequency), C = 1.7*109 1/s2 (which is close to experimental values equal 1.4*109 1/s2 [ref. 37]) and the correlation times for the rotation of methyl groups derived from the simulations, the value of T1 for bulk crystalline nimodipine was calculated. In the case of the reorientation of methyl groups, the fraction n/N is equal to 3/26 (3 hydrogen atoms in a methyl group, and 26 in the nimodipine molecule). The results are shown in Figure 9 and we can see that the calculated T1 is in very good agreement with experiment. However, an analogous calculation for confined nimodipine using only methyl group reorientations gave T1 values very similar to those of the crystal, clearly overestimating the experimental ones. Therefore, the motions of the nimodipine identified and described in the first part of this work, i.e. fragments of its chains were included in the calculation. For the


18

Page 19 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


reorientation around the O15b-C17 axis, the ratio is equal to n/N = 7/14, while for the reorientation around the C18-C19 axis n/N is equal to 5/26. We see (Figure 9) that taking these motions into account (with proportional constant equal C = 2.6*1010 1/s2) reduces significantly the relaxation time T1 and also the activation energy (and therefore the slope of lnT1(1000/T)), bringing the T1 extracted from MD simulations in quite good agreement with NMR.

Fig. 9. T1 relaxation time observed experimentally and calculated from MD simulations. ■ – experimental, nimodipine bulk; ● – experimental, nimodipine confined; □ – MD simulations, nimodipine bulk, - MD simulations – nimodipine confined, considering only reorientation of methyl groups;  - MD simulations – nimodipine confined, considering reorientation of methyl groups and fragments of its chains. Error bars are of the size of the symbols.

4. Conclusions We have performed MD simulations to study the behaviour of nimodipine confined in SBA-15 and compare it to that in the crystalline form, aiming to provide a molecular explanation to the differences observed experimentally between both systems by NMR. The main conclusions can be summarized as follows:


19


Page 20 of 26

1. At the pore fillings corresponding to experimental conditions, nimodipine molecules form one incomplete layer at the silica surface and cluster into small groups, each containing only a few molecules. Apart from a few exceptions, the distance of one molecule to the nearest one is significantly larger than in the crystal, so under these conditions nimodipine interacts mainly with the silica surface and not with other nimodipine molecules. 2. Nimodipine molecules are connected with the silica surface by hydrogen bonds, which are formed, broken and rebuilt continuously during the simulation. Due to these strong interactions translational diffusion of confined nimodipine molecules is completely inhibited. 3. Molecular reorientations are also limited due to the formation of these hydrogen bonds, but they are not completely inhibited. MD simulations show that confined nimodipine molecules can still perform partial molecular rotations in a limited angular range, as probed by the reorientation of the 1,4-dihydropyridinium ring. Additionally the isopropyl and methoxyethyl ester chains also exhibit a degree of flexibility that is absent in the crystal. 4. In the crystal, the motion of the five methyl groups is enough to reproduce the T1 values measured by NMR. Using the correlation times derived from the simulation for each methyl group, the computed T1 are in excellent agreement with experiment, as shown in Fig. 13. This supports the validity of the force field and conditions employed in our simulations. 5. For confined nimodipine, the same approach fails. But when the reorientation of the entire molecule in a limited angular range and of the isopropyl and methoxyethyl ester chains around specific axis is included, we recover again a very good agreement between the values of T1 extracted from the MD simulations and the experimental ones. This indicates that such motions are clearly responsible of the differences observed between the crystal and the confined system.

Acknowledgements This research was supported in part by PL-Grid Infrastructure. The work has been partially financed by the National Science Centre of Poland, Grant No. 2015/17/B/ST5/00104. A.P. would like to thank the Institut Laue Langevin (France) for an invitation as a visiting scientist.


20

Page 21 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60



21


Page 22 of 26

References 1. Davis, M.E. Ordered Porous Materials for Emerging Applications, Nature 2002, 417, 813-821. 2. Li, W.; Zhao, D. An Overwiew of the Synthesis of ordered Mesoporous Materia;s, Chem. Commun. 2013, 49, 943-946. 3. Corma, A. From Microporous to Mesoporous Molecular Sieve Materials and Their Use in Catalysis. Chem. Rev. 1997, 97, 2373-2419. 4. Zhao, D.; Feng, J.; Huo Q.; Melosh, N.; Fredrickson, G. H.; Chmelka, B. F.; Stucky, G. D. Triblock Copolymer Syntheses of Mesoporous Silica with Periodic 50 to 300 Angstrom Pores. Science, 1998, 279, 548–552. 5. Xia, Y.; Dosseh, G.; Morineau, D.; Alba-Simionesco, Ch.; Phase Diagram and Glass Transition of Confined Benzene. J. Phys. Chem. B, 2006, 110, 19735–19744. 6. Alba-Simionesco, Ch.; Dosseh, G.; Dumont, E.; Frick, B.; Geil, B.; Morineau, D.; Teboul, V.; Xia, Y. Confinement of Molecular Liquids: Consequences on Thermodynamic, Static and Dynamical Properties of Benzene and Toluene. Eur. Phys. J. E, 2003, 12,19–28. 7. Schonhals, A.; Goering, H.; Schick, Ch.; Frick, B.; Zorn, R. Glassy Dynamics of Polymers Confined to Nanoporous Glasses Revealed by Relaxational and Scattering Experiments. Eur. Phys. J. E, 2003, 12, 173–178. 8. Frick, B.; Alba-Simionesco, Ch.; Dosseh, G.; Quellec, C. Le.; Moreno, A.J.; Colmenero, J.; Schonhals ,A.; Zorn, R.; Chrissopoulou, K.; Anastasiadis, S.H.; Dalnoki-Veress, K. Inelastic Neutron Scattering for Investigating the Dynamics of Confined Glass-Forming Liquids. J. NonCrystal. Solids, 2005, 351, 2657 – 2667. 9. Findenegg ,G. H.; Jahnert, S.; Akcakayiran, D.; Schreiber, A. Freezing and Melting of Water Confined in Silica Nanopores. ChemPhysChem, 2008, 9, 2651–2659. 10. Gallo, P.; Rovere, M.; Chen, S.-H. Dynamic Crossover in Supercooled Confined Water: Understanding Bulk Properties Through Confinement. J. Phys. Chem. Lett., 2010, 1, 729–733. 11. de Souza, N. R.; Kolesnikov, A. I.; Burnham, C. J.; Loong, C-K. Structure and Dynamics of Water Confined in Single-Wall Carbon Nanotubes. J. Phys.: Condens. Matter, 2006, 18, S2321. 12. Mallamace, F.; Broccio, M.; Corsaro, C.; Faraone, A.; Liu, A.; Mou, Ch-Y.; Chen, Ch.-Y.; Dynamical Properties of Confined Supercooled Water: an NMR Study. J. Phys.: Condens. Matter, 2006, 18, S2285.


22

Page 23 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


13. Vallet-Regi, M.; Ramila ,A.; del Real, R. P.; Perez-Pariente, J. A New Property of MCM-41: Drug Delivery System. Chem. Mater., 2000, 13, 308–311. 14. Manzano, M.; Vallet-Regi ,M. New Developments in Ordered Mesoporous Materials for Drug Delivery. J. Mater. Chem., 2010, 20, 5593–5604. 15. A. Paudel, M. Geppi, G. van de Mooter, J. Pharm. Sci. 2014, 105, 2635-2662. 16. Azais, T.; Tourné-Péteilh, C.; Aussenac, F.; Baccile, N.; Coelho, C.; Devoisselle ,J.M.; Babonneau, F. Chem. Mater. 2006, 18, 6382-6390. 16. Yang ,X.; Ong, T.C.; Michaelis, V.K.; Heng, S.; Huang, J.; Griffin, R.G.; Myerson, A.S. CrystEngComm 2014, 16, 9345-9352. 18. Mäkilä, E.; Ferreira, M.P.A.; Kivelä ,H.; Niemi, S.M.; Correia, A.; Shahbazi, M.A.; Kauppila, J.; Hirvonen, J.; Santos, H.A.; Salonen, J. Langmuir 2014, 30, 2196-2205. 19. Kiwilsza, A.; Milanowski, B.; Drużbicki, K.; Jenczyk, J.; Jarek, M.; Mielcarek, J.; Lulek, J.; Pajzderska, A.; Wasicki, J. Molecular Dynamics and the Dissolution Rate of Nifedipine Encapsulated in Mesoporous Silica. Microp. Mesop. Materials, 2017, 250, 186-194. 20. Kiwilsza, A.; Milanowski, B.; Drużbicki, K.; Emmerson ,C.; Grzeszkowiak, M.; Jarek, M.; Mielcarek, J.; Lulek, J.; Pajzderska, A.; Wasicki, J. Mesoporous Drug Carrier Systems for Enhanced Delivery Rate of Poorly Water-Soluble Drug: Nimodipine. J. Porous Mater. 2015, 22, 817-829. 21. Inzitari, D.; Poggesi ,A. Calcium Channel Blockers and Stroke. Aging Clin. Exp. Res.

2005, 17, 16–30. 22. St-Onge, M.; Dube, P.A.; Gosselin, S.; Guimont, C.; Godwin, J.; Archambault, P.M.; Chauny, J.M.; Frenette, A.J.; Darveau, M.: Le Sage, N.; Poitras, J., Provencher, J., Juurlink, D.N.; Blais, R. Treatment for Calcium Channel Blocker Poisoning: a Systematic Review. Clin. Toxicol. 2014, 52, 926–944. 23. Kiwilsza, A.; Pajzderska, A.; Mielcarek, J.; Jenczyk, J.; Wasicki, J. Dynamical Properties of Nimodipine Molecules Confined in SBA-15 matrix, Chem. Phys. 2016, 475, 126-130. 24. Todorov, T.; Smith, W.; Trachenko, K.; Dove, M.T. Dl_POLY_3: New Dimensions in Molecular Dynamics Simulations via Massive Parallelism. J. Mater. Chem. 2006, 16, 1911−1918.


23


Page 24 of 26

25. Langs, D.A.; Strong, P.D.; Triggle, D.J. Receptor Model for the Molecular Basis of Tissue Selectivity of 1,4-Dihydropyridine Calcium Channel Drugs. J.Comput.-Aided Mol. Des. 1990, 4, 215-230. 26. Pajzderska, A.; Gonzalez, M. A.; Mielcarek, J. Wasicki, J. Water Behaviour in MCM-41 as a Function of Pore Filling and Temperature Studied by NMR and Molecular Dynamics Simulations, J. Phys. Chem. C, 2014, 118, 23701 – 23710. 27. Dassault Systèmes BIOVIA, Materials Studio, Release 2016, San Diego: Dassault Systèmes, 2016. 28. Berendsen, H. J. C.; Postma, J. P. M.; Gunsteren, W. F.; DiNola, A.; Haak, J. R. Molecular Dynamics with Coupling to an External Bath. J. Chem. Phys. 1984, 81, 3684−3690. 29. Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. Development and Testing of a General AMBER Force Field”. J. Comput. Chem, 2004, 25, 1157-1174. 30. Brodka, A.; Zerda, T. W. Molecular Dynamics of SF6 in Porous Silica. The Journal of Chemical Physics, 1991, 95, 3710–3718. 31. Milischuk, A. A.; Krewald, V.; Ladanyi, B. Water dynamics in silica nanopores: The selfintermediate function, J. Chem. Phys. 2012, 136, 224704. 32. Lerbret, A.; Lelong, G.; Mason, P. E.; Saboungi, M.-L.; Brady, J. W. Molecular Dynamics and Neutron Scattering Study of Glucose Solutions Confined in MCM-41, J. Phys. Chem. B,

2011, 115, 910-918. 33. Hinsen, K.; Pellegrini, E.; Stachura, S.; Kneller, G. R. nMoldyn 3: Using Task Farming for a Parallel Spectroscopy-Oriented Analysis of Molecular Dynamics Simulations. J. Comput. Chem.

2012, 33, 2043−2048. 34. Durrant, J. D.; McCammon, J. A. HBonanza: A Computer Algorithm for MolecularDynamics-Trajectory Hydrogen-Bond Analysis, J. Mol. Graphics Modell., 2011, 31, 5–9. 35. Desiraju, G.R.; Steiner ,Th. The weak hydrogen bond in structural chemistry and biology, Oxford University Press, 2001. 36. Delle Piane, M.; Corno, M.; Ugliengo, P. Propionic acid derivatives confined in mesoporous silica: monomeric or dimers? The case of ibuprofen investigated by static and dynamic ab abitio simulations. Theor Chem Acc 2016, 135, 53. 37. Drużbicki, K.; Pajzderska, A.; Kiwilsza, A.; Jenczyk, J.; Chudoba, D.; Jarek, M.; Mielcarek, J.; Wasicki, J. In Search of the Mutual Relationship Between the Structure, Solid-State


24

Page 25 of 26 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Spectroscopy and Molecular Dynamics in Selected Calcium Channel Blockers Eur. J. Pharm. Scien., 2016, 85, 68-83. 38. J.-P. Korb, S.Xu, F. Cros, L. Malier, J. Jonas, Quenched molecular reorientation and angular momentum for liquids confined to nanopores of silica glasses, J. Chem. Phys. 1997, 107, 40444050. 39. Korb, J.-P.; Xu ,S.; Jonas, J. Confinement effects on dipolar relaxation by translational dynamice of liquids in porous silica glasses, J. Chem. Phys. 1993, 98, 2411-2422. 40. Descressain ,R.; Mansare, T.; Cors, C.; Cochin, E.; Carpentier, L. Dynamic 1H NMR Investigation of the Liquid Crystal MBBA Confined in Porous Silica, Mol. Cryst. Liq. Cryst.,

2006, 460, 43-61. 41. Fatkullin, N.; Kimmich, R.; Fischer, E.; Mattea, C.; Beginn, U.; Kroutieva, M. The Confinedto-Bulk Dynamic Transition of Polymer Melts in Nanoscopic Pores of Solid Matrices with Varying Pore Diamater, New J. Phys. 2004, 6, 46. 42. Abragam, A. The principles of nuclear magnetism; Oxford University Press: Amen House London, 1961.


25


Page 26 of 26

TOC Graphic


26

Molecular Dynamics Simulations of Nimodipine ... - ACS Publications

Recommend Documents