Simulation of a single freely-jointed chain

import unittest
import math
import numpy as np
import pandas as pd

params = {"num_monomers": 10000}
params["length"] = 3*params["num_monomers"]

script="""
MonteCarlo
Configuration cubic_side_length {length} particle_type0 chain.fstprt.swp
Potential Model IdealGas
ThermoParams beta 1 chemical_potential0 1
Metropolis
TrialAdd particle_type 0
Run until_num_particles 1
RemoveTrial name TrialAdd
TrialGrowFile grow_file grow_chain.fstprt.swp
set_variable trials_per 1
Log trials_per_write trials_per output_file chain.txt
Movie trials_per_write trials_per output_file chain.xyz
EndToEndDistance output_file end_to_end.txt
RadiusOfGyration output_file rg.txt
Run num_trials 1e2
""".format(**params)

def linear_single_site(n, filename):
    # write fstprt file
    with open (filename, "w") as f:
        f.write(
"# LAMMPS-inspired data file\n\n" + \
str(n) + " sites\n" + \
str(n-1) + " bonds\n\
\n\
1 site types\n\
1 bond types\n\
\n\
Site Properties\n\
\n\
0 epsilon 1 sigma 1.0 cutoff 3.0\n\
\n\
Bond Properties\n\
\n\
0 RigidBond length 1.0 delta 0.00001\n\
\n\
Sites\n\
\n\
")

        for i in range(n):
            f.write(str(i) + " 0 " + str(i) + ".0 0.0 0.0\n")
        f.write("\nBonds\n\n")
        for i in range(n - 1):
            f.write(str(i) + " 0 " + str(i) + " " + str(i + 1) + "\n")
        f.write("\n")

    # write TrialGrowFile
    with open('grow_' + filename, 'w') as f:
        f.write("TrialGrowFile\n\nparticle_type 0 regrow true site 0\n")
        for i in range(1, n):
            f.write("bond true mobile_site "+str(i)+" anchor_site "+str(i-1)+"\n")

linear_single_site(params["num_monomers"], "chain.fstprt.swp")

with open('script.txt', 'w') as file: file.write(script)
import subprocess
syscode = subprocess.call("../../../build/bin/fst < script.txt > script.log", shell=True, executable='/bin/bash')
with open('script.log', 'r') as file: print(file.read(), '\n', 'exit:', syscode)

# FEASST version: v0.24.2-6-gf3cf33172b-dirty-hwh/whatnow
MonteCarlo
Configuration cubic_side_length 30000 particle_type0 chain.fstprt.swp
Potential Model IdealGas
ThermoParams beta 1 chemical_potential0 1
Metropolis
TrialAdd particle_type 0
Run until_num_particles 1
# initializing random number generator with seed: 1702503390
RemoveTrial name TrialAdd
TrialGrowFile grow_file grow_chain.fstprt.swp
Log output_file chain.txt trials_per_write 1
Movie output_file chain.xyz trials_per_write 1
EndToEndDistance output_file end_to_end.txt
RadiusOfGyration output_file rg.txt
Run num_trials 1e2

 exit: 0

class TestFreelyJointedIdealChain(unittest.TestCase):
    def test(self):
        end_to_end=pd.read_csv('end_to_end.txt')
        self.assertAlmostEqual(end_to_end['average'][0], math.sqrt(params['num_monomers']), delta=20)
        rg=pd.read_csv('rg.txt')
        self.assertAlmostEqual(rg['average'][0], np.sqrt(params['num_monomers']/6), delta=200)
unittest.main(argv=[''], verbosity=2, exit=False)

test (__main__.TestFreelyJointedIdealChain) ... ok

----------------------------------------------------------------------
Ran 1 test in 0.003s

OK

<unittest.main.TestProgram at 0x7fc914528820>

Did this tutorial work as expected? Did you find any inconsistencies or have any comments? Please contact us. Any feedback is appreciated!