Molecular dynamics simulations of rupture in lipid bilayers

Abstract

Magnetic fluid hyperthermia is a promising cancer therapy in which magnetic nanoparticles are acted upon by a high-frequency oscillating magnetic field. While the accepted mechanism is localized hyperthermia, it is plausible that shear stresses due to nanoparticles rotating near a cell membrane may induce rupture, enhancing the effectiveness of the treatment. With the goal of understanding this further, molecular dynamics simulations were carried out on a model cell membrane. A bilayer composed of dipalmitoylphosphatidylcholine lipids was subjected to an incremental tension as well as an incremental shear stress. In both cases, it was found that the bilayer could withstand a surface tension of approximately 90 mN/m prior to rupture. Under tension, the bilayer ruptured at double its initial area, whereas under shear, the bilayer ruptured at 1.8 times its initial area. The results show that both incremental tension and incremental shearing are able to produce bilayer rupture, with shear being more injurious, yielding a larger surface tension for a smaller deformation. This information allows for comparison between the estimated energy required to rupture a cell membrane and the energy that a magnetic nanoparticle would be able to generate while rotating in a cellular environment. Our estimates indicate that magnetically blocked nanoparticles with diameters larger than 50 nm may result in rupture due to shear.

Keywords

molecular dynamics simulation lipid bilayer shear stress magnetic fluid hyperthermia

Introduction

Current treatments for cancer are either highly invasive, such as surgery, or very taxing on the patient, as in radiation and chemotherapy. In some cases, the side-effects associated with treatment result in patients foregoing such treatment. Also, for many types of cancer, established treatments have very low success rates. These limitations stimulate the development of alternative treatments, and recently, the application of nanotechnology and nanoparticles in the early stage detection and treatment of cancer.¹ One promising alternative is magnetic fluid hyperthermia (MFH), in which magnetic nanoparticles are delivered to a cancer tumor and a high-frequency oscillating magnetic field is applied. When subjected to an oscillating magnetic field, the nanoparticles' non-equilibrium response results in energy dissipation and a localized rise in temperature, usually to a target temperature above 46°C. Because cancer cells are more susceptible to heat damage than healthy cells, the treatment has been shown to result in tumor suppression. Jordan et al. ^2,3 examined iron oxide nanoparticles functionalized with either dextran or silane for the treatment of breast cancer cells in vitro. When the cells were subjected to a high-frequency alternating magnetic field, a local temperature increase was observed. It was determined that maintaining a tumor temperature of 47°C resulted in cell death. A similar study was conducted by Johannsen et al. ⁴ to determine the effect of rotating nanoparticles in vivo to treat prostate cancer in mice. Although the authors were able to shrink the tumors, substantial peripheral heating of healthy tissue surrounding the tumor was also observed. Using the heat generated from magnetic nanoparticles to kill cells resulted in the authors coining the term MFH.

Although energy dissipation is widely accepted as the dominant mechanism leading to cell death, it is possible that shear due to rotation of nanoparticles near the cell membrane may lead to membrane rupture, especially for larger nanoparticles with magnetic dipoles that are locked into the nanoparticle crystal structure (i.e. magnetically blocked nanoparticles). There are reports in the literature which indicate that a rise in temperature is not required to induce cell death. In a study conducted by Halbreich et al. ⁵ magnetic nanoparticles in an oscillating magnetic field were used to damage mouse liver cells in vivo. Cell death was found to occur in the absence of any observable temperature increase, leading to the possibility that mechanical means were involved. In this paper, we aim to estimate the magnitude of the energy required to induce membrane rupture due to shear stresses, and gain understanding of how nanoparticle physical parameters may be engineered to improve the effectiveness of MFH. Molecular dynamics (MD) simulations are performed in the absence of heating and have the ultimate goal of predicting whether a rotating particle interacting with the cell membrane could induce membrane rupture.

