Charge density for finite molecules

General discussion of the Cambridge quantum Monte Carlo code CASINO; how to install and setup; how to use it; what it does; applications.
Post Reply
Yasmine_AlHamdani
Posts: 3
Joined: Tue Feb 04, 2014 4:29 pm

Charge density for finite molecules

Post by Yasmine_AlHamdani »

Hi,
I notice in the manual it is said that for the charge density expectation value, "This can only be done for periodic systems and finite single atoms at the moment; the finite molecule case has not been implemented."
Any ideas whether this will be implemented in the near future?
Thanks.
Mike Towler
Posts: 239
Joined: Thu May 30, 2013 11:03 pm
Location: Florence
Contact:

Re: Charge density for finite molecules

Post by Mike Towler »

Hi Yasmine!

Sadly, there are no plans to do this anytime soon that I am aware of. The main problem is, well, most of our developers are (understandably in this modern life) too busy to actually develop anything. I'm working in the background to improve this situation, but it may be some time before things change.

If you look at this page:

http://vallico.net/casinoqmc/things-to-do/

we note three things:

(1) Your issue is at number two in the list.

(2) Nobody has their name next to it in green ink.

(3) You're listed as a developer! ;)

I'm working pretty flat out at the moment to produce a major new version of the code (2.14). I'm currently expecting this will be completed in a few months time, then I have a whole bunch of other things to do. I could conceivably take a look at it in Vallico over the Easter holidays if you really want me to.. Let me know.

M.
Yasmine_AlHamdani
Posts: 3
Joined: Tue Feb 04, 2014 4:29 pm

Re: Charge density for finite molecules

Post by Yasmine_AlHamdani »

Hi Mike,

Forgive my naivety, but it feels like there's a catch to doing this that I'm missing...? If the charge density can be accumulated for a single atom in a finite system, I would have thought a molecule is not much more beyond that. If there is a known complication it would be great to be made aware of this...

Yasmine

PS. I'm not saying I'll do it, so no need to put my name in green ink just yet!
Mike Towler
Posts: 239
Joined: Thu May 30, 2013 11:03 pm
Location: Florence
Contact:

Re: Charge density for finite molecules

Post by Mike Towler »

Hi Yasmine,
If the charge density can be accumulated for a single atom in a finite system, I would have thought a molecule is not much more beyond that. If there is a known complication it would be great to be made aware of this...
OK, you asked.. :D

What's a charge density? Well, forget the charge bit, as the charge is 1 in atomic units. So it's a probability density associated with some point r. To make a probability density into a probability, we need to multiply it by a volume, normally the volume dr of an infinitesimal little box centred on r. So if you define the probability density as the absolute square of the many-body wave function (with the coordinates of electrons 2 to N integrated over all space), multiplied by N, that will give the probability for finding an electron (any electron) in that little box.

How to calculate this in CASINO? Well - as you saw - the manual says it uses two different methods, a reciprocal space method for periodic systems, and a real space method for finite systems (in practice restricted to atoms).

Let's get periodic systems out of the way first. CASINO just represents the density as a Fourier series i.e. as a linear combination of plane-waves
exp(iG.r) of different wavelength. The coefficient of each term is essentially just the overlap integral of the given plane wave and the function you're trying to expand (the density). In a QMC calculation the density is represented by the distribution in space and time of the electrons (i.e on the VMC random walk, or the set of configurations that are being propagated in imaginary time in a DMC calculation). So, if you think about it, that means that every time you move an electron to point r, you work out e(iG.r). You sum that over the all the different r after each move, then divide by the number of moves at the end to give the average value. That's just the overlap integral, i.e. the Fourier coefficient, and there you have your answer.

You're interested in finite systems though, where we use the real space method. Despite the massive overcomplication of the relevant section in the manual, this just means:

"Divide space into little boxes. As the electrons move around, count how many times an electron lands in each little box. Draw a histogram of the results at the end."

The only possible complications are that I'm assuming for the moment all boxes have the same volume, but they might not have, and I might need to include the *weight* of a given DMC configuration, rather than just counting 1 per electron visit."

OK, let's say we have a molecule like benzene or something. I start by assuming the density will be zero outside a box 10 Angstrom cubed (or whatever) with the molecule centred in the middle of it. I think 5000 points per side will give me an adequate resolution for the density.. This gives me 125 billion little boxes. For each box, I have to store an integer (the number of electron visits to that box) which takes 4 bytes per box i.e. 500 billion bytes or 500 Gb of information just to store the charge density.

This is in general too much for a computer to cope with, and obviously most boxes just won't get visited in a normal million move simulation..

