zadorlab / sella

A Python software package for saddle point optimization and minimization of atomic systems.
https://www.ecc-project.org/
Other
64 stars 19 forks source link

constrains in IRC #4

Open GengSS opened 4 years ago

GengSS commented 4 years ago

Hi,

I hope I find the right place to ask a question. I would like to use the IRC in Sella.

The transition state is located by another algorithm with a simple EMT force field, which is simply an Al atom migration on Al surface and some of the Al atoms in the slab are fixed.

The coordinates of TS are: Al 1.0000000000000000 8.1227999999999998 0.0000000000000000 0.0000000000000000 0.0000000000000000 5.7436869622220881 0.0000000000000000 0.0000000000000000 0.0000000000000000 13.7436869622220872 Al 17 Selective dynamics Cartesian 0.0000000000000000 0.0000000000000000 5.4359217405555222 F F F 2.0306999999999999 1.4359217405555220 4.0000000000000000 F F F 0.0777722201139265 -0.0018913176189019 8.2945732741720111 T T T 2.0307032040535633 1.4648076410935980 6.8798630655734492 T T T 0.0000000000000000 2.8718434811110440 5.4359217405555222 F F F 2.0306999999999999 4.3077652216665658 4.0000000000000000 F F F -0.3262762671500816 2.8700209572327453 8.2319835045797269 T T T 2.0307032101652940 4.2757822283465670 6.8803484683802045 T T T 4.0613999999999999 0.0000000000000000 5.4359217405555222 F F F 6.0921000000000003 1.4359217405555220 4.0000000000000000 F F F 3.9836224186428559 -0.0018913482524109 8.2945752309502883 T T T 6.0921191920564599 1.3672281815067406 6.7658911068951042 T T T 4.0613999999999999 2.8718434811110440 5.4359217405555222 F F F 6.0921000000000003 4.3077652216665658 4.0000000000000000 F F F 4.3876211260504459 2.8700209967848611 8.2320230173653748 T T T 6.0921190498270379 4.3742831571633722 6.7661580694779211 T T T 2.0307321321540455 2.8672631794380150 9.3982316719762213 T T T

Here is my script for IRC with Sella,

from ase.calculators.emt import EMT from ase.io import read from sella.irc import IRC from sella.peswrapper import PESWrapper atoms=read("POSCAR")

fixed_indices = [i for i in atoms.constraints[0].get_indices()] a = atoms.copy() a.set_constraint() a.set_calculator(EMT())

irc=IRC(a, trajectory="irc.traj",constraints=dict(fix=fixed_indices)) irc.run()

While, unfortunately, I can not run it smoothly. Here is the error information.

Traceback (most recent call last): File "run_irc.py", line 13, in irc.run() File "/home/gengsun/Downloads/test_sella/sella-master/selle_build/lib.linux-x86_64-3.7/sella/irc.py", line 165, in run for converged in self.irun(fmax, steps, direction): File "/home/gengsun/Downloads/test_sella/sella-master/selle_build/lib.linux-x86_64-3.7/sella/irc.py", line 147, in irun self._calc_v0ts() File "/home/gengsun/Downloads/test_sella/sella-master/selle_build/lib.linux-x86_64-3.7/sella/irc.py", line 134, in _calc_v0ts Hw = self.H0 / np.outer(self.sqrtm, self.sqrtm) ValueError: operands could not be broadcast together with shapes (27,27) (51,51)

I wonder did I impose the constraints correctly? or this is a bug in the code?

Thank you very much.

Best Wishes,

Geng

ehermes commented 4 years ago

Unfortunately, I did not consider constraints when writing the IRC module. Currently it is not possible to perform constrained IRC. I will fix this eventually, but there are other higher priority tasks that need to be completed first.

ehermes commented 4 years ago

On the internal branch, I've added the ability to perform IRC calculations with constraints. Fixing atoms should work, but I have no idea whether other kinds of constraints (e.g. fixing bond distances or angles) will work with IRC. I'll be merging the internal branch into master eventually, after I have verified that there are no regressions.

If you wish to test the internal branch, please modify your script to look like the following:

from ase.calculators.emt import EMT
from ase.io import read
from sella import IRC

atoms = read("POSCAR")

atoms.calc = EMT()

irc = IRC(atoms, trajectory="irc.traj")
irc.run()

The constraints will be automatically constructed from the ASE FixAtoms object attached to atoms, which should be determined automatically from the Selective Dynamics flags in your POSCAR.

Let me know if you run into any problems.

ehermes commented 2 years ago

I have updated the code and hopefully fixed many of the IRC bugs we had before (also @PattanaikL). Currently it appears possible to perform IRC calculations on systems with fixed-atom Constraints, though other constraints (e.g. fixed bond distances, fixed bending angles, fixed dihedral angles) still do not work correctly with IRC. Please try your examples again with the current git master branch and let me know whether or not your issue has been solved. I will probably be making a new release with these fixes soon.