As a first approximation, it is possible to model the cell membrane as a lipid bilayer in water. This system has been previously studied using MD simulations, and researchers have been able to reproduce experimental results for many equilibrium thermodynamic properties such as density, area and volume per lipid, bilayer repeat distance, etc. (see Tieleman et al. ⁶ and Feller⁷ for reviews). MD has also been used to study pore formation and stability of lipid bilayers subject to stress.^8,9 Tieleman et al. ¹⁰ induced pore formation using an applied mechanical stress. They found that when tension was applied in the form of a lateral pressure >−200 bar, a pore would form in the bilayer. This pore would grow rapidly in size until the bilayer became destabilized. Leontiadou et al. ¹¹ performed a similar study. They began with a bilayer containing a pre-existing pore and subjected the bilayer to a tension in the form of a lateral pressure. They found that at varying levels of stress, up to a surface tension of 38 mN/m, the pores were stable at a particular size that was correlated to the applied stress. If the applied stress was >38 mN/m, the pores would grow indefinitely, resulting in bilayer rupture. The speed of rupture depended on the magnitude of the applied stress. Additionally, in comparison to a bilayer without a preformed pore, Leontiadou et al. ¹¹ found results consistent with Tieleman et al., in that a very high stress of approximately 90 mN/m (−200 bar lateral pressure) is required to form an unstable pore. This led to the conclusion that pore formation and expansion, and hence bilayer strength, is dependent on the loading rate of the stress. This is consistent with the theory of Evans and Heinrich,¹² who hypothesized that at low loading rates, pore expansion becomes rate limiting, whereas at high loading rates, pore formation is rate limiting. Groot and Rabone¹³ took a more in-depth look at the rupture of a lipid bilayer. Using dissipative particle dynamics (DPD)¹⁴ they looked at bilayer rupture in the presence of surfactants. To do so, they subjected the bilayer to normal strain increments while monitoring the surface tension of the bilayer. When the surface tension reached a peak and subsequently sharply decreased, this was taken as signifying rupture of the bilayer. Their simulations showed that a lipid bilayer could withstand a surface tension of 67 mN/m prior to rupture.

Similar systems have also been studied experimentally. The most common procedure for examining properties such as area expansion, area compressibility modulus and rupture tension is through the micropipette aspiration technique.^15–19 Here, a small portion of a vesicle is sucked through the tip of a micropipette under a defined pressure drop, and the vesicle deformation and rupture tension can be measured using a high-speed camera. Generally, lipid vesicles rupture following an area expansion of ∼5% and an external tension on the order of 10–30 mN/m. Even very tough bilayers made from high molecular weight diblock co-polymers cannot sustain more than 20% area expansion prior to rupture.²⁰ However, values for the rupture tension and area expansion are both dependent on the loading rate.^15–19

None of the aforementioned studies focused on the rupture of lipid bilayers under shear, which is what may occur in MFH. The existing literature pertaining to simulations of lipid bilayers under shear is sparse. Blood et al. ²¹ simulated the response of a solvent shear flow producing a force on the bilayer of 2.9×10⁻³ pN, finding that the flow increased the order of the lipids in the bilayer. Shkulipa et al. ²² used Lees–Edwards sliding boundary conditions²³ to obtain the surface shear viscosity and intermolecular friction of amphiphilic bilayers containing lipids of varying tail lengths. To the best of our knowledge, simulations of lipid bilayer rupture under shear stress are non-existent. In addition, given that intermolecular interactions of a bilayer under tension are different from those occurring under shear, it is not a straightforward expectation that a bilayer would have the same strength for rupture in these two cases.

The paper is organized as follows: the next section describes the system we are using and the methodology employed to simulate it. The following subsections describe the simulation methodology of a bilayer under tension and explains the methodology for subjecting a bilayer to a shear stress. The subsequent section gives the results for each method and compares the two, discussing the results. And finally we have the conclusions.

Materials and methods

All MD simulations are performed with the software package GROMACS, version 3.3.1.²⁴ The force field used was developed by Berger et al. ²⁵ and is based on the united-atom version of the optimized potential for liquid simulations force field.²⁶ In this approach, all the CH_x groups (with x = 1, 2 or 3) are modeled as single units, while other heavy atoms are modeled as exact atoms. For water, the Simple Point Charge (SPC) model of Berendsen et al. ²⁷ is used. We chose a cut-off for Lennard–Jones interactions of 1.4 nm and electrostatics were treated using a 1.0 nm cut-off and the particle mesh Ewald algorithm²⁸ for long-range interactions.

