backflow to the g-orbitals

[Vladimir_Konjkov]

Hello, CASINO developers.

It is a well-known that as a trial wave function for the CASINO you can use HF or DFT wfn, the latter in many cases give a better WFN approximation.
I also tried to use orbital-optimized MP2 trial wave function but it's almost never better than DFT, at least it does not give significant improvements.
After that it's time to use multi-determinant trial wave function, this approach can improve DMC energy further, but the number of determinants that must be taken to obtain the exact energy is frightening.
I suppose that after taking into account the necessary number of determinants (in some cases, only one), when we get the correct topology of the nodal surface, we can not add new determinants, but use the backflow to improve the quality of the nodal surface. So I plan to use the backflow in my calculations next year.
In most cases, the generation of a multi-determinant wave function is performed using CI/CC methods, I can recall only a few simple systems such as Be, B, C atoms for which multi-determinant wave function can be generated from CASSCF or written by hands.
CI/CC methods requires the use of a high angular momentum basis (like cc-pVQZ, cc-pV5Z, cc-pV6Z) because of nature of methods. So applying the backflow after this also requires support as high orbital moment as possible.
Unfortunately, I noticed that backflow supports only up to f-orbitals.

I would like to know what plans for the new functionality of backflow next year.

I found use_gbackflow key word

CASINO HELP SYSTEM
==================

Keyword : use_gbackflow
Title : DEV: use gbackflow
Type : Boolean
Level : Expert

DESCRIPTION
-----------
RESERVED KEYWORD FOR DEVELOPMENT, NO EFFECT.

Is still no effect of this?

How difficult is to support of g-orbitals in backflow, for me if I would start programming in FORTRAN?
How slower will the backflow calculation be with g-orbitals?
In the preliminary calculations, I found that it is sufficient to use only ETA and MU terms in backflow to improve MD WFN energy and I also making full electron calculations. How much does this simplify or complicate the task?

Best Vladimir.

[Pablo_Lopez_Rios]

Hi Vladimir,

Implementing backflow for g orbitals requires coding up the evaluation of the matrix of second derivatives of g Gaussians, i.e., d2/dx2, d2/dy2, d2/dz2, d2/dxdy, d2/dxdz, and d2/dydz of phi(r), in gauss_mol_bf.f90. This is regardless of which backflow terms you wish to use. Module gauss_mol.f90 can handle the evaluation of the Laplacian of g Gaussians, i.e., d2/dx2 + d2/dy2 + d2/dz2, which is the trace of the second derivative matrix; this might be helpful in the implementation.

BTW, the (missing) functionality associated with the "use_gbackflow" keyword has nothing to do with g orbitals - it is intended to be analogous to "use_gjastrow", triggering the use of a different implementation of the evaluation of backflow transformations.

Best,
Pablo

[Vladimir_Konjkov]

Hello Pablo

I've read the code in gauss_mol_bf.f90 and found that.
pgra array serves to store gradients in a form (d orb(m=0)/dx, d orb(m=0)/dy, d orb(m=0)/dz, ....) and so on where m in (0, 1, -1, 2, -2, ...)
psderiv array serves to store second derivatives in a form (d2 orb(m=0)/dx2, d2 orb(m=0)/dy2, d2 orb(m=0)/dz2, d2 orb(m=0)/dxdy, d2 orb(m=0)/dxdz, d2 orb(m=0)/dydz) and so on where m in (0, 1, -1, 2, -2, ...)

I do not fully understand the normalization coefficients, so I wrote a python script which using sympy http://www.sympy.org makes symbolic differentiation of gaussian orbitals.

And I try to compare code from CASINO source and generated by the script.

I'm sorry but today it's too late to continue working.
Best Vladimir

[Vladimir_Konjkov]

Hello Pablo

I've found that G-orbitals should be in the following format in order to obtain pgra from the source code

Code: Select all

    (35.0 * z*z*z*z - 30.0 * z*z*r2 + 3.0 * r4) / 8.0 * exp(-alpha*r2),     # (4,  0)
    2.5 * (7.0 * z*z*z*x - 3.0 * x*z*r2) * exp(-alpha*r2),                  # (4,  1)
    2.5 * (7.0 * z*z*z*y - 3.0 * y*z*r2) * exp(-alpha*r2),                  # (4, -1)
    7.5 * (7.0 * (x*x*z*z - y*y*z*z) - (x*x*r2 - y*y*r2)) * exp(-alpha*r2), # (4,  2)
    15 * (7.0 * z*z*x*y - x*y*r2) * exp(-alpha*r2),                         # (4, -2)
    105 * (x*x*x*z - 3.0 * y*y*x*z) * exp(-alpha*r2),                       # (4,  3)
    105 * (3.0 * x*x*y*z - y*y*y*z) * exp(-alpha*r2),                       # (4, -3)
    105 * (x*x*x*x - 6.0 * x*x*y*y + y*y*y*y) * exp(-alpha*r2),             # (4,  4)
    420 * (x*x*x*y - y*y*y*x) * exp(-alpha*r2),                             # (4, -4)

