1Computational/Mathematical Physics 3MA010 2Three Subjects Partial differential equations and numerical methods (PDE) (Prof. dr. Federico Toschi, Dr. Gianluca Di Staso, MSc. Alessandro Corbetta, MSc Ilian Pihlajamaa) Brownian dynamics (BD) (Dr. Alexey Lyulin, MSc. I. Pihlajamaa, MSc. Maarten Boomstra, Msc Nikos Sigalas) Monte-Carlo (MC) (Dr. Henk Huinink, MSc. I. Pihlajamaa, MSc. Maarten Boomstra, Msc Nikos Sigalas) 3● Check Canvas and OnCourse pages for the detailed planning, and all course material ● PDE: 6 lectures/instructions ● BD: 4 lectures/instructions ● MC: 4 lectures/instructions Software: MatLab 4PDE: lectures, instructions, home work. OnCourse testing will be used. BD and MC parts: students will work on theoretical and modelling exercises. – Written reports are expected, with findings from the simulations Deadlines for reports: November 6 retake 5 Evaluation: testing for PDE part, 2 reports for BD and MC parts, NO final exam Final mark: 40% – PDE 30% – BD report 30% – MC report Lecture 1: Random numbers and distributions 6 – Stochastic modelling – Random number generators – Distribution functions and its momenta – Central Limit Theorem – Non-uniform distributions – Markov processes – Random walks and applications 7 John von NeumannNicholas Metropolis 8 Andrey Markov 9 Random walks in polymer physics t0=0 t1=t Idea: the collective activity of many randomly moving participants can be effectively predictable, even if their individual motions are not 10 Polymers as long molecular chains Number of monomer units in a chain N>>1 (for DNA N~109 – 1010) 11 12 Freely-jointed flexibility mechanism The flexibility is located in the freely rotating junction points. This mechanism is normally not characteristic for real chains, but it is used for model theoretical calculations. 13 Portrait of a polymeric coil Typical conformations of few PE-like chains simulated on a computer Chain conformation – trajectory of a Brownian particle 14 Let us consider ideal N-segment freely-jointed chain ( each segment of length l ) 15 Thus, the conformation of an ideal chain is far from the rectilinear one. The chain conformation is equivalent to the trajectory of a Brownian particle. 16 17 – Stochastic modelling – Random number generators – Distribution functions and its momenta – Central Limit Theorem – Non-uniform distributions – Markov processes – Random walks and applications Lecture 1 recap 18 Lecture 2: Correlation functions and diffusion equations – Time series analysis – ACF computation in practice – Brownian motion – Diffusion equation – Langevin equation – Ermak-McCammon algorithm 19 20 Autocorrelation functions time correlation function: if A=B – autocorrelation function: ( ) (0) ( )AAC t A A t= ( ) (0) ( )ABC t A B t= (0) ( ) ( ) (0)A B t A t B= 2(0)AAC A= normalization: 2 22 (0) ( ) ( )AA A A t AC t A A = (0) 1AAC = ( ) 0AAC t = ∞ = 21 Time correlation function may be evaluated as a time average, assuming the system is ergodic – the phase space average is equal to a time average 0 1( ) lim ( ) ( ) T t AA T C t d A A t T t τ τ τ →∞ = + ∫ simulations: phase space trajectory is determined at discrete time steps, the integral is expressed as a sum + = = = ∑1 1( ) ( ) ( ), 0,1,2…, N k AA j k j c j C k t A t A t k N N k N – total number of time steps t – time step Nc << N why 22 time AC F 1 0 exp(-t/τ) τ Reality, i.e. simulations Nc << N why 23 Practical realization tNti+1tit1 t3t2 time a(tN)a(t1) a(ti) ACF=t 21 1 2 2( ) ( ) ( ) ( ) ... ( ) ( )N Na t a t a t a t a t a t a N + + + = 1 2 2 3 1( ) ( ) ( ) ( ) ... ( ) ( ) .. (0) ( ) 1 i ia t a t a t a t a t a t a a t N ++ + + + = t2 1 3 2 4 2) 3 ( ) ( ) ( ) ( ) ... ( ) ( ) ( ) 2 N Na t a t a t a t a t a t ACF t N + + + = t3 ….. tN equidistant time intervals How many terms exist in order to calculate this point Important example - white noise Brownian motion t0=0 t1=t 25 ( ) ( ) ( ) ( )Ri i i it V t m tζ + =r r F r Langevin equation 27 if /L m mτ τ ζ<< = 0 ( ) ( ) 2 ( ) 2 ( ) R i R R i j ij B i ijt t Z t t k T t tδ δ ζ δ δ = ′ ′ ′= = F F F Random force FR 28 Drag force Fd mζ >>v a d i i iζ= F v 2~ , ~ , ~t t te e eω ω ωω ωr v a /m mζ ω ω ζ>> → << /BD v mτ τ ζ>> = ( ) ( ) ( ) 0Ri i it V tζ + =r r F Langevin equation in non-inertial approximation /BD v mτ τ ζ>> = 29 30 Ermak-McCammon algorithm 1 1pot R pot R i i i i i i i i i v F F v F Fζ ζ ζ = + = + t trttr tv iii + = )()( )( 1 1( ) ( ) ( ) ( ) i i pot R i i i ir t t r t F t t F tζ ζ+ = + + B( ) ( ) 2 R R i j i ijF t F t k T tζ δ = t = 31 Time scales
σ length ζ friction T temperature m mass 3 2 3 10 20 3 23 6 6 3.14 10 10 4 102~ ~ ~ ~ 20 1.38 10 300BD B B Pa s m ps k T k T J σπηζστ /cor r BDmτ τ ζ τ<< = << 2, ~ 10BD BDt tτ τ << Lecture 3: Random force, Brownian dynamics of a polymer chain 32 - Random force, spectral density, magnitude - Langevin equation (continued) - Fokker-Planck equation - Einstein-Smolukhovski equation - Equilibration of the initial configuration - Applications for polymer modelling, hydrodynamic interactions Lecture 4: Brownian dynamics and diffusion Brownian dynamics of polymers 35 Ηψ = Εψ F = ma exp(- E/kT) domain quantum chemistry Brownian dynamics, MD Monte Carlo mesoscale continuum What and Where: Scales in Simulations Length scale 10-10 M 10-8 M 10-6 M 10-4 M 10-12 S 10-8 S 10-6 S 36 Historical Background Brownian motion Robert Brown (1773-1858) 1827: 1905-1966: Einstein, Smolukhovski, Langevin, Fokker, Planck 1978: Ermak, McCammon, BD simulations of dimers and trimers in solution 37 Brownian Dynamics rij v(rij) i j 0tot pot d Ri i i i m= + + =F F F F a Langevin equation 38 Fpot, the Force Field intermolecular interactions intramolecular nonbonded torsional bond stretch valence angle bend 39 Vstretch = 1 2 kij(rij – r 0 ij)2 Vbend = 1 2 kijk(Θijk – Θ 0 ijk)2 Vtor= 1 2 kijkl(1- cosnφijkl) VLJ = Aijrij-12 – Cijrij-6 VCoulomb = qiqj rij V(r) = Vstretch + Vbend + Vtor + VLJ + VCoulomb ( )pot V= F r HYBRANE from DSM g=1 g=2 g=3 g=4 g=0 NN O O O OO O OOC HO O O N OOC OOC O O N O O O O N O O O O N O O O O N OH O O O N O C11H23 O O O O C11H23 C11H23 C11H23 C11H23 C11H23 C11H23 C11H23 fatty acid ester 41 v v l Hyperbranched Polymers synthesized in a much more economical one-pot reaction 42 ( ) ( )ttvFkTtrr 0i0i0j j 0 ij 0 ii Φ+ + += ∑ D/ - velocity (the shear force, ) - Lennard-Jones + SHAKE for rigid bonds (+Coulomb) - random force - diffusion tensor, hydrodynamic interactions via Rotne-Prager-Yamakawa tensor γ 0 iv 0 jF 0 iΦ 0 ijD l - bond length a - hydrodynamic bead radius Brownian dynamics Langevin equation; Ermak-McCammon algorithm Random force Potential force Diffusion tensor 000 2)()( ijji ttt D = Φ Φ 0 0 0 1 ( ) N k j k jk j UF rr σλ = = ∑ r 2 1/ 2 2 2 2 32( / 3) (3 / 4 )( / ) ( ) ( ) 3 ij ij ij ij ij ij ij ij ij R R R RaD h kT l R R R R α β α β αβ αβ αβπ ζ δ δ = + + + = ij ijijij ij R RR aa R kTD βα αβ αβ δζ ) 32 3() 32 9 1()/( 2ijR a≥ 2ijR a< 44 Problems (neutral) dendrimers in shear flow structure of charged dendrimers complexes with linear polyelectrolytes 2 max 3 3ln ( 1) ln(1 ) ( 1) 2 Q B B FF g g s k T k T φ α α φε φ = + + + + Viscosity, shear thinning shear viscosity dimensionless intrinsic viscosity low shear data were fit using xyτη γ = * ( )[ ] s Bk T η η η λ = [ ] * 20 1[ ] aη η γ = Mark-Houwink: [ ]0 aKMη = 47 Explanation of [η] Maximum 0 1 2 3 4 5 6 1 10 [η ] 0* /N generation dendrimers with ( ) and without ( ) spacers [η] ~ [molecular density]-1 ~ V/M ; M ~ N Linear Polymer Dendrimer V ~ (g+1)3 M ~ 2(g+1) [ ] 2~ ~R Nη [ ] ( )( ) 3 1 1 2 g g ~η + + INITIALLY GROWS FASTER EVENTUALLY TAKES OVER Mark-Houwink behaviour 3 ~ V R 0.5~R N dendrimer 48 Charged Dendrimers Interplay between excluded-volume and Coulomb interactions influence of perturbations of the regular dendrimer structure v v + + + + + + + + + + + + 49 ,4)( 612 = rr rU LJLJ σσε εLJ = 0.3kbT σ = 0.8l Lennard-Jones: / 2 2( ) , 4 Dr r C b B B b e eU r k Tz r k T λ λ πε = = Debye-Hückel: 0 < z2 < 5 λB= l Debye radius rD=100l Force Field 50 2 max 3 3ln ( 1) ln(1 ) ( 1) 2 Q LJ B B FF g g s k T k T φ α α φε φ = + + + + 1. stretching and confinement entropy 2. excluded volume and van der Waals attraction 3. Coulomb Mean-Field Theory 1 2 single linear branch in a g=6 dendrimer - linear expansion factor; ( ) g g R R ideal α = 3 ( ) 2 2 / 2~ 1 D Q R rB D B F Q r e k T R λ - mean segment densityφ 51 Neutral dendrimer Linear-expansion factor vs. solvent quality BD and the mean-field theory are in agreement 52 Layering of monomers for charged dendrimers Distribution of Monomers 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 0.5 0.6 ρ(r) r / l gint = 0 gint = 1 gint = 2 gint = 3 gint = 4 gint = 5 0 1 2 3 4 5 6 0.0 0.2 0.4 0.6 0.8 ρ(r) r / l gint = 0 gint = 1 gint = 2 gint = 3 gint = 4 gint = 5 g=5 dendrimer uncharged charged 53 Perspectives Drug delivery with charged HBPs low salt, good solvent high salt, poor solvent