Help with DFT and Cosmo
Clicked A Few Times
2:06:17 PM PDT - Tue, Apr 10th 2012
I'm trying to get started with DFT optimization using Cosmo in NWChem v6.1. I've used NWChem's DFT in the past, but not with Cosmo.

I'm currently trying to run the charge neutral and -1 states of ethylene carbonate, in a dielectric constant of around 40.

Here's an example of what's going on for EC-:
http://imgur.com/3qNr4
The energy drops in a few discrete jumps, but after the jumps, the oscillations never damp out and converge.

I'm probably doing something stupid, since, like I said, I haven't used Cosmo in NWChem before. I've run a lot of Jaguar cases using their Poisson-Boltzmann solvation, and there the geometry converges after 20 or so geometry steps.

Here is the input file I'm using. Can anyone give me some pointers to making this work better?

Thanks!

Rick


start tmp_ec-

title "EC-"
charge -1

geometry units angstroms print xyz
   symmetry group c1

   O    0.0000000000    0.0000000000   -2.2090362900

   C    0.0000000000    0.0000000000   -0.9640410700

   O    0.0001589200    1.1411701400   -0.2829492300

   C   -0.0012770200    0.7024679000    1.5854052800

   C    0.0012770200   -0.7024679000    1.5854052800

   O   -0.0001589200   -1.1411701400   -0.2829492300

   H   -0.9178669200    1.2411461900    1.8072578800

   H    0.9127606700    1.2448048800    1.8088623600

   H   -0.9127606700   -1.2448048800    1.8088623600

   H    0.9178669200   -1.2411461900    1.8072578800
end

basis
 * library 6-31G**
end

cosmo
dielec 40.0

rsolv 0.5
end

dft
 xc b3lyp

 mult 2

 iterations 100
end

task dft optimize

Forum Regular
3:14:07 PM PDT - Thu, Apr 12th 2012
Hi Rick,

I have tried running your case myself and for me it does not converge either. Now COSMO in NWChem does not have analytical gradients and the gradients for the geometry optimization are obtained from finite differences. This can introduce some numerical issues. So I want to test a few things to see if I can come up with a recipe that works.

Best wishes,

    Huub

Forum Regular
10:32:41 AM PDT - Wed, Apr 18th 2012
Some progress
Hi Rick,

Initially I thought that perhaps various accuracy aspects of the energy expressions and numerical differentiation play into this problem. So I tried somethings the crank up the accuracy of the energy (finer DFT grids) as well as the gradient evaluation (using a 5 point expression for the numerical gradients) but all to no effect. In particular I noticed that the geometry kept changing rather a lot from point to point.

It seems that this problem is caused by the fact that the system you are studying essentially consists of 2 non-bonded fragments. The Hessian along the coordinates connecting the 2 fragments is virtually 0 causing the optimizer to take rather large steps. So I am currently looking at ways to reduce the coordinate updates in the hope that that will allow the calculation to converge.

Best wishes,

Huub

Forum Regular
11:16:45 AM PDT - Fri, Apr 20th 2012
Problem resolution
Hi Rick,

In order to address this problem I have changed the geometry optimizer a little bit. The problem was that for this system the are some small eigenvalues in the Hessian. On the first point the smallest eigenvalue of the Hessian is 0.4e-3 and along the optimization this declines further to 0.3e-5. Obviously in this context using a Newton-Raphson update of your coordinates is not a good idea as -g/H becomes large due to the developing singularity. So I changed the optimizer to allow me to specify a minimum Hessian eigenvalue, below which the algorithm would switch to using steepest decent for the coordinates corresponding to such eigenvalues. This approach resulted in a much more stable optimization behavior (not exactly fast though). This led to a structure where one of the carbonate oxygens attaches to one of the ethene carbons. At that point it seems that something happens to the electronic structure as the SCF runs into trouble and fails to converge. I have also tried our newly developed second order SCF for DFT calculations but this does not seem to be properly integrated with COSMO yet.

So I think we now know what the issues are. However, at the present time there does not seem to be a way to automatically complete this calculation successfully. If you want I can send you the results of the calculations I have done, just let me know.

Best wishes,

Huub

Just Got Here
7:09:20 AM PDT - Sun, Apr 22nd 2012
Gradient Module
I'm doing a DFT geometry optimization with a neutral molecule containing 51 atoms. Without cosmo, the driver converge in 11 steps. The output gives for each step the energy after the symbol @. But before, there is a table called DFT energy gradient. When the cosmo block is added:

basis spherical
  *  library  6-31G*
end

cosmo
  dielec 32.6
  rsolv 0.50
end
set cosmo:map cosmo.par

dft
  xc  b3lyp
end

set driver:linopt 0
driver
  maxiter 50
end

task dft optimize


The DFT gradient energy section never appear. The energy profile for each step (which is given without the @ symbol, but with the word "step") look like the one gave by Rpmuller. Does these problem are related?

Clicked A Few Times
7:48:40 AM PDT - Mon, Apr 23rd 2012
Quote:Huub Apr 20th 6:16 pm
Hi Rick,

So I think we now know what the issues are. However, at the present time there does not seem to be a way to automatically complete this calculation successfully. If you want I can send you the results of the calculations I have done, just let me know.

Best wishes,

Huub


Huub,

Thanks very much for looking into this. I need to think about this a little bit more to see whether I'm making bad assumptions about the system, or giving it a bad geometry or something. Can't imagine why anything out of the ordinary would be happening here.

(Pause as he looks at what he sent it.)

I sent you the harder of the two geometries I was looking at with solvation, the anion of ethylene carbonate, which does (at some point) decompose after accepting an electron. Jaguar has this being weakly bound (and converge-able) after accepting the electron, but when the decomposition occurs doesn't seem clear to me, and I was hoping to hit this with a higher level of theory in NWChem.

For the problems you looked at, it seems like the unstable geometry makes this an unfair case. However, I get something very similar with neutral ethylene carbonate, oddly enough, which is perfectly stable in ethylene carbonate solution.

FYI, this is the electrolyte in Li-ion batteries, and the decomposition here is what forms the solid-electrolyte interphase, which passivates the anode against thermal runaway.


Forum >> NWChem's corner >> General Topics