The cell membrane in our simulations is modeled as a lipid bilayer, taken from the website of Dr P Tieleman,²⁹ consisting of 128 dipalmitoylphosphatidylcholine (DPPC) lipids (64 lipids/leaflet) and 3655 SPC water molecules, yielding a fully solvated bilayer. All simulations began by performing an energy minimization on the bilayer, followed by a 5 ns equilibration period in an NPT (fixed pressure P, temperature T and number of atoms N) ensemble (Berendsen temperature coupling at 310°K and isotropic Berendsen pressure coupling at 1 bar³⁰). Figure 1 shows the lipid bilayer following energy minimization and equilibration with a resulting simulation cell size of L _x × L _y × L _z = (6.4 nm)³, the initial configuration for all rupture runs.

Figure 1

Dipalmitoylphosphatidylcholine (DPPC) lipid bilayer following energy minimization and 5 ns of equilibration time. (a) Side view of the bilayer. (b) Top-down view of the bilayer, with water molecules omitted for clarity

Assessment of bilayer rupture is performed in three ways. First, we perform direct examination of the simulation snapshots for a qualitative view of the rupture. Second, we compute variations in the surface tension as a function of bilayer area. Finally, we also monitor the diffusion coefficient of water molecules through the bilayer.

The surface tension of a bilayer with a surface normal in the z-direction can be found from the pressure tensor P through the relation³¹:

with γ being the surface tension, P _N the pressure normal to the bilayer and P _T the tangential pressure. In an MD simulation, this equation becomes³²

where L _Z is the size of the simulation cell in the z-direction and P _i is the pressure in the i = x-, y-, and z-direction, respectively. We compute the surface tension as a function of bilayer area, noting that a sharp decrease in surface tension indicates that the lipid bilayer has ruptured.

To determine the diffusion coefficient, we calculate the mean-squared displacement (MSD) of the water molecules in the z-direction during each relaxation period. A least-squares fit is used to fit the MSD to a straight line, and the Einstein relation can then be used to calculate the diffusion coefficient D,

where n is the number of dimensions, t is time and r _i is the position of molecule i.

Bilayer rupture by incremental tension

Following equilibration, incremental stretches are performed on the bilayer. The x- and y-coordinates of all atoms are incremented by a stretch factor of λ, while the z-coordinates are divided by a stretch factor of λ ². This transformation conserves the total volume of the bilayer. A stretch is applied instantaneously and followed by an NP_NAT relaxation period, which keeps the area of the bilayer fixed and a constant normal pressure of 1 bar. Performing the simulation in this manner allows the bilayer to respond to the applied stress, maintaining quasi-equilibrium. We have analyzed three different values of λ: λ = 1.03, 1.04 and 1.05. For every value of λ, each stretch is followed by a 3 ns relaxation period and an average measure of the surface tension of the bilayer is calculated over this relaxation period to obtain an estimate of the surface tension associated with that particular value of the bilayer area.

Bilayer rupture by incremental shearing

The rotation of a nanoparticle responding to a magnetic field interacting in a cellular environment would likely result in a shear stress imparted on the cell surface. In our simulations, an acceleration is applied in the x-direction over the whole simulation box with the form

where A is the amplitude of the acceleration, z is the z-coordinate of a particular atom and L _z is the height of the simulation box in the z-direction. In this way, shear acceleration is applied at the bilayer surface and smoothly shifted to zero at the center, as one moves away from the surface. The acceleration is applied for 10 ps in an NPT ensemble followed by a 3 ns relaxation period in an NP_NAT ensemble with the same parameters as described above. The chosen value A = 0.40 nm/ps² for the amplitude of the acceleration in each increment of shear yields a similar area increase as obtained in each increment of the tension simulations described above.

Results

We have looked at bilayer rupture using two methods, an in-plane stretch of the bilayer as in reference,¹³ and shearing the surface of the bilayer, a process more akin to the forces that a rotating nanoparticle would likely create.

Bilayer rupture by incremental tension

