?a???s?as? t?? PowerPoint

1
Simulations of Soft Matter under Equilibrium and
Non-equilibrium Conditions
VAGELIS HARMANDARIS International Conference on
Applied Mathematics Heraklion, 16/09/2013
2
Outline

• Introduction Motivation, Length-Time Scales,
  Simulation Methods.

• Multi-scale Particle Approaches Atomistic and
  systematic coarse-grained simulations of
  polymers.

• Application Equilibrium polymeric systems.

• Application Non-equilibrium (flowing) polymer
  melts.

• Conclusions Open Questions.

3

• COMPLEX SYSTEMS TIME - LENGTH SCALES
• A wide spread of characteristic times (15 20
  orders of magnitude!)
• -- bond vibrations 10-15 sec
• -- dihedral rotations 10-12 sec
  
• -- segmental relaxation 10-9 - 10-12 sec
• -- maximum relaxation time, t1 1 sec (for ? lt
  ?m)?

4
Modeling of Complex Systems Molecular Dynamics

• Classical mechanics solve classical equations of
  motion in phase space, GG(r, p).
• In microcanonical (NVE) ensemble

The evolution of system from time t0 to time t
is given by
Liouville operator 
5
Modeling of Complex Systems Molecular Dynamics

• Various methods for dynamical simulations in
  different ensembles.
• In canonical (NVT) ensemble
• -- Langevin (stochastic) Thermostat
• -- Nose-Hoover thermostat Nosé 1984 Hoover,
  1985 add one more degree of freedom ?.

6
Molecular Interaction Potential (Force Field)
Molecular model Information for the functions
describing the molecular interactions between
atoms.
-- Potential parameters are obtained from more
detailed simulations or fitting to experimental
data.
7
MULTI-SCALE DYNAMIC MODELING OF COMPLEX SYSTEMS
Atomistic MD Simulations Quantitative
predictions of the dynamics in soft matter.
Limits of Atomistic MD Simulations (with usual
computer power) -- Length scale few Å O(10
nm)? -- Time scale few fs - O(1 µs) (10-15
10-6 sec) 107 109 time steps -- Molecular
Length scale (concerning the global
dynamics) up to a few Me for simple
polymers like PE, PB much below Me for more
complicated polymers (like PS)?
Need - Simulations in larger length time
scales. - Application in molecular weights
relevant to polymer processing. - Quantitative
predictions. Proposed method - Coarse-grained
particle models obtained directly from the
chemistry.
8
Systematic Coarse-Graining Overall Procedure
1. Choice of the proper CG description.
-- Microscopic (N particles)
-- Mesoscopic (M super particles)
-- Usually T is a linear operator (number of
particles that correspond to a super-particle
9
Systematic Coarse-Graining Overall Procedure
2. Perform microscopic (atomistic) simulations
of short chains (oligomers) (in vacuum) for short
times.
3. Develop the effective CG force field using the
atomistic data-configurations.
4. CG simulations (MD or MC) with the new
coarse-grained model.
Re-introduction (back-mapping) of the atomistic
detail if needed.
10
Effective (Mesoscopic) CG Interaction Potential
(Force Field)
CG Potential In principle UCG is a function of
all CG degrees of freedom in the system and of
temperature (free energy)

• CG Hamiltonian Renormalization Group Map

• Remember

• Assumption 1

11
Bonded CG Interaction Potential

• Bonded Potential
• Degrees of freedom bond lengths (r), bond
  angles (?), dihedral angles (?)

• Procedure
• From the microscopic simulations we calculate
  the distribution functions of the degrees of
  freedom in the mesoscopic representation,
  PCG(r,?,?).
• PCG(r,?, ?) follow a Boltzmann distribution

• Assumption 2

• Finally

12
Non-bonded CG Interaction Potential Reversible
Work

• Assumption 3 Pair-wise additivity

• Reversible work method McCoy and Curro,
  Macromolecules, 31, 9362 (1998)
• By calculating the reversible work (potential of
  mean force) between the centers of mass of two
  isolated molecules as a function of distance

• Average lt gt over all degrees of freedom G that
  are integrated out (here orientational ) keeping
  the two center-of-masses fixed at distance r.

13
CG MD DEVELOPMENT OF CG MODELS DIRECTLY FROM THE
CHEMISTRY APPLICATION POLYSTYRENE
(PS)? Harmandaris, et al. Macromolecules, 39,
6708 (2006) Macromol. Chem. Phys. 208, 2109
(2007) Macromolecules 40, 7026 (2007) Fritz et
al. Macromolecules 42, 7579 (2009)
1) CHOICE OF THE PROPER COARSE-GRAINED MODEL

