Dear Charles,
thanks for your swift replies!
The maxDisplacement did wonders for the calculation, although it needs 35
Angstrom for the set of Tb-induced PCSs. Is there any way to optimize the
tensor such that
all Lanthanides are at the same position? Right now I use following input:
lnAtom = AtomSel('resid 200 and name TB+3')[0]
def init_tensor(name, pos, ln):
t = create_VarTensor(name)
t.setDa(abs(get_X_components(ln)['dXax']) * 1e30)
t.setRh(abs(get_X_components(ln)['dXrh']) * 1e30)
return t
def fmt_tensor(t):
dXax = "dXax: %7.4e" % t.Da()
dXrh = "dXrh: %7.4e" % t.Rh()
xyz = "(x,y,z): %s" % lnAtom.pos()
return ', '.join([dXax, dXrh, xyz])
LANTHANIDES = ['Tb', 'Ho', 'Yb', 'Tm']
pcsSpecs = [(f'{ln}PCS', f'expdata/{ln}.tbl',
init_tensor(f'{ln}PCS', lnAtom, ln)) for ln in LANTHANIDES]
pcsPots = PotList('pcs')
for (name, pcs_list, tensor) in pcsSpecs:
pcsPot = rdcPotTools.create_RDCPot(name, pcs_list, oTensor=tensor)
pcsPot.setUseDistance(True)
pcsPot.setScale(1)
pcsPot.setShowAllRestraints(True)
pcsPot.setVerbose(True)
tensors[name] = tensor
pcsPots.append(pcsPot)
for pot in pcsPots:
print_header(f"{pot.instanceName()} tensor calculation")
calcXTensor(pot.oTensor, maxDisplacement=35)
print(fmt_tensor(pot.oTensor))
giving me the following output:
=========================== TbPCS tensor calculation
===========================
Current tensor paramaters are: [7.85, 6.337, -13.026, 0.42100000000000004,
0.047151999999999986, 360.0, 360.0, 360.0]
Opt failure. Trying with different metal coords
> maxDisplacement. Trying with different metal coords
Opt failure. Trying with different metal coords
Opt failure. Trying with different metal coords
Optimized tensor paramaters: [20.048839814743694, -14.238790696623267,
12.497844758778207, 6.440972097357804, 3.2500711156689386, 36.70794394974678,
61.71489965526524, 150.59060317636602]
dXax: 6.4410e+00, dXrh: 5.0459e-01, (x,y,z): [ 20.049, -14.239, 12.498]
=========================== HoPCS tensor calculation ===========================
Current tensor paramaters are: [20.048839814743694, -14.238790696623267,
12.497844758778207, 0.18500000000000003, 0.01073000000000001, 360.0, 360.0,
360.0]
Optimized tensor paramaters: [12.20767939491271, -3.0928574866850096,
10.45786248651365, -0.2613893990138944, -0.15388824548624594,
143.81489717326673, 69.85647207953299, 40.5484721773534]
dXax: -2.6139e-01, dXrh: 5.8873e-01, (x,y,z): [ 12.208, -3.093, 10.458]
=========================== YbPCS tensor calculation
===========================
Current tensor paramaters are: [12.20767939491271, -3.0928574866850096,
10.45786248651365, 0.083, 0.004814000000000005, 360.0, 360.0, 360.0]
Optimized tensor paramaters: [13.05731929409086, -1.4633119142431874,
9.026098800357362, -0.15565782500677988, -0.0536901499956898,
85.84530366261492, 133.59986034415567, 56.780559845445055]
dXax: -1.5566e-01, dXrh: 3.4492e-01, (x,y,z): [ 13.057, -1.463, 9.026]
=========================== TmPCS tensor calculation
===========================
Current tensor paramaters are: [13.05731929409086, -1.4633119142431874,
9.026098800357362, 0.21899999999999997, 0.04401899999999999, 360.0, 360.0,
360.0]
Opt failure. Trying with different metal coords
Optimized tensor paramaters: [18.55987618612005, -8.735037801651153,
6.7649054404504625, 0.4318129276333319, 0.2578006080978475, 151.9353062239399,
75.5710173982344, 26.013551062651004]
dXax: 4.3181e-01, dXrh: 5.9702e-01, (x,y,z): [ 18.560, -8.735, 6.765]
I also wonder why the axial/rhombic components of the tensor are not always the
same as the ones reported by the calculation.
Now to the real problem, the simulated annealing. My annealing function starts
like this:
def anneal(loopInfo):
global ramps
InitialParams(ramps)
# ---------------- ramp up experimental restraints ---------------- #
for i in range(0, expramp_output_steps):
scale = scalepcs_for_time(expramp_duration * i /
expramp_output_steps)
print_header("PCS now scaled at %5.3f (%i/%i)" %
(scale, i, expramp_output_steps))
protocol.initDynamics(torsionDynamics,
potList=hightempPotList,
bathTemp=init_temp,
initVelocities=1,
finalTime=hightemp_timestep,
printInterval=100)
torsionDynamics.setETolerance(init_temp / 100)
torsionDynamics.run()
Which gives the output:
======================== PCS now scaled at 0.000 (0/10)
========================
*--- Dynamics ---- step= 0 ---- time= 0 ---- delta_t=
0.001 --*
| E(kin)+E(poten)= 12220.585 E(kin)= 894.087 temperature=
952.222 |
| E(poten)= 11326.4981465 grad= 52.0421334 ANGL=
181.3980536 |
| BOND= 47.9098246 CDIH= 0.0000000 IMPR=
13.7081550 |
| TbCoord= 11083.4821133 pcs= 0.0000000 xrdRMSD=
0.0000000 |
*------------------------------------------------------------------------------*
calcD_G: singular D matrix: { { 0 } }
H matrix: { { 0.749263, 0.556428, 0.359155, 0, 0, 0 } }
node level: 2
number of children: 1
atoms:
2545 ANI 5000 Z : { 0.727034, 0.536782, 0.339489 }
InternalDynamics::step: caught assertion failure:
calcD_G: singular D matrix. Bad topology?
Traceback (most recent call last):
File "/opt/xplor-nih-3.0.3/python/simulationTools.py", line 574, in run
s.structLoopAction(s)
File "./run.py", line 380, in anneal
torsionDynamics.run()
File "/opt/xplor-nih-3.0.3/python/ivm.py", line 569, in run
(done,s.stepsize_) = PublicIVM.step(s,s.stepsize());
File "/opt/xplor-nih-3.0.3/python/wrappers/publicIVM.py", line 202, in
step
return _publicIVM.PublicIVM_step(self, *args, **kwargs)
_publicIVM.IVMError: IVM error: calcD_G: singular D matrix. Bad topology?
Error calculating structure 1: IVM error: calcD_G: singular D matrix. Bad
topology?
Thanks again for your help, and I wish a good start into the new year!
All the best,
Till
########################################################################
To unsubscribe from the XPLOR-NIH list, click the following link:
http://list.nih.gov/cgi-bin/wa.exe?SUBED1=XPLOR-NIH&A=1
