Generate elastic network

The genelastic tool generates elastic network bonds for coarse-grained simulation topologies (.itp files). Elastic networks are widely used in CG/HPS/Go-like models to enforce:

  • tertiary structure stability

  • domain rigidity

  • elastic network modeling (ENM)

  • local structural preservation across residues

genelastic reads:

  • a PDB structure (for atom coordinates),

  • an input ITP topology (for atom list & existing bonds),

  • an elastic residue cluster file (ASCII list of atom groups),

and generates new Bond entries between atom pairs whose distances fall within user-specified bounds.

This tool is implemented in genelastic.py.

Overview

The program builds an elastic network by:

  1. Loading atomic coordinates from a PDB file

  2. Loading a topology from an ITP file

  3. Reading an ASCII file specifying clusters of residue IDs

  4. For each cluster: - consider all atom pairs - skip pairs already bonded in topology - skip pairs already added earlier in this run - compute the distance from PDB coordinates - if the distance is within the cutoff window

    elastic_lower r elastic_upper → add an elastic bond

  5. Write updated topology to an output ITP file

Elastic bonds are harmonic:

\[V(r) = \frac{1}{2} k(r - r_0)^2\]

where k is the elastic-force constant and r_0 is the native contact distance.

Usage

Minimal example:

dps genelastic \
    -f structure.pdb \
    -p input.itp \
    -o output.itp \
    -er clusters.txt

Arguments

Required arguments

-f, --structure FILE

Input PDB structure. Coordinates are used to determine native distances.

-p, --topology FILE

Input ITP file whose topology is to be extended. Must end with .itp.

-o, --output FILE

Output ITP filename. Extension .itp is added automatically if missing.

-er, --elastic-residues FILE

ASCII file where each line lists a cluster of atoms (1-based indices). Example:

3  5  8 10
2  4  7

Elastic-network parameters

-ef, --elastic-force-constant FLOAT

Force constant (k) in kJ / mol / nm². Default: 5000.

-el, --elastic-lower FLOAT

Lower distance cutoff in nm. Only pairs with r > elastic_lower are considered. Default: 0.5.

-eu, --elastic-upper FLOAT

Upper distance cutoff in nm. Only pairs with r < elastic_upper are considered. Default: 0.9.

Input File Formats

### 1. PDB File (structure)

Used only for coordinates (x, y, z) and box vectors.

### 2. ITP File (topology)

Must contain:

  • [ atoms ]

  • [ bonds ]

The elastic network is added to the [ bonds ] section using:

Bond(atom1, atom2, bond_length, force_constant)

### 3. Elastic Cluster File (ASCII)

Each non-empty line:

idx1 idx2 idx3 ...

These indexes are 1-based atom indices, not residue numbers.

A cluster means every atom in the line can potentially form an elastic bond with every other atom on that line, subject to distance criteria.

Program Workflow

  1. Read PDB structure

atoms, box = read_pdb(args.structure)

If reading fails:

## An exception occurred when trying to open structure file...
  1. Read ITP topology

topology = read_itp(args.topology)

If input ITP file does not end in .itp → immediate error.

  1. Read elastic clusters

elastic_cluster_list = [
    [int(i)-1 for i in line.split()]
    for line in open(args.elastic_residues)
]

Each atom index is converted to 0-based.

  1. Prepare existing bond lookup

exist_bonds = [[bond.a1, bond.a2] for bond in topology.bonds]

Used to avoid overwriting or duplicating bonds.

  1. Loop over each cluster & atom pair

For every pair (i, j) in a cluster:

  • skip if i == j

  • skip if pair appears in exist_bonds

  • skip if pair appears in previously added_bonds

  • compute distance from PDB coordinates:

    distance = np.linalg.norm(coor_1 - coor_2)
    
  • check:

    elastic_lower ≤ distance ≤ elastic_upper
    
  • if passed → add an elastic bond with:

    • bond_length = distance (in nm)

    • k = elastic-force-constant (converted to OpenMM units)

Then:

topology.bonds.append(Bond(atomid_1, atomid_2, bond_length, bond_k))

Program message:

## Adding bond (k=...) between atom X and Y with distance Z.
  1. Save output

write_itp(output_file_name, topology)

Final message:

## A total number of N bonds have been added.
## Output topology file written to <output>.itp.

Example

Example elastic residue file:

5  6  7
10 20 30 35

Run:

dps genelastic \
    -f model.pdb \
    -p protein.itp \
    -o protein_elastic.itp \
    -er elastic_groups.txt \
    -el 0.4 -eu 1.0 -ef 6000

Example output:

## Open structure file model.pdb which contains 414 atoms.
## Open topology file protein.itp which contains 414 atoms.
## Reference structure contains 2 clusters.
## Will generate elastic bonds between residue within 0.40 ~ 1.00 nm.
## Elastic bond constant will be 6000.00 kJ / nm ^ 2
## Adding bond (k=...) between atom 5 and 6 with distance 0.87.
## Adding bond (k=...) between atom 5 and 7 with distance 0.71.
...
## A total number of 48 bonds have been added.
## Output topology file written to protein_elastic.itp.

Error Messages

“ERROR: Only itp file can be treated as input.” Input topology must be .itp format.

“An exception occurred when trying to open structure file …” PDB could not be loaded.

“An exception occurred when trying to open angle file …” Actually refers to elastic-residue file loading error.

Duplicate or conflicting bonds Automatically skipped silently (existing bonds take priority).

Summary

dps genelastic is a powerful generator of elastic networks for coarse-grained simulations. It provides:

  • cluster-based ENM construction

  • distance-based filtering

  • duplicate-bond protection

  • clean integration into existing ITP files

  • precise unit handling (nm, kJ/mol/nm²)

It is ideal for stabilizing tertiary structures or imposing distance constraints in CGPS / DROPPS simulations.