In the atomic case, CASINO uses a simplification i.e. the spherical symmetry. If atoms have a spherical charge density, then I don't need to partition the space into lots of little boxes as before; I partition it into lots of little spherical shells centred on the nucleus, and I count 1 every time an electron lands between a radius r and a radius r+dr. Now I can use 5000 shells, and that will give me a nice resolution - and the density can be store in 20 Kb. Note that I have to worry about the volume here, as not all the shells have the same volume (see the manual for how to do that).

This is why CASINO can do atoms and not molecules.

OK now, I don't want to be seen to give you projects (I'm not your supervisor after all!) but feel free to think about what might be a better algorithm
to accumulate a molecular charge density... (there's only something like 137,000 people subscribed to this forum now, and I know that each one of them will be watching this thread with great interest to see your solution - the whole of your future career could depend on impressing them!).

If that isn't motivation enough, note that we're so desperate for active developers these days that there is even money, and sweets, and beer and so on available as a reward.. I know at least you like beer.

Mike

PS: Is an atom spherical? Is an atom like boron 1s2 2s2 2p1 with '1 p orbital occupied' spherical? Well, as Bill Clinton would say, that depends on what the meaning of the word 'is' is - a common problem in quantum mechanics.

The obvious - and indeed quite correct - answer is 'yes' because of the rotational symmetry. Take your nonspherical boron, rotate it an infinitesimal amount, and you get a rotated density with a wave function that must be degenerate with the unrotated atom. Taking the whole set of such wave functions, you can make linear combinations and end up with the 2L+1 functions that form an irreducible representation of the rotation group etc. etc. Applying any symmetry operator to the wave function should give back the original wave function, except for an irrelevant phase factor. Thus, an exact wave function of the boron atom has to be spherically symmetric.

That said, the actual answer as far as we're concerned is often 'no' because UHF and DFT calculations often give 'broken symmetry solutions' with e.g. particular p orbitals occupied, and these are just not spherically symmetric.. In principle one therefore ought to use some sort of multireference procedure allowing equal weighting of degenerate occupations, e.g. where the HF wave function is a linear combination of three determinants each with a different p orbital occupied. This is in fact how typical atomic SCF codes work. With general purpose HF/DFT codes, it isn't, and the usual way of dealing with that problem is to ignore it! In practice, the energy differences involved are usually much smaller than other sources of error in the calculation.

Note that in such a case, the CASINO procedure above is giving you the 'spherical average' of the non-spherical density (this issue is almost
certainly not mentioned in the manual..).

Now of course, coming back to President Clinton's question, it's amusing to note that using interpretations of quantum mechanics that are clear about what they believe exists, it is perfectly possible for a boron atom to be non-spherical.

In the usual interpretation, people are very careful not to make statements about what 'is', just about the results of measurements. Quantum mechanics is not a dynamical theory of the motion of particles, but a statistical theory of measurement. In between measurements, what happens is 'meaningless', or in the more extreme cases, the particles 'do not exist'. You might - possibly quite rightly, think that too (but then you should consider very carefully the ethics of doing molecular dynamics calculations next time you fire up the CASTEP program or whatever..) :-)

If you calculate your density by an ensemble of thousands of position measurements, it may actually be the case that you are measuring the positions of electrons associated with a non-spherical atom, but that your state preparation procedure does not fix the orientation of the atom in space, and the spherical density is thus an artefact of 'ensemble averaging' of the orientation.

For an interpretation which is clear about its ontology ('what exists') consider the de Broglie-Bohm one, which reproduces all the experimental
predictions of QMC, but in which electrons - shock, horror - exist and have trajectories, and the Schroedinger wave function is the mathematical
representation of an objectively existing 'energy field'. If an atom genuinely 'is' in a state with a '1 p orbital occupied' Schroedinger wave - and of course this is a solution of the Schroedinger equation - then the resultant electron trajectory looks like this (the result of an actual calculation):

Image

Not very spherical..

Note that deBB is in fact a dynamical theory of the motion of particles (which nevertheless reproduces all the results about measurement) so all the molecular dynamics people are in fact secret de Broglie-Bohmists.

All this is of course irrelevant to your question. It's just interesting to note that when one says definitely 'an isolated atom IS spherical', the
universe may in fact be so arranged that it is not, and their apparent sphericity is just an artefact of the measurement process. So, as I said, it
all depends on what the meaning of the word 'IS' is.

Frankly, my dear, I don't give a damn, I hear you say. You can never prove it one way or the other anyway. Quite right too!
Yasmine_AlHamdani
Posts: 3
Joined: Tue Feb 04, 2014 4:29 pm

Re: Charge density for finite molecules

Post by Yasmine_AlHamdani »

Hi Mike,

Thanks for that thorough response and explaining the issue.
It's generated some discussion amongst my peers so it's good to put it out here and I, amongst others I'm sure, will keep it in mind...

Cheers,

Yasmine
Post Reply