I also found that the pgra(27) has a wrong sign, should be pgra(27)=840.d0*alpha*xyz*(yy-xx)

Code: Select all

    
    case(6) ! g (harmonic)
    ....
    pgra(27)=840.d0*alpha*xyz*(xx-yy)

Best, Vladimir

[Mike Towler]

Hi Vladimir,

There's more information about normalization of d/f/g and related issues (including various stupid conventions we use..) in my notes at CASINO/examples/generic/gauss_dfg/README

Cheers,
Mike

[Vladimir_Konjkov]

Hi Vladimir,

There's more information about normalization of d/f/g and related issues (including various stupid conventions we use..) in my notes at CASINO/examples/generic/gauss_dfg/README

Cheers,
Mike

Thank you Mike for all your advices. I checked my formulas are identical to there in the CASINO/examples/generic/gauss_dfg/README (also for d,f,g-orbitals)

It seems to me that using symbolic mathematics from SymPy one can easily generate expressions for the second derivatives of G-orbitals.

I got the following results (script in the attachments):

Code: Select all

sderiv(1)=alpha**2*x**2*(17.5*z**4 - 15.0*z**2*(x**2 + y**2 + z**2) + 1.5*(x**2 + y**2 + z**2)**2) - 4*alpha*x**2*(1.5*x**2 + 1.5*y**2 - 6.0*z**2) - 2*alpha*(4.375*z**4 - 3.75*z**2*(x**2 + y**2 + z**2) + 0.375*(x**2 + y**2 + z**2)**2) + 4.5*x**2 + 1.5*y**2 - 6.0*z**2
sderiv(2)=alpha**2*y**2*(17.5*z**4 - 15.0*z**2*(x**2 + y**2 + z**2) + 1.5*(x**2 + y**2 + z**2)**2) - 4*alpha*y**2*(1.5*x**2 + 1.5*y**2 - 6.0*z**2) - 2*alpha*(4.375*z**4 - 3.75*z**2*(x**2 + y**2 + z**2) + 0.375*(x**2 + y**2 + z**2)**2) + 1.5*x**2 + 4.5*y**2 - 6.0*z**2
sderiv(3)=-15.0*alpha**2*x**2*z**4 - 15.0*alpha**2*y**2*z**4 + 2.5*alpha**2*z**6 + 1.5*alpha**2*z**2*(x**2 + y**2 + z**2)**2 + 31.5*alpha*x**2*z**2 + 31.5*alpha*y**2*z**2 - 17.25*alpha*z**4 - 0.75*alpha*(x**2 + y**2 + z**2)**2 - 6.0*x**2 - 6.0*y**2 + 12.0*z**2
sderiv(4)=x*y*(alpha**2*(17.5*z**4 - 15.0*z**2*(x**2 + y**2 + z**2) + 1.5*(x**2 + y**2 + z**2)**2) - 4*alpha*(1.5*x**2 + 1.5*y**2 - 6.0*z**2) + 3.0)
sderiv(5)=x*z*(alpha**2*(17.5*z**4 - 15.0*z**2*(x**2 + y**2 + z**2) + 1.5*(x**2 + y**2 + z**2)**2) - 2*alpha*(1.5*x**2 + 1.5*y**2 - 6.0*z**2) + 2*alpha*(6.0*x**2 + 6.0*y**2 - 4.0*z**2) - 12.0)
sderiv(6)=y*z*(alpha**2*(17.5*z**4 - 15.0*z**2*(x**2 + y**2 + z**2) + 1.5*(x**2 + y**2 + z**2)**2) - 2*alpha*(1.5*x**2 + 1.5*y**2 - 6.0*z**2) + 2*alpha*(6.0*x**2 + 6.0*y**2 - 4.0*z**2) - 12.0)
sderiv(7)=x*z*(-30.0*alpha**2*x**4 - 30.0*alpha**2*x**2*y**2 + 40.0*alpha**2*x**2*z**2 + 105.0*alpha*x**2 + 45.0*alpha*y**2 - 60.0*alpha*z**2 - 45.0)
sderiv(8)=-x*z*(alpha**2*y**2*(30.0*x**2 + 30.0*y**2 - 40.0*z**2) - 60.0*alpha*y**2 - 2*alpha*(7.5*x**2 + 7.5*y**2 - 10.0*z**2) + 15.0)
sderiv(9)=x*z*(-30.0*alpha**2*x**2*z**2 - 30.0*alpha**2*y**2*z**2 + 40.0*alpha**2*z**4 + 45.0*alpha*x**2 + 45.0*alpha*y**2 - 140.0*alpha*z**2 + 60.0)
sderiv(10)=-y*z*(alpha**2*x**2*(30.0*x**2 + 30.0*y**2 - 40.0*z**2) - 30.0*alpha*x**2 - 2*alpha*(22.5*x**2 + 7.5*y**2 - 10.0*z**2) + 15.0)
sderiv(11)=alpha**2*x**2*z**2*(-30.0*x**2 - 30.0*y**2 + 40.0*z**2) + 2*alpha*x**2*(7.5*x**2 + 7.5*y**2 - 30.0*z**2) + 2*alpha*z**2*(22.5*x**2 + 7.5*y**2 - 10.0*z**2) - 22.5*x**2 - 7.5*y**2 + 30.0*z**2
sderiv(12)=-x*y*(alpha**2*z**2*(30.0*x**2 + 30.0*y**2 - 40.0*z**2) - 30.0*alpha*z**2 - 2*alpha*(7.5*x**2 + 7.5*y**2 - 30.0*z**2) + 15.0)
sderiv(13)=-y*z*(alpha**2*x**2*(30.0*x**2 + 30.0*y**2 - 40.0*z**2) - 60.0*alpha*x**2 - 2*alpha*(7.5*x**2 + 7.5*y**2 - 10.0*z**2) + 15.0)
sderiv(14)=y*z*(-30.0*alpha**2*x**2*y**2 - 30.0*alpha**2*y**4 + 40.0*alpha**2*y**2*z**2 + 45.0*alpha*x**2 + 105.0*alpha*y**2 - 60.0*alpha*z**2 - 45.0)
sderiv(15)=y*z*(-30.0*alpha**2*x**2*z**2 - 30.0*alpha**2*y**2*z**2 + 40.0*alpha**2*z**4 + 45.0*alpha*x**2 + 45.0*alpha*y**2 - 140.0*alpha*z**2 + 60.0)
sderiv(16)=-x*z*(alpha**2*y**2*(30.0*x**2 + 30.0*y**2 - 40.0*z**2) - 30.0*alpha*y**2 - 2*alpha*(7.5*x**2 + 22.5*y**2 - 10.0*z**2) + 15.0)
sderiv(17)=-x*y*(alpha**2*z**2*(30.0*x**2 + 30.0*y**2 - 40.0*z**2) - 30.0*alpha*z**2 - 2*alpha*(7.5*x**2 + 7.5*y**2 - 30.0*z**2) + 15.0)
sderiv(18)=alpha**2*y**2*z**2*(-30.0*x**2 - 30.0*y**2 + 40.0*z**2) + 2*alpha*y**2*(7.5*x**2 + 7.5*y**2 - 30.0*z**2) + 2*alpha*z**2*(7.5*x**2 + 22.5*y**2 - 10.0*z**2) - 7.5*x**2 - 22.5*y**2 + 30.0*z**2
sderiv(19)=-30.0*alpha**2*x**6 + 180.0*alpha**2*x**4*z**2 + 30.0*alpha**2*x**2*y**4 - 180.0*alpha**2*x**2*y**2*z**2 + 135.0*alpha*x**4 - 450.0*alpha*x**2*z**2 - 15.0*alpha*y**4 + 90.0*alpha*y**2*z**2 - 90.0*x**2 + 90.0*z**2
sderiv(20)=-30.0*alpha**2*x**4*y**2 + 180.0*alpha**2*x**2*y**2*z**2 + 30.0*alpha**2*y**6 - 180.0*alpha**2*y**4*z**2 + 15.0*alpha*x**4 - 90.0*alpha*x**2*z**2 - 135.0*alpha*y**4 + 450.0*alpha*y**2*z**2 + 90.0*y**2 - 90.0*z**2
sderiv(21)=-30.0*alpha**2*x**4*z**2 + 180.0*alpha**2*x**2*z**4 + 30.0*alpha**2*y**4*z**2 - 180.0*alpha**2*y**2*z**4 + 15.0*alpha*x**4 - 450.0*alpha*x**2*z**2 - 15.0*alpha*y**4 + 450.0*alpha*y**2*z**2 + 90.0*x**2 - 90.0*y**2
sderiv(22)=alpha*x*y*(-30.0*alpha*x**4 + 180.0*alpha*x**2*z**2 + 30.0*alpha*y**4 - 180.0*alpha*y**2*z**2 + 60.0*x**2 - 60.0*y**2)
sderiv(23)=x*z*(-30.0*alpha**2*x**4 + 180.0*alpha**2*x**2*z**2 + 30.0*alpha**2*y**4 - 180.0*alpha**2*y**2*z**2 - 120.0*alpha*x**2 + 180.0*alpha*y**2 - 180.0*alpha*z**2 + 180.0)
sderiv(24)=y*z*(-30.0*alpha**2*x**4 + 180.0*alpha**2*x**2*z**2 + 30.0*alpha**2*y**4 - 180.0*alpha**2*y**2*z**2 - 180.0*alpha*x**2 + 120.0*alpha*y**2 + 180.0*alpha*z**2 - 180.0)
sderiv(25)=x*y*(-60.0*alpha**2*x**4 - 60.0*alpha**2*x**2*y**2 + 360.0*alpha**2*x**2*z**2 + 210.0*alpha*x**2 + 90.0*alpha*y**2 - 540.0*alpha*z**2 - 90.0)
sderiv(26)=x*y*(-60.0*alpha**2*x**2*y**2 - 60.0*alpha**2*y**4 + 360.0*alpha**2*y**2*z**2 + 90.0*alpha*x**2 + 210.0*alpha*y**2 - 540.0*alpha*z**2 - 90.0)
sderiv(27)=-x*y*(alpha**2*z**2*(60*x**2 + 60*y**2 - 360.0*z**2) + 720.0*alpha*z**2 - 2*alpha*(15*x**2 + 15*y**2 - 90.0*z**2) - 180.0)
sderiv(28)=alpha**2*x**2*y**2*(-60*x**2 - 60*y**2 + 360.0*z**2) + 2*alpha*x**2*(15*x**2 + 45*y**2 - 90.0*z**2) + 2*alpha*y**2*(45*x**2 + 15*y**2 - 90.0*z**2) - 45*x**2 - 45*y**2 + 90.0*z**2
sderiv(29)=-y*z*(alpha**2*x**2*(60*x**2 + 60*y**2 - 360.0*z**2) + 360.0*alpha*x**2 - 2*alpha*(45*x**2 + 15*y**2 - 90.0*z**2) - 180.0)
sderiv(30)=-x*z*(alpha**2*y**2*(60*x**2 + 60*y**2 - 360.0*z**2) + 360.0*alpha*y**2 - 2*alpha*(15*x**2 + 45*y**2 - 90.0*z**2) - 180.0)
sderiv(31)=x*z*(420.0*alpha**2*x**4 - 1260.0*alpha**2*x**2*y**2 - 1470.0*alpha*x**2 + 1890.0*alpha*y**2 + 630.0)
sderiv(32)=x*z*(420.0*alpha**2*x**2*y**2 - 1260.0*alpha**2*y**4 - 210.0*alpha*x**2 + 3150.0*alpha*y**2 - 630.0)
sderiv(33)=2*alpha*x*z*(105*x**2 - 315.0*y**2)*(2*alpha*z**2 - 3)
sderiv(34)=y*z*(420.0*alpha**2*x**4 - 1260.0*alpha**2*x**2*y**2 + 630.0*alpha*x**2 + 630.0*alpha*y**2 - 630.0)
sderiv(35)=alpha**2*x**2*z**2*(420*x**2 - 1260.0*y**2) - 2*alpha*x**2*(105*x**2 - 315.0*y**2) - 2*alpha*z**2*(315*x**2 - 315.0*y**2) + 315*x**2 - 315.0*y**2
sderiv(36)=x*y*(alpha**2*z**2*(420*x**2 - 1260.0*y**2) + 1260.0*alpha*z**2 - 2*alpha*(105*x**2 - 315.0*y**2) - 630.0)
sderiv(37)=y*z*(1260.0*alpha**2*x**4 - 420.0*alpha**2*x**2*y**2 - 3150.0*alpha*x**2 + 210.0*alpha*y**2 + 630.0)
sderiv(38)=y*z*(1260.0*alpha**2*x**2*y**2 - 420.0*alpha**2*y**4 - 1890.0*alpha*x**2 + 1470.0*alpha*y**2 - 630.0)
sderiv(39)=2*alpha*y*z*(315.0*x**2 - 105*y**2)*(2*alpha*z**2 - 3)
sderiv(40)=x*z*(1260.0*alpha**2*x**2*y**2 - 420.0*alpha**2*y**4 - 630.0*alpha*x**2 - 630.0*alpha*y**2 + 630.0)
sderiv(41)=x*y*(alpha**2*z**2*(1260.0*x**2 - 420*y**2) - 1260.0*alpha*z**2 - 2*alpha*(315.0*x**2 - 105*y**2) + 630.0)
sderiv(42)=alpha**2*y**2*z**2*(1260.0*x**2 - 420*y**2) - 2*alpha*y**2*(315.0*x**2 - 105*y**2) - 2*alpha*z**2*(315.0*x**2 - 315*y**2) + 315.0*x**2 - 315*y**2
sderiv(43)=420.0*alpha**2*x**6 - 2520.0*alpha**2*x**4*y**2 + 420.0*alpha**2*x**2*y**4 - 1890.0*alpha*x**4 + 6300.0*alpha*x**2*y**2 - 210.0*alpha*y**4 + 1260.0*x**2 - 1260.0*y**2
sderiv(44)=420.0*alpha**2*x**4*y**2 - 2520.0*alpha**2*x**2*y**4 + 420.0*alpha**2*y**6 - 210.0*alpha*x**4 + 6300.0*alpha*x**2*y**2 - 1890.0*alpha*y**4 - 1260.0*x**2 + 1260.0*y**2
sderiv(45)=2*alpha*(2*alpha*z**2 - 1)*(105*x**4 - 630.0*x**2*y**2 + 105*y**4)
sderiv(46)=x*y*(420.0*alpha**2*x**4 - 2520.0*alpha**2*x**2*y**2 + 420.0*alpha**2*y**4 + 1680.0*alpha*x**2 + 1680.0*alpha*y**2 - 2520.0)
sderiv(47)=2*alpha*x*z*(2*alpha*(105*x**4 - 630.0*x**2*y**2 + 105*y**4) - 420*x**2 + 1260.0*y**2)
sderiv(48)=2*alpha*y*z*(2*alpha*(105*x**4 - 630.0*x**2*y**2 + 105*y**4) + 1260.0*x**2 - 420*y**2)
sderiv(49)=840*x*y*(2*alpha**2*x**4 - 2*alpha**2*x**2*y**2 - 7*alpha*x**2 + 3*alpha*y**2 + 3)
sderiv(50)=840*x*y*(2*alpha**2*x**2*y**2 - 2*alpha**2*y**4 - 3*alpha*x**2 + 7*alpha*y**2 - 3)
sderiv(51)=840*alpha*x*y*(x**2 - y**2)*(2*alpha*z**2 - 1)
sderiv(52)=1680*alpha**2*x**4*y**2 - 1680*alpha**2*x**2*y**4 - 840*alpha*x**4 + 840*alpha*y**4 + 1260*x**2 - 1260*y**2
sderiv(53)=840*alpha*y*z*(2*alpha*x**2*(x**2 - y**2) - 3*x**2 + y**2)
sderiv(54)=840*alpha*x*z*(2*alpha*y**2*(x**2 - y**2) - x**2 + 3*y**2)

