Using the 'or' syntax can evaluate to False (if parameter is set to 0 for example). A more correct implementation for default values is to use kwargs.get('keyword', default_value) in place of kwargs.get('keyword') or default_value. This is found all throughout the code currently.

Dear Pysimm Developer,
Thank you so much for your contributions to Pysimm, which is a very useful tool in the MD simulation area.
I have one question about the atom-type order of data that generated by pysimm. Actually, I have many frames of my system in the Materials Studio platform, and I randomly choose two frames in one trajectory and save them as .mol2 file. But when I convert these two .mol2 file to lammps data file with pysimm, the atom type sequence is different in the generated two data files, like the following picture.

And my command is just like this:
from pysimm import system, lmps, forcefield
p = system.read_mol2('frame1.mol2')
p.apply_forcefield(forcefield.dreiding.Dreiding(),charges='gasteiger' )
p.write_lammps('frame2.data')

Do you know how to solve this problem? Many thanks in advance.
Best,
Yang

When updating box dimensions in System.read_lammpstrj(), there is no check to see if the current frame is correct. This results in the final dimensions always being the dimensions of the last frame in the file.

The random_walk module is great, but when building a polymer composed of over 1000 monomer, it is getting slower and slower. If there is another rude approach with constant last_backbone_vector, so I can chose the random or rude way as I need, that would be great :)

However, simply repeat the random_walk process already helps me out, thank you!

