Hello Charles,
I have written a script for ensemble refinement based on
dna_ensembe/ensemble.py and gb1_rdc/refine.py. The script is along with the
mail below.
When running the simulation, one of the child process stops working and in
the log files I found this line :
barrier: error reading from process 1
The simulation proceeds with computing only one structure.
While the killed process has :
StructureLoop: this process has no work. Exiting...
Can you tell me where am I going wrong ? The below script is executed this
way:
xplor -smp 2 -py -o refine.out 5050.py
xplor.requireVersion("2.17")
numberOfStructures=1
ensembleSize=2 # number of ensemble members
outFilename = "SCRIPT_%d_STRUCTURE_MEMBER.sa" % ensembleSize
import protocol
protocol.initRandomSeed()
command = xplor.command
protocol.initParams("protein")
#load pdb and covalently minimize the protein
protocol.loadPDB("model.pdb")
xplor.simulation.deleteAtoms("not known")
protocol.fixupCovalentGeom(maxIters=1000,useVDW=1)
#initilize ensemble simulations
from ensembleSimulation import EnsembleSimulation
esim = EnsembleSimulation("ensemble",ensembleSize)
#
# a PotList contains a list of potential terms. This is used to specify
which
# terms are active during refinement.
#
from potList import PotList
potList = PotList()
# parameters to ramp up during the simulated annealing protocol
#
from simulationTools import MultRamp, StaticRamp, InitialParams
rampedParams=[]
highTempParams=[]
# orientation Tensor - used with the dipolar coupling term
# one for each medium
# For each medium, specify a name, and initial values of Da, Rh.
#
from varTensorTools import create_VarTensor
media={}
expdatadir = 'xplor5050/'
import os,re
mediaList = open(expdatadir+'mediaList.txt').readlines()
#media[medium] = create_VarTensor(medium)
for mediaName in mediaList:
mediaName = mediaName.strip()
if not mediaName in media.keys():
media[mediaName] = create_VarTensor(mediaName)
from rdcPotTools import create_RDCPot, scale_toNH,correctGyromagneticSigns
correctGyromagneticSigns()
bondTypes = ['NH','CAC','CAHA','CN','CHN']
#xplorBondTypes = {'NH':'NH','CN':'NCO','HNC':'CHN','CAHA':'CAHA',
'CAC':'CACO'}
#bondTypeScales = {'NH':15, 'CAHA':3, 'CN':100, 'CAC':100, 'CHN':3}
bondTypeScales = {'NH':1, 'CAHA':1, 'CN':1, 'CAC':1, 'CHN':1}
rdcs = PotList('rdc')
for medium in media.keys():
#rdc = create_RDCPot("%s_%s"%(medium,expt),file,media[medium])
for btype in bondTypes:
if os.path.isfile(expdatadir+medium+'_'+btype+'.tbl'):
rdc =
create_RDCPot("%s_%s"%(medium,btype),expdatadir+medium+'_'+btype+'.tbl',media[medium])
rdc.setScale(bondTypeScales[btype])
rdc.setAveType("sum") # its set to be average by default
rdc.setShowAllRestraints(1)
rdcs.append(rdc)
potList.append(rdcs)
rampedParams.append( MultRamp(0.05,5.0, "rdcs.setScale( VALUE )") )
# calc. initial tensor orientation
# and setup tensor calculation during simulated annealing
#
from varTensorTools import calcTensorOrientation, calcTensor
for medium in media.keys():
calcTensor(media[medium])
rampedParams.append( StaticRamp("calcTensor(media['%s'])" % medium) )
pass
#add other things about NOE later
from avePot import AvePot
from xplorPot import XplorPot
potList.append( AvePot(XplorPot,'VDW') )
rampedParams.append( StaticRamp("protocol.initNBond()") )
rampedParams.append( MultRamp(0.9,0.8,
"command('param nbonds repel VALUE end
end')") )
rampedParams.append( MultRamp(.004,4,
"command('param nbonds rcon VALUE end end')")
)
potList.append( AvePot(XplorPot,"BOND") )
potList.append( AvePot(XplorPot,"ANGL") )
potList['ANGL'].setThreshold( 5 )
rampedParams.append( MultRamp(0.4,1,"potList['ANGL'].setScale(VALUE)") )
potList.append( AvePot(XplorPot,"IMPR") )
potList['IMPR'].setThreshold( 5 )
rampedParams.append( MultRamp(0.1,1,"potList['IMPR'].setScale(VALUE)") )
# Give atoms uniform weights, except for the anisotropy axis
#
protocol.massSetup()
# IVM setup
# the IVM is used for performing dynamics and minimization in
torsion-angle
# space, and in Cartesian space.
#
from ivm import IVM
dyn = IVM()
dyn.reset()
for m in media.values():
m.setFreedom("varyDa, varyRh") #vary tensor Rh, Da, vary
orientation
protocol.torsionTopology(dyn)
#for ending with cartesian minimzation
minc = IVM()
protocol.initMinimize(minc)
for m in media.values():
m.setFreedom("varyDa, varyRh") #allow all tensor parameters float
here
pass
protocol.cartesianTopology(minc)
# object which performs simulated annealing
#
from simulationTools import AnnealIVM
init_t = 3000. # Need high temp and slow annealing to converge
cool = AnnealIVM(initTemp =init_t,
finalTemp=25,
tempStep =12.5,
ivm=dyn,
rampedParams = rampedParams)
#criteria to accept a structure
def accept(potList):
"""
return True if current structure meets acceptance criteria
"""
#if potList['noe'].violations()>0:
# return False
if potList['rdc'].rms()>1.2: #this might be tightened some
return False
#if potList['CDIH'].violations()>0:
# return False
if potList['BOND'].violations()>0:
return False
if potList['ANGL'].violations()>0:
return False
if potList['IMPR'].violations()>1:
return False
return True
def calcOneStructure(loopInfo):
""" this function calculates a single structure, performs analysis on
the
structure, and then writes out a pdb file, with remarks.
"""
# initialize parameters for high temp dynamics.
InitialParams( rampedParams )
# high-temp dynamics setup - only need to specify parameters which
# differfrom initial values in rampedParams
InitialParams( highTempParams )
# high temp dynamics
#
protocol.initDynamics(dyn,
potList=potList, # potential terms to use
bathTemp=init_t,
initVelocities=1,
finalTime=10, # stops at 10ps or 5000 steps
numSteps=5000, # whichever comes first
printInterval=100)
dyn.setETolerance( init_t/100 ) #used to det. stepsize. default:
t/1000
dyn.run()
# initialize parameters for cooling loop
InitialParams( rampedParams )
# initialize integrator for simulated annealing
#
protocol.initDynamics(dyn,
potList=potList,
numSteps=100, #at each temp: 100 steps or
finalTime=.2 , # .2ps, whichever is less
printInterval=100)
# perform simulated annealing
#
cool.run()
# final torsion angle minimization
#
protocol.initMinimize(dyn,
printInterval=50)
dyn.run()
# final all- atom minimization
#
protocol.initMinimize(minc,
potList=potList,
dEPred=10)
minc.run()
#do analysis and write structure
loopInfo.writeStructure(potList)
pass
from simulationTools import StructureLoop, FinalParams
StructureLoop(numStructures=numberOfStructures,
pdbTemplate=outFilename,
structLoopAction=calcOneStructure,
genViolationStats=1,
averagePotList=potList,
averageTopFraction=0.5, #report only on best 50% of structs
averageAccept=accept, #only use structures which pass
accept()
averageContext=FinalParams(rampedParams),
averageFilename="SCRIPT_ave.pdb", #generate regularized
ave structure
averageFitSel="name CA",
averageCompSel="not resname ANI and not name H*" ).run()
Thanks for your support,
Santhosh
-------------- next part --------------
An HTML attachment was scrubbed...
URL:
<http://dcb.cit.nih.gov/pipermail/xplor-nih/attachments/20130513/834c7127/attachment.html>
