Hi,
I have been using a script to refine a NMR structure using NOEs, RDCs, SAXS and
other data without problems until recently. Now, the same script results in a
Segmentation fault (core dumped) error. After doing some trouble shooting, I
noticed that the program runs fine when I remove the RDC part in the script,
but I don't know what the issue is. I use Xplor-NIH in NMRBOX. Here is the code
and the last lines in the output.
Thanks for the help,
Cristian
###########################################################################
import psfGen
outfilename = 'SCRIPT_STRUCTURE.pdb'
import protocol
protocol.initRandomSeed(3422) # by specific seed number
# Load paramaters.
protocol.initParams(files=['nucleic'])
psfGen.pdbToPSF(pdb_filename)
protocol.initCoords(pdb_filename)
xplor.simulation.deleteAtoms("not known")
from potList import PotList
potList = PotList()
crossTerms = PotList()
from simulationTools import StaticRamp, MultRamp, InitialParams, AnnealIVM
highTempParams = []
rampedParams = []
#
# Set up distance restraint potential
#
import noePotTools
noe = PotList('noe')
for (name, scale, table) in [('all', 1, noe_tbl),
]:
pot = noePotTools.create_NOEPot(name, table)
pot.setScale(scale)
noe.append(pot)
potList.append(noe)
rampedParams.append(MultRamp(2, 50, "noe.setScale( VALUE )"))
import xplorPot
# Set up dihedral angles
protocol.initDihedrals(dihedral_tbl,
#useDefaults=False # by default, symmetric sidechain
# restraints are included
)
potList.append( xplorPot.XplorPot('CDIH') )
highTempParams.append( StaticRamp("potList['CDIH'].setScale(10)") )
rampedParams.append( StaticRamp("potList['CDIH'].setScale(200)") )
#
# Set up potential for base-pair planarity restraints.
#
protocol.initPlanarity(plane_tbl)
potList.append(xplorPot.XplorPot('PLAN'))
# Set up statistical torsion angle potential (torsionDB).
#
import torsionDBPotTools
torsiondb = torsionDBPotTools.create_TorsionDBPot(name='torsiondb',system='rna')
potList.append(torsiondb)
rampedParams.append(MultRamp(0.5, 4, "torsiondb.setScale(VALUE)"))
# Setup parameters for atom-atom repulsive term (van der Waals-like term).
#
from repelPotTools import create_RepelPot,initRepel
repel = create_RepelPot('repel')
potList.append(repel)
rampedParams.append( StaticRamp("initRepel(repel,use14=False)") )
rampedParams.append( MultRamp(.004, 4, "repel.setScale( VALUE)") )
# nonbonded interaction only between C1' atoms
highTempParams.append( StaticRamp("""initRepel(repel,
use14=True,
scale=0.004,
repel=1.2,
moveTol=45,
interactingAtoms="name C1'"
)""") )
# Selected 1-4 interactions.
import torsionDBPotTools
repel14 = torsionDBPotTools.create_Terminal14Pot('repel14')
potList.append(repel14)
highTempParams.append(StaticRamp("repel14.setScale(0)"))
rampedParams.append(MultRamp(0.004, 4, "repel14.setScale(VALUE)"))
# Set up bond length potential.
potList.append(xplorPot.XplorPot('BOND'))
# (The setup of this term remains unchanged throughout; no need to involve
# highTempParams and/or rampedParams.)
#
# Set up bond angle potential.
potList.append(xplorPot.XplorPot('ANGL'))
rampedParams.append(MultRamp(0.4, 1.0, "potList['ANGL'].setScale(VALUE)"))
#
# Set up improper dihedral angle potential.
potList.append(xplorPot.XplorPot('IMPR'))
rampedParams.append(MultRamp(0.1, 1.0, "potList['IMPR'].setScale(VALUE)"))
###RDC
from varTensorTools import create_VarTensor
from rdcPotTools import create_RDCPot, scale_toNH
from varTensorTools import calcTensorOrientation, calcTensor
media={}
for medium, Da, Rh in tensor_list:
#print(medium,Da, Rh, tensor_list)
oTensor = create_VarTensor(medium)
oTensor.setDa(Da)
oTensor.setRh(Rh)
oTensor.setFreedom ("varyDa, varyRh" )
#oTensor.setFreedom ("fixDa, fixRh" )
media[medium] = oTensor
highTempParams.append(StaticRamp("""
for medium in media.values():
calcTensorOrientation(medium)
""") )
rdcs = PotList('rdc')
for medium, expt, rdc_file, scale in rdc_exp_list:
#print(experiments_list, medium, expt, rdc_file, scale, media_dictionary)
rdc = create_RDCPot("%s_%s"%(medium,expt), file=rdc_file,
oTensor=media[medium])
#rdc.setScale(scale)
scale_toNH(rdc)
rdc.setShowAllRestraints(1) #all restraints are printed during analysis
rdc.setThreshold(1.5) # in Hz
rdcs.append(rdc)
potList.append(rdcs)
rampedParams.append( MultRamp(0.01,1, "rdcs.setScale( VALUE )") )
### SAXS term ###########
#from solnXRayPotTools import create_solnXRayPot
import solnXRayPotTools
xray=solnXRayPotTools.create_solnXRayPot('xray',
experiment=saxs_data, #data file with columns q, I, err
numPoints=100, #specifies the number of datapoints to
sample
normalizeIndex=-3, #specifies which grid point to use
to normalize data. -3 specifies normalization which minimizes the Chi^2 value.
preweighted=False)
xrayCorrect=solnXRayPotTools.create_solnXRayPot('xray-c',
experiment=saxs_data,
numPoints=100, #should be the same as 'xray'
normalizeIndex=-3,
preweighted=False)
solnXRayPotTools.useGlobs(xray) #uses atom globbing aproximation
xray.setNumAngles(100) #number of angles in solid angle averaging
xrayCorrect.setNumAngles(500) # use a large number for correction term
xray.setScale(50) #restrint weight value
xray.setCmpType("plain") #'plain'
potList.append(xray) #add to potlist energy terms
crossTerms.append(xrayCorrect)
from solnScatPotTools import fitParams
rampedParams.append(
StaticRamp(
"fitParams(xrayCorrect);xray.calcGlobCorrect(xrayCorrect.calcd())",
stride=100))
# Give atoms uniform weights, except for anisotropy axes (if any).
#
protocol.massSetup()
#
# Set up IVM object(s).
#
# IVM object for torsion-angle dynamics/minimization.
import ivm
dyn = ivm.IVM()
protocol.torsionTopology(dyn, flexRiboseRing='resid 1:24')
### Optional IVM object for final Cartesian minimization.
minc = ivm.IVM()
protocol.cartesianTopology(minc)
#
# Temperature set up.
#
temp_ini = 3000.0 # initial temperature
temp_fin = 25.0 # final temperature
def calcOneStructure(loopInfo):
"""Calculate a structure.
"""
# Fix up covalent geometry.
# (The torsion restraints may include ring torsions and distort geometry.)
while True:
try:
protocol.fixupCovalentGeom(maxIters=100, useVDW=1)
break
except protocol.CovalentViolation:
pass
#
# High Temperature Dynamics Stage.
#
# Initialize parameters for high temperature dynamics.
InitialParams(rampedParams)
InitialParams(highTempParams) # purposedly overides some
# setups in rampedParams
# Set up IVM object and run.
protocol.initDynamics(dyn,
potList=potList,
bathTemp=temp_ini,
initVelocities=1,
finalTime=15, # run for finalTime or
numSteps=15001, # numSteps * 0.001, whichever is less
printInterval=100)
dyn.setETolerance(temp_ini/100) # used to find stepsize (default: temp/1000)
dyn.run()
#
# Simulated Annealing Stage.
#
# Initialize parameters for annealing.
InitialParams(rampedParams)
# Set up IVM object for annealing.
protocol.initDynamics(dyn,
potList=potList,
finalTime=0.63, # run for finalTime or
numSteps=631, # numSteps * 0.001, whichever is less
printInterval=100)
# Set up cooling loop and run.
AnnealIVM(initTemp=temp_ini,
finalTemp=temp_fin,
tempStep=12.5,
ivm=dyn,
rampedParams=rampedParams).run()
#
# Torsion angle minimization.
#
protocol.initMinimize(dyn,
potList=potList,
printInterval=50)
dyn.run()
#
# Cartesian minimization (optional).
#
protocol.initMinimize(minc,
potList=potList,
dEPred=10)
minc.run()
from simulationTools import StructureLoop
StructureLoop(numStructures=nstructures,
#pdbFilesIn=infilename,
pdbTemplate=outfilename,
doWriteStructures=True,
structLoopAction=calcOneStructure,
# Arguments for generating structure statistics:
genViolationStats=True,
averageSortPots=[potList['noe'], # terms for structure sorting
potList['PLAN'],
potList['xray'],
potList['CDIH'],
potList['rdc']],
averageTopFraction=FRACTION, # top fraction of structs. to
report on
averageFilename="SCRIPT_ave.pdb",
averagePotList=potList, # terms analyzed
averageFitSel='not (name H* or resname ANI)', # selection to fit
).run() # to average struct.
# and report precision
#################################################################
Log:
GU12_ref_RDC_SAXS.py(373): StructureLoop(numStructures=nstructures,
StructureLoop: calculating structure 0
GU12_ref_RDC_SAXS.py(306): try:
GU12_ref_RDC_SAXS.py(307): protocol.fixupCovalentGeom(maxIters=100,
useVDW=1)
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
GU12_ref_RDC_SAXS.py(308): break
GU12_ref_RDC_SAXS.py(317): InitialParams(rampedParams)
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
[stdin](1): xplor.execfile('GU12_ref_RDC_SAXS.py')
Segmentation fault (core dumped)
########################################################################
To unsubscribe from the XPLOR-NIH list, click the following link:
http://list.nih.gov/cgi-bin/wa.exe?SUBED1=XPLOR-NIH&A=1