I think it also would be useful to refactor the existing code.

Best, Vladimir

[Mike Towler]

Hi Vladimir,

(1) Sure - it's not that difficult to implement - just a bit boring. You just need to add the relevant g mixed derivatives to the routines gm2_bf and gm3_bf in gauss_mol_bf.f90. There's no need to do this in the periodic Gaussian routines like gauss_per_bf.f90 because the only serious periodic Gaussian DFT code (CRYSTAL) doesn't support g functions. If you want to have a go at this then analyze very carefully how the f function second derivatives are handled in the same routines then just redo that in an analogous way for the g.

The mixed second derivatives weren't coded up at the time backflow was implemented because I didn't have immediate access to a DFT code that did g functions, and there was some question over whether the g implementation in CASINO worked at all in the normal (non-backflow) case - the relevant converter utilities either explicitly didn't support g functions or were written in languages I didn't understand or the person responsible had left. I always meant to come back to it but never did.

(2) As for this:

I also found that the pgra(27) has a wrong sign, should be pgra(27)=840.d0*alpha*xyz*(yy-xx)

Can you triple-check that before I have a look? Everything was supposed to have been numerically tested with the g-ane examples in the examples directory a few years back..

M.

[Vladimir_Konjkov]

Hi Vladimir,

I also found that the pgra(27) has a wrong sign, should be pgra(27)=840.d0*alpha*xyz*(yy-xx)

