Slab, cylindrical, spherical and mixed confinementΒΆ

[1]:

import unittest
import feasst as fst

def shape(confinement):
    if confinement == "slab":
        shape = fst.MakeSlab(fst.args({"dimension": "2", "bound0": "-1", "bound1": "1"}))
    elif confinement == "cylinder":
        shape = fst.MakeCylinder(
            fst.args({"radius": "4"}),
            fst.Position(fst.args({"x": "0", "y": "0", "z": "0"})),
            fst.Position(fst.args({"x": "0", "y": "0", "z": "1"})))
    elif confinement == "sphere":
        shape = fst.MakeSphere(
            fst.args({"radius": "4"}),
            fst.Position(fst.args({"x": "0", "y": "0", "z": "0"})))
    elif confinement == "union":
        shape = fst.MakeShapeUnion(
            fst.MakeSphere(
                fst.args({"radius": "2.5"}),
                fst.Position(fst.args({"x": "0", "y": "0", "z": "0"}))),
            fst.MakeSlab(fst.args({"dimension": "2", "bound0": "-1", "bound1": "1"})))
    elif confinement == "network":
        shape = fst.MakeSphere(
            fst.args({"radius": "5"}),
            fst.Position(fst.args({"x": "0", "y": "0", "z": "0"})))
        shape = fst.MakeShapeUnion(shape, fst.MakeCylinder(
            fst.args({"radius": "2"}),
            fst.Position(fst.args({"x": "0", "y": "0", "z": "0"})),
            fst.Position(fst.args({"x": "1", "y": "0", "z": "0"}))))
        shape = fst.MakeShapeUnion(shape, fst.MakeCylinder(
            fst.args({"radius": "2"}),
            fst.Position(fst.args({"x": "0", "y": "0", "z": "0"})),
            fst.Position(fst.args({"x": "0", "y": "1", "z": "0"}))))
        shape = fst.MakeShapeUnion(shape, fst.MakeCylinder(
            fst.args({"radius": "2"}),
            fst.Position(fst.args({"x": "0", "y": "0", "z": "0"})),
            fst.Position(fst.args({"x": "0", "y": "0", "z": "1"}))))
    else:
        assert False # unrecognized
    return shape

def run(confinement, cubic_box_length, num_particles):
    monte_carlo = fst.MonteCarlo()
    monte_carlo.set(fst.lennard_jones(fst.args({"cubic_box_length": str(cubic_box_length)})))
    monte_carlo.add(fst.MakePotential(fst.MakeModelHardShape(shape(confinement))))
    monte_carlo.set(fst.MakeThermoParams(fst.args(
        {"beta": "1.5", "chemical_potential": "1."})))
    monte_carlo.set(fst.MakeMetropolis())
    monte_carlo.add(fst.MakeTrialTranslate(fst.args(
        {"weight": "1.", "tunable_param": "2."})))
    fst.SeekNumParticles(num_particles).with_trial_add().run(monte_carlo)
    monte_carlo.add(fst.MakeLogAndMovie(fst.args({"steps_per": str(int(1e4)), "file_name": confinement})))
    monte_carlo.add(fst.MakeCheckEnergyAndTune(fst.args({"steps_per": str(int(1e4))})))
    monte_carlo.attempt(int(1e5))

class TestConfinement1LJ(unittest.TestCase):
    def test(self):
        for confinement in ["slab", "cylinder", "sphere", "union", "network"]:
            cubic_box_length = 8
            num_particles = 50
            if confinement == 'network':
                cubic_box_length=20
                num_particles=500
            run(confinement=confinement,
                cubic_box_length=cubic_box_length,
                num_particles=num_particles)

[2]:

%%time
unittest.main(argv=[''], verbosity=2, exit=False)

test (__main__.TestConfinement1LJ) ...

CPU times: user 14.2 s, sys: 16 ms, total: 14.2 s
Wall time: 14.2 s

ok

----------------------------------------------------------------------
Ran 1 test in 14.167s

OK

[2]:

<unittest.main.TestProgram at 0x7f2530f4a390>

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