Figure 2 presents the response of the lipid bilayer surface tension to the imposed strain as a function of time for the case λ = 1.03. The value plotted is the running average of the surface tension over each strain increment. The instantaneous surface tension directly following each stretch is quite large; hence, the axes of Figure 2 have been truncated in order to show the physically relevant region in which the bilayer surface tension has reached its equilibrium value after each stretch. The equilibrium value is used in Figure 3, which shows the surface tension of the bilayer as a function of the normalized area of the bilayer.

Figure 2

Running average of the surface tension following each bilayer stretch for λ = 1.03, corresponding to an area increase of 6% (A color version of this figure is available in the online journal)

Figure 3

Surface tension versus stretch area for a dipalmitoylphosphatidylcholine (DPPC) bilayer stretched in its principal directions for stretch values λ = 1.03, 1.04 and 1.05 (A color version of this figure is available in the online journal)

Figure 3 shows the change in surface tension of the bilayer with respect to the area of the membrane for each magnitude of strain examined. The abscissa is normalized to the area of an equilibrated bilayer, taken as the result of the equilibration run. As the area of the bilayer increases, the surface tension also increases up to a peak value. Further stretching of the bilayer beyond the peak value results in a decrease in surface tension. Interestingly, the bilayer ruptured at roughly the same tension, independent of the magnitude of the strain. For each strain value, the bilayer could be stretched to approximately double its initial area, with a corresponding surface tension of approximately 90 mN/m.

Figure 4 shows typical snapshots of the evolution of bilayer structure during a simulation. For the first several increments, stretching the bilayer does not seem to induce any noticeable pore formation. With subsequent stretches, small pores in the bilayer begin to develop. These pores are transient and reseal within tens to hundreds of picoseconds. Eventually, the pores in the bilayer become very large, as indicated in Figure 4c. Pores of this size are stable throughout the entire relaxation period of 3 ns and were observed under conditions that led to bilayer rupture.

Figure 4

Representative snapshot of pore formation in the bilayer during incremental strain with λ = 1.04. (a) Initial configuration of the bilayer, (b) following five stretches (15 ns) and (c) following 9 stretches (27 ns). In all pictures, water molecules are omitted for clarity

Figure 5 shows the transverse diffusion coefficient of water molecules as a function of bilayer area. When stretched to double the initial area, the diffusion coefficient perpendicular to the bilayer increased by roughly an order of magnitude, indicating water molecules flowing through the pores in the bilayer, and thus rupture. Initially, the diffusion of water perpendicular to the bilayer is quite small, ∼3 × 10⁻⁷ cm²/s, compared with the diffusion coefficient of bulk SPC water, 2.3 × 10⁻⁵ cm²/s.²⁷ The small but noticeable downward trend present in the figure is a simulation artifact resulting from periodic boundary conditions and the progressive shrinking of the water layer as the bilayer area expands.

Figure 5

Diffusion coefficient perpendicular to the bilayer surface as a function of normalized bilayer area (A color version of this figure is available in the online journal)

Our results are consistent with those found in reference,¹³ in that the bilayer can be stretched to double its initial area before stable pores begin to develop. This is not entirely unexpected because even though different approaches were used, the volumes of each lipid molecule are roughly the same. Stretching the bilayer causes an increase in the average distance between lipid molecules, leading to an increase in the potential energy of the system. Eventually, the increase in intermolecular lipid distance becomes large enough that it is energetically favorable to form a stable pore. The value for the bilayer in our system to rupture is slightly higher than the one found in reference,¹³ largely due to the fact that we are working with different lipids. In addition, the authors in reference¹³ used DPD, while our simulations were carried out with MD. Systems studied with MD are generally smaller and inherently more stable, and will therefore require larger forces for rupture.

Bilayer rupture by incremental shearing

We have looked at rupture of a lipid bilayer under strain to determine the maximum surface tension it can withstand in order to obtain an estimate of the energy that a rotating nanoparticle would need to produce to rupture a cell. However, a nanoparticle rotating in a cellular milieu would likely induce a shear force, rather than a tension on the lipid bilayer. Thus, it is of interest to study the response of a lipid bilayer to a shear stress.