Can you triple-check that before I have a look?
M.

OK let's look at the following piece if code:

Code: Select all

     
     t16=xx*(3.d0+twoalphay2)
     t17=2.d0*alphax4+yy
     pgra(25)=420.d0*y*(t16-t17)
     t14=2.d0*alphay4-3.d0*yy
     t15=xx*(twoalphay2-1.d0)
     pgra(26)=420.d0*x*(t14-t15)
     pgra(27)=840.d0*alpha*xyz*(xx-yy)

where pgra(25) is d(420 * x^3y - 420 * xy^3)*exp(-alpha*r2)/dx Wolfram algebra solution
where pgra(26) is d(420 * x^3y - 420 * xy^3)*exp(-alpha*r2)/dy Wolfram algebra solution
where pgra(27) is d(420 * x^3y - 420 * xy^3)*exp(-alpha*r2)/dz Wolfram algebra solution

simplify expressions:

Code: Select all

     
     pgra(25)=420.d0*y*(xx*(3.d0+2.d0*alpha*yy)-2.d0*alpha*x4-yy)
     pgra(26)=420.d0*x*(2.d0*alpha*y4-3.d0*yy-xx*(2.d0*alpha*y2-1.d0))
     pgra(27)=840.d0*alpha*xyz*(xx-yy)