• 21 model Each chemical repeat unit replaced by
  two CG spherical beads (PS 16 atoms or 8 united
  atoms replaced by 2 beads).

• CG operator T from CHx to A and B
  description.
• Each CG bead corresponds to O(10) atoms.

s? 4.25 Å sB 4.80 Å

• Chain tacticity is described through the
  effective bonded potentials.
• Relatively easy to re-introduce atomistic detail
  if needed.

2) ATOMISTIC SIMULATIONS OF ISOLATED PS RANDOM
WALKS
14
CG MD Simulations Structure in the Atomistic
Level after Re-introducing the Atomistic Detail
in CG Configurations.

• Simulation data atomistic configurations of
  polystyrene obtained by reinserting atomistic
  detail in the CG ones.
• Wide-angle X-Ray diffraction measurements
  Londono et al., J. Polym. Sc. B, 1996.

grem total g(r) excluding correlations between
first and second neighbors.
15
CG Polymer Dynamics is Faster than the Real
Dynamics
PS, 1kDa, T463K
Free Energy Landscape

• -- CG effective interactions are softer than the
  real-atomistic ones due to lost degrees of
  freedom (lost forces).
• This results into a smoother energy landscape.
• CG MD We do not include friction forces.

16
CG Polymer Dynamics Quantitative Predictions
CG dynamics is faster than the real dynamics.
Time Mapping (semi-empirical) method

• Find the proper time in the CG description by
  moving the raw data in time. Choose a reference
  system. Scaling parameter, tx, corresponds to the
  ratio between the two friction coefficients.

Time Mapping using the mean-square displacement
of the chain center of mass

• Check transferability of tx for different
  systems, conditions (?, T, P, ).

17
Polymer Melts through CG MD Simulations Self
Diffusion Coefficient

• Correct raw CG diffusion data using a time
  mapping approach.
• V. Harmandaris and K. Kremer, Soft Matter, 5,
  3920 (2009)
• Crossover regime from Rouse to reptation
  dynamics. Include the chain end (free volume)
  effect.

-- Rouse D M-1 -- Reptation D M-2
Crossover region -- CG MD Me 28.000-33.000
gr/mol -- Exp. Me 30.0000-35.000 gr/mol

-- Exp. Data NMR Sillescu et al. Makromol.
Chem., 188, 2317 (1987)
18
CG Simulations Application Non-Equilibrium
Polymer Melts

• Non-equilibrium molecular dynamics (NEMD)
  modeling of systems out of equilibrium - flowing
  conditions.

• NEMD Equations of motion (pSLLOD)

simple shear flow

• In canonical ensemble (Nose-Hoover) C. Baig et
  al., J. Chem. Phys., 122, 11403, 2005

19
CG Simulations Application Non-Equilibrium
Polymer Melts

• NEMD equations of motion are not enough we
  need the proper periodic boundary conditions.
• Steady shear flow

Lees-Edwards Boundary Conditions
simple shear flow
20
CG Polymer Simulations Non-Equilibrium Systems

• CG NEMD - Remember CG interaction potentials
  are calculated as potential of mean force (they
  include entropy). In principle UCG(x,T) should be
  obtained at each state point, at each flow field.

Important question How well polymer systems
under non-equilibrium (flowing) conditions can be
described by CG models developed at equilibrium?
Method C. Baig and V. Harmandaris,
Macromolecules, 43, 3156 (2010)
Use of existing equilibrium CG polystyrene (PS)
model.

• Direct comparison between atomistic and CG NEMD
  simulations for various flow fields. Strength of
  flow (Weissenberg number, Wi 0.3 - 200)

• Study short atactic PS melts (M2kDa, 20
  monomers) by both atomistic and CG NEMD
  simulations.

21
CG Non-Equilibrium Polymers Conformations

• Properties as a function of strength of flow
  (Weissenberg number)

• Conformation tensor

R

• Atomistic cxx asymptotic behavior at high Wi
  because of (a) finite chain extensibility, (b)
  chain rotation during shear flow.
• CG cxx allows for larger maximum chain
  extension at high Wi because of the softer
  interaction potentials.

22
CG Non-Equilibrium Polymers Conformation Tensor

• cyy, czz excellent agreement between atomistic
  and CG configurations.

23
CG Non-Equilibrium Polymers Dynamics

• Is the time mapping factor similar for different
  flow fields?
• C. Baig and V. Harmandaris, Macromolecules, 43,
  3156 (2010)

Translational motion

