Skip to content

Commit f08bb02

Browse files
V-FrancoisV-Francois
authored
docs: otp example (#66)
**Summary** This PR adds a new example, to show how to recover results from Romain's paper on OTP.
1 parent 892f5b4 commit f08bb02

12 files changed

Lines changed: 499 additions & 0 deletions

File tree

‎examples/ortho-terphenyl/1-create-config-right-density/create-initial-config.py‎

Lines changed: 69 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,69 @@
1+
# /// script
2+
# requires-python = ">=3.13"
3+
# dependencies = [
4+
# "numpy>=2.4.1",
5+
# ]
6+
# ///
7+
import itertools
8+
9+
import numpy as np
10+
11+
M = 1000
12+
N = M * 3
13+
number_density = 0.2 # This is per particle, not per molecule
14+
# rho = N/V
15+
V = N / number_density
16+
L = V ** (1 / 3)
17+
18+
sigmas = [0.9, 1, 1.1]
19+
20+
21+
def main(output_file_name: str) -> None:
22+
"""Create a regular lattice of molecules. These planar molecules are in the x-y plane."""
23+
24+
number_in_each_direction = round(M ** (1 / 3))
25+
dxdydz = L / number_in_each_direction
26+
assert number_in_each_direction**3 == M, "M is not an int to power 3"
27+
28+
r_ab = (sigmas[0] + sigmas[1]) / 2
29+
r_ac = (sigmas[0] + sigmas[2]) / 2
30+
cos_alpha = np.cos(60 / 180)
31+
sin_alpha = np.sin(60 / 180)
32+
33+
f = open(output_file_name, "w")
34+
35+
f.write(f"{N}\n")
36+
f.write(f"columns:molecule,species,position cell:{L},{L},{L}\n")
37+
38+
counter = 1
39+
for i, j, k in itertools.product(
40+
range(number_in_each_direction),
41+
range(number_in_each_direction),
42+
range(number_in_each_direction),
43+
):
44+
r_a = (i * dxdydz, j * dxdydz, k * dxdydz)
45+
r_b = (i * dxdydz, j * dxdydz + r_ab, k * dxdydz)
46+
r_c = (i * dxdydz + r_ac * cos_alpha, j * dxdydz + r_ac * sin_alpha, k * dxdydz)
47+
48+
f.write(f"{counter} 1 {r_a[0]} {r_a[1]} {r_a[2]}\n")
49+
f.write(f"{counter} 2 {r_b[0]} {r_b[1]} {r_b[2]}\n")
50+
f.write(f"{counter} 3 {r_c[0]} {r_c[1]} {r_c[2]}\n")
51+
52+
counter += 1
53+
54+
number_of_bonds = M * 3
55+
f.write(f"{number_of_bonds}\n")
56+
f.write("columns:bond\n")
57+
for i in range(M):
58+
index_of_a = 1 + i * 3
59+
# AB
60+
f.write(f"{index_of_a} {index_of_a + 1}\n")
61+
# BC
62+
f.write(f"{index_of_a + 1} {index_of_a + 2}\n")
63+
# AC
64+
f.write(f"{index_of_a} {index_of_a + 2}\n")
65+
66+
67+
if __name__ == "__main__":
68+
output_file_name = "inputframe.xyz"
69+
main(output_file_name)

‎examples/ortho-terphenyl/1-create-config-right-density/params-template.toml‎

Lines changed: 82 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,82 @@
1+
[system]
2+
config = "./inputframe.xyz"
3+
temperature = 2.0
4+
density = DENSITY
5+
list_type = "LinkedList"
6+
7+
[model]
8+
9+
[model."1-1"]
10+
name = "GeneralKG"
11+
epsilon = 1.0
12+
sigma = 0.9
13+
k = 37.03703703703703
14+
r0 = 1.35
15+
16+
[model."1-2"]
17+
name = "GeneralKG"
18+
epsilon = 1.0
19+
sigma = 0.95
20+
k = 33.24099722991689
21+
r0 = 1.425
22+
23+
[model."1-3"]
24+
name = "GeneralKG"
25+
epsilon = 1.0
26+
sigma = 1.0
27+
k = 30.0
28+
r0 = 1.5
29+
30+
[model."2-2"]
31+
name = "GeneralKG"
32+
epsilon = 1.0
33+
sigma = 1.0
34+
k = 30.0
35+
r0 = 1.5
36+
37+
[model."2-3"]
38+
name = "GeneralKG"
39+
epsilon = 1.0
40+
sigma = 1.05
41+
k = 27.2108843537415
42+
r0 = 1.575
43+
44+
[model."3-3"]
45+
name = "GeneralKG"
46+
epsilon = 1.0
47+
sigma = 1.1
48+
k = 24.79338842975207
49+
r0 = 1.65
50+
51+
52+
53+
[simulation]
54+
type = "Metropolis"
55+
steps = 10000
56+
seed = 10
57+
parallel = false
58+
59+
[[simulation.move]]
60+
action = "Displacement"
61+
probability = 0.8
62+
policy = "SimpleGaussian"
63+
parameters = {sigma = 0.05}
64+
65+
[[simulation.move]]
66+
action = "MoleculeFlip"
67+
probability = 0.2
68+
policy = "DoubleUniform"
69+
70+
[[simulation.output]]
71+
algorithm = "StoreCallbacks"
72+
callbacks = ["energy", "acceptance"]
73+
scheduler_params = {linear_interval = 1}
74+
75+
[[simulation.output]]
76+
algorithm = "StoreLastFrames"
77+
scheduler_params = {linear_interval = 10000}
78+
fmt = "XYZ"
79+
80+
[[simulation.output]]
81+
algorithm = "PrintTimeSteps"
82+
scheduler_params = {linear_interval = 1}

‎examples/ortho-terphenyl/1-create-config-right-density/run.sh‎

Lines changed: 16 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,16 @@
1+
#!/bin/bash
2+
3+
uv run --script create-initial-config.py
4+
5+
density=0.2
6+
7+
while (( $(echo "$density < 1.2" | bc -l ) ))
8+
do
9+
echo "Density" $density
10+
11+
sed "s/DENSITY/$density/" params-template.toml > params.toml
12+
particlesmc params.toml
13+
cp trajectories/1/lastframe.xyz inputframe.xyz
14+
15+
density=$(echo "$density" | awk '{printf "%f", $1 * 1.1}')
16+
done

‎examples/ortho-terphenyl/2-equilibrate-at-different-temperatures/create-folders.sh‎

Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,7 @@
1+
#!/bin/bash
2+
3+
for temperature in 3.0 2.0 1.6 1.4 1.25 1.2 1.15 1.1 1.05 1.0
4+
do
5+
mkdir -p $temperature
6+
sed "s/TEMPERATURE/$temperature/" params-template.toml > $temperature/params.toml
7+
done

‎examples/ortho-terphenyl/2-equilibrate-at-different-temperatures/params-template.toml‎

Lines changed: 87 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,87 @@
1+
[system]
2+
config = "../../1-create-config-right-density/inputframe.xyz"
3+
temperature = TEMPERATURE
4+
density = 1.2
5+
list_type = "LinkedList"
6+
7+
[model]
8+
9+
[model."1-1"]
10+
name = "GeneralKG"
11+
epsilon = 1.0
12+
sigma = 0.9
13+
k = 37.03703703703703
14+
r0 = 1.35
15+
16+
[model."1-2"]
17+
name = "GeneralKG"
18+
epsilon = 1.0
19+
sigma = 0.95
20+
k = 33.24099722991689
21+
r0 = 1.425
22+
23+
[model."1-3"]
24+
name = "GeneralKG"
25+
epsilon = 1.0
26+
sigma = 1.0
27+
k = 30.0
28+
r0 = 1.5
29+
30+
[model."2-2"]
31+
name = "GeneralKG"
32+
epsilon = 1.0
33+
sigma = 1.0
34+
k = 30.0
35+
r0 = 1.5
36+
37+
[model."2-3"]
38+
name = "GeneralKG"
39+
epsilon = 1.0
40+
sigma = 1.05
41+
k = 27.2108843537415
42+
r0 = 1.575
43+
44+
[model."3-3"]
45+
name = "GeneralKG"
46+
epsilon = 1.0
47+
sigma = 1.1
48+
k = 24.79338842975207
49+
r0 = 1.65
50+
51+
52+
53+
[simulation]
54+
type = "Metropolis"
55+
steps = 5000000
56+
seed = 10
57+
parallel = false
58+
59+
[[simulation.move]]
60+
action = "Displacement"
61+
probability = 0.8
62+
policy = "SimpleGaussian"
63+
parameters = {sigma = 0.05}
64+
65+
[[simulation.move]]
66+
action = "MoleculeFlip"
67+
probability = 0.2
68+
policy = "DoubleUniform"
69+
70+
[[simulation.output]]
71+
algorithm = "StoreCallbacks"
72+
callbacks = ["energy", "acceptance"]
73+
scheduler_params = {linear_interval = 1000}
74+
75+
[[simulation.output]]
76+
algorithm = "StoreLastFrames"
77+
scheduler_params = {linear_interval = 10000}
78+
fmt = "XYZ"
79+
80+
[[simulation.output]]
81+
algorithm = "StoreTrajectories"
82+
scheduler_params = {linear_interval = 1000}
83+
fmt = "XYZ"
84+
85+
[[simulation.output]]
86+
algorithm = "PrintTimeSteps"
87+
scheduler_params = {linear_interval = 1000}

‎examples/ortho-terphenyl/3-run-production/create-folders.sh‎

Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,7 @@
1+
#!/bin/bash
2+
3+
for temperature in 3.0 2.0 1.6 1.4 1.25 1.2 1.15 1.1 1.05 1.0
4+
do
5+
mkdir -p $temperature
6+
sed "s/TEMPERATURE/$temperature/" params-template.toml > $temperature/params.toml
7+
done

‎examples/ortho-terphenyl/3-run-production/params-template.toml‎

Lines changed: 87 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,87 @@
1+
[system]
2+
config = "../../2-equilibrate-at-different-temperatures/TEMPERATURE/trajectories/1/lastframe.xyz"
3+
temperature = TEMPERATURE
4+
density = 1.2
5+
list_type = "LinkedList"
6+
7+
[model]
8+
9+
[model."1-1"]
10+
name = "GeneralKG"
11+
epsilon = 1.0
12+
sigma = 0.9
13+
k = 37.03703703703703
14+
r0 = 1.35
15+
16+
[model."1-2"]
17+
name = "GeneralKG"
18+
epsilon = 1.0
19+
sigma = 0.95
20+
k = 33.24099722991689
21+
r0 = 1.425
22+
23+
[model."1-3"]
24+
name = "GeneralKG"
25+
epsilon = 1.0
26+
sigma = 1.0
27+
k = 30.0
28+
r0 = 1.5
29+
30+
[model."2-2"]
31+
name = "GeneralKG"
32+
epsilon = 1.0
33+
sigma = 1.0
34+
k = 30.0
35+
r0 = 1.5
36+
37+
[model."2-3"]
38+
name = "GeneralKG"
39+
epsilon = 1.0
40+
sigma = 1.05
41+
k = 27.2108843537415
42+
r0 = 1.575
43+
44+
[model."3-3"]
45+
name = "GeneralKG"
46+
epsilon = 1.0
47+
sigma = 1.1
48+
k = 24.79338842975207
49+
r0 = 1.65
50+
51+
52+
53+
[simulation]
54+
type = "Metropolis"
55+
steps = 5000000
56+
seed = 10
57+
parallel = false
58+
59+
[[simulation.move]]
60+
action = "Displacement"
61+
probability = 0.8
62+
policy = "SimpleGaussian"
63+
parameters = {sigma = 0.05}
64+
65+
[[simulation.move]]
66+
action = "MoleculeFlip"
67+
probability = 0.2
68+
policy = "DoubleUniform"
69+
70+
[[simulation.output]]
71+
algorithm = "StoreCallbacks"
72+
callbacks = ["energy", "acceptance"]
73+
scheduler_params = {linear_interval = 1000}
74+
75+
[[simulation.output]]
76+
algorithm = "StoreLastFrames"
77+
scheduler_params = {linear_interval = 10000}
78+
fmt = "XYZ"
79+
80+
[[simulation.output]]
81+
algorithm = "StoreTrajectories"
82+
scheduler_params = {linear_interval = 4096, log_base=2.0}
83+
fmt = "XYZ"
84+
85+
[[simulation.output]]
86+
algorithm = "PrintTimeSteps"
87+
scheduler_params = {linear_interval = 1000}

0 commit comments

Comments
 (0)