PEO-LiTFSI with AmberTools¶

After reading this page you will be able to parameterize molecules with the Amber workflow (antechamber, parmchk2, tleap), build polymer chains with AmberPolymerBuilder, and assemble a multi-component electrolyte system.

!!! warning "External dependencies" This guide requires AmberTools (via conda), RDKit, and Packmol. All three must be installed and accessible. Without AmberTools, no code on this page will run.

??? note "Setting up AmberTools" Install AmberTools in a dedicated conda environment:

```bash
conda create -n AmberTools25 -c conda-forge ambertools=25
conda activate AmberTools25
# Verify installation
which antechamber   # should print a path
which tleap         # should print a path
```

MolPy's wrapper classes activate the conda environment automatically when running commands, so you do not need to keep it active in your shell. The `env="AmberTools25"` parameter in the code below tells the wrapper which environment to activate.

If you use a different environment name, replace `"AmberTools25"` throughout this guide.

Six stages take a molecule from SMILES to a LAMMPS-ready electrolyte¶

The workflow begins with the hardest piece — parameterizing TFSI, the anion, using the standard Amber small-molecule pipeline of antechamber, parmchk2, and tleap. Li⁺ is simpler: its nonbond parameters are taken directly from Åqvist (1990) and written into an frcmod by hand. With both ions parameterized, PEO chains are built using the AmberPolymerBuilder, which wraps prepgen and tleap internally. The three component force fields are then merged, Packmol places the molecules at realistic density, and the final system is exported to LAMMPS.

Antechamber assigns GAFF types and BCC charges to TFSI¶

The Amber workflow for small molecules is: antechamber (assign types + charges) → parmchk2 (missing parameters) → tleap (topology + coordinates).

In [1]:

from pathlib import Path
import molpy as mp
from molpy.adapter import RDKitAdapter
from molpy.tool import Generate3D
from molpy.io.writers import write_pdb
from molpy.wrapper import AntechamberWrapper, Parmchk2Wrapper, TLeapWrapper

output_dir = Path("07_output")
ions_dir = output_dir / "ions"
ions_dir.mkdir(parents=True, exist_ok=True)

# Create TFSI from SMILES and generate 3D coordinates
tfsi = mp.parser.parse_molecule("O=S(=O)(C(F)(F)F)[N-]S(=O)(=O)C(F)(F)F")
adapter = RDKitAdapter(internal=tfsi)
adapter = Generate3D(
    add_hydrogens=False,
    embed=True,
    optimize=True,
    update_internal=True,
)(adapter)
tfsi = adapter.get_internal()

# Write PDB for antechamber input
write_pdb(ions_dir / "tfsi.pdb", tfsi.to_frame())

from pathlib import Path import molpy as mp from molpy.adapter import RDKitAdapter from molpy.tool import Generate3D from molpy.io.writers import write_pdb from molpy.wrapper import AntechamberWrapper, Parmchk2Wrapper, TLeapWrapper output_dir = Path("07_output") ions_dir = output_dir / "ions" ions_dir.mkdir(parents=True, exist_ok=True) # Create TFSI from SMILES and generate 3D coordinates tfsi = mp.parser.parse_molecule("O=S(=O)(C(F)(F)F)[N-]S(=O)(=O)C(F)(F)F") adapter = RDKitAdapter(internal=tfsi) adapter = Generate3D( add_hydrogens=False, embed=True, optimize=True, update_internal=True, )(adapter) tfsi = adapter.get_internal() # Write PDB for antechamber input write_pdb(ions_dir / "tfsi.pdb", tfsi.to_frame())

In [2]:

conda_env = "AmberTools25"

# Step 1: antechamber — assign GAFF types and BCC charges
ac = AntechamberWrapper(
    name="antechamber", workdir=ions_dir, env=conda_env, env_manager="conda"
)
ac.atomtype_assign(
    input_file=(ions_dir / "tfsi.pdb").absolute(),
    output_file=(ions_dir / "tfsi.mol2").absolute(),
    input_format="pdb",
    output_format="mol2",
    charge_method="bcc",
    atom_type="gaff2",
    net_charge=-1,
)