Figure 6 shows the running averages of the surface tension for each shear–relax iteration. In total, 20 shear–relax iterations were performed. The plot shows an oscillation in the surface tension during the shearing phase followed by a relaxation period of 3 ns, where the surface tension decays. The equilibrium value of the surface tension is taken as the final average value during the relaxation period prior to a subsequent shear increment. Thus, each of these averages corresponds to a point value in surface tension in Figure 7. Initial trials were run with different values of the relaxation time, and it was determined that 3 ns was sufficient for the surface tension to decay to a constant value following a particular shear increment.

Figure 6

Running average of the surface tension as a function of time for each shear increment (A color version of this figure is available in the online journal)

Figure 7

Surface tension versus area for a dipalmitoylphosphatidylcholine (DPPC) bilayer subjected to incremental shearing (A color version of this figure is available in the online journal)

Figure 7 shows a plot of the surface tension of the bilayer versus its area. As can be seen, the surface tension increases with increasing area until approximately 84 mN/m, whereupon it begins to sharply decrease down to 52 mN/m. This decrease is indicative of bilayer rupture. At this point, a large number of water molecules are able to flow through the pores in the bilayer, neutralizing the tension that was imparted by the shearing. The surface tension does not decrease to its initial value at the start of the simulation due to the structure that is still present in the bilayer.

Figure 8 shows the evolution in structure of the bilayer during the shearing simulations. Figure 8a shows the initial configuration of the bilayer prior to any external stimulus. Figure 8b shows the bilayer following 10 shear increments. Evident in the snapshots is that the bilayer is still much intact, although at times one does see small transient pores starting to develop during the 3 ns relaxation period. Figure 8c shows the bilayer following 20 shear increments, with a large pore that is stable over the entire 3 ns relaxation period. This provides visual confirmation that bilayer rupture has in fact occurred. The transverse water diffusion coefficient is plotted in Figure 9, obtaining similar results to the strain simulations.

Figure 8

Representative snapshot of pore formation in the bilayer during incremental shearing. (a) Initial configuration of the bilayer following equilibration. (b) The bilayer following 10 shear increments (30.1 ns). (c) The bilayer following 20 shear increments (60.2 ns). The surface tension at this point is 52 dyn/cm. In all pictures, water molecules are omitted for clarity

Figure 9

Diffusion coefficient perpendicular to the bilayer surface as a function of normalized bilayer area for the case of a bilayer under shear (A color version of this figure is available in the online journal)

We have found that we can produce bilayer rupture through incremental shearing as evidenced by the surface tension and water diffusion measurements. In comparing our shearing simulations to the incremental stretching simulations, it is seen that the two different methods of producing rupture do not behave in the same way. The strain and shearing simulations showed bilayer rupture at slightly similar values of the surface tension, 84 and 90 mN/m, respectively; however, the bilayer under shear ruptured at approximately 1.8 times its initial area, lower than that for the strain simulation. This suggests that a shear stress is more injurious to a lipid bilayer than in-plane stretching, i.e. lipid bilayers rupturing under shear sooner than in tension indicates that the rotating magnetic nanoparticles will require less energy to produce cell membrane rupture.

Comparison with experimental results

The small system size and fast loading rate present in our simulations have a large influence on the system. In our simulations we are examining a very small portion of a lipid bilayer, not an entire vesicle or cell. In a laboratory setting, the area expansion and rupture tension are generally measured using a micropipette aspiration technique.^15–19 In these studies the area expansion is rarely >5%. However, this value is relative to the entire vesicle rather than the small portion in contact with the micropipette and is dependent on loading rate.¹⁸

Additionally, the micropipette aspiration studies^15–19 report the surface tension at which rupture occurs in the order of 5–10 times lower than the value reported here. This is also as a result of the applied loading rate in our system, which is 7–9 orders of magnitude higher than in the experiments. In fact, using the theory of Evans et al. ¹⁷ for the dependence of rupture tension on loading rate, we can estimate a theoretical rupture tension for our applied loading rate (1.4 mN/m/ns). For the high loading rate regimen, Evans et al. derive the relation:

