Skip to content

Commit

Permalink
ZnDraw Visualization added
Browse files Browse the repository at this point in the history
  • Loading branch information
@keerthirk1995
keerthirk1995 committed Aug 12, 2024
1 parent 82c962d commit 35e0617
Show file tree
Hide file tree
Showing 4 changed files with 693 additions and 268 deletions.
618 changes: 556 additions & 62 deletions doc/tutorials/charged_system/charged_system.ipynb
View file

Large diffs are not rendered by default.

245 changes: 52 additions & 193 deletions doc/tutorials/lattice_boltzmann/lattice_boltzmann_sedimentation.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -79,6 +79,7 @@
"import espressomd.shapes\n",
"import espressomd.observables\n",
"import espressomd.accumulators\n",
"import espressomd.zn\n",
"\n",
"espressomd.assert_features([\"LENNARD_JONES\", \"WALBERLA\"])\n",
"\n",
Expand Down Expand Up @@ -266,21 +267,47 @@
"system.thermostat.set_langevin(kT=0., gamma=15., seed=12)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "02a4a0c2-e468-4723-b285-78868bd4bde7",
"metadata": {},
"outputs": [],
"source": [
"parts = system.part.add(pos=initial_positions, ext_force=[f_gravity] * n_parts)\n"
]
},
{
"cell_type": "markdown",
"id": "d90aa40c-a974-41fe-a069-18939ff0c6e1",
"metadata": {},
"source": [
"### Visualization: Langevin Dynamics"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "4fd712ff",
"metadata": {},
"outputs": [],
"source": []
"source": [
"color = {0: \"#00f0f0\"}\n",
"radii = {0: lj_sigma/2.0}\n",
"\n",
"print(\"Particle dynamics in the absence of hydrodynamics\")\n",
"vis = espressomd.zn.Visualizer(system, colors=color, radii=radii, folded=True)\n",
"\n",
"# Visualize constraints\n",
"#vis.draw_constraints([wall_shape_b, wall_shape_t])"
]
},
{
"cell_type": "markdown",
"id": "f093b544",
"metadata": {},
"source": [
"We can now sample the particle positions as a function of time.\n",
"We will use a particle observable and an accumulator to record the trajectory."
"We can now sample the particle positions as a function of time.\n"
]
},
{
Expand All @@ -290,17 +317,9 @@
"metadata": {},
"outputs": [],
"source": [
"parts = system.part.add(pos=initial_positions, ext_force=[f_gravity] * n_parts)\n",
"obs_particle_pos = espressomd.observables.ParticlePositions(ids=list(range(n_parts)))\n",
"acc_particle_pos = espressomd.accumulators.TimeSeries(obs=obs_particle_pos, delta_N=1)\n",
"\n",
"system.integrator.run(0)\n",
"\n",
"for step in tqdm.tqdm(range(sampling_steps)):\n",
" acc_particle_pos.update()\n",
" system.integrator.run(25)\n",
"\n",
"data_ld = np.remainder(np.reshape(acc_particle_pos.time_series(), (-1, n_parts, 3)), system.box_l)"
" vis.update()"
]
},
{
Expand Down Expand Up @@ -396,71 +415,41 @@
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "b05ea198",
"cell_type": "markdown",
"id": "8db1ba81-15a0-4a73-b0ac-9910a307256e",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"id": "b2a4deaf",
"id": "4e73ce63-0d94-43db-8395-ea5264ad1530",
"metadata": {},
"source": [
"We will plot the fluid flow field in the final video using 2D vectors.\n",
"To this end, we need to record the fluid trajectory with the same frequency as the particle positions.\n",
"\n",
"**Exercise:**\n",
"* create a LB velocity profile in Cartesian coordinates\n",
" ([user guide](https://espressomd.github.io/doc/espressomd.html#espressomd.observables.LBVelocityProfile))\n",
"* register that observable in a time series accumulator named ``acc_lb_vel``\n",
" ([user guide](https://espressomd.github.io/doc/analysis.html#time-series))\n",
"### Visualization with hydrodynamics\n",
"\n",
"**Hints:**\n",
"* the velocity observable takes parameters in MD units, not LB units\n",
"* there is no fluid inside the top an bottom boundaries, therefore the number of bins for the *y*-axis\n",
" is smaller than ``n_height`` by 2 units (parameters ``min_y`` and ``max_y`` are also affected)\n",
"* use ``n_z_bins=1`` to average the velocity along the *z*-direction\n",
"* for ``sampling_delta_*`` and ``sampling_offset_*``, use ``spacing`` and ``0.5 * spacing`` respectively"
"Here, we will visualize the particle and the fluid dynamics with hydrodynamics on. For that we will use the inbuilt ZnDraw feature in Espresso. The fluid is visualized using a 2D vector flow field. "
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "a18ee853",
"id": "8f4c3326",
"metadata": {},
"outputs": [],
"source": [
"# SOLUTION CELL\n",
"obs_lb_vel = espressomd.observables.LBVelocityProfile(\n",
" n_x_bins=n_width,\n",
" n_y_bins=n_height - 2, # skip data inside the LB boundaries (top and bottom walls)\n",
" n_z_bins=1, # averaged velocity along the z-direction\n",
" min_x=0.0,\n",
" min_y=spacing,\n",
" min_z=0.0,\n",
" max_x=system.box_l[0],\n",
" max_y=system.box_l[1] - spacing,\n",
" max_z=system.box_l[2],\n",
" sampling_delta_x=spacing,\n",
" sampling_delta_y=spacing,\n",
" sampling_delta_z=spacing,\n",
" sampling_offset_x=0.5 * spacing,\n",
" sampling_offset_y=0.5 * spacing,\n",
" sampling_offset_z=0.5 * spacing,\n",
" allow_empty_bins=True)\n",
"acc_lb_vel = espressomd.accumulators.TimeSeries(obs=obs_lb_vel, delta_N=1)"
"print(\"Particle dynamics in the presence of hydrodynamics\")\n",
"color = {0: \"#00f0f0\"}\n",
"radii = {0: lj_sigma/2.0}\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=1)\n",
"vis1 = espressomd.zn.Visualizer(system, colors=color, radii=radii, folded=True, vector_field=lbfield)\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "8f4c3326",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"id": "31e5cc87",
Expand All @@ -476,146 +465,16 @@
"metadata": {},
"outputs": [],
"source": [
"obs_particle_pos = espressomd.observables.ParticlePositions(ids=list(range(n_parts)))\n",
"acc_particle_pos = espressomd.accumulators.TimeSeries(obs=obs_particle_pos, delta_N=1)\n",
"\n",
"for step in tqdm.tqdm(range(sampling_steps)):\n",
" acc_lb_vel.update()\n",
" acc_particle_pos.update()\n",
" system.integrator.run(25)\n",
"\n",
"data_lb = np.remainder(np.reshape(acc_particle_pos.time_series(), (-1, n_parts, 3)), system.box_l)\n",
"data_flowfield = acc_lb_vel.time_series()[:, :, :, 0, 0:2]"
]
},
{
"cell_type": "markdown",
"id": "f225cd58",
"metadata": {},
"source": [
"## Visualization\n",
"\n",
"Now let's visualize the data. First some imports and definitions for inline visualization."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "db083422",
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.animation as animation\n",
"import matplotlib.quiver\n",
"import tempfile\n",
"import base64\n",
"\n",
"VIDEO_TAG = \"\"\"<video controls>\n",
" <source src=\"data:video/x-m4v;base64,{0}\" type=\"video/mp4\">\n",
" Your browser does not support the video tag.\n",
"</video>\"\"\"\n",
"\n",
"# set ignore 'divide' and 'invalid' errors\n",
"# these occur when plotting the flowfield containing a zero velocity\n",
"np.seterr(divide='ignore', invalid='ignore')\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",
"animation.Animation._repr_html_ = anim_to_html"
]
},
{
"cell_type": "markdown",
"id": "7fefebcb",
"metadata": {},
"source": [
"And now the actual visualization code.\n",
"\n",
"Note: if Jupyter encapsulates the video in a scrollable area, click on the ``Out[]:`` text to the left of\n",
"the output cell to toggle auto-scrolling off. This can also be achieved by hitting the key combination\n",
"Shift+O after highlighting the output cell."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "d67242dd",
"metadata": {
"scrolled": false
},
"outputs": [],
"source": [
"# setup figure and prepare axes\n",
"fig = plt.figure(figsize=(2 * 5, 5 / box_width * box_height))\n",
"gs = fig.add_gridspec(1, 2, wspace=0.1)\n",
"ax1 = plt.subplot(gs[0])\n",
"ax2 = plt.subplot(gs[1], sharey=ax1)\n",
"\n",
"ax1.set_title(\"Langevin\")\n",
"ax1.set_xlim((0, box_width))\n",
"ax1.set_ylim((0, box_height))\n",
"\n",
"ax2.set_title(\"LB Fluid\")\n",
"ax2.set_xlim((0, box_width))\n",
"ax2.set_ylim((0, box_height))\n",
"\n",
"# draw walls\n",
"for ax in [ax1, ax2]:\n",
" ax.hlines((1, box_height-1), 0, box_width, color=\"gray\")\n",
"\n",
"# create meshgrid for quiver plot\n",
"xs = np.array([x for x in range(n_width)]) * spacing\n",
"ys = np.array([y for y in range(1, n_height-1)]) * spacing\n",
"X, Y = np.meshgrid(xs, ys)\n",
"\n",
"# create a transposed flow field for quiver plot\n",
"data_flowfield_t = np.transpose(data_flowfield, axes=(0, 2, 1, 3))\n",
"\n",
"# set quiver scale (fraction of the highest velocity in the XY plane)\n",
"lb_vel_max = np.sum(np.square(data_flowfield), axis=-1)\n",
"quiver_scale = np.sqrt(np.max(lb_vel_max))\n",
"\n",
"def plot_lb_vel(ax, X, Y, flowfield, t, scale):\n",
" return ax.quiver(X, Y,\n",
" flowfield[t, :, :, 0],\n",
" flowfield[t, :, :, 1],\n",
" scale_units=\"xy\", scale=scale)\n",
"\n",
"# initialize plot objects\n",
"lb_ff = plot_lb_vel(ax2, X, Y, data_flowfield_t, 0, quiver_scale)\n",
"lb_particles, = ax2.plot([], [], 'o')\n",
"ld_particles, = ax1.plot([], [], 'o')\n",
"\n",
"def draw_frame(t):\n",
" # manually remove Quivers from ax2\n",
" for artist in ax2.get_children():\n",
" if isinstance(artist, matplotlib.quiver.Quiver):\n",
" artist.remove()\n",
"\n",
" # draw new quivers\n",
" lb_ff = plot_lb_vel(ax2, X, Y, data_flowfield_t, t, quiver_scale)\n",
"\n",
" # draw particles\n",
" ld_particles.set_data(data_ld[t, :, 0], data_ld[t, :, 1])\n",
" lb_particles.set_data(data_lb[t, :, 0], data_lb[t, :, 1])\n",
"\n",
" return [ld_particles, lb_particles, lb_ff]\n",
"\n",
"animation.FuncAnimation(fig, draw_frame, frames=sampling_steps, blit=True, interval=0, repeat=False)"
" vis1.update()\n",
"\n"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
Expand All @@ -629,7 +488,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.12"
"version": "3.10.12"
}
},
"nbformat": 4,
Expand Down
24 changes: 21 additions & 3 deletions doc/tutorials/lennard_jones/lennard_jones.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -210,6 +210,8 @@
"import espressomd.observables\n",
"import espressomd.accumulators\n",
"import espressomd.analyze\n",
"import espressomd.zn \n",
"\n",
"required_features = [\"LENNARD_JONES\"]\n",
"espressomd.assert_features(required_features)"
]
Expand Down Expand Up @@ -632,13 +634,28 @@
"system.thermostat.set_langevin(kT=TEMPERATURE, gamma=GAMMA, seed=42)"
]
},
{
"cell_type": "markdown",
"id": "10510b16-bc19-4c16-a881-a3b88cc298f6",
"metadata": {},
"source": [
"### System visualization "
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "42c573df",
"metadata": {},
"outputs": [],
"source": []
"source": [
"# Visualizer Parameters\n",
"color = {0: \"#00f0f0\"} # Particle color by type\n",
"radii = {0: LJ_SIG/2.0} # Particle size by type\n",
"\n",
"# Initializing Visualizer\n",
"vis = espressomd.zn.Visualizer(system, colors=color, radii=radii)"
]
},
{
"cell_type": "markdown",
Expand Down Expand Up @@ -696,6 +713,7 @@
" e_total[i] = energy['total']\n",
" e_kin[i] = energy['kinetic']\n",
" system.integrator.run(STEPS_PER_SAMPLE)\n",
" vis.update()\n",
"T_inst = 2. / 3. * e_kin / N_PART"
]
},
Expand Down Expand Up @@ -1235,7 +1253,7 @@
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
Expand All @@ -1249,7 +1267,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.10"
"version": "3.10.12"
}
},
"nbformat": 4,
Expand Down
Loading

0 comments on commit 35e0617

Please sign in to comment.