Split LB tutorial and create polymers and Langevin simulation tutorial

doc/tutorials/langevin_dynamics/langevin_dynamics.ipynb

@@ -46,7 +46,7 @@
     "where one can see that the strength of the fluctuation force depends on the friction coefficient and the system temperature $T$.\n",
     "(The Boltzmann constant is denoted as $k_\\text{B}$ and the two $\\delta$s are Dirac delta functions with respect to particle id and time, respectively.)\n",
     "\n",
-    "The Langevin equation obviously provides a very straightforward approach to represent both Stokes friction and Brownian motion acting on a particle. However, please be aware that due to its simplicity this technique also has its limitations, i.e., the drag force and diffusion to adjacent particles are uncorrelted and hydrodynamic interactions with the particle's surroundings are neglected. Langevin dynamics should therefore be used with caution, particularly in systems with high particle densities or strong hydrodynamic coupling.\n",
+    "The Langevin equation obviously provides a very straightforward approach to represent both Stokes friction and Brownian motion acting on a particle. However, please be aware that due to its simplicity this technique also has its limitations, i.e., the drag force and diffusion to adjacent particles are uncorrelated and hydrodynamic interactions with the particle's surroundings are neglected. Langevin dynamics should therefore be used with caution, particularly in systems with high particle densities or strong hydrodynamic coupling.\n",
     "\n",
     "In the Langevin equation, only ensemble averaged properties of the random ${\\bf f}$ are specified. Consequently, it doesn't  make sense to look at a single deterministic solution of Eq. (1). Instead, one should measure ensemble averaged quantities that characterize the dynamics of the spherical Brownian particle. The simplest quantity is the so-called mean square displacement (MSD) after time $\\tau$\n",
     "\n",

doc/tutorials/lattice_boltzmann/NotesForTutor.md

@@ -1,19 +1,28 @@
 # Part 1: The Lattice-Boltzmann method
 
-## Physical learning goals
+## Physics learning goals
 
-* Learning basic concepts of the LBM
+After the tutorial, students should be able to:
 
-## Espresso learning goals
+* Explain the basic idea of the Lattice-Boltzmann method
 
-* Learning the LB interface in ESPResSo
+## ESPResSo learning goals
+
+In the course of this tutorial, students should learn to:
+
+* Set up a Lattice-Boltzmann fluid
 
 # Part 2: Planar Poiseuille flow
 
-## Physical learning goals
+## Physics learning goals
+
+After the tutorial, students should be able to:
+
+* Describe the geometry and expected flow field of a Poiseuille flow
 
-* Modeling a laminar fluid flow between two walls
+## ESPResSo learning goals
 
-## Espresso learning goals
+In the course of this tutorial, students should learn to:
 
-* Setting up an LB fluid and LB boundaries
+* Set up an LB fluid with LB boundaries
+* Retrieve data from the LB grid

doc/tutorials/lattice_boltzmann/lattice_boltzmann_poiseuille_flow.ipynb