where σ is the tension at rupture, σ _δ is a tension scale for defect formation, ℜ _σ is the loading rate and ν _0δ is the frequency for spontaneous nucleation of defects. Using the values of Evans et al. for DOPC (similar to DPPC), where σ _δ = 4.0 mN/m and ν _0δ = 0.22 s⁻¹ and a loading rate of 1.4 mN/m/ns, the calculated rupture tension is 85 mN/m, consistent with the 84 mN/m determined from our simulations. Thus, although in a laboratory setting the area expansion and rupture tension are lower than reported here, these values are highly dependent on the loading rate and the portion of the bilayer being examined. Our results are consistent with the experiments, confirming that the value for the surface tension obtained in our simulations is very close to 85 mN/m for a 1.4 mN/m/ns loading rate.

Estimate of the energy required to rupture a cell membrane

Using the results of our simulations, we can estimate the energy required to produce rupture in a lipid bilayer. The surface tension, γ, is defined as the amount of energy needed to increase the area of a surface by a value δ _A. Thus, multiplying the surface tension at which rupture occurs by the normalized bilayer area increase results in the energy needed for rupture. In our shearing simulations, δγ = 81.2 mN/m and δ _A = 32.2 nm² yielding an energy for rupture of 6.3 × 10⁻² J/m². In our simulations, this energy can be thought of as an upper bound. Our system consists of a bilayer film in a simulation cell with cross-sectional area of ∼40 nm². As such it is unable to undergo thermal fluctuations and undulations with wavelength larger than the length of the simulation cell. A much larger 3D bilayer on the other hand undergoes periodic long wavelength undulations and fluctuations that could enhance the formation of pore–nucleation sites. In our stretching simulations with λ = 1.04, the energy for rupture is calculated as 7.8 × 10⁻² J/m². We have performed additional simulations with a bilayer half the size (64 lipids) and confirmed that the force increases with decreasing bilayer size. Our simulations on the smaller bilayer result in an energy for rupture of 1.1 × 10⁻¹ J/m². This is in accordance with the findings of Marrink and Mark,³³ who found that as a result of the size dependence of the area compressibility modulus, the surface tension of stressed bilayers decreases with increasing system size.

Our calculated value of the energy needed to rupture a lipid bilayer under shear can be compared with the energy imparted on a rotating nanoparticle due to an applied oscillating field. From this comparison, an estimate of the size of magnetic nanoparticles needed to induce membrane rupture can be obtained. When a magnetic nanoparticle with a rigidly locked dipole is placed in a magnetic field, the magnetostatic energy is given by³⁴

where μ ₀ is the permeability of free space, equal to 4π × 10⁻⁷ Wb/A/m, M is the domain magnetization, V is the volume and H is the applied magnetic field. Magnetic nanoparticles commonly used in MFH are generally composed of a magnetite core. Magnetite has a domain magnetization of 430 kA/m.³⁵ Johannsen et al. ⁴ worked with nanoparticles composed of magnetite with an average particle core of 15 nm. The particles were coated with ‘an aminosilan-type shell’. The authors used variable field strengths of 0–18 kA/m, with the most common field being 12.6 kA/m. Using Eq. (6) and our value of the energy at 18 kA/m, we obtained a particle diameter sufficient to produce mechanical rupture approximately 80 nm, and at 12.6 kA/m, it is 90 nm, 5–6 times larger than that used in reference 4. These estimates indicate that it is unlikely that cell death observed in reference 4 was augmented by mechanical rupture. Still, our estimates indicate that nanoparticles smaller than 100 nm, a widely accepted upper limit for biomedical applications of nanoparticles, may enhance the effectiveness of MFH by mechanical disruption of the cell membrane.

Conclusions

