[ase-users] Parallelizing over NEB images
Gaël Donval
G.Donval at bath.ac.uk
Tue Nov 21 17:05:30 CET 2017
Hi Divya,
> Dear ASE users,
>
> I am trying to perform my first CI-NEB calculation where I have 4
> images between the initial and final images and the optimization of
> the images is parallelized (script below). I would greatly appreciate
> your help on the below 2 questions:
>
>
> 1. What does the fmax for the optimizer specify exactly? Is
> this taken as the curvature at the initial and final states of the
> NEB interpolation? How does NEB determine the convergence criterion
> for the each of the images on the minimum energy path and where can I
> see this?
The `fmax` keyword is always a maximum force of some sort. To go into
more details, your LBFGS optimizer is going to call the `get_forces`
method to whatever object is provided (in your case a NEB object),
calculate their norm and check if none of them is greater than `fmax`.
For an NEB object, `get_forces` is going to gather all the forces from
all the inner images (that is, all images except the end-points --that
are untouched).
So the `fmax` you provide *is* the convergence criterion for all
images: the images are interdependent; they are coupled; "converging"
one by itself makes little sense.
>
>
> 2. When I read the log file, it seems like instead of logging
> the energy and max force for each image, it is logging the potential
> energy and max force for the image with the highest energy. I want to
> confirm that despite this, each of the images are separately
> optimized towards the convergence criteria and that this logging does
> not affect the optimization itself. Also, how can I make the log the
> forces and energies of each image separately for it to make more
> sense?
The energy in your logfile is indeed the maximum energy. What NEB is
usually used for is determining an energy barrier providing known end-
points. In that context, providing the maximum energy point on the NEB
path as a summary of the advance of the algorithm is the only thing
that really makes sense. Though this is just a log file, a summary: all
the other inner images are constantly optimised too.
Actually, you should really think of all the images as a single system.
Regarding local information, you can extract it either from the
trajectory file generated by LBGFS (provided you gave a value to its
`trajectory` argument); you could extract that from the images
themselves or you could use (and study if you need something more)
`ase.neb.NEBTools`.
Cheers,
Gaël
>
> Thanks in advance.
>
> Best regards,
> Divya
>
> from __future__ import print_function
> from ase.constraints import FixAtoms
> from ase.optimize import BFGS, QuasiNewton, LBFGS, MDMin
> from ase.lattice.cubic import FaceCenteredCubic
> from ase.build import fcc110, add_adsorbate
> from gpaw import GPAW, PW, FermiDirac, Mixer
> from ase.io import write, read,Trajectory
> from ase import Atoms, Atom
> from ase.neb import NEB
> from ase.optimize.fire import FIRE
> from ase.parallel import rank, size
> import numpy as np
>
> calc = GPAW('Ag110_2x3x6_BEEF.gpw')
> CO2_slab=read('CINEB_Ag110_2x3x6_CO2-H_reference.traj')
> COOH_slab=read('Ag110_2x3x6_COOHupLbridge_rot_BEEF.traj')
>
> initial = CO2_slab
> final = COOH_slab
>
> constraint = FixAtoms(mask=[atom.tag > 2 for atom in initial])
>
> n = size // 4 # number of cpu's per image
> j = 1 + rank // n # my image number
> assert 4 * n == size
>
> images = [initial]
>
> for i in range(4):
> ranks = range(i * n, (i + 1) * n)
> image = initial.copy()
>
> if rank in ranks:
> calc=GPAW(mode=PW(450), xc='BEEF-vdW', basis='dzp',
> gpts=(48,48,160), kpts=(3,3,1), txt='CINEB_Ag110_2x3x6_CO2-
> COOH_Tafel_neb%d.txt' % j, occupations=FermiDirac(0.1),
> mixer=Mixer(beta=0.1, nmaxold=5, weight=50.0), communicator=ranks)
> image.set_calculator(calc)
>
> image.set_constraint(constraint)
> images += [image]
>
> images += [final]
>
> neb = NEB(images, parallel=True, climb=True)
> neb.interpolate('idpp')
> #Run NEB calculation
> opt = LBFGS(neb, logfile='CINEB_Ag110_2x3x6_CO2-COOH_Tafel.log')
> opt.run(fmax=0.05)
>
> for i in range(1, 5):
> opt.attach(io.write('CINEB_Ag110_2x3x6_CO2-COOH_Tafel-%d.cif' %
> i, 'w', images[i]))
> opt.attach(io.Trajectory('CINEB_Ag110_2x3x6_CO2-COOH_Tafel-
> %d.traj' % i, 'w', images[i]))
>
>
> #Logfile
>
> LBFGS: 0 18:17:31 -5379.937600 6.8668
> LBFGS: 1 18:44:44 -5380.096707 6.1825
> LBFGS: 2 19:11:00 -5380.155402 4.4230
> LBFGS: 3 19:39:07 -5380.249108 3.3944
> LBFGS: 4 20:07:05 -5380.368769 3.2199
> LBFGS: 5 20:35:12 -5380.521739 3.4305
> LBFGS: 6 21:03:11 -5380.727266 3.9711
> LBFGS: 7 21:31:13 -5380.953959 4.5234
> LBFGS: 8 21:57:42 -5381.208836 5.0993
> LBFGS: 9 22:24:52 -5381.490995 5.6717
> LBFGS: 10 22:51:58 -5381.801246 6.2453
> LBFGS: 11 23:19:50 -5382.139225 6.7877
> LBFGS: 12 23:47:46 -5382.503905 7.2781
> LBFGS: 13 00:14:56 -5382.896339 7.7335
> LBFGS: 14 00:43:06 -5383.314196 8.1039
> LBFGS: 15 01:12:07 -5383.755902 8.3441
> LBFGS: 16 01:42:02 -5384.209731 8.3273
> LBFGS: 17 02:12:50 -5384.660808 7.9158
> LBFGS: 18 02:44:25 -5385.012915 7.1416
> LBFGS: 19 03:14:20 -5385.247200 6.9705
> LBFGS: 20 03:45:01 -5385.395184 7.2958
> LBFGS: 21 04:15:33 -5385.502550 7.5201
> LBFGS: 22 04:46:22 -5385.588210 7.6800
> LBFGS: 23 05:17:00 -5385.661367 7.7952
> LBFGS: 24 05:47:43 -5385.726172 7.8750
> LBFGS: 25 06:18:26 -5385.784934 7.9238
> LBFGS: 26 06:49:07 -5385.839004 7.9446
> LBFGS: 27 07:20:34 -5385.889114 7.9381
> LBFGS: 28 07:52:51 -5385.935424 7.9049
> LBFGS: 29 08:25:06 -5385.974499 7.8480
> LBFGS: 30 08:57:26 -5386.006282 7.7645
> LBFGS: 31 09:29:42 -5386.029139 7.6479
> LBFGS: 32 10:02:02 -5386.039673 7.4886
> LBFGS: 33 10:35:16 -5386.032930 7.2809
> LBFGS: 34 11:07:37 -5386.006100 7.0377
> LBFGS: 35 11:39:01 -5385.960810 6.7818
> LBFGS: 36 12:11:23 -5385.900440 6.5308
> LBFGS: 37 12:44:37 -5385.829873 6.2943
> LBFGS: 38 13:17:49 -5385.753279 6.0716
> LBFGS: 39 13:50:59 -5385.675096 6.1741
> LBFGS: 40 14:23:24 -5385.599475 6.4399
> LBFGS: 41 14:54:54 -5385.530148 6.6489
> LBFGS: 42 15:28:03 -5385.469514 6.8101
>
>
> _______________________________________________
> ase-users mailing list
> ase-users at listserv.fysik.dtu.dk
> https://listserv.fysik.dtu.dk/mailman/listinfo/ase-users
More information about the ase-users
mailing list