Skip to content

Commit

Permalink
Add Zndraw in tutorials
Browse files Browse the repository at this point in the history
  • Loading branch information
RebeccaStephan authored and jngrad committed Aug 23, 2024
1 parent 3d5ca31 commit cfa136a
Show file tree
Hide file tree
Showing 7 changed files with 277 additions and 263 deletions.
44 changes: 41 additions & 3 deletions doc/tutorials/constant_pH/constant_pH.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -153,6 +153,7 @@
"import espressomd.electrostatics\n",
"import espressomd.reaction_methods\n",
"import espressomd.polymer\n",
"import espressomd.zn\n",
"from espressomd.interactions import HarmonicBond"
]
},
Expand Down Expand Up @@ -569,7 +570,7 @@
"\n",
"# add thermostat and short integration to let the system relax further\n",
"system.thermostat.set_langevin(kT=KT_REDUCED, gamma=1.0, seed=7)\n",
"system.integrator.run(steps=1000)\n",
"system.integrator.run(1000)\n",
"\n",
"if USE_ELECTROSTATICS:\n",
" COULOMB_PREFACTOR=BJERRUM_LENGTH_REDUCED * KT_REDUCED\n",
Expand Down Expand Up @@ -739,6 +740,40 @@
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"id": "cd417295-bab3-46a5-a574-fc60d06371a5",
"metadata": {},
"source": [
"We will initialize the zndraw visualizer to visualize the simulation for increasing pH value:"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "04dd303c",
"metadata": {
"scrolled": true
},
"outputs": [],
"source": [
"color = {TYPES[\"HA\"]: \"#7fc454\", #green\n",
" TYPES[\"A\"]: \"#225204\", #dark green\n",
" TYPES[\"B\"]: \"#fca000\", #orange\n",
" TYPES[\"Na\"]: \"#ff0000\", #red\n",
" TYPES[\"Cl\"]: \"#030ffc\" #blue\n",
" }\n",
"radii = {TYPES[\"HA\"]: 2,\n",
" TYPES[\"A\"]: 2,\n",
" TYPES[\"B\"]: 2,\n",
" TYPES[\"Na\"]: 2,\n",
" TYPES[\"Cl\"]: 2\n",
" }\n",
"\n",
"vis = espressomd.zn.Visualizer(system, colors=color, radii=radii)\n",
"vis.update()"
]
},
{
"cell_type": "markdown",
"id": "07daa583",
Expand Down Expand Up @@ -778,11 +813,14 @@
"outputs": [],
"source": [
"# SOLUTION CELL\n",
"def perform_sampling(type_A, num_samples, num_As:np.ndarray, reaction_steps, \n",
"def perform_sampling(type_A, num_samples, num_As:np.ndarray, reaction_steps,\n",
" prob_integration=0.5, integration_steps=1000):\n",
" for i in range(num_samples):\n",
" if USE_WCA and np.random.random() < prob_integration:\n",
" system.integrator.run(integration_steps)\n",
" for _ in range(integration_steps):\n",
" system.integrator.run(1)\n",
" global vis\n",
" vis.update()\n",
" # we should do at least one reaction attempt per reactive particle\n",
" RE.reaction(steps=reaction_steps)\n",
" num_As[i] = system.number_of_particles(type=type_A)"
Expand Down
66 changes: 52 additions & 14 deletions doc/tutorials/electrodes/electrodes_part2.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -214,6 +214,7 @@
"import espressomd.observables\n",
"import espressomd.accumulators\n",
"import espressomd.shapes\n",
"import espressomd.zn\n",
"\n",
"espressomd.assert_features(['WCA', 'ELECTROSTATICS'])\n",
"rng = np.random.default_rng(42)\n",
Expand Down Expand Up @@ -476,20 +477,20 @@
"source": [
"# SOLUTION CELL\n",
"offset = LJ_SIGMA # avoid unfavorable overlap at close distance to the walls\n",
"init_part_btw_z1 = offset \n",
"init_part_btw_z1 = offset\n",
"init_part_btw_z2 = box_l_z - offset\n",
"ion_pos = np.empty((3), dtype=float)\n",
"ion_pos = np.empty((3,), dtype=float)\n",
"\n",
"for i in range (N_IONPAIRS):\n",
" ion_pos[0] = rng.random(1) * system.box_l[0]\n",
" ion_pos[1] = rng.random(1) * system.box_l[1]\n",
" ion_pos[2] = rng.random(1) * (init_part_btw_z2 - init_part_btw_z1) + init_part_btw_z1\n",
"for i in range(N_IONPAIRS):\n",
" ion_pos[0] = rng.random(1)[0] * system.box_l[0]\n",
" ion_pos[1] = rng.random(1)[0] * system.box_l[1]\n",
" ion_pos[2] = rng.random(1)[0] * (init_part_btw_z2 - init_part_btw_z1) + init_part_btw_z1\n",
" system.part.add(pos=ion_pos, type=types[\"Cation\"], q=charges[\"Cation\"])\n",
" \n",
"for i in range (N_IONPAIRS):\n",
" ion_pos[0] = rng.random(1) * system.box_l[0]\n",
" ion_pos[1] = rng.random(1) * system.box_l[1]\n",
" ion_pos[2] = rng.random(1) * (init_part_btw_z2 - init_part_btw_z1) + init_part_btw_z1\n",
"\n",
"for i in range(N_IONPAIRS):\n",
" ion_pos[0] = rng.random(1)[0] * system.box_l[0]\n",
" ion_pos[1] = rng.random(1)[0] * system.box_l[1]\n",
" ion_pos[2] = rng.random(1)[0] * (init_part_btw_z2 - init_part_btw_z1) + init_part_btw_z1\n",
" system.part.add(pos=ion_pos, type=types[\"Anion\"], q=charges[\"Anion\"])"
]
},
Expand Down Expand Up @@ -1008,7 +1009,40 @@
"\n",
"With the above knowledge, we can now assess the \n",
"differential capacitance of the system, by changing the applied voltage\n",
"difference and determining the corresponding surface charge density."
"difference and determining the corresponding surface charge density.\n",
"We will use the ZnDraw visualizer to visualize our system:"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "ef5d908f-a7f0-4845-9555-34822128c141",
"metadata": {},
"outputs": [],
"source": [
"color = {\n",
" types[\"Cation\"]: \"#ff0000\", #red\n",
" types[\"Anion\"]: \"#030ffc\" #blue\n",
" }\n",
" \n",
"radii = {\n",
" types[\"Cation\"]: LJ_SIGMA,\n",
" types[\"Anion\"]: LJ_SIGMA\n",
" }\n",
"\n",
"vis = espressomd.zn.Visualizer(system, colors=color, radii=radii)\n",
"#vis.draw_constraints([floor, ceiling])\n",
"\n",
"#note: The system displayed is the enlarged system as you can see the ELG-GAP along the nonperiodic direction. The particles are only in the smaller subsystem.\n",
"#note: The particles are shown bigger for visualization purpose."
]
},
{
"cell_type": "markdown",
"id": "21da58e8-e666-4d06-8538-a8d6af521148",
"metadata": {},
"source": [
"Do the sampling from high to low potential:"
]
},
{
Expand All @@ -1029,7 +1063,9 @@
"for potential_diff in tqdm.tqdm(np.linspace(MIN_PHI, MAX_PHI, N_PHI)[::-1]):\n",
" system.auto_update_accumulators.clear()\n",
" system.electrostatics.solver = setup_electrostatic_solver(potential_diff)\n",
" system.integrator.run(N_SAMPLES_EQUIL_CAP * STEPS_PER_SAMPLE)\n",
" for i in range(N_SAMPLES_EQUIL_CAP * STEPS_PER_SAMPLE//50):\n",
" system.integrator.run(50)\n",
" vis.update()\n",
" sigmas = []\n",
"\n",
" for tm in range(N_SAMPLES_CAP):\n",
Expand All @@ -1039,7 +1075,9 @@
" system.auto_update_accumulators.clear()\n",
" system.auto_update_accumulators.add(density_accumulator_cation)\n",
" system.auto_update_accumulators.add(density_accumulator_anion)\n",
" system.integrator.run(STEPS_PER_SAMPLE)\n",
" for j in range(int(STEPS_PER_SAMPLE/50)):\n",
" system.integrator.run(50)\n",
" vis.update()\n",
"\n",
" cation_profile_mean = density_accumulator_cation.mean()[0, 0, :]\n",
" anion_profile_mean = density_accumulator_anion.mean()[0, 0, :]\n",
Expand Down
Loading

0 comments on commit cfa136a

Please sign in to comment.