molden2qmc

General discussion of the Cambridge quantum Monte Carlo code CASINO; how to install and setup; how to use it; what it does; applications.
Vladimir_Konjkov
Posts: 165
Joined: Wed Apr 15, 2015 3:14 pm

molden2qmc

Post by Vladimir_Konjkov »

Dear all

The last few days I have been refactoring script molden2qmc.py, trying to understand the algorithm of conversion.
Normalization algorithm of MO-coefficients, and Cartesian->Spherical transformation remains unclear to me.

Because the script is written in python, I did not break the tradition, therefore I also used python.
The current version of the script, is in the repo https://github.com/Konjkov/molden2qmc

Could you explain how the conversion coefficient were calculated for the d and f orbitals.
I use them in the next functions :

Code: Select all

    def d_normalize(self, coefficient):
        """
        The following order of D functions is expected:
            5D: D 0, D+1, D-1, D+2, D-2
        """
        return (coefficient[0] / sqrt(3),
                coefficient[1] * 2.0,
                coefficient[2] * 2.0,
                coefficient[3],
                coefficient[4] * 2.0)

    def f_normalize(self, coefficient):
        """
        The following order of F functions is expected:
            7F: F 0, F+1, F-1, F+2, F-2, F+3, F-3
        """
        return (coefficient[0] * sqrt(8.0/15.0),
                coefficient[1] * 2.0 / sqrt(45),
                coefficient[2] * 2.0 / sqrt(45),
                coefficient[3] * sqrt(2) / 15.0,
                coefficient[4] * sqrt(2) / 15.0,
                coefficient[5] / sqrt(3) / 15.0,
                coefficient[6] / sqrt(3) / 15.0)
In the initial code, they look like this:

Code: Select all

d5_vectors.append(d5_vectors2[0] * one_over_rt_three)
d5_vectors.append(d5_vectors2[1] * 2.)
d5_vectors.append(d5_vectors2[2] * 2.).
d5_vectors.append(d5_vectors2[3] ).
d5_vectors.append(d5_vectors2[4] * 2.)

harm_f_1 = (8./15.0)**0.5
harm_f_2 = 2.0 * (1/(5**0.5))/3.0
harm_f_3 = rt2 / 15.0
harm_f_4 = (2.0**1.5) / 30.0
harm_f_5 = one_over_rt_three/15.0
So also could you tell how to calculate the similar coefficients for the normalization of g- orbitals.

My second question is about the transformation of Cartesian-> Spherical.

I found an article in wiki http://en.wikipedia.org/wiki/Table_of_s ... _harmonics
but I did not notice the correspondence in coefficients in the article and in the script.

It seemed to me that in the conversion of d-orbitals omitted (1/4)sqrt(15/pi) multiplier.

However, I failed to find a pattern in the following block of code:

Code: Select all

        f7_vectors.append(zero  *  0.331990278)
        f7_vectors.append(plus_1 * 0.059894236)
        f7_vectors.append(minus_1* 0.064894235)
        f7_vectors.append(plus_2*  0.059227155)
        f7_vectors.append(minus_2* 0.050227159)
        f7_vectors.append(plus_3 * 0.010915709)
        f7_vectors.append(minus_3* 0.010915709)
Is there any analytical expression for this coefficient ?

I would also like to know the considerations about convertation of cartesian representation of g-orbital to spherical.

Best regard, Vladimir.
In Soviet Russia Casino plays you.
Katharina Doblhoff
Posts: 84
Joined: Tue Jun 17, 2014 6:50 am

Re: molden2qmc

Post by Katharina Doblhoff »

]Hi!
Maybe this helps in understanding:
(I just copied in something I noted down some time ago. It is just my personal notes, so it may not be well explained. Feel free to ask.)

CASINO does not internally take care of the normalization of the Gaussian orbitals. A
normalization factor given by
eq1.PNG
eq1.PNG (4.64 KiB) Viewed 80059 times
thus has to be multiplied into the contraction coefficients of the primitive Gaussians. In
this equation, alpha denotes the exponent of the primitive Gaussian and l denotes the angular
momentum quantum number. (Note that part a cannot be multiplied into the orbital
coefficients, since these are contracted and thus cannot take the different exponents
of the primitive Gaussians into account any more, while part b could alternatively be
multiplied into the orbital coefficients.)
The expression given aboves, assumes the real solid harmonics R^m_l building the primitive
Gaussian type orbitals
eq3.PNG
eq3.PNG (2.78 KiB) Viewed 80059 times
to follow Racah’s normalization
eq2.PNG
eq2.PNG (5.66 KiB) Viewed 80059 times
Katharina Doblhoff
Posts: 84
Joined: Tue Jun 17, 2014 6:50 am

Re: molden2qmc

Post by Katharina Doblhoff »