PySIMM: reading lammps data file "pim_unit.lmps"
PySIMM: read masses for 11 ParticleTypes
PySIMM: read "lj/cut" nonbonded parameters for 11 ParticleTypes
PySIMM: read "harmonic" bond parameters for 14 BondTypes
PySIMM: read "harmonic" angle parameters for 25 AngleTypes
PySIMM: read "harmonic" dihedral parameters for 38 DihedralTypes
PySIMM: read 35 particles
PySIMM: read 40 bonds
PySIMM: read 66 angles
PySIMM: read 97 dihedrals
PySIMM: Molecule 1 inserted
02:59:25: 1/5 monomers added
02:59:25: 2/5 monomers added
02:59:25: starting relax_002 LAMMPS simulation
02:59:36: relax_002 simulation using LAMMPS successful
02:59:36: 3/5 monomers added
02:59:36: starting relax_003 LAMMPS simulation
02:59:49: relax_003 simulation using LAMMPS successful
02:59:49: 4/5 monomers added
02:59:50: starting relax_004 LAMMPS simulation
Traceback (most recent call last):
File "/mnt/c/Users/sptiw/Desktop/github_desktop_downloads/ubuntu_softwares/pysimm/pysimm/lmps.py", line 955, in run
call_lammps(self, np, nanohub, prefix=prefix)
File "/mnt/c/Users/sptiw/Desktop/github_desktop_downloads/ubuntu_softwares/pysimm/pysimm/lmps.py", line 1045, in call_lammps
simulation.system.read_lammps_dump('pysimm.dump.tmp')
File "/mnt/c/Users/sptiw/Desktop/github_desktop_downloads/ubuntu_softwares/pysimm/pysimm/system.py", line 1741, in read_lammps_dump
with open(fname) as f:
FileNotFoundError: [Errno 2] No such file or directory: 'pysimm.dump.tmp'

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
File "run.py", line 116, in
sngl_chain = random_walk(mloc, 5, density=0.01,
File "/mnt/c/Users/sptiw/Desktop/github_desktop_downloads/ubuntu_softwares/pysimm/pysimm/apps/random_walk.py", line 436, in random_walk
sim.run(np=settings.get('np'))
File "/mnt/c/Users/sptiw/Desktop/github_desktop_downloads/ubuntu_softwares/pysimm/pysimm/lmps.py", line 957, in run
raise PysimmError('There was a problem calling LAMMPS with {}'.format(prefix))
pysimm.utils.PysimmError: There was a problem calling LAMMPS with mpiexec

hi, i was able to install pysimm but in the methanol example in the lmps.quick_min () part I get an error this PysimmError: A pair_style must be defined during initialization, could it help me

Right now check_lmps_exec() is called upon import of the lmps module. The module can be used without a working LAMMPS executable in theory, so this check should be optional, as opposed to mandatory.
Additionally, Popen causes problem when using mpi4py.

For some reason whenever I try to run the first example it throws the following error: PysimmError: A pair_style must be defined during initialization

from pysimm import system, lmps, forcefield

s = system.System()

m = s.molecules.add(system.Molecule())

f = forcefield.Dreiding()
#f = forcefield.Gaff2()

dreiding_C_ = s.particle_types.add(f.particle_types.get('C_3')[0].copy())
dreiding_H_ = s.particle_types.add(f.particle_types.get('H_')[0].copy())

#gaff_c3 = s.particle_types.add(f.particle_types.get('c3')[0].copy())
#gaff_hc = s.particle_types.add(f.particle_types.get('hc')[0].copy())

#c1 = s.particles.add(system.Particle(type=gaff_c3, x=0, y=0, z=0, charge=0, molecule=m))
c1 = s.particles.add(system.Particle(type=dreiding_C_, x=0, y=0, z=0, charge=0, molecule=m))

#h1 = s.add_particle_bonded_to(system.Particle(type=gaff_hc, charge=0, molecule=m), c1, f)
#h2 = s.add_particle_bonded_to(system.Particle(type=gaff_hc, charge=0, molecule=m), c1, f)
#h3 = s.add_particle_bonded_to(system.Particle(type=gaff_hc, charge=0, molecule=m), c1, f)
#h4 = s.add_particle_bonded_to(system.Particle(type=gaff_hc, charge=0, molecule=m), c1, f)

h1 = s.add_particle_bonded_to(system.Particle(type=dreiding_H_, charge=0, molecule=m), c1, f)
h2 = s.add_particle_bonded_to(system.Particle(type=dreiding_H_, charge=0, molecule=m), c1, f)
h3 = s.add_particle_bonded_to(system.Particle(type=dreiding_H_, charge=0, molecule=m), c1, f)
h4 = s.add_particle_bonded_to(system.Particle(type=dreiding_H_, charge=0, molecule=m), c1, f)

s.apply_charges(f, charges='gasteiger')

s.set_box(padding=10)

s.bond_style='lj'
s.pair_style='harmonic'
s.angle_style='harmonic'

lmps.quick_min(s, min_style='fire', name='fire_min')

#lmps.quick_min(s, min_style='sd', name='min_sd')
#lmps.quick_min(s, min_style='cg', name='min_cg')

s.write_pdb('methane.pdb')

Do you have any idea what could be wrong?
Thanks in advance!

How to use pysimm to cross-link custom molecules？

Some examples and creation of model monomers don't work without internet connection, because they make requests to the PubChem API. We can cache the response from PubChem and encapsulate in try/except in case the users cannot connect to the internet.

Docstrings (used for Sphinx documentation) should use correct syntax for referencing other pysimm classes/functions/methods to make documentation easier to travel through.
pysimm.system.System should be :class:`~pysimm.system.System`

Dear Developers,
I would like to build a polymer (Like PIM-1) with two heads and two tails. I can find examples with single head and tail and procedures to assign them. However, I am finding it difficult to define two heads in the code. Could you please provide some examples?

If you can please edit the following code to assign two heads that would be great.

def monomer():

    s = system.read_mol('My.mol'))
    m = s.molecules[1]
    f = forcefield.Gaff2()
    
    for b in s.bonds:
        if b.a.bonds.count == 3 and b.b.bonds.count == 3:
            b.order = 4
    
    s.apply_forcefield(f)
    
    c1 = s.particles[1] #head 1
    c2 = s.particles[2]  #head 2

    c5 = s.particles[5] #tail 1
    c5 = s.particles[6]  #tail 2

    
    for b in c1.bonds:
        if b.a.elem == 'H' or b.b.elem == 'H':
            pb = b.a if b.b is c1 else b.b
            s.particles.remove(pb.tag, update=False)
            break
        
    for b in c5.bonds:
        if b.a.elem == 'H' or b.b.elem == 'H':
            pb = b.a if b.b is c1 else b.b
            s.particles.remove(pb.tag, update=False)
            break
            
    s.remove_spare_bonding()

    c1.linker = 'head'
    #c2.linker= 'head'
    c5.linker = 'tail'
    #c6.linker= 'tail'
    
    lmps.quick_min(s, min_style='fire')
    
    s.add_particle_bonding()
    
    return s

The current pysimm version has the Dreiding FF with the harmonic improper style. The original paper suggests using the umbrella functional form for impropers. The results from both styles are similar but can lead to slight differences in structure. Is it possible to change the Dreiding force constants to be consistent with the umbrella style functional form from the original paper and then modify the lammps module to generate umbrella style inputs?

Thank you.

Hi,
I am new in Pysimm and tried to run some examples recently. When I run the example 2_Methanol script, the Dreiding force field example run successful, but the Gaff2 return below error. Could anyone give me any suggestion about this? This situation also occur in example 8. Thank you very much!