@@ -203,13 +203,12 @@
    "outputs": [],
    "source": [
     "logging.info(\"Extract fluid velocities along the x-axis\")\n",
-    "fluid_positions = np.zeros(lbf.shape[0])\n",
-    "fluid_velocities = np.zeros(lbf.shape[0])\n",
-    "for x in range(lbf.shape[0]):\n",
-    "    # Average over the nodes in y direction\n",
-    "    v = [lbf[x, y, 0].velocity[1] for y in range(lbf.shape[1])]\n",
-    "    fluid_velocities[x] = np.average(v)\n",
-    "    fluid_positions[x] = (x + 0.5) * AGRID\n",
+    "\n",
+    "fluid_positions = (np.arange(lbf.shape[0]) + 0.5) * AGRID\n",
+    "# get all velocities as Numpy array and extract y components only\n",
+    "fluid_velocities = (lbf[:,:,:].velocity)[:,:,:,1]\n",
+    "# average velocities in y and z directions (perpendicular to the walls)\n",
+    "fluid_velocities = np.average(fluid_velocities, axis=(1,2))\n",
     "\n",
     "\n",
     "def poiseuille_flow(x, force_density, dynamic_viscosity, height):\n",

doc/tutorials/lattice_boltzmann/lattice_boltzmann_theory.ipynb

@@ -39,15 +39,15 @@
     "### Notes on the **ESPResSo** version you will need\n",
     "\n",
     "**ESPResSo** offers both CPU and GPU support for LB. We recommend using the GPU code,\n",
-    "as it is much (100x) faster than the CPU code.\n",
+    "as it is much (10x) faster than the CPU code. It should be noted, that the GPU code only uses single precision floating point numbers whereas the CPU code uses double precision. For most applications, however, single precision is sufficient.\n",
     "\n",
     "For the tutorial you will have to compile in the following features:\n",
     "```c++\n",
     "#define LB_BOUNDARIES\n",
     "#define LB_BOUNDARIES_GPU\n",
     "#define LENNARD_JONES\n",
     "```\n",
-    "Please uncomment them in the <tt>myconfig.hpp</tt> and compile **ESPResSo** using this <tt>myconfig.hpp</tt>."
+    "Please uncomment the features in the <tt>myconfig.hpp</tt> and compile **ESPResSo** using this <tt>myconfig.hpp</tt>. This is not necessary if you do not use a custom <tt>myconfig.hpp</tt>, since the features are activated by default. For more information on configuring **ESPResSo** and how to activate CUDA (for GPU computation), refer to the [documentation](https://espressomd.github.io/doc/installation.html). "
    ]
   },
   {
@@ -298,7 +298,9 @@
     "+ <tt>pressure_tensor_neq</tt>: the nonequilibrium part of the pressure tensor (3x3 matrix).\n",
     "+ <tt>population</tt>: the 19 populations of the D3Q19 lattice.\n",
     "+ <tt>boundary</tt>: the boundary flag.\n",
-    "+ <tt>density</tt>: the local density."
+    "+ <tt>density</tt>: the local density.\n",
+    "\n",
+    "Slicing is supported, e.g. to obtain all velocity vectors in the LB fluid as a Numpy array, use <tt>lbf[:,:,:].velocity</tt>."
    ]
   },
   {

doc/tutorials/polymers/NotesForTutor.md

@@ -1,16 +1,20 @@
 # Polymer theory and Brownian motion of a polymer in a fluid
 
-## Physical learning goals
+## Physics learning goals
 
 * Getting a short introduction to the Flory theory of polymers
 * Modeling the Brownian motion of a polymer
 * Learning about the hydrodynamic radius, radius of gyration, end-to-end distance
 * Reproducing theoretical results of the Flory theory of polymers
 * Using the Green-Kubo relation to calculate the diffusion coefficient of a polymer
 * Estimating the correlation-corrected standard error of the mean of a time series
-* Observing the effect of implicit hydrodynamics vs. LB hydrodynamics on polymer diffusion
 
-## Espresso learning goals
+After the tutorial, students should be able to:
+
+* Explain how hydrodynamic interactions qualitatively affect the diffusion of a polymer
+* Name the three important regimes for polymer diffusion and the corresponding scaling laws
+
+## ESPResSo learning goals
 
 * Creating linear polymers with ESPResSo functions
 * Using tools from the analysis module not based on the Observable framework

doc/tutorials/polymers/polymers.ipynb

@@ -25,7 +25,7 @@
     "are present, one expects a scaling law $D \\propto N ^{- 1}$ and if they are present, a scaling law\n",
     "$D \\propto N^{- \\nu}$ is expected. Here $\\nu$ is the Flory exponent that plays a very prominent\n",
     "role in polymer physics. It has a value of $\\sim 3/5$ in good solvent conditions in 3D.\n",
-    "Discussions on these scaling laws can be found in polymer physics textbooks like [1–3].\n",
+    "Discussions on these scaling laws can be found in polymer physics textbooks like <a href='#[1]'>[1]</a>, <a href='#[2]'>[2]</a>, and <a href='#[3]'>[3, chapter 8]</a>.\n",
     "\n",
     "The reason for the different scaling law is the following: when being transported, every monomer\n",
     "creates a flow field that follows the direction of its motion. This flow field makes it easier for\n",
@@ -38,7 +38,55 @@
     "\\end{equation}\n",
     "\n",
     "This hydrodynamic radius exhibits the scaling law  $R_h \\propto N^{\\nu}$\n",
-    "and the diffusion coefficient of a long polymer is proportional to its inverse $R_h$.\n",
+    "and the diffusion coefficient of a long polymer is proportional to its inverse $R_h$."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Polymer models\n",
+    "\n",
+    "The diffusion coefficient $D$ of a spherical particle in a carrier fluid experiencing drag is\n",
+    "related to the friction coefficient $\\zeta$ via the Einstein relation:\n",
+    "\n",
+    "\\begin{equation}\n",
+    "  D = \\frac{k_{\\mathrm{B}}T}{\\zeta},\n",
+    "\\end{equation}\n",
+    "\n",
+    "with $k_{\\mathrm{B}}$ the Boltzmann constant and $T$ the temperature. For a sphere of radius $R$\n",
+    "moving in a fluid of viscosity $\\eta$, the friction coefficient is obtained via the Stokes law:\n",
+    "\n",
+    "\\begin{equation}\n",
+    "  \\zeta = 6\\pi\\eta R.\n",
+    "\\end{equation}\n",
+    "\n",
+    "Combining both equations yields the Stokes–Einstein relation:\n",
+    "\n",
+    "\\begin{equation}\n",
+    "  D = \\frac{k_{\\mathrm{B}}T}{6\\pi\\eta R}.\n",
+    "\\end{equation}\n",
+    "\n",
+    "The simplest description of a polymer is the Rouse model, where beads are connected by springs.\n",
+    "All beads experience a drag from the solvent, and the friction coefficient $\\gamma$ is identical\n",
+    "for all beads. The solvent flows freely between beads and hydrodynamic interactions are neglected.\n",
+    "The diffusion coefficient takes the following form:\n",
+    "\n",
+    "\\begin{equation}\n",
+    "  D_{\\mathrm{R}} = \\frac{D_0}{N} = \\frac{k_{\\mathrm{B}}T}{\\gamma N},\n",
+    "\\end{equation}\n",
+    "\n",
+    "where $D_0$ is the diffusion coefficient of a single bead.\n",
+    "\n",
+    "To account for hydrodynamic interactions mediated by the solvent, i.e. the transport of solvent\n",
+    "in contact with the beads and the correlation in the motion of beads due to the carried solvent,\n",
+    "the Zimm model was created. For an ideal chain, it takes the following form:\n",
+    "\n",
+    "\\begin{equation}\n",
+    "  D_{\\mathrm{Z}} = \\frac{8}{3\\sqrt{6\\pi^3}}\\frac{k_B T}{\\eta R} \\simeq 0.196\\frac{k_B T}{\\eta b N^{\\nu}},\n",
+    "\\end{equation}\n",
+    "\n",
+    "with $R$ the radius of the polymer and $b$ the length of the spring connecting the beads.\n",
     "\n",
     "For shorter polymers there is a transition region. It can be described\n",
     "by the Kirkwood–Zimm model:\n",
@@ -58,9 +106,11 @@
     "  D=\\frac{D_0}{N} + \\frac{k_B T}{6 \\pi \\eta } \\left\\langle \\frac{1}{R_h} \\right\\rangle \\left( 1- \\left\\langle\\frac{R_h}{L} \\right\\rangle \\right)\n",
     "\\end{equation}\n",
     "\n",
-    "It is quite difficult to prove this formula computationally with good accuracy.\n",
-    "It will need quite some computational effort and a careful analysis. So please don't be\n",
-    "too disappointed if you don't manage to do so.\n",
+    "It is quite difficult to fit this analytical expression to simulation data with good accuracy.\n",
+    "It will need a LB fluid, long simulation times and a careful analysis. For this tutorial, we\n",
+    "will use an implicit solvent and short polymer lengths to keep the runtime short. If you want\n",
+    "to collect data suitable for the Zimm model, simply set the global variable `POLYMER_MODEL` to\n",
+    "`'Zimm'`, and increase the box size and number of beads in the polymer.\n",
     "\n",
     "We want to determine the long-time self diffusion coefficient from the mean square\n",
     "displacement of the center-of-mass of a single polymer. For large $t$ the mean square displacement is\n",
@@ -71,7 +121,7 @@
     "    D = \\lim_{t\\to\\infty}\\left[ \\frac{1}{6t} \\left\\langle \\left(\\vec{r}(t) - \\vec{r}(0)\\right)^2 \\right\\rangle \\right].\n",
     "\\end{equation}\n",
     "\n",
-    "This equation can be found in virtually any simulation textbook, like [4]. We will set up a\n",
+    "This equation can be found in virtually any simulation textbook, like <a href='#[4]'>[4]</a>. We will set up a\n",
     "polymer in an implicit solvent, simulate for an appropriate amount of time, calculate the mean square\n",
     "displacement as a function of time and obtain the diffusion coefficient from a linear\n",
     "fit. However we will have a couple of steps inbetween and divide the full problem into\n",
@@ -161,11 +211,12 @@
     "def build_polymer(system, n_monomers, polymer_params, fene):\n",
     "    positions = espressomd.polymer.linear_polymer_positions(\n",
     "        beads_per_chain=n_monomers, **polymer_params)\n",
+    "    p_previous = None\n",
     "    for i, pos in enumerate(positions[0]):\n",
-    "        pid = len(system.part)\n",
-    "        system.part.add(id=pid, pos=pos)\n",
-    "        if i > 0:\n",
-    "            system.part[pid].add_bond((fene, pid - 1))\n",
+    "        p = system.part.add(pos=pos)\n",
+    "        if p_previous is not None:\n",
+    "            p.add_bond((fene, p_previous))\n",
+    "        p_previous = p\n",
     "```"
    ]
   },