CASINO uses a different normalization for the real solid harmonics. Normalization
factors given by
eq4.PNG
eq4.PNG (3.69 KiB) Viewed 80059 times
thus have to be multiplied in to the orbital coefficients. (Note that these factors cannot be
multiplied into the m-independent contraction coefficients, since they explicitly depend
on the rotational quantum number m.)
Furthermore, the d-orbitals (l = 2) are treated specially and have their own, nonstandard,
additional normalization factor in CASINO. The additional normalization
eq5.PNG
eq5.PNG (4.1 KiB) Viewed 80059 times
.
These additional factors also have to be multiplied into the orbital coefficients.


In the Gaussian to CASINO interface these factors are quite illogically multiplied into
the contraction and orbital coefficients. At the moment the converter does the following:
• l = 0 and l = 1 (s and p orbitals): a and b are multiplied into the contraction
coefficients, d = 1 in all cases and does not have to be considered.
• l = 2 (d orbitals): a is multiplied into the contraction coefficient, b, c and f are
multiplied together into the orbital coefficients, where a·b·f is given by 1 p3 , 2, 2, 1, 2
for m = 0, 1,−1, 2,−2.
• l = 3, 4 (f and g orbitals): a is multiplied into the contraction coefficients, b and c
are multiplied into the orbital coefficients.
Kevin_Gasperich
Posts: 7
Joined: Wed Mar 18, 2015 7:46 am

Re: molden2qmc

Post by Kevin_Gasperich »

Hello Vladimir,

I recently spent some time working on the Cartesian->spherical transformation to produce gwfn files from GAMESS output.

The d-shell transformation is still not completely clear to me, but the one that is used in the gamess2qmc perl script seems to be correct, so that is what is used in the internal GAMESS converter (initially implemented by Albert DeFusco; updated version to be included with the next GAMESS release).

I summarized some of my findings here (the explanation on that page could be made more thorough, so let me know if you have any questions).

Regarding the form of the spherical harmonics found on Wikipedia: those formulas are "correct" (i.e. the formulas you should be using) up to an l-dependent scale factor. I believe you should be using the regular solid harmonics instead.


Kevin
Vladimir_Konjkov
Posts: 165
Joined: Wed Apr 15, 2015 3:14 pm

Re: molden2qmc

Post by Vladimir_Konjkov »

Sorry....
I found some usefull definition in CASINO/examples/generic/gauss_dfg/README
d functions
-----------
l,m Orange book (and CASINO)
2,0 0.5 * (3z^2 - r^2)
2,1 3xz
2,-1 3yz
2,2 3(x^2-y^2)
2,-2 6xy

f functions
-----------
l,m Orange book (and CASINO)
3,0 z^3 - 3/2 * (x^2z + y^2z)
3,1 6 * xz^2 - 3/2 * (x^3 + xy^2)
3,-1 6 * yz^2 - 3/2 * (x^2y + y^3)
3,2 15 * (x^2z - y^2z)
3,-2 30 * xyz
3,3 15 * x^3 - 45 * xy^2
3,-3 45 * x^2y - 15 * y^3

g functions (not actually implemented in CRYSTAL)
-----------
l,m Orange book (and CASINO).......
4,4 105 * x^4 - 630 * x^2y^2 + 105 * y^4
4,3 105 * x^3z - 315 * xy^2z
4,2 -15/2 * x^4 + 45 * x^2z^2 + 15/2 * y^4 - 45 * y^2z^2
4,1 -15/2 * x^3z - 15/2 * xy^2z + 10 * xz^3
4,0 3/8 * x^4 + 3/4 * x^2y^2 - 3 * x^2z^2 + 3/8 * y^4 - 3 * y^2z^2 + z^4
4,-1 -15/2 * x^2yz - 15/2 * y^3z + 10 * yz^3.
4,-2 -15 * x^3y - 15 * xy^3 + 90 xyz^2
4,-3 315 * x^2yz - 105 * y^3z
4,-4 420 * x^3y - 420 * xy^3
I think this will answer some of my questions.
In Soviet Russia Casino plays you.
Katharina Doblhoff
Posts: 84
Joined: Tue Jun 17, 2014 6:50 am

Re: molden2qmc

Post by Katharina Doblhoff »

Is your problem the conversion from cartesian to spherical solid harmonics or is it the normalization i.e. is your molden file in a 5d7f format or in 6d10f?

Your first question was: where do the factors come from. If you read my answer above, you will see, that for the d-orbitals b·c·f is given by 1 /sqrt(3) , 2, 2, 1, 2 for m = 0, 1,−1, 2,−2, which is exactly the coefficients you see. Factor a will be in the contraction coefficients, I assume.
Last edited by Katharina Doblhoff on Thu Apr 30, 2015 6:39 am, edited 1 time in total.
Mike Towler
Posts: 239
Joined: Thu May 30, 2013 11:03 pm
Location: Florence
Contact:

Re: molden2qmc

Post by Mike Towler »

Hi all,
Furthermore, the d-orbitals (l = 2) are treated specially and have their own, nonstandard,
additional normalization factor in CASINO. The additional normalization..
I'm not sure I would describe this as an 'additional normalization factor' as Katharina does. This makes it sound like something mathematically profound is being done, whereas it's in reality just an efficiency issue.

Just to be clear and simple, let's say we have a single orbital expanded in a Gaussian basis consisting of a single uncontracted d shell with exponent a. For the harmonic case this shell defines 5 separate basis functions g_i * exp(-a|r|^2) where the g_i are polynomials in x,y,z. The orbital is then given by Phi(r)=(sum_i c_i g_i ) * exp(-a|r|^2) where the c_i are the suitably normalized orbital coefficients..

For a d shell, as Vladimir has correctly divined, the real solid harmonics g_i that we use are:

Code: Select all

g_1 = 0.5 * (3z^2 - r^2) 
g_2 = 3xz 
g_3 = 3yz 
g_4 = 3(x^2-y^2)
g_5 = 6xy 
Looking at e.g. g_5 = 6xy we see there is a constant factor 6 and since the expression for the orbital involves products c_i * g_i this factor can either (1) be premultiplied into the orbital coefficients c_i, or (2) multiplied into each new value of xy 'on the fly. Option (1) would be done by the Gaussian/GAMESS/Orca/MOLDEN (or whatever..) converter, Option (2) by CASINO itself. From the point of view of 'tidiness' (2) is preferable since we maintain a nice separation between the basis functions as defined in our standard reference, and the orbital coefficients as generated by the external program. However if we look at efficiency we see that (1) requires a single multiplication, and (2) requires several billion repetitions of the same multiplication over the course of a long simulation, so (1) is clearly preferable and that's the convention I chose when I originally developed this.

To understand where we are today, let's just look at the historical development:

I wrote the first version of the CASINO Gaussian orbital evaluator long ago in 1997. At that time the only Gaussian basis set code I had access to was the CRYSTAL program (which - despite the name - does atoms and molecules as well as periodic systems and so it was all I needed). However, CRYSTAL only implemented s,p,d and so that's what I did originally for CASINO. At that time, the convention above was clear: 'all constant numerical factors in the Gaussian basis functions must be premultiplied into the orbital coefficients by the converter program'. The s and p functions have no such numerical factors - only the d - and so everything was consistent at that point.