• Purely convective contributions from the applied
  strain rate are excluded.

• Very good qualitative agreement between
  atomistic and CG (raw) data at low and
  intermediate flow fields.

24
CG Non-Equilibrium Polymers Dynamics
Orientational motion

• Rotational relaxation time small variations at
  low strain rates, large decrease at high flow
  fields.
• Good agreement between atomistic and CG at low
  and intermediate flow fields.

25
CG Non-Equilibrium Polymers Dynamics

• Time mapping parameter as a function of the
  strength of flow.

• Strong flow fields smaller time mapping
  parameter ? effective CG bead friction decreases
  less than the atomistic one.
• Reason CG chains become more deformed than the
  atomistic ones.

26
Conclusions

• Hierarchical systematic CG models, developed
  from isolated atomistic chains, correctly predict
  polymer structure and dimensions.

• Time mapping using dynamical information from
  atomistic description allow for quantitative
  dynamical predictions from the CG simulations,
  for many cases.

• Overall speed up of the CG MD simulations,
  compared with the atomistic MD, is 3-5 orders
  of magnitude.

• System at non-equilibrium conditions can be
  accurately studied by CG NEMD simulations at low
  and medium flow fields.

• Deviations between atomistic and CG NEMD data at
  high flow fields due to softer CG interaction
  potentials.

27
Challenges Current Work

• Estimation of CG interaction potential (free
  energies) Check improve all assumptions
• Ongoing work with M. Katsoulakis, D.
  Tsagarogiannis, A. Tsourtis

• Quantitative prediction of dynamics based on
  statistical mechanics
• e.g. Mori-Zwanzig formalism (Talk by Rafael
  Delgado-Buscalioni)

• Parameterizing CG models under non-equilibrium
  conditions
• e.g. Information-theoretic tools (Talk by Petr
  Plechac)

• Application of the whole procedure in more
  complex systems
• e.g. Multi-component biomolecular systems,
  hybrid polymer based nanocomposites
• Ongoing work with A. Rissanou

28
ACKNOWLEDGMENTS
Prof. C. Baig School of Nano-Bioscience and
Chemical Engineering, UNIST University, Korea
Funding ACMAC UOC Regional Potential Grant
FP7 DFG SPP 1369 Interphases and Interfaces
, Germany Graphene Research
Center, FORTH Greece
29
EXTRA SLIDES
30
APPLICATION PRIMITIVE PATHS OF LONG POLYSTYRENE
MELTS

• Describe the systems in the levels of primitive
  paths
• V. Harmandaris and K. Kremer, Macromolecules,
  42, 791, (2009)

• Entanglement Analysis using the Primitive Path
  Analysis (PPA) method
• Evereaers et al., Science 2004, 303, 823.

PP PS configuration (50kDa)
CG PS configuration (50kDa)

• Calculate directly PP contour length Lpp,, tube
  diameter

-- PP CG PS Ne 180 20 monomers
31
CALCULATION OF Me in PS Comparison Between
Different Methods

• Several methods to calculate Me broad spread
  of different estimates
• V. Harmandaris and K. Kremer, Macromolecules,
  42, 791 (2009)

Method T(K) Ne (mers) Reference
Rheology 423 140 15 Liu et al., Polymer, 47, 4461 (2006)
Self-diffusion coefficient 458 280-320 Antonieti et al., Makrom. Chem., 188, 2317 (1984)
Self-diffusion coefficient 463 240-300 This work
Segmental dynamics 463 110 30 This work
Entanglement analysis 463 180 20 This work
Entanglement analysis 413 124 Spyriouni et al., Macromolecules, 40, 3876 (2007)
32
MESOSCOPIC BOND ANGLE POTENTIAL OF PS
Distribution function PCG(?,T)?
CG Bending potential UCG(?,T)?
33
CG Simulations Applications Equilibrium
Polymer Melts

• Systems Studied Atactic PS melts with molecular
  weight from 1kDa (10 monomers) up to 50kDa (1kDa
  1000 gr/mol).

• NVT Ensemble.
• Langevin thermostat (T463K).
• Periodic boundary conditions.

34
STATIC PROPERTIES Radius of Gyration
RG
35
SMOOTHENING OF THE ENERGY LANDSCAPE
Qualitative prediction due to lost degrees of
freedom (lost forces) in the local level ?
Local friction coefficient in CG mesoscopic
description is smaller than in the
microscopic-atomistic one
Rouse
Reptation
CG diffusion coefficient is larger than the
atomistic one
36
Time Mapping Parameter Translational vs
Orientational Dynamics