@@ -200,7 +251,7 @@
     "def correlator_msd(pids_monomers, tau_max):\n",
     "    com_pos = espressomd.observables.ComPosition(ids=pids_monomers)\n",
     "    com_pos_cor = espressomd.accumulators.Correlator(\n",
-    "        obs1=com_pos, tau_lin=16, tau_max=tau_max, delta_N=1,\n",
+    "        obs1=com_pos, tau_lin=16, tau_max=tau_max, delta_N=5,\n",
     "        corr_operation=\"square_distance_componentwise\", compress1=\"discard1\")\n",
     "    return com_pos_cor\n",
     "\n",
@@ -408,7 +459,7 @@
     "### 3. Data analysis\n",
     "\n",
     "We will calculate the means of time series with error bars obtained from\n",
-    "the correlation-corrected standard error of the mean [5,6]."
+    "the correlation-corrected standard error of the mean [<a href='#[5]'>5</a>,<a href='#[6]'>6</a>]."
    ]
   },
   {
@@ -433,7 +484,7 @@
     "    '''\n",
     "    Calculate the mean and the correlation-corrected standard error\n",
     "    of the mean of time series by integrating the autocorrelation\n",
-    "    function. See Janke 2002 [8] and Weigel, Janke 2010 [9].\n",
+    "    function. See Janke 2002 [5] and Weigel, Janke 2010 [6].\n",
     "\n",
     "    Due to the short simulation length, it is not possible to fit an\n",
     "    exponential to the long-time tail. Instead, return a percentile.\n",
@@ -864,12 +915,12 @@
    "source": [
     "## References\n",
     "\n",
-    "[1] P. G. de Gennes. *Scaling Concepts in Polymer Physics*. Cornell University Press, Ithaca, NY, 1979.  \n",
-    "[2] M. Doi. *Introduction to Polymer Physics.* Clarendon Press, Oxford, 1996.  \n",
-    "[3] Michael Rubinstein and Ralph H. Colby. *Polymer Physics.* Oxford University Press, Oxford, UK, 2003.  \n",
-    "[4] Daan Frenkel and Berend Smit. *Understanding Molecular Simulation*, section 4.4.1. Academic Press, San Diego, second edition, 2002.  \n",
-    "[5] W. Janke, Statistical analysis of simulations: Data correlations and error estimation, *Quantum Simulations of Complex Many-Body Systems: From Theory to Algorithms*, Lecture Notes, J. Grotendorst, D. Marx, A. Muramatsu (Eds.), John von Neumann Institute for Computing, 10:423–445, 2002. https://www.physik.uni-leipzig.de/~janke/Paper/nic10_423_2002.pdf  \n",
-    "[6] M. Weigel, W. Janke, Error estimation and reduction with cross correlations, *Phys. Rev. E*, 81:066701, 2010, doi:[10.1103/PhysRevE.81.066701](https://doi.org/10.1103/PhysRevE.81.066701); Erratum-ibid 81:069902, 2010, doi:[10.1103/PhysRevE.81.069902](https://doi.org/10.1103/PhysRevE.81.069902).  "
+    "<a id='[1]'></a>[1] P. G. de Gennes. *Scaling Concepts in Polymer Physics*. Cornell University Press, Ithaca, NY, 1979.  \n",
+    "<a id='[2]'></a>[2] M. Doi. *Introduction to Polymer Physics.* Clarendon Press, Oxford, 1996.  \n",
+    "<a id='[3]'></a>[3] Michael Rubinstein and Ralph H. Colby. *Polymer Physics.* Oxford University Press, Oxford, UK, 2003. ISBN: 978-0-19-852059-7  \n",
+    "<a id='[4]'></a>[4] Daan Frenkel and Berend Smit. *Understanding Molecular Simulation*, section 4.4.1. Academic Press, San Diego, second edition, 2002.  \n",
+    "<a id='[5]'></a>[5] W. Janke, Statistical analysis of simulations: Data correlations and error estimation, *Quantum Simulations of Complex Many-Body Systems: From Theory to Algorithms*, Lecture Notes, J. Grotendorst, D. Marx, A. Muramatsu (Eds.), John von Neumann Institute for Computing, 10:423–445, 2002. https://www.physik.uni-leipzig.de/~janke/Paper/nic10_423_2002.pdf  \n",
+    "<a id='[6]'></a>[6] M. Weigel, W. Janke, Error estimation and reduction with cross correlations, *Phys. Rev. E*, 81:066701, 2010, doi:[10.1103/PhysRevE.81.066701](https://doi.org/10.1103/PhysRevE.81.066701); Erratum-ibid 81:069902, 2010, doi:[10.1103/PhysRevE.81.069902](https://doi.org/10.1103/PhysRevE.81.069902).  "
    ]
   }
  ],