Several years later the CRYSTAL people implemented f functions, and consequently I did the same in CASINO. They had told me they were in the process of doing g functions, and so I coded that up in CASINO too. However the CRYSTAL g functions never appeared (and still haven't even today) and so the g part of the code was left unused and untested for many years. My original g implementation in CASINO in fact, somewhat miraculously, turned out to be correct - but as Katharina showed relatively recently the only converter ever to try to deal with g functions - Andrew Porter's gaussian2qmc - was not implemented correctly. We think it now is.

Now the main point is - when I did the f and g implementation I abandoned the convention about constant numerical factors being premultiplied into the orbital coefficients by the converter. In retrospect this was probably a mistake, and I can't quite remember why I did that. It was probably a mixture of (1) not caring about efficiency thinking no-one would ever use f and g, the implementation only being for 'completeness', (2) the numerical factors being much more 'messy' than in the d case, (3) the complication of the expressions particularly when you get to the analytic second derivatives meaning the premultiplication gives less of a saving - not that it buys you all that much even in the d case, but this is QMC and every little helps, (4) being very young, and (5) I may even have forgotten about this convention given the several years gap between s,p,d and f,g.

Now once we've adopted a convention for this it becomes very difficult to change it without breaking all pre-existing gwfn.data files, and that probably wouldn't be a very good idea. Once I became aware of the issue I took the view that in the documentation we should just write down the solid harmonics as defined in the original standard reference that I used (the 'Orange book' referred to in Vladimir's post - see the relevant pages at http://www.tcm.phy.cam.ac.uk/~mdt26/cry ... monics.pdf ) along with the formulae for the normalization coefficients and then accompany them with a statement of the convention: 'all constant factors in our real solid harmonics up to L=2 (effectively d only) must be premultiplied into the orbital coefficients, but not for L>2.' I appreciate this is a little inelegant, but it's clear and it shouldn't cause any problems for the writers of CASINO interfaces to external programs.

Never having had to deal with a Cartesian to spherical conversion myself (CRYSTAL uses 5d7f spherical solid harmonics natively..) I allowed converter authors to have their own adventures doing this for their own codes, and hoped they would be able to work it out...

I've been threatening for several weeks to synthesize all relevant information on this in a new appendix to the CASINO manual. All the stuff currently in the gauss_dfg/REAME file will be put in there, and hopefully some example Cartesian --> spherical stuff of the kind that Kevin and our Pittsburgh friends have been working on. However, as so often this year, life keeps intruding. Hopefully soon.

Mike
Vladimir_Konjkov
Posts: 165
Joined: Wed Apr 15, 2015 3:14 pm

Re: molden2qmc

Post by Vladimir_Konjkov »

Katharina Doblhoff wrote:Is your problem the conversion from cartesian to spherical solid harmonics or is it the normalization i.e. is your molden file in a 5d7f format or in 6d10f? Could you please specify your question better! Your initial question was: where do the factors come from. If you read my answer above, you will see, that for the d-orbitals b·c·f is given by 1 /sqrt(3) , 2, 2, 1, 2 for m = 0, 1,−1, 2,−2, which is exactly the coefficients you see. Factor a will be in the contraction coefficients, I assume.
Hello, Katharina.

I used to write code in an object oriented style and I used to write code based on detailed specification.
and because format specification of molden and gwfn.data files exist and is quite detailed, I wrote a Molden() base class in molden2qmc.py,
which can read a molden file and store ALL available data in structures atom_list and mo_matrix whose format I have documented in detail.

Molden() class also has a method which can write atom_list, mo_matrix data to gwfn.data file.

It must be noted that the Molden class itself reads and writes data in atom_list, mo_matrix structures without any normalization and conversion
that makes impossible its direct use. Molden class only "helps" to work with data in a more convenient format.

So another thing I have to do is implement conversion and normalization procedures.
Since these procedures may be specific to different programs (ORCA, CFour, PSI4, etc), I implement subclasses for each program separately.
in which I must implement next three methods:

1. data conversion form cartesian to spherical - applying to [MO] section of molden data
2. m-independent normalisation - applying to [GTO] section of molden data
3. m-dependent normalisation - applying to [MO] section of molden data

Since these procedures are not well documented, I asked the relevant questions in this topic.

To answer your question, I can say that cartesian -> spherical conversion procedures will be applied or not,
depending on the particular Molden file format (5d7f or 6d10f or other possible)
and appropriate normalisation procedures will be applied depending on which program created the Molden file.

So "my problem" is that I must implement 'em all.

Vladimir.
Last edited by Vladimir_Konjkov on Thu Apr 30, 2015 4:11 pm, edited 1 time in total.
In Soviet Russia Casino plays you.
Katharina Doblhoff
Posts: 84
Joined: Tue Jun 17, 2014 6:50 am

Re: molden2qmc

Post by Katharina Doblhoff »

Hi Vladimir!

Believe me, I know that this is a hell of a work, and I know that normalizations (also in the molden format) are not well documented. Otherwise I would not even try to help you. I myself never worked with a cartesian file format. I think it is dangerous, since one would have to test whether the additional functions were eliminated during the MO optimization (which should actually be straight forward: if there are less MOS than cartesian contracted GTOs then everything is fine, if there are as many MOs as cartesian contracted GTOs then there is no way this can be transferred to casino, unless you add "effective" basis functions).

But lets begin step by step: The normalization I gave in my first post should be OK for spherical GTOs. Why don't you go ahead and test it. Then you can take the second step and try the cartesians. Have you looked at what Kevin wrote? I do not think (although I am not sure since I did not test) that you can directly use the transformation matrices you wrote down from the casino documentation. I believe they have a different normalization of the cartesians, but I may be mistaken. Did you have a look at Kevin's documentation?
Vladimir_Konjkov
Posts: 165
Joined: Wed Apr 15, 2015 3:14 pm

Re: molden2qmc

Post by Vladimir_Konjkov »

Hello, Katharina.

I can not quite understand what justified the use of ratios (1/sqrt(3), 2, 2, 1, 2) for the d-orbitals?
I think they do not needed.

I calculated the cubane molecule in basis dev2-TZVP (up to f-orbitals) in ORCA, set print options and output the overlap matrix and the MO coefficients when converged.

Code: Select all

%output
 Print[P_Mos] 1
 Print[P_Overlap] 1
end
Then with a script I processed the data and carried out the following checks:

1. Compute Density Matrix: D = Cocc*Cocc_transpose from eigenvectors C contains the MO coefficients.
2. Check Density Matrix for symmetry, trace, and idempotency conditions.

Density Matrix symmetrical, idempotency conditions is not satisfied, trace condition 2*Tr(DS) = 55.9999935755 agreed with the exact number of electrons in the molecule (56).

MO coefficient is okay, any ideas? Why idempotency conditions is not satisfied?

Vladimir.
In Soviet Russia Casino plays you.
Post Reply