Home

CASINO manual

1. r but with the constant value 39 chosen so that the effective one electron local energy is continuous at Te Hydrogen is treated as a special case as the 1s orbital of the isolated atom is only half filled and we use Hideal p By 0 5 We wish to choose 0 so that E r is as close as possible to EX r for 0 lt r lt re i e the effective one electron local energy is required to follow the ideal curve as closely as possible In our current implementation we find the best 0 by minimizing the maximum square deviation from the ideal energy E r Eideal 7 2 within this range Beginning with 4 0 0 we first bracket the minimum then refine 0 using a simple golden section search In principle we are more interested in Ej r being close to Fideel r near re than near zero because the probability of an electron being near re is normally much greater than it being near the nucleus One might therefore consider using a 78 weighting factor and minimizing e g r E r Eis r In practical calculations this was found not to improve the result in general and weighting factors were not used in our final implementation It is clearly important to find an automatic procedure for choosing appropriate values of the cusp correction radii Although the final quality of the wave function in QMC calculations is expected to have only a relatively weak dependence on its precise value the optimal cusp correct
2. 0000000000044777 0000000000000000 0000000000000002 0089256999176878 53 Output files 7 9 vmc hist file This file records information about each block of a VMC run The first line gives the average of the local energy over the block and its standard error in Hartrees per primitive cell per electron for HEG The second line contains six integers and one energy i 2 Di 7 Number of moves in the block at which the energy is evaluated NMOVE Number of ions Number of up spin electrons Number of down spin electrons Periodicity the number of dimensions in which the system is periodic 0 3 Number of energy terms in file to average needed by the REBLOCK utility There are fifteen energy terms at present see below The ion ion Coulomb interaction energy per ion Below these for each of the moves in the block rows containing fifteen energy terms are given 10 11 12 13 14 15 o ON Don pc Local energy calculated using the 1 r or Ewald Coulomb interaction energy Total potential energy KE 1 Expectation of the kinetic energy operator see Section 10 17 Ti See Section 10 17 F 12 Square modulus of the drift vector see Section 10 17 Electron electron interaction energy 1 r Ewald or MPC Electron ion interaction energy Non local energy MPC energy Short range part of MPC Long range part of MPC Electron ion CPP term Electron pa
3. 14 1 5 Linear scaling QMC Tf the numbers of orbitals and basis functions to be evaluated are both O 1 then the CPU time for a configuration move scales as O N This is what is meant by linear scaling QMC Note however that the number of configuration moves required to achieve a given error bar scales as O N hence in practice linear scaling QMC calculations scale as O N and it is better to refer to them as quadratic scaling calculations This has not however been the practice in the literature 14 2 Using CASINO to carry out linear scaling QMC calculations 14 2 1 Generation of Bloch orbitals At present the Bloch orbitals must be represented in a plane wave basis A plane wave DFT or HF code should be used to generate a pwfn data file as though an ordinary QMC calculation with a plane wave basis were to be carried out Note that the following restrictions apply 1 The simulation cell must be orthorhombic 2 The same orbitals must be available for spin up and spin down electrons 3 Only one k point may be used the T point The first two restrictions may be lifted in a later release 14 2 2 Transformation to Wannier orbitals The xwannier utility written by R Q Hood can be used to generate maximally localized Wannier functions using the method of Berghold et al 11 This program requires a pwfn data file along with an input wannier file of format 4 No init states included in Wannier tr
4. none PERIODICITY system dimension 0 1 2 3 gt molecule polymer slab solid SPIN_UNRESTRICTED true or false i e different orbitals for different spins EIONION nuclear nuclear repulsion energy au atom NUM_ELECTRONS number of electrons per primitive cell GEOMETRY NUM_ATOMS number of atoms per primitive cell BASIS 3 NUM_ATOMS atomic positions au ATNO NUM_ATOMS atomic numbers for each atom VALENCE_CHARGE NUM_ATOMS valence charges for each atom bhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhh PERIODIC INSERTA ALALLA Ahh hhhhhhhhhhhhhhhhhh LATTICE_VECTORS 3 3 primitive lattice vectors au K SPACE NET 41 NUM_K no of k points in BZ NUM_REAL_K no of real k points all components of real k points are either zero or half a reciprocal lattice vector KVEC 3 NUM_K k point coordinates a u NB coordinates of real k points must occupy the first num_real_k positions in kvec bhhhhhhhhhhhhhhhhhhhhhhhhhhhhhh END PERIODIC INSERT bbhhhhhhhhhhhhhbhhhhh BASIS SET NUM_CENTRES number of centres with associated Gaussian functions i e number of atoms number of non atom Gaussian sites per primitive cell NUM_SHELLS number of shells per primitive cell NUM_AO number of basis functions AO per primitive cell NUM_PRIMS number of Gaussian primitives per primitive cell HIGHEST_ANG_MOM highest angular momentum of shells SHELL_AM NUM_SHELLS c
5. solution with V v a can be made position dependent in a fashion that roughly distinguishes between the drift velocity being large due to proximity to a node and proximity to a nucleus This is described in Section 10 7 However because the dependence on a is very weak we can simply choose a constant value of a 0 1 see Section 10 6 2 if bare nuclei are not present Note that the limited velocity is only substantially different from the unlimited velocity when the magnitude of the latter is large compared to 1 2aT Also the magnitude of the limited velocity is always less than the unlimited velocity So Equation 37 is indeed a limiting scheme for the single electron drift velocity By proposing electron moves using the limited drift velocity we may improve upon the short time approximation since the variation of the drift velocity over moves close to nodes is taken into account Furthermore although the unlimited velocity diverges as the node is approached the limited velocity remains bounded it tends to 2 ar Thus the various pathologies caused by the divergent drift velocity are ameliorated 10 6 4 Theoretical background to the UNR limiting scheme for the local energy Although CASINO uses an electron by electron DMC algorithm the derivation below is given in the context of an entire configuration move When close to the nodal surface the limited energy is calculated as the average of the local energy over a configura
6. Eq 138 contains 5 parameters for each ion a 7 7or 711 727 and fr whose values are entered at the end of the xx_pp data file see Section 7 3 Suitable values of the parameters are given in the paper by Shirley and Martin 32 If fo T r Far Fr the second term in Eq 138 is zero and it is not calculated The second term in Eq 138 is short ranged because f x 1 at large x This term is calculated in real space The second term in Eq 138 is added to the pseudopotential and the core radii LCUTOFFTOL and NLCUTOFFTOL are determined from the resulting potential The electric field evaluation is 86 activated by the presence of core polarization terms in the pseudopotential files they are not calculated by default since they may be expensive especially when periodic boundary conditions are used In periodic boundary conditions the electric fields are evaluated directly from the analytic first derivatives of the Ewald potential see section 10 20 Calculations using CPPs may be 5 10 percent slower than ones without CPPs in periodic systems NB the first derivatives of the periodic potential in 1D polymers has not yet been implemented and thus the core polarization energy cannot be evaluated in such systems 10 20 Evaluation of infinite Coulomb sums 10 20 1 3D Ewald interaction In three dimensionally periodic systems the periodic potential of a neutralized lattice of point charges may be evaluated using the Ewald method 1 2 Cons
7. cally gives the information about the pair correlation function When PCF_SPH_MODE has a non zero value the accumulation on such radial bins is performed In order to get the spherically averaged pair correlation function the radial distribution of 1 n r should be normalized by the radial distribution of unity i e the bin volume So users should first accumulate the normalized distribution with PCF_SPH MODE 1 and then accumulate 1 n r with PCF SPH_MODE 2 Theoretical background is given in the paper PRB 68 165103 2003 PERMIT_DEN_SYMM Logical If this flag is set then 1 it will be assumed that if the SCF charge density expansion coefficients are real then the QMC ones are too 2 the symmetry of the SCF charge density will be imposed upon the QMC charge density data for use in the MPC interaction It is possible that DMC will break the symmetry of the SCF calculation in this case the user should turn off PERMIT_DEN_SYMM POPRENORM Logical If poprenorm is T then in DMC after every block the number of con figs walkers on each node will be renormalized to NCONFIG This is achieved in subroutine RENORM using the following method 1 each node tells node 0 how many configs it now has 2 Node 0 instructs random configs on random nodes to be deleted or copied until the total number equals NCONFIG NNODES Note that it is not normally required to do this and you should turn it on only when your run is showing abnormally large population flu
8. J Needs and G Rajagopal Phys Rev B 59 12344 1999 95 C Filippi and C J Umrigar J Chem Phys 105 213 1996 95 S Manten and A L chow J Chem Phys 115 5362 2001 76 K E Schmidt and J W Moskowitz J Chem Phys 93 4172 1990 95 N Metropolis A W Rosenbluth M N Rosenbluth A M Teller and E Teller J Chem Phys 21 1087 1953 62 S Fahy X W Wang and S G Louie Phys Rev B 42 3503 1990 62 85 D M Ceperley and M H Kalos in Monte Carlo Methods in Statistical Physics edited by K Binder Springer Berlin 1979 K E Schmidt and M H Kalos in Monte Carlo Methods in Statistical Physics I edited K Binder Springer Berlin 1984 M Dewing J Chem Phys 113 5123 2000 63 B L Hammond W A Lester Jr and P J Reynolds Monte Carlo methods in ab initio quantum Chemistry World Scientific Singapore 1994 64 73 W M C Foulkes L Mitas R J Needs and G Rajagopal Rev Mod Phys 73 33 2001 62 63 64 67 V R Saunders R Dovesi C Roetti M Caus N M Harrison R Orlando and C M Zicovich Wilson CRYSTAL98 User s Manual University of Torino Torino 1998 4 V R Saunders R Dovesi C Roetti R Orlando C M Zicovich Wilson N M Harrison K Doll B Civalleri I Bush Ph D Arco M Llunell CRYSTAL2003 User s Manual University of Torino Torino 2003 4 128 26 27 28 29 30 31 32 33 34 35 36 37 38 39
9. Ti 94 B Gat where the Gz are the reciprocal lattice vectors of the primitive unit cell belonging to the Bth star of vectors that are equivalent under the space group symmetry of the crystal and the means that if Gg is included in the sum G p is excluded 10 14 4 Old Jastrow factor The old redundant Jastrow factor jasfun data used by CASINO before Easter 2004 is described in this section as it is still supported can be written as follows We denote the distance of electron i from ion I by riz and the distance between electrons and j by Tij J Y uolrig Sry YOYO Salrir a I i gt j j y 5 Ss rir rjr Sa riz rij Sa rgr rij S5 rir r5r rij 95 i gt jj I The one electron term S2 is allowed to depend on the electron spin while the other terms which all involve two electrons are allowed to depend on the relative spins At the moment the terms for two up spin electrons are the same as for two down spin electrons In three dimensions the function uy is taken to be A 7 7 uo riz 1 exp J lexp 2 96 1 en F where F is chosen so that the cusp conditions are obeyed i e Fy V2A and Fi VA Lois chosen so that uo is effectively zero at the shortest distance to the surface of the Wigner Seitz cell Lws We normally take Lo 0 3Lws and A is either set to the inverse of the plasma energy or obtained by an optimization procedure In two dimensional electron and el
10. Voth P Salvador J J Dannenberg V G Zakrzewski S Dapprich A D Daniels M C Strain O Farkas D K Malick A D Rabuck K Raghavachari J B Foresman J V Ortiz Q Cui A G Baboul S Clifford J Cioslowski B B Stefanov G Liu A Liashenko P Piskorz I Komaromi R L Martin D J Fox T Keith M A Al Laham C Y Peng A Nanayakkara M Challacombe P M W Gill B Johnson W Chen M W Wong C Gonzalez and J A Pople Gaussian Inc Pittsburgh PA 2003 4 M D Segall P L D Lindan M J Probert C J Pickard P J Hasnip S J Clark and M C Payne J Phys Cond Matt 14 11 pp 2717 2743 2002 4 8 First principles computation of material properties the ABINIT software project X Gonze J M Beuken R Caracas F Detraux M Fuchs G M Rignanese L Sindic M Verstraete G Zerah F Jollet M Torrent A Roy M Mikami Ph Ghosez J Y Raty D C Allan Computational Materials Science 25 478 492 2002 4 8 S Baroni A Dal Corso S de Gironcoli and P Giannozzi http www pwscf org 4 8 W Miller J Flesch and W Meyer J Chem Phys 80 3297 1984 86 E L Shirley and R M Martin Phys Rev B 47 15413 1993 85 86 D M Ceperley M H Kalos and J L Lebowitz Computer simulation of properties of a polymer chain Macromolecules 14 1472 1981 66 J B Foresman and M J Frisch Exploring Chemistry with Electronic Structure Methods Gaus sian Inc Pittsb
11. however both must be isotropic There should only be one set of u terms and the spherical harmonic l and m values should therefore be set to 0 e The expansion order of u determines the number of parameters Typically it is given a value between 4 and 8 e The spin dep params line allows the user to specify whether the same u parameters are to be used for parallel and antiparallel spin electron pairs If the value is set to 0 then the same parameters are used for parallel and antiparallel pairs this option should not be used in general 34 if the value is set to 1 then different parameters are used for parallel and antiparallel pairs if it is set to 2 then different parameters are used for up spin down spin and opposite spin pairs this is useful for spin polarized systems The cutoff length is given a default value if it is set to 0 Note that it is possible to optimize the cutoff length using variance minimization This can be useful because the choice of cutoff length has a significant effect on the optimized energy and variance Unfortunately optimizing the cutoff length is numerically difficult variance minimization will take many more iterations to converge if the cutoff is optimizable Setting the truncation order to 3 can help somewhat It is possible to specify exactly which of the expansion coefficients parameter values are optimizable and which are not One does not need to specify all of the expansion co
12. is activated each orbital expressed as a sum of s type Gaussians is replaced within r re by ISIGNxezxp p r where p r is the polynomial p a agr asr aar asr such that Pre p re p re are continuous p 0 Znuc EL 0 is set by a smoothness criterion This procedure can be expected to greatly reduce fluctuations in the local energy in all electron Gaussian calculations See Section 10 13 for more details CUSP INFO Logical If CUSP_CORRECTION is set to true for an all electron Gaussian basis set calculation then CASINO will alter the orbitals inside a small radius around each nucleus in such a way that they obey the electron nuclear cusp condition IF CUSP_INFO is set to true then information about precisely how this is being done will be printed to the output file Be aware that in large systems this may produce a lot of output Furthermore if you create a file called orbitals in containing an integer triplet specifying which orbital ion spin you want the code will additionally print graphs of the specified orbital radial gradient Laplacian and one electron local energy to the files orbitals dat gradients dat laplacians dat and local_energy dat These graphs may be viewed using xmgr grace or similar plotting programs See Section 10 13 for more details DBARRC Jnteger Maximum number of blocks between recalculation of DBAR matrices see Sec tion 10 15 DENFT_THRESHOLD Real When calculatin
13. k k pair the presence of the other is assumed To make a many body Bloch wave function satisfying the condition that it is multiplied by a phase factor under the replacement ri gt r R where R is a translation vector of the simulation cell the offset must satisfy ks 1 2G where G is a reciprocal lattice vector of the simulation cell 12 see Section 10 12 Finite systems are treated as if they had a single k point 74 10 11 2 Plane wave basis set Soon 10 11 3 Blip function basis set Soon 10 11 4 Atomic wave functions on a grid Soon 10 11 5 Wannier orbital evaluation Soon 10 12 Constructing real orbitals The Bloch orbitals at an arbitrary point in k space are complex If the set of wavevectors consists of k pairs then one can always construct a set of real orbitals spanning the same space as the original complex set A necessary and sufficient condition for the mesh of Eq 61 to consist of k pairs is that k G 2 where G is a reciprocal lattice vector the simulation cell lattice It is four times more efficient to use real orbitals than complex ones because it takes four multiplications to evaluate the product of two complex numbers but only one multiplication for the product of two real numbers An orbital satisfying Bloch s theorem can be written as dx r u_ r e T y 62 where ux has the periodicity of the primitive lattice The function is a Bloch function with wavevector k
14. sians per angular momentum channel e G98 s ECP integral package crashes when one attempts to do a large both in terms of basis and ECP expansion calculation with symmetry switched on The solution to this is to switch symmetry off using Nosymm e One cannot use basis functions containing g and higher basis functions with a pseudopotential in G94 e Obviously gaussiantoqmc needs the molecular orbitals MOs produced by the calculation It gets these from the formatted checkpoint file which is produced by putting Form check Basis MO in the route section of the Gaussian job Alternatively it may be obtained from the binary checkpoint file chk using the formchk utility see the Gaussian manual e Many but not all of the IOp s mentioned here are described on Gaussian s website at http www gaussian com iops htm 13 2 3 CIS e Performing a HF test for an excited state It is possible to get Gaussian to output a breakdown of the energy of a CIS excited state which may be compared with the results of a determinant only VMC run The key to this is the density used to perform the population analysis and other post SCF calculations By default Gaussian uses the density produced by the original SCF run To get the kinetic nuclear nuclear potential and electron nuclear potential energies you must tell it to use the one particle CI density via density RhoCI You must also specify the excited state that you ar
15. will turn off the collection of subroutine timings You might want to do this as the the timing routines can adversely affect system performance on certain computers such as Alpha or PC clusters especially for small systems TPDMC Integer TPDMC T is the number of timesteps for which the effects of changes in the theoretically constant reference energy should be undone in order to estimate the DMC energy at a given point It is assumed that the best estimates of the DMC energy separated 29 by an amount greater than this are not correlated by fluctuations in the reference energy Thus Tp should exceed the timescale of fluctuations in the reference energy Umrigar suggests using T 10 7 where 7 is the timestep If you set it to 9999 in the input then the code will automatically use this value if you set it to 0 then the reweighting scheme for population control biasing will not be used The latter is the default See Section 10 8 TRIP_POPN Integer In the course of a DMC simulation it is possible for a configuration pop ulation explosion to occur If TRIP_POPN is set to 0 then nothing will be done about this If TRIP_POPN is greater than 0 then it will attempt to restart the block if the single node popu lation exceeds TRIP_POPN A general suggestion for its value would be two times NCONFIG See the discussion in Section 10 10 for more details USE_COEFF_FILE Logical Use a coefficient file with the New Optimization Method
16. 0983 40 078 44 95591 47 867 50 9415 51 9961 54 93805 55 845 58 93320 58 6934 63 546 65 39 69 723 72 61 74 92160 78 96 79 904 83 80 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 Rb Sr Y Zr Nb Mo Tc Ru Rh Pd Ag Cd In Sn Sb Te I Xe 85 4678 87 62 88 90585 91 224 92 90638 95 94 98 0 101 07 102 90550 106 42 107 8682 112 411 114 818 118 710 121 760 127 60 126 90447 131 29 55 56 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 Cs Ba Lu Hf Ta W Re Os Ir Pt Au Hg Tl Pb Bi Po At Rn 132 90545 137 327 174 967 178 49 180 9479 183 84 186 207 190 23 192 217 195 078 196 96655 200 59 204 3833 207 2 208 98038 209 0 210 0 222 0 87 88 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 Fr Ra Lr Rf Db Sg Bh Hs Mt Ds Uuu Uub Uut Uuq Uup Uuh Uus Uuo 223 0 226 0 262 0 261 0 262 0 263 0 264 0 265 0 268 0 57 58 59 60 61 62 63 64 65 66 67 68 69 70 Lanthanoids La Ce Pr Nd Pm Sm Eu Gd Tb Dy Ho Er Tm Yb 138 9055 140 116 140 90765 144 24 145 0 150 36 151 964 157 25 158 92534 162 50 164 93032 167 26 168 93421 173 04 89 90 91 92 93 94 95 96 97 98 99 100 101 102 Actinoids Ac Th Pa U Np Pu Am Cm Bk Cf Es Fm Md No 227 0 232 0381 231 03588 238 0289 237 0 244 0 243 0 247 0 247 0 251 0 252 0 257 0 258 0 259 0 Table 2 The default nuclear masses used in CASINO 101 11 Making movies with CASINO 11 1 How to make movies In the input
17. 2 there is an initial period of equilibration which should be excluded from the averaging We use the blocking method 9 in which adjacent data points are averaged to form block averages 1 a 5 t2i 1 t 172 This procedure is carried out recursively so that after each blocking transformation the number of data points is reduced by one half For simplicity we assume that the number of data points M is a power of two The mean value of the block means is unchanged by the block transformations but the variance of the block means tends to increase with the number of block transformations 1 Mp 2 X 2 52 j l where Map is the number of blocks Zp is the mean value of the jth block and Z is the overall mean value The error in o is estimated from 2 2 l 174 M M 1 ae The value of o increases with block length Mp until a limiting value is reached which corresponds to the true standard deviation of the mean For very large block sizes the error in op becomes large and if M is not large enough a limiting value of o is not reached but because it is simple to calculate the error in oy such behaviour is easy to identify These operations are carried out by the utility reblock which reads data from a dmc hist or vmc hist file 10 22 1 Estimate of the correlation time given by CASINO The correlation time of the energy is computed and shown for every block in a VMC calculation and also when using
18. 3 3 6 4 5 12 5 5 18 6 7 26 7 11 50 Within a VMC calculation it is often possible to use a low order quadrature rule because the error cancels over the run but higher accuracy is required for wave function optimization and DMC calculations which are biased by errors in the non local integration In principle the non local energy should be summed over all the ionic cores and all electrons in the system However since the non local potential of each ion is short ranged one need only sum over the few atoms nearest to each electron Exact sampling of the non local energy with DMC is problematic and we use the localization approximation in which the non local operator acts on the trial wave function in exactly the same way as in VMC The error introduced by this approximation is proportional to Y Yo 3 where Vo is the exact wave function 10 19 The core polarization potential energy Core polarization potentials CPPs account for the polarization of the pseudo ion cores by the fields of the other charged particles in the system The polarization of the pseudo ion cores by the fields of the valence electrons is a many body effect which includes some of the core valence correlation energy 32 In the CPP approximation the polarization of a particular core is determined by the electric field at 85 the nucleus The electric field acting on a given ion core at Ry due to the other ion cores at Ry and the electrons at r is R D Ae
19. AAU dme bist and dime Dista les oa aa we as e A A 8 Specifying the Slater determinants 8 1 Gaussian basis SEUS is oca a a A Aa A e 8 2 Plane wave Dasis ciu a eee eee eae ed eee eee eed So B pand spline basis c paa ai ea Aa Aa eS RG BA ee ee ROE 9 Pseudocode 91 Subroutme VMG ss s cs a ag ad A QoS ee a a a aeons Bee Ba acaba 92 Subroutme DMC sesei dada aa oS a SS aE SRO a a ew eB ad ww Ss 10 Theoretical background 10 1 The trial wave function s sa s 224 soaa eee eee a RA 10 2 The variational Monte Carlo method 0 000000 10 3 The diffusion Monte Carlo method 2 6 cee ae ass we PE e AAA 104 Dif and USOS lt y 2 6 4 pace a ew eR Ea ae CeO Ee ee Ee A 105 Branching and population Control lt coc 3 ce a oR PA AR eRe GS A 10 6 Modifications to the Green s Function 0 0202002020005 10 7 Modifications to the DMC Green s function at bare nuclei 10 8 Evaluating expectation values of observables 0 0 0 20000000008 10 9 Growth estimator of the energy aa E 10 10Automatic block resetting o 10 11Evaluation of orbitals in the determinant part of the wave function 10 12 Constructina real orbitals sak ke e ea a Se a ee BEA 10 13Cusp corrections for Gaussian orbitals s s s sa saaa ave o ee es w uw uw w N A 10 14 The Jastrow factor ise a a ee ees 1015 Waxe function Update s4 raaa aada eS AO RG a oe a Yee
20. CRYSTAL run Including the qmc flag as an argument to the run script means these files will be grabbed and renamed as silicon f12 silicon f10 and silicon f30 or whatever They will be kept in the scratch directory since they can become very large Sit in the directory where these three files live and type crysgen This will run the crystaltogmc program which will ask you some questions then generate the gwfn data file crystaltogmc will ask you for the size of the Monkhorst Pack k point net in the CRYSTAL cal culation and the desired size of supercell in the QMC calculation These need not be the same if not you are plucking to use the local vernacular but the former must be divisible by the latter For example 12x12x12 MP net in CRYSTAL will allow you to generate gwfn data for 1x1x1 2x2x2 3x3x3 4x4x4 and 6x6x6 supercell cases Note that it is desirable to carry out the calculation on a higher density k point grid than you actually need so that the orbital coefficients are calculated accurately 104 Note that for polymer and slab calculations the last one and two numbers in the MP net and supercell specifications should be 1 to reflect the fact that the system is not periodic in those dimensions For example polymer 12x1x1 1x1x1 2x1x1 3x1x1 4x1x1 6x1x1 slab 12x1x1 1x1x1 2x2x1 3x3x1 4x4x1 6x6x1 For molecular calculations just imagine you have a 1x1x1 MP net So to summarize to produce gwfn data from th
21. H 178 CG t 1 H where the numerator is the average of the measured values of its arguments over the configurations indexed by 1 lt t lt N t and 0 is the variance of these measures One possibility for setting the cut off is to check against the self consistent inequality tmar lt 3T tmax while computing the sum and truncate it as soon as it stops being true This allows for an estimate of the error in above expression to be calculated 2 2tmaz 1 179 Er az P 10 22 2 Further considerations for DMC Each iteration is equally weighted in a VMC calculation however for DMC each iteration is weighted by the total weight of the configurations multiplied by the Il weight In either case let the iteration weights be w the total number of data points be M and the energy from iteration i be e Consider the bth reblocking transformation in which the block length is By 2 The data range may be divided into M blocks the last of which is usually incomplete For each block j the block weight is Wo wi 180 tej and the corresponding block energy is 181 The reblocked energy is 2 Fj Wo eS j Woj Fetal F 182 Wi which is therefore independent of reblocking transformation On the other hand the reblocked variance is E 2 2j We Hej Ey On Sw 2j Woj Sw 183 which does depend on the reblocking transformation number 93 The number of blocks
22. Open e Click on Display and untick the box for Bonds e Click on Extras Animate An animation tool bar will appear To start the movie click on the play symbol 103 12 Using CASINO with the CRYSTAL program see CASINO utils wfn_converters crystal_9x and crystal_03 12 1 The CRYSTAL program CRYSTAL is a commercially available quantum mechanical electronic structure package which is able to calculate the electronic structure both of molecules and of systems with periodic boundary condi tions in 1 2 or 3 dimensions polymers slabs and crystals using either Hartree Fock or density func tional theory The latest version of the program was written by Roberto Dovesi Vic Saunders Carla Roetti Roberto Orlando Nic Harrison Claudio Zicovich Wilson Klaus Doll Bartolomeo Civalleri lan Bush Philippe D Arco and Miquel Llunell See www tcm phy cam ac uk mdt26 crystal html and the various links thereon for more information CASINO supposedly supports the official releases CRYSTAL95 CRYSTAL98 and CRYSTALO3 but not 88 or 92 Later versions of CRYSTALO3 contain a routine which will write out a CASINO gwfn data file firectly If you want to use CRYSTAL95 or CRYSTAL98 then you must use a separate utility crystaltogmc which transforms the output of CRYSTAL into a gwfn data file 12 2 Generating gwfn data files with CRYSTAL95 98 and crystaltoqme The utility crystaltogmc written by MDT together with the driver script crysgen i
23. Therefore we can make two real orbitals from fx and fy as follows 1 Py r a ox 1 Ha 1 u r dia 63 The orbitals 4 and _ are orthogonal if fx and dj f_x are orthogonal which is true unless k k Gy i e their wavevectors differ by a reciprocal lattice vector of the primitive lattice In this case 6 and _k are linearly dependent and we must use only one of them Therefore the scheme is Casel Ifk F use Py and g_ Case2 Ifk use d or g_ 64 IK 2 In the second case if one of py or p_ is zero then obviously one must use the other one It may happen that we have in our k point grid the vectors k and k which are not related by a reciprocal lattice vector of the primitive lattice We may wish to occupy only one of the orbitals from these k points We can then form a real orbital as cos 4 sin 0 _ where 0 is a phase angle between zero and 27 If both k and k orbitals are supplied CASINO chooses one which in general is sufficient to generate all linearly independent orbitals 10 13 Cusp corrections for Gaussian orbitals One of the problems with Gaussian basis sets is that they are unable to describe the cusps in the single particle orbitals at the nuclei that would be present in the exact HF orbitals because the Gaussian basis functions have zero gradient at the nuclei on which they are centered This can lead to considerable difficulties in QMC simulations In bot
24. VMC energies have equilibrated The user may choose whether to optimize the Jastrow function determinant expansion coefficients or the pairing parameter or Pad coefficients in a Wigner crystal by setting the vm_opt_jasfun vm_opt_detcoeff and vm_opt_pairing flags as appropriate For most applications it is only nec essary to optimize the Jastrow function Exactly which parameters in the Jastrow function are to be optimized are indicated in the jastrow data file See Section 10 14 When optimizing a Jastrow function if non linear parameters such as the cutoffs are not being optimized it is possible to accelerate the variance minimization process by setting vm_mode to linear CASINO will then exploit the fact that the Jastrow function has a linear dependence on all the parameters appearing in it enormously speeding up the process of recomputing configuration energies when the parameters are changed This is expensive in memory however If memory requirements prevent the use of linear mode or cutoffs are to be included in the optimization then vm_mode should be set to direct and the local energies of each of the configurations will be calculated from scratch when the parameters are changed NOTE THIS STILL ONLY WORKS FOR THE OLD JASTROW FACTOR TO BE CHANGED SOON MDT 6 2004 When optimizing Pad coefficients the condition that the coefficients are the same for electrons of either spin can be enforced by setting the opti
25. We are of course always happy to discuss such requests 127 References 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 P P Ewald Ann Phys 64 25 1921 87 M P Tosi in Solid State Physics Vol 16 edited by H Ehrenreich and D Turnbull Academic New York 1964 p 1 87 L Mit s E L Shirley and D M Ceperley J Chem Phys 95 3467 1991 85 P J Reynolds D M Ceperley B J Alder and W A Lester Jr J Chem Phys 77 5593 1982 66 68 C J Umrigar M P Nightingale and K J Runge J Chem Phys 99 2865 1993 65 68 69 12 T3 C J Everett and E D Cashwell 4 Third Monte Carlo Sampler Los Alamos Technical Report No LA 9721 MS 1983 72 C J Umrigar and C Filippi unpublished 65 68 M F Depasquale S M Rothstein and J Vrbik J Chem Phys 89 3629 1988 69 H Flyvbjerg and H G Petersen J Chem Phys 91 461 1989 92 L M Fraser W M C Foulkes G Rajagopal R J Needs S D Kenny and A J Williamson Phys Rev B 53 1814 1996 89 A J Williamson G Rajagopal R J Needs L M Fraser W M C Foulkes Y Wang and M Y Chou Phys Rev B Rapid Communications 55 4851 1997 89 P R C Kent R Q Hood A J Williamson R J Needs W M C Foulkes and G Rajagopal Phys Rev B 59 1917 1999 74 89 C J Umrigar K G Wilson and J W Wilkins Phys Rev Lett 60 1719 1988 P R C Kent R
26. a random number drawn from a uniform distribution on the interval 0 1 In weighted DMC each configuration carries a weight that is simply multiplied by M a m after each move only if the weight of a configuration goes outside certain bounds currently above 2 or below 0 5 is it allowed to branch or be combined with another configuration Throughout this section we denote the best estimate of the ground state energy at time step a by Evest m At the start of a DMC run we set Ebest 0 Er 0 Ev where Ey is the average local energy of the initial configurations During equilibration Ebest is updated after each block of moves as 3 1 Exest m q Polock q Pest m a 1 25 At the end of each block various book keeping tasks are performed i e in VMC data is written to file at the end of each block and in DMC data is written out the population may be renormalized POPRENORM and Epest is updated during equilibration Er is updated after every time step as 1 Er m 1 Ebest m 2 log E 26 67 where g min 1 7cg cp is a constant which must be set in the INPUT file CEREFDMC but is usually set equal to one My NCONFIG x NNODES is the target number of configurations and Noontig M Mirom 5 walm 27 a 1 where Noonfig m is the number of configurations and wa is the weight of configuration a Note that Evest m is the best energy at time step m while Er m 1 is the trial energy to be used in the n
27. at the origin is 52 Zac 65 Or r 0 where W is the spherical average of the many body wave function about r 0 For a determinant of orbitals to obey the Kato cusp condition at the nuclei it is sufficient for every orbital to obey Eq 65 at every nucleus We need only correct the orbitals which are non zero at a particular nucleus because the others already obey Eq 65 This is sufficient to guarantee that the local energy is finite at the nucleus provided at least one orbital is non zero there In the unlikely case that all of the orbitals are zero at the nucleus then the probability of an electron being at the nucleus is zero and it is not important whether Y obeys the cusp condition An orbital 4 expanded in a Gaussian basis set can be written as p o n 66 where is the part of the orbital arising from the s type Gaussian functions centered on the nucleus in question which for convenience is at r 0 and 77 is the rest of the orbital The spherical average of Y about r 0 is given by Y n 67 In our scheme we seek a corrected orbital 4 which differs from 4 only in the part arising from the s type Gaussian functions centered on the nucleus i e p or n 68 76 The correction 4 4 is therefore spherically symmetric about the nucleus We now demand that Y obeys the cusp condition at r 0 dla 42 Zio 69 To Note that 7 is cusp less because it arises from the Gaussian basis
28. data file is unformatted There exists a utility called format_spline that reads in a swfn data file and produces a formatted swfn data formatted file and vice versa This is useful if swfn data is to be transferred from one platform to another It is also possible to insert a title in the formatted file if desired 14 2 4 Running a QMC calculation The swfn data file can now be used by CASINO The usual input files input jastrow data and pseudopotentials should also be supplied The btype input parameter should be set to 6 Please note the following e The isperiodic input parameter may be set to either true or false In the former case the orbitals are evaluated using the minimum image of the distance to the Wannier centre and a periodic electron electron interaction is used In the latter case the contents of the simulation cell are treated as a finite system Although the original orbitals were represented in a periodic plane wave basis once the localized orbitals have been constructed represented by splines and truncated we may dispense with the periodicity e An explicit list of the occupied orbitals must be provided using statelist_up and statelist_down blocks in the input file For example block statelist_up 1234 Zendblock statelist_up means that orbitals 1 2 3 and 4 are occupied by spin up electrons An excited state could be specified by replacing the 4 with a 5 This differs from the way in which excitations are
29. file input contains all the parameters needed to control the QMC calculation A complete list of the input parameters is given below Further details including default values may be found by using the casinohelp utility Type casinohelp all to get a list of all keywords that CASINO knows about or casinohelp keyword for detailed help on a particular keyword Type casinohelp search text to search for the string text in all the keyword descriptions Create a blank file containing the keyword input_example and type runvmc this will create a sample input file containing all valid input keywords and their default values in the correct format The input file is meant to be very flexible and the list of understood keywords can vary with time without breaking anything See the comments at the top of the esdf f90 module for more information e Each line is of the form keyword value There must be a space either side of the colon e The parameters are divided into types string integer single double physical boolean Vari ables of physical type must be supplied with a unit such as hartree eV rydberg or Joules Megaparsec All reasonable physical units are understood e The parameter names are case insensitive e g EnerGyCuToFF is equivalent to energycutoff and punctuation insensitive energycutoff is equivalent to energy_cutoff is equiva lent to energy cutoff Punctuation chara
30. functions centered on the origin with non zero angular momentum whose spherical averages are zero and the tails of the Gaussian basis functions centered on other sites which must be cusp less at the nucleus in question We therefore obtain do x a Z 90 n 0 70 We use Eq 70 as the basis of our scheme for constructing cusp corrected orbitals We have found that it is important to include 0 when correcting molecular orbitals expanded in Gaussian basis functions We suspect this will also be so for Slater type functions and our cusp correction scheme could also be applied in that case to impose the cusp conditions exactly 10 13 2 Cusp correction algorithm One could conceive of correcting the orbitals either by adding a function to the Gaussian orbital inside some reasonably small radius multiplying by a function e g using the Jastrow factor as mentioned earlier or by replacing the orbital near the nucleus by a function which obeys the cusp condition However as the local energy obtained from Gaussian orbitals shows wild oscillations close to the nucleus the best option seems to be the latter one replacement of the orbital inside some small radius by a well behaved form We apply a cusp correction to each orbital at each nucleus at which it is non zero Inside some cusp correction radius re we replace the part of the orbital arising from s type Gaussian functions centered on the nucleus in question by C
31. is normally taken to be a VMC Slater Jastrow wave function with the parameters in the Jastrow factor optimized with a variance minimization procedure 14 So we begin by assuming we have an input file an xwfn data file containing the determinantal part of the wave function and an optimized jastrow data file containing the Jastrow factor A DMC calculation then consists of three basic steps e VMC config generation e DMC equilibration e DMC statistics accumulation The wave function in DMC is not represented in terms of a basis set but by the time dependent distribution in configuration space of a set of configurations or walkers The shape of the many electron wave function in configuration space is built up by killing or duplicating individual walkers according to certain rules and thus the population of walkers fluctuates The first step in DMC then is to generate these configurations in their initial distribution i e accord ing to the VMC trial wave function This is done through a VMC calculation in exactly the same way as we generated configurations for variance minimization in the previous section Therefore one should again set nwrcon to be the desired number of configs node and nmove to be gt nwrcon Setting the correlation period corper to a higher value than in VMC is less important than in variance minimization as the DMC configs are time dependent Note this is a common point of confusion that the irun flag
32. liable to lead to parallel inefficiency however especially on small numbers of nodes 26 NEIGHPRINT Integer NEIGHPRINT n will generate a printout of the first n stars of neigh bours of each atom in the primitive cell with the relevant interatomic distances given in both Angstrom and a u If n 0 or if you are an atom free electron or electron hole fluid phase then the keyword has no effect NEQUIL Integer Number of Metropolis equilibration steps Note that the correlation period is not accounted for i e NEQUIL configuration move attempts are made NEWOPT_METHOD Text NEWOPT_METHOD specifies the method used to minimize the quartic LSF NEWOPT_METHOD should be one of CG conjugate gradients MC Monte Carlo LM line minimization SD steepest descents BFGS Broyden Fletcher Goldfarb Shanno BFGS_MC BFGS and Monte Carlo CG_MC conjugate gradients and Monte Carlo or GN Gauss Newton See section 10 24 NEU NED ntegers Number of up spin and down spin electrons NHU NHD Integers Number of up spin and down spin holes NEWOPT ITERATIONS Integer NEWOPTITERATIONS specifies the maximum num ber of conjugate gradients steepest descent or BFGS iterations to be performed if NEWOPT_METHOD is CG SD BFGS CG_MC or BFGS_MC If NEWOPT_METHOD is MC LM CG_MC or BFGS_MC then it specifies the number of line minimizations to be perform
33. list of supported architectures with definitions and library requirements The instructions below require you to set environment variables The procedure for doing this depends on what Unix shell you use Roughly speaking if you use the t csh shell then environment variables should be set in your cshrc file using commands like setenv QMC_ARCH sun setenv QMC_TMPDIR temp user set path path HOME CASINO_bin_qmc utils QMC_ARCH If you use the bash shell then add something like the following to your bashrc file export QMC_ARCH sun export QMC_TMPDIR temp USER export PATH PATH HOME CASINO bin_qmc utils QMC_ARCH After modifying the appropriate file you should type source cshre or source bashrc Make sure these files are read automatically during the login procedure that you use if not include the source line above in your login or bash_profile file If you use a shell other than these three then you re on your own Set the environment variable QMC_ARCH to be machine_type where machine_type is one of the supported architectures listed in the CASINO ARCH file This variable identifies what kind of computer you are running on so that CASINO knows what Fortran compiler to use what libraries are required how jobs are run and so on If you are not sure which one you should use take a look at the list in this file and or the files in the CASINO src zmakes directory 2 Add CASINO bin_qmc ut
34. look at its contents then use the format _spline utility to generate a formatted swfn data formatted file This is also useful if you want to transfer a swfn data file from one platform to 17 4 Carrying out a QMC calculation Now finally we can carry out QMC calculations using the localized orbitals For example type runumc to carry out a HF VMC calculation The usual CASINO inputs are required along with the swfn data file 6 3 How to run the code run scripts As you may have already gathered although CASINO can be run by sitting in a directory containing input files and typing casino gt amp out it is recommended that one run it using the shell scripts provided in the CASINO scripts directory There are three basic scripts meant to control VMC DMC and variance minimization calculations Predictably their names are runumc rundmc and runvarmin These may be different for different machine types but set up should be automatic following definition of the QMC_ARCH environment variable i e the relevant scripts for your machine type should be copied into your path on typing make in CASINO utils The scripts are designed to reduce the effort of doing QMC calculations to just one command They are most useful when using parallel machines with a batch queueing system They are written to detect all common errors that a user may make in setting up a calculation so that the fault is detected immediately rather than when the batch job sta
35. number of bins in the heg data file All data about g r for all particle and spin combinations is then output at the end of the simulation in the vmc corr and dmc corr files The utility plot_corr must be used in order to obtain g r for the desired pair of particle types and to apply correct normalization to the data 10 27 Relativistic corrections to atomic energies Relativistic corrections to the ground state energies of closed shell atoms can be calculated using first order perturbation theory which gives results accurate to the order 1 c where c is the velocity of light This method works well for atoms of low nuclear charge Z when the relativistic corrections are small but is not satisfactory when Z is large In CASINO the relativistic corrections can be calculated for VMC method 1 and DMC by setting the relativistic flag in the input file to I By default the relativistic corrections are not calculated First we consider the mass polarization term 1 which accounts for the correction due to the finite mass of the nucleus to order 1 M where M is the nuclear mass in amu This term is given by 1 E1 yhen 206 where v R is the drift vector of electron i v R Y R V R By default CASINO uses nuclear masses averaged over isotopes which are listed in Table 2 If a nuclear mass for a specific isotope is required the default setting can be overwritten by the isotope mass keyword in the input file The relativistic ter
36. of Bloch orbitals for a periodic system There exist efficient algorithms for computing the transformation from Bloch orbitals to so called maximally localized Wannier functions 40 41 The transformation from Bloch to Wannier functions is unitary so the orthogonality of the orbitals is preserved Wannier orbitals are spatially localized so they can be truncated to zero when sufficiently far from their centres Therefore when an electron is moved only a few orbitals need to be evaluated the others are zero because the electron is outside their truncation radii So the number of orbitals to be computed is O 1 Although the Wannier functions are localized they are not rigorously zero at their truncation radii Hence their truncation results in small discontinuities in the Slater determinant These are potentially serious for QMC because they result in the presence of Dirac delta functions in the kinetic energy integrand The delta functions cannot be sampled so their contribution to the total energy is lost However in practice the resulting bias is extremely small provided the truncation radii are sufficiently large The discontinuities can be avoided if the orbitals are truncated smoothly over spherical shells by multiplying them by polynomials that are equal to one at inner truncation radii and zero at outer truncation radii 14 1 4 Localized bases Suppose the basis functions are zero outside fixed radii about their centres
37. one of various schemes to evaluate corresponding limited vector Use limited drift velocities to calculate Metropolis acceptance probability for proposed electron move Decide whether or not move is accepted Evaluate components of limited electron energy at whichever position old or new that it ends up at End loop over electrons in configuration Evaluate local energy of the new configuration and the limited energy Calculate branching factor for new configuration and so determine the number of configurations that proceed to next generation at point in phase space occupied by this configuration Use the limited energy for evaluating the branching factor End loop over configurations Update best estimate of energy Update reference energy End loop over moves in a block End loop over blocks of time steps 61 10 Theoretical background 10 1 The trial wave function Throughout this section we assume the trial wave function is of the Slater Jastrow form vr R exp J R D r dea TRI a wes tN 1 where R ri ry is a point in the configuration space of the N electron system J is the Jastrow function and D and D are Slater determinants of orbitals of spin up and spin down electrons respectively See reference 23 for further information 10 2 The variational Monte Carlo method 10 2 1 Evaluating expectation values The expectation value of the Hamiltonian H with respect to the trial wave function wr can be writ
38. out the number of electrons to which neu and ned have to be set 8 2 4 Phases Why phases We know that for many but not all orbitals we can extract two real orbitals the real and the imaginary part of a complex orbital that we can use to build a determinant In some cases we might only use one of the two We could use the real or the imaginary part however we could also use some linear combination of the two It is possible to achieve this by using a multiple determinant calculation one involving the real part and one involving the imaginary part in conjunction with an appropriate weighting of the determinants Since one determinant is more efficient than several CASINO allows the specification of phases to choose the linear combination of real and imaginary part see Section 10 12 GSEPC and DSPEC modify the phases for a whole determinant ASINGLE modifies the phase of a single orbital ASINGLE is followed by a phase between zero and 27 followed by a spin a determinant a k point and a band index i e ASINGLE 0 5 2 3 2 2 sets the phase to 0 5 for the orbital at k point 2 band 2 in the 3rd determinant for the down spin electrons This is independent of whether this orbital is actually used The phases are speci fied in multiples of 7 The units for the phases can be altered using UNITA followed by a real number UNITA 1 0 Now setting a phase to ir as in the example above has to be down using ASINGLE 1 570796327 2 3 2 2 Instead of s
39. out via a correlated sampling approach in which a set of configura tions distributed according to W ag is generated where ao is an initial set of parameter values The variance 07 a is then evaluated as 186 _J V ao w a Ej a Ey a dR 2 1 opla J 2 a0 w a dR 87 where fS Vao wla Er a dR E 1 vio Jla w a dR G and the integrals contain a weighting factor w a given by V a The parameters a are adjusted until o a is minimized The advantage of the correlated sampling approach is that one does not have to generate a new set of configurations every time the parameter values are changed In order to generate a set of configurations distributed according to W ao a VMC configura tion generation run must be carried out first The subsequent variance minimization using these configurations is handled by the subroutine VARMIN The minimization itself is carried out by the routine NL2SNO in module NL2SOL which performs an unconstrained minimization without re quiring derivatives of a sum of m squares of functions which contain n variables where m gt n Information on the minimization routine can be found in Reference 39 Having generated a new set of parameters we then carry out a configuration generation run with these parameters and if necessary a further variance minimization run and so forth Before carrying out this process the user must decide how they wish to parameteri
40. output files Brings each of the CAS configurations into maximum co incidence with the reference configuration and optionally calls resum_cas Outputs the CAS wave function i e the determinant ex pansion to the gwfn data file MODULE holds data defining the CIS and other multi determinant wave functions Multiplies the common part of shell normalization factors into the contraction coefficients and adjusts their storage for improved accessibility Echoes a string and then kills the program MODULE holds the data about the type of Gaussian run as well as the data defining the MOs etc Main driver unit Evaluates a primitive d type Gaussian basis function at a specified location in 3D space Evaluates a primitive s type Gaussian basis function at a specified location in 3D space Reads the Gaussian output file and identifies whether it is from Gaussian 94 or Gaussian 98 MODULE holds parameters defining granularity of plot ting and integration grids as well as which MO to plot test Brings a CIS configuration into maximum coincidence with a specified reference determinant Tests the normalization of a specified MO Normalizes the CIS expansion Not necessary for QMC but keeps things tidy and output gives an idea of how complete the expansion is Sort an array into descending numeric order and option ally keep track of reordering As for numsrt but has additional argument to specify as cending or descending order As f
41. runvarmin n 3 16 runvarmin n 3 16 night Run 3 iterations of a config generation variance minimization cycle on 16 nodes Again some machines might require you to specify test short medium long etc Note if you make any changes to the scripts in the appropriate CASINO scripts SQMC_ARCH directory typing make in the utils directory will copy the modified scripts to the appropriate bin_qme directory which should be in your path There is also a runvarmin1 script on some machines which is designed for machines where you can waste weeks sitting around in batch queues such as Cambridge HPCF Origin 2000 machines Basically what it does is to generate complicated batch scripts which will do lots of separate CASINO jobs within a single allocated time slot in some queue See the top of the script for more information Note that ultimately the scripts for the different machine types wil be merged into two one for batch queuing systems one for interactive use but as this is deeply tedious I keep putting it off 20 7 Files used by CASINO CASINO uses a variety of files during the course of its operation to store things in Here is a list of them Input Files input Parameters controlling the QMC calculation See Section 7 1 xx_pp data Pseudopotential if required xx element symbol in lower case e g si c See Section 7 3 jastrow data Jastrow factor if required See Section 7 2 xwfn data x g p b a Trial wave function in Gaussian
42. the drawing method the bond resolution to be 1 and the sphere resolution to be 15 Click Apply On VMD Main click on Extensions gt vmdmovie For version 1 8 3 click on Extensions Visualisation Movie Maker A VMD Movie Generator will pop up e In Movie Settings choose Trajectory The movie can be saved in the AVI or MPEG format Choose by clicking on Format and tick the preferred format Then check whether the the name of the temporary directory suggested is right this is where the RGB files are created Note that this directory should be free of RGB files belonging to other users If this has to be changed then click on the Set working directory button and browse for the directory Type in the name of the movie in the box provided Click on the Make Movie button e The movie will be displayed in the Open GL Display screen The mpg or avi movie file will be produced in the working directory being specified They have to be viewed with other viewers for example mpeg_play for mpg files For an example movie made with VMD a CASINO VMC simulation of cyclohexane see www tcm phy cam ac uk mdt26 downloads cyclohexane2 mpg 11 2 2 Jmol Jmol is a free open source molecule viewer It supports computers running Windows Mac OS X and Linux Unix systems Jmol can be downloaded from the following website http jmol sourceforge net e Type jmol Click on File Open Browse for the file movie out and click
43. the energy is not calculated and can thus be implemented into the equilibration process A second remark is that DTVMC can be regarded as a length it measures the width of the probability distribution that is sampled to move the particles at every step and as such is connected with a characteristic length of the system which can be influenced by its size density etc For example for low densities it is found that the optimized DTVMC is much larger than those that one gets used to seeing what happens is that the distribution flattens into a constant or almost so so that the electrons can move to a random place in the simulation box not only within a small environment of their previous positions hence better significance decorrelation in the sampling of the desired observables is achieved The last comment is that the procedure of optimization has not been proven infallible for all initial values or extreme conditions such as systems in which the small number of particles might prevent gathering of good statistics for the interpolation Nor is it the case that there is an analytic connection of the acceptance ratio being 50 with the minimization of the correlation time However test cases show that the technique works and is useful specially when tackling a new system for the first time 10 3 The diffusion Monte Carlo method See e g references 22 and 23 for general information about DMC Here we present some more detailed inf
44. the reblock utility The correlation time measures the average number of Monte Carlo steps between two uncorrelated values of the energy and should be unity for optimal statistics Note that by Monte Carlo Step what is meant is every step for which an energy is stored in this context For example in a VMC calculation every CORPERxNVMCAVE steps are merged into one energy Increasing any of the two variables would decrease the calculated correlation time meaning that the values that have been stored are more uncorrelated The definition of the correlation time 7 of an observable H is ica A t dt i Hv Hin ieee MOVE y 175 00 On T where A t is the value of the autocorrelation function at an interval of t 0 Hy Heye e is the variance of the expectation values and the latter are taken with respect to their subscript which 92 we shall remove for the sake of clarity For a discrete set of values equally spaced by an amount At 1 2 On T A09 1 2Y A0142 eS BOTs A 176 t oco t and for a finite set of length N N 1 pi 2 y He NN Heys H ns t 1 CH Numerically the problem with this expression is that if averages are used instead of proper ex pectation values which are of course unknown great fluctuations will appear at the tail of the autocorrelation function This problem is solved by introducing a cut off in the summation T tmax 14 255 Hy FO Te
45. the user about excited states and resumming Plot an MO debugging option Called by wfn_test Asks user about plotting and normalization testing Op tionally calls wfn_construct and normalization_check 111 14 Use of localized orbitals and bases in CASINO 14 1 Theoretical background 14 1 1 Introduction The rate determining step in practical QMC calculations is the evaluation of the orbitals in the Slater part of the trial wavefunction after electron moves We explain how the CPU time for evaluating the orbitals can be made essentially independent of system size This makes it feasible to simulate systems of hundreds of electrons with sufficient accuracy to resolve their optical gaps 14 1 2 Standard QMC Let us assume that in standard QMC the orbitals are HF or DFT eigenfunctions extending over the entire system and that the basis functions also extend over the entire system as is the case with a plane wave basis or even with a Gaussian basis if the Gaussians are not truncated Let N be the system size number of electrons In standard QMC after each electron move we must update O N orbitals expanded in O N basis functions Hence the time taken for a configuration move of all the electrons scales as O N3 14 1 3 Localized orbitals It is easy to show that a nonsingular linear transformation of a set of orbitals can only change the normalization of the Slater determinant of those orbitals Consider a set
46. to derive the relevant formula is to take the infinite separation limit of the 3D sum for a periodic stack of finite width slabs CASINO uses this algorithm when treating two dimensional slabs of atoms with local Gaussian basis sets useful in modelling surfaces and also when treating 2D electron and electron hole phases either as strict 2D planes 2D slabs with finite thickness or strict 2D bilayers In the case of periodic arrays of slabs separated by a finite vacuum gap necessary when using plane wave basis sets then the regular 3D algorithm is used Consider a finite width slab with a charge density periodic in two dimensions consisting of a unit point charge at rj in every simulation cell plus a uniform cancelling background 1 pitt Y sr 8 3 149 R where R now denotes the 2D lattice translation vectors in the ry plane and A is the area of the simulation cell in that plane The Parry formula for the periodic potential corresponding to this charge density is adh erfe 72 r rj R Wen Le eR 4 expire eR xp 2G eric G 3 ex 2G er ss 2 Ade a E p 2G erf 5 p 2G erfe 5 o T Eea exp y2 a fe f zy a a 150 88 where G denotes the reciprocal lattice translation vectors in the ry plane and z is the z component of the r r vector In two dimensions CASINO sets the screening parameter y to 2 4 A2 which again should approximately minimize the cost over different Bravais lattices The fu
47. to zero The atomic mass unit amu in this sense means the ratio of the average mass per atom of the element to 1 12 of the mass of 12C ISPERIODIC Logical TRUE if and only if system is periodic in 1 2 or 3 dimensions ITERAC Integer Determines type of electron electron interaction 1 Ewald interaction see Section 10 20 or 1 r for finite systems 2 MPC interaction see Section 10 20 4 3 Ewald and MPC Ewald used in DMC propagation 4 Ewald and MPC MPC used in DMC propagation JASBUF Logical If JASBUF is TRUE then the chi term in the Jastrow function for each electron in each configuration is buffered saves time at the expense of memory Clearly this will have no effect in systems without one body terms in the Jastrow such as the HEG JASTROW_PLOT Block This utility allows you to plot u ri x r and f ri rj rij terms in the new Jastrow factor The first line is a flag for whether the Jastrow factor is to be plotted 0 NO 1 YES the second line holds the spin of particle i the third line holds the spin of particle 7 the fourth line holds the x y z position of particle j for plots of f the fifth line holds a vector with the direction in which is moved for plots of f and the sixth line holds the position vector of a point on the straight line along which electron moves for plots of f Note that the nucleus is assumed to lie at the origin If f is plotted the jastrow_value_f_ out files contain the
48. want to diffuse in this fashion when the electron is likely to cross the nucleus Let II be the plane with normal e that contains the nucleus For the usual Gaussian diffusion process the probability that an electron drifts assuming that nuclear overshoot is permitted and diffuses across II is 1 Z ULT q 1 p erfc 47 q P 5e B 47 So with probability p we sample w from 3 2 w gi w 277 exp Tr J 48 and set the new electron position to be r r w otherwise we sample w from CF ga w E exp 2 w1 49 and set r Rz w We have defined by 1 4 22 50 T which reduces to Z for large timesteps giving the desired cusp however this choice of causes the second moments of g and ga to be equal to O r Hence the Green s function remains correct to O r The single electron Green s function for the move from r to r is given by gir r pgi r r Ggo r R2 51 In order to calculate the Green s function for the reverse move need to perform all of the steps above apart from the random diffusion starting at point r and ending up at r 10 7 4 Using the modifications in CASINO These three modifications to the DMC Green s function are applied if the NUCLEUS_GF_MODS keyword is set to T Note that they can only be used if bare nuclei are actually present 10 8 Evaluating expectation values of observables The reference energy Er is varied to
49. zero variance because each time an electron energy is evaluated we will if the electron move is accepted change the energies that the preceding electrons should have in the resulting configuration In summary better trial wave functions favour the standard Metropolis method In addition if the electron energies are very expensive to evaluate then it may be preferable to evaluate them once per electron at the end of the configuration move rather than twice per electron as required for use with Equation 4 We refer to the standard Metropolis procedure as method 1 and the method in which quan tities are calculated during the sweep as method 2 Both methods are available in CASINO the VMC_METHOD input parameter selects which is to be used The local energy does not have to be calculated every configuration move In particular energies are not calculated during the first NEQUIL moves of a VMC simulation when the Metropolis algorithm has yet to reach its equilibrium Furthermore because the configurations are serially correlated the expense of calculating the energy at every configuration move is unjustified it is more efficient to calculate the energy once every CORPER moves where typically the input parameter CORPER might be 4 or 5 Note that it is not necessary to write out every local energy calculated NVMCAVE successive energies can be averaged over before being written out to vmc hist this helps ensure that the vmc hist file is not exces
50. 0 3 GSPEC IP 0 000 1 endblock wavefunction GSPEC sets the default values for all determinants If it is omitted CASINO uses built in defaults DPSEC can then be used instead of GSPEC if the user wishes to set flags only for a specific determinant The syntax is equivalent to GSPEC except for two integers in the place of the two letter token These are ISPIN and IDET ISPIN may be 1 or 2 referring to up and down spin determinants and IDET can take an integer value from one to the maximum set by MD i e DSPEC 120 0001 sets the phase DIR SPR and SPRBND for the up electron sub determinant in determinant number 2 In addition the value for POL can be specified by using the line POL token ALL or POL token number where token can be chosen from IP Pl P2 and UP Using ALL results in POL being set for all determinants Using a number instead means POL will only be set for the determinant given by the number 58 The occupation of orbitals can be further modified by using the the specifier EPOL EPOL 0 32800768 ALL or POL 0 32800768 number EPOL temporarily lowers the energy eigenvalues of the down electron Kohn Sham eigenvalues by an amount given by the real number that follows the token In effect it biases the occupation scheme and can be used to generate strongly spin polarized systems If the use of EPOL or any other flag leads to the number of spin up and spin down electrons not being equal to neu and ned CASINO stops and prints
51. 00000 0 306218280918 0 306218280918 0 306218280918 0 306218280918 0 306218280918 0 306218280918 2 143527966427 2 755964528263 1 531091404591 1 531091404591 3 368401090099 0 918654842754 Blip grid 20 20 20 46 WAVE Numbe 1 FUNCTION r of k points k point of bands up spin down spin k point coords au 32 O 0 0000000000000000 0 0000000000000000 0 0000000000000000 1 Band spin eigenvalue au 1 Blip Band spin eigenvalue au 2 Blip 7 6 4 1 0 015443213499 coefficients 0 347534087942 0 517756216863 0 528529751776 0 450993755155 0 024823593723 0 022992059309 1 0 009057759900 coefficients 0 097602672632 0 176031268912 0 188226226118 for 32 bands swfn data file specification Orbitals expanded in spline functions Soon 7 6 5 awfn data file specification Atomic orbitals given on a radial grid Below is an example of an awfn data file for beryllium Atomic Be wave function in real space Atomi 4 Total 2 c number number of orbitals The 1s 2 2s 2 1S state electronic configuration Number of up down spin electrons 2 2 State NENE NENE S oooo ooo 0 Radial grid a u 30 DO OQ OG OOOO 1 000000000000000E 00 457890972218354E 02 477372816593466E 02 497683553179352E 02 518858448775392E 02 540934270674177E 02 563949350502860E 02 label of spin up electron quantum number n 1 m label of spin dow
52. 14136 etc for all 8 k points with 4 bands per k point in this case OTHER STUFF Anything you like can be written here intended for copy of input output files from the PW DFT program 7 6 3 bwfn data file specification Orbitals expanded in blip functions An example of a bwfn data file for bulk silicon is given below bulk Si BASIC INFO 45 Generated by PWSCF Method DFT DFT Functional unknown Pseudopotential unknown Plane wave cutoff au 7 5 Spin polarized F Total energy au per primitive cell 62 76695352625196 Kinetic energy au per primitive cell 0 E 0 Local potential energy au per primitive cell 0 E 0 Non local potential energy au per primitive cell 0 E 0 Electron electron energy au per primitive cell 0 E 0 Ion ion energy au per primitive cell 67 208144397949283 Number of electrons per primitive cell 64 GEOMETRY Number of atoms per primitive cell 16 Atomic number and position of the atoms au 14 1 282415538750 1 282415538750 1 282415538750 14 1 282415538750 3 847246616250 3 847246616250 14 11 541739848750 6 412077693750 6 412077693750 14 11 541739848750 11 541739848750 11 541739848750 Primitive lattice vectors au 10 259324310000 10 259324310000 0 000000000000 0 000000000000 10 259324310000 10 259324310000 10 259324310000 0 000000000000 10 259324310000 G VECTORS Number of G vectors 2085 Gx Gy Gz au 0 000000000000 0 000000000000 0 0000000
53. 155389173550 1 2824155389173550 14 1 2824155389173550 1 2824155389173550 1 2824155389173550 Primitive lattice vectors au 0 000000000000000E 000 5 12966215566942 5 12966215566942 5 12966215566942 0 000000000000000E 000 5 12966215566942 5 12966215566942 5 12966215566942 0 000000000000000E 000 G VECTORS Number of G vectors 341 Gx Gy Gz au 0 000000000000000E 000 0 000000000000000E 000 0 000000000000000E 000 lt snip gt 2 44974624702536 2 44974624702536 2 44974624702536 44 WAVE FUNCTION Number of k points 8 k point of bands up spin down spin k point coords au 1400 0 0 0 0 0 Band spin eigenvalue au 1 1 2 432668987064920E 002 Eigenvector coeficients 0 953978956601082 lt snip gt 0 000000000000000E 000 Band spin eigenvalue au 2 1 0 417639684278199 Eigenvector coeficients 1 576227515378431E 012 lt snip gt 0 000000000000000E 000 Band spin eigenvalue au 3 1 0 417639711691157 Eigenvector coeficients 6 348405363267598E 008 lt snip gt 0 000000000000000E 000 Band spin eigenvalue au 4 1 0 417639745897256 Eigenvector coeficients 2 030042411030856E 009 lt snip gt 0 000000000000000E 000 k point of bands up spin down spin k point coords au 2 4 0 0 3062182808781698 0 3062182808781698 0 3062182808781698 Band spin eigenvalue au 1 1 6 261478573034787E 002 Eigenvector coeficients 0 683590846620266 lt snip gt 0 000000000000000E 000 Band spin eigenvalue au 2 1 0 1568569314
54. 334347 2 79252680319093 0 005315599024948 The lines belong at the end of the file just before END JASTROW The Jastrow factor will then include a term J 5 C k sin kz k 122 20 Appendix 2 Programming Guide for CASINO You are not allowed to make modifications to CASINO without explicit written permission from the Cambridge group Should you be in such a position please read the following 20 1 STYLE CASINO has a FORTRANO9O style which should be adhered to when writing code either for the main source or for utilities This is because we think it is desirable that the package has a relatively homogeneous look and feel Everybody has their own style Yours is different and may even be better but we ve decided on one for CASINO and here it is If you don t write your code like this the likelihood is that somebody else will reformat it for you and they will probably accidentally delete a crucial minus sign while correcting your routine an error which may take two weeks of your life which could be much better spent doing other things to track down You can get most of this just by looking at the code but let s emphasize the main points e Don t use more than 1 module per physical file as Hitachi compilers won t allow it e Capitalization outside of character context Use upper case letters for keywords whose use is primarily at compile time statements that delimit program and subprogram boundaries declaration
55. 40 41 Gaussian 98 Revision A 7 M J Frisch G W Trucks H B Schlegel G E Scuseria M A Robb J R Cheeseman V G Zakrzewski J A Montgomery R E Stratmann J C Burant S Dapprich J M Millam A D Daniels K N Kudin M C Strain O Farkas J Tomasi V Barone M Cossi R Cammi B Mennucci C Pomelli C Adamo S Clifford J Ochterski G A Petersson P Y Ayala Q Cui K Morokuma D K Malick A D Rabuck K Raghavachari J B Foresman J Cioslowski J V Ortiz B B Stefanov G Liu A Liashenko P Piskorz I Komaromi R Gomperts R L Martin D J Fox T Keith M A Al Laham C Y Peng A Nanayakkara C Gonzalez M Challacombe P M W Gill B G Johnson W Chen M W Wong J L Andres M Head Gordon E S Replogle and J A Pople Gaussian Inc Pittsburgh PA 1998 4 Gaussian 03 Revision B 03 M J Frisch G W Trucks H B Schlegel G E Scuseria M A Robb J R Cheeseman J A Montgomery Jr T Vreven K N Kudin J C Burant J M Millam S S Iyengar J Tomasi V Barone B Mennucci M Cossi G Scalmani N Rega G A Petersson H Nakatsuji M Hada M Ehara K Toyota R Fukuda J Hasegawa M Ishida T Nakajima Y Honda O Kitao H Nakai M Klene X Li J E Knox H P Hratchian J B Cross C Adamo J Jaramillo R Gomperts R E Stratmann O Yazyev A J Austin R Cammi C Pomelli J W Ochterski P Y Ayala K Morokuma G A
56. AL_K NUM_AO NUM_AO NUM_COMPLEX_K NUM_AO NUM_AO block as follows spin spin polarized calcs only k for solids bands solids _or MOs molecules A0 basis functions grouped by shell Complex coefficients are given as 2 adjacent real numbers real imaginary Ordering of d orbital coefficients Z2 XZ yZ x2 y2 xy Ordering of f g h OTHER STUFF m 0 1 1424 24 3 354 ita Anything you like can be written here intended for copy of input output files from DFT HF etc program 43 7 6 2 pwfn data file specification Orbitals expanded in plane waves pwfn data file example should be self explanatory Si diamond BASIC INFO Generated by CASTEP Method DFT DFT Functional LDA Pseudopotential LDA Trouiller Martin 1551 coeff Plane wave cutoff au 7 5 Spin polarized F Total energy au per primitive cell 7 84635517057440 Kinetic energy au per primitive cell 3 24780615624099 Local potential energy au per primitive cell 1 03508540683458 Non local potential energy au per primitive cell O 1448987 14904598 Electron electron energy au per primitive cell 1 80295658630558 Ion ion energy au per primitive cell 8 40101804857984 Number of electrons per primitive cell 8 GEOMETRY Number of atoms per primitive cell 2 Atomic numbers and positions of atoms au 14 1 2824155389173550 1 2824
57. CASINO User s Guide Version 1 7 7 Richard Needs Mike Towler Neil Drummond and Paul Kent Theory of Condensed Matter Group Cavendish Laboratory Madingley Road Cambridge CB3 0HE UK January 2005 Contents 1 Introduction 2 The quantum Monte Carlo method 3 Miscellaneous issues dl ISHPPO sencilla dada 3 2 Legal agreements and being nice 3 3 Getting the latest version of the code o o iey t eee 4 Functionality of CASINO 5 Installation 51 Supported architechIter s netos a a a de aa 5 2 Unsupported architecture i sa a a eee aa a as e a A a 5 3 Further installation noted ocre pe 320303900 a nd 6 How to use CASINO 6 1 Getting started Summary a a e a roe a E a a a 6 2 Getting started i details o o soa a a a e AA 6 3 How to run the code trun scripts lt 6 es c ia tte dde Se ee de de e dd de a a 7 Files used by CASINO tal Basi mput ilet OPUS ee A e aa se A 2 Jastrow tactor file jastrowsdata ossos se so al dae a ad T gt Pseudopotential fle xxppsdate s conn oe u e RE ee ke OD we 7 4 Charge density file density data s es saa t i ia a a a a kaea eee ee To MPC file eepotdata ico sordas aate aa ad A a aa a 7 6 Trial wave function files awfn data bwfn data gwfn data pwfn data swfn data TE adata nea maa p a a 44 rra A ale a SE OSE a 7 8 Jellium slabs heg data file with ETYPE 6 TO ac Met ble a a RA A aa a y a
58. Dye a x 134 JAI where Rj R Rz and riz r Rz The CPP energy is then 1 Vopp 7 Dk F 135 where a is the dipole polarizability of core I Eq 134 assumes a classical description which is valid when the valence electrons are far from the core When a valence electron penetrates the core the classical result is a very poor approximation diverging at the nucleus To remove this unphysical behaviour each contribution to the electric field in Eq 134 is multiplied by a cutoff function f rir Fr which tends to unity at large rir A further possible modification is to allow the one electron term in Eq 135 which takes the form a7 2r1 to depend on the angular momentum component l so that 77 in the cutoff function is replaced by T r With these modifications the CPP energy operator becomes Do DLE ALD es E 2 i jAi Ere 2 ri Raz rir Ryr Aa a Er 0 i JAI JAI where is the projector onto the lth angular momentum component of the ith electron with respect to the Ith ion We use the cutoff function 31 32 f x 1 cry 137 For efficient evaluation Eq 136 is written as Vor 3 abia Do SO ay TU where pt pl GE JAI and the maximum angular momentum is mx 2 In our approach the cutoff parameter for all angular momenta l gt 2 is rz which is slightly different from Shirley and Martin 32 who use Fay 10 19 1 Implementation of CPPs in molecules and solids
59. File vmc out 3 Total energy 0 611320944695 0 001157863567 Time 20 9800 seconds We see that optimizing the parameters in this way lowers the error bar and total energy significantly and also that most of the improvement is in the first iteration One could try additional things such as optimizing the cutoff or using more parameters or having different functional forms between different spin types to optimize the Jastrow still further Feel free to try this Attempting to do this more carefully for me gave HF 0 644241947395 0 000416713194 VMC 0 601693183856 0 000110549234 DMC 0 599042511676 0 000185847226 where it can be seen that VMC gives around 94 of the correlation energy if we assume the DMC result is the correct answer If we were optimizing a real system with atoms then we would have included the atom centred electron nucleus x term one for each type of atom and possibly the electron electron nucleus f term again one for each atom type There are two types of additional term p and q but these are rarely used and then only for periodic systems 6 2 8 How to do a DMC calculation Diffusion Monte Carlo calculations are the main point of doing QMC They are generally ex tremely accurate comparable or better than the best quantum chemistry correlated wave function techniques and yet remain applicable to very large systems However they require an accu rate trial wave function to be efficient This
60. Go 10 16Evaluating the local energy ae ee giiia aa Re ee ee ed 10 1 7Evaluatine the kinetioenerfy sa a saad 6008 a OED ER AR a a 10 18Evaluating the non local pseudopotential energy o 10 19The core polarization potential energy e e 10 20Evaluation of infinite Coulomb sums e a 10 21Estimating equilibration times and correlation periods 10 22 Statistical analysis of data sa ao fea keh eS a a a a da A 10 23Wave function optimization standard method o 10 24Wave function optimization new method o e e e 10 25Further considerations in electron hole systems o 10 26Pair correlation function calculation for electron electron hole systems 10 27Relativistic corrections to atomic energies o o e 11 Making movies with CASINO Till Howto make Movies io 4d 44a aw e Edd HS 11 2 Visualisation 2 144 524 dera 44 88 oe aa a A 12 Using CASINO with the CRYSTAL program 121 The CRYSTAL program is ha ee eee he ba ad 12 2 Generating gwfin data files with CRYSTAL95 98 and crystaltogmc 12 3 Generating gwfn data files with CRYSTALO3 0 13 Using CASINO with the Gaussian94 98 03 programs 13 1 How to Use Geussiontogme lt di ed Oe beeen a Pe bd ee di de es 13 2 Other features Of ga ussiantogme a s s a a aoa dor a a ee eA 13 3 S
61. MODS is set to false otherwise a will be calculated as described in Section 10 7 In the UNR 5 scheme the modified local energy S a m is given by S a m Epest m Evest m Ez a mp 30 68 where V a m 1 Vw Note that we define S slightly differently from Ref 5 S is used only in the branching factor and when evaluating the average energy we sum the unlimited local energies Ez Finally a very simple option is simply to cut off the local energy drift vector when their magnitudes becomes large using the method of Depasquale et al 8 S R Ey sign 2 y7 E R Ey for EL R Ey gt 2 Y7 v R sign 1 7 vj for w gt 1 7 31 In the paper of Depasquale et al 8 they recommend using the variational energy for Ey but we use the best estimate of the energy We prefer the UNR scheme 10 6 3 Theoretical background to the UNR limiting scheme for the drift velocity Recall the approximation that the drift velocity remains constant over a configuration move Unfor tunately the drift velocity varies rapidly as it diverges near the nodes UNR suggest therefore that the Green s function for the drift process should be calculated on the assumption that Vyr rather than V Y 1V pr is constant over a configuration move when close to the nodal surface 5 Consider the drift of a single electron in a single configuration in the electron by electron DMC algorithm Let the electron have position vecto
62. Monte Carlo QMC electronic structure calculations for finite and periodic systems Its development was inspired by a Fortran77 development code known simply as the QMC code written in the early 1990s in Cambridge by Richard Needs and Guna Rajagopal assisted by many helpful discussions with Matthew Foulkes This was later extended by Andrew Williamson up to 1995 and then by Paul Kent and Mike Towler up to 1998 Various different versions of this were able to treat fcc solids single atoms and the homogeneous electron gas By the late 1990s it was clear that a modern general code capable of treating arbitrary systems e g at least atoms molecules polymers slabs crystals and electron phases was required not only for the use of the Cambridge QMC group but for public distribution At that time a user friendly general publically available code did not exist at least for periodic systems and it was felt to be a good thing to create one to allow people around the world to join in the fun So beginning in 1999 a new Fortran90 code CASINO was gradually developed in the group of Richard Needs largely by Mike Towler considerably assisted from 2002 by tall Ph D student Neil Drummond some routines from the old code were retained translated and reused although most were gradually replaced The main aims of the new code were generality speed ease of use and transferability over a wide range of computational hardware It is hoped that these o
63. Must be set to true at the moment default as newopt can only cope with a varmin_coeffs data file for the time being See section 10 24 USE_NEWJAS Logical If this input parameter is set to true then CASINO will use the new form of Jastrow factor introduced around Easter 2004 If it is set to false then the old Jastrow factor employed prior to that will be used The former requires a jastrow data file while the latter requires a jasfun data file The new Jastrow almost always gives a significantly lower VMC energy and variance than the old one USE_NEWOPT Logical If USE NEWOPT is T then the coefficients of the quartic LSF will be accumulated in VMC and the New Optimisation Method will be used in a varmin run Note that only linear Jastrow parameters can be optimized See section 10 24 USE JASTROW Logical Use a wave function of the Slater Jastrow form where the Jastrow factor exp J is an optimizable object that multiplies the determinant part in order to intoduce correlations in the system There are currently two forms of the Jastrow factor available in CASINO selectable via the flag USE_NEWJAS For historical reasons the default value for USE_JASTROW is T and is overriden when IRUN is 1 which is the old synonym for IRUN 2 VMC and USE_JASTROW F i e for a Hartree Fock VMC HF VMC calculation VMC_METHOD Integer VMC_METHOD selects which version of VMC to use 1 Evaluate configuration energies at the end of the configurati
64. NPCELL Block Vector of length three giving the number of primitive cells in each dimension that make up the simulation cell N B for the 1D periodic case NCELL 2 and NCELL 3 must be 1 for 2D slab case NCELL 3 must be 1 NUCLEUS GF MODS Logical This keyword is the switch for enabling the use of the modifica tions to the DMC Green s function for the presence of bare nuclei suggested in J Chem Phys 99 2865 1993 expected to reduce timestep errors in all electron calculations NVMCAVE Integer VMC only Instead of writing out local energies etc every time they are calculated i e every CORPERth configuration move we average over the NVMCAVE energies before writing to the vmc hist file This saves both space and time Note that the total number of moves of all electrons in a block is given by NMOVExCORPERxNVMCAVE so you will need to reduce NMOVE proportionately if you increase NVMCAVE this will be changed soon NWRCON Integer Number of configurations to be written out VMC or read in VARMIN 27 OPT_DTVMC Integer The purpose of OPT_DTVMC is to optimize the value of the time step DTVMC in VMC calculations in order to minimize serial correlation This is done during the Metropolis equilibration phase The keyword may take two possible values OPT_DTVMC 0 default turns off the optimization OPT DTVMC 1 means DTVMC is optimized in order to achieve an acceptance ratio of roughly 50 Note further that attempts have been
65. RM Spin dep 0 gt u d 1 gt u d 0 Number of primitive cell G vectors NB cannot have both G amp G 10 G vector in terms of rec latt vects label 0 t 2 2 1 1 0 2 10d 2 1 1 1 1 1 0 0 1 0O 0 1 1 o 1 0 1 1 1 2 2 2 1 1 2 1 2 1 2 Parameter value Optimizable 0 N0 1 YES END Q TERM END JASTROW Notes e There are 5 types of term available u isotropic electron electron terms x isotropic electron nucleus terms f isotropic electron electron nucleus terms p plane wave expansions in electron electron separation and q plane wave expansions in electron position All of these terms are optional e g omitting all the lines from START Q TERM to END Q TERM will give a Jastrow factor with no g terms In many cases particularly in the presence of pseudopotentials only u and x terms are needed e The truncation order should be either 2 or 3 If it is of value C then then the Jastrow factor is C times differentiable everywhere it must therefore be at least 2 for the kinetic energy integrand to be well defined If it is 3 then the local energy is continuous in configuration space assuming the orbitals are smooth This is not strictly required and leads to the loss of some variational freedom but it makes the numerical optimization of cutoff lengths easier e In a future release the u and x terms may be made anisotropic there are placeholders for this in jastrow data At present
66. ROMO o ae 10 X vel 5 ln 6 0 C p 0 Xs 77 Aa X5 01 X4 X X X X Xs XP A a ea 2 ae r2 2 as aR a XS pat y ghey 22 i r3 r2 To r2 r3 re O e 236 O UA e se 78 Our procedure is to solve Eq 78 using an initial value of 0 0 We then vary 0 so that the effective one electron local energy mo 5 jo tajo i Seno A m is well behaved Here the effective nuclear charge Zeg is given by n 0 Za Z 1 ey aa 80 which ensures that Ej 0 is finite when the cusp condition of Eq 76 is satisfied We studied the effective one electron local energies obtained using Eq 79 with Zef Z for the ls and 2s all electron Hartree Fock orbitals of neutral atoms calculated by numerical integration on fine radial grids for atoms up to Z 82 We noticed that the quantity Ej r Z is only weakly dependent on Z in the range r lt 1 5 Z We therefore chose an ideal effective one electron local energy curve given by Ejideal La Bo Bir Bar Bar Gar Bar fg6r Grr 81 The values chosen for the coefficients were 61 3 25819 Gg 15 0126 33 33 7308 G4 42 8705 Ps 31 2276 66 12 1316 G7 1 94692 obtained by fitting to the data for the 1s orbital of the carbon atom The value of 84 depends on the particular atom and its environment The ideal effective one electron local energy for a particular orbital is chosen to have the functional form of E
67. The move of the th electron of a configuration is accepted with probability min 1 e ae E meer eee E Dice o oe Ue iy E dg T min 1 exp r r Elvin tea Vel ora vifne otin rh ory Pudi tata thy o Eh P I A A rN r1 ri ri E rh Eip o y 19 This leads to the single electron detailed balance condition 1 2 f Silti Ei ri ri iyi SAS ri Ti rN Z 2 A A r SH ri Tipi o A AA Es 20 where A A A f E Silao Di ri ri q Py ee Piya o py 1 1 ie sd 21 is the effective single electron transition probability density once the accept reject step has been introduced Hence we find that the effective transition probability density for the entire configuration move satisfies 66 47441 N S R R isidro aa A EN i 1 N 2 Vello THT Ep Si T1 65 TET KH ri Tha r ve II i 1 PICAR AS 47 R And so detailed balance in configuration space is satisfied It is more efficient to use an electron by electron algorithm than the perhaps more straightforward configuration by configuration algorithm in which moves of entire configurations are proposed and then accepted or rejected This is because for a given timestep a configuration will travel further on average if the accept reject step is carried out for each electron in turn For example it is clear that it is very unlikely for a configuration not to be moved at all in an electron by electr
68. The slab width is the extent in the z direction of the region containing the positive background charge The z dependency of the wave function is described using Fourier sine series within a region determined by the wave function width outside this region it is assumed to be zero Sub bands are allowed and a set of Fourier coefficients must be supplied for each only one set of wavevectors is used The bands may be filled automatically by setting the Occupy by energy option Alternatively the band occupation for each spin may be specified directly The energy width of each band is given next It is used when occupying bands automatically and also as a consistency check The set of wave vectors k is given next followed by the set of Fourier coefficients C k for each band The final z dependence of the wave function is the same as within the Jastrow factor shown above QMC data file for jellium slab 2 00000000000000 Slab width wave function width 20 0000000000000 102 400000000000 Number of bands number of coefficients for Fourier sine series 29 29 25 25 21 21 21 21 13 13 9 1 1 Width of each subband 0 519117510000000 0 488868040000000 0 436931490000000 0 364359310000000 0 280391730000000 0 184220470000000 8 550147999999999E 002 Wavevectors in z direction 52 31 31 Coefficients o etc 0306796157577128 0613592315154256 3852469201402166 4159265358979294 0260839121091309 0000000000000000
69. Then the only functions that have to be evaluated when an electron is moved are those with the electron inside their radii So the number of basis functions to be computed simply depends on the local environment of the electron and is therefore independent of system size Gaussian basis functions can be regarded as localized if they are truncated to zero outside a certain radius This is done by default in CASINO Gaussian functions exp ar are assumed to be zero when exp ar lt 10 7 where Gr is the gautol input parameter Gaussian basis sets cannot yet be used in conjunction with localized orbitals however 15The time spent updating the cofactor matrices will in principle dominate in the limit of large system size but this limit is not reached in practice 16In CASINO the Gaussian basis functions are truncated 17The Slater part of the wavefunction is written as the product of Slater determinants for spin up and spin down electrons After the move of a spin up electron the number of orbitals to be evaluated is equal to the number of spin up electrons after the move of a spin down electron the number to be evaluated is the number of spin down electrons 112 On the other hand orbitals can be represented numerically using splines on a grid If the orbitals are localized then the memory requirements are greatly reduced because we only have to store the spline coefficients needed to evaluate the orbital within its truncation radius
70. Use of Slater Jastrow wave functions where the Slater part may consist of multiple determinants of spin orbitals Trial wave functions expanded in Gaussian basis sets s sp p d f or g functions centred on atoms or elsewhere produced by the following programs using DFT HF or various multideter minant methods CRYSTAL9X 03 24 25 GAUSSIAN9X 03 26 27 and TURBOMOLE Trial wave functions expanded in plane waves generated from PWSCF 30 ABINIT 29 GP CASTEP 28 and k207 Trial wave functions expanded in blip functions generated by post processing plane wave DFT solutions Trial wave functions in spline functions generated by post processing plane wave DFT solutions Numerical atomic calculations with the orbitals interpolated from a radial grid Improved scaling with system size through use of Wannier orbitals and localized basis functions Computation of excitation energies corresponding to either promotion or addition subtraction of electrons Computation of distribution functions such as the pair correlation function and density matrices electron and electron hole systems only for the moment 2D 3D electron phases in fluid or crystal wave functions with arbitrary cell shape spin polar ization density including excited state capability 2D layer separation between opposite spin electrons possible Pair correlation function 2D 3D electron hole phases with fluid crystal pairing wave functions with arbitrary cell sha
71. able since plane waves are the absolute worst basis set you can choose for QMC every basis function contributes at every point in space adding a factor of N to the scaling with system size Splines and blips are essentially the same thing written by two different groups However splines have now been flagged as obsolescent and users of this functionality should move over to using blip functions the splines will be removed in a future version of CASINO To make the CASINO calculation scale independently of system size energy per atom properties or quadratically energy per cell properties then the delocalized orbitals generated by most DFT HF programs need to undergo a unitary transformation to localized Wannier form This can be done for plane waves orbitals only at the moment either by using the CASINO zwannier utility splines or localizer utility blips The orbitals should then be re expanded in splines using the generate_spline utility or blips using the blipl utility Generating trial wave functions currently requires you to have access to one of the following codes This list may increase in the future contact Mike Towler if you want your code to be supported Gaussian basis sets CRYSTAL95 98 03 www tcm phy cam ac uk mdt26 crystal html and links from there GAUSSIAN94 98 03 www gaussian com TURBOMOLE www chemie uni karlsruhe de PC TheoChem turbomole intro en html Plane wave basis sets ABINIT www abinit o
72. ad For optimization Header not read S_2 to depend on spin Header not read B of S_2 beta_1 of S_2 beta_2 of S_2 beta_3 of S_2 beta_4 of S_2 Denotes end of set 1 Denotes end of file If the LO parameter is set to zero then it takes on a default value For periodic cases the default is 0 3Lwyys where Lws is the shortest distance to the boundary of the Wigner Seitz cell so that the Gaussian factor in Eq 96 is small 1 5x10 at rj Lws The default LO for a finite system atom or molecule is just a very large number 1 x 101 so the function is not cut off at all Note that any value given to the LC parameter other than zero will override the LO setting LC is a parameter outside of which the uy function of Eq 96 is set to zero If LC lt 0 then LC is interpreted as the value of uy below which it is set to zero if LC 0 no cutoff is applied and if LC gt 0 then the value given is interpreted as the cutoff radius LC gives a hard cutoff of uo whereas LO is a soft cutoff We normally choose LC to be zero and use LO to apply a soft cutoff If the L parameter is set to zero then it takes on a default value For periodic cases the default is Lws i e the shortest distance to the boundary of the Wigner Seitz cell For molecular systems the default is the molecular size currently defined as the maximum interatomic distance times 1 25 For atoms the default is 2 5 if the L parameter in an
73. ad an swfn data file and plot selected states plot_vmc_energy Simple program to read a vmc hist file and plot the total energy as a function of move number plot_vmc_hist Simple utility to plot any of the records from a vmc hist file as a function of move number pretty_printer mpr pretty printing script for source code listings which saves trees as suming a2ps has been set up correctly for your system ptm Manipulate pseudopotentials on radial grids or as Gaussian expansions quad fit Program for carrying out a quadratic fit to a set of data so as to find a local extremum Can be used e g for graphing DMC or VMC energies against a parameter on which they depend e g Gaussian exponent in orbitals of Wigner crystal quickblock Simple reblocking utility intended for analyzing very large dmc hist files since it only looks at one column at a time Alternative to reblock reblock Perform statistical analysis and reblocking of QMC data Use this rmstore Remove old vmc hist files produced with the big option from TMPDIR STORE on all machines in TCM Cambridge interest only simple_qmc_codes Some very simple QMC codes freely available on the web for teaching purposes trimdmc Trims dmc hist and dmc hist2 files to the end of the last complete block This is useful when restarting DMC runs that have been halted e g by running out of time on a parallel machine with a batch queue system ve Simple utili
74. al form of CASINO s new Jastrow factor CASINO s new Jastrow factor is the sum of homogeneous isotropic electron electron terms u isotropic electron nucleus terms x centred on the nuclei isotropic electron electron nucleus terms f also cen tred on the nuclei and in periodic systems plane wave expansions of electron electron separation and electron position p and q The form is N 1 N Nions N Nions N 1 N dt Ea DO Do L fi rir jr Tij i 1 j i 1 I 1 i 1 T 1 i l j i 1 N 1 N N T 5 5 tis Xale 83 i 1 j i 1 i l where N is the number of electrons Nions is the number of ions rij ri rj rir ri ry r is the position of electron and r is the position of nucleus J In periodic systems the electron electron 79 and electron ion separations rij and riz are evaluated under the minimum image convention Note that u x f p and q may also depend on the spins of electrons 7 and j The plane wave term p will describe similar sorts of correlation to the u term In periodic systems the u term must be cut off at a distance less than or equal to the Wigner Seitz radius of the simulation cell and therefore the u function includes electron pairs over less than three quarters of the simulation cell The p term adds variational freedom in the corners of the simulation cell which could be important in small cells The p term can also describe anisotropic correlations such as might be encountered in a layered comp
75. alue gradient and Laplacian of Jastrow factor are consistent Recommended to be set TRUE 30 VM_MODE Character May take values linear or direct see Section 10 14 linear for optimizing parameters which occur linearly in J direct for optimizing parameters which occur non linearly in J i e A VM_OPT_DET_COEFF Logical Optimize the coefficients of the determinants in variance mini mization VM_OPT_JASFUN Logical Optimize the Jastrow factor in variance minimization VM_OPT_PAIRING Logical Optimize the Rex parameter used in the electron hole excitonic insulator phase in variance minimization VM_REWEIGHT Logical If TRUE then use reweighting in VARMIN If FALSE set all weights to unity VM SMOOTH LIMITS Logical When set to T the optimizing routine used in variance mini mization is sent a smoothed version of the set of parameters This only affects those which are to remain bounded such as Jastrow cutoffs The result is a set of parameters which can vary in the range co 00 which can be more convenient than ignoring out of range values without the minimizer knowing A suitable hyperbolic function is used for mapping limited values into extended ones and viceversa VM_USE_E_GUESS Logical Setting this flag to true will cause the variance in varmin to be evaluated using VM_E_GUESS in place of the average energy of the config set in an attempt to combine energy minimisation wit
76. and we can give you some support in using it but please understand that this may well be rather limited in extent Having noted the above feel free to mail Richard Needs rn11 at cam ac uk or Mike Towler mdt26 at cam ac uk with any questions you may have For more general questions you can use our mailing list casino_users at phy cam ac uk 3 2 Legal agreements and being nice We are making CASINO available to a number of groups now We are not currently asking for any payment for academic use of the code but we ask users to sign an agreement concerning use of CASINO The practical upshot of this is that users may not redistribute the code they may not incorporate any part of it into any other program system nor may they modify it in any way whatsoever without prior agreement of the Cambridge group Please read this agreement carefully and abide by what it says You may not use or even retain a copy of CASINO if you have not signed this agreement A copy of the agreement can be found in the CASINO distribution in file CASINO docs consent_form ps New versions of CASINO are produced on a regular basis nightly builds for the beta version As a courtesy to the authors it would also be appreciated if users refrain from publishing scientific articles using recently introduced facilities in the code if the authors of those facilities have not yet published their work at least without first checking with them Thanks for your understanding on
77. and Frisch s book 34 getting CASSCF to converge for a singlet state is difficult The following procedure normally works 1 Run a ROHF calculation for the lowest triplet state of the system and save the checkpoint file 2 Run CASSCF for the second triplet state Nroot 2 taking the initial guess from the checkpoint file Gaussian calculates the first and second triplets but converges on the second 3 Run CASSCF for the first triplet Nroot 1 taking the initial guess from the checkpoint file 4 Run CASSCF for the singlet excited state Nroot 2 taking the initial guess from the triplet checkpoint file 5 Finally run CASSCF for the singlet ground state Nroot 1 By default Gaussian uses spin configurations combinations of Slater determinants in a CASSCF calculation It is best to converge the CASSCF state that you want using this option However for input to the QMC code the wave function must be in terms of Slater determinants In principle the SlaterDet option to CASSCF will do this but I never succeeded in getting it to work Instead specifying IOp 4 21 10 does the trick as does IOp 4 46 3 For large CASSCF calculations which in Silane equates to an active space of more than six electrons and eight orbitals the Cray T3E Turing at Manchester does not have enough memory per node and simply cannot be used For large CASSCF calculations on the alpha cluster Columbus the diagonalizatio
78. and N is the number of particles of type 91 similarly r is the position and N is the number of particles of type 02 the sums on and j run over all N and N particles respectively and Q denotes the volume of the simulation cell When 0 a2 the sum is replaced by Dus j r In practice g r is evaluated during the simulation by accumulating in radial bins that go up to the Wigner Seitz radius Lws The quantity that CASINO actually evaluates is Q Kn 200 N N V 200 90 02 Tn where rn is the radial position corresponding to the nth bin Vn is the volume of the nth bin and Kn is an average over the simulation of the number of pairs of particles whose distance falls in the nth bin In practice this means that at each step of the simulation N N3 records of distances are made Nj No 1 when 01 02 The radial position and volume of the nth bin in 3 dimensions is given by 47 3_ 3 2 n 5n n 1 nA 201 r a 201 3 2 1 Vn 47 A n n 3 202 Similarly in 2 dimensions 2_ 1 asi es 203 O 2 1 Vn 2TA n 2 204 The quantity A is the width of each bin given by A 4 5 where npin is the total number of bins Nbin From Eqs 199 and 200 it follows that the normalization of g r is such that 1 ane r 4rr dr gt X Vn Joro2 Tn Q N fa 205 01 n 99 To tell CASINO to evaluate g r set the pair correlation function flag to true and specify the total
79. and antiparallel components 5 625899453763526E 004 1 122777200374748E 003 3 487554343373239E 004 5 611467605982052E 004 9 941108935490427E 005 2 579395477618149E 004 1 175244704897280E 004 2 103573951529612E 004 2 542237859441641E 005 4 327076908591314E 005 No of sets of atom centred functions No of powers in ri rirj ririj 2 rirjrij 2 terms 1 SET 1 No of ions in set 1 Ions in set 1 500000 L a u 0 0 ri dependent terms Enable optimization T Enable spin polarization F Components i 8 773613686917481E 003 810651282199839E 003 410044958605923E 004 466197634637714E 004 119083682175356E 005 SET 1 JASTROW Notes 1 The file is read by SUBROUTINE READ _JASFUN 121 must be JASTROW not not Header read read Header Header read u_0 and cutoff read Header not A L_O of Header not A not for optimization Header not read L of S_1 Header not read Header not read No of powers in S_1 Header not read For optimization Header not read S_1 same for and anti Header not read B of spin 1 and spin alpha_1 of spin 1 and alpha_2 of spin 1 and alpha_3 of spin 1 and alpha_4 of spin 1 and Header not read NNNNN Number of sets of atoms A set has same S_1 S_5 Header not read Number of ions in set 1 Header not read Identifies ions in set 1 Header not read No of powers in S_2 S_5 Header not read L of S_2 S_5 Header not read Header not re
80. ansform TRUE Compute Wannier functions FALSE Print additional files If false cannot restart 1d 12 Tol of converged Wannier functions Default is 1d 12 1 The first line allows the user to specify the orbitals to be included in the Wannier transforma tion these are the first orbital in pwfn data up to the number given If one is interested in ground state calculations then all of the orbitals up to the highest occupied orbital should be transformed If on the other hand one is interested in calculating optical gaps then it is essen tial that both the highest occupied and lowest unoccupied orbitals are not transformed since the trial wavefunction for the first excited state is generated by replacing the highest occupied orbital with the lowest unoccupied one 2 The second line should always be set to true 3 It is possible to restart xwannier if the third line is set to true 4 The fourth line contains the localization tolerance The default of 1071 is usually adequate The output of xwannier is a pwfn data wannier file of the same format as pwfn data holding the Wannier orbitals expanded in a plane wave basis A second file called wannier_centers dat is also produced This holds the Cartesian coordinates of the Wannier centres 113 14 2 3 Representing orbitals with splines The input for generate_spline consists of a pwfn data file with the restrictions listed above an optional wannier_centers dat file and a
81. at the bth reblocking transformation is N M B Note that N is not necessarily an integer because the last block may be incomplete The standard error in the energy estimate at reblocking transformation number b is 184 and the error in 6 may be estimated as b 2 N 1 The reblock utility produces a table of 5 and e against b the user can then look for a region in which the standard error 6 has plateaued as a function of b and choose an appropriate reblocking transformation number which is then used to calculate the error bars on all the different components of energy 185 10 23 Wave function optimization standard method Optimization of the trial wave function is a crucial part of a VMC or DMC calculation CASINO allows optimization of the parameters in the Jastrow factor see Section 7 2 of the coefficients of the determinants of a multi determinant wave function of the parameters in the backflow n function of pairing parameters in electron hole gases and Pad coefficients in orbitals in Wigner crystals Wave function optimization in CASINO is achieved by minimizing the variance of the energy o a J Y a EL a Ev a dR JAR where a is a set of parameters and Ey is the variational energy The most important ground for preferring variance minimization to say energy minimization is that it shows better numerical stability particularly in large systems Minimization of 0 is carried
82. ation average of KE_i 2x Ti F2 Configuration average of F_i Configuration average of local part of potential energy Configuration average of short range part of MPC energy Configuration average of long range part of MPC energy Configuration average of non local energy Configuration average of local energy Configuration average of electron ion CPP energy Configuration average of electron CPP energy Configuration average of electron electron CPP energy Average weight of configuration for branching DMC Minimum weight of configuration for branching DMC Maximum weight of configuration for branching DMC Average age of the configurations number of time steps Maximum age of a configuration number of time steps 17 Effective timestep All energies are quoted in Hartrees per primitive cell per electron for electron phases per particle for electron hole phases The reblock utility should be used to analyze the dmc hist and dmc hist2 files Information about the different components of the energy and how they are calculated may be found in Sections 10 16 10 20 4 55 8 Specifying the Slater determinants 8 1 Gaussian basis sets The file gwfn data contains the information about the orbitals in the determinant s The gwfn data file contains the orbitals which will have been generated by GAUSSIAN94 98 03 or CRYSTAL95 98 03 and other information is added by the utilit
83. be set to T 2 There is no need to write out any configurations so the user can set the nwrcon keyword to zero 3 If desired the user may change the method used to minimize the LSF with respect to the set of parameters by using the newopt_method keyword This can take the values CG conjugate gradients SD steepest descents GN Gauss Newton MC Monte Carlo 96 line minimization LM simple line minimization CG_MC alternate conjugate gradi ents and Monte Carlo line minimization BFGS Broyden Fletcher Goldfarb Shanno or BFGS_MC alternate BFGS and Monte Carlo line minimization Tf the newopt_method keyword is omitted then the BFGS method will be used by default If you experience difficulty optimizing a large set of parameters then the Gauss Newton method is worth trying The BFGS method seems to be the most efficient method in general however 4 If desired the user can change both the maximum number of iterations and the number of line minimizations to be performed by means of the newopt_iterations keyword If this keyword is omitted or it is given a negative value then a default number of iterations will be performed The cutoff lengths in the Jastrow factor are important variational parameters and some attempt to optimize them should always be made It is recommended that a relatively cheap calculation using the standard variance minimization method should be carri
84. bjectives have been largely attained but the code continues to be actively developed Over the years additional contributions have been made by Andrew Porter Randy Hood Andrew Williamson Dario Alf Gavin Brown Chris Pickard Rene Gaudoin Ben Wood Zoltan Radnai Andrea Ma Pablo Lopez Rios Ryo Maezono John Trail and possibly others for which we are grateful Finally note that users need to sign the relevant piece of paper available from Mike Towler or Richard Needs before using the code and that the following citation quoted in full is required in any publication describing results obtained with CASINO R J Needs M D Towler N D Drummond and P R C Kent CASINO version 1 7 User Manual University of Cambridge Cambridge 2004 Further public information and resources are available at the CASINO web page www tcm phy cam ac uk mdt26 casino html j The CASINO pseudopotential library and download facility can be found on the private CASINO users page www tcm phy cam ac uk mdt26 casino_users html j 2 The quantum Monte Carlo method The correlated motion of electrons plays a crucial role in the aggregation of atoms into molecules and solids in electronic transport properties and in many other important physical phenomena Ab initio electronic structure calculations in which the properties of such correlated electron systems are computed from first principles are a vital tool in modern condensed matter p
85. c energy is discussed in Section 10 17 The evaluation of the non local energy is discussed in Section 10 18 while the core polarization potential energy is discussed in Section 10 19 The local part of the external potential energy is divided into a short range part around each ion which is evaluated directly and a long range Coulomb part which is evaluated using the Ewald potential in periodic systems see Section 10 20 or simply as a sum of 1 r potentials in finite systems The electron electron interaction energy is evaluated either using the Ewald interaction or the MPC interaction see Section 10 20 4 in periodic systems or simply as a sum of 1 r potentials in finite systems 83 10 17 Evaluating the kinetic energy The kinetic part of the local energy K can be expressed as a sum of contributions from each electron N N K Ki Y UR VYR 120 i 1 Because of the exponential form of the Jastrow factor it is convenient to re express K in terms of the logarithm of Y We define 1 1V2U 1 V Uy T V In Y y g f 121 pe 37 121 and the drift vector F 1 1 Viv F V In 122 Wp In W aU 122 Therefore K 2T F 123 In VMC an integration by parts shows that K IF 7 124 where the angle brackets denote averages over the variational distribution W R Eq 124 provides a useful consistency check for VMC calculations but note that it does not hold exactly
86. ccupying the up and down electrons separately until the neu lowest up electrons are occupied and the ned lowest down electrons too TP can be replaced by P1 P2 and UP In these cases CASINO first evaluates the total number of electrons in the system ntot neu ned It continues to occupy the lowest lying orbitals In the case of degeneracies it aims to P1 occupy as many up electrons as possible P2 occupy as many down electrons as possible UP align the number of occupied electrons as closely as possible with neu and ned The phase specified by GSPEC is the default value used instead of the built in default The three following flags are DIR If set to one the orbitals will be occupied from first to last used k point as they appear in pwfn data and reverse if set to zero SPR If set to one CASINO attempts to spread the occupied orbitals over several k points and over as few as possible if set to zero SPRBND Same as SPR except that now the focus is on different bands at the same k point These flags obviously only take effect if the highest occupied energy levels are not all occupied i e if there are degeneracies In that case GSPEC can be used efficiently to generate different occupation schemes without using any promotions see below 8 2 3 Multi determinant calculations In a multiple determinant calculation i e using MD instead of SD followed by the number of determinants and their relative weights block wavefunction MD
87. ce that the optimization of difficult parameters such as cutoff lengths works better when vm_reweight is TRUE It is possible to choose how much information CASINO will provide during the minimization process Setting vm_info to 1 will provide no information during the minimization setting it to 2 will provide a list of the parameters the mean energy and variance at each iteration 3 will provide the same information as 2 but the weights will also be computed and their mean and variance displayed note that the weights are not actually used unless vm_reweight is set 4 provides an enormous amount of detail and is only likely to be of use for development or debugging purposes Some things to check when carrying out variance minimization 95 e Does the VMC energy of successive configuration generation runs decrease by a statistically significant amount Once it has ceased to fall there is little point in continuing the minimization process Usually for a successful run we find that a single configuration generation variance minimization loop is sufficient but this should always be checked If the energy fails to settle down over successive loops then there may be too much variational freedom in the trial wave function e Does the reported variance at each iteration in the variance minimization runs fall If not then there is a problem with the minimizer NL2SOL itself e Does the reported mean energy at each iteratio
88. cell and noting that the Fourier transform of the Ewald interaction is Gaz we have He e gt fri r 9 fal pace Sm o aie opaa 0 C 163 i gt j i GX0 oG where Q2 2 Ar Q C gt gt ES fe PA aPA G gt fa 0Pa G 0PA G 0 164 G40 The calculation of fa is achieved using the following scheme developed by Randy Hood The integrand in Eq 161 diverges at the origin and we separate out the divergent behaviour by writing f r g r h r 165 where gr y r r lt L r gt L 166 and hr G vo r lt L 0 r gt L 167 where L is the radius of the largest sphere which is contained within the Wigner Seitz cell and y r is chosen to be y 3 y r 3 an so that both g and h have continuous first derivatives at r L The Fourier transform of h r is calculated analytically as r 3 iGrcos h tGTrcos G 5 af L os xr Qnr7e dcos6 dr 168 127 sinGL ae tape fos Po ea from which the G 0 value can be extracted as Qn L cas 170 5Q The Fourier transform of g is calculated numerically with an FFT added to the analytic Fourier transform of h and stored by CASINO in the file eepot data Note that for historical reasons the quantity stored in eepot data is actually O fa 10 21 Estimating equilibration times and correlation periods The rms distance diffused by a particle in a period T of imaginary time is V2Np DAT where A is the move acceptance ratio which is usually clo
89. ch orbitals may only be generated by post processing a pwfn data file a localized orbital facility for Gaussian basis sets will be added soon Additional limitations mean that currently only cuboidal simulation cells with a single k point k 0 for spin unpolarized systems may be treated A detailed discussion on the use of localized orbitals in QMC will appear in the Theoretical Background section fairly soon Let s now look at the procedure for setting up a calculation with localized orbitals using the silane molecule SiH as an example The required input files have been set up appropriately in CASINO examples molecule splines silane_localized It is assumed that the user is trying out this example on a serial workstation The generate_spline utility can be run in parallel an appropriate job submission script has to be used 1 Carrying out the transformation from Bloch to Wannier orbitals The first stage requires the pwfn data and input _wannier files Type rwannier to run the wannier con version utility The output consists of a new plane wave file pwfn data wannier and a file wannier_centers dat that contains a list of the Wannier function centres 2 Type mv pwfn data pwfn data bloch and then type mv pwfn data wannier pwfn data From now on we will work with localized orbitals 3 Generation of a spline representation of the orbitals We could run CASINO using the localized orbitals represented in the plane wave basis however it
90. conds with Jastrow Total energy 29 221776570372 0 002571501258 Time 12 2700 seconds After any CASINO calculation you can normally delete all the output that CASINO produces and return the directory back to the state it was in at the start of the calculation by typing clearup 6 2 7 How to do a variance minimization calculation NOTE there is now a new method of optimization used in CASINO for linear parameters which is about thirty times faster than the old one It is not yet documented in this introductory section but see section 10 24 The values for the parameters in the Jastrow factor used in the last part of the previous section and other adjustable parameters such as coefficients in a multideterminant expansion are optimized through a variance minimization calculation This is probably the most difficult part of QMC for beginners and perseverance will help The user might first like to take a look at section 10 23 which contains a detailed summary of the theory and best practice Here we simply summarize 12 We wish to minimize the variance the value of which is given by the first integral expression in section 10 23 This integral is approximated by summing over a fixed set of configurations i e values of the vector R giving the positions of all the electrons and associated energies These configurations must clearly be distributed according to the square of the trial wave function and the first part of a
91. cters are and e Some of the parameters are of block type which means that multiple parameters must be supplied which may be spread over several lines See e g wavefunction or lineplot Here is the current list of the input parameters in alphabetical order Note that the casinohelp facility is always up to date it runs CASINO to find out what it knows about but this manual may not be Input parameters ALIMIT Real Parameter required when LIMDMC 2 or 3 A value of 0 25 was suggested by Umrigar et al for all electron calculations but a value of 1 0 may be more appropriate for pseudopotential calculations The answers are insensitive to the precise value ALIMIT is not required if NUCLEUS_GF_MODS is set to true See Section 10 6 BTYPE Integer Determines type of basis set Possible values 1 Plane waves input from CASTEP PWSCF ABINIT GP k207 Gaussians from CRYSTAL95 98 03 GAUSSIAN94 98 03 TURBOMOLE Numerical atoms only Blip function B spline basis from transformation of pwfn data using the blip converter Gaussian Wannier orbitals non functional Spline basis from transformation of pwfn data using the generate_spline converter Note that the definition of BT YPE was changed in December 2003 This affected only blip and spline inputs and obviously the various electron hole phases Doe w bv CALC_VARIANCE Logical If CALC_VARIANCE is true then the variance of
92. ctuations Note also that in versions of CASINO prior to sometime in early 2001 this procedure was carried out automatically but the DMC population control algorithm was then changed and the renormal ization procedure became unnecessary Note also that using this will introduce a large bias into the results It should not be used for any serious calculations PRINTGSCREENING Logical Before doing a periodic Gaussian calculation CASINO prepares lists of potentially significant primitive cells and sites in each such cell which could contain Gaussians having a non zero value in a reference primitive cell centred on the origin Zero is defined as 10 G4UTOL Turning on the PRINTGSCREENING flag prints out the important information about this screening 28 QMC_DENSITY_MPC Logical If this flag is set then the QMC charge data at the end of the density data file will be used for the MPC interaction RANLUXLEVEL Integer To generate the parallel streams of pseudo random numbers for its stochastic algorithms CASINO uses an implementation of the RANLUX algorithm This is an advanced pseudo random number generator based on the RCARRY algorithm proposed in 1991 by Marsaglia and Zaman RCARRY used a subtract and borrow algorithm with a period on the order of 10171 but still had detectable correlations between numbers Martin Luescher proposed the RANLUX algorithm in 1993 RANLUX generates pseudo random numbers using RCARRY but throws away numbers to destr
93. do two iterations of the config generation optimization cycle the v flag means when you ve finished these two iteration do a final VMC calculation with the optimized Jastrow factor so we can check that it really does lower the variance sometimes it doesn t Note if you just type runvarmin then the default on workstations is equivalent to runvarmin n 8 v The runvarmin script will create a directory io in which all the input and output files from the different stages of the calculation will be placed with a numeric suffix indicating which iteration of config generation optimization they are associated with The most important things to look at are vmc out 1 vmc out 2 and vmc out 3 These are the VMC output files with respectively the unoptimized Jastrow factor the optimized Jastrow factor after the first iteration and the final optimized Jastrow factor obviously there will be even more of these if you increase the argument to the n flag If you type ve vmc out etc then you can quite clearly see whether the optimization has worked at each stage Be careful It won t always work and you should choose the best of the Jastrow factors from the different iterations from jasfun out 1 jasfun out 2 etc in the o directory In this case we see File vmc out 1 Total energy 0 643382740589 0 002557624213 Time 13 8200 seconds File vmc out 2 Total energy 0 613263097995 0 001254910540 Time 21 1500 seconds
94. e crystal input file dna run qmc dna the qmc flag invokes generation of relevant QMC files in temp cd temp mdt or whatever the temp directory is called in your CRYSTAL run script 7 crysgen then answer the questions mv gwin data rm dna then run CASINO 12 3 Generating gwfn data files with CRYSTALO3 The gwfn data file can be generated directly using the CRYSTALO3 properties program Note that this facility only exists in the official version of CRYSTAL in binaries produced after December 2003 and that some early versions of CRYSTALO3 had a broken pseudopotential evaluator which should now be fixed The CRYSTAL run script included with CASINO will actually force the production of gwfn data automatically if you invoke it with the qmc flag when running the calculations i e all you need to is type run qmc input_filename If the calculation is periodic the script will ask you how many different supercell sizes N you wish to generate and then to input N sets of integer triplets indicating the supercell sizes these must be a subset of the MP shrinking factors in the CRYSTAL input file For example run qmc h Number of different QMC supercell sizes to calculate Maximum 5 1 Size of cell 1 e g 2 2 2 22 2 Put the script in the background with Ctrl Z then bg If you want to do this by hand rather than letting the run script do it for you then the relevant part of the CRYSTAL properties inpu
95. e density pa either from an independent particle calculation or a VMC calculation We can avoid the need to know p exactly by constructing a new interaction operator which involves only pa fee DIRE PAD lele r Fx de i gt j 1 gt I palpa e vel r f r r dear 158 Because of the presence of the third term on the right hand side which is a constant there is no double counting and the interaction energy becomes Ese W H W J plr pa r veler x dr dr pa r pa z up r r dr dr ws ws S WPS Im Med f pepate r dede 1 f 1 1 7 e aa ds 159 Noting that Eee ES F e r pato o r paz ve r f r r dr dr 160 WS we see that the error due to this approximation is second order in p pa and in addition the operator vg f becomes very small as the size of the simulation cell goes to infinity The error term is therefore usually small and is neglected although it could be calculated after the simulation We use the MPC expressions of Eqs 158 and 159 in both VMC and DMC calculations The first term of the Hamiltonian of Eq 158 is evaluated in real space and the second term in Fourier space The third term is a constant which is evaluated in reciprocal space at the start of the calculation Introducing the Fourier transformed quantities 1 fa qf Oid 161 1 iG r pa a fee Crar 162 90 where 2 is the volume of the
96. e here in Cambridge CASINO reads pseudopotential data from a file called xx_pp data where xx is the chemical symbol of the element in lower case letters The file contains a logarithmic radial grid and the values of the different angular momentum components of the pseudopotential at each grid point A library of pseuopotentials suitable for use with CASINO is available at www tcm phy cam ac uk mdt26 casino_users html i An important point is that exactly the same pseudopotential should be used in the DFT HF calculation that generates the trial wave function and in CASINO Other programs do not in general understand the CASINO pseudopotential format and so the information must somehow be transformed so that they do For programs using Gaussian basis sets such as GAUSSIAN9X and CRYSTAL the pseudopotential must be reexpanded in Gaussian functions times powers of r If you are using the Cambridge pseudopotentials these expansions in formats suitable for these two programs are included in the on line library If you are using some other pseudopotential then there is a CASINO utility ppfit which can be persuaded to do the expansion for you There is an additional utility ptm which can manipulate pseudopotential files on grids and their Gaussian expansions in various useful ways As for plane wave programs CASTEP understands the CASINO grid format and can read such files directly Other plane wave programs require conversion util
97. e interested in via Root N in the CIS options With all of this done properly Gaussian produces some output like N N 6 9546274D 00 E N 2 3781323D 01 KE 3 3771062D 00 units of Ha which is hidden in the density analysis right at the end of the output e Some trouble has been encountered with the Add N option to CIS which reads converged ex cited states from the checkpoint file and then calculates N more The IOp alternative which does work is IOp 9 49 2 use guess vectors from the checkpoint file combined with IOp 9 39 N make N additional guesses to those present 107 Using the 50 50 option to CIS to calculate singlet and triplet excitations simultaneously can cause problems It appears best to do the singlet singlets and triplet triplets calculations separately For QMC we want the complete CIS expansion Gaussian may be persuaded to output all excited states with coefficients gt 107 by using IOp 9 40 N Typically N 5 is good enough gaussiantoqme outputs the sum of the square of the coefficients so that the user can see how complete the wave function is Standard Gaussian output has the coefficients normalized so that their sum of squares for a complete expansion should be unity 13 2 4 CISD Although gaussiantogmc cannot read a CISD wave function it might be worth mentioning that IOp 9 6 N is equivalent to MAXCYCLE N for such a calculation 13 2 5 CASSCF As described in Foresman
98. e ug is chosen to obey the cusp electron electron conditions S1 S3 54 and S are chosen to be cusp less L is the range of the short distance correlation functions which is normally chosen to be of order a bond length and L is the range of the long distance correlation function which is normally chosen to be the Wigner Seitz radius for a periodic system or approximately the size of the system for finite systems Functions which depend on the distances to ions can be chosen to be different for different sets of ions The most important terms are uy and S which give a homogeneous u function and S2 which amounts to a x function 10 15 Wave function updating Consider a Slater Jastrow wave function v R e Dl ri suas Tn D rn 41 ena Ty gt 111 where J is the Jastrow factor and D and D are Slater determinants for the spin up and spin down electrons respectively We will need to calculate the ratio of the new wave function to the old when 82 new for example the ith spin up electron is moved from r to r Y The wave function ratio can be written as j RY _ d e RAW 112 Y R eT Rola where i DV 91 7a 2 25 TIM T q Tiro 2 N 113 D I r1 Pajas r yt EN The contribution from the Jastrow factor is easily evaluated A direct calculation of the determi nants in q at every move by LU decomposition is time consuming and instead we use an updating method We define the Slater matr
99. ectron hole gases uy is set to 5 A Tij Tij Taa uo Tiz E exp 2 4 3 exp 97 a Tij 2F F L 81 where cusp conditions make F 424 and Fy 4 3 The functions S to Ss are l1 L S Tij LY ri 5 aTilT B rij TA T 5 2 rs 98 1 0 12 L Sa ri Dri AE Br L G ri 99 lm3 lm3 S3 ra LP ry LP XO YO imil Tm 5 100 l 0 m 0 l4 m4 Sa ri a L rij m 2L 5 rr 2 Ne 5 Elm Til FT m Tiz 101 l 0 m 0 lm5 lm5 nm5 Ss r LP rj Drririra Y YO mn TFT PT a Fis 102 l 0 m 0 n 0 where 7 denotes the lth Chebyshev polynomial and the expansion range is rescaled to the orthogo nality interval of the Chebyshev polynomials 1 1 so that 2rij L Tij ag 103 2r L 7 cy 104 r 104 ij L We assume that terms are zero according to S 0 when rij gt L 106 So 0 when r gt L 107 S 0 when r gt L or r gt L 108 S4 0 when r gt L or rij gt 2L 109 S5 0 when r gt L or rj gt L 110 To maintain the even symmetry of the Jastrow factor under interchange of the positions of two electrons Yim and wimn are forced to be symmetric under interchange of l and m The parameters in S to Ss are obtained from variance minimization The computational cost of including terms other than uo S and S2 can be large S is chosen to be cusp less as we assume that the ionic potentials are finite at the ion centre Becaus
100. ed See section 10 24 NLCUTOFFTOL Real This is used to define the cutoff radius for the non local parts of the pseudopotential It is defined as the maximum deviation of the non local potentials from the local potential at the non local cutoff radius NLRULE1 Integer Rule for non local integration in production VMC or DMC run see Sec tion 10 18 Currently assumed to be the same for all atoms if it is controlled through this keyword you can provide an override value of NURULE1 for particular atoms at the top of the corresponding pseudopotential file NLRULE2 Integer Rule for non local integration in configuration generation for DMC or wave function optimization see Section 10 18 Currently assumed the same for all atoms if it is con trolled through this keyword you can provide an override value of NLRULE2 for particular atoms at the top of the corresponding pseudopotential file It is usual to set NLRULE2 greater than NLRULE1 as the accuracy of the integration is more critical in config generation errors in the non local energy bias wave function optimization and DMC calculations but note that set ting it to something different than NLRULE1 can lead to a significant performance disadvantage when VMC_METHOD 1 is used as all the energies have to be computed again NMOVE Integer For DMC NMOVE is the number of moves of all electrons in a block For VMC the total number of moves of all electrons in a block is NUOVExCORPERxNVMCAVE
101. ed out in order to optimize the cutoff lengths followed by an accurate optimization of the linear parameters using the new method For some systems good values of the cutoff lengths can be supplied immediately for example in the fluid phase of an electron gas the cutoff length La should be set equal to the Wigner Seitz radius of the simulation cell and one can make use of the new optimization method straight away The quartic LSF coefficients are stored in a file called varmin_coeffs data This file is gen erated in a VMC simulation when the use_newopt keyword is set to T If a VMC simulation is performed in a directory in which a varmin_coeffs data file is already present the existing data in the file will be added to it is therefore possible to refine the quality of quartic coefficient data by extending a VMC calculation Note that varmin_coeffs data is an unformatted file how ever a formatting unformatting utility called FORMAT_VARMIN_COEFFS exists This utility reads varmin_coeffs data and produces a formatted varmin_coeffs data_formatted file and vice versa varmin_coeffs data is deleted by the RUNVARMIN script at the end of each optimization calculation so that one starts afresh at every VMC coefficient generation stage The new optimization method can only be used in conjunction with CASINO s new Jastrow factor The method is not yet implemented for electron hole systems 10 24 3 Further developments to the new optimi
102. eeds to know about it will ask you about these if it can t find input If all files are present reblock will read these files and then starting asking you questions What units do you want the answer in Whatever you want At which line of dmc hist do you want to start computing statistics The move number where the green line in the graphit output becomes approximately constant e g 501 in the NiO case Average over how many lines Usually input 0 meaning use the rest of the now equilibrated data in the file from lines 501 to the end Reblock will then show you a reblocking analysis and ask you to choose a block length Choose a value where the Std err of mean column starts to plateau again the plot_reblock utility can help you with this A value of 64 might be typical Reblock will then print out the final DMC energies and error bars together with an analysis of the population fluctuations effective time steps and acceptance ratios In the case of molecular hydrogen the final energy is 1 1744731 0 00056 The exact energy is 1 1744757 a u so this result is pretty good without paying particular attention to getting it absolutely right which is promising 6 2 9 How to generate localized orbitals For large systems a significant speed up may be gained by using localized orbitals in CASINO instead of the canonical delocalized orbitals generated by the DFT HF code Currently this facility is not fully implemented and su
103. efficients any that are not listed in jastrow data are assumed to be zero Furthermore if too many are given then the extra parameters will be ignored If all of the parameter values in all of the Jastrow terms are blank as is the case in the example given above then only the Slater wave function will be used for the first configuration generation run when performing variance minimization Different x functions are used for different species of ion the number of sets should be chosen to be equal to the number of chemically distinct species The ions in each set are specified by giving a list of the numbers that label them these are the same as the labels used in the xwfn data file that specifies the geometry of the system It is possible to make the Jastrow factor enforce the electron nucleus cusp condition This should only be done if the y set contains bare nuclei and the orbitals do not satisfy the cusp condition With a Gaussian basis set it is much better to use the in built cusp correction algorithm activated with the input keyword CUSP_CORRECTION rather than use the Jastrow There are two spin dependence options for x if it is 0 then the same parameters are used for up and down spin electrons if it is 1 then different parameters are used Similar comments to those made for u apply to the spherical harmonic labels the cutoff length and the parameter values of x The f function contains terms that approximatel
104. egenerate orbital These files detail all excitations out of the specified orbitals along with a percentage giving their contribution to the CIS expansion as a whole The sum of the percentages final column from each of the fromi_j dat should be 100 if everything is working OK 106 13 2 Other features of gaussiantoqme The code also contains some crude normalization and plotting routines that are really just debugging aids By setting the flag test true in gaussiantoqmc f90 and recompiling the user is given the option of plotting and testing the normalization of individual molecular orbitals The axis along which the plotting is done is set in wfn_construct f 90 and this must be hacked if the user wishes to change things 13 2 1 Getting Gaussian to do what you want In principle Gaussian can do an awful lot of things In fact some of these things seem to require magical incantations These will be described in this section broken down into the different calcu lations to which they apply The comments on G98 refer to revision A9 and may depend on which version is used 13 2 2 General bits and pieces Some points to note e Both G94 and G98 appear to have have formatting errors when printing out the Gaussians used in the ECP expansion large exponent values are replaced with stars This does not affect the subsequent calculation e G94 s ECP pseudopotential package will not accept expansions containing more than 13 Gaus
105. er crystal Dimensionality 1 0 1 0 true Flag to separate layers z axis separation of layers au relevant only in 2D 1 0 Include all cells contributing to wave function true Number of gaussians Ng ndet 2 only Ng lines coefficients a_i Gaussian exponents b_i 1 1 0 2 0 Fit Gaussian parameters to exp r true or allow to vary freely false false Accumulate pair correlation function number of bins false 5000 Accumulate 1 particle density matrix Type of crystal 3D C ubic B cc F cc 2D R ectangular H exagonal T riangular manual M 50 F Ferromagnetic F or if not frustrated antiferromagnetic A F Offset x y between lattices in each layer relevant only for 2D systems 0 d0 0 d0 optimize all parameters true Pade coefficients k1 k2 for in order up elecs down elecs up holes down holes 0 354388900376931 0 412126133902296 0 706867779971532 0 287471090722791 Number of localized orbitals 1 Positions of localized wave functions and spin of associated particle ci alt c2 a2 c3 a3 spin 0 d0 0 d0 0 d0 1 5l 7 8 Jellium slabs heg data file with ETYPE 6 When ETYPE 6 the heg data file format is entirely different from that of the various homogeneous electron or hole systems available The simulation cell size in the x and y directions is determined by a combination of the density given here by the r s parameter the number of electrons and the slab width
106. etting the phases for single orbitals it is possible to set phases for an entire block using ABL ABL 0 5231122 This command has the same effect as the ASINGLE instruction above except that the the last four numbers refer to the minimum k point and minimum band index and the maximum k point and maximum band index of the block in this order 8 2 5 Promotions GSPEC DPSEC POL and EPOL are very flexible however using a few flags is unlikely to cover all possible occupations of the given orbitals so once the number of electrons has been set promotions can be used to generate any possible occupation scheme The syntax is DET idet ispin PR bfrom kfrom bto kto 59 ispin and idet determining which determinant and spin are being dealt with bfrom and kfrom refer to the band and k point of the electron that is to be moved bto and kto are the positions the electron is to be moved to Promotions my be combined with setting phases by using DET idet ispin PRA bfrom kfrom phfrom bto kto phto If applicable phfrom is the phase to be used for the electron remaining at bfrom and kfrom and phto is the phase of the target orbital e g DET 1 3 PRA 4 3 0 2 4 4 0 2 8 2 6 Defaults No values except the number of up and down electrons have to be specified The CASINO defaults for all determinant are POL 1 DIR SPR SPRBND 1 ie true UNITA 7 EPOL 0 0 All phases are set to 0 0 8 3 Blip and spline basis Note that the setup for t
107. ext time step Terr is the current best estimate over all configurations and time steps of the mean effective time step calculated using Eq 54 with A Tela m Note that g is the time scale in terms of time steps over which the population attempts to return to Mo During accumulation Ebest is set equal to the current value of the mixed estimator of the energy given by Eq 54 with A Ez a m To further limit the size of population fluctuations the population can be renormalized to Mo at the start of each block POPRENORM This should not be necessary however and should be avoided if at all possible as it can lead to large time step errors 10 6 Modifications to the Green s Function 10 6 1 The effective time step Time step errors can be reduced and the stability of the DMC algorithm improved by modifying the Green s function An important modification is to introduce an effective time step Teg into the branching factor 4 When the accept reject step is included the mean distance diffused by each electron each move which should go as the square root of the time step is reduced because some moves are rejected When calculating branching factors it is therefore more accurate to use a time step appropriate for the actual distance diffused Umrigar and Filippi 7 have suggested using an effective time step for each configuration at each time step This helps to eliminate the problem of persistent configurations which are low ener
108. f all keywords that CASINO knows about or casinohelp keyword for detailed help on a particular keyword Type casinohelp search text to search for the string text in all the keyword descriptions Create a blank file containing the keyword input_example and type runvmc this will create a sample input file containing all valid input keywords and their default values in the correct format Although there are many keywords the beginning user can play around by changing only a few of them Here s a very rough guide Advice on good values to use will be given in the subsequent sections explaining how to do VMC DMC and variance minimization calculations General system dependent keywords you should always change in addition to the title given in the input file neu ned number of electrons of up and down spin isperiodic whether the system is periodic or not npcells number of primitive cells making up the simulation cell not required for finite systems Important VMC keywords nequil number of Metropolis equilibration steps nmove number of VMC moves dtvmc VMC time step corper VMC correlation period Important DMC keywords nmove number of VMC moves to perform in the config generation step nwrcon number of configs per node to generate in the config generation step nconfig number of configs per node to use in DMC nmove_dmc_equil number of moves to use in DMC equilibration nblock_dmc_equil number of blocks to use
109. fe side and the last one in theory tells it to just print the configurations although in actual fact it still proceeds to do the calculation as well Restarting a CASSCF calculation from a previously converged run fails when symmetry is switched on Use Nosymm to avoid this problem 13 2 6 TD HF and TD DFT As far as converting the resulting wave function for use in CASINO is concerned these two methods are no different to CIS apart from the issue of normalization In CIS the default output is normalized such that the sum of the squares of the coefficients is equal to unity In a TD HF or TD DFT calculation which involves solving a non Hermitian eigenvalue problem a different scheme is used which essentially means that the sum of the squares of the coefficients is arbitrary It should also be noted that Gaussian cannot do gradients within TD DFT yet and so cannot relax excited states 109 13 3 Summary of routines used in gausstantoqmc Routine Purpose analyse_cis_state awk_like cas_win cas_write cis_data con_coefts fatal g94_wave function gaussiantoqmc g d_type g s type get_gauss_version integ_params max_coincidence normalisation_check normalise_ci numsrt numsrt_2way numsrt_signchange pack_evcoeffs paramfile psi Break the selected CIS TD DFT state down into excita tions from each distinct occ MO Module and subroutines to give AWK like functionality Used in parsing the Gaussian
110. file the following block entry controls the movie making process MOVIES makemovie T Make movie Boolean movieplot 1 Plot every moves Inte movienode 0 Node to plot Integer moviecells F Plot nn cells Boolean Set the keyword makemovie to T to enable the movie making facility You could set movieplot to an integer greater than 1 so that the particle positions are only written out every movieplot moves The node which plots the particle positions are chosen by the keyword movienode which has to be zero or a positive integer less than the total number of nodes If moviecells is set to F then the unit cell will be plotted if set to T then nearest neighbour cells in the x y plane will also be written Type runvmc and an output file movie out will be produced This is a file in the standard xyz molecular format And example movie out file will look like 4 Input geometry 0 000000 0 000000 1 407651 1 138185 054434 1 000000 894354 0 554315 263301 1 000000 0 528074 1 081535 0 755823 1 000000 4 Input geometry C 0 000000 0 000000 000000 5 000000 H 1 407651 1 138185 054434 1 000000 H 0 894354 0 554315 263301 1 000000 H 000000 5 000000 TETTA o e OO e OO 0 528074 1 081535 0 755823 1 000000 Line 1 indicates the total number n of ions and particles In this case we have 4 particles Line 2 is a line of comment The following n lines consist of 5 columns Column 1 specifies the t
111. g the Fourier transform of the charge density and writing the result to the density data file which CASINO will do if you set IRUN 6 the program only writes data for those G vectors where the absolute value of either the real part or the imaginary part of the Fourier coefficient is larger than a threshold This is that threshold The default of 107 should be good enough DTVMC Real Time step for VMC run atomic units DTVMC2 Real Time step for VMC run atomic units In an electron hole calculation DTVMC is used to move the electrons while DT VMC2 is used to move the holes DTDMC Real Time step for DMC run atomic units EDIST_BY_ION EDIST BY IONTYPE Block The optional EDIST BY ION block allows fine control of the initial distribution of the electrons before equilibration starts The standard al gorithm shares out the electrons amongst the various ions weighted by the pseudo charge atomic number of the ion Each electron is placed randomly on the surface of a sphere surrounding its parent ion There are certain situations for example a simple crystal with a very large lattice constant where the standard algorithm in the POINTS routine may give a bad initial distribu tion which cannot be undone by equilibrating for a reasonable amount of time This keyword allows a user defined set of electron ion associations to be supplied The syntax is to supply N_ION lines within the block which look like e g 1 4 4 where the three numbers a
112. g the non local pseudopotential energy The action of the non local pseudopotential on the wave function can be written as a sum of contribu tions from each electron and each angular momentum channel The contribution to the local energy made by the non local pseudopotential is Va YY y Vy wa YO Vasa 132 where for simplicity we consider the case of a single atom placed at the origin Va may be written as 19 ai DES Vas vel SE f Pileos e 1 Wir Di Ep Pity EN V r i Y 1 1 5 1 541 T y dy 133 where P denotes a Legendre polynomial CASINO currently performs the non local projections for l 0 1 2 only The integral over the surface of the sphere in Eq 133 is evaluated numerically The r dependence of the many body wave function is expected to have predominantly the angular momentum character of the orbitals in the Slater part of the wave function A suitable integration scheme is therefore to use a quadrature rule that integrates products of spherical harmonics exactly up to some maximum value Imax The quadrature grids currently available are listed in Table 1 To avoid bias the orientation of the axes is chosen randomly each time such an integral is evaluated Table 1 Quadrature grids for the non local integration NLRULE is the label for the rule see Section 7 1 max is the maximum value of l which is integrated exactly and Np is the number of points in the grid NLRULE max Np 1 0 1 2 2 4
113. generate dat file The pwfn data file should usually be the renamed pwfn data wannier generated by xwannier The corresponding wannier_centers dat file should also be supplied However it is possible to bypass the Wannier transformation and generate a spline representation of Bloch orbitals In this case no wannier_centers dat file is required and the pwfn data file with the Bloch orbitals should be used The format of the generate dat file is shown below false 2 d0 1 1 10 d0 0 0 05d0 2 Plot out the truncated splined orbitals Multiplication factor of the FFT grid 1 gt spherical cutoff 2 gt cuboidal cutoff 1 gt Use norm 2 gt Use fixed diameter Fixed cutoff diameter au if set to 2 above No of states for which the whole box is to be used Minimum shell width au Verbosity level 0 997d0 Cutoff criterion if set to 1 above If the first line is set to true then plots of the splined orbitals in the z y and z directions through the Wannier centres or the origin for Bloch orbitals are produced The second line controls the fineness of the real space spline grid The pwfn data file holds plane wave coefficients on a cuboidal array of reciprocal lattice vectors The orbital is evaluated at a grid of points in real space using an inverse fast Fourier transform The real space grid spacing in a particular direction is given by 27 over the length of the cuboid in the correspo
114. gy configurations for which all moves are rejected and which can multiply and bias the calculation The effective time step is given by i Pi Ara Dard a Teala m where the averages are over all attempted moves of the electrons 7 in configuration a at time step m The Ara are the diffusive displacements i e the distance travelled by the electrons without the drift displacement see Eq 15 and p is the acceptance probability of the electron move see Eq 19 The values averaged over the current run are written in the output file We calculate Terr m using Eq 54 with A Ter Qt m 10 6 2 Drift vector and local energy limiting The drift vector diverges at the nodal surface and a configuration which approaches a node can exhibit a very large drift resulting in an excessively large move in the configuration space One can improve the Green s function by cutting off the drift vector when its magnitude becomes large The total drift vector is defined in Eq 14 We use the smoothly cut off drift vector suggested by UNR 5 For each electron with drift vector vi we define the smoothly cut off drift vector v where 1 1 2 il ae 29 alv 7 where a is a constant which can be chosen to minimize the bias in UNR 5 a value of a 1 4 was suggested for all electron calculations but a 1 may be more appropriate for pseudopotential calculations The value of a ALIMIT must be entered by the user if NUCLEUS_GF_
115. h VMC and DMC methods the energy is calculated as the average over many points in the electron configuration space of the local energy 75 E V AW where is the Hamiltonian and W is the many electron trial wave function When an electron approaches a nucleus of charge Z the potential energy contribution to E diverges as Z r where r is the distance from the nucleus The kinetic energy operator acting on the cusps in the wave function must therefore supply an equal and opposite divergence in the local kinetic energy because the local energy is constant everywhere in the configuration space if Y is an eigenstate of the Hamiltonian Unfortunately when using orbitals expanded in a Gaussian basis set the kinetic energy is finite at the nucleus and therefore Ey diverges In practice one finds that the local energy has wild oscillations close to the nucleus which give rise to a large variance in the energy This is undesirable in VMC but within DMC it can lead to severe bias and even to catastrophic numerical instabilities It might seem that using local basis functions which can be made individually to obey the cusp conditions at the nucleus on which they are centered such as a Slater type basis of cusped exponential functions would solve the problem This is only partially true however One issue is the fact that the codes available for performing calculations with Slater type functions are more limited and less widely used than those which e
116. h varmin Otherwise the least squares function will simply be the variance of the configuration local energies VM_W_MIN Real Minimum value that a configuration weight may take during weighted variance minimization This parameter should have a value between zero and one Note that the limiting is not applied if VM_W_MAX 0 VM_W_MAX Real Maximum value that a configuration weight may take during weighted variance minimization Set this to zero if you do not wish to limit the weights otherwise it should be greater than 1 WAVEFUNCTION Block The wave function block specifies the nature of the Slater part of the many electron wave function The xwfn data file specifies a reference configuration which may consist of one or more determinants The input file allows use of either the reference configuration GS or to specify excitations additions subtractions from GS If the excita tions additions subtractions are made to a single determinant GS configuration then we have to specify SD in the input file and if they are made to a multi determinant GS configuration we have to specify MD Excitations additions subtractions are made using the keywords PR PL and MI respectively Several changes to the GS configuration can be made at once WDMCMIN WDMCMAX Reals Weighted DMC only minimum and maximum weights See entry for LWDMC 31 7 2 Jastrow factor file jastrow data The file that holds the parameter values for the Jastro
117. hat because of the stochastic nature of the method the energy has an associated error bar users of DFT etc may find this disconcerting The true Hartree Fock energy of the hydrogen atom is exactly 0 5 a u so we have seemingly correctly performed out first VMC calculation Next look in the file out We see a complete report of the calculation and at the end we see the total energy and its components Just before the total energy we see the acceptance ratio This should be around 50 If it isn t then you should play with the VMC timestep dtvmc either manually or by activating the opt_dtvmc flag with nequil set to at least 500 whereupon the timestep will be optimized automatically If the acceptance ratio is very large then the electrons are only attempting to move very short distances at each step of the random walk They are thus likely to be bumbling about in one corner of the configuration of space and the serial correlation is likely to be large increase dtvmc to prevent this If the acceptance ratio is small then the electrons are attempting very large jumps which are likely to be to regions of lower probability and are thus likely to be rejected by the Metropolis algorithm Thus the electrons don t move very much and again the serial correlation is large Reduce dtvmc to prevent this Next look at the file vmc hist see section 7 9 This contains mainly energy data produced at each step of the random walk you can average
118. he billy distribution but can be used on its own e density_plotter Three programs denconvert d2rs plot_density The denconvert script runs d2rs which converts the density data file density expanded in plane waves into a real space density on a grid This can be converted into a data format readable by gnuplot using the plot_density utility e extrapolate_N Extrapolate HEG energies to infinite electron numbers using a chi squared fit to Ceperley s extrapolation formula 36 Documentation available in the relevant directory e extrapolate_N_eh As above for electron hole systems e extrapolate_tau Program to extrapolate DMC energies to zero stepsize by fitting a polyno mial form to the DMC energy as a function of stepsize Documentation in the relevant utils directory e format_configs Format an unformatted config in file config in_formatted or vice versa e format_pos Format an unformatted vmc posin file gt vmc posin_formatted or vice versa e format_spline Utility for formatting unformatting swfn data files e format_varmin coeffs Format an unformatted varmin_coeffs data file gt varmin_coeffs data formatted or vice versa e generate spline Utility to construct a spline interpolation in 3D e graphit Plot energy vs move data from a DMC run i e the numbers in dmc hist in an xamgr plot e help Simple script to invoke the CASINO help system casinohelp lt casino_keyword gt tells
119. he blip function module and the spline function module has not yet been fully merged with the plane wave setup routines so that one cannot currently use the fancy orbital occupation schemes excited states or multi determinant wave functions with blip or spline orbitals This functionality will almost certainly be included shortly or even sooner if someone specifically requests it 60 9 Pseudocode 9 1 Subroutine VMC Loop over blocks of time steps Loop over time steps Loop over configurations Loop over electron in configuration Calculate local energy at Roa if required Propose move of electron Calculate ratio U Rnyew U Roia Calculate local energy at Rnew if required Perform Metropolis accept reject step If move is accepted update cofactor matrix and determinant Accumulate local energy etc if required End loop over electrons in configuration End loop over configurations End loop over time steps End loop over blocks of time steps 9 2 Subroutine DMC Loop over blocks of time steps Loop over moves in a block Loop over configurations Loop over electrons in configuration Evaluate drift vector at old position Use one of various schemes to evaluate corresponding limited vector Evaluate diffusive component of the electron displacement Generate new position of electron using limited drift vector If crossed node then reject move of this electron Evaluate drift vector at new position Use
120. he reference energy which acts as a renormalization factor see Section 10 5 Ey is the local energy E R vw Av 13 where H is the Hamiltonian and V is the drift vector V R v vv 14 Note that V vi vw where v is the drift vector of electron i 10 3 2 The ensemble of configurations The f distribution is represented by an ensemble of electron configurations which are propagated according to rules derived from the Green s function of Eq 10 Gp represents a drift diffusion process while Gg represents a branching process The branching process leads to fluctuations in the population of configurations and or fluctuations in their weights We will introduce labels for the different configurations present at each time step m From now on R represents the electronic positions of a particular configuration in the ensemble and i labels a particular electron In CASINO electrons are moved one at a time This is much more efficient than making whole configuration moves for large systems This is standard procedure for large systems but most of the algorithms described in the literature are for moving all electrons at once 10 4 Drift and diffusion We now discuss the practical implementation of the drift diffusion process using the electron by electron algorithm in which electrons are moved one at a time as used in CASINO To implement the drift diffusion step each electron 7 in each configuration a is moved from r a t
121. hes Relevant only if the TRIP_POPN keyword is set to a non zero value See the discussion in Section 10 10 for more details MOVIECELLS Logical If FALSE then MAKEMOVIE will plot the unit cell if TRUE then nearest neighbour cells in the x y plane will also be written see Section 11 MOVIENODE Integer Plot the particle positions on node MOVIENODE see Section 11 MOVIEPLOT Integer Plot the particle positions every MOVIEPLOT moves see Section 11 NBLOCK Integer The total number of blocks of NMOVE moves in the run N B you should use the REBLOCK utility in the utils directory to see the effect of varying the block size on the variance after the calculation has completed In VMC NBLOCK just determines how often block averaged quantities are written to the output file In DMC the configuration population is redistributed among the nodes and may in extremis be renormalized at the end of each block NCONFIG Integer DMC only Initial number of configurations per node This can be different from NWRCON number of configs written out by preliminary VMC run due to the possibility of the config generation and DMC simulation being run on different numbers of nodes on a parallel machine this might be desirable as the VMC part is often very quick and can be run on relatively few nodes It is possible to carry out a DMC run with a target weight of NCONFIG 1 configurations per node provided you are running on more than one node This is
122. hould never cross the nucleus rather they should end up on top of it In order to impose this condition at finite timesteps we work in cylindrical polar coordinates with the z axis lying along the line from the nucleus to the electron Let the position of the closest nucleus be R GF The position of the electron relative to the nucleus is r Rz 2 e 43 while the limited drift velocity can be resolved as V Uz z UpEp 44 where e is a unit vector orthogonal to e The new z coordinate after drifting for one timestep is z max z 0 7 O 45 which cannot lie beyond the nucleus The drift in the radial direction over one timestep is 20 72 mo P 46 The new radial coordinate is approximately 0 7 when far from the nucleus but it is forced to go to zero as the nucleus is approached Hence if the electron attempts to overshoot the nucleus it will end up on top of it So timestep errors caused by drifting across nuclei are eliminated Let the electron position at the end of the drift process be r ze pep 71 10 7 3 Diffusion close to a bare nucleus Close to a nucleus f is proportional to the square of the hydrogenic 1s orbital assuming the trial wavefunction has the correct behaviour This cusp cannot be reproduced by Gaussian diffusion at finite timesteps In fact starting from the nucleus we would like our electron to take a random step w distributed according to exp 2Z w However we only
123. hysics and molecular quantum chemistry Currently the most popular way to include the effects of electron correlation in these calculations is density functional theory This method is in principle exact in reality fast and often very accurate but does have a certain number of well known limitations In particular with only limited knowledge available concerning the exact mathematical form of the so called exchange correlation functional the accuracy of the approximate form of the theory is non uniform and non universal and there are important classes of materials for which it gives qualitatively wrong answers An important and complementary alternative for situations where accuracy is paramount is the quantum Monte Carlo QMC method which has many attractive features for probing the electronic structure of real systems It is an explicitly many body method applicable to both finite and periodic systems which takes electron correlation into account from the outset It gives consistent highly accurate results while at the same time exhibiting favourable scaling of computational cost with system size This is in sharp contrast to the accurate methods used in mainstream quantum chemistry such as configuration interaction or coupled cluster theory which are impractical for anything other than small molecules and cannot generally be applied to condensed matter problems The use of quantum Monte Carlo has been greatly hampered over the last two decades b
124. ider the periodic charge density consisting of a unit point charge at r in every simulation cell plus a uniform cancelling background 1 pile 001 8 5 140 R where R denotes the lattice translation vectors and Q is the volume of the simulation cell The Ewald formula for the periodic potential corresponding to this charge density is erfe y2 r rj R ve r rj yz r r R y Ey RES a map aa G 0 where G denotes the reciprocal lattice translation vectors The value of up r r is in principle independent of the screening parameter y and also in practice providing enough vectors are included in the sums The larger the value of y the more rapidly convergent is the real space sum but the more slowly convergent is the reciprocal space sum A compromise is required to minimize the overall computational cost In CASINO this parameter is set to 2 8 01 3 which approximately minimizes the cost for a wide variety of Bravais lattices 37 The full periodic potential of the simulation cell is obtained by adding the potentials of all the N charges and their cancelling backgrounds which sum to zero because the cell has no net charge N v r 5 jun r rj 142 j 1 The potential acting on the charge at r is therefore N v r 5 qjue ri rj q5 143 i i where E lim werd 144 R r2 Qy MIncreasing the input parameter EWALD CONTROL will increase the number of reciprocal vector
125. ies provided in the directory utils wfn_converters The gwfn data specifies a reference configuration which may consist of one or more determi nants The input file allows use of either the reference configuration GS or to specify excita tions additions subtractions from GS If the excitations additions subtractions are made to a single determinant GS configuration then we have to specify SD in the INPUT file and if they are made to a multi determinant GS configuration we have to specify MD The simplest case is when we want to use exactly the GS configuration specified in gwfn data The INPUT file could then contain 3232 NEU NED 32 spin up and down electrons GS Ground state calculation as specified in gwfn data Excitations additions subtractions are using the keywords PR PL and MI respectively Sev eral changes to the GS configuration can be made at once For example 32 32 SD DET11PR2151 NEU NED 32 spin up and down electrons Single determinant calculation Promote PR electron in determinant 1 spin 1 up band 2 k point 1 to band 5 k point 1 in determinant 1 spin 1 up Add electron spin 2 down to det 1 band 6 k point 1 Remove the added electron DET 12 PL61 DET 12MI61 Here is an example for a multi determinant case 32 32 NEU NED 32 spin up and down electrons MD Multi determinant calculation 3 3 determinants 1d0 Determinant 1 prefactor 1d0 Determinant 1 p
126. ils SQMC_ARCH to your path using shell dependent commands see above for examples 3 Type make in the CASINO sre directory 4 Type make in the CASINO utils directory 5 When compilation has finished type rehash if using csh or tcsh All required programs will now be in your path 6 That s it If e g the default location of the MPI library doesn t work on your machine you can personalize it by setting an environment variable QMCID your machine name Create a file called CASINO src zmakes users QMC_ARCH QMC ID inc containing override definitions for the library locations If you send such files to Mike mdt26 at cam ac uk then they can be perma nently included in the CASINO distribution if you wish The QMC_ID variable need not be set if you are happy with the defaults Active QMC_IDs are listed at the bottom of the CASINO ARCH file Note that the MPI library is not required to compile CASINO on single processor machines No other external libraries are required 5 2 Unsupported architecture 1 Write a CASINO src zmakes machine_nameJ inc file containing compiler name and appropriate flags for optimized debugged profiled code using the existing include files as a model Copy it with the evident modifications into CASINO utils zmakes machine nameJ inc If your ma chine uses non standard library locations or compiler flags or whatever you may also create a CASINO src zmakes machine type machine_name inc fi
127. in DMC equilibration nmove_dmc stats number of moves to use in DMC statistics accumulation nblock_dmc_stats number of moves to use in DMC stats accumulation dtdmc DMC time step Important variance minimization keywords nmove number of VMC moves to perform in the config generation step nwrcon number of configs per node to generate in the config generation step corper VMC correlation period 6 2 4 Jastrow factor file If you run a CASINO VMC calculation using a trial wave function consisting only of a single determinant of orbitals then the result will be the Hartree Fock energy If the orbitals were generated using a Hartree Fock calculation then the VMC HF energy should agree with the HF energy from the other code This is a good check that everything is being done correctly Obviously if the determinant is made up of Kohn Sham orbitals from a DFT calculation then the energies will not agree since the DFT program adds an exchange correlation energy deduced from the self consistent density The full Slater Jastrow trial function normally used in QMC requires the determinantal part of the wave function stored in xwfn data to be multiplied by a separate Jastrow factor which defines the functional form of explicit interparticle correlations In a typical VMC calculation one might recover 60 80 per cent of the correlation energy in this way This is not really enough to be generally useful and in practice the main use of VMC i
128. individual single electron drift velocities The ratio of drift velocities is therefore calculated as VR of R 34 R an VR V R l The limited single electron drift velocities are calculated at the same time as the local energy of the configuration after all of the electron positions in the configuration have been updated 10 7 Modifications to the DMC Green s function at bare nuclei 10 7 1 Modifications to the limiting of the drift velocity The limiting of the drift velocity is intended to remove the divergence at the nodal surface interference with the cusps at bare nuclei is an undesirable side effect which may introduce bias In order to distinguish between nodes and nuclei UNR make the a parameter in their limiting scheme dependent on electron position Immediately before the limited drift velocity of an electron at r is calculated a is evaluated as 1 222 Lx UA a ray a where v is the unit vector in the direction of the unlimited drift velocity e is the unit vector from the closest bare nucleus to the electron z is the distance of the electron from the nucleus and Z is the atomic number This formula makes a small and hence the limiting weak if the electron is both close to the nearest nucleus and drifting towards it 10 7 2 Preventing electrons from overshooting nuclei In its immediate vicinity the single electron drift velocity is always directed towards a bare nucleus Therefore drifting particles s
129. ion radius re for each orbital and nucleus should depend on the quality of the basis set and on the shape of the orbital in question In particular one would expect the cusp correction radii to become smaller as the quality of the basis set is improved Although clearly many other schemes are possible we choose the re in our implementation as follows The maximum possible cusp correction radius is taken to be fremas 1 Z The actual value of re is then determined by a universal parameter ce CUSP_CONTROL in the CASINO input file for which a default value of 50 was found to be reasonable The cusp correction radius re for each orbital and nucleus is set equal to the largest radius less than remax at which the deviation of the effective one electron local energy calculated with from the ideal curve has a magnitude greater than Z ce Appropriate polynomial coefficients a and the resulting maximum deviation of the effective one electron local energy from the ideal curve are then calculated for this re As a final refinement one might then allow the code to vary re over a relatively small range centered on the initial value recomputing the optimal polynomial cusp correction at each radius in order to optimize further the behavior of the effective one electron local energy This is done by default in the implementation When a Gaussian orbital can be readily identified as for example a 1s orbital it generally does not have a node within remax In man
130. ions are not accepted and can lead to moody lack of cooperation lasting years 20 3 PERFORMANCE CASINO prints a timing analysis at the end of the output file when it has finished a calculation This information can be useful to the programmer in deciding what optimizations to do to make the program run faster The total time is divided into User and System time the sum of which is printed as CPU time These timings are broken down according to specific tasks although by no means all activities are timed only those thought likely to contribute significantly For very small systems such as atoms a significant amount of CPU time will remain unaccounted for since the traditional heavyweight tasks are very fast For medium and large systems often the calculation of the orbitals and their derivatives will dominate closely followed by Jastrow factors and Ewald interactions In systems with pseudopotentials the non local integration is very expensive indeed since it calls the orbital evaluation routine a lot and this expense will increase the larger the non local cutoff radius lower nicutofftol in the input to see the effect of this The User time is time spent executing the instructions that make up the user code including e g calls to subroutine libraries linked to the code 126 The System time is time spent in kernel mode processing I O requests for example or other situations which require the interventi
131. ions e g INTEGER ialloc except where required e g INTEGER INTENT in n e In lists of declared variables adhere to the following order INTENTed dummy arguments first order in inout out INTEGER INTENT in INTEGER ITENT inout etc INTEGER INTENT out etc REAL dp INTENT in COMPLEX dp INTENT in CHARACTER 12 INTENT in LOGICAL INTENT in followed by things which are not arguments INTEGER REAL dp COMPLEX dp CHARACTER 12 LOGICAL Within each class put standard variables on the first line followed by things with attributes like ALLOCATABLE PARAMETER etc on subsequent lines in whatever order seems aes thetically pleasing Note the CHARACTER 12 not CHARACTER 12 which is not in the Fortran90 standard e Don t use tab characters anywhere in the code e All units for reading and writing to be allocated unit numbers using the standard call open_units procedure which you can figure out by looking at other examples in the source e Initial letter of comments to be in capitals Comments to be spelled correctly Comments to be useful 124 This is a legitimate comment this is an illegitimate comment X tihs one is evn more illetigimate as i cant spell XX allocate variables XXX really allocate a 1 e If you add a module USE statement anywhere don t forget to change the dependency list in the Makefile e module USE statements in alphabetical order fo
132. is faster and more memory efficient to represent the localized orbitals numerically The next stage requires a pwfn data file and corresponding wannier_centers dat file if this exists and a generate dat file Type generate_spline to run the utility for producing the splined orbitals The output of this code is a swfn data file 1The only parameter the user is ever likely to want to change in input_wannier is the number of states included in the Wannier transformation all the states in the pwfn data file from 1 up to this number are transformed any remaining states are unchanged In the present case we are including all the states in the transformation 2For this example we shall use localized grids for all the splined states if however some states were non localized then it would be advisable to represent those states over the whole simulation box The multiplication factor in the generate dat file controls the fineness of the spline grid Generally this needs to be set to at least about 2 before the kinetic energy of the splined orbitals is the same as that of the orbitals in the plane wave basis to an high degree of accuracy The orbitals have a spherical truncation surface such that 99 7 of the square of the orbital norm lies inside this radius A further shell of width at least 0 05 a u is saved from truncation the orbital can be brought smoothly to zero over this shell region 3The swfn data file is unformatted If you wish to
133. is is given in section 10 23 Basically you should use as many configurations as possible i e greater than 10000 in general remembering that nwrcon refers to the number of configs per node So if you have a ten node parallel computer then you can use 10000 configurations by setting nwrcon to 1000 Conversely you should use as few parameters as you can get away with This might be surprisingly few Also don t forget to make sure that the config generating run has equilibrated before you start to write out configurations i e nequil should be suitable high i e at least 1000 and probably more Use the plot_umc_energy utility and zmgrace to check this visually The parameters in the Jastrow factor see section 7 2 are of multiple types the u term the x term and the f term together with their associated cutoffs are the ones normally varied If the non linear cutoff parameters are not being optimized it is possible to greatly accelerate the optimization process by turning on linear mode This is done by setting the vm_mode keyword to linear This is expensive in memory however and if this is a problem or if the cutoff parameters are to be included in the optimization then vm_mode should be set to direct NOTE NOT IMPLEMENTED FOR NEW JASTROW FACTOR ONLY FOR THE OLD jasfun data FORM MDT 6 2004 So let s play To see a particular example go and look in examples electron_phases 3D_fluid This is a homogeneous elec
134. is set to go for nmove 250 moves in order to produce nwrcon 250 configs For DMC equilibration we run 100 blocks nblock_dmc equil of 10 moves nmove_dmc_equil For the statistics accumulation we run 200 blocks nblock_dmc stats of 200 moves nmove_dmc stats The DMC nconfig parameter is set to 250 the same as nwrcon The only reason nconfig and nwrcon would ever differ is on a parallel machine config generation is very fast and can be done quickly on far fewer nodes than the DMC calculation Running config generation on batch queues with small numbers of nodes thus usually reduces time spent waiting in queues increase nwrcon in proportion to the reduction in the number of nodes to maintains the same total config number The dtdmc parameter is the DMC timestep note that it is much smaller than the VMC timestep it needs to be small because the DMC Green s function is only exact in the limit as the timestep goes to zero A typical value in a DMC simulation might be 0 003 Note that for accurate work you need to consider extrapolating to zero timestep see the utility extrapolate_tau Note finally that irun is set to 2 indicating that we do a VMC calculation first Now type rundmc and wait till the file 16 DMC_FINISHED appears Amuse yourself by looking at the output files in the o directory Now type reblock in the directory where the dmc hist and dmc hist2 files are reblock also reads input to determine some input parameters it n
135. ities which may be included in the CASINO utils pseudo_converters directory For atomic calculations on radial grids awfn data files generated with the Cambridge pseudopoten tials are available in the on line library 6 2 2 Trial wave functions emergencies If you don t have access to any of these codes and have a specific system in mind then Mike Towler is known for being able to generate trial wave functions in record time on payment of a suitable fee just to get you started he likes in no particular order bright shiny things books about romance yummy things to eat poems authorship on papers for which he hasn t really done any work cute cuddly toys especially bears and money 6 2 3 Input file Having prepared a trial wave function file and perhaps a pseudopotential file you need to tell CASINO exactly what to do with them CASINO takes its instructions from a file called input which contains a flexible list of keywords which control the behaviour of the calculation switch on and off various options and so on Take a moment to examine the various input files in CASINO examples to get a feel for what they look like A complete list of input keywords together with their definitions is contained in section 7 1 Further details including default values may be found by using the casinohelp utility This tends to be more up to date than the manual since it interrogates CASINO directly Type casinohelp all to get a list o
136. ition for Coulomb interactions is 1 dv _ 2059 Hij Y dr r 0 d 1 ae where q and q are the charges in units of the charge of the electron pij m mj m my is the reduced mass and d is the dimensionality The minus sign is used for distinguishable particles e g anti parallel spin electrons or electron and holes and the plus sign for indistinguishable particles e g parallel spin electrons When combined with Eq 96 cusp conditions force Fi WAi 2ui for distinguishable paricles and Fi y A 5 puj for the indistinguishable case For 2D Eq 97 these become Fi Aiz 6u e and Fij Aging respectively The kinetic energy term in the local energy is modified to include the mass N N 1 K Y Kk Y W R VYR 196 i 1 2m i 1 98 Similarly 1 1 VY 1 V w Pi in eh i i 197 i V a SY oe and for the drift vector F 1 VY aan Vi In Wl Ey 198 e The time steps within VMC for the electrons and holes are set independently This can be used to improve the efficiency 10 26 Pair correlation function calculation for electron electron hole sys tems For electron electron hole systems i e when etype is specified CASINO allows the evaluation of the spherically averaged pair correlation function defined by O 2 Joris U R 6 r r r dR Ary2 No No y Joro T 199 where a label o denotes the type electron or hole and spin of a particle r is the position
137. ix D via De Vee 114 where pj is the jth one electron orbital of the spin up Slater determinant and r is the position of the kth spin up electron The transpose of the inverse of D which we call D may be expressed in terms of the cofactors and determinant of D E cof D je det D n The move of electron i changes only the ith column of D and so does not affect any of the cofactors associated with this column The new Slater determinant may be expanded in terms of these cofactors and the result divided by the old Slater determinant to obtain det DT pio neen rew q det det Dtold Za p r Dji 116 If the D matrix is known one can compute q in a time proportional to N Evaluating D using LU decomposition takes of order N operations but once the initial ou matrix has been calculated it can be updated at a cost proportional to N using the formulae T new 1 1 old Di gPa 117 in the case where k i and H1 new old T old new old Dy De Dy X tele Die 118 when k 4 1 10 16 Evaluating the local energy The local energy is given by Noa N N N 2 s E R gt 34 R V Y R 2 V R 23 RIVA UR Vorp R Deel 119 where the terms are the kinetic energy the local part of the external potential energy the non local part of the potential energy the core polarization potential energy if present and the electron electron interaction energy The evaluation of the kineti
138. izable 0 N0 1 YES END SET 2 END CHI TERM START F TERM Number of sets 2 START SET 1 Number of atoms in set 1 Labels of the atoms in this set 1 Prevent duplication of u term 0 N0 1 YES 1 Prevent duplication of chi term 0 N0 1 YES 1 Electron nucleus expansion order 2 Electron electron expansion order 2 Spin dep params 0 gt uu dd ud 1 gt uu dd ud 2 gt uu dd ud 0 Cutoff a u Optimizable 0 NO 1 YES 0 d0 0 Parameter values Optimizable 0 N0 1 YES END SET 1 START SET 2 Number of atoms in set 4 Labels of the atoms in this set 2345 Prevent duplication of u term 0 N0 1 YES 1 Prevent duplication of chi term 0 N0 1 YES 1 Electron nucleus expansion order 2 Electron electron expansion order 2 Spin dep params 0 gt uu dd ud 1 gt uu dd ud 2 gt uu dd ud 0 Cutoff a u Optimizable 0 NO 1 YES 0 d0 0 Parameter values Optimizable 0 N0 1 YES END SET 2 END F TERM END JASTROW 33 There are two additional terms which may be of benefit in periodic systems but which in practice are very rarely used insert these between END F TERM and END JASTROW above START P TERM Spin dep 0 gt uu dd ud 1 gt uu dd ud 2 gt uu dd ud 1 Number of simulation cell G vectors NB cannot have both G amp G 6 G vector in terms of rec latt vects label 0 0 1 1 0 1 0 1 1 0 0 1 1 1 1 1 1 1 0 2 1 0 1 2 Parameter value Optimizable 0 NO 1 YES END P TERM START Q TE
139. l wave function given in the xJwfn data file by using runvmc to run CASINO with IRUN set to 5 in the input file Note that the density data file is also required for this The format of the file should be obvious from looking at some of the examples NB it is not necessary to regenerate this file for different supercell sizes provided the lattice vectors in eepot data are obtainable from the lattice vectors in x wfn data through scaling by a single parameter 39 7 6 Trial wave function files awfn data bwfn data gwfn data pwfn data swfn data These files contain the data defining the determinant part of the trial wave function The data can be given in a Gaussian basis set gwfn data in a plane wave basis set pwfn data in a blip function basis set bwfn data or in splines swfn data The awfn data file contains a trial wave function for an atom with the orbitals given explicitly on a radial grid Without going into details of specific formats the files basically contain the following information Basic info about the trial wave function generating calculation e g DFT HF etc including total energy and components Geometry of the system Details of the k space net used by the program that generated the trial wave function only if the systems is periodic Details of the basis set Exponents contraction coefficients shell types Gaussians G vectors of the plane waves Blip grid Multi dete
140. le containing personalized settings which can be accessed through the QMCID environment variable see note about this in Supported architectures above 2 Now you need some run scripts to run the code If your machine is a workstation multiprocessor PC cluster which does not use a batch queue system then you can just use the standard runumc runvarmin and rundmc from the CASINO scripts workstation directory You will need to add a few lines to the CASINO utils Makefile so that these scripts are copied from this directory if QMC_ARCH is set to the new machine type If your machine is a parallel machine using a batch queue system then you will probably need to adapt existing run scripts Obviously you can use the ones in the origin t3e hitachi ibm_sp3 etc directories as models Running CASINO without these scripts is possible but not really recommended since the scripts have been carefully tuned to check for most common errors to run CASINO in as efficient a way as possible and to tidy up after it 3 You may also need to do a little further hacking around in the CASINO utils Makefile in order to get all the utilities working properly Note that many of these utilities require you to be able run inter active jobs which may not be impossible on certain batch based machines so make sure the Makefile doesn t bother to compile and set them up Have a look at the CASINO utils help casinohelp script too to make sure it behaves appropriately
141. le xx_pp data where xx is the symbol for the element in question in lower case This file contains the different angular momentum components of the pseudopotential given on a radial grid and some auxiliary information in the following format LSDA Pseudopotential in real space for Si Atomic number and pseudo charge 14 4d0 Energy units rydberg hartree ev rydberg Angular momentum of local component 0 s 1 p 2 d 1 NLRULE override 1 VMC DMC 2 config gen 0 gt input default value 00 Number of grid points 2476 R i in atomic units 0 000000000000000E 000 7 178690774405650E 012 39 6567798705125 40 0553371342383 r potential L 0 in Ry 0 000000000000000E 00 0 131424063731862E 10 8 00000000000000 8 00000000000000 r potential L 1 in Ry 0 000000000000000E 00 0 574393968349227E 10 8 00000000000000 8 00000000000000 r potential L 2 in Ry 0 000000000000000E 00 0 128711823920867E 09 8 00000000000000 8 00000000000000 Core polarization terms 0 1650 0 5446 alpha a u rbaree a u usually set to 0 5 rbars rbarp 0 5216 0 5676 0 7172 rbar for s p d a u The NRULE parameters control the grid on which the angular integration is performed see Section 10 18 and are normally specified in the file input The values given in input are the default for all atoms in the system but they can be overridden for particular elements by setting the parameters in this file to non zero values Unless you e
142. lectron systems coefficients for spin up pairs coefficients for antiparallel pairs if spin dependence is 1 or 2 coefficients for spin down pairs if spin dependence is 2 For each spin pair the number of parameters is equal to the expansion order For electron hole systems a spin dependence value of 0 means that the same parameters are used for all spins and particles 1 means that different parameters are used for electrons and holes 8For electron systems the parameter values are given for spin up electrons then if the spin dependence is 1 for spin down electrons For each spin type the number of parameters is equal to the expansion order 35 e For p a list of simulation cell G vectors must be provided These are specified in terms of the reciprocal lattice vectors Only one out of each G and G should be specified The same parameter value is used for G vectors with the same label e For q a list of primitive cell G vectors must be provided Again G vectors with the same label have the same parameter value It is possible to specify a negative relationship between parameter values by using a negative label for the appropriate G vectors 36 7 3 Pseudopotential file xx_pp data CASINO can carry out all electron or pseudopotential calculations but it is normally advantageous to replace the core electrons by a pseudopotential CASINO will automatically treat an atom as all electron unless there exists a pseudopotential fi
143. lip_input file In the blip_input file we need to specify the multiplication factor which defines the fineness of the blip grid we suggest using values greater than one Multiplication factors of at least two may be needed in some cases in order to get accurate VMC energies but this is less of an issue in a DMC calculation The second entry in the input file is a flag which if set to true allows some testing of the quality of the blip representation see Section 15 2 for more details The third entry is the number of localized states If this is lower than the total number of states than the remainder will be treated in the usual non localized way and the blip grid extends over the whole system for those states The fourth entry is a starting guess for the radius of the spherical localization region and the fifth is the percentage of the norm of the orbital contained in the localization region Finally the sixth entry is the thickness of the shell region over which the orbital is smoothly truncated to zero Type blipl to run the utility to produce a blwfn data file 116 16 Using CASINO with other supported programs 16 1 TURBOMOLE See the README file in CASINO utils wfn_converters turbomole 16 2 CASTEP 2002 See the README file in CASINO utils wfn_converters castep 16 3 PWSCF See the README file in CASINO utils wfn_converters pwscf 16 4 ABINIT See the README file in CASINO utils wfn_converters abinit 17 Using CASINO with unsuppor
144. ll periodic potential of the simulation cell is obtained by following a procedure analogous to that described for the 3D case with the self term given by IN II B vete a 151 II 152 R0 G0 The first derivatives of the 2D Ewald potential required for the evaluation of the core polarization energy in 2D slabs are different in directions parallel and perpendicular to the plane of the slab The x and y derivatives are given by OPCs y R oe 10 A rF r R F E R exp y r rj wp 153 Sni 2y3 7 R 2 A en lexp C erte G nt exp 2G erfe where z or y and the z derivative is given by vP r r _ z erfe y2 r r R I 271 ee Oz 7 2 r rj R lr 1 R t E exp y r rj R i 1 2 cae oe lexp Gente 1 exp 2G erfc i aat G40 23 T am k a A Note finally that in things such as 2D bilayer systems electrons in one layer holes in the other say there is an additional capacitor term due to interaction of the backgrounds about which I shall write more here in a minute CASINO does not evaluate this term 10 20 3 1D Coulomb interaction Coulomb sums in systems which are periodic in one dimension the x direction according to CASINO are performed using an algorithm based on the Euler Maclaurin summation formula See for example Eq 4 8 in Comp Phys Commun 84 156 1994 More details to appear here la
145. llowed by USE ONLY s in alphabetical order USE a USE b USE c USE ai ONLY x USE b1 ONLY x USE c1 ONLY x e Don t use return at the end of each routine just use END SUBROUTINE or END FUNCTION etc 20 2 CONTENT Important policy statement Following advice from ABINIT chief developer Xavier Gonze and in order to reduce the workload on Mike and Neil the policy of the Cambridge group regarding inclusion of routines written for CASINO is now as follows Such modifications will not be accepted under any circumstances unless the modifications are made to the most recent version of the code email Mike when commencing the merge so that the most recent version does not change too much during this process The reason for adopting this policy is simply that merging complex code written by someone else into a program as large and complex as CASINO is a very time consuming procedure and is prone to introducing errors This process can and has taken months with previous submissions Even if modifications are made to the latest version you may find that they are still not incorporated into the public release If you want this to happen there are a variety of ways to lower the energy barrier for the inclusion of your routines The general theme here is to try to reduce the amount of effort that Mike has to perform to validate and check your code If your submission includes lots of things from the Good things list and no
146. ls that are not truncated must be listed on the next Ng lines If Ng lt 0 then no orbitals will be truncated The eighth line holds the minimum permitted thickness in a u of the shell region between the inner and outer truncation radii for the smooth truncation of localized orbitals 18 Assuming Np 0 114 9 The ninth line controls the amount of output from 1 minimal to 4 for debugging Tf it is set to 2 or above then useful information about the kinetic energy of the truncated and untruncated orbitals is provided The difference between the total kinetic energies of the truncated and untruncated orbitals occupied appropriately gives the order of magnitude of the likely biasing due to the discontinuities in the trial wavefunction resulting from abrupt truncation It should be ensured that this bias is much less than the desired QMC error bar The splined truncated orbitals can be generated independently Correspondingly generate_spline can be run in parallel if desired It uses the MPI library The output of generate_spline is a swfn data file containing the splined orbitals and geometry information If the verbosity level is set to 3 or above then information about the location of the truncation radii will be placed in output files for each processor If requested the plots of the orbitals are placed in files with names such as plotx 003 containing a plot along the z direction through the centre of orbital 3 The swfn
147. luctuates as the simulation progresses There is a feedback process in operation to limit the population fluctuations so this should just oscillate around the total initial number of configurations that we chose 1280 in this case If you were smart and printed this manual out in colour then in the bottom panel you will see a red line a green line and a wobbly black line The black line is the instantaneous value of the local energy the thing which gets averaged over the red line is the reference energy Er which is adjusted to control the feedback process that keeps the population in check the green line is the best estimate of the DMC energy as the simulation progresses You should note that in the picture above the best estimate of the energy falls from its initial value from the VMC config generation run to a much lower constant value around 55 72 a u as the wave function evolves to the ground state You should look at what point the energy becomes roughly constant i e at around 500 moves in this case This splits the graph quite neatly into an equilibration phase and a statistics accumulation phase When we average energies to produce the final energy and error bar we should only include those moves between 500 and the end of the run To see DMC in action go to examples molecules h2 RHF dmc and let s calculate the DMC energy of a hydrogen molecule The input files should already be set up correctly We see that an equilibrated VMC run
148. made to design an algorithm for optimizing the timestep by directly minimizing the correlation time In the end however it was found that this required too many equiibration steps to get a clear min imum and there was not much benefit over the simpler procedure of choosing a 50 acceptance ratio OPT_MAXEVAL Logical Maximum number of evaluations during optimization default 200 OPT_MAXITER nteger Largest permitted number of variance minimizer iterations default 150 ORBBUF Logical Setting ORBBUF TRUE turns on DMC orbital buffering This is an efficiency device in which buffered copies of orbitals gradients Laplacians are kept for later reuse This has a high memory cost Orbital buffering should always be used unless you start running out of memory hence the ability to turn it off ORB_NORM Real Allows user to change normalization of orbitals by multiplying all of them by this constant Of course this should have no effect on the energy but it can be useful if the Slater determinant starts going singular as it might for example for some very dilute low density systems PAIR _ CORR Logical Set PAIR CORR to T to activate calculation of the pair correlation function Currently restricted to periodic systems Note you also need to give values for the FIXED E POS block Spherical averaging may be performed see the PCF_SPH_MODE keyword PCF_SPH_MODE Integer Real space accumulation on radial bins around a fixed electron basi
149. maintain a reasonably steady population However this procedure can result in a bias in the estimate of expectation values especially for small populations We evaluate expectation values using the method of UNR 5 which attempts to eliminate bias due to population control Using the label m for time step Eq 10 becomes f R m J Cono R T f R m 1 dR 52 Clearly in the absence of the accept reject step the effect of including the time step dependent reference energy Er m in Gpuc can be undone by multiplying the right hand side of Eq 52 by exp TEr m In a similar fashion the effect of including the reference energy from the previous time step can be eliminated by multiplying by exp TEr m 1 When the accept reject step is 13In order to sample w from g2 w we sample the polar angle uniformly on 0 7 the azimuthal angle uniformly on 0 27 and the magnitude w from 4 w exp 2 w This is achieved by sampling 71 r2 and r3 uniformly on 0 1 and setting w log rir2r3 2C see reference 6 for further information 72 present we can approximately undo the effect of the reference energy by using our best estimate of the effective time step Tarp See Section 10 6 in the undoing factors Continuing this process we may eliminate the effect of changing the reference energy from f m by multiplying it by TI m II exp TerrE7r m m 53 m 0 where in principle the product run
150. meters By contrast in the new optimization method the quartic expansion coefficients of the unreweighted variance are accumulated directly in VMC there is no need to write out a set of configurations Furthermore when the unreweighted variance referred to as the least squares function LSF is evaluated during the subsequent opti mization stage there is no need to repeatedly sum over a set of configurations the quartic LSF can be evaluated directly The fact that the LSF can be evaluated as a quartic function of the parameters gives two significant advantages over the standard variance minimization algorithm i the LSF can be evaluated extremely rapidly typically thousands of times per second furthermore the CPU time required is independent of the system size and ii the LSF along any line in parameter space is a simple quartic polynomial so that the exact global minimum of the LSF along that line can be computed analytically The method has two drawbacks i only linear Jastrow parameters can be optimized in this fashion and ii the number of quartic coefficients to be evaluated and stored in memory grows as the fourth power of the number of parameters to be optimized 10 24 2 Using the new optimization method A CASINO variance minimization calculation using the new optimization method is carried out in exactly the same way as an ordinary variance minimization calculation except that 1 The use_newopt keyword in input should
151. mize all parameters flag in heg data to FALSE The second pair of Pad coefficients in the file are then ignored If Pad coefficients are to be optimized direct mode must be used The process of variance minimization requires careful monitoring by the user Unless proper action is taken the correlated sampling procedure may become numerically unstable particularly for large system sizes The characteristic of these instabilities is that during the minimization procedure a few configurations often only one acquire a very large weight w a The estimate of the variance is then reduced almost to zero by a set of parameters which are found to give extremely poor results in a subsequent QMC calculation One can overcome this instability by using more configurations Various alternative ways of dealing with this instability have been devised One method is to limit the upper value of the weights 15 or to set the weights equal to unity 17 14 In our calculations for large systems we normally set the weights to unity which is achieved by setting vm_reweight to FALSE It can be verified that the effect of this approximation is in general negligible by regenerating configurations and carrying out a new variance minimization process On the other hand if the initial trial wave function is poor then some of the configuration genera tion variance minimization cycles can be bypassed by using the weights There is also some eviden
152. mploy Gaussian basis sets The more important point is that the use of a Slater basis in which each basis function obeys the cusp conditions at the nucleus on which it is centered does not automatically lead to a full solution of the electron nucleus cusp problem This is due to the contributions from the tails of basis functions centered on other nuclei which in general prevent molecular orbitals expanded in that basis from obeying the electron nucleus cusp condition Manten and Liichow 16 have developed a scheme for applying cusp corrections in QMC calculations but as it similarly relies on correcting individual atom centered basis functions it is clearly not a full solution An alternative solution to the cusp problem might be to enforce the electron nucleus cusp condition using the Jastrow factor This is feasible and we have implemented it but we found this to be unsatisfactory because a very large number of variable parameters are required in order to obtain a good trial wave function The solution we have finally adopted in CASINO involves the direct modification of the molecular orbitals so that each of them obeys the cusp condition at each nucleus This ensures that the local energy remains finite whenever an electron is in the vicinity of a nucleus although it generally has a discontinuity at the nucleus 10 13 1 Electron nucleus cusp corrections The Kato cusp condition 42 applied to an electron at r and a nucleus of charge Z
153. ms can be written as a sum of the mass velocity term Darwin terms and the retardation term The mass velocity term 2 arises from the relativistic variation of mass with velocity and is written as 1 aea NV gt vi v 207 i This term is proportional to the square of the nonrelativistic kinetic energy The spread of electronic charge is described by the electron nucleus and electron electron Darwin terms 3 4 which are expressed together as S Vivi 2v E3 E4 3 4 Ig x y 5 208 E e Ey Ta with Z as the nuclear charge The last term known as the retardation term es arises from the interaction between spin magnetic moments which are not mutually penetrating It is given by the expression i lt j ij Tij where rj r rj Calculations for the beryllium atom 46 show that the total relativistic correction to the energy is approximately 0 00239Ha with the mass velocity term having the greatest contribution of 0 0145Ha followed by the Darwin terms of 0 0119Ha 100 I 2 H He 1 00794 4 00260 3 4 5 6 7 8 9 10 Li Be B C N O F Ne 6 941 9 012187 10 811 12 0107 14 00674 15 9994 18 99840 20 1797 11 12 13 14 15 16 17 18 Na Mg Al Si P S Cl Ar 22 98977 24 3050 26 98154 28 0855 30 97376 32 066 35 4527 39 948 19 20 21 22 23 34 25 26 27 28 29 30 31 32 33 34 35 36 K Ca Sc Ti V Cr Mn Fe Co Ni Cu Zn Ga Ge As Se Br Kr 39
154. n electron n 1 m Number of radial grid points are given Distance from centre of atom r 47 O O O 108431746952108E 04 113045182349640E 04 117854905151603E 04 Orbital 1 1s 010 000000000000000E 00 659528705936585E 01 687054582044431E 01 115705591983797E 01 oooo O O O 000000000000000E 00 000000000000000E 00 000000000000000E 00 Orbital 2 2s 020 000000000000000E 00 120187930416515E 01 125203784849961E 01 130424628552784E 01 Ooooo oooo 000000000000000E 00 000000000000000E 00 000000000000000E 00 000000000000000E 00 spin type O unpolarized 1 up 2 down n 1 r Value of trial wave function at point r 48 7 7 heg data file The heg data file is used to define the trial wave function used in calculations of electron and electron hole systems It is the equivalent of the x wfn data file used for calculations of molecular crystalline systems but it is entirely self contained and does not need to be generated by an external program See CASINO examples electron_phases and CASINO examples electron hole phases for practical examples with different wave function types 2D 3D fluid crystal excitonic insulator NB in the case of the jellium slab system ETYPE 6 the heg data file takes on quite a different form see the next section The basic data required are the dimensionality 2D or 3D the Wigner Seitz radius parameter rs and the translation vectors of
155. n fall Since energy and variance minima do not necessarily coincide it is possible for it to rise slightly however if it rises substantially then it suggests there may be a problem such as the use of too many parameters e Do the parameters themselves change values substantially in successive optimizations when the energy is not being lowered If yes then they may well be redundant 10 24 Wave function optimization new method 10 24 1 Background Consider the linear parameters in CASINO s new Jastrow factor the expansion coefficients a Gm Yimn aa and bg in the u x f p and q terms respectively The local energy of a single configuration can be shown to be a quadratic function of the linear parameters hence the variance of the local energies of a fixed finite set of configurations is a quartic function of the parameters But this is precisely the quantity that is minimized in an unreweighted variance minimization calculation The process of variance minimization can therefore be greatly simplified and accelerated if only linear Jastrow parameters are to be optimized In an ordinary variance minimization calculation the VMC method is used to generate a set of configurations distributed according to the initial trial wave function During the optimization process the variance of the local energies of this configuration set is computed for different sets of parameters and the variance is minimized with respect to the para
156. n method must be changed by specifying IOp 5 51 1 In large CASSCF calculations where very many determinants are involved Gaussian currently only prints the first fifty determinants in the expansion those with the largest coefficients The significance of the truncation may be judged by looking at the sum of the squares of the coefficients that gaussiantogmc outputs when it reads the wave function Unfortunately this truncation prevents a full HF test 7 e running the wave function in VMC without a Jastrow factor and checking that the result agrees with that of Gaussian but the energy returned by such a test should be above that reported by Gaussian As well as this truncation to fifty determinants Gaussian has a formatting error which means that if the index number of a determinant is greater than 99 999 then it is replaced by stars and is thus useless gaussiantogmc deals with this by simply throwing away such configurations which further truncates the expansion 108 e Gaussian switches to direct mode for large CASSCF calculations and in doing so automatically stops printing the definitions of the Slater determinants used in the calculation In order to reconstruct the wave function we do of course need to know what the Slater determinants are Gaussian may be persuaded to print them by using IOp 4 46 3 IOp 4 21 10 IOp 4 21 100 The first two of these both tell Gaussian to use Slater determinants I specified both to be on the sa
157. nding direction in reciprocal space So the easiest way to increase the fineness of the real space grid is to enlarge the cuboid by including extra reciprocal lattice vectors with zero coefficients The parameter on the second line is the factor by which the cuboid is enlarged in each direction hence for example doubling it has the effect of doubling the fineness of the real space grid in each direction giving eight times as many grid points The third line specifies whether the truncation surfaces of the localized orbitals are spherical 1 or cubic 2 Note that smooth truncation can only be applied if the truncation surfaces are spherical The fourth line allows the user to choose the method for determining the truncation radii If this line is set to 1 then the truncation radius for each orbital will be chosen so that a certain percentage of the square of the orbital norm lies inside it if set to 2 then the user specifies a fixed radius for all the orbitals The fifth line holds the desired fraction of the norm squared to be contained within the inner truncation radius It is ignored if the previous line is set to 2 The sixth line holds the fixed inner truncation diameter in a u if the fourth line is set to 2 If the fourth line is set to 1 then the sixth line is used as an initial guess at the inner truncation diameter The seventh line holds the number of orbitals Ng that are not to be truncated If Ng gt 0 then the orbita
158. ne at all from the Bad things list then inclusion is semi automatic and will normally take place within days We are aware that it requires a lot of work to do the things on the Good list so sulk about it but we remind the reader that we will have to do it if you don t 20 2 1 Good things e Well written absolutely standard Fortran90 in the standard CASINO format Do not use for tran95 or later constructs or non standard but useful extensions like flush e The results of extensive testing of your routine with a range of examples with before and after numbers including timings e A detailed description of whatever your code is supposed to do I shouldn t have to go through each line of your code at the same detailed level you did or I might as well have written it myself 125 A bit of TeX for inclusion in the manual if the code requires the user to know something over and above what the current manual describes e Some example input output which demonstrate the new capabilities that Mike can just run and look at without having to invent his own examples e Evidence of testing on both single processor and parallel machines Remember only the master node should be writing to the output file Also CASINO is supposed to compile and work on single processor machines without an installed MPI library Use the fake comms_serial f90 module supplied both in the utils directory and with the main code e You should no
159. ng positron problems contact MDT if you want this to be implemented In this section the changes required to the CASINO code and to the basic equations in the presence of holes are discussed These largely stem from the possibility of having a variable mass ratio between the positively and negatively charged particles The basic differences are 97 The diffusion Green s function Eq 11 becomes 1 Re R 27D V R Gp R R 7 5 e ene p R R 7 4r Der 3Ne 2 oo 4D 1 Ra Ri 27D Va R y 1 A Dpr 3 N2 P 4DpT oy where e and h denote electron and hole quantities Ne and N are the numbers of elec trons and holes the diffusion constants are defined as De st and D 27 where me and mz are the electron and hole masses in atomic units i e in units of the mass of the electron When particle is moved Eq 15 becomes r r x 2D 7v i R 191 where x is a 3 dimensional vector of normally distributed numbers with variance 2D 7 and zero mean The probability of accepting this move Eq 19 is then Y R pi min 1 exp TDR AA vR FS 192 The effective time step Eq 28 is given by Zi MipiAra i Teff Q m mi r 193 The drift vector limiting Eq 29 takes the form 14 yTH4Dalw Pr AAA y 194 vi 2a v Dit Separate Jastrow factors must be defined for the electron electron hole hole and electron hole interactions The general form of the cusp cond
160. nite systems 0 Density not written out INEW Integer Determines whether this is a new run or a continuation of a previous one 1 VMC Start a new run DMC Read old VMC configurations and recalculate the current best estimate of the ground state energy as the mean of the configuration energies 0 VMC Continuation of old run Read in old configurations DMC Continuation of old run Do not recompute the current best estimate of the energy INPUT_EXAMPLE Logical If INPUT_EXAMPLE is TRUE then write out an example of a QMC input file with all currently known keywords and their default values A modified version of this can be used as an input file in future runs IRUN Integer Determines type of calculation Possible values 1 VMC run with wave function of determinant s only VMC run with a Slater Jastrow wave function DMC run Variance minimization run for optimizing a trial wave function Generate eepot data file for geometry in xwfn data x p g for MPC interac tion See Section 10 20 4 6 Generate density data file for geometry given in xwfn data x p g ote ow ISOTOPE_MASS Real Nuclear mass override for the relativistic corrections activated by RELA TIVISTIC This keyword can be used to define the nuclear mass in amu if you need to override the default value which is averaged over isotopes according to their abundances The default given in the table in Section 10 27 is used if ISOTOPE_MASS is set
161. ns of the code compiled with different compiler options which are automatically stored by the Makefile in different bin directories Note that on a batch machine the script will create a batch queue file for the first job in the list show it to you and ask you whether it should submit it You should answer y or n Then if you have requested that the script perform some action in the future when the job it has just submitted has finished then the script will stay alive after submitting the job and demand to be put in the background Do this by holding down the control key and pressing z to suspend the script followed by the command bg to put it in the background The script will lurk and monitor the queues every so often to check when the current job has finished then perform the appropriate actions to submit the next job in the sequence and so on The following notes try to demonstrate typical usage of these scripts on various types of machine another since the unformatted file can usually only be read by the machine that generated it Run format_spline in the presence of swfn data formatted file to produce the corresponding unformatted file 4Note the following in the input file 1 We have set isperiodic to false because we are simulating a molecule Now that we have generated truncated localized orbitals we can dispense with periodicity 2 We have provided lists of the occupied states for spin up and spin down electr
162. nteraction because CASINO can evaluate it much more quickly than the Ewald interaction and also because it gives rise to smaller finite size effects More details are given in section 10 20 4 10 6 2 6 How to do a VMC calculation Let s begin by calculating the Hartree Fock energy of a hydrogen atom Change directory to CASINO examples atom hydrogen You will see a gwfn data file generated by a GAUSSIAN94 calculation see the bottom of gwfn data for the GAUSSIAN output and a CASINO input file No pseudopotential file is supplied so CASINO will assume you wish to do an all electron calculation No jastrow data file is supplied because the correlation energy in a one electron atom zero is not difficult to calculate without one Look in the input file You will see that neu and ned have been given the correct values one spin up electron present isperiodic is false and as this is a finite system the npcells block is not required It is not necessary to equilibrate the electron distribution since there is only one electron but nequil is set to 1000 to remind you that this should normally be done The nmove parameter is set to 100000 The corper parameter is set to 3 to reduce serial correlation Type runumc Three files out vmc hist and vmc posout should be produced First type ve This is a quick way to pull the VMC result out of the output file It says Total energy 0 500009481495 0 000046640057 Time 2 9000 seconds Note t
163. o r a in turn according to 65 ri r X TVs 11 Ei 1 1 15 where y is a three dimensional vector of normally distributed numbers with variance 7 and zero mean v R denotes those components of the total drift vector V R due to electron i Hence each electron is moved from r to r with a transition probability density of 1 ri r TV T1 T 1 T Ty Eli Ti ri S rh hippe Ny SHE exp Ln a oe aN 16 For a complete sweep through the set of electrons the transition probability density for a move from R r r y to R r1 ry is simply the probability that each electron i moves from r to r So the transition probability density for the configuration move is N T R R J tn wee Te ri HU Tay Ty 17 i l In the limit of small timesteps the drift velocity V is constant over the small configuration move Evaluating the product in this case we find that the transition probability density is T R R Gp R lt R 7 18 so that the drift diffusion process is described by the drift diffusion Green s function Gp At finite timesteps however the approximation that the drift velocity is constant leads to the violation of the detailed balance condition We may enforce the detailed balance condition on the DMC Green s function by means of a Metropolis style accept reject step introduced by Ceperley et al 33 This has been shown to greatly reduce timestep errors 1
164. ode for shell type s 1 sp 2 p 3 d 4 f 5 etc harmonic representation d 4 f 5 cartesian representation NUMPRIMS_IN_SHELL NUM_SHELLS Number of primitive Gaussians in each shell FIRST_SHELL NUM_CENTRES 1 Sequence number of first shell on each Gaussian centre Allows e g do n 1 num_centres do shell first_shell n first_shell n 1 1 blah enddo enddo to loop over shells on each centre Note dimension EXPONENT NUM_PRIMS exponents of Gaussian primitives C_PRIM NUM_PRIMS correctly normalized contraction coefficients C_PRIM2 NUM_PRIMS must be omitted if no sp shells in basis correctly normalized 2nd contraction coefficients i e p coefficient for sp shells zero otherwise SHELL_POS 3 NUM_SHELLS positions of shells not necessarily atom centred MULTIDETERMINANT INFORMATION GS Ground state calculation 42 SD Single determinant inc Example SD DET11PR2151 or DET 1 2 PL6 1 or DET 1 2 MI 6 1 or MD Multiple determinants MD 3 1 d0 2 d0 1 d0 DET31PR2151 ORBITAL COEFFICIENTS excitations additions subtractions Single det calculation Promote electron in band 2 kpoint 1 to band 5 kpoint 1 in determinant 1 spin 1 up Add electron spin 2 down to det 1 band 6 kpoint 1 Remove the added electron Multideterminant 3 determinants Determinant 1 prefactor Determinant 2 prefactor Determinant 3 prefactor Promotion as examples above amp CK NUM_RE
165. oing research with CASINO is to generate a trial wave function using for example a DFT or Hartree Fock calculation multide terminant quantum chemistry methods can also be used This must be done using an external program which must in the past have been persuaded to write out the wave function in a format that CASINO understands either all by itself or through the transformation of the standard output of the program using a separate utility which will probably be included in the directory CASINO utils wfn_converters Note that writing out the wave function basically means writing out the geometry the basis set and the orbital coefficients The information defining the trial wave function generated by the external program lives in files whose name depends on the basis set in which the orbitals in the determinantal part of the wave function are expanded They are called gwfn data Gaussians pwfn data plane waves bwfn data blip functions swfn data splines or awfn data atomic orbitals given explicitly on a radial grid These files will often be referred to generically with the name xwfn data Choice of basis set has been found to depend largely on personal prejudice though some consideration should be given to issues of computational efficiency Gaussian and plane wave functions are generated directly Blip and spline wave functions are generated by post processing a plane wave pwfn data file using a CASINO utility this is desir
166. on algorithm Hence the sampling of configuration space in an electron by electron algorithm is more efficient After each move of each electron we check whether the configuration has crossed the nodal surface by checking the sign of the Slater part of the trial wave function If it has then the move is rejected This has been found to be the least biased method of imposing the fixed node approximation 23 10 5 Branching and population control The branching Green s function can be implemented by altering the population of configurations and or their weights At the start of the calculation one chooses a target population Mo and the actual population Mix m see Eq 27 is controlled so that it does not deviate too much from Mo The population control is principally exerted by altering the reference energy Er m Large changes in Er can lead to a bias and therefore it is varied smoothly over the simulation For each move of all the electrons in configuration a the branching factor is calculated as Mulas exp 5 Som Sau Erl rosa 23 where Teg is the effective time step for configuration a at time step m see Section 10 6 S is the local energy we denote it by S because it is usually a modified version of Ey see Section 10 6 Unless weighted DMC is used i e unless LWDMC TRUE the number of copies of this configuration that continue to the next time step is given by M a m INT n Mila m 24 where 7 is
167. on move 2 Evaluate the energy during the move by averaging electron energies over the old and new positions or 3 Like 1 but propose entire configuration moves instead of single particle ones at the accept reject stage Method 1 is the default VM_DERIV_BUFFER Logical Setting this flag to T will speed up variance minimization by buffering the kinetic energy terms associated with the different factors present in the wave func tion This is useful because NL2SOL spends most function evaluations in computing numerical derivatives hence changing only one parameter at a time which affects only a portion of the wave function It is currently possible to buffer all five terms in the new Jastrow function and the Slater determinant For the old Jastrow this feature is not yet implemented VM_E_GUESS Real If VM_USE_E_GUESS is true then VM_E_GUESS should be supplied as an estimate of the ground state energy VM_FIXNL Logical Fix nonlocal in VARMIN Currently only TRUE is allowed VM_FORGIVING Logical Do not whinge about calculated energies not agreeing with those read in Recommended to be set to FALSE VM INFO Integer Controls amount of information displayed during variance minimization Takes values 1 Display no information 2 Display energies at each function evaluation 3 As 2 but calculate weights as well 4 Write out configs and their energies etc as they are read in VM_JASCHECK Logical Numerical check of whether v
168. on of the operating system Together these comprise the total CPU time Note that the User time is usually significantly greater than the System time Indeed a high ratio of system to user CPU time may indicate a problem Exceptions such as page faults misaligned memory references and floating point exceptions consume a large amount of system time Note that time spent doing things like waiting for input from the terminal seeking on a disk or rewinding a tape does not appear in the CPU time these situations do not require the CPU and thus it can and generally does work on something else while waiting The elapsed wallclock or real time is simply the total real time that passed during the execution of the program This will typically be greater than the total CPU time due primarily to sharing the CPU with other programs In addition programs that perform a large number of I O operations require more memory bandwidth than is available on the machine i e the CPU spends significant time waiting for data to arrive from memory or is paged or swapped will show a relatively lower ratio of CPU time used to elapsed time 20 4 BUG REPORTS We cannot of course guarantee that CASINO is free of bugs and if you find one please tell us If you don t know what the problem is yourself than please include as much detailed information you can about the nature of the error in order to ensure a quick fix 20 5 REQUESTS FOR NEW FEATURES
169. ons the statelist_up and statelist_down blocks 18 Workstations parallel machines without batch queues runvmc runvmc big runvmc nnodes 4 Do a VMC calculation The big option will create the potentially large vmc hist files in a pre defined scratch space and create a link to them There is not usually a n option to chain multiple jobs because there is no time limit on a typical workstation The nnodes flag says how many nodes to use on parallel machines runvarmin n 2 v runvarmin n 2 v nnodes 16 Run 2 iterations of a config generation variance minimization cycle 4 CASINO jobs al together The v flag tells it to run a post optimization VMC calculation with the final optimized Jastrow The optional nnodes flag says how many nodes to use on parallel machines rundmc rundmc c e s rundmc s rundmc nnodes 32 Run a DMC calculation Default is to do one config generation one DMC equilibration one DMC stats accumulation though this may be changed by supplying the c e or s flags explicitly Let s say you have completed a DMC calculation but the standard error is too large You can accumulate more statistics on top of the ones you already have by typing rundmc s Remember to do a DMC calculation with the rundmc script you need to setup your input file to do a VMC config generation All parameter changes to run DMC equilibration and stats accumulation will then be performed automatically by the script The op
170. ons of atoms Cartesians in atomic units 57 0 8425817327 0 8425817327 0 8425817327 Atomic and position 6 0 8425817327 0 8425817327 0 8425817327 Atomic and position Number of G vectors 6861 Number of G vectors G vectors Cartesians in atomic units 0 000000000000 0 000000000000 0 000000000000 List of G vectors 13 049860754934 0 000000000000 13 049860754934 Charge density from SC calculation 8 00000000000000 List of rho G 2 020855523387648E 006 BOOLEAN for QMC charge density to follow TRUE BOOLEAN for QMC charge density 8 of primitive cells in simulation cell 2 99883837683579 of moves 32 Number of up spin electrons 32 Number of down spin electrons 6861 Number of G vectors 191 925656117491 List of QMC rho G To get the correct correct normalization divide by the of moves and number of primitive cells 9 16791786818292 Notes 1 The utility d2rs can be used to convert the charge density in reciprocal space into a real space density on a grid This can be converted into a format readable by gnuplot using the plot_density utility 38 7 5 MPC file eepot data This contains the Fourier transform of the function f r defined in Section 10 20 4 It is required when using the Modified Periodic Coulomb MPC interaction to compute the electron electron interactions Before doing a QMC calculation with the MPC interaction this file should be generated from the tria
171. or not the move is accepted we can evaluate the contribution to the average from this particular electron as E 1 Pacc old PaccEnew 4 Calculating quantities using Equation 4 allows the contribution of configurations with high energies but low probabilities such as when two electrons come close together to be more readily included in averages fluctuations in the system energy caused by the occasional acceptance of such unlikely configurations are reduced decreasing the error bars on the estimated quantities 19 On the other hand in the limit of perfect importance sampling that is where the trial wave function is exact and so the local energy is constant in configuration space the standard Metropolis 9By the energy of an electron we mean the sum of its kinetic energy electron ion potential energy and half the sum of its electron electron interaction energy Summing these over all electrons we get the total energy of the configuration 10Note that an electron energy in the existing configuration is not the same as the electron s energy at the end of the previous sweep some of the other electrons will almost certainly have been moved in the intervening period Only the electron ion component of energy ignoring core polarization effects in pseudoatoms will be the same 62 algorithm gives the correct result with zero variance Evaluating the electron energies as we sweep through the set gives a result with non
172. or numsrt but returns an associated sign change given by multiplying by 1 for each exchange Stores the alpha and beta eigenvector coefficients separately and multiply in remaining normalization factors which dif fer between dyz dz2_y2 etc MODULE contains define pi and also defines conversion factors for Ha to eV and Bohr to Angstrom Evaluates an MO of a given spin at a specified point in space 110 Routine Purpose qmc_write resum read_G9xout read_fchk rejig resum_cas set_parameter_values shell_centres sum_degen_excite user_control win_construct win_test Writes the gwfn data file Resums a CIS TD DFT wave function Reads the output file produced by the Gaussian job Gets the nuclear nuclear potential energy CIS TD DFT CAS expansion if present and HF eigenvalues Reads the formatted checkpoint file produced by the Gaus sian job Gets the MOs etc MODULE contains numsrt and hence prototypes it which is necessary because it has an optional argument Partially resums a CAS expansion Sets the value of Pi and related constants Identifies the positions of the distinct shell centres and store the first shell index corresponding to each Called by analyse_cis Loops over excitations out of a given range i j of degenerate occ MOs and sums those that correspond to degenerate final virtual MOs Outputs the results to a from i_j dat file Calls the major reading routines and asks
173. or xmgr grace in the file lineplot dat The syntax of the block is as follows LINE 1 what_to_plot This can be 1 orbitals 2 grad_x of orb 3 grad_y of orb 4 grad_z of orb 5 laplacian of orb 6 local e i potential 7 local energy with N 1 electrons fixed in random positions or specified positions see later Then for what_to_plot 1 5 add lines 2 7 for what_to_plot 6 add lines 5 7 for what_to_plot 7 add lines 2a 3a 5 7 LINE 2 number of orbitals norb LINE 3 norb orbital sequence numbers identifying the orbitals to be plotted LINE 4 spins 1 or 1 for each of the orbitals in the previous line LINE 5 the number of points to be plotted along the line LINE 6 xyz coordinates of point A in a u LINE 7 xyz coordinates of point B in a u When plotting the local energy the positions of up to N 1 electrons may be fixed e g to investigate e e coalescences To do this add the following lines in place of lines 2 4 LINE 2a number of electrons to fix 0 if you don t want to fix any LINE 3a position of each electron you want to fix in a u one per line LWDMC Logical Weighted DMC only allow configuration weight to vary between WDMCMIN and WDMCMAX before killing or branching MAKEMOVIE Logical Plot the particle positions every MOVIEPLOT moves see Section 11 MAX_REC_ATTEMPTS Integer This is the maximum number of times DMC will attempt to restart a block if it continues to encounter catastrop
174. ormation about the implementation of DMC within CASINO Our algorithm follows that of 111t is better to use the Slater wave function to compute the first level acceptance probability because p is much lower in general than the corresponding acceptance probability for the Jastrow part p2 hence we are less likely to need to compute the second level acceptance probability than if the situation were reversed 64 Umrigar Nightingale and Runge 5 referred to simply as UNR with some additional developments due to Umrigar and Filippi 7 10 3 1 Imaginary time propagation Let R be a point in the configuration space of an N electron system The importance sampled DMC method propagates the distribution f R t U R R t in imaginary time t where is the trial wave function and is the DMC wave function The importance sampled imaginary time Schrodinger equation may be written in integral form SRi GR OR tt RA aR 8 where the Green s function G R R t t satisfies the initial condition G R R 0 R R 9 The Green s function used in DMC is an approximation to the exact form which is accurate for short time steps 7 t t DTDMO Gpmc R R 7 Gp R R T Gp R lt R 7 10 where E Gp R R r E exp R R a 11 is the drift diffusion Green s function and Gp R R 7 exp Z E R EL R 2Er 12 is the branching Green s function Er is t
175. ound It is expected that the u term will be considerably more important than the p term which cannot describe the electron electron cusps and is therefore best limited to describing longer ranged correlations The q term will describe similar electron nucleus correlations to the xz terms 10 14 2 The u x and f terms in the new Jastrow factor The u term consists of a complete power expansion in rj up to order a which satisfies the Kato cusp conditions at rj 0 goes to zero at the cutoff length ri Lu and has C 1 continuous derivatives at Lu N Tis C lt u rij rij Lu x O Lu riz x o A er ag sc Tij Yan 84 where O is the Heaviside function and T 1 2 if electrons i and j have opposite spins and T 1 4 if they have the same spin In this expression C determines the behaviour at the cutoff length If C 2 the gradient of u is continuous but the second derivative and hence the local energy is discontinuous and if C 3 then both the gradient of u and the local energy are continuous The form of x is N Z C x xa rir rir Lyn x O Lyr rit x Bord I g1 Por rit 5 Pmirir 85 Lx1 Lyr m 2 It may be assumed that Bm Bm z where I and J are equivalent ions The term involving the ionic charge Zz enforces the electron nucleus cusp condition The expression for f is the most general expansion of a function of r rir and rjr that does not interfere with the Kato cusp condi
176. over n successive steps using the nvmcave keyword to prevent this file getting too large The vmc hist file contains all the information necessary to perform a statistical analysis of the data If you want an accurate estimate of the error bar you should then you must analyze the data in vinc hist using the utility reblock which analyzes the serial correlation this utility will also let you change the units of the answer Type reblock in the directory containing the vmc hist file It will first ask you what units you want Then it will print out a reblocking analysis see section 10 22 The error bar will generally be too small for small values of the block length and then should hopefully go up to a roughly constant value for higher block lengths and for very high block lengths it will oscillate as the error bar on the error bar is very high This constant value in the middle is the accurate error bar you want Note that in general VMC calculations will reach this plateau more quickly than DMC calculations since they generally use a much larger time step and there is thus less serial correlation A file reblock plot is normally produced as part of the output of the reblock procedure If you have the programs xzmgr or xmgrace set up on your system then you can visualize the results of reblock by typing plot_reblock you should be able to see a plateau in the variance versus block length curve like in the figure below Note that xmgrace i
177. oy correlations RANLUX trades execution speed for quality through the choice of a luxury level given in CASINO by the RANLUXLEVEL input keyword By choosing a larger luxury setting one gets better random numbers slower By the tests available at the time it was proposed RANLUX at its higher settings appears to give a significant advance in quality over previous generators The luxury setting must be in the range 0 4 Level 0 equivalent to the original RCARRY of Marsaglia and Zaman very long period but fails many tests Level 1 considerable improvement in quality over level 0 now passes the gap test but still fails spectral test Level 2 passes all known tests but theoretically still defective Level 3 DEFAULT any theoretically possible correlations have very small chance of being observed Level 4 highest possible luxury all 24 bits chaotic RANPRINT Integer Setting this keyword to a value greater than zero will cause the first RAN PRINT numbers generated by the CASINO random number generator to be printed to a file random log On parallel machines the numbers generated on all nodes are printed The run scripts should pick out random log files from different stages of a calculation e g VMC config gen DMC equil DMC stats accumulation and rename them appropriately REDIST PERIOD Integer This is the number of DMC moves between configuration redistribu tions on a parallel machine where the number of configurations on each p
178. pe spin polarization density including excited state capability Variable electron hole mass ratio 2D layer separation between holes and electrons possible Pair correlation function s Calculation of electron electron interactions using either Ewald and or our modified periodic Coulomb interaction which is faster and has smaller Coulomb finite size effects Parallelized using MPI tested in parallel on Hitachi SR2201 and SR11000 Cray T3E SGI Origin 2000 SGI Altix IBM SP3 and IBM Pseries Fujitsu Primepower Alpha servers and SunFire Galaxy Linux PC clusters Opteron clusters Also set up for workstation use on DEC Alphas SGI Octane and 02 Linux PC with various compilers MPI libraries not required on single processor machines Keyword input handled using highly flexible ESDF electronic structure data format Self documenting help system 5 Installation Unpack the CASINO_vxxx tar gz distribution in your home directory then go to 5 1 or 5 2 depending on whether your architecture is supported or not 5 1 Supported architecture This currently includes at least Workstations DEC Alpha SGI Octane and O2 Linux PC ifc pgf90 nag compaq lahey 1f95 compilers Parallel Machines Hitachi SR2201 Cray T3E SGI Origin 2000 SGI Altix IBM SP3 IBM P series Fujitsu Primepower SunFire Galaxy Alpha multiprocessors Linux PC clusters Linux PC multiprocessors Opteron clusters See the file CASINO ARCH for the definitive
179. pts which will have been placed in your path Change the units and reblock the answers using the utility reblock The behaviour of the QMC calculation is determined by a list of keywords in the file input To access the internal CASINO help system which tells you about this type casinohelp all to get a list of all keywords that the program knows about or casinohelp keyword for detailed help on a particular keyword Create a blank file containing the keyword input_example and type runumc this will create a sample input file containing all valid input keywords and their default values in the correct format 6 2 Getting started details This section gives basic practical details for running QMC calculations with CASINO and is intended for new users It assumes you already know something about the theory of VMC DMC and variance minimization If you don t then please see the standard references for general details about the methods and the Theoretical Background section in this manual for details about how the algorithms are implemented in CASINO Following the instructions in this section will not necessarily lead to publication quality results but should at least allow you to play with the code and get some feel for how it works 6 2 1 Trial wave functions and pseudopotentials Unless you re interested in electron and electron hole phases in the absence of an external potential in which case you can start straight away the first hurdle to d
180. r r t as it drifts with all other electrons fixed at positions r j 4 Consider the case where the configuration is close to the nodal surface and V Wr R is nonzero at the surface Then r t must be close to a point s such that wr i Ti 1 8i i 1 n 0 So we can write wr ri Ti 1 Ti t tigi rN x Vidr ti Vi 1 i Vigi N tilt si Ary t 32 where A Viwr ri i 1 8i Tit1 Pw is a constant over the move and r t is the component of r t s in the direction of Viwr ri Vi 1 i Titi TN The equation of motion for the single electron drift process is dr t dl Vilt 1 1 15 0 Fi41 T 0 vi t L Vibr ES gt 33 br ri t Hence dr t 1 i t E S e 4 vit g ri t 34 Integrating this over one timestep from time t to t 7 we find that ter ri E y PO 27 ri t altr 35 where the magnitude of the limited drift velocity is given by v rE 27 ri t 14 1427 40 _ 14 4142007 36 T r t vi t T a t So finally the limited drift velocity is written as If Viwr R 0 at the nodal surface then the drift velocity of electron i will be constant close to the nodal surface rather than divergent 69 z t 1 1 Ole 37 av t r where the parameter a has been introduced such that a 1 corresponds to the solution close to the node and the limit a 0 corresponds to the normal
181. radius for the reciprocal space sum in the Ewald interaction used for calculating electrostatic interac tions between particles in periodic systems Its default value is zero Increasing this will cause more vectors to be included in the sum the effect of which is to increase the range of the Ewald gamma parameter over which the energy is constant the default gamma should lie somewhere in the middle of this range This need only be done in exceptional circumstances and the default should be fine for the general user FIXED_E_POS Block This is the position r at which to fix an electron during calculation of the pair correlation function g r r activated with the PAIR CORR keyword To be given in atomic units GAUTOL Real Tolerance for Gaussian orbital evaluation Neglect contribution if exp ar lt 1Q GAUTOL GROWTH_ESTIMATOR Logical Turn on calculation of the growth estimator of the total energy in DMC calculations The difference between mixed estimator and the growth estimator for the energy provides a rough estimate of the time step bias IACCUMULATE Logical TRUE for accumulation of averages Normally turned off during equilibration 24 IBRAN Logical If set to TRUE then enables weighting and branching in DMC IBRAN FALSE may be used to check the DMC algorithm IDEN Jnteger Controls writing of Fourier space charge density to file 1 Density appended to file density data or written to mol_density data in fi
182. rage acceptance probability converges very rapidly to its final value hence by making a couple of extremely short trial VMC runs it is easy to find an appropriate value for DTVMC If you use at least 500 steps in the VMC equilibration phase then CASINO will automatically optimize the VMC time step if you set the keyword OPT_DTVMC to 1 10 2 5 Automatic DTVMC optimization Optimization of DTVMC to achieve an acceptance ratio of roughly 50 is performed in CASINO by means of linear inter extrapolation of two consecutive trial pairs of values log DTV MC Acc Ratio As we are dealing with random electronic configurations these numbers are subjected to fluctuations due to which divergences i e vertical lines arising from the interpolation are bound to appear near the solution In order to avoid these for any two consecutive DIT VMC values that are too close from each other the next trial point is chosen depending only on the last one Also for acceptance ratios too far away from 50 the function becomes irregular so a similar method is used There are some remarks that must be made regarding the optimum value of DTVMC the most important one being that the appropriate way to determine it is by minimizing the correlation time of the energy of the electronic configurations Nevertheless it is found in practice that taking the acceptance ratio to 50 yields almost the same results but with the advantage of being much faster to compute because
183. re the ion sequence number the number of up spin electrons associated with this ion the number of down spin electrons associated with this ion Alternatively one may use the EDIST_ BY IONTYPE keyword block where you replace the ion sequence number with the ion type sequence number and the information is supplied only for each particular type of ion 23 EH_SAFEMODE Logical In some cases the electron hole distribution resulting from a variance minimization run can yield completely wrong VMC results It is found that the problem is introduced by a Jastrow factor which pulls electrons and holes together so strongly that the wave function collapses This problem can be solved by setting EH_SAFEMODE to T which prevents the e h component of the Jastrow factor from cancelling the other terms It is recom mended to set VM_SMOOTH_LIMITS T OPT_MAXITER 3 and running many around 50 variance minimization iterations when using this feature the runvarmin script will allow this if EH SAFEMODEST After this it is a good idea to set EH SAFEMODE F and run a normal variance minimization from the output of the previous one ENERGY_CUTOFF Physical This is the energy cutoff for G vectors used in a the FFT of the MPC interaction required when generating eepot data IRUN 5 and b the FFT of the one particle density required when generating density data IRUN 6 The program will suggest a value for ENERGY_CUTOFF if the existing value is unsuitable or if
184. refactor 1d0 Determinant 1 prefactor DET 31PR2151 Promotion as examples above Note that fancier orbital occupation schemes allowing finer control of which orbitals are used will be introduced shortly along the lines of the one now used in the case of plane wave basis sets see next section 56 8 2 Plane wave basis The file pwfn data specifies the plane wave orbitals used in a CASINO plane wave calculation These orbitals are usually generated using a density functional computation and the data file contains information on the geometry the k points the reciprocal lattice vectors and the expansion coefficients of Kohn Sham orbitals for several bands as well as the corresponding Kohn Sham energies The pwfn data file typically also contains orbitals that are not or only partially occupied in the density functional computation CASINO being a Monte Carlo program dealing with real electrons is based on Slater determinants of orbitals containing zero or one electron per spin The file input specifies which orbitals of pwfn data make up the Slater determinant s used by CASINO 8 2 1 The ground state In many cases one wishes to use a single Slater determinant occupying the lowest lying Kohn Sham orbitals In this case the file input might contain neu 32 ned 32 specifying the number of electrons and block wavefunction GS endblock wavefunction GS is not used with any other input token and instructs CASINO
185. rg See reference 29 Contact Mike Towler mdt26 at cam ac uk to discuss the current status of the interface CASTEP The modern CASTEP 28 is an entirely new Fortran90 version of the venerable Cambridge plane wave program which it is designed to replace It is distributed by Accelrys see www accelrys com mstudio castep html Mail Chris Pickard cjp20 phy cam ac uk to discuss CASTEP and how to get hold of a copy k207 A basic plane wave program apparently only used by Richard Needs GP If you don t work at Lawrence Livermore you re not allowed to use this as it is a Government Secret If you do then obviously you already know about it PWSCE www pwscf org See reference 30 Atoms using a radial grid TCM atomic code mail Richard Needs rn11 at cam ac uk Pseudopotentials CASINO is capable of running all electron calculations where core and valence electrons are explicitly included in the simulation or pseudopotential calculations where the core electrons are replaced by an effective potential The latter is normally advantageous since the computer time required for all electron calculations scales rather badly with atomic number Z this scaling is improved by using pseudopotentials However recent improvements to the algorithms used by CASINO mean that all electron calculations may be performed for systems containing much heavier atoms than is normally believed to be sensible test calculations up to Z 54 have been don
186. rgy with a low enough error bar This is the statistics accumulation part During its operation the rundme script creates a directory called io containing all the input out and config in files from the various stages of the calculation renamed according to an obvious scheme e g out congen out equil and out stats The important files left in the original working directory are dmc hist and dmc hist2 which contain energy and other data as a function of move number see sections 7 10 for precise definitions of what they contain note that dmc hist contains the most important data and dmc hist2 the less important The progress of a DMC simulation is easily visualized by means of the utility graphit which reads the dmc hist file and calls the plotting programs zmgr or xmgrace to show you the re sulting pretty picture This requires you to have xmgr or xmgrace set up on your machine see plasma gate weizmann ac il Grace if you haven t Typing graphit will produce something like this 15 1500 T T T 1400 1300 1200 E POPULATION 1100 4 1 000 5 500 1000 1500 55 4 55 5 Local energy Ha J Reference energy 55 6 H Best estimate 55 7 55 8 0 500 1000 1500 Number of moves This picture shows the results of the simulation of an antiferromagnetic NiO crystal with 10 configs per node on 128 nodes of a parallel computer The upper panel shows how the population f
187. rminant properties of the trial wave function if any Specification of the orbitals Gaussian coefficients Plane wave coefficients Blip coefficients Eigenvalue spectrum used to work out which orbitals to occupy and crudely whether the thing is a metal or an insulator The output file from the program that generated the trial wave function not read by CASINO for reference only These files are generated automatically by various utilities available for different electronic structure programs see Sections 12 13 15 16 and 16 The format of the different files should be clear from looking at the various examples for different dimensionalities and basis sets 40 7 6 1 gwfn data file specification Dimensions of allocated arrays are shown For the larger fields 1 use the following formatting strings title code method functional may be up to 80 characters long free format FORMAT 3 1PE20 13 for reals naturally grouped in triples e g coords amp FORMAT 4 1PE20 13 other reals e g gaussian exponents FORMAT 8110 lengthy lists of integers use e g write 10 format array i i 1 size First line of file is TITLE field MDT 1997 TITLE the title BASIC INFO CODE name of code producing this file i e CRYSTAL95 CRYSTAL98 GAUSSIAN94 GAUSSIAN98 GAUSSIANO3 TURBOMOLE METHOD Comment RHF ROHF UHF DFT S DFT CI etc FUNCTIONAL DFT functional name If method not DFT then
188. rocessor is equalised by transferring configurations between processors The default is 1 but changing it to something larger could produce a performance benefit on some machines RELATIVISTIC Logical If RELATIVISTIC is T then calculate relativistic corrections to the energy using perturbation theory Note that for the moment this can only be done for closed shell atoms See Section 10 27 for further details SPARSE Logical CASINO is capable of using sparse matrix algebra in some algorithms for ef ficiency purposes For systems which are definitely not sparse orbitals not well localized then attempting to use sparse algorithms might actually slow things down Thus until we work out a better way you can toggle this behaviour with the SPARSE flag In these algo rithms a matrix element is considered to be zero if it less than the value of the input keyword SPARSE_THRESHOLD SPARSE_THRESHOLD Real CASINO sometimes uses sparse matrix algebra for efficiency pur poses This keyword defines a threshold such that a matrix element is taken to be zero if it is less than this threshold Changing this quantity might be used to trade speed for accuracy in some cases SPIN_DENSITY Logical Setting SPIN DENSITY to T will activate the accumulation of separate up and down spin densities in the density data file TESTRUN Logical If TRUE then read input files print information and stop TIMING INFO Logical Setting TIMING_INFO to F the default
189. rt of CPP term Electron electron part of CPP Average of the local energy calculated with the MPC All energies are quoted in Hartrees per primitive cell per electron for electron phases per particle for electron hole phases The REBLOCK utility should be used to analyse vmc hist files and change energy units if desired Information about the different components of the energy and how they are calculated may be found in Sections 10 16 10 20 4 54 7 10 dmc hist and dmc hist2 files The file dmc hist contains the most important information from a DMC run while the dmc hist2 file contains less important information There are 13 columns of data in the dmc hist file 1 o 0 Nn DO a A Ww N e pp yN e oO Time step counting from the beginning of the run Number of configurations at this time step Local energy averaged over the configurations Reference energy Best estimate of ground state energy Acceptance ratio at this time step Number of the step within the current block Block number Ratio of the number of configurations to the target number Number of ions Periodicity the number of dimensions in which the system is periodic 0 3 Growth estimator of electron energy 13 Total weight of all configurations NEW 1 2003 The dmc hist2 file contains seventeen columns of data o 0 ND o FF W N KF e e e e Rep Q 0 gt wo N e O Configuration average of Ti Configur
190. rts which may be many hours later If the job would take longer than the maximum time slot on such a machine the scripts allow you to chain jobs so that for example typing runumc n 5 16 will run 5 successive VMC jobs on a 16 node queue with say an 8 hour time limit on each job with the script handling all the required i o decisions The result is 40 hours worth of statistics DMC calculations and variance minimization calculations are more complicated since they require you to submit different types of job within a single calculation DMC config genera tion equilibration statistics accumulation VARMIN n iterations of config generation variance minimization These are handled automatically with all changes to the input files being performed by the scripts For example the first step in a DMC calculation is to perform a VMC calculation to generate a bunch of configurations distributed according to the square modulus of the trial wave function The next step is to use DMC to equilibrate the configurations Amongst other things the rundme script will automatically change the IRUN flag in the input file from 2 to 3 telling CASINO to switch to a DMC calculation See the usage notes at the top of each script for precise information about how they work or on parallel batch machines at least just type the name of the script to get a usage note With all scripts there are normally debug prof and dev flags to run the versio
191. s plane waves blips or on a grid heg data Control file for electron and electron hole systems Mandatory input plus either an xwfn data file or a heg data file Others optional In Output Files eepot data Fourier transform of 1 r required for the MPC interaction Generated by CASINO with IRUN 5 and used in subsequent runs density data Charge density in Fourier space required for the MPC interaction Generated by CASINO with IRUN 6 and used in subsequent runs vmc posin vmc posout List of final electron coordinates on all nodes config in config out List of positions energies etc of a set of configurations required by DMC and varmin Output Files out CASINO output file sent to standard out named out by the run scripts vmc hist VMC history file See Section 7 9 dmc hist DMC history file most important info See Section 7 10 dmc hist2 DMC history file less important info See Section 7 10 Gaussians input CRYSTAL95 98 03 zwfn data GAUSSIAN94 98 TURBOMOLE Splines swfn data A Plane waves ABINIT m pwfn data out CASTEP K270 y JEEP bwfn data PWSCF Blips Numerical orbitals pseudopotential file TCM atomic code 21 Input Files The contents and format of the various input and output files will now be described 7 1 Basic input file input The
192. s a very useful thing to have and it is used by various CASINO utilities If you don t have it set up on your system 11 then you can download it for free from the following website plasma gate weizmann ac il Grace 8e 05 6e 05 4e 05 Standard error of the mean 2e 05 0 L L L L 0 5 10 15 20 Reblock transformation number The file vmc posout contains the final positions of the electrons and the current state of the random number generator so that the VMC run may be continued if desired For example if the error bar is still too large after a certain number of moves then you can make it smaller by running the simulation for longer To do this set the input keyword inew to 0 and rename the vmc posout file to vmc posin in fact the runmvc script will normally ask if it can do this for you then run the calculation for another nmove moves All the extra data will be put onto the end of vmc hist and the error bar from the reblocking analysis will be smaller and can in principle be made as small as desired Finally go and find a system with more than one electron and a Jastrow factor such as the beryllium dimer in molecules be2 in the examples directory Experiment with switching irun between 1 VMC HF and 2 VMC with Jastrow and see the energy and error bar lowering effect of the Jastrow factor without Jastrow Total energy 29 102495953655 0 012473252263 Time 9 3100 se
193. s included in the sum the effect of which is to increase the range of y over which the energy is constant The default value of y should normally lie in the middle of the constant energy region for the default number of reciprocal vectors so playing with EWALD_CONTROL is normally unnecessary 87 4r exp G 47 E 145 G 0 is the potential acting on the charge at r due to its own images and cancelling background The full Ewald potential energy appearing in the QMC Hamiltonian is therefore N U ri ry qiqjve ri rj oo 146 Me s z gt Il amp SS Qe YK Il dj vp ri rj gt 147 i Me Mz gt Il an gt Sa Gu YK II where we have used the charge neutrality condition q gt gt i i G The interaction up r rj approaches 1 r r as r gt r and is independent of the choice of the zero of potential The gradient of the 3D Ewald potential Eq 141 is required for evaluation of the core polarization contribution to the total energy Section 10 19 It is given by amp r rj R ferfe y3lr r R 2 l 2 Vug r rj 2 r FR Po T exp y r r R I Ar exp G 4y 5 Y G Es ERTE 148 G0 10 20 2 2D Ewald interaction Infinite Coulomb sums in systems which are periodic in two dimensions the xy plane according to CASINO are performed using the standard 2D Ewald method originally developed by Parry 38 One way
194. s over all previous time steps In practice it is sufficient to include Tp TPDMC terms in the product provided that T is greater than the number of iterations over which the DMC data are correlated by fluctuations in the reference energy T NINT 10 7 is generally sufficient Let m Tp fs exp TErr m Er m m Then the mixed estimator of the expectation value of a local operator A may be written as VAS sf W R A R R dR Ud f U R R dR m Neontig m A ya _ m T Diem wy mi Ala m m Noontfig m Na ma DO mi where wa m is the weight of configuration a at the end of time step m For unweighted DMC wa is simply the branching factor Note that if we choose A a m Y Ra m AY Ra m then Eq 54 gives us our mixed estimator of the ground state energy at time step m This is used as our best estimate of the ground state energy Epest see Section 10 5 Note that the terms in the II weights are exponential functions of the reference energy hence the II weights are potentially very large or small However it can be seen that any constant contribution to the reference energy will cancel in Eq 54 Therefore in practice we evaluate the Il weights as 54 Tp 1 I m T exp err Ev reer m Er m m 55 0 m where Ey is the variational energy This is necessary in order to avoid floating point errors 10 9 Growth estimator of the energy The total weight of a DMC sim
195. s to prepare an accurate trial wave function given as input to a diffusion Monte Carlo DMC calculation The DMC energy does not in principle depend on the Jastrow factor since it is positive definite and the DMC energy depends only on the nodal surface the set of points in configuration space where the many electron wave function is zero However it makes the calculation vastly more efficient and in general Jastrow factors should always be used The Jastrow factor is stored in a file called jastrow data Again you should look at the examples to see what these look like The various parameters in the files are defined in section 7 2 The adjustable parameters in the file must be optimized for a specific system and this is the purpose of the variance minimization procedure 6 2 5 Other files If for a two or three dimensionally periodic system you want to use the modified periodic Coulomb MPC interaction to calculate the electron electron energies instead of or as well as the standard Ewald interaction then you need to generate two more files for the given geometry before you start doing VMC DMC calculations These are eepot data containing the Fourier transform of the 1 r interaction and density data containing the Fourier transform of the charge density These should be prepared for a given trial wave function and geometry by setting the input keyword irun to 5 and 6 respectively and typing runumc One might consider using the MPC i
196. s used to convert the output of CRYSTAL98 into a gwfn data file readable by CASINO I will assume you know how to run the CRYSTAL program I am aware this is a big assumption You need to run the program using the run script which lives in CASINO wfn_converters crystal9x run_script using the qmc flag as an argument this should be set up automatically for you during utilities compilation The publicly available run script on MDT s CRYSTAL web site does not contain this flag so you need to use the script included with CASINO instead You will need to change some environment variable definitions in the run and crysgen scripts to get them to work properly on your system If you want to use your own copy of CRYSTAL95 or CRYSTAL98 to generate QMC wave func tions you will need to make some minor modifications to the source code see the accompanying CASINO wfn_converters crystal_to_qmc README_CHANGES file NOTE 4 2003 the above changes involve turning off the use of symmetry in k space for some metallic calculations this seems to cause a problem Gilat net unless you also manually turn off all symmetry in the input file with the keyword SYMMREMO This might be fixed one day CRYSTAL will produce three binary files in the scratch directory you define in the run script namely fort 12 basis set geometry common variables fort 10 orbital coefficients and fort 30 eigenvalues These files would normally be deleted at the end of a
197. se to 1 in DMC and 1 2 in VMC Np is the dimensionality of the system which is usually 3 unless a strict 2D system is being studied and D 1 2m is the diffusion constant where m is the particle mass NB D 1 2 for electrons We expect that correlation effects will disappear when the particles have diffused through distances in excess of the physically relevant length scale A For example in an electron gas the density parameter r is the relevant length scale Let T Nmove X T where Ninove is the number of moves and 7 is the timestep Then the number of moves over which we expect correlation effects to be present is 91 2 2NpDTA The number of equilibration moves should be substantially larger than the above estimate of the correlation period in order to ensure that all of the transient effects due to the initial distribution die away The required equilibration period is often considerably greater than one might expect by simply examining the variation of the total energy with time In practical QMC calculations with sensible choices of timestep we often find the VMC correlation period to be about 5 configuration moves and the DMC correlation period to be about 1000 moves Nmove 171 10 22 Statistical analysis of data The Monte Carlo data must be analysed to obtain the required mean values along with a measure of their statistical errors Two factors complicate the analysis 1 the data is normally serially correlated and
198. sgn 0 explp r C R r 71 In this expression sgn 0 is 1 reflecting the sign of at the nucleus and C is a shift chosen so that C is of one sign within re This shift is necessary since the uncorrected s part of the orbital q may have a node where it changes sign inside the cusp correction radius and we wish to replace by an exponential function which is necessarily of one sign everywhere The polynomial p is given by p ao ar azr asr aar 72 and we determine ag 41 Q2 a3 and a4 by imposing five constraints on We demand that the value and the first and second derivatives of db match those of the s part of the Gaussian orbital at r re We also require that the cusp condition is satisfied at r 0 We use the final degree of freedom to optimize the behavior of the local energy in a manner to be described below However if we impose such a constraint directly the equations satisfied by the a cannot be solved analytically This is inconvenient and we found that a superior algorithm was obtained by imposing a fifth constraint which allows the equations to be solved analytically and then searching over the value of the fifth constraint for a good solution To this end we chose to constrain the value of 0 With these constraints we have 1 ln l re C plr X1 73 2 1 d Raya Pre Xi 74 dl 1 dy A lA Roa a p rc p re X3 75 77 A e
199. should be set to 2 VMC with Jastrow at the start of a DMC calculation not as might be thought to 3 DMC When it is required to actually carry out the DMC part of the calcula tion after the configs have been generated the rundmc script will automatically flick irun from 2 to 3 An appropriate number of configs to use for DMC might be 640 10 configs per node on a 64 node parallel machine or more You can usually get away with using 200 configs on a single processor workstation for small systems this will clearly start getting too expensive for a single workstation quite quickly as you increase the system size Parallel machines are in general a good thing Following config generation the distribution of walkers is then allowed to propagate in imaginary time according to the usual DMC rules During a certain period of equilibration the distribution will change until the walkers are distributed according to the ground state wave function of the system subject to the constraint that the wave function has the same nodes as the VMC trial function that we started with This part of the process is called DMC equilibration The best estimate of the energy will fall from the initial VMC value to around the correct ground state energy during this process After equilibration the best estimate of the energy will be roughly correct and we must now propagate for a long period of imaginary time to accumulate enough energy data to estimate the DMC ene
200. sively large when accurate calculations are carried out Furthermore on parallel machines the local energies computed on each processor node are averaged over before being written out The input parameter NMOVE gives the number of energies written out to the vmc hist file in one block of moves Thus for a single processor and single block of moves the total number of configuration moves made is in fact NMOVExNVMCAVExCORPER 10 2 3 Two level sampling The Metropolis acceptance probability for a move from R to R in the standard algorithm is rr RS DARIO PEI MRR min LE f cn 5 It is straightforward to show that if detailed balance 23 in configuration space is satisfied then the resulting ensemble of configurations will be distributed according to the square of the trial wave function However CASINO employs a two level sampling algorithm which has been shown to be considerably more efficient 21 Let us define the first level acceptance probability p R R min f a a 2 HTA DARIDAR and the second level acceptance probability pQXR 7 In the two level algorithm we accept trial moves from R to R with probability p R R po R R It can be shown that provided detailed balance in configuration space is satisfied this procedure also results in an ensemble of configurations distributed according to the square of the trial wave function The two level approach is computationall
201. specified for plane wave and Gaussian orbitals The contents of the wavefunction block are ignored e There are two ways of storing the spline coefficients of localized orbitals They can be stored in pointer arrays which are specially allocated for each occupied orbital or they can be stored in a single array with an index for the orbital The former method permits the memory requirements to be tailored to the appropriate size for each orbital however pointer arrays are relatively slow to look up The latter method is faster but not so memory efficient because the size of the array is determined by the size of the largest orbital e The bsmooth input parameter is used to specify whether localized orbitals with spherical trun cation surfaces should be abruptly or smoothly truncated If it is set to false then the abrupt truncation occurs at the outer truncation radius e For large systems it is advisable to set the sparse input parameter to true in order to exploit the fact that most of the orbital values are zero when updating cofactor matrices 19Tn fact two arrays are used one for the truncated orbitals and one for the untruncated orbitals Since only one or two orbitals are left untruncated using separate arrays for large and small orbitals allows a big memory saving to be made 115 15 Using CASINO with Blip functions 15 1 Blip functions Blip functions were devised by Mike Gillan and implemented in CASINO b
202. statements of variables Use lower case letters for everything else including the bulk of run time code Description of routines to be enclosed in a little box with a description of what it does and an author Later changes to be documented at the bottom of this box For example SUBROUTINE rubbish Routine not to do anything at all MDT 8 2000 3 2001 MDT added capability to write Hello 5 2001 MDT added additional capability to write Everyday Spanish Q Pp E 09 O un USE dsp USE parallel IMPLICIT NONE INTEGER i j k some_integer REAL dp a CHARACTER 11ength my_name i 0 do j 1 10 write 6 Hello do k 1 1000000 i i 1 write 6 Everyday Spanish enddo enddo etc END SUBROUTINE rubbish e Contents of if do and case blocks to be indented by one space 123 e In general don t put spaces between words e g if something_is_true a btc not if something_is_true a b c e Use f90 versions of eq ne It le etc i e lt lt e Maximum of 80 characters per line If you go over this use amp continuation characters at the end of the line and at the beginning of the next one e Two blank lines between subroutines in a given physical file e Error conditions to be handled using the errstop and errwarn routines in the utilities f 90 mod ule e Use endif enddo not end if end do e No double colons in simple variable declarat
203. strophes due to the occurrence of persistent electrons have been investigated Of these the one that seems to perform best in practice involves returning to an earlier point in the simulation and changing the random number sequence If the TRIP_POPN input variable is set to a non zero value then a config backup file will be created This contains the data in the config out file from the previous block If the population on a single node exceeds TRIP_POPN then the data from config backup will be read in the last block of lines will be erased from the dmc hist and dmc hist2 files and the random number generator will be called a few times so that the configurations go off on new random walks hopefully avoiding the catastrophe that led to TRIP_POPN being exceeded If TRIP_POPN is exceeded in the first or second blocks then config in will be read in instead of config backup If a block has to be reset more than MAX_REC_ATTEMPTS times then the program will abort with an error Great care should be taken when choosing a value for TRIP_POPN It should be sufficiently large that it cannot interfere with normal population fluctuations this would lead to population control biasing Note that the population often grows rapidly at the start of equilibration again it must be ensured that automatic block resetting does not interfere with this natural process On the other hand TRIP_POPN should be sufficiently small that persistence is dealt with q
204. t file looks like this for periodic systems QMC 2 want to generate 2 supercell sizes i e 222 a 2x2x2 one 333 and a 3x3x3 one END For molecules only the keyword QMC is required with no additional input 105 13 Using CASINO with the Gaussian94 98 03 programs see CASINO utils wfn_converters gaussian9x 03 Gaussian is an extremely large and widely used commercially available quantum chemistry package More details are available from www gaussian com gaussiantoqme is a utility to read the wave function from the output of a Gaussian94 98 03 G94 G98 G03 calculation and output it in a form compatible with CASINO gaussiantogmc was originally written by Andrew Porter 2000 Note that there is an associated utility called egaussian which extracts the SCF energy and components from a Gaussian output file The gaussiantoqme code requires the existence of a formatted checkpoint file produced by putting FormCheck MO Basis in the route section of the Gaussian job file for G94 G98 produced automatically by G03 It expects this file to have a Fchk suffix The output file of the Gaus sian job is also required It is assumed that this has a out suffix If the original Gaussian job file is present then it will be appended to the end of the QMC input file as is the Gaussian output file 13 1 How to use gausstantoqme If you have a Gaussian job file called dna say and run it to produce dna out and Test FChk
205. t use any calls to external libraries apart from MPI This includes BLAS LA PACK NAG etc e It uses the standard CASINO facilities for handling errors e The routine implements something that we desperately want or is in the TODO list 20 2 2 Bad things e We do not sacrifice speed for additional functionality Think of another way to write it if it slows the code down It s slow enough as it is e CASINO does not like to do things unless they are completely general The development of general complex electronic structure codes can be set back years by people choosing to implement functionality which works only for their current project Non general algorithms tend to get completely thrown away and re implemented which wastes the time of both parties It is not usually much harder to write a general program than it is to write a specific one we understand this is not always true e When writing to the output file the new routines should not write out lots of weird arrays whose meaning is understood only by the author in an incredibly scruffy format full of spelling mistakes In general be as beautiful and informative as you can when writing to output or Mike will just have to make it so which takes ages If there are all sorts of special cases which require different output then so be it select case and if blocks are very useful in this regard e Prior agreement for modification of the code was not obtained Unsolicited submiss
206. ted programs If the program in question uses Gaussian plane wave or blip basis sets then write your own converter to write the output of the program in the appropriate x wfn data format and send it to us mdt26 at cam ac uk for inclusion in future releases It may be a stand alone program or a routine to be inserted into the electronic structure package in question If your program uses some other basis set then this could be a major project Please ask 117 18 Utilities provided with the CASINO distribution A large variety of little programs which do useful things are provided in the CASINO utils directory Here is a reasonably current list of them No other documentation is provided here there may be some in the appropriate directory somewhere in utils e abinit_to_casino_pp Converts pseudopotentials for the ABINIT program into CASINO for mat e billy Shell script for optimizing basis sets and geometrical parameters with CRYS TAL95 98 03 Reasonably vital for developing decent trial wave functions for CASINO with these programs e casinokill Utility for killing CASINO cleanly on workstations sorting out the scratch disk etc e casino_to_abinit_pp Converts CASINO pseudopotentials into the correct format for the ABINIT program e clearup Script for cleaning up after CASINO by removing stray output and indicator files e dfit Draw a polynomial through a set of energy points and find the minimum lives in t
207. ten as H J Ex R y7 R dR 2 SY RAR where E R 47 R R Yr R is the local energy We can evaluate this expectation value by using the Metropolis algorithm 18 to generate a sequence of configurations R distributed according to y R and averaging the corresponding local energies 10 2 2 The electron by electron algorithm The implementation of VMC in CASINO involves repeatedly sweeping through the set of electrons making a trial move of each in turn and accepting or rejecting the move in accordance with the Metropolis algorithm The Metropolis transition probability density is Gaussian of variance 7 where T is the VMC timestep DTVMC In the standard Metropolis algorithm quantities to be averaged are evaluated at the end of each configuration move that is after all the electrons have attempted to move However quantities may be evaluated after each individual electron move so long as care is taken to ensure that for a randomly chosen step the probability of arriving at a particular configuration R having just attempted to move a particular electron i is pi R vp R N 3 Calculating quantities in this way has the advantage that we can average over those trial electron moves that are rejected Suppose for a particular electron move the probability of acceptance is Pacc and the energies of the electron in the existing and proposed configurations are iq and Enew respectively Then irrespective of whether
208. ter 10 20 4 Model Periodic Coulomb Interaction The Model Periodic Coulomb MPC interaction 10 11 12 is used to reduce finite size effects in periodic calculations The exact MPC interaction operator is Hees Y fein Y o0 lets r Ha de 155 where f r is the 1 r Coulomb interaction treated within the nearest image convention which corresponds to reducing the vector r into the Wigner Seitz WS cell of the simulation cell vg is the Ewald potential and p is the electronic charge density from the many electron wave function Y The 89 electron electron contribution to the total energy is evaluated as the expectation value of H exact with Y minus a double counting term EZ TA b I pr p x ver r f r r dr dr 156 Evaluating the expectation value gives 1 go F p r p r up r r dr dr 2 Jws PX Arpa f AI 157 i gt j WS where the first term on the right hand side is the Hartree energy and the term in brackets is the exchange correlation energy We can see that the Hartree energy is calculated with the Ewald interac tion while the exchange correlation energy expressed as the difference between a full Coulomb term and a Hartree term is calculated with the cutoff interaction f In a DMC calculation we require the local energy at every step but we only know the DMC charge density p at the end of the run Normally we have a good approximation to the charg
209. the local energy is calculated in each block of VMC If several block are used then the block variances will be averaged and the standard error in the variance will be calculated The variance calculation assumes that the individual energies are uncorrelated so CORPER should be given a large value Note that VMC_METHOD must be 1 if the variance is to be calculated CEREFDMC Real Constant used in updating the reference energy in DMC algorithm See Sec tion 10 5 22 CHECKWEN Logical Enable a numerical check of the analytic orbital derivatives coded in the various routines such as gaussians_periodic gaussians_mol bwfdet pwfdet etc CONFIGS_VERBOSE Logical Setting CONFIGS_VERBOSE to true will cause the printout of condata_xxx files on each node where xxx is the id number of the node These files contains extra information about the energies of generated configs which is useful in debugging and so on but is not normally required CON_LOC Character Directory to which configurations written CORPER Integer VMC only Local energy calculated only every CORPER moves CUSP_CORRECTION Logical When expanded in a basis set of Gaussian functions the electron nuclear cusp present in all electron calculations is not represented correctly since the gradient of a Gaussian is necessarily zero there Clearly this only matters for s type GTFs since all functions of higher angular momentum have a node at the nucleus When the CUSP_CORRECTION flag
210. the simulation cell CASINO works out the scaling of the translation vectors from the number of electrons and the value of r and builds the determinant part of the ground state wave function accordingly Running through the parameters in the remainder of the heg data file e Exciton units may be selected for calculations of electron hole phases e The mass values only affect electron hole phases and the logical flag indicates whether to use particle symmetries when optimizing the Jastrow factor Only set this true when the masses are equal In 2D you may separate up and down spin electrons for electron only phases or electrons and holes for electron hole phases This is equivalent to having two 2D planes dz apart e For excitonic insulator phases of the electron hole system you may set the orbital parameter Rex in exp r Rox Alternatively you may choose to use a linear combination of Gaussian orbitals If Ng gt 0 the Gaussian form will be selected and should be set to one unless you wish to fit the coefficients to the exponential form Note by default you should include all contributions of orbitals from other cells into the zero cell when dealing with pairing wave functions e For crystalline phases you may select from a set of predefined crystal structures or define your own by selecting manual If there is a separation along the z axis in 2D as mentioned earlier the crystal structure in the two planes must be the same tho
211. the user inputs a value of zero E_OFFSET Physical This keyword gives a constant shift in the total energy such that the the final result will be E Ecale Eos fset The default is zero ESUPERCELL Logical By default total energies and their components are printed as energies per primitive cell Switching this flag to T forces printing of energies per simulation cell in the output file ETYPE Integer Determines type of electron hole phase assuming btype 0 Possible values 1 Homogeneous electron gas fluid phase plane wave basis Homogeneous electron gas crystal phase Pad function basis Homogeneous electron hole gas fluid phase plane wave basis Homogeneous electron gas excitonic insulator phase exponential pairing basis Homogeneous electron hole gas crystal phase Pad function 6 Quasi 2d jellium slab fluid phase with hard walls oo BR WwW bd EWALD_CHECK Logical CASINO and the wave function generating program should be able to calculate the same value for the nuclear repulsion energy given the same crystal structure By default CASINO therefore computes the Ewald interaction and compares it with the value given in the wave function file If they differ by more than 1075 then CASINO will stop and complain If you have a justifiable reason for doing so you may turn off this check by setting EWALD_CHECK to F EWALD_CONTROL Real This is the percentage increase from the default of the cutoff
212. then you must gt mv Test FChk dna Fchk run gaussiantogmc and follow the prompts gt mv dna qmc gwfin data run CASINO The code should automatically detect what sort of Gaussian job it is and give you the opportunity to construct an excited state wave function if applicable It can deal with the following sorts of calculation e HF and DFT ground states open and closed shell e CIS excited states open and closed shell e CASSCF states Getting Gaussian to output these can be problematic for large calculations It is possible that gaussiantogmc will get confused if you use some combination of IOps other than those described in Section 13 2 5 e Time dependent HF TD HF or DFT TD DFT excited states If the user chooses to output a CIS or TD DFT wave function they are given the option of resumming it The wave function must also be resummed if the user wishes to analyse its composition As distributed gaussiantogmc will not do this analysis but if you wish to switch this option on then the flag analyse_cis in cis data f90 must be set to true and gaussiantogmc recompiled With this flag set the code evaluates the percentage contribution of each single excitation to the CIS expansion Degenerate virtual and occupied orbitals are identified and their contributions summed The final output takes the form of files called from1_j dat where i j indicates a range of degenerate occupied orbitals if i j then i is a non d
213. this point It would be greatly appreciated if you could forward a copy of any article published using the results of CASINO calculations to us both for our interest and so we can add references to the CASINO web pages 3 3 Getting the latest version of the code Registered users can download recent versions of CASINO from our secure web site www tcm phy cam ac uk mdt26 casino_users html It should ask you for a username and a password These are normally assigned to particular institutions if you don t have one then please contact Mike Towler who will be pleased to help The site will normally contain the most recent stable release and a beta version The latter is a nightly build released most evenings with the latest changes Clearly it isn t guaranteed not to fall over and you use it at your own risk To help you decide which version you need a link to the CASINO DIARY file is provided which contains a comprehensive list of all modifications to the code Users in the Cambridge TCM group who belong to the Unix group casino may additionally copy the latest source from the directory casino current_beta 4 Functionality of CASINO The CASINO program is actively being developed and many improvements and revision are envisaged but in its present state its capabilities are as follows Variational Monte Carlo including variance minimization of wave functions Diffusion Monte Carlo branching DMC and pure DMC
214. tion drift starting from the configuration s final position This gives a limited local energy that remains bounded as the node is approached and is physically sensible for computing the branching factor albeit approximate Consider the average value of the local energy over a configuration drift to R t 7 from R t R close to a node ER f TRE de z 1 t r Broa f ELE Brest de 38 t where Epest is the current best DMC estimate of the ground state energy The divergence in the local energy close to a node can be written as A B E R Yr R AYr R Ebest Ry Ebest BV R 39 where B is a constant and Ry is the distance of point R from the nearby nodal surface Note that by the same arguments as in the single electron case the full drift velocity is approximately given by V R 1 R So we find that B t r EL R Ebest V R t dt Brow 2 Ru t 7 Ra 8 Ebest BV R _ V R gt Ebest uN Ez R Ebest V R 40 E reduces to Er when the drift velocity is small away from nodes and nuclei In addition Ez is always closer to yest than Ez Hence Equation 40 is a formula for limiting the local energy 70 Since both the drift velocity and the local energy diverge as 1 R1 the limited local energy Ez remains finite as the node is approached In practice even though the local energy is calculated at the end of the configuration move the limiting scheme of Section 10 6 3 is applied to the
215. tional nnodes flag says how many nodes to use on parallel machines Parallel machines with batch queues Depends entirely on the computer architecture In all cases just typing the name of the script will present you with a usage note detailing exactly what options queues nodes and time limits are available runvmc n 2 16 runvmc n 2 16 night Chain 2 successive VMC jobs on 16 nodes Some machines require an extra parameter such as test short medium long or day night weekend or some such to indicate which 16 node queue you want to submit to test might denote a 10 minute queue weekend a 24 hour one etc rundmc c e s 4 16 rundmc c 4 test e s 4 16 long Run 1 config generation 1 DMC equilibration followed by 4 chained DMC statistics accu mulation jobs on 16 nodes Again some machines might require test short medium long Note that some scripts will allow you to run DMC config generation jobs which are gen erally very quick on a small number of nodes say 4 with a short time limit then do the actual DMC on lots of nodes with a long time limit which is what the second example above does Note that it is not currently possible to use the scripts to chain multiple DMC equilibration jobs only 1 is allowed This behaviour could be simulated through multiple rundmce e jobs but the scripts will be changed to do this properly eventually it s only important for incredibly large systems anyway 19
216. tions and goes smoothly to zero when either rj or rj reach cutoff lengths NSF NGT NGS frlrir Tir Tij riur Lp rjr Lpr Olp riz O L pr rjr 5 5 5 Yn yrr 86 l 0 m 0 n 0 Various restrictions are placed on Ymnr To ensure the Jastrow factor is symmetric under electron exchanges it is demanded that Yimnr Ym in VI m l n If ions I and J are equivalent then it is demanded that Yimnr Yimns The condition that the f term has no electron electron cusps is of a Tij 0 de ee TEI TRT which implies that NN 5 5 mirri rir Lyr 0 88 1 0 m 0 for all riz Hence Vk 0 2Ner 5 Yim r 9 89 lym l4m k 80 The condition that the f term has no electron nucleus cusps is of 34 rir 0 N 20 Ee Ree which gives NSD NSS 5 X Cromnr Lim EL rr Lyr 0 91 m 0 n 0 for all r r It is therefore required that Vk 0 NEP Ni y CYomnt LyrYimnt 0 92 m n m n k 10 14 3 The p and q terms in the new Jastrow factor The p term takes the cuspless form p r ys aA S cos G4 Yij 93 A Gi where the Gy are the reciprocal lattice vectors of the simulation cell belonging to the Ath star of vectors that are equivalent under the full symmetry group of the Bravais lattice and means that if Ga is included in the sum G a is excluded For systems with inversion symmetry the q term takes the cuspless form q ri 5 bg X cos Gp
217. to set up a single ground state determinant according to its in build defaults see below The neu and ned parameters can be varied independently in units of one from 0 to the maximum number that is consistent with the supplied Kohn Sham orbitals For example neu 32 ned 31 and block wavefunction GS endblock wavefunction yields a 63 electron system This system is similar to the one above except for an electron hole assuming that the underlying neutral system is actually a 64 electron system Often the GS keyword will be all you need However if the Kohn Sham energy that this hole corresponds to appears at several k points or bands the ground state occupation is not unique We have several degenerate possibilities and CASINO chooses one using its default scheme Note that such degeneracies can also appear for a neutral ground state 8 2 2 Choosing between degenerate ground states Since degeneracies are not uncommon there are several options for altering CASINO s defaults to control the population of orbitals The token GS allows no further tokens We have to use block wavefunction SD GSPEC IP 0 0001 endblock wavefunction 57 SD on its own is identical to GS However now we can use GSPEC to specify how the ground state should be populated GSPEC GSPEC is followed by a two letter token IP P1 P2 and UP determining the internal variable POL a phase between zero and 27 and three flags zero or one IP refers to o
218. tron gas calculation with r 1 and 54 electrons per cell Temporarily remove the old Jastrow factor and copy examples generic new_jastrow jastrow data blank u chi to this directory rename it to jastrow data The variable parameters are not given explicitly in this example file so CASINO assumes they are all zero and generates the initial configurations from the Hartree Fock wave function Edit the blank jasfun data file First change the title from Insert title here to Homogeneous electron gas rs 1 0 Then delete all the lines from START CHI TERM to END CHI TERM we have just deleted the atom centred x parameters which is sensible because we don t have any atoms in a homogeneous electron gas For this first example we don t want to optimize the non linear cutoff parameter because it can take a long time to optimize so let s just accept the default To do this go to the point after the line where it says Cutoff Optimizable 0 NO 1 YES set the Optimizable flag to 0 How many parameters do we want to use Not many since we re just playing so where it says Expansion order set the value to 4 13 Now edit the input file Let s do 1000 configurations nmove nwrcon 1000 Note this is in general not enough you should normally have at least 10000 but we re just playing Let s use corper 10 to reduce serial correlation Make these changes now Now type runvarmin n 2 v The n flag means
219. ty to extract energy standard error and time taken from VMC standard output file single block runs only wfn_converters GAUSSIAN9X 03 gt gwfn data CRYSTAL9X 03 gt gwfn data PW gt bwfn data blips PWSCF gt pwfn data TURBOMOLE gt gwfn data ABINIT gt pwfn data CASTEP gt pwfn data 119 Utilities to read data from GAUSSIAN 9X 03 CRYSTAL9X 03 PWSCF CASTEP ABINIT and TURBOMOLE calculations and convert to CASINO x wfn data format e xwannier Given real PW wave functions in a pwfn data file this finds the Wannier centres and if desired the Wannier functions of the specified orbitals 120 19 Appendix 1 Old Jastrow factor file jasfun data In Easter 2004 the new Jastrow factor jastrow data described in Section 19 was introduced It was designed to replace an old unpublished form of Jastrow factor jasfun data which had been used in CASINO since 1999 The old Jastrow is now redundant but CASINO still fully supports it for purposes of backwards compatibility hence it is documented here The jasfun data file contains details of a Jastrow factor of the form described in Section 10 14 4 JASTROW Silicon atom s 2p 2 ground state Non linear rij term A LO LC 0 44 0 0 0 0 Enable optimization L gt F a u 0 0 Homogeneous rij term No of parameters per spin 5 Enable optimization T Parallel and antiparallel components present T Spin parallel
220. ugh an offset in the x and y direction may be defined The localized orbitals are defined as Pad functions k lr R r exp _ _ __ i r el Irek R where the ith orbital is centred at R The positive Pad coefficients k and k2 are defined next in the input file and you may select to optimize them all or constrain them due to particle symmetry The final section of the file defines the primitive crystal lattice For the same structure as the supercell you would select one basis function and centre it on the origin To define a different structure say BCC primitive with a cubic supercell you would define the appropriate number and position of basis functions The positions written in terms of fractions 0 1 of the supercell vectors The spin variable defines a relative spin orientation between sites assuming you have an antiferromagnetic phase Depending on which lattice you select you will be allowed to simulate different numbers of electrons one per site e Note that the one particle density matrix and relativistic options exist only in test programs external to CASINO and have not yet been re implemented in the main source The program will stop if you turn on these flags e Some excited state calculations can be performed control these through the wave function block in the input file 49 QMC DATA FILE FOR HOMOGENEOUS ELECTRON and ELECTRON HOLE PHASES Electron gas 3D ferromagnetic fcc Wign
221. uickly and that there is insufficient time for a population of configurations containing a persistent electron to stabilise Choosing larger block lengths allows the program to return to an earlier point in the simulation increasing the likelihood that the catastrophe will be avoided 10 11 Evaluation of orbitals in the determinant part of the wave function 10 11 1 Gaussian basis set CASINO can handle Gaussian basis sets up to and including angular momentum l 4 s p sp d f and g functions Suppose we have N Gaussians gy 1 v located at positions r in the primitive cell Let R label the other primitive cells There are therefore copies of the basic set of Gaussian functions located at positions r R We want to evaluate the Bloch orbitals at some point r The orbital is labelled by band v and k point pult Y Cox Y gult ru R explik R 60 w R CASINO implicitly assumes the following about the k points when in Gaussian mode e The k points form a grid l m n kimn b by bs k 61 q1 q2 43 where the q are integers the b are the primitive reciprocal lattice vectors and 0 lt 1 lt q 1 0 lt m lt q2 1 0 lt n lt q3 1 In CRYSTAL the offset k is always zero but CASINO does not assume this e CASINO deals only with real orbitals which can be formed by making linear combinations of the states at k and k The list of k points in the file gwfn data must contain only one of each
222. ulation at time t mz is given by Neonfig m W t l FR idRR XO wa m Mia m 56 a 1 Assuming the DMC simulation to be equilibrated so the ground state is the only remaining component of the wave function the time dependence of the total weight is simply given by 22 W t 7 exp Eo Er r W t 57 where Eo is the fixed node ground state energy Hence the growth estimator of the ground state energy for a single iteration is given by AR A igs 1 7 Mios m 2 l T Mio m By evaluating the expectation value of the argument of the logarithm using the method of UNR 5 and using our estimate of the effective timestep we obtain a much better estimate of the ground state 58 m j exp Er m 1 rerr m Mio m 1 1 mat HM Tp Meios m Mio m Egrow e l m growth m TEFF di 2m 1 H Tp Mot m Nm Him 1 Tp 1 Meor m 2 59 TErF m D1 Ur Tp Mot m 73 Equation 59 is used to evaluate the growth estimator of the energy in CASINO if the GROWTH_ESTIMATOR flag is set to TRUE in the input file Note that when evaluating the growth estimator using Equation 59 one of the constant multi plicative terms in the I weights of Equation 55 exp rerr 1 Ey does not cancel out Hence we need to include one such term in the numerator of the argument of the logarithm 10 10 Automatic block resetting Numerous schemes for preventing population control cata
223. ummary of routines used in gaussiantogmc e e e sas ee 14 Use of localized orbitals and bases in CASINO 14 1 Theoretical backeroutid 4 ae each se cid ew OR GE ROR OE eS RS 14 2 Using CASINO to carry out linear scaling QMC calculations 15 Using CASINO with Blip functions 15 1 Blip ORCOS basis eee RARE DA a we A da A 15 2 The blip conversion utility gt e sagada cee ee eRe aoe Oo Oh doe oe ea os 153 The blip conversion utility orgias a bea bE eh A 16 Using CASINO with other supported programs 16 1 TURBOMOLE 2 4225 4244444 oe sal de e e tacto dl de 4444404 162 CASTED 2002 244544684244 44866 ee BRIDA de Rd A RAR A de A 16 3 PWSOE E SSS EE EERE a addo Stee Yay ae Se 164 OAR INI sara ee Pe a Ree eae es 17 Using CASINO with unsupported programs 18 Utilities provided with the CASINO distribution 19 Appendix 1 Old Jastrow factor file jasfun data 20 Appendix 2 Programming Guide for CASINO CU AMES so ae ew bg REA SER RRO org ee OE ON Se Eee oe a a Se ae 202 CONTENT sico E05 ERE Ea Gale bed ad he ee 203 PERFORMANCE 2 442404 4284 seme SHR Pee eee ete A EE BO eS 204 BUG REPORES gt o 4 4 5 a aes eA Rae eae eh ne OS ee eee ee ee A 20 5 REQUESTS FOR NEW FEATURES Bibliography 106 106 107 110 112 112 113 116 116 116 116 117 117 117 117 117 117 118 121 123 123 125 126 127 127 128 1 Introduction CASINO is a code for performing quantum
224. urgh PA 2nd edition 1996 108 D Alf and M J Gillan Phys Rev Lett 70 161101 2004 116 D M Ceperley and B J Alder Phys Rev B 45 566 1980 118 V R Saunders C Freyria Fava R Dovesi L Salasco and C Roetti Mol Phys 77 629 1992 87 See e g D E Parry Surf Sci 49 433 1975 and erratum Surf Sci 54 195 1976 88 J E Dennis D M Gay and R E Welsch Algorithm 573 NL2SOL An Adaptive Nonlinear Least Squares Algorithm ACM Transactions on Mathematical Software 7 369 1981 94 N Marzari and D Vanderbilt Maximally localized generalized Wannier functions for composite energy bands Phys Rev B 56 12847 1997 112 G Berghold C J Mundy A H Romero J Hutter and M Parrinello General and efficient algorithms for obtaining maximally localized Wannier functions Phys Rev B 61 10040 2000 112 113 129 42 43 44 45 46 47 T Kato Commun Pure Appl Math 10 151 1957 76 Y Kwon D M Ceperley and R M Martin Phys Rev B 48 12037 1993 Y Kwon D M Ceperley and R M Martin Phys Rev B 58 6800 1998 M Holzmann D M Ceperley C Pierleoni and K Esler Phys Rev E 68 046707 2003 Relativistic corrections to atomic energies from quantum Monte Carlo calculations S D Kenny G Rajagopal and R J Needs Phys Rev A 51 1898 1995 100 L M Fraser PhD thesis Imperial College of Science Technology and Medicine 1995 99 130
225. value of f against the distance from the point given in line 6 KWARN Logical The KWARN flag is relevant only in calculations using a plane wave basis set If the flag is set to TRUE then the routine PWFDET_SETUP will issue a warning whenever the kinetic energy calculated from the supplied orbitals differs from the DFT kinetic energy given in the pwfn data file by more than an internal tolerance usually set to 10 If the flag is FALSE then CASINO will stop with an error message on detecting this condition Note 25 that in cases where the DFT calculation which generated the orbitals used fractional occupation numbers the kinetic energy mismatch is very likely to occur since QMC deals in principle only with integer occupation numbers hence the existence of this flag LCUTOFFTOL Real This is used to define the cutoff radius for the local part of the pseudopo tential It is the maximum deviation of the local potential from z r at the local cutoff radius LIMDMC Integer Set modifications to Green s function see Section 10 6 May take values 0 no modifications applied 1 Depasquale et al limits to drift vector and energy 2 Umrigar et al modifications need to supply ALIMIT parameter LINEPLOT Block This utility allows you to plot the value of certain quantities such as orbitals and their derivatives or various local energy terms along a line from point A to point B The data will be plotted in a format suitable f
226. variance minimization calculation is to generate them This is done through an ordinary VMC run with the keyword irun set to 2 and nwrcon set to the number of configurations to generate clearly nmove must be greater than or equal to nwrcon you might want it to be greater if you need a longer run to get the error bar acceptable this is useful for judging the success of an optimization at each stage Also corper should be set to a high value gt 10 in order that the configurations are completely uncorrelated The configs generated are written to a file called config in These config urations are then read and fed into the variance minimizer through a CASINO run with irun set to 4 The optimizer then adjusts the value of the parameters in such a way that the variance is minimized At this stage the configurations are distributed according to the original unoptimized wave function Thus it is often a good idea to recalculate the configurations using the new optimized wave function to generate a new config in and then perform the minimization again This iterated procedure can typically be carried on as long as desired but doing more than two or three such iterations is normally of little benefit The whole variance minimization procedure is actually automated by the runvarmin script so the two main things to worry about are how many configurations to use and how many parameters to include in the Jastrow factor Good advice about th
227. w factor used by CASINO is called jastrow data Note this replaces an older form of Jastrow factor called jasfun data which was used until Easter 2004 CASINO still supports this and it is therefore described in Section 10 14 4 and Appendix 19 The format of the jastrow data file is as follows This represents the state of the file before optimization since the variable parameters are not given explicitly and are therefore assumed by CASINO to be zero START JASTROW Title Title of system goes here Truncation order 3 START U TERM Number of sets 1 START SET 1 Spherical harmonic 1 m 00 Expansion order 6 Spin dep params 0 gt uu dd ud 1 gt uu dd ud 2 gt uu dd ud 1 Cutoff a u Optimizable 0 N0 1 YES 0 d0 1 Parameter values Optimizable 0 N0 1 YES END SET 1 END U TERM START CHI TERM Number of sets 2 START SET 1 Spherical harmonic 1 m 00 Number of atoms in set 1 Labels of the atoms in this set 1 Impose electron nucleus cusp 0 N0 1 YES 0 Expansion order 6 Spin dep params 0 gt u d 1 gt u d 0 Cutoff a u Optimizable 0 NO 1 YES 0 d0 0 Parameter values Optimizable 0 N0 1 YES END SET 1 START SET 2 Spherical harmonic 1 m 00 Number of atoms in set 4 Labels of the atoms in this set 2345 Impose electron nucleus cusp 0 N0 1 YES 0 Expansion order 32 6 Spin dep params 0 gt u d 1 gt u d 0 Cutoff a u Optimizable 0 NO 1 YES 0 d0 0 Parameter values Optim
228. with your new machine type 4 Send 1 2 and 3 to Mike Towler at mdt26 at cam ac uk 5 Go to supported architecture list above 5 3 Further installation notes e By default all Fortran90 files will be compiled with full optimization If you want to compile with debugging or profiling compiler options instead of the optimizing ones or to compile a separate development version then type make debug make prof or make dev The Makefile will keep separate copies of all object files and binaries for the different compilation levels e Note that on a multi machine environment like the Cambridge High Performance Computing Centre where many different machines access the same filespace you can use the same CASINO source distribution for all machines The Makefile keeps separate compilation records for each machine and source type opt debug prof dev thus you can compile the same source on more than one different machine at the same time 6 How to use CASINO 6 1 Getting started summary You should first learn how to use CASINO using the examples provided Example input files for finite systems atoms molecules and for systems with periodicity in one two or three dimensions polymers slabs and crystalline solids plus various electron and electron hole systems can be found in the CASINO examples directory If you have set up the code as described in the installation section you can run the code with the runumc rundmc and runvarmin scri
229. within DMC except in the limit of perfect importance sampling In VMC the kinetic energy may be evaluated using any of the three estimators in Eq 124 CASINO automatically uses XK for the evaluation of the total energy because this normally leads to the lowest variance However the lowest variance of the kinetic energy itself is often obtained from T In DMC the three estimators of Eq 124 are not exactly equivalent and XK should always be used For the Slater Jastrow wave function of Eq 111 we have VD Vit pa tV 125 V2D V DY 2 Anly i i T 12 V ln w po pa tY 126 The terms involving Slater determinants may be evaluated by expanding D in terms of the cofactors of the ith column of the Slater matrix D i If electron has spin up for example the required expansion is Dt det D X py ri cof D 127 j Since all the cofactors appearing in this equation are independent of r we obtain V iD zl pr D Vt Dis 128 j v2D t pr D Vive Di 129 J When moving electron i from r 4 to r it is useful to be able to evaluate the kinetic energy at the new position before updating the D matrix Since the cofactors in Eq 127 are independent of r Eqs 128 and 129 become yo 1 new 7 old DT new q 3 Vip ri Dji 130 Jj yo 1 2 new z old DT new ql y Vi W t Dii gt 131 J where q Dtnew Dt cld 84 10 18 Evaluatin
230. xplicitly wish to use a core polarization potential see Section 10 19 the few lines at the bottom relating to this should be omitted CASINO pseudopotential library www tcm phy cam ac uk mdt26 casino_users html 37 7 4 Charge density file density data This contains the Fourier transform of the charge density and it is used in the evaluation of the Modified Periodic Coulomb MPC interaction Before doing a QMC calculation with the MPC interaction this file should be generated from the trial wave function given in the x wfn data file by using runvmc to run CASINO with IRUN set to 6 in the input file note that the eepot data file is also required for this Note that if the IDEN parameter is set to 1 in the input file then the QMC density will be accumulated at the end of the file following the description of the charge density of the trial wave function Plutonium doped lanthanum nickelate Charge density in reciprocal space Real space translation vectors Cartesians in atomic units 0 000000000000 3 370326931161 3 370326931161 al 3 370326931161 0 000000000000 3 370326931161 a2 3 370326931161 3 370326931161 0 000000000000 a3 Reciprocal space translation vectors Cartesians in atomic units 0 932132911066 0 932132911066 0 932132911066 bi 0 932132911066 0 932132911066 0 932132911066 b2 0 932132911066 0 932132911066 0 932132911066 b3 Number of atoms in the unit cell 2 of atoms in cell Positi
231. y Dario Alf 35 15 2 The blip conversion utility The blip utility is used to convert pwfn data files generated by a plane wave DFT package into bwfn data files where the wave function is represented in a basis set of localized blip functions on a grid This should make the code go faster and scale better with system size than plane waves which are horrible for QMC calculations It can require a lot of memory and disk though See CASINO utils wfn_converters pw_to_blips The quality of the blip expansion can be increased with the input parameter mul the converter will ask you about this when it is run This results in an increasing number of blip coefficients so a larger memory occupancy in CASINO but the CPU time should be the same To test the quality of the blip expansion set ltest t when the converter asks you about it The converter then samples the wave function the Laplacian and the gradient at 1000 points in the simulation cell and computes the quantity BW PW V BW BW PW PW The closer to 1 this quantity is the better the representation Through increasing mul one should always be able to get alpha as close to one as wanted The gain in speed with respect to plane waves should be of the order of nPW 64 15 3 The blipl conversion utility Localized orbitals see Section 14 can also be represented in terms of blips In order to do this we require a pwfn data file a wannier_centers dat file and a b
232. y a com bination of insufficient computer power and a lack of available efficient QMC computer programs general enough to treat a similar range of problems to regular DFT codes Fast parallel computers are now widespread and you are now in possession of just such a computer program that is capable of carrying out QMC calculations for a wide range of interesting chemical and physical problems on a variety of computational hardware CASINO has been designed to make the power of the quantum Monte Carlo method available to everyone and we hope you enjoy using it 3 Miscellaneous issues 3 1 Support CASINO has documentation and examples etc but it is a research code and learning how to use it is a significant task This is particularly the case if the user does not have relevant experience such as familiarity with VMC and DMC calculations and knowledge of DFT HF methods for atoms molecules and solids We are finding that supporting users can take a large amount of our time and so we are having to limit the number of such groups that we can work with directly We are finding that most people need quite a lot of help and the project turns into a collaboration but of course we cannot enter into too many projects of this type as our time is limited We do have people visiting here to learn about the codes and how to do calculations and this seems to work well In general you are welcome to have a copy of CASINO if you sign and adhere to our agreement
233. y advantageous because the Metropolis accept reject step can be carried out in two stages if the first level is accepted with probability p then we compute the Jastrow functions of the new and old configurations and determine whether the second level is po R R min 1 rane 63 accepted with probability p2 Thus if a move is rejected at the first level then we do not have to compute the Jastrow functions for the new and old configuration It is unclear how two level sampling should be applied in the context of VMC method 2 since we do not know the full acceptance probability for the trial electron move if the move is rejected at the first level Therefore when VMC_METHOD is 2 single level sampling is used for those moves on which electron energies are calculated with p p2 being the acceptance probability two level sampling is used elsewhere 10 2 4 The optimal value of the VMC timestep The VMC timestep DTVMC should be chosen such that the correlation period as determined by reblocking analysis of the resulting configuration local energies is minimized For large timesteps the move rejection probability is high which obviously leads to serial correlation On the other hand for low timesteps the configuration does not move very far at a given step and so serial correlation is again large A rule of thumb for choosing an appropriate DTVMC is that the average acceptance probability should be about 1 2 The ave
234. y cases however some of the molecular orbitals have small s components which may have nodes close to the nucleus The possible presence of nodes inside the cusp correction radius complicates the procedure because the effective one electron local energy diverges there One could simply force the cusp correction radius to be less than the radius of the node closest to the nucleus but in practice nodes can be very close to the nucleus and such a constraint severely restricts the flexibility of the algorithm In practice we define small regions around each node where the effective one electron local energies are not taken into account during the minimization and from which the cusp correction radius is excluded Gaussian cusp correction is activated through the input keyword CUSP_CORRECTION this has an effect only if the system contains at least one all electron atom In periodic systems it has proved difficult to implement this scheme efficiently and while it works perfectly well the performance of the code is significantly affected Check the timings More information about the cusp correction of each orbital at each nucleus can be produced with the keyword CUSPINFO which can be useful in fine tuning Clearly this can produce a lot of output so beware 10 14 The Jastrow factor The Slater Jastrow wave function is v R exp J X cn D D 82 where the Dn are determinants of up and down spin orbitals and J is the Jastrow factor 10 14 1 Gener
235. y duplicate the u and x terms additional constraints can be placed on f to remove these terms if desired This is usually a good idea if the cutoff length for f is about half that of u and about the same as that of the corresponding x set If the cutoff length of f is very different then duplication of u and x should be permitted The number of f parameters grows very rapidly with the electron electron and electron nucleus expansion orders These should normally be either 2 or 3 The spin dependence options for f are exactly the same as for u The p and q terms can only be present in periodic systems In practice neither appears to be particularly useful Note that q should only be used if the origin is a centre of inversion symmetry of the charge density The spin dependence options for p and q are the same as those of u and x respectively 5For electron hole systems a spin dep params value of 0 means that the same coefficients are used for all spin particle types 1 means that different parameters are used for electron electron electron hole and hole hole pairs 2 means that different parameters are used for parallel and antiparallel electron and hole pairs but that electron hole pairs all have the same parameters 3 means that different parameters are used for parallel and antiparallel spin pair types and different parameters are used for each type of particle pair 6The u parameters are listed in the following order for e
236. y of the sets of atom centred functions is set to zero then it is given a default value The default is equal to the longest of the nearest neighbour distances for the atoms in the set CASINO will perform a neighbour analysis to work this out For single atoms which don t have any neighbours it is currently given the value of 2 0 Where functions are spin dependent the first column 1 is for the spin parallel or spin up case and the second column 2 is for the antiparallel or spin down case For the No of powers in ri rirj ririj 2 rirjrij 2 terms the ririj and rirjrij terms require two integers for the l4 and m4 of Eq 101 and the lm5 and nm5 of Eq 102 At the moment only the terms in S3 with l m have been coded The extra terms will be included in a later release Modifications for jellium slab system The jellium slab although homogeneous in the plane of the slab x and y directions is inhomogeneous in the z direction The introduction of the usual two electron term in the Jastrow factor modifies the density away from the desired z dependence a one electron term depending only on z is required to restore the correct form This term is only available for ETYPE 6 and is introduced by including the following lines in jasfun data Enable one body z dependent term T Enable optimization T No of Fourier components present 2 Fourier components frequencies and coefficients 0 930842267730309 0 335794237
237. you the definition and type of keyword casinohelp search lt text gt finds lt text gt in descriptions of all keywords casinohelp all lists all possible keywords casinohelp basic lists all basic level keywords casinohelp inter lists all intermediate level keywords casinohelp expert lists all expert level keywords 118 ion_dist Program for automatic generation of the edist_by_ion block in the input file for CASINO Currently works only for antiferromagnetic Wigner crystals that have been generated by the CRYSTAL program since the numbers are generated from an analysis of gwfn data modify_inputs Simultaneously modify all input files in subdirectories off the current directory e g CASINO examples by adding or deleting keywords or by changing the value of a keyword nstring Generate integer number sequence from a to b step c Useful for strings required in jastrow data files plot_corr Plot calculated pair correlation function ee eh hh as calculated by CASINO and dumped in a vmc corr or dmc corr file plot_jastrow Convert functions defined in the old jasfun data Jastrow factor file into a form suitable for plotting by mgr or whatever With the new jastrow data form this functionality has been subsumed into CASINO itself see the jastrow_plot input keyword plot_hist General program to read a vmc hist dmc hist dmc hist2 files and plot selected quantities as a function of move number plot_spline Code to re
238. ype of particle H electron O hole and C other atoms of course this labelling is silly for our purposes but we didn t invent the format Columns 2 3 4 are the x y z coordinates of the particle Column 5 specifies the charge of the particle This information is then repeated for the number of times specified in the CASINO input file For example if we run a vmc calculation for nequil 2 nblock 1 nmoves 2 and movieplot 1 then 5 sets of geometry will be created in the movie out file 11 2 Visualisation Having generated the movie out file we are now able to visualise the results with VMD or Jmol 11 2 1 VMD VMD Visual Molecular Dynamics is a molecular visualization program It supports computers running MacOS X Unix or Windows is distributed free of charge and includes source code VMD can be downloaded from the following website http www ks uiuc edu Research vmd e Type umd a VMD console and a VMD Display window will appear e In the VMD console window type menu main on An extra VMD Main menu bar will appear e On the VMD Main menu bar click on File New Molecule A Molecule File Browser will appear 102 e In the Molecule File Browser browse for the file movie out and choose the file type to be xyz Click Load to open the file e On VMD Main click on Graphics Representations A Graphical Representations menu bar will pop up Choose CPK as
239. zation method Extensive tests are currently being carried out to establish the importance of locating the global minimum of the LSF and to determine the best optimization strategy The following developments to the new optimization method will be made in the near future 1 The default behaviour of the new optimization method will be improved 2 The ability to use the new optimization method for electron hole systems will be added 3 The new optimization method will be given the ability to read in a config in file containing a set of configurations and construct the quartic coefficients for that configuration set 4 The new optimization method may be combined with the standard variance minimization al gorithm in such a fashion that the new method is used to optimize linear Jastrow parameters while the nonlinear parameters are optimized by the standard method Each time the nonlin ear parameters are adjusted the quartic coefficients will be constructed and the LSF will be minimized with respect to the linear parameters using the new method 10 25 Further considerations in electron hole systems CASINO has the ability to include positively charged particles of variable mass holes in the simu lation in addition to electrons Currently these may only be used in electron hole phases without an external potential but the code needs only a few trivial changes for these things to be able to wander around inside real crystals useful for studyi
240. ze the trial wave function and with how many parameters They must also decide on the number of configurations to 94 be used in the optimization These choices are system specific and depend on the level of accuracy to which the user wishes to work A systematic approach to deciding on an appropriate number of variational parameters is to start by optimizing a few parameters then to add more and re optimize and so on until the decrease in energy resulting from the inclusion of additional parameters is comparable with the error bars in the VMC energy Once this stage has been reached adding further parameters can be of no benefit indeed the extra variational freedom will then become problematic leading to poorer trial wave functions as parameters are optimized for specific configuration sets If further parameters are to be optimized the user must increase the number of configurations It is clearly desirable for the VMC generated configurations to be completely uncorrelated This can be achieved by giving corper a large value Reblocking VMC energies in a preliminary VMC run will allow the user to determine the correlation period for VMC energies which in turn sug gests a suitable value for corper It is also essential that the VMC configuration generation run is fully equilibrated Since VMC equilibration is usually computationally inexpensive this should be straightforward enough The utility plot_umc_energy can be used to verify that the

CASINO manual

Contents

Download Pdf Manuals

Related Search

Related Contents