Traceback (most recent call last):
File "/home/kun/Desktop/Support/pysimm/Examples/02_methanol/gaff2/create.py", line 24, in
run()
File "/home/kun/Desktop/Support/pysimm/Examples/02_methanol/gaff2/create.py", line 15, in run
lmps.quick_min(s,min_style='fire',settings={'np':1})
File "/home/kun/pysimm/pysimm/lmps.py", line 1121, in quick_min
sim.run(np, nanohub)
File "/home/kun/pysimm/pysimm/lmps.py", line 957, in run
raise PysimmError('There was a problem calling LAMMPS with {}'.format(prefix))

Dear developers of Pysimm,

I've encountered the following output for the dihedral coefficients when running the 04_polyethylene example with the GAFF2, which also happened when I tried to run it on the PAMPSA polymer I created on my own. These files are for the input in lammps.


I'm wondering why are there multiple sets of data for the same dihedrals, and which ones should I use to successfully run the model with lammps. What dihedral style should I use in lammps? Thank you very much.

Hi,

I'm trying to read Charmm parameter file to construct pcbaa polymer. My code is,

`from pysimm import lmps
from pysimm import system
from pysimm import forcefield

from pysimm.apps.random_walk import random_walk, random_walk_tacticity
from pysimm.apps.random_walk import check_tacticity
from pysimm.apps.random_walk import random_walk, random_walk_tacticity, check_tacticity
import matplotlib.pyplot as mplp

s = system.System()
f = forcefield.Charmm()
m = s.molecules.add(system.Molecule())`

then I constructed molecule with,

