Geometry optimization with constrained DFT

Does anyone have any tips for doing geometry optimizations with constrained DFT? It seems that trying to constrain some molecules to oxidized states on a carbon surface model and relaxing geometry is difficult. The total energy seems to not track smoothly to a lower energy and frequently goes up. Is this just a matter of needing higher thresholds all the way around than for a non constrained optimization or is there a simpler fix? I know this post is a little vague but any tips on geometry optimizations with constrained DFT are greatly appreciated.

This may be a bug in the constrained DFT code. We have reported it as a ticket (#2347). For a test job, analytic gradients do not agree with finite-difference results (IDERIV=0).

Just an update to say that we are looking into this. Errors between analytic and finite-diff gradient are not huge (~ 1e-4 in a.u.), but larger than I’m accustomed to seeing, and we see the behavior that you describe, i.e., energy does not go down consistently in some cases. Bottom line is that I think you can continue using the gradients for now, but you are not imagining things or doing something wrong, so far as I can tell.

I’m no longer 100% convinced this is a bug (though still investigating). As a consequence of looking into this, however, Evgeny @ Q-Chem discovered 3 undocumented job control variables, which we will add to the manual for the 5.3.2 bug-fix release but I can tell you here. They are:
CDFT_THRESH ==> threshold for convergence of CDFT weights
CDFT_MAXITER ==> max number of iterations for converging the CDFT weight
CDFT_PRINT ==> enable printing of CDFT iterations

The most important one in the present context is the first one, which defaults to 5 and that was giving erratic behavior. A better choice may be to set it to 8 for geometry optimizations, which is the default value for SCF_CONVERGENCE in the case of a geometry optimization. For the example that we were looking at, this makes things look a lot better. Please try setting CDFT_THRESH = 8 in your $rem section and (if you don’t mind) reporting back whether or not this solves your problem.

-John

Hi John,

I tried using the higher threshold but still the total energy climbs quite quickly for my test case. I will try another test case and see what that does.

thanks,
Jason

Hi John,

You might find it useful that it seems that everything works okay if doing a MD run but it is when doing an OPT run the energy keeps rising.

thanks,
Jason

Is your molecule top secret or do you mind posting an input file that I could submit as a bug ticket? That way the Q-Chem scientific staff (of which I am not) could take a look. And probably I would look too.

Hi John,

Do you know what would cause the Becke populations not to add up to the constrained value?

thanks,
Jason

Poor quality grid. Note there is an XC_GRID for c-DFT even if you use Hartree-Fock theory. Try something like XC_GRID = 000099000590, i.e., (99,590). That is “high quality” so I would hope is good enough. The “kitchen sink” grid that we use as a stand-in for grid limit is (250,974).