We have modeled two plausible mechanisms that may lead to cell membrane rupture during MFH, with the goal of estimating the energies needed to produce rupture of the lipid bilayer and relating these to particle size. Our simulations of membrane rupture via incremental shear stress yield a value for the surface tension of the bilayer equal to 80–90 mN/m, with energy for rupture of 6.3 × 10⁻² J/m². These results are consistent with previous results obtained for similar systems in the literature.¹³ Under conditions representative of experimental MFH studies, the estimated energy required to produce rupture obtained from our simulations yielded a required nanoparticle diameter of 80 nm. This is 5–10 times larger than typical particle diameters used in MFH. However, a particle size of 80 nm is plausible for biomedical applications, in which the widely accepted upper limit is between 100 and 200 nm. We also note that these estimates are based on the value of the energy required to rupture a lipid bilayer of approximately 36 nm². This energy can be considered an upper bound, since at the nanoscale size (40 nm²) the bilayer is more stable resulting from reduced thermal fluctuations as compared with a real bilayer; therefore, the force and the energy necessary to rupture real cells in vivo might be lower.

Our next step is to increase the complexity of the system to better model a realistic natural environment. Cell membranes are much more complicated than the simple lipid bilayer system simulated here. For example, cellular membranes in vivo contain a significant amount of cholesterol, which acts to order the lipids resulting in a more solid-like membrane. Similarly, cellular membranes also contain several types of lipids, proteins and surface attached carbohydrates. All of these factors could contribute to the stability of a membrane in response to an applied stress.

Future simulations will involve the addition of ions to the system to simulate the chemical and electrical potential gradient that is present across all cells, for it is the loss of an ion gradient that will ultimately lead to the death of a cell. Following this, it would be beneficial to study the interaction of the magnetic nanoparticle with the cell membrane. As the nanoparticles will likely include a polymeric coating such as dextran, silane, polyethylene glycol, polyvinyl alcohol, etc., simulating the interaction of the polymers coating the nanoparticle with a lipid bilayer would provide useful insights for the study of MFH.

Footnotes

Acknowledgements

The authors would like to thank the NSF NIRT Grant CBET-0609117 and NSF IGERT Nanopharmaceutical Engineering and Science Grant #0504497 for funding this research.

References

Ferrari

. Cancer nanotechnology: opportunities and challenges. Nat Rev Cancer 2005;5:161–71

Jordan

, Scholtz

, Wust

, Fähling

, Felix

. Magnetic fluid hyperthermia (MFH): cancer treatment with AC magnetic field induced excitation of biocompatible superparamagnetic nanoparticles. J Magn Magn Mater 1999;201:413–9

Jordan

, Scholtz

, Wust

, Schirra

, Schiestel

, Schmidt

, Felix

. Endocytosis of dextran and silan-coated magnetite nanoparticles and the effect of intracellular hyperthermia on human mammary carcinoma cells in vitro . J Magn Magn Mater 1999;194:185–96

Johannsen

, Thiesen

, Jordan

, Taymoorian

, Gneveckow

, Waldöfner

, Scholtz

, Kock

, Lein

, Jung

, Loening

. Magnetic fluid hyperthermia (MFH) reduces prostate cancer growth in the orthotopic dunning R3327 rat model. Prostate 2005;64:283–92

Halbreich

, Groman

, Raison

, Bouchaud

, Paturance

. Damage to the protein synthesizing apparatus in mouse liver in vivo by magnetocytolysis in the presence of hepatospecific magnetic nanoparticles. J Magn Magn Mater 2002;248:276–85

Tieleman

, Marrink

, Berendsen

HJC

. A computer perspective of membranes: molecular dynamics studies of lipid bilayer systems. Biochim Biophys Acta 1997;1331:235–70

Feller

. Molecular dynamics simulation of lipid bilayers. Curr Opin Colloid Interface Sci 2000;5:217–23

Wang

, Frenkel

. Pore nucleation in mechanically stretched bilayer membranes. J Chem Phys 2005;122:234711

Tolpekina

, den Otter

, Briels

. Nucleation free energy of pore formation in an amphiphilic bilayer studied by molecular dynamics simulations. J Chem Phys 2004;121:12060

10.

Tieleman

, Leontiadou

, Mark

, Marrink

. Simulation of pore formation in lipid bilayers by mechanical stress and electric fields. J Am Chem Soc 2002;125:6382–3

11.

Leontiadou

, Mark

, Marrink

. Molecular dynamics simulations of hydrophilic pores in lipid bilayers. Biophys J 2004;86:2156–64

12.

Evans

, Heinrich

. Dynamics strength of fluid membranes. C R Phys 2003;4:265–74

