Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Not only querying, but also modifying system parameters #18

Merged
merged 14 commits into from
May 1, 2023
Merged
Original file line number Diff line number Diff line change
Expand Up @@ -5,14 +5,14 @@
"id": "3a4fdeed",
"metadata": {},
"source": [
"## Querying Charges and Other Parameters\n",
"## Querying and Modifying Charges and Other Parameters\n",
"\n",
"Sometimes you want to inspect the charges or other parameters of the particles or bonds in a System. Force field parameters are stored in the Force objects added to a System. As an example, let's load a PDB file and model it using the Amber14 force field."
]
},
{
"cell_type": "code",
"execution_count": 1,
"execution_count": 2,
"id": "9323aaae",
"metadata": {},
"outputs": [],
Expand All @@ -25,6 +25,28 @@
"system = forcefield.createSystem(pdb.topology)"
]
},
{
"attachments": {},
"cell_type": "markdown",
"id": "6d26c481",
"metadata": {},
"source": [
"Let's also create a `Simulation` with this system."
]
},
{
"cell_type": "code",
"execution_count": 11,
"id": "552481f3",
"metadata": {},
"outputs": [],
"source": [
"integrator = LangevinMiddleIntegrator(\n",
" 300 * unit.kelvin, 1 / unit.picosecond, 0.004 * unit.picoseconds)\n",
"simulation = Simulation(pdb.topology, system, integrator)\n",
"simulation.context.setPositions(pdb.positions)"
yuanqing-wang marked this conversation as resolved.
Show resolved Hide resolved
]
},
{
"cell_type": "markdown",
"id": "3c6c73f4",
Expand All @@ -35,7 +57,7 @@
},
{
"cell_type": "code",
"execution_count": 2,
"execution_count": 3,
"id": "1734190f",
"metadata": {},
"outputs": [],
Expand All @@ -53,7 +75,7 @@
},
{
"cell_type": "code",
"execution_count": 3,
"execution_count": 4,
"id": "bf38fc9a",
"metadata": {},
"outputs": [],
Expand All @@ -76,7 +98,7 @@
},
{
"cell_type": "code",
"execution_count": 4,
"execution_count": 6,
"id": "933928f9",
"metadata": {},
"outputs": [
Expand All @@ -98,12 +120,98 @@
" if particle1 == 0 or particle2 == 0:\n",
" print(f'Particles ({particle1}, {particle2}), length = {length}, k = {k}')"
]
},
{
"cell_type": "markdown",
"id": "50f5157c",
"metadata": {},
"source": [
"In some cases you may want to modify those parameters to make a bond behave differently from what the `ForceField` assigned.\n",
"We show you here that to modify these parameters are also easy."
]
},
{
"attachments": {},
"cell_type": "markdown",
"id": "cd5396c5",
"metadata": {},
"source": [
"Let's say we would want to modify the force constant of the bond connecting particle 1 and 0 to a new number."
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "086c7629",
"metadata": {},
"outputs": [],
"source": [
"for i in range(bonded.getNumBonds()):\n",
" particle1, particle2, length, k = bonded.getBondParameters(i)\n",
" if particle1 == 1 and particle2 == 0:\n",
" bonded.setBondParameters(\n",
" i, particle1, particle2, length, \n",
" 2666 * unit.kilojoules_per_mole/unit.nanometer**2\n",
" )"
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

No need to split this onto four lines. It makes the code less readable, and it can easily fit on a single line.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sure! If we don't worry about the 80 char limit.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If someone wants to view the entry on a VT52 terminal, I admit it will be inconvenient. But I think we can say VT100 is the earliest terminal we support. :)

]
},
{
"cell_type": "markdown",
"id": "8f9728ad",
"metadata": {},
"source": [
"Now if you query again you can see new parameters."
]
},
{
"cell_type": "code",
"execution_count": 9,
"id": "68498ff5",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Particles (4, 0), length = 0.1471 nm, k = 307105.5999999999 kJ/(nm**2 mol)\n",
"Particles (1, 0), length = 0.101 nm, k = 2666.0 kJ/(nm**2 mol)\n",
"Particles (2, 0), length = 0.101 nm, k = 363171.19999999995 kJ/(nm**2 mol)\n",
"Particles (3, 0), length = 0.101 nm, k = 363171.19999999995 kJ/(nm**2 mol)\n"
]
}
],
"source": [
"bonded = [f for f in system.getForces() if isinstance(f, HarmonicBondForce)][0]\n",
yuanqing-wang marked this conversation as resolved.
Show resolved Hide resolved
"for i in range(bonded.getNumBonds()):\n",
" particle1, particle2, length, k = bonded.getBondParameters(i)\n",
" if particle1 == 0 or particle2 == 0:\n",
" print(f'Particles ({particle1}, {particle2}), length = {length}, k = {k}')"
]
},
{
"attachments": {},
"cell_type": "markdown",
"id": "47729e2e",
"metadata": {},
"source": [
"Modifying `Force` objects will affect any new `Simulations` or `Contexts` you create, but it will have no effect on ones that already exist. \n",
"If you want your modifications to apply to an existing `Simulation`, you can copy them over by calling"
]
},
{
"cell_type": "code",
"execution_count": 12,
"id": "77075eed",
"metadata": {},
"outputs": [],
"source": [
"bonded.updateParametersInContext(simulation.context)"
]
}
],
"metadata": {
"tags": ["force field", "inspection", "forces"],
"kernelspec": {
"display_name": "Python 3",
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
Expand All @@ -117,8 +225,13 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.1"
}
"version": "3.11.3"
},
"tags": [
"force field",
"inspection",
"forces"
]
},
"nbformat": 4,
"nbformat_minor": 5
Expand Down