# Step 2: parmchk2 — generate missing parameters
parmchk2 = Parmchk2Wrapper(
    name="parmchk2", workdir=ions_dir, env=conda_env, env_manager="conda"
)
parmchk2.run(args=["-i", "tfsi.mol2", "-o", "tfsi.frcmod", "-f", "mol2", "-s", "gaff2"])

# Step 3: tleap — generate prmtop and inpcrd
leap_script = """source leaprc.gaff2
TFSI = loadmol2 tfsi.mol2
loadamberparams tfsi.frcmod
saveamberparm TFSI tfsi.prmtop tfsi.inpcrd
quit
"""
(ions_dir / "tfsi_leap.in").write_text(leap_script)

tleap = TLeapWrapper(name="tleap", workdir=ions_dir, env=conda_env, env_manager="conda")
tleap.run(args=["-f", "tfsi_leap.in"])

conda_env = "AmberTools25" # Step 1: antechamber — assign GAFF types and BCC charges ac = AntechamberWrapper( name="antechamber", workdir=ions_dir, env=conda_env, env_manager="conda" ) ac.atomtype_assign( input_file=(ions_dir / "tfsi.pdb").absolute(), output_file=(ions_dir / "tfsi.mol2").absolute(), input_format="pdb", output_format="mol2", charge_method="bcc", atom_type="gaff2", net_charge=-1, ) # Step 2: parmchk2 — generate missing parameters parmchk2 = Parmchk2Wrapper( name="parmchk2", workdir=ions_dir, env=conda_env, env_manager="conda" ) parmchk2.run(args=["-i", "tfsi.mol2", "-o", "tfsi.frcmod", "-f", "mol2", "-s", "gaff2"]) # Step 3: tleap — generate prmtop and inpcrd leap_script = """source leaprc.gaff2 TFSI = loadmol2 tfsi.mol2 loadamberparams tfsi.frcmod saveamberparm TFSI tfsi.prmtop tfsi.inpcrd quit """ (ions_dir / "tfsi_leap.in").write_text(leap_script) tleap = TLeapWrapper(name="tleap", workdir=ions_dir, env=conda_env, env_manager="conda") tleap.run(args=["-f", "tfsi_leap.in"])

---------------------------------------------------------------------------
FileNotFoundError                         Traceback (most recent call last)
Cell In[2], line 7
      3 # Step 1: antechamber — assign GAFF types and BCC charges
      4 ac = AntechamberWrapper(
      5     name="antechamber", workdir=ions_dir, env=conda_env, env_manager="conda"
      6 )
