Day 9 Mesoscale: DPD Tutorial Exercise 2: Applying DPD to molecular systems

Following on from the first tutorial, we are now going to extend our DPD calculations to include additional interactions between some of our particles and form mesoscopic representations of molecules. This vastly extends the range of systems DPD can model, which is well-suited to finding larger scale structures relatively quickly (certainly compared with atomistic MD).

The principal learning objectives for this tutorial are:

  1. to demonstrate the simplicity of extending DPD simulations to molecular systems,

  2. to show DPD can simulate the formation of large-scale structures,

  3. to see the effect molecular topology (e.g. bond straightness) has on those structures, and

  4. to illustrate how you can determine mesophases using order parameters.

You will be working exclusively with DL_MESO_DPD in this tutorial and the exercises in Jupyter notebooks can be run in your work environment (STFC Cloud Environment), although the second exercise will benefit from using SCARF to accelerate its DPD calculations: please make sure you have setup access to SCARF before starting this one.

As you should have already done with the first tutorial, access the working directory for this exercise:

~/WORKSHOP/Day_9Meso

and then open the notebook Day9DPDTutorial2Ex1.ipynb, and follow the instructions within it for the first exercise. The second exercise will follow with the notebook Day9DPDTutorial2Ex2.ipynb. While these notebooks are fully contained, we include some additional background information here for you to look at while DL_MESO_DPD is running. (The simulations for each exercise will take about half an hour to run.)

Ex1. DPD mesophases

The simplest molecular structure we can create in a DPD simulation is a dimer, a molecule consisting of two beads joined together with a bond (e.g. a harmonic spring). If the dimer is amphiphilic, the two beads interact differently with other beads representing a solvent; the hydrophobic bead repels solvent beads more strongly, while the hydrophilic bead has a stronger affinity for the solvent (or repels it less strongly).

The concentration of amphiphilic molecules in solution (and temperature) have an effect on the structures that the molecules can form. Our dimers can produce one of three distinctive phases: isotropic phases with spherical ‘blobs’ of material, hexagonal phases of long tubes laid out in a hexagonal pattern, and lamellar phases of parallel ‘sheets’ (lamellae).

image0

image1

image2

Isotropic \((L_{1})\) phase

Hexagonal \((H_{1})\) phase

Lamellar \((L_{\alpha})\) phase

While it might be possible to see these patterns using standard trajectories, it is often easier to make them out by rendering particular parts of the molecules (e.g. the hydrophobic tails) as isovolumes on a grid. These density maps can then be used to calculate order parameters to help identify particular mesophases without needing to visualise the simulations, which makes it possible to automate the process of finding phase diagrams (e.g. discovering concentrations and temperatures where one phase transforms into another).

Ex2. Lipid bilayers, micelles and vesicles

Following on from amphiphilic dimers, we can construct longer, more realistic amphiphilic molecules - often with smaller hydrophilic head groups and longer hydrophobic tails - such as surfactants and lipids. Lipids are naturally occurring molecules (e.g. fats, waxes, sterols, glycerides, fat-soluble vitamins) that are often either hydrophobic or amphiphilic and can form mesoscopic structures.

We know from the first exercise that DPD can quickly find large-scale mesoscopic structures, which are functions of concentration, temperature, molecular topology etc. In terms of topology, we can add further interactions between the beads: one example includes controlling the angles between adjacent bonds in each molecule.

In this exercise, we will use a DPD model for a lipid molecule consisting of a hydrophilic bead and six hydrophobic beads joined together with harmonic bonds. We will also include a cosine bond angle potential:

(1)\[U_{ijk} = A\lbrack 1 + \cos{\left( m\theta_{ijk} - \theta_{0} \right)\rbrack\ }\]

that will control the angle between pairs of bonds. The angles will have an effect on the structures that form out of the lipids, which include:

  • micelles - spherical structures with hydrophobic tails pointing inwards

  • bilayers - two planar layers of lipids with hydrophobic tails pointing towards each other

  • vesicles (liposomes) - spherical form of bilayer with hydrophilic head groups and solvent inside