ONIOM run through xTB and CREST behaves differently

[tmachonkin]

I have been trying to set up and run ONION geometry-optimization calcs, and have been running into some problems.

If I run an open-shell system through xTB, it ignores my attempts to specify the number of unpaired electrons. I tried it using the input file ($spin 3) and I tried it in the command line (--uhf 3), and both times in the output file it gave

...................................................
: SETUP :
:.................................................:

      :  unpaired electrons                  0          :
      ...................................................

CoTpOPhCl2.zip

Now, if I run the same calc through CREST, it seems to work. At least, it converges to a different geometry, one that is much closer than to what I got when I ran the whole molecule at the gfn2 level. And the output file says:

---

Calculation info

Automatic ONIOM setup calculation level 1:
: GFN-FF calculation via GFNFF lib
: Molecular charge : 0
: UHF parameter : 3

Which seems to imply that it knows that the molecule has 3 unpaired electrons.

CoTpOPhCl2_7c.zip

So, I tried to set up and run some larger calcs on on open-shell systems through CREST. I had big problems getting any geometry optimization to converge. As a test case, I tried a larger diamagnetic system, and tried the exact same ONIOM calc using both xTB and CREST. The xTB calc converged just fine, and gave what looks like a perfectly reasonable structure. The CREST calc, however, failed after only 24 cycles, and the geometry it finished with looked very very strange.

CpRd_2.zip

I'm running xTB ver. 6.7.0, and CREST ver. 3.0.2

At this point, I'm a bit stuck. Running an ONIOM calc through xTB seems to work great if the system is diamagnetic, but I can't get it work on paramagnetic systems. Running an ONION calc through CREST seems to work on a small-molecule paramagnetic system, but the geometry-optimization has failed on several different metalloproteins.

Any help would be much appreciated.

Thanks.

[pprcht]

Hi,
unfortunately I can't tell you much about the xtb ONIOM implementation, theirs and the one in CREST were developed independently. However, I have an idea what's going on with the metalloprotein: The charge is botched during the automatic setup. You can see this in the calculation info printout

 ----------------
 Calculation info
 ----------------
> Automatic ONIOM setup calculation level 1:
 : GFN-FF calculation via GFNFF lib
 :   Molecular charge    : -14
 :   Solvation model     : alpb
 :   Solvent             : water
 :   ONIOM frag ("root") : 1
> Automatic ONIOM setup calculation level 2:
 : xTB calculation via tblite lib
 : GFN2-xTB level
 :   Molecular charge    : -1
 :   Fermi temperature   : 300.00000
 :   Accuracy            : 1.00000
 :   max SCC cycles      : 500
 :   ONIOM frag ("high") : 2
> Automatic ONIOM setup calculation level 3:
 : GFN-FF calculation via GFNFF lib
 :   Molecular charge    : -14
 :   Solvation model     : alpb
 :   Solvent             : water
 :    ONIOM frag ("low") : 2

where the "low-level" lower ONIOM GFN-FF replica of the Ga metal region gets a molecular charge of -14 (because the calculator was copied from the "high-level/root" region of the upper layer). This is of course way too much charge for the small fragment....

I don't remember why we are not simply copying the charge from the GFN2-xTB region (in this case), but there was a reason. Might be because there can be some confusion with the arbitrary number of ONIOM layers that could be defined. We instead opted for an additional definition of charges in the [lwoniom] block. In your case it would require something like

# Fragment specific charges
chrg.1 = -14
chrg.2 = -1

Which sets up the charges correctly:

 ----------------
 Calculation info
 ----------------
> Automatic ONIOM setup calculation level 1:
 : GFN-FF calculation via GFNFF lib
 :   Molecular charge    : -14
 :   ONIOM frag ("root") : 1
> Automatic ONIOM setup calculation level 2:
 : xTB calculation via tblite lib
 : GFN2-xTB level
 :   Molecular charge    : -1
 :   Fermi temperature   : 300.00000
 :   Accuracy            : 1.00000
 :   max SCC cycles      : 500
 :   ONIOM frag ("high") : 2
> Automatic ONIOM setup calculation level 3:
 : GFN-FF calculation via GFNFF lib
 :   Molecular charge    : -1
 :    ONIOM frag ("low") : 2

And converges in 540-ish opt. cycles. You can even remove the chrg= options from the [[calculation.level]] setups in that case if you want to have a cleaner input.

I hope that helps.

By the way, uhf won't do anything for GFN-FF, the method has no electronic description, even when it was printed as recognized keyword. It will only affect the xTB levels.

[tmachonkin]

Thanks. That seemed to work.
A question about how it handles solvation, though. In the .toml file, in the [[calculation.level.]] block, I had specified:

alpb='water'

In the output file, where it sets up the initial calculation, it says:

Calculation info

Automatic ONIOM setup calculation level 1:
: GFN-FF calculation via GFNFF lib
: Molecular charge : -14
: Solvation model : alpb
: Solvent : water
: ONIOM frag ("root") : 1
Automatic ONIOM setup calculation level 2:
: xTB calculation via tblite lib
: GFN2-xTB level
: Molecular charge : -1
: Fermi temperature : 300.00000
: Accuracy : 1.00000
: max SCC cycles : 500
: ONIOM frag ("high") : 2
Automatic ONIOM setup calculation level 3:
: GFN-FF calculation via GFNFF lib
: Molecular charge : -1
: Solvation model : alpb
: Solvent : water
: ONIOM frag ("low") : 2

Does that mean that when it runs the inner region at the low level that it is solvating it with water? I would assume that the point of the ONIOM calc is that for the inner region, when run at both the high and low levels, the "solvent" is the rest of the protein.

Thanks.