----> 7 ac.atomtype_assign(
      8     input_file=(ions_dir / "tfsi.pdb").absolute(),
      9     output_file=(ions_dir / "tfsi.mol2").absolute(),
     10     input_format="pdb",

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/wrapper/antechamber.py:98, in AntechamberWrapper.atomtype_assign(self, input_file, output_file, input_format, output_format, charge_method, atom_type, net_charge, formal_charges)
     95 if formal_charges:
     96     args.extend(["-cf", "y"])
---> 98 return self.run_raw(args=args)

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/wrapper/antechamber.py:45, in AntechamberWrapper.run_raw(self, args, input_text)
     30 def run_raw(
     31     self,
     32     args: list[str],
     33     *,
     34     input_text: str | None = None,
     35 ) -> subprocess.CompletedProcess[str]:
     36     """Execute antechamber with raw arguments.
     37 
     38     Args:
   (...)     43         The completed process result.
     44     """
---> 45     return self.run(args=args, input_text=input_text)

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/wrapper/base.py:241, in Wrapper.run(self, args, input_text, capture_output, check)
    238 if real_cwd is not None:
    239     real_cwd.mkdir(parents=True, exist_ok=True)
--> 241 return subprocess.run(
    242     final_args,
    243     cwd=str(real_cwd) if real_cwd is not None else None,
    244     input=input_text,
    245     capture_output=capture_output,
    246     text=True,
    247     env=self._merged_env(),
    248     check=check,
    249 )

File ~/.asdf/installs/python/3.13.3/lib/python3.13/subprocess.py:554, in run(input, capture_output, timeout, check, *popenargs, **kwargs)
    551     kwargs['stdout'] = PIPE
    552     kwargs['stderr'] = PIPE
--> 554 with Popen(*popenargs, **kwargs) as process:
    555     try:
    556         stdout, stderr = process.communicate(input, timeout=timeout)

File ~/.asdf/installs/python/3.13.3/lib/python3.13/subprocess.py:1039, in Popen.__init__(self, args, bufsize, executable, stdin, stdout, stderr, preexec_fn, close_fds, shell, cwd, env, universal_newlines, startupinfo, creationflags, restore_signals, start_new_session, pass_fds, user, group, extra_groups, encoding, errors, text, umask, pipesize, process_group)
   1035         if self.text_mode:
   1036             self.stderr = io.TextIOWrapper(self.stderr,
   1037                     encoding=encoding, errors=errors)
-> 1039     self._execute_child(args, executable, preexec_fn, close_fds,
   1040                         pass_fds, cwd, env,
   1041                         startupinfo, creationflags, shell,
   1042                         p2cread, p2cwrite,
   1043                         c2pread, c2pwrite,
   1044                         errread, errwrite,
   1045                         restore_signals,
   1046                         gid, gids, uid, umask,
   1047                         start_new_session, process_group)
   1048 except:
   1049     # Cleanup if the child failed starting.
   1050     for f in filter(None, (self.stdin, self.stdout, self.stderr)):

File ~/.asdf/installs/python/3.13.3/lib/python3.13/subprocess.py:1969, in Popen._execute_child(self, args, executable, preexec_fn, close_fds, pass_fds, cwd, env, startupinfo, creationflags, shell, p2cread, p2cwrite, c2pread, c2pwrite, errread, errwrite, restore_signals, gid, gids, uid, umask, start_new_session, process_group)
   1967     err_msg = os.strerror(errno_num)
   1968 if err_filename is not None:
-> 1969     raise child_exception_type(errno_num, err_msg, err_filename)
   1970 else:
   1971     raise child_exception_type(errno_num, err_msg)

FileNotFoundError: [Errno 2] No such file or directory: 'conda'

Li⁺ needs no charge calculation — literature parameters go directly into an frcmod file¶

Li⁺ has no bonded terms and no partial charges to compute, so antechamber is not needed. Instead, write the nonbond parameters from Åqvist (1990) directly into an frcmod file and create the prmtop with tleap.

Li⁺ nonbond parameters — Åqvist (1990), J. Phys. Chem. 94, 8021–8024, DOI: 10.1021/j100384a009. These were fitted to hydration free energies and are the standard choice for polymer electrolyte simulations with GAFF.

Parameter	Value
Rmin/2	1.137 Å
ε	0.0183 kcal/mol

In [3]:

from molpy.io import read_amber

li_dir = output_dir / "li"
li_dir.mkdir(parents=True, exist_ok=True)

# Write Åqvist (1990) frcmod — NONBON uses Rmin/2 and epsilon
li_frcmod = """Li+ Aqvist 1990 parameters
MASS
LI    6.941               0.0000000

BOND

ANGLE

DIHE

IMPROPER

NONBON
  LI        1.137        0.0183

"""
(li_dir / "li.frcmod").write_text(li_frcmod)

# Minimal mol2 for a single Li+ atom (net charge = +1)
li_mol2 = """@<TRIPOS>MOLECULE
LIT
 1 0 0 0 0
SMALL
USER_CHARGES

@<TRIPOS>ATOM
      1 LI          0.0000    0.0000    0.0000 LI    1  LIT      1.000000
@<TRIPOS>BOND
"""
(li_dir / "li.mol2").write_text(li_mol2)

# tleap: generate prmtop for Li+
li_leap = """source leaprc.gaff2
loadamberparams li.frcmod
LIT = loadmol2 li.mol2
saveamberparm LIT li.prmtop li.inpcrd
quit
"""
(li_dir / "li_leap.in").write_text(li_leap)

tleap_li = TLeapWrapper(
    name="tleap", workdir=li_dir, env=conda_env, env_manager="conda"
)
tleap_li.run(args=["-f", "li_leap.in"])

li_frame, li_ff = read_amber(li_dir / "li.prmtop", li_dir / "li.inpcrd")
print(
    f"Li+: {li_frame['atoms'].nrows} atom, charge={li_frame['atoms']['charge'][0]:.1f}"
)

from molpy.io import read_amber li_dir = output_dir / "li" li_dir.mkdir(parents=True, exist_ok=True) # Write Åqvist (1990) frcmod — NONBON uses Rmin/2 and epsilon li_frcmod = """Li+ Aqvist 1990 parameters MASS LI 6.941 0.0000000 BOND ANGLE DIHE IMPROPER NONBON LI 1.137 0.0183 """ (li_dir / "li.frcmod").write_text(li_frcmod) # Minimal mol2 for a single Li+ atom (net charge = +1) li_mol2 = """@<TRIPOS>MOLECULE LIT 1 0 0 0 0 SMALL USER_CHARGES @<TRIPOS>ATOM 1 LI 0.0000 0.0000 0.0000 LI 1 LIT 1.000000 @<TRIPOS>BOND """ (li_dir / "li.mol2").write_text(li_mol2) # tleap: generate prmtop for Li+ li_leap = """source leaprc.gaff2 loadamberparams li.frcmod LIT = loadmol2 li.mol2 saveamberparm LIT li.prmtop li.inpcrd quit """ (li_dir / "li_leap.in").write_text(li_leap) tleap_li = TLeapWrapper( name="tleap", workdir=li_dir, env=conda_env, env_manager="conda" ) tleap_li.run(args=["-f", "li_leap.in"]) li_frame, li_ff = read_amber(li_dir / "li.prmtop", li_dir / "li.inpcrd") print( f"Li+: {li_frame['atoms'].nrows} atom, charge={li_frame['atoms']['charge'][0]:.1f}" )

---------------------------------------------------------------------------
FileNotFoundError                         Traceback (most recent call last)
Cell In[3], line 50
     46 
     47 tleap_li = TLeapWrapper(
     48     name="tleap", workdir=li_dir, env=conda_env, env_manager="conda"
     49 )
---> 50 tleap_li.run(args=["-f", "li_leap.in"])
     51 
     52 li_frame, li_ff = read_amber(li_dir / "li.prmtop", li_dir / "li.inpcrd")
     53 print(

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/wrapper/base.py:241, in Wrapper.run(self, args, input_text, capture_output, check)
    238 if real_cwd is not None:
    239     real_cwd.mkdir(parents=True, exist_ok=True)
--> 241 return subprocess.run(
    242     final_args,
    243     cwd=str(real_cwd) if real_cwd is not None else None,
    244     input=input_text,
    245     capture_output=capture_output,
    246     text=True,
    247     env=self._merged_env(),
    248     check=check,
    249 )

File ~/.asdf/installs/python/3.13.3/lib/python3.13/subprocess.py:554, in run(input, capture_output, timeout, check, *popenargs, **kwargs)
    551     kwargs['stdout'] = PIPE
    552     kwargs['stderr'] = PIPE
--> 554 with Popen(*popenargs, **kwargs) as process:
    555     try:
    556         stdout, stderr = process.communicate(input, timeout=timeout)

File ~/.asdf/installs/python/3.13.3/lib/python3.13/subprocess.py:1039, in Popen.__init__(self, args, bufsize, executable, stdin, stdout, stderr, preexec_fn, close_fds, shell, cwd, env, universal_newlines, startupinfo, creationflags, restore_signals, start_new_session, pass_fds, user, group, extra_groups, encoding, errors, text, umask, pipesize, process_group)
   1035         if self.text_mode:
   1036             self.stderr = io.TextIOWrapper(self.stderr,
   1037                     encoding=encoding, errors=errors)
-> 1039     self._execute_child(args, executable, preexec_fn, close_fds,
   1040                         pass_fds, cwd, env,
   1041                         startupinfo, creationflags, shell,
   1042                         p2cread, p2cwrite,
   1043                         c2pread, c2pwrite,
   1044                         errread, errwrite,
   1045                         restore_signals,
   1046                         gid, gids, uid, umask,
   1047                         start_new_session, process_group)
   1048 except:
   1049     # Cleanup if the child failed starting.
   1050     for f in filter(None, (self.stdin, self.stdout, self.stderr)):

File ~/.asdf/installs/python/3.13.3/lib/python3.13/subprocess.py:1969, in Popen._execute_child(self, args, executable, preexec_fn, close_fds, pass_fds, cwd, env, startupinfo, creationflags, shell, p2cread, p2cwrite, c2pread, c2pwrite, errread, errwrite, restore_signals, gid, gids, uid, umask, start_new_session, process_group)
   1967     err_msg = os.strerror(errno_num)
   1968 if err_filename is not None:
-> 1969     raise child_exception_type(errno_num, err_msg, err_filename)
   1970 else:
   1971     raise child_exception_type(errno_num, err_msg)

FileNotFoundError: [Errno 2] No such file or directory: 'conda'

Three monomer variants define the chain start, interior, and end¶

Port markers in BigSMILES ([>] and [<]) define connection points. A monomer with both ports is an interior repeat unit; one with a single port is an end cap. The builder uses the port annotations to decide which prepgen variant (HEAD / CHAIN / TAIL) to generate.

In [4]:

def parse_monomer_3d(bigsmiles):
    mol = mp.parser.parse_monomer(bigsmiles)
    adapter = RDKitAdapter(internal=mol)
    adapter = Generate3D(
        add_hydrogens=True,
        embed=True,
        optimize=True,
        update_internal=True,
    )(adapter)
    return adapter.get_internal()


# Head cap:  only < port → start of chain
me_head = parse_monomer_3d("{[][<]C[]}")

# Chain unit: both < and > → interior repeat
eo_chain = parse_monomer_3d("{[][<]COC[>][]}")

# Tail cap:  only > port → end of chain
me_tail = parse_monomer_3d("{[]C[>][]}")

library = {"MeH": me_head, "EO": eo_chain, "MeT": me_tail}

def parse_monomer_3d(bigsmiles): mol = mp.parser.parse_monomer(bigsmiles) adapter = RDKitAdapter(internal=mol) adapter = Generate3D( add_hydrogens=True, embed=True, optimize=True, update_internal=True, )(adapter) return adapter.get_internal() # Head cap: only < port → start of chain me_head = parse_monomer_3d("{[][<]C[]}") # Chain unit: both < and > → interior repeat eo_chain = parse_monomer_3d("{[][<]COC[>][]}") # Tail cap: only > port → end of chain me_tail = parse_monomer_3d("{[]C[>][]}") library = {"MeH": me_head, "EO": eo_chain, "MeT": me_tail}

AmberPolymerBuilder runs the full Amber pipeline internally¶

AmberPolymerBuilder wraps the monomer library, connector rules, and Amber tool chain (prepgen + tleap) into one builder that produces fully parameterized chains. Each unique chain length writes its Amber intermediate files into its own subdirectory under work_dir to prevent file conflicts.

In [5]:

from molpy.builder.polymer.ambertools import AmberPolymerBuilder

polymer_dir = output_dir / "polymer"
polymer_dir.mkdir(exist_ok=True)

builder = AmberPolymerBuilder(
    library=library,
    force_field="gaff2",
    charge_method="bcc",
    env="AmberTools25",
    env_manager="conda",
    work_dir=polymer_dir,
)

result = builder.build("{[#MeH][#EO]|10[#MeT]}")

from molpy.builder.polymer.ambertools import AmberPolymerBuilder polymer_dir = output_dir / "polymer" polymer_dir.mkdir(exist_ok=True) builder = AmberPolymerBuilder( library=library, force_field="gaff2", charge_method="bcc", env="AmberTools25", env_manager="conda", work_dir=polymer_dir, ) result = builder.build("{[#MeH][#EO]|10[#MeT]}")

---------------------------------------------------------------------------
FileNotFoundError                         Traceback (most recent call last)
Cell In[5], line 15
     11     env_manager="conda",
     12     work_dir=polymer_dir,
     13 )
     14 
---> 15 result = builder.build("{[#MeH][#EO]|10[#MeT]}")

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/builder/polymer/ambertools/amber_builder.py:106, in AmberPolymerBuilder.build(self, cgsmiles)
    103 self._validate_ir(ir)
    105 # Prepare all monomers (antechamber → parmchk2 → prepgen)
--> 106 self._prepare_monomers(ir.base_graph)
    108 # Generate and run tleap
    109 result = self._build_with_tleap(ir.base_graph, output_prefix="polymer")

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/builder/polymer/ambertools/amber_builder.py:215, in AmberPolymerBuilder._prepare_monomers(self, graph)
    207 mol2_file = monomer_dir / f"{label}.mol2"
    209 antechamber = AntechamberWrapper(
    210     name="antechamber",
    211     workdir=monomer_dir,
    212     env=self.env,
    213     env_manager=self.env_manager,
    214 )
--> 215 r = antechamber.atomtype_assign(
    216     input_file=input_pdb,
    217     output_file=mol2_file,
    218     input_format="pdb",
    219     output_format="mol2",
    220     charge_method=self.charge_method,
    221     atom_type=self.force_field,
    222 )
    223 if r.returncode != 0:
    224     raise RuntimeError(
    225         f"antechamber (mol2) failed for monomer '{label}':\n"
    226         f"{r.stderr or r.stdout}"
    227     )

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/wrapper/antechamber.py:98, in AntechamberWrapper.atomtype_assign(self, input_file, output_file, input_format, output_format, charge_method, atom_type, net_charge, formal_charges)
     95 if formal_charges:
     96     args.extend(["-cf", "y"])
---> 98 return self.run_raw(args=args)

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/wrapper/antechamber.py:45, in AntechamberWrapper.run_raw(self, args, input_text)
     30 def run_raw(
     31     self,
     32     args: list[str],
     33     *,
     34     input_text: str | None = None,
     35 ) -> subprocess.CompletedProcess[str]:
     36     """Execute antechamber with raw arguments.
     37 
     38     Args:
   (...)     43         The completed process result.
     44     """
---> 45     return self.run(args=args, input_text=input_text)

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/wrapper/base.py:241, in Wrapper.run(self, args, input_text, capture_output, check)
    238 if real_cwd is not None:
    239     real_cwd.mkdir(parents=True, exist_ok=True)
--> 241 return subprocess.run(
    242     final_args,
    243     cwd=str(real_cwd) if real_cwd is not None else None,
    244     input=input_text,
    245     capture_output=capture_output,
    246     text=True,
    247     env=self._merged_env(),
    248     check=check,
    249 )

File ~/.asdf/installs/python/3.13.3/lib/python3.13/subprocess.py:554, in run(input, capture_output, timeout, check, *popenargs, **kwargs)
    551     kwargs['stdout'] = PIPE
    552     kwargs['stderr'] = PIPE
--> 554 with Popen(*popenargs, **kwargs) as process:
    555     try:
    556         stdout, stderr = process.communicate(input, timeout=timeout)

File ~/.asdf/installs/python/3.13.3/lib/python3.13/subprocess.py:1039, in Popen.__init__(self, args, bufsize, executable, stdin, stdout, stderr, preexec_fn, close_fds, shell, cwd, env, universal_newlines, startupinfo, creationflags, restore_signals, start_new_session, pass_fds, user, group, extra_groups, encoding, errors, text, umask, pipesize, process_group)
   1035         if self.text_mode:
   1036             self.stderr = io.TextIOWrapper(self.stderr,
   1037                     encoding=encoding, errors=errors)
-> 1039     self._execute_child(args, executable, preexec_fn, close_fds,
   1040                         pass_fds, cwd, env,
   1041                         startupinfo, creationflags, shell,
   1042                         p2cread, p2cwrite,
   1043                         c2pread, c2pwrite,
   1044                         errread, errwrite,
   1045                         restore_signals,
   1046                         gid, gids, uid, umask,
   1047                         start_new_session, process_group)
   1048 except:
   1049     # Cleanup if the child failed starting.
   1050     for f in filter(None, (self.stdin, self.stdout, self.stderr)):

File ~/.asdf/installs/python/3.13.3/lib/python3.13/subprocess.py:1969, in Popen._execute_child(self, args, executable, preexec_fn, close_fds, pass_fds, cwd, env, startupinfo, creationflags, shell, p2cread, p2cwrite, c2pread, c2pwrite, errread, errwrite, restore_signals, gid, gids, uid, umask, start_new_session, process_group)
   1967     err_msg = os.strerror(errno_num)
   1968 if err_filename is not None:
-> 1969     raise child_exception_type(errno_num, err_msg, err_filename)
   1970 else:
   1971     raise child_exception_type(errno_num, err_msg)

FileNotFoundError: [Errno 2] No such file or directory: 'conda'

AmberPolymerBuilder.build() internally runs antechamber, parmchk2, prepgen, and tleap. The result carries the polymer Frame, ForceField, and paths to the intermediate Amber files.

In [6]:

peo_frame = result.frame
peo_ff = result.forcefield
print(f"PEO 10-mer: {peo_frame['atoms'].nrows} atoms")

peo_frame = result.frame peo_ff = result.forcefield print(f"PEO 10-mer: {peo_frame['atoms'].nrows} atoms")

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
Cell In[6], line 1
----> 1 peo_frame = result.frame
      2 peo_ff = result.forcefield
      3 print(f"PEO 10-mer: {peo_frame['atoms'].nrows} atoms")

NameError: name 'result' is not defined

Merging three force fields before packing prevents type conflicts¶

Merging is done before packing rather than after because Packmol operates on coordinates only — it has no awareness of force field types. If two components share an atom type name with different parameters, a post-packing merge would silently overwrite one of them. Merging first makes any type name collision an error before coordinates are generated.

In [7]:

import numpy as np
from molpy.io import read_amber
from molpy.pack import Molpack, InsideBoxConstraint

# Read TFSI from Amber files generated in Stage 1
tfsi_frame, tfsi_ff = read_amber(
    ions_dir / "tfsi.prmtop",
    ions_dir / "tfsi.inpcrd",
)

# Merge all three force fields: PEO + TFSI + Li+
combined_ff = peo_ff.merge(tfsi_ff).merge(li_ff)

# Pack system
box_size = 60.0
packer = Molpack(workdir=output_dir / "packmol")
constraint = InsideBoxConstraint(length=[box_size] * 3, origin=[0.0] * 3)
packer.add_target(peo_frame, number=3, constraint=constraint)
packer.add_target(li_frame, number=10, constraint=constraint)
packer.add_target(tfsi_frame, number=10, constraint=constraint)

system = packer.optimize(max_steps=20000, seed=12345)
system.box = mp.Box.cubic(box_size)

import numpy as np from molpy.io import read_amber from molpy.pack import Molpack, InsideBoxConstraint # Read TFSI from Amber files generated in Stage 1 tfsi_frame, tfsi_ff = read_amber( ions_dir / "tfsi.prmtop", ions_dir / "tfsi.inpcrd", ) # Merge all three force fields: PEO + TFSI + Li+ combined_ff = peo_ff.merge(tfsi_ff).merge(li_ff) # Pack system box_size = 60.0 packer = Molpack(workdir=output_dir / "packmol") constraint = InsideBoxConstraint(length=[box_size] * 3, origin=[0.0] * 3) packer.add_target(peo_frame, number=3, constraint=constraint) packer.add_target(li_frame, number=10, constraint=constraint) packer.add_target(tfsi_frame, number=10, constraint=constraint) system = packer.optimize(max_steps=20000, seed=12345) system.box = mp.Box.cubic(box_size)

---------------------------------------------------------------------------
FileNotFoundError                         Traceback (most recent call last)
Cell In[7], line 6
      2 from molpy.io import read_amber
      3 from molpy.pack import Molpack, InsideBoxConstraint
      4 
      5 # Read TFSI from Amber files generated in Stage 1
----> 6 tfsi_frame, tfsi_ff = read_amber(
      7     ions_dir / "tfsi.prmtop",
      8     ions_dir / "tfsi.inpcrd",
      9 )

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/io/readers.py:288, in read_amber_prmtop(prmtop, inpcrd, frame)
    286 prmtop_path = Path(prmtop)
    287 reader = AmberPrmtopReader(prmtop_path)
--> 288 frame, ff = reader.read(frame)
    290 if inpcrd is not None:
    291     from .data.amber import AmberInpcrdReader

File ~/.asdf/installs/python/3.13.3/lib/python3.13/site-packages/molpy/io/forcefield/amber.py:34, in AmberPrmtopReader.read(self, frame)
     33 def read(self, frame: Frame):
---> 34     with open(self.file) as f:
     35         lines = filter(
     36             lambda line: line, map(AmberPrmtopReader.sanitizer, f.readlines())
     37         )
     39     # read file and split into sections

FileNotFoundError: [Errno 2] No such file or directory: '07_output/ions/tfsi.prmtop'

Exporting skips pair_style because long-range electrostatics need it in the script¶

In [8]:

from molpy.io.writers import write_lammps_data, write_lammps_forcefield

lammps_dir = output_dir / "lammps"
lammps_dir.mkdir(exist_ok=True)
write_lammps_data(lammps_dir / "system.data", system, atom_style="full")
write_lammps_forcefield(lammps_dir / "system.ff", combined_ff, skip_pair_style=True)

from molpy.io.writers import write_lammps_data, write_lammps_forcefield lammps_dir = output_dir / "lammps" lammps_dir.mkdir(exist_ok=True) write_lammps_data(lammps_dir / "system.data", system, atom_style="full") write_lammps_forcefield(lammps_dir / "system.ff", combined_ff, skip_pair_style=True)

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
Cell In[8], line 5
      1 from molpy.io.writers import write_lammps_data, write_lammps_forcefield
      2 
      3 lammps_dir = output_dir / "lammps"
      4 lammps_dir.mkdir(exist_ok=True)
----> 5 write_lammps_data(lammps_dir / "system.data", system, atom_style="full")
      6 write_lammps_forcefield(lammps_dir / "system.ff", combined_ff, skip_pair_style=True)

NameError: name 'system' is not defined

skip_pair_style=True omits the pair_style line from the force-field file. This is required when using kspace (long-range electrostatics), because the pair_style must be set by the simulation input script rather than the force-field file.

Troubleshooting¶

Symptom	Check
Antechamber fails	Verify PDB has correct atom names and no duplicate IDs
TFSI charge wrong	Use charge_method="bcc" and verify -nc -1
tleap fails for Li⁺	Confirm the mol2 atom type (LI) matches the frcmod NONBON entry
Polymer build fails	Check port markers in monomer SMILES
Force field merge conflict	Inspect atom type names for collisions between PEO and TFSI
Packing fails	Increase box size or reduce molecule count

See also: Force Field Typification, Wrapper and Adapter.