search

leeping / forcebalance

Systematic force field optimization.
Other
146 stars 75 forks source link

crash in CopyPeriodicTorsionParameters #168

Closed mkrompiec closed 4 years ago

mkrompiec commented 4 years ago

Hi, I'm trying to fit some torsions using AbInitio_OpenMM target, and I'm getting this error consistently: Traceback (most recent call last): File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/bin/ForceBalance.py", line 45, in Run_ForceBalance optimizer.Run() File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/optimizer.py", line 319, in Run xk = self.OptTab[self.jobtype]() File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/optimizer.py", line 923, in NewtonRaphson return self.MainOptimizer(b_BFGS=0) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/optimizer.py", line 492, in MainOptimizer data = self.Objective.Full(xk,Ord,verbose=True) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/objective.py", line 294, in Full Objective = self.Target_Terms(vals, Order, verbose, customdir) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/objective.py", line 235, in Target_Terms Ans = Funcs[Order](mvals, customdir=customdir) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/target.py", line 319, in get_H Ans = self.meta_get(mvals,1,1,customdir=customdir) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/target.py", line 584, in meta_get Answer = self.get(mvals, AGrad, AHess) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/abinitio.py", line 1092, in get Answer_EF = self.get_energy_force(mvals, AGrad, AHess) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/abinitio.py", line 698, in get_energy_force dM_all[:,p,:], ddM_all[:,p,:] = f12d3p(fdwrap(callM, mvals, p), h = self.h, f0 = M_all) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/finite_difference.py", line 109, in f12d3p fm1, f1 = [f(ih) for i in [-1, 1]] File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/finite_difference.py", line 109, in fm1, f1 = [f(ih) for i in [-1, 1]] File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/finite_difference.py", line 160, in func1 return func(mvals,*kwargs) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/abinitio.py", line 696, in callM return self.energy_force_transform() File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/abinitio.py", line 516, in energy_force_transform M = self.energy_force_all() File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/abinitio.py", line 509, in energy_force_all return self.engine.energy_force() File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/openmmio.py", line 1000, in energyforce Result = self.evaluate(force=True, traj=True) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/openmmio.py", line 965, in evaluate_ self.update_simulation() File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/openmmio.py", line 887, in update_simulation UpdateSimulationParameters(self.system, self.simulation) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/openmmio.py", line 325, in UpdateSimulationParameters CopySystemParameters(src_system, dest_simulation.system) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/openmmio.py", line 320, in CopySystemParameters Copiersnm File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/forcebalance/openmmio.py", line 274, in CopyPeriodicTorsionParameters dest.setTorsionParameters(i,src.getTorsionParameters(i)) File "/sw-pv/sdk/anaconda-python/2019.03/envs/QUBEKit/lib/python3.6/site-packages/simtk/openmm/openmm.py", line 15072, in setTorsionParameters return _openmm.PeriodicTorsionForce_setTorsionParameters(self, index, particle1, particle2, particle3, particle4, periodicity, phase, k) Exception: Assertion failure at PeriodicTorsionForce.cpp:60. Index out of range

It seems that the torsions do not exist in the "dest" system. A toy example is attached. torsion.zip

leeping commented 4 years ago

Hello Michal,

I see the source of your problem. OpenMM does not create torsion terms in System.PeriodicTorsionForce when the k parameter is equal to zero, presumably because it is cheaper to skip the interaction altogether. As a result, there is an unequal number of terms in the System when a small displacement is applied to k for computing the numerical derivative.

As a workaround, please try setting your initial value of k to a small value such as 1e-6. Because this behavior is confusing, I should look into whether it's possible to have OpenMM create all of the terms even when the force constant is set equal to zero.

Thanks,

mkrompiec commented 4 years ago

Hi Lee-Ping, Thank you, this makes perfect sense! And changing the k values to nonzero fixes the issue.