13.

Groot

, Rabone

. Mesoscopic simulation of cell membrane damage, morphology change, and rupture by nonionic surfactants. Biophys J 2001;81:725–36

14.

Groot

, Warren

. Dissipative particle dynamics: bridging the gap between atomistic and mesoscopic simulation. J Chem Phys 1997;107:4423–35

15.

Needham

, Nunn

. Elastic deformation and failure of lipid bilayer membranes containing cholesterol. Biophys J 1990;58:997–1009

16.

Olbrich

, Rawicz

, Needham

, Evans

. Water permeability and mechanical strength of polyunsaturated lipid bilayers. Biophys J 2000;79:321–7

17.

Evans

, Heinrich

, Lugwig

, Rawicz

. Dynamic tension spectroscopy and strength of biomembranes. Biophys J 2003;85:2342–50

18.

Heinrich

, Rawicz

. Automated, high-resolution micropipette aspiration reveals new insight into the physical properties of fluid membranes. Langmuir 2005;21:1962–71

19.

Rawicz

, Smith

, McIntosh

, Simon

, Evans

. Elasticity, strength, and water permeability of bilayers that contain raft microdomain-forming lipids. Biophys J 2008;94:4725–36

20.

Discher

, Won

, Ege

, Lee

JCM

, Bates

, Discher

, Hammer

. Polymersomes: tough vesicles made from diblock copolymers. Science 1999;284:1143–6

21.

Blood

, Ayton

, Voth

. Probing the molecular-scale lipid bilayer response to shear flow using nonequilibrium molecular dynamics. J Phys Chem 2005;109:18673–9

22.

Shkulipa

, den Otter

, Briels

. Surface viscosity, diffusion, and intermonolayer friction: simulating sheared amphiphilic bilayers. Biophys J 2005;89:823–9

23.

Allen

, Tildesley

. Computer Simulation of Liquids. Oxford, UK: Oxford University Press, 1987

24.

Lindahl

, Hess

, van der Spoel

. GROMACS 3.0: a package for molecular simulation and trajectory analysis. J Mol Model 2001;7:306–17

25.

Berger

, Edholm

, Jähnig

. Molecular dynamics simulation of a fluid bilayer of dipalmitoylphosphatidylcholine at full hydration, constant pressure, and constant temperature. Biophys J 1997;72:2002–13

26.

Jorgensen

, Madura

, Swenson

. Optimized intermolecular potential functions for liquid hydrocarbons. J Am Chem Soc 1984;106:6638–46

27.

Berendsen

HJC

, Postma

JPM

, van Gunsteren

, Hermans

. Interaction models for water in relation to protein hydration. In: Pullman

, Dordrecht

, eds. Intermolecular Forces. Dordrecht, Holland: Reidel, 1981:331–42

28.

Essmann

, Perera

, Berkowitz

, Darden

, Lee

, Pedersen

LGJ

. A smooth particle mesh Ewald routine. J Chem Phys 1995;103:8577–93

29.

Tieleman

. Biocomputing at the University of Calgary, 2006. See http://moose.bio.ucalgary.ca (last accessed 17 August 2009)

30.

Berendsen

HJC

, Postma

JPM

, van Gunsteren

, Di Nola

, Haak

. Molecular dynamics with coupling to an external bath. J Chem Phys 1984;81:3684–90

31.

Rowlinson

, Windom

. Molecular Theory of Capillarity. Oxford, UK: Oxford University Press, 1982

32.

Chiu

, Clark

, Balaji

, Subramaniam

, Scott

, Jakobsson

. Incorporation of surface tension into molecular dynamics simulation of an interface: a fluid phase lipid bilayer membrane. Biophys J 1995;69:1230–45

33.

Marrink

, Mark

. Effect of undulations on surface tension in simulation bilayers. J Chem Phys 2001;105:6122–7

34.

Jakubovics

. Magnetism and Magnetic Materials. London: Institute of Metals, 1987

35.

Kaneko

, Tokoro

, Nakabayashi

, Adachi

, Fujii

. Magnetic beads with high saturation magnetization. IEEE Trans Magn 2007;43:2424–6