`gaff_CG2DC3 = s.particle_types.add(f.particle_types.get('CG2DC3')[0].copy())
gaff_CG2DC1 = s.particle_types.add(f.particle_types.get('CG2DC1')[0].copy())
gaff_CG2O1 = s.particle_types.add(f.particle_types.get('CG2O1')[0].copy())
gaff_CG321 = s.particle_types.add(f.particle_types.get('CG321')[0].copy())
gaff_CG324 = s.particle_types.add(f.particle_types.get('CG324')[0].copy())
gaff_CG2O3 = s.particle_types.add(f.particle_types.get('CG2O3')[0].copy())
gaff_HGA5 = s.particle_types.add(f.particle_types.get('HGA5')[0].copy())
gaff_HGA4 = s.particle_types.add(f.particle_types.get('HGA4')[0].copy())
gaff_HGA2 = s.particle_types.add(f.particle_types.get('HGA2')[0].copy())
gaff_HGP2 = s.particle_types.add(f.particle_types.get('HGP2')[0].copy())
gaff_OG2D1 = s.particle_types.add(f.particle_types.get('OG2D1')[0].copy())
gaff_OG2D2 = s.particle_types.add(f.particle_types.get('OG2D2')[0].copy())
gaff_NG3P2 = s.particle_types.add(f.particle_types.get('NG3P2')[0].copy())

C1 = s.particles.add(system.Particle(type=gaff_CG2DC3, x=0, y=0, z=0, charge=-0.417000, molecule=m))
C2 = s.add_particle_bonded_to(system.Particle(type=gaff_CG2DC1, x=0, y=0, z=1.5419, charge=-0.097000, molecule=m), C1, f)
H3 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA5, x=1.0218 , y=0, z=-0.3878, charge=0.210000, molecule=m), C1, f)
H4 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA5, x=-0.4291 , y=-0.9533, z=-0.3291, charge=0.210000, molecule=m), C1, f)
C5 = s.add_particle_bonded_to(system.Particle(type=gaff_CG2O1, x=0.4192, y=1.3138 , z=2.2238, charge=0.479000, molecule=m), C2, f)
H6 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA4, x=-1.0564 , y=-0.1017, z=1.8266, charge=0.150000, molecule=m), C2, f)
O7 = s.add_particle_bonded_to(system.Particle(type=gaff_OG2D1, x=-0.3907, y=1.8257, z=3.0075, charge=-0.531000, molecule=m), C5, f)
N8 = s.add_particle_bonded_to(system.Particle(type=gaff_NG2S1, x= 1.6654, y=1.8812, z=2.0153, charge=-0.32900, molecule=m), C5, f )
C9 = s.add_particle_bonded_to(system.Particle(type=gaff_CG321, x=2.6544, y=1.4433, z=1.0462, charge=-0.083000, molecule=m),N8, f )
H10 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA2, x=2.6622, y=0.3620, z=0.9367, charge=0.090000, molecule=m), C9, f )
H11 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA2, x=3.6519, y= 1.7467, z=1.3793, charge=0.090000, molecule=m), C9, f)
H12 = s.add_particle_bonded_to(system.Particle(type=gaff_HGP1, x=2.1280, y=1.7515, z=2.8923, charge=0.274000, molecule=m), N8, f)
C13 = s.add_particle_bonded_to(system.Particle(type=gaff_CG321, x=2.3494, y=2.1061, z=-0.3100 charge=-0.153000, molecule=m), C9, f)
H14 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA2, x=2.3839, y=1.3676, z=-1.0836, charge=0.090000, molecule=m), C13, f)
H15 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA2, x=1.3747, y=2.5465, z=-0.2806, charge=0.090000, molecule=m), C13, f)
C16 = s.add_particle_bonded_to(system.Particle(type=gaff_CG324, x=3.3975, y=3.1976, z=-0.5953 charge=0.129000, molecule=m),C16, f)
H17 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA2, x=4.3790, y=2.7784, z=-0.5187, charge=0.090000, molecule=m), C16, f)
H18 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA2, x=3.2912, y=3.9890, z=0.1170, charge=0.090000, molecule=m), C16, f)
N19 = s.add_particle_bonded_to(system.Particle(type=gaff_NG3P2, x=3.1958, y=3.7283, z=-1.9512, charge=-0.3910000, molecule=m),C16,f)
H20 = s.add_particle_bonded_to(system.Particle(type=gaff_HGP2, x=3.9652, y=3.4612, z=-2.5315, charge=0.320000, molecule=m),N19,f )
H21 = s.padd_particle_bonded_to(system.Particle(type=gaff_HGP2, x=2.3476, y=3.3589, z=-2.3310, charge=0.320000, molecule=m),N19,f)
C22 = s.add_particle_bonded_to(system.Particle(type=gaff_CG324, x=3.1101, y=5.1948, z=-1.8957 charge=0.173000, molecule=m),N19,f)
H23 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA2, x=3.9215, y=5.5777, z=-1.3127, charge=0.090000, molecule=m),C22,f)
H24 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA2, x=2.1820, y=5.4811, z=-1.4467 charge=0.090000, molecule=m),C22,f)
C25 = s.add_particle_bonded_to(system.Particle(type=gaff_CG321, x=3.1881, y=5.7679, z=-3.3230 charge=-0.242000, molecule=m),C22,f)
H26 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA2, x=4.1187, y=5.4867, z=-3.7699 charge=0.090000, molecule=m),C25,f)
H27 = s.add_particle_bonded_to(system.Particle(type=gaff_HGA2, x=2.3801, y=5.3804, z=-3.9076 charge=0.090000, molecule=m),C25,f)
C28 = s.add_particle_bonded_to(system.Particle(type=gaff_CG2O3 , x=3.0898, y=7.3037, z=-3.2655, charge=0.598000, molecule=m),C25,f)
O29 = s.add_particle_bonded_to(system.Particle(type=gaff_OG2D2, x=2.7764, y=7.9519, z=-2.0299, charge=-0.760000, molecule=m),C28,f)
O30 = s.add_particle_bonded_to(system.Particle(type=gaff_OG2D2, x=3.2842, y=7.9902, z=-4.3021, charge=-0.760000, molecule=m),C28,f)

`
with that I tried,

`s.apply_forcefield(f, charges='none')
lnkr_atoms = {'head': 1, 'tail': 2}

for itm, nm in lnkr_atoms.items():
s.particles[nm].linker = itm

p = forcefield.Charmm()
p.parse_add_file('pcbaa.prm')`

It gave me the error,

KeyError Traceback (most recent call last)
/scratch/local/44140053/ipykernel_189726/397443237.py in
1 p = forcefield.Charmm()
----> 2 p.parse_add_file('pcbaa.prm')

/orange/colina/dkodituwakku/peggox/pysimm/pysimm/forcefield/charmm.py in parse_add_file(self, file)
590 tmp = getattr(self, 'nbfix_types')._dict
591 for t in tmp.keys():
--> 592 outp_obj['nbfix_types'].append(tmp[t].dict)
593
594 DATA_PATH = os.path.join(os.path.dirname(os.path.realpath(file)), os.pardir, 'data', 'forcefields', 'charmm')

KeyError: 'nbfix_types'

Could you let me know what I'm doing wrong?

Thank you

lmps.md(), lmps.minimize(), and lmps.relax() are being deprecated. Their replacement functions (lmps.quick_md/min) should be used in apps.polymatic in their place.