pgra(25) and pgra(26) has a right sign, but pgra(27) doesn't.

Best, Vladimir.

[Mike Towler]

Hi Vladimir,

OK - I've now had a look at this. You're correct that there is a problem in the analytic first derivatives of g functions.

One can see that - for example - the 'g-ane' example in examples/generic/gauss_dfg/gaussian/g-ane fails the numerical check of the analytic derivatives enabled by setting checkwfn : T in input.

The trouble is that, if I correct the sign of pgra(27) in the way that you suggest, it still fails this test. So there must be something else wrong.

Something seems to have been overlooked here. The tests a few years ago that I referred to seem to have pre-assumed that the analytic derivatives were correct (there is a statement to that effect in DIARY 2.13.350 and it's hardly very difficult to check). The recent tests of g-ane etc. were supposed to check whether e.g. GAUSSIAN and CASINO gave the same answer for the Hartree-Fock energy, i.e. verifying that all the wave function/basis set/geometry information was correctly passed between the two codes by the converter. Funnily enough they did give the same answer within error bars according to examples/generic/gauss_dfg/RESULTS - perhaps the error bars weren't small enough.

I can't immediately see any other errors in the analytic derivatives.. Can you?

M.

[Vladimir_Konjkov]

Hi Vladimir,

OK - I've now had a look at this. You're correct that there is a problem in the analytic first derivatives of g functions.

One can see that - for example - the 'g-ane' example in examples/generic/gauss_dfg/gaussian/g-ane fails the numerical check of the analytic derivatives enabled by setting checkwfn : T in input.

The trouble is that, if I correct the sign of pgra(27) in the way that you suggest, it still fails this test. So there must be something else wrong.

Something seems to have been overlooked here. The tests a few years ago that I referred to seem to have pre-assumed that the analytic derivatives were correct (there is a statement to that effect in DIARY 2.13.350 and it's hardly very difficult to check). The recent tests of g-ane etc. were supposed to check whether e.g. GAUSSIAN and CASINO gave the same answer for the Hartree-Fock energy, i.e. verifying that all the wave function/basis set/geometry information was correctly passed between the two codes by the converter. Funnily enough they did give the same answer within error bars according to examples/generic/gauss_dfg/RESULTS - perhaps the error bars weren't small enough.

I can't immediately see any other errors in the analytic derivatives.. Can you?

M.

Hello Mike.

Because I carry out all my calculations in the cc-pVQZ basis I'm interested in getting the correct results. I checked selectively only one-third of all G-terms (one per every m), tomorrow I will check them all.
could you tell me what the analytical gradient affects? Kinetic energy, VMC/DMC variance or total energy, etc?

In the gm2_bf SUBROUTINE evaluates the elements of the prim array in CASE of g-orbitals, it seems to me that in this function such calculations are not required.

Best, Vladimir.