\n",
- "Figure 2: Real Ferrofluid exposed to an external magnetic field (neodymium magnet) [4]\n",
- "
\n",
+ "\n",
+ "
Figure 2: Real ferrofluid exposed to an external magnetic field (neodymium magnet) [4]
\n",
""
]
},
@@ -79,7 +73,7 @@
"id": "095339c1",
"metadata": {},
"source": [
- "## The Model"
+ "## The model"
]
},
{
@@ -124,7 +118,7 @@
" \\lambda = \\frac{\\tilde{u}_{\\text{DD}}}{u_T} = \\gamma \\frac{\\mu^2}{k_{\\text{B}}T\\sigma^3}\n",
"\\end{equation}\n",
"\n",
- "where $u_\\mathrm{T} = k_{\\text{B}}T$ is the thermal energy and $\\tilde{u}_{DD}$ is the absolute value of the dipole-dipole interaction energy at close contact (cc) and head-to-tail configuration (htt) (see figure 4) per particle, i.e. in formulas it reads\n",
+ "where $u_\\mathrm{T} = k_{\\text{B}}T$ is the thermal energy and $\\tilde{u}_{DD}$ is the absolute value of the dipole-dipole interaction energy at close contact (cc) and head-to-tail configuration (htt) (see figure 3) per particle, i.e. in formulas it reads\n",
"\n",
"\\begin{equation}\n",
" \\tilde{u}_{\\text{DD}} = \\frac{ \\left| u_{\\text{DD}}^{\\text{htt, cc}} \\right| }{2}\n",
@@ -143,11 +137,9 @@
"id": "62b95d9c",
"metadata": {},
"source": [
- "\n",
- "\n",
- "
\n",
- "Figure 4: Schematic representation of the head-to-tail configuration of two magnetic particles at close contact.\n",
- "
\n",
+ "\n",
+ "\n",
+ "
Figure 3: Schematic representation of the head-to-tail configuration of two magnetic particles at close contact.
\n",
""
]
},
@@ -182,7 +174,7 @@
"id": "b6d5a65e",
"metadata": {},
"source": [
- "## Compiling ESPResSo for this Tutorial"
+ "## Compiling ESPResSo for this tutorial"
]
},
{
@@ -218,7 +210,7 @@
"id": "a0bffbeb",
"metadata": {},
"source": [
- "## A Monolayer-Ferrofluid System in ESPResSo"
+ "## A monolayer-ferrofluid system in ESPResSo"
]
},
{
@@ -257,6 +249,7 @@
"import espressomd.magnetostatics\n",
"import espressomd.cluster_analysis\n",
"import espressomd.pair_criteria\n",
+ "import espressomd.zn\n",
"\n",
"espressomd.assert_features(['DIPOLES', 'DP3M', 'LENNARD_JONES'])\n",
"\n",
@@ -323,7 +316,8 @@
"id": "136aa645",
"metadata": {},
"source": [
- "## Exercise:\n",
+ "**Exercise:**\n",
+ "\n",
"How large does `BOX_SIZE` have to be for a system of `N_PART` particles with a volume (area) fraction `PHI`?\n",
"Define `BOX_SIZE`."
]
@@ -383,7 +377,8 @@
"outputs": [],
"source": [
"# Lennard-Jones interaction\n",
- "system.non_bonded_inter[0, 0].lennard_jones.set_params(epsilon=LJ_EPSILON, sigma=LJ_SIGMA, cutoff=LJ_CUT, shift=\"auto\")"
+ "system.non_bonded_inter[0, 0].lennard_jones.set_params(\n",
+ " epsilon=LJ_EPSILON, sigma=LJ_SIGMA, cutoff=LJ_CUT, shift=\"auto\")"
]
},
{
@@ -393,7 +388,7 @@
"source": [
"Now we generate random positions and orientations of the particles and their dipole moments. \n",
"\n",
- "**Hint:**\n",
+ "Hint:\n",
"It should be noted that we seed the random number generator of numpy. Thus the initial configuration of our system is the same every time this script will be executed. You can change it to another one to simulate with a different initial configuration."
]
},
@@ -402,9 +397,11 @@
"id": "c49729f7",
"metadata": {},
"source": [
- "## Exercise:\n",
+ "**Exercise:**\n",
+ "\n",
"How does one set up randomly oriented dipole moments?\n",
- "*Hint*: Think of the way that different methods could introduce a bias in the distribution of the orientations.\n",
+ "\n",
+ "Hint: Think of the way that different methods could introduce a bias in the distribution of the orientations.\n",
"\n",
"Create a variable `dip` as a `N_PART x 3` numpy array, which contains the randomly distributed dipole moments."
]
@@ -522,7 +519,7 @@
"For the simulation of our system we choose the velocity Verlet integrator.\n",
"After that we set up the thermostat which is, in our case, a Langevin thermostat to simulate in an NVT ensemble.\n",
"\n",
- "**Hint:**\n",
+ "Hint:\n",
"It should be noted that we seed the Langevin thermostat, thus the time evolution of the system is partly predefined.\n",
"Partly because of the numeric accuracy and the automatic tuning algorithms of Dipolar P3M and DLC where the resulting parameters are slightly different every time.\n",
"You can change the seed to get a guaranteed different time evolution."
@@ -561,7 +558,7 @@
"execution_count": null,
"id": "fd12b1f6",
"metadata": {
- "scrolled": true
+ "scrolled": false
},
"outputs": [],
"source": [
@@ -591,7 +588,7 @@
"execution_count": null,
"id": "95c4a8b6",
"metadata": {
- "scrolled": true
+ "scrolled": false
},
"outputs": [],
"source": [
@@ -607,7 +604,7 @@
"id": "a6c4c46b",
"metadata": {},
"source": [
- "## Sampling"
+ "### Sampling and cluster analysis"
]
},
{
@@ -615,6 +612,8 @@
"id": "af0906f5",
"metadata": {},
"source": [
+ "To quantify the number of clusters and their respective sizes, we now want to perform a cluster analysis.\n",
+ "For that we can use ESPResSo's [cluster analysis class](https://espressomd.github.io/doc/analysis.html#cluster-analysis).\n",
"The system will be sampled over 100 loops."
]
},
@@ -628,165 +627,12 @@
"LOOPS = 100"
]
},
- {
- "cell_type": "markdown",
- "id": "516c2f33",
- "metadata": {},
- "source": [
- "As the system is two dimensional, we can simply do a scatter plot to get a visual representation of a system state. To get a better insight of how a ferrofluid system develops during time we will create a video of the development of our system during the sampling. If you only want to sample the system simply go to [Sampling without animation](#Sampling-without-animation)"
- ]
- },
- {
- "cell_type": "markdown",
- "id": "60dc1653",
- "metadata": {},
- "source": [
- "### Sampling with animation"
- ]
- },
- {
- "cell_type": "markdown",
- "id": "5f7fe3a7",
- "metadata": {},
- "source": [
- "To get an animation of the system development we have to create a function which will save the video and embed it in an html string."
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "id": "ed0ee691",
- "metadata": {},
- "outputs": [],
- "source": [
- "import matplotlib.animation as animation\n",
- "import tempfile\n",
- "import base64\n",
- "\n",
- "VIDEO_TAG = \"\"\"\"\"\"\n",
- "\n",
- "\n",
- "def anim_to_html(anim):\n",
- " if not hasattr(anim, '_encoded_video'):\n",
- " with tempfile.NamedTemporaryFile(suffix='.mp4') as f:\n",
- " anim.save(f.name, fps=20, extra_args=['-vcodec', 'libx264'])\n",
- " with open(f.name, \"rb\") as g:\n",
- " video = g.read()\n",
- " anim._encoded_video = base64.b64encode(video).decode('ascii')\n",
- " plt.close(anim._fig)\n",
- " return VIDEO_TAG.format(anim._encoded_video)\n",
- "\n",
- "\n",
- "animation.Animation._repr_html_ = anim_to_html\n",
- "\n",
- "\n",
- "def init():\n",
- " # Set x and y range\n",
- " ax.set_ylim(0, BOX_SIZE)\n",
- " ax.set_xlim(0, BOX_SIZE)\n",
- " x_data, y_data = [], []\n",
- " part.set_data(x_data, y_data)\n",
- " return part,"
- ]
- },
- {
- "cell_type": "markdown",
- "id": "89d9844a",
- "metadata": {},
- "source": [
- "## Exercise:\n",
- "\n",
- "In the following an animation loop is defined, however it is incomplete.\n",
- "Extend the code such that in every loop the system is integrated for 100 steps.\n",
- "Afterwards `x_data` and `y_data` have to be populated by the folded $x$- and $y$- positions of the particles.\n",
- "\n",
- "(You may copy and paste the incomplete code template to the empty cell below.)\n",
- "\n",
- "```python\n",
- "def run(i):\n",
- " # < excercise >\n",
- " \n",
- " # Save current system state as a plot\n",
- " x_data, y_data = # < excercise >\n",
- " ax.figure.canvas.draw()\n",
- " part.set_data(x_data, y_data)\n",
- " print(f'progress: {(i + 1) * 100. / LOOPS:3.0f}%', end='\\r')\n",
- " return part,\n",
- "```"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "id": "91004d62",
- "metadata": {},
- "outputs": [],
- "source": [
- "# SOLUTION CELL\n",
- "def run(i):\n",
- " system.integrator.run(100)\n",
- "\n",
- " # Save current system state as a plot\n",
- " x_data, y_data = system.part.all().pos_folded[:, 0], system.part.all().pos_folded[:, 1]\n",
- " ax.figure.canvas.draw()\n",
- " part.set_data(x_data, y_data)\n",
- " print(f'progress: {(i + 1) * 100. / LOOPS:3.0f}%', end='\\r')\n",
- " return part,"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "id": "95caca0c",
- "metadata": {},
- "outputs": [],
- "source": []
- },
- {
- "cell_type": "markdown",
- "id": "0f7e9093",
- "metadata": {},
- "source": [
- "Now we use the animation class of matplotlib to save snapshots of the system as frames of a video which is then displayed after the sampling is finished. Between two frames are 100 integration steps.\n",
- "\n",
- "In the video chain-like and ring-like clusters should be visible, as well as some isolated monomers."
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "id": "3d11e830",
- "metadata": {
- "scrolled": true
- },
- "outputs": [],
- "source": [
- "fig, ax = plt.subplots(figsize=(10, 10))\n",
- "part, = ax.plot([], [], 'o')\n",
- "\n",
- "animation.FuncAnimation(fig, run, frames=LOOPS, blit=True, interval=0, repeat=False, init_func=init)"
- ]
- },
- {
- "cell_type": "markdown",
- "id": "a7ab793c",
- "metadata": {},
- "source": [
- "## Cluster analysis\n",
- "\n",
- "To quantify the number of clusters and their respective sizes, we now want to perform a cluster analysis.\n",
- "For that we can use ESPREsSo's [cluster analysis class](https://espressomd.github.io/doc/analysis.html?highlight=cluster%20analysis#cluster-analysis)."
- ]
- },
{
"cell_type": "markdown",
"id": "2f6d201a",
"metadata": {},
"source": [
- "## Exercise:\n",
+ "**Exercise:**\n",
"\n",
"Setup a cluster analysis object (`ClusterStructure` class) and assign its instance to the variable `cluster_structure`.\n",
"As criterion for the cluster analysis use a distance criterion where particles are assumed to be\n",
@@ -847,15 +693,19 @@
"id": "6b8d8e7f",
"metadata": {},
"source": [
- "### Sampling without animation"
+ "We initialize an instance of the ZnDraw visualizer to animate the sampling. As our system is constrained to 2D, the particles are on a single face of the cubic system box."
]
},
{
- "cell_type": "markdown",
+ "cell_type": "code",
"id": "2098e033",
+ "execution_count": null,
"metadata": {},
+ "outputs": [],
"source": [
- "The following code just samples the system and does a cluster analysis every loops (100 by default) simulation steps."
+ "color = {0: \"#7fc454\"} \n",
+ "radii = {0: LJ_SIGMA/2}\n",
+ "vis = espressomd.zn.Visualizer(system, colors=color, radii=radii)"
]
},
{
@@ -863,7 +713,7 @@
"id": "98b99dd6",
"metadata": {},
"source": [
- "## Exercise:\n",
+ "**Exercise:**\n",
"\n",
"Write an integration loop which runs a cluster analysis on the system, saving the number of clusters `n_clusters` and the size distribution `cluster_sizes`.\n",
"Take the following as a starting point:\n",
@@ -878,6 +728,7 @@
" for c in cluster_structure.clusters:\n",
" cluster_sizes.append(# < excercise >)\n",
" system.integrator.run(100)\n",
+ " vis.update()\n",
"```"
]
},
@@ -899,7 +750,8 @@
" n_clusters.append(len(cluster_structure.clusters))\n",
" for c in cluster_structure.clusters:\n",
" cluster_sizes.append(c[1].size())\n",
- " system.integrator.run(100)"
+ " system.integrator.run(100)\n",
+ " vis.update()"
]
},
{
@@ -915,7 +767,7 @@
"id": "0042827e",
"metadata": {},
"source": [
- "You may want to get a visualization of the current state of the system. For that we plot the particle positions folded to the simulation box using matplotlib."
+ "You may want to get a 2D visualization of the current state of the system. For that we plot the particle positions folded to the simulation box using matplotlib."
]
},
{
@@ -947,7 +799,7 @@
"id": "4c0828ea",
"metadata": {},
"source": [
- "## Cluster distribution"
+ "### Cluster distribution"
]
},
{
@@ -963,12 +815,12 @@
"id": "1c6be57c",
"metadata": {},
"source": [
- "## Exercise:\n",
+ "**Exercise:**\n",
"\n",
"Use `numpy` to calculate a histogram of the cluster sizes and assign it to the variable `size_dist`.\n",
"Take only clusters up to a size of 19 particles into account.\n",
"\n",
- "*Hint*: In order not to count clusters with size 20 or more, one may include an additional bin containing these.\n",
+ "Hint: In order not to count clusters with size 20 or more, one may include an additional bin containing these.\n",
"The reason for that is that `numpy` defines the histogram bins as half-open intervals with the open border at the upper bin edge.\n",
"Consequently clusters of larger sizes are attributed to the last bin.\n",
"By not using the last bin in the plot below, these clusters can effectively be neglected."
diff --git a/doc/tutorials/ferrofluid/ferrofluid_part2.ipynb b/doc/tutorials/ferrofluid/ferrofluid_part2.ipynb
index 453d5f9e97..80dcc881bb 100644
--- a/doc/tutorials/ferrofluid/ferrofluid_part2.ipynb
+++ b/doc/tutorials/ferrofluid/ferrofluid_part2.ipynb
@@ -59,6 +59,7 @@
"source": [
"import espressomd\n",
"import espressomd.magnetostatics\n",
+ "import espressomd.zn\n",
"\n",
"espressomd.assert_features(['DIPOLES', 'DP3M', 'LENNARD_JONES'])\n",
"\n",
@@ -166,7 +167,8 @@
"pos = box_size * np.hstack((np.random.random((N_PART, 2)), np.zeros((N_PART, 1))))\n",
"\n",
"# Add particles\n",
- "particles = system.part.add(pos=pos, rotation=N_PART * [(True, True, True)], dip=dip, fix=N_PART * [(False, False, True)])\n",
+ "particles = system.part.add(pos=pos, rotation=N_PART * [(True, True, True)],\n",
+ " dip=dip, fix=N_PART * [(False, False, True)])\n",
"\n",
"# Remove overlap between particles by means of the steepest descent method\n",
"MASS = 1.0\n",
@@ -269,7 +271,7 @@
"execution_count": null,
"id": "7750295c",
"metadata": {
- "scrolled": true
+ "scrolled": false
},
"outputs": [],
"source": [
@@ -305,95 +307,27 @@
"plt.show()"
]
},
- {
- "cell_type": "markdown",
- "id": "2f52d782",
- "metadata": {},
- "source": [
- "## Video of the development of the system"
- ]
- },
{
"cell_type": "markdown",
"id": "eb06ab55",
"metadata": {},
"source": [
- "You may want to get an insight of how the system develops in time. Thus we now create a function which will save a video and embed it in an html string to create a video of the systems development "
+ "You may want to get an insight of how the system develops in time. We use the ZnDraw visualiser:"
]
},
{
"cell_type": "code",
"execution_count": null,
- "id": "e25ba03b",
+ "id": "a6807ebd-2f97-476c-94e1-ea5fe1d383a4",
"metadata": {},
"outputs": [],
"source": [
- "import matplotlib.pyplot as plt\n",
- "import matplotlib.animation as animation\n",
- "import tempfile\n",
- "import base64\n",
- "\n",
- "VIDEO_TAG = \"\"\"\"\"\"\n",
- "\n",
- "\n",
- "def anim_to_html(anim):\n",
- " if not hasattr(anim, '_encoded_video'):\n",
- " with tempfile.NamedTemporaryFile(suffix='.mp4') as f:\n",
- " anim.save(f.name, fps=20, extra_args=['-vcodec', 'libx264'])\n",
- " with open(f.name, \"rb\") as g:\n",
- " video = g.read()\n",
- " anim._encoded_video = base64.b64encode(video).decode('ascii')\n",
- " plt.close(anim._fig)\n",
- " return VIDEO_TAG.format(anim._encoded_video)\n",
- "\n",
- "\n",
- "animation.Animation._repr_html_ = anim_to_html\n",
- "\n",
- "\n",
- "def init():\n",
- " # Set x and y range\n",
- " ax.set_ylim(0, box_size)\n",
- " ax.set_xlim(0, box_size)\n",
- " xdata, ydata = [], []\n",
- " part.set_data(xdata, ydata)\n",
- " return part,\n",
- "\n",
- "\n",
- "def run(i):\n",
+ "color = {0: \"#7fc454\"} \n",
+ "radii = {0: LJ_SIGMA/2}\n",
+ "vis = espressomd.zn.Visualizer(system, colors=color, radii=radii)\n",
+ "for _ in range(100):\n",
" system.integrator.run(50)\n",
- "\n",
- " # Save current system state as a plot\n",
- " xdata, ydata = particles.pos_folded[:, 0], particles.pos_folded[:, 1]\n",
- " ax.figure.canvas.draw()\n",
- " part.set_data(xdata, ydata)\n",
- " print(f'progress: {(i + 1):3.0f}%', end='\\r')\n",
- " return part,"
- ]
- },
- {
- "cell_type": "markdown",
- "id": "25b6af88",
- "metadata": {},
- "source": [
- "We now can start the sampling over the animation class of matplotlib"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "id": "8c9fe146",
- "metadata": {
- "scrolled": true
- },
- "outputs": [],
- "source": [
- "fig, ax = plt.subplots(figsize=(10, 10))\n",
- "part, = ax.plot([], [], 'o')\n",
- "\n",
- "animation.FuncAnimation(fig, run, frames=100, blit=True, interval=0, repeat=False, init_func=init)"
+ " vis.update()"
]
},
{
@@ -401,7 +335,7 @@
"id": "278cfc7f",
"metadata": {},
"source": [
- "In the visualization video we can see that the single chains break and connect to each other during time. Also some monomers are present which break from and connect to chains. If you want to have some more frames, i.e. a longer video, just adjust the frames parameter in FuncAnimation."
+ "In the visualization video we can see that the single chains break and connect to each other during time. Also some monomers are present which break from and connect to chains."
]
},
{
@@ -538,7 +472,8 @@
"pos = box_size * np.hstack((np.random.random((N_PART, 2)), np.zeros((N_PART, 1))))\n",
"\n",
"# Add particles\n",
- "particles = system.part.add(pos=pos, rotation=N_PART * [(True, True, True)], dip=dip, fix=N_PART * [(False, False, True)])\n",
+ "particles = system.part.add(pos=pos, rotation=N_PART * [(True, True, True)],\n",
+ " dip=dip, fix=N_PART * [(False, False, True)])\n",
"\n",
"# Remove overlap between particles by means of the steepest descent method\n",
"system.integrator.set_steepest_descent(f_max=0, gamma=0.1, max_displacement=0.05)\n",
@@ -870,7 +805,7 @@
"outputs": [],
"source": [
"def dL_dB(alpha):\n",
- " return (1. / (alpha**2.) - 1. / ((np.sinh(alpha))**2.)) * dipm / (KT)"
+ " return (1. / (alpha**2) - 1. / np.sinh(alpha)**2) * dipm / (KT)"
]
},
{
@@ -951,8 +886,8 @@
"plt.xlabel(r'$\\alpha$', fontsize=20)\n",
"plt.ylabel(r'$M^*$', fontsize=20)\n",
"plt.plot(x, L(x), label='Langevin function')\n",
- "plt.plot(x, magnetization_approx_perp(x), label='q2D approximation $\\perp$')\n",
- "plt.plot(x, magnetization_approx_para(x), label='q2D approximation $\\parallel$')\n",
+ "plt.plot(x, magnetization_approx_perp(x), label=r'q2D approximation $\\perp$')\n",
+ "plt.plot(x, magnetization_approx_para(x), label=r'q2D approximation $\\parallel$')\n",
"# < exercise >\n",
"plt.legend(fontsize=20)\n",
"plt.show()\n",
@@ -974,10 +909,10 @@
"plt.xlabel(r'$\\alpha$', fontsize=20)\n",
"plt.ylabel(r'$M^*$', fontsize=20)\n",
"plt.plot(x, L(x), label='Langevin function')\n",
- "plt.plot(x, magnetization_approx_perp(x), label='q2D approximation $\\perp$')\n",
- "plt.plot(x, magnetization_approx_para(x), label='q2D approximation $\\parallel$')\n",
- "plt.plot(alphas, magnetization_perp / N_PART, 'o', label='simulation results $\\perp$')\n",
- "plt.plot(alphas, magnetization_para / N_PART, 'o', label='simulation results $\\parallel$')\n",
+ "plt.plot(x, magnetization_approx_perp(x), label=r'q2D approximation $\\perp$')\n",
+ "plt.plot(x, magnetization_approx_para(x), label=r'q2D approximation $\\parallel$')\n",
+ "plt.plot(alphas, magnetization_perp / N_PART, 'o', label=r'simulation results $\\perp$')\n",
+ "plt.plot(alphas, magnetization_para / N_PART, 'o', label=r'simulation results $\\parallel$')\n",
"plt.legend(fontsize=20)\n",
"plt.show()"
]
diff --git a/doc/tutorials/ferrofluid/ferrofluid_part3.ipynb b/doc/tutorials/ferrofluid/ferrofluid_part3.ipynb
index b7f896f46a..b3be7c3c32 100644
--- a/doc/tutorials/ferrofluid/ferrofluid_part3.ipynb
+++ b/doc/tutorials/ferrofluid/ferrofluid_part3.ipynb
@@ -134,6 +134,7 @@
"source": [
"import espressomd\n",
"import espressomd.magnetostatics\n",
+ "import espressomd.zn\n",
"\n",
"espressomd.assert_features(['DIPOLES', 'DP3M', 'LENNARD_JONES'])\n",
"\n",
@@ -272,7 +273,7 @@
"execution_count": null,
"id": "89966f9b",
"metadata": {
- "scrolled": true
+ "scrolled": false
},
"outputs": [],
"source": [
@@ -333,11 +334,11 @@
"execution_count": null,
"id": "bd24f8eb",
"metadata": {
- "scrolled": true
+ "scrolled": false
},
"outputs": [],
"source": [
- "# initialize array for hold the sampled dipole moments\n",
+ "# initialize array that will hold the sampled dipole moments\n",
"dipms = np.full((loops, 3), np.nan)\n",
"\n",
"# sample dipole moment\n",
@@ -536,6 +537,14 @@
"As in **part II** we use the same system for every value of the Langevin parameter $\\alpha$. Thus we use that the system is already pre-equilibrated from the previous run so we save some equilibration time. For scientific purposes one would use a new system for every value for the Langevin parameter to ensure that the systems are independent and no correlation effects are measured. Also one would perform more than just one simulation for each value of $\\alpha$ to increase the precision of the results."
]
},
+ {
+ "cell_type": "markdown",
+ "id": "395a4fc0-baaf-48d6-8b88-b8ba42408bfe",
+ "metadata": {},
+ "source": [
+ "Additionally, we will visualize the sampling of the magnetization curve using ZnDraw. Note, that the magnetic field is increasing during the video."
+ ]
+ },
{
"cell_type": "code",
"execution_count": null,
@@ -579,7 +588,8 @@
" magnetizations[ndx] = np.mean(magn_temp)\n",
"\n",
" # remove constraint\n",
- " system.constraints.clear()"
+ " if not alpha == alphas[-1]:\n",
+ " system.constraints.clear()"
]
},
{
@@ -668,6 +678,39 @@
"At sufficiently high volume fractions $\\phi$ and dipolar interaction parameters $\\lambda$ these clusters can be so rigid, that simulations with normal methods are impossible as the relaxation times exceeds normal simulation times by far, resulting in strongly correlated configurations and thus measurements."
]
},
+ {
+ "cell_type": "markdown",
+ "id": "d2333fe6-c427-48b6-92f0-d946f4084b30",
+ "metadata": {},
+ "source": [
+ "## Visualization"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "7fd720f4-858e-4d2e-8d4a-23644ab996e3",
+ "metadata": {},
+ "source": [
+ "Finally, we use ZnDraw to visualize the evolving system for high $\\alpha$. We just need to continue running our system as it is already warmed up for high magnetic field."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "8e1d8ade-efa9-4d73-a328-27b0e63993b4",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# initializing zndraw visualizer\n",
+ "color = {0: \"#7fc454\"} \n",
+ "radii = {0: lj_sigma/2}\n",
+ "vis = espressomd.zn.Visualizer(system, colors=color, radii=radii)\n",
+ "\n",
+ "for _ in range(400):\n",
+ " system.integrator.run(5)\n",
+ " vis.update()"
+ ]
+ },
{
"cell_type": "markdown",
"id": "81e06b9b",
diff --git a/doc/tutorials/langevin_dynamics/langevin_dynamics.ipynb b/doc/tutorials/langevin_dynamics/langevin_dynamics.ipynb
index 804725d179..6c602341e1 100644
--- a/doc/tutorials/langevin_dynamics/langevin_dynamics.ipynb
+++ b/doc/tutorials/langevin_dynamics/langevin_dynamics.ipynb
@@ -159,14 +159,6 @@
"## 2. Simulating Brownian motion"
]
},
- {
- "cell_type": "markdown",
- "id": "c38456fa",
- "metadata": {},
- "source": [
- "We will simulate the diffusion of a single particle that is coupled to an implicit solvent."
- ]
- },
{
"cell_type": "code",
"execution_count": null,
@@ -181,6 +173,7 @@
"import espressomd\n",
"import espressomd.accumulators\n",
"import espressomd.observables\n",
+ "import espressomd.zn\n",
"\n",
"logging.basicConfig(level=logging.INFO, stream=sys.stdout)\n",
"\n",
@@ -188,15 +181,76 @@
"KT = 1.1\n",
"STEPS = 1000000\n",
"\n",
+ "gammas = [1.0, 2.0, 4.0, 10.0]"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "c38456fa",
+ "metadata": {},
+ "source": [
+ "We will simulate the diffusion of a single particle that is coupled to an implicit solvent. First, we have to setup our system."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "45a6ad79-fae4-4bd7-b041-012b193b83bf",
+ "metadata": {},
+ "outputs": [],
+ "source": [
"# System setup\n",
"system = espressomd.System(box_l=[16] * 3)\n",
"system.time_step = 0.01\n",
"system.cell_system.skin = 0.4\n",
"\n",
- "particle = system.part.add(pos=[0, 0, 0])\n",
+ "particle = system.part.add(pos=[0, 0, 0])"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "3d594176-bd64-4ce8-9e25-f0efefbf996c",
+ "metadata": {},
+ "source": [
+ "In order to get an idea of the Brownian motion, we will visualize our system evolving for a few steps with the ZnDraw visualizer:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "80788a7f-4bf2-4fa6-9292-5de8b2a926ff",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "system.thermostat.set_langevin(kT=KT, gamma=gammas[0], seed=42)\n",
"\n",
+ "# Setup ZnDraw visualizer\n",
+ "color = {0: \"#7fc454\"} \n",
+ "radii = {0: 1}\n",
+ "vis = espressomd.zn.Visualizer(system, colors=color, radii=radii)\n",
+ "\n",
+ "system.integrator.run(1000)\n",
+ "for _ in range(200):\n",
+ " system.integrator.run(50)\n",
+ " vis.update()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "933cea77-4bca-47a6-8a82-cc13d07ec673",
+ "metadata": {},
+ "source": [
+ "Now we can do the actual production run where we simulate for different friction coefficients and sample our observables:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "bf997a0a",
+ "metadata": {},
+ "outputs": [],
+ "source": [
"# Run for different friction coefficients\n",
- "gammas = [1.0, 2.0, 4.0, 10.0]\n",
"tau_results = []\n",
"msd_results = []\n",
"vacf_results = []\n",
diff --git a/doc/tutorials/raspberry_electrophoresis/raspberry_electrophoresis.ipynb b/doc/tutorials/raspberry_electrophoresis/raspberry_electrophoresis.ipynb
index 1d91441855..a72c35294e 100644
--- a/doc/tutorials/raspberry_electrophoresis/raspberry_electrophoresis.ipynb
+++ b/doc/tutorials/raspberry_electrophoresis/raspberry_electrophoresis.ipynb
@@ -71,6 +71,7 @@
"import espressomd.interactions\n",
"import espressomd.electrostatics\n",
"import espressomd.lb\n",
+ "import espressomd.zn\n",
"\n",
"import sys\n",
"import tqdm\n",
@@ -683,6 +684,56 @@
"system.thermostat.set_lb(LB_fluid=lb, seed=123, gamma=20.0)"
]
},
+ {
+ "cell_type": "markdown",
+ "id": "a7388457-8c1b-4839-bfc6-e34778ae9511",
+ "metadata": {},
+ "source": [
+ "## Visualization"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "f3deb299-de31-4901-8f65-0ff06112232c",
+ "metadata": {},
+ "source": [
+ "We use the ZnDraw visualizer and run our simulation for a few steps in order to get an idea what is happening in our system:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "c24f4ee4-2a6a-46ff-ba89-4f6c7f7fb227",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Initializing ZnDraw visualizer\n",
+ "colors = {\n",
+ " TYPE_CENTRAL: \"#fca000\", #orange\n",
+ " TYPE_SURFACE: \"#7fc454\", #green\n",
+ " TYPE_CATIONS: \"#ff0000\", #red\n",
+ " TYPE_ANIONS: \"#030ffc\" #blue\n",
+ " }\n",
+ "radii = {\n",
+ " TYPE_CENTRAL: sig_ss/2,\n",
+ " TYPE_SURFACE: sig_ss/2,\n",
+ " TYPE_CATIONS: sig_ss/2,\n",
+ " TYPE_ANIONS: sig_ss/2\n",
+ " }\n",
+ "arrows_config = {'colormap': [[-0.5, 0.9, 0.5], [-0.0, 0.9, 0.5]],\n",
+ " 'normalize': True,\n",
+ " 'colorrange': [0, 1],\n",
+ " 'scale_vector_thickness': True,\n",
+ " 'opacity': 1.0}\n",
+ "\n",
+ "lbfield = espressomd.zn.LBField(system, step_x=2, step_y=2, step_z=5, scale=0.3)\n",
+ "vis = espressomd.zn.Visualizer(system, colors=colors, radii=radii, folded=True, vector_field=lbfield)\n",
+ "\n",
+ "for _ in range(50):\n",
+ " system.integrator.run(10)\n",
+ " vis.update()"
+ ]
+ },
{
"cell_type": "markdown",
"id": "f848c501",
@@ -709,6 +760,7 @@
"with open('posVsTime.dat', 'w') as f: # file where the raspberry trajectory will be written to\n",
" for _ in tqdm.tqdm(range(num_iterations)):\n",
" system.integrator.run(num_steps_per_iteration)\n",
+ " vis.update()\n",
" pos = central_part.pos - initial_pos\n",
" f.write(f\"{system.time:.2f} {pos[0]:.4f} {pos[1]:.4f} {pos[2]:.4f}\\n\")\n",
"\n",
diff --git a/doc/tutorials/widom_insertion/widom_insertion.ipynb b/doc/tutorials/widom_insertion/widom_insertion.ipynb
index 04e9140afd..cf06a7c507 100644
--- a/doc/tutorials/widom_insertion/widom_insertion.ipynb
+++ b/doc/tutorials/widom_insertion/widom_insertion.ipynb
@@ -195,6 +195,7 @@
"Now we are ready to set up the system. Because we will need to rescale the simulation box, we will initially only set up the short-range interactions.\n",
"\n",
"**Exercise:**\n",
+ "\n",
"* Set up a system with box length $10\\,\\sigma$, integrator time step $\\Delta t=0.01 \\,\\tau$ ($\\tau$ is the Lennard-Jones timescale, i.e. equal to 1 in the employed unit system) and Verlet skin width of $\\delta = 0.4\\sigma$\n",
"* Add a total of ``N_ION_PAIRS`` ion pairs at random positions to the simulation box\n",
"* Add a WCA interaction (purely repulsive Lennard-Jones interaction) between the particles\n",
@@ -244,10 +245,12 @@
"metadata": {},
"source": [
"**Exercise:**\n",
+ "\n",
"* Implement a function ``system_setup(c_salt_SI)`` that takes the desired salt concentration ``c_salt_SI`` in SI-units (mol/l) as an argument and rescales the simulation box accordingly\n",
"* Afterwards, the function should also add the P3M actor to account for electrostatic interactions\n",
"\n",
"**Hints:**\n",
+ "\n",
"* Use the prefactor ``PREF`` defined above to convert the concentration from SI-units to the unit system employed in the simulation\n",
"* An accuracy of $10^{-3}$ for the P$^3$M should be enough for this exercise"
]
@@ -290,9 +293,11 @@
"metadata": {},
"source": [
"**Exercise:**\n",
+ "\n",
"* Implement a function ``warmup()`` that removes potential overlaps between particles with 10000 steps of the steepest descent integrator and performs a warmup integration with 100 loops of 100 simulation steps\n",
"\n",
"**Hints:**\n",
+ "\n",
"* keep in mind that after the overlaps have been removed, the integrator should be switched\n",
" back to Velocity-Verlet with the Langevin thermostat (using ``GAMMA`` and ``KT``)\n",
" in order to perform a warmup integration\n",
@@ -340,10 +345,11 @@
"metadata": {},
"source": [
"**Exercise:**\n",
- "* Add the functionality to perform Widom insertion moves. Make sure that you always insert an ion pair in order to conserve the system electroneutrality.\n",
"\n",
+ "* Add the functionality to perform Widom insertion moves. Make sure that you always insert an ion pair in order to conserve the system electroneutrality.\n",
"\n",
"**Hints:**\n",
+ "\n",
"* Refer to the documentation to find out how to set up Widom particle insertion ([Widom Insertion](https://espressomd.github.io/doc/reaction_methods.html#widom-insertion-for-homogeneous-systems))"
]
},
@@ -377,6 +383,7 @@
"metadata": {},
"source": [
"**Exercise:**\n",
+ "\n",
"* Implement a function ``calculate_excess_chemical_potential(n_md_steps, n_mc_steps, sample_size)``\n",
" that performs ``n_mc_steps`` Widom particle insertions every ``n_md_steps`` MD steps in a loop\n",
" that repeats ``sample_size`` times and returns the chemical potential and its error as a tuple via\n",
@@ -471,6 +478,7 @@
"a logarithmically scaled x-axis.\n",
"\n",
"**Exercise:**\n",
+ "\n",
"* Explain the observed behaviour qualitatively\n",
"* How do you expect the excess chemical potential to behave in the limit $c_{\\mathrm{salt}}\\rightarrow 0$ mol/l?"
]
@@ -524,9 +532,11 @@
"$B \\simeq 1.0$ and $C = 0.2$ in reduced units for NaCl in an aqueous solution at 25°C.\n",
"\n",
"**Exercise:**\n",
+ "\n",
"* Compare your simulation results with the Debye-Hückel theory and the Davies equation by plotting all three curves in a single plot\n",
"\n",
"**Hint:**\n",
+ "\n",
"* re-use the code from the previous figure\n",
"* create a logarithmic concentration range from $10^{-4}$ to $10^{0}$ mol/l with ``np.logspace(-4, 0.0, num=500, base=10)`` to plot the analytical solutions"
]
![](\"figures/Electro-Steric_Stabilization.jpg\")
\n",
- "
\n",
- "
\n",
- "