acpype.topol.ACTopol
- class acpype.topol.ACTopol(inputFile, binaries={'ac_bin': 'antechamber', 'obabel_bin': 'obabel'}, chargeType='bcc', chargeVal=None, multiplicity='1', atomType='gaff2', force=False, basename=None, debug=False, outTopol='all', allhdg=False, timeTol=10800, qprog='sqm', ekFlag=None, verbose=True, gmx4=False, merge=False, direct=False, is_sorted=False, chiral=False, amb2gmx=False, level=20)
Bases:
acpype.topol.AbstractTopolClass to build the AC topologies (Antechamber AmberTools)
- __init__(inputFile, binaries={'ac_bin': 'antechamber', 'obabel_bin': 'obabel'}, chargeType='bcc', chargeVal=None, multiplicity='1', atomType='gaff2', force=False, basename=None, debug=False, outTopol='all', allhdg=False, timeTol=10800, qprog='sqm', ekFlag=None, verbose=True, gmx4=False, merge=False, direct=False, is_sorted=False, chiral=False, amb2gmx=False, level=20)
Methods
__init__(inputFile[, binaries, chargeType, ...])balanceCharges(chargeList[, FirstNonSoluteId])Note that python is very annoying about floating points.
Check FRCMOD file
checkLeapLog(log)Check Leap log
checkSmiles()Check XYZ and TOP files
Convert PDB to MOL2 by using obabel
Convert Smiles to MOL2 by using obabel
If successful, Amber Top and Xyz files will be generated
Create molTop obj
Delete temporary output files
execAntechamber([chargeType, atomType])To call Antechamber and execute it
Execute obabel
Execute parmchk
Execute tleap
Get non-bonded coefficients
Get Angles
getAtoms()Set a list with all atoms objects build from dat in acFileTop Set also if molTopol atom type system is gaff or amber Set also list atomTypes Set also resid Set also molTopol total charge
getBonds()Get Bonds
Get chiral atoms, its 4 neighbours and improper dihedral angle to store non-planar improper dihedrals for CNS (and CNS only!)
For a given acFileXyz file, return a list of coords as: [[x1,y1,z1],[x2,y2,z2], etc.]
Get dihedrals (proper and imp), condensed list of prop dih and atomPairs
getFlagData(flag)For a given acFileTop flag, return a list of the data related
Get a 3 capital letters code from acFileTop Returns a list.
Guess the charge of a system based on antechamber Returns None in case of error
locateDat(aFile)locate a file pertinent to $AMBERHOME/dat/leap/parm/
makeDir()Make Dir
Example
printDebug([text])Debug log level
printDebugQuoted([text])Print quoted messages
printError([text])Error log level
printErrorQuoted([text])Print quoted messages
printMess([text])Info log level
printWarn([text])Warn log level
readMol2TotalCharge(mol2File)Reads the charges in given mol2 file and returns the total
search([name, alist])returns a list with all atomName matching 'name' or just the first case
Atom types names in Gromacs TOP file are not case sensitive; this routine will append a '_' to lower case atom type.
It takes self.condensedProperDihedrals and returns self.properDihedralsCoefRB, a reduced list of quartet atoms + RB.
Set a 3 letter residue name and check coords for issues like duplication, atoms too close or too sparse
signal_handler(_signum, _frame)Signal handler
Re-sort atoms for gromacs, which expects all hydrogens to immediately follow the heavy atom they are bonded to and belong to the same charge group.
Write CHARMM topology files
Write CNS topology files
Write GRO files
Write GMX topology file
Write GMX topology Files
Write MDP for test with GROMACS
writePdb(afile)Write a new PDB file with the atom names defined by Antechamber
writePosreFile([fc])Write file with positional restraints for heavy atoms http://www.mdtutorials.com/gmx/complex/06_equil.html
- balanceCharges(chargeList, FirstNonSoluteId=None)
Note that python is very annoying about floating points. Even after balance, there will always be some residue of order e-12 to e-16, which is believed to vanished once one writes a topology file, say, for CNS or GMX, where floats are represented with 4 or 5 maximum decimals.
- checkFrcmod()
Check FRCMOD file
- checkLeapLog(log)
Check Leap log
- checkXyzAndTopFiles()
Check XYZ and TOP files
- convertPdbToMol2()
Convert PDB to MOL2 by using obabel
- convertSmilesToMol2()
Convert Smiles to MOL2 by using obabel
- createACTopol()
If successful, Amber Top and Xyz files will be generated
- createMolTopol()
Create molTop obj
- delOutputFiles()
Delete temporary output files
- execAntechamber(chargeType=None, atomType=None) bool
To call Antechamber and execute it
- Parameters
chargeType ([str], optional) – bcc, gas or user. Defaults to None/bcc.
atomType ([str], optional) – gaff, amber, gaff2, amber2. Defaults to None/gaff2.
- Returns
True if failed.
- Return type
bool
Welcome to antechamber 21.0: molecular input file processor. Usage: antechamber -i input file name -fi input file format -o output file name -fo output file format -c charge method -cf charge file name -nc net molecular charge (int) -a additional file name -fa additional file format -ao additional file operation crd : only read in coordinate crg : only read in charge radius: only read in radius name : only read in atom name type : only read in atom type bond : only read in bond type -m multiplicity (2S+1), default is 1 -rn residue name, overrides input file, default is MOL -rf residue topology file name in prep input file, default is molecule.res -ch check file name for gaussian, default is 'molecule' -ek mopac or sqm keyword, inside quotes; overwrites previous ones -gk gaussian job keyword, inside quotes, is ignored when both -gopt and -gsp are used -gopt gaussian job keyword for optimization, inside quotes -gsp gaussian job keyword for single point calculation, inside quotes -gm gaussian memory keyword, inside quotes, such as "%mem=1000MB" -gn gaussian number of processors keyword, inside quotes, such as "%nproc=8" -gdsk gaussian maximum disk usage keyword, inside quotes, such as "%maxdisk=50GB" -gv add keyword to generate gesp file (for Gaussian 09 only) 1 : yes 0 : no, the default -ge gaussian esp file generated by iop(6/50=1), default is g09.gesp -tor torsional angle list, inside a pair of quotes, such as "1-2-3-4:0,5-6-7-8" ':1' or ':0' indicates the torsional angle is frozen or not -df am1-bcc precharge flag, 2 - use sqm(default); 0 - use mopac -at atom type gaff : the default gaff2: for gaff2 (beta-version) amber: for PARM94/99/99SB bcc : bcc sybyl: sybyl -du fix duplicate atom names: yes(y)[default] or no(n) -bk component/block Id, for ccif -an adjust atom names: yes(y) or no(n) the default is 'y' for 'mol2' and 'ac' and 'n' for the other formats -j atom type and bond type prediction index, default is 4 0 : no assignment 1 : atom type 2 : full bond types 3 : part bond types 4 : atom and full bond type 5 : atom and part bond type -s status information: 0(brief), 1(default) or 2(verbose) -eq equalizing atomic charge, default is 1 for '-c resp' and '-c bcc' and 0 for the other charge methods 0 : no use 1 : by atomic paths 2 : by atomic paths and structural information, i.e. E/Z configurations -pf remove intermediate files: yes(y) or no(n)[default] -pl maximum path length to determin equivalence of atomic charges for resp and bcc, the smaller the value, the faster the algorithm, default is -1 (use full length), set this parameter to 10 to 30 if your molecule is big (# atoms >= 100) -seq atomic sequence order changable: yes(y)[default] or no(n) -dr acdoctor mode: yes(y)[default] or no(n) -i -o -fi and -fo must appear; others are optional Use 'antechamber -L' to list the supported file formats and charge methods List of the File Formats file format type abbre. index | file format type abbre. index -------------------------------------------------------------- Antechamber ac 1 | Sybyl Mol2 mol2 2 PDB pdb 3 | Modified PDB mpdb 4 AMBER PREP (int) prepi 5 | AMBER PREP (car) prepc 6 Gaussian Z-Matrix gzmat 7 | Gaussian Cartesian gcrt 8 Mopac Internal mopint 9 | Mopac Cartesian mopcrt 10 Gaussian Output gout 11 | Mopac Output mopout 12 Alchemy alc 13 | CSD csd 14 MDL mdl 15 | Hyper hin 16 AMBER Restart rst 17 | Jaguar Cartesian jcrt 18 Jaguar Z-Matrix jzmat 19 | Jaguar Output jout 20 Divcon Input divcrt 21 | Divcon Output divout 22 SQM Input sqmcrt 23 | SQM Output sqmout 24 Charmm charmm 25 | Gaussian ESP gesp 26 Component cif ccif 27 | GAMESS dat gamess 28 Orca input orcinp 29 | Orca output orcout 30 -------------------------------------------------------------- AMBER restart file can only be read in as additional file. List of the Charge Methods charge method abbre. index | charge method abbre. index -------------------------------------------------------------- RESP resp 1 | AM1-BCC bcc 2 CM1 cm1 3 | CM2 cm2 4 ESP (Kollman) esp 5 | Mulliken mul 6 Gasteiger gas 7 | Read in charge rc 8 Write out charge wc 9 | Delete Charge dc 10 --------------------------------------------------------------
- execObabel()
Execute obabel
- execParmchk()
Execute parmchk
- execTleap()
Execute tleap
- getABCOEFs()
Get non-bonded coefficients
- getAngles()
Get Angles
- getAtoms()
Set a list with all atoms objects build from dat in acFileTop Set also if molTopol atom type system is gaff or amber Set also list atomTypes Set also resid Set also molTopol total charge
- getBonds()
Get Bonds
- getChirals()
Get chiral atoms, its 4 neighbours and improper dihedral angle to store non-planar improper dihedrals for CNS (and CNS only!)
- getCoords()
For a given acFileXyz file, return a list of coords as: [[x1,y1,z1],[x2,y2,z2], etc.]
- getDihedrals()
Get dihedrals (proper and imp), condensed list of prop dih and atomPairs
- getFlagData(flag)
For a given acFileTop flag, return a list of the data related
- getResidueLabel()
Get a 3 capital letters code from acFileTop Returns a list.
- guessCharge()
Guess the charge of a system based on antechamber Returns None in case of error
- locateDat(aFile)
locate a file pertinent to $AMBERHOME/dat/leap/parm/
- makeDir()
Make Dir
- pickleSave()
Example
to restore:
from acpype import * # import cPickle as pickle import pickle mol = pickle.load(open('DDD.pkl','rb'))
- printDebug(text='')
Debug log level
- printDebugQuoted(text='')
Print quoted messages
- printError(text='')
Error log level
- printErrorQuoted(text='')
Print quoted messages
- printMess(text='')
Info log level
- printWarn(text='')
Warn log level
- readMol2TotalCharge(mol2File)
Reads the charges in given mol2 file and returns the total
- search(name=None, alist=False)
returns a list with all atomName matching ‘name’ or just the first case
- setAtomType4Gromacs()
Atom types names in Gromacs TOP file are not case sensitive; this routine will append a ‘_’ to lower case atom type.
Example
>>> CA and ca -> CA and ca_
- setProperDihedralsCoef()
It takes self.condensedProperDihedrals and returns self.properDihedralsCoefRB, a reduced list of quartet atoms + RB. Coefficients ready for GMX (multiplied by 4.184)
self.properDihedralsCoefRB = [ [atom1,…, atom4], C[0:5] ]
For proper dihedrals: a quartet of atoms may appear with more than one set of parameters and to convert to GMX they are treated as RBs.
The resulting coefficients calculated here may look slighted different from the ones calculated by amb2gmx.pl because python is taken full float number from prmtop and not rounded numbers from rdparm.out as amb2gmx.pl does.
- setResNameCheckCoords()
Set a 3 letter residue name and check coords for issues like duplication, atoms too close or too sparse
- signal_handler(_signum, _frame)
Signal handler
- sortAtomsForGromacs()
Re-sort atoms for gromacs, which expects all hydrogens to immediately follow the heavy atom they are bonded to and belong to the same charge group.
Currently, atom mass < 1.2 is taken to denote a proton. This behaviour may be changed by modifying the ‘is_hydrogen’ function within.
JDC 2011-02-03
- writeCharmmTopolFiles()
Write CHARMM topology files
- writeCnsTopolFiles()
Write CNS topology files
- writeGroFile()
Write GRO files
- writeGromacsTop()
Write GMX topology file
- writeGromacsTopolFiles()
Write GMX topology Files
# from ~/Programmes/amber10/dat/leap/parm/gaff.dat #atom type atomic mass atomic polarizability comments ca 12.01 0.360 Sp2 C in pure aromatic systems ha 1.008 0.135 H bonded to aromatic carbon #bonded atoms harmonic force kcal/mol/A^2 eq. dist. Ang. comments ca-ha 344.3* 1.087** SOURCE3 1496 0.0024 0.0045 * for gmx: 344.3 * 4.184 * 100 * 2 = 288110 kJ/mol/nm^2 (why factor 2?) ** convert Ang to nm ( div by 10) for gmx: 1.087 A = 0.1087 nm # CA HA 1 0.10800 307105.6 ; ged from 340. bsd on C6H6 nmodes; PHE,TRP,TYR (from ffamber99bon.itp) # CA-HA 367.0 1.080 changed from 340. bsd on C6H6 nmodes; PHE,TRP,TYR (from parm99.dat) # angle HF kcal/mol/rad^2 eq angle degrees comments ca-ca-ha 48.5* 120.01 SOURCE3 2980 0.1509 0.2511 * to convert to gmx: 48.5 * 4.184 * 2 = 405.848 kJ/mol/rad^2 (why factor 2?) # CA CA HA 1 120.000 418.400 ; new99 (from ffamber99bon.itp) # CA-CA-HA 50.0 120.00 (from parm99.dat) # dihedral idivf barrier hight/2 kcal/mol phase degrees periodicity comments X -ca-ca-X 4 14.500* 180.000 2.000 intrpol.bsd.on C6H6 *convert 2 gmx: 14.5/4 * 4.184 * 2 (?) (yes in amb2gmx, not in topolbuild, why?) = 30.334 or 15.167 kJ/mol # X -CA-CA-X 4 14.50 180.0 2. intrpol.bsd.on C6H6 (from parm99.dat) # X CA CA X 3 30.334 0.000 -30.33400 0.000 0.000 0.000 ; intrpol.bsd.on C6H6 ;propers treated as RBs in GMX to use combine multiple AMBER torsions per quartet (from ffamber99bon.itp) # impr. dihedral barrier hight/2 phase degrees periodicity comments X -X -ca-ha 1.1* 180. 2. bsd.on C6H6 nmodes * to convert to gmx: 1.1 * 4.184 = 4.6024 kJ/mol/rad^2 # X -X -CA-HA 1.1 180. 2. bsd.on C6H6 nmodes (from parm99.dat) # X X CA HA 1 180.00 4.60240 2 ; bsd.on C6H6 nmodes ;impropers treated as propers in GROMACS to use correct AMBER analytical function (from ffamber99bon.itp) # 6-12 parms sigma = 2 * r * 2^(-1/6) epsilon # atomtype radius Ang. pot. well depth kcal/mol comments ha 1.4590* 0.0150** Spellmeyer ca 1.9080 0.0860 OPLS * to convert to gmx: sigma = 1.4590 * 2^(-1/6) * 2 = 2 * 1.29982 Ang. = 2 * 0.129982 nm = 1.4590 * 2^(5/6)/10 = 0.259964 nm ** to convert to gmx: 0.0150 * 4.184 = 0.06276 kJ/mol # amber99_3 CA 0.0000 0.0000 A 3.39967e-01 3.59824e-01 (from ffamber99nb.itp) # amber99_22 HA 0.0000 0.0000 A 2.59964e-01 6.27600e-02 (from ffamber99nb.itp) # C* 1.9080 0.0860 Spellmeyer # HA 1.4590 0.0150 Spellmeyer (from parm99.dat) # to convert r and epsilon to ACOEF and BCOEF # ACOEF = sqrt(e1*e2) * (r1 + r2)^12 ; BCOEF = 2 * sqrt(e1*e2) * (r1 + r2)^6 = 2 * ACOEF/(r1+r2)^6 # to convert ACOEF and BCOEF to r and epsilon # r = 0.5 * (2*ACOEF/BCOEF)^(1/6); ep = BCOEF^2/(4*ACOEF) # to convert ACOEF and BCOEF to sigma and epsilon (GMX) # sigma = (ACOEF/BCOEF)^(1/6) * 0.1 ; epsilon = 4.184 * BCOEF^2/(4*ACOEF) # ca ca 819971.66 531.10 # ca ha 76245.15 104.66 # ha ha 5716.30 18.52For proper dihedrals: a quartet of atoms may appear with more than one set of parameters and to convert to GMX they are treated as RBs; use the algorithm:
for(my $j=$i;$j<=$lines;$j++){ my $period = $pn{$j}; if($pk{$j}>0) { $V[$period] = 2*$pk{$j}*$cal; } # assign V values to C values as predefined # if($period==1){ $C[0]+=0.5*$V[$period]; if($phase{$j}==0){ $C[1]-=0.5*$V[$period]; }else{ $C[1]+=0.5*$V[$period]; } }elsif($period==2){ if(($phase{$j}==180)||($phase{$j}==3.14)){ $C[0]+=$V[$period]; $C[2]-=$V[$period]; }else{ $C[2]+=$V[$period]; } }elsif($period==3){ $C[0]+=0.5*$V[$period]; if($phase{$j}==0){ $C[1]+=1.5*$V[$period]; $C[3]-=2*$V[$period]; }else{ $C[1]-=1.5*$V[$period]; $C[3]+=2*$V[$period]; } }elsif($period==4){ if(($phase{$j}==180)||($phase{$j}==3.14)){ $C[2]+=4*$V[$period]; $C[4]-=4*$V[$period]; }else{ $C[0]+=$V[$period]; $C[2]-=4*$V[$period]; $C[4]+=4*$V[$period]; } } }
- writeMdpFiles()
Write MDP for test with GROMACS
- writePdb(afile)
Write a new PDB file with the atom names defined by Antechamber
The format generated here use is slightly different from:
respected to atom name. Using GAFF2 atom types. CU/Cu Copper, CL/cl Chlorine, BR/br Bromine
- Parameters
afile ([str]) – file path name
- writePosreFile(fc=1000)
Write file with positional restraints for heavy atoms http://www.mdtutorials.com/gmx/complex/06_equil.html