I am not an expert for pseudopotentials, so please correct me if I am wrong ...
Actually, i think the requirement is that the integrated charge of each wavefunction inside the cutoff is equal for the pseudopotential and the all electron calculations. This is the norm-conversation condition in the norm-conserving pseudopotentials and is ensured by the generation procedure. It obviously does not apply for the core wavefunctions, since there are none for the pseudopotential. The PS core charge given by the atom program is probably just the charge used for the partial core correction, it is zero if no core correction is used. By the way, there is a nice introduction to pseudopotential generation by Paolo Giannozzi (http://www.fisica.uniud.it/~giannozz/software.html). Best regards, Wolfgang On Tuesday 03 July 2007 22:52, Sterling Paramore wrote: > Thanks! Actually, that seems to have made the DFT calculations stable. > > However, I must not be understanding what some features of the atom > program mean. When I had the core correction radius at the very small > value of 0.075, I got excellent agreement between the AE core charge > and PS core charge (as given by gnuplot charge.gplot) at all values of > r. When I increase the value to 0.66, I get good agreement beyond > 0.66, but very poor before it (the PS charge curves are usually much > smaller than AE). In the atom manual, it says that the 4th column of > the PSCHARGE file contains the core charge density, multiplied by > 4*pi*r**2 (so integrating this function up to some value r, gives the > total charge enclosed in r). With such large differences in the > charge density, the total charge contained within the cutoff is very > different between the PS and AE results. I read in Troullier and > Martins (1991) that a good pseudopotential encloses the same charge > density at the cutoff as the AE result. So as far as I'm > understanding it, by increasing the core correction radius to 0.66, > I'm getting stable DFT results, but no longer have "good" > pseudopotentials. What am I missing here? > > I have attached a figure of gnuplot charge.gplot when the core > correction radius is set to 0.66. > > Sincerely, > Sterling Paramore > > On 7/3/07, Wolfgang Voegeli <[EMAIL PROTECTED]> wrote: > > Dear Sterling, > > > > Maybe your problem is caused by the very small radius for the core > > correction in the pseudopotential input file. This leads to a very hard > > pseudopotential, which is difficult to describe with the real-space grid. > > Increasing the meshcutoff reduces the error, but the forces for a > > meshcutoff of 2000 Ry are still 10 times higher than for the other atoms. > > i think a rule of thumb for the radius for the partial core correction > > is the radius where the valence and core charge are about equal. I > > recommend you to try to increase the radius to something like that. The > > core radii in your input file are also quite small, the O pseudopotential > > on the siesta webpage uses 1.15 bohr. > > For good results you will still have to increase the meshcutoff to > > 200-300 Ry. > > > > Best regards, > > Wolfgang > > > > On Tuesday 03 July 2007 07:32, Sterling Paramore wrote: > > > Dear SIESTA users, > > > > > > I'm very new to SIESTA, and I'm having a few problems generating > > > pseudopotentials (PP) for Oxygen for use in simulations of > > > orthosilicic acid (SiO4H4). I'm trying to follow the work of Ma, > > > Foster, and Nieminen (2005,JCP 122, p144709) carefully, so I'm using > > > GGA with the PBE functional (I've also been checking it with RPBE, > > > which supposedly gives better agreement with experiment). I generated > > > the PPs for O using the atom program with the following input file: > > > > > > pe Oxygen > > > tm2 2.0 # PS flavor, logder R > > > n=O c=rp # Symbol, XC flavor,{ |r|s} > > > 0.0 > > > 1 4 # norbs_core, norbs_valence > > > 2 0 2.00 0.00 # 2s2 > > > 2 1 4.00 0.00 # 2p2 > > > 3 2 0.00 0.00 # 3d0 > > > 4 3 0.00 0.00 # 4f0 > > > 0.60 0.60 0.60 0.60 0.00 0.075 > > > > > > With the above values, I get good agreement between the all-electron > > > and pseudopotential orbital eigenvalues. I generated PP for Si and H > > > using very similar input files as well. The problem I'm having is > > > that when I go to run a SIESTA calculation (e.g, verlet MD or CG > > > minimization), the forces on the oxygen atoms are 3-4 ORDERS of > > > MAGNITUDE larger than the forces on the Si or H atoms. Obviously, > > > this is not correct. > > > > > > The strangest issue I'm seeing is that if I generate the O PP with CA > > > instead of RP (or PB), and then run SIESTA using GGA and PBE, the > > > forces on the Os are comparable to those on Si and H (more like what I > > > would expect). This suggests that there's some sort of instability > > > that is being generated in the PP when using PB. However, I have > > > tried varying just about everything I can in the PP generation with no > > > luck. I have tried changing the cutoff from 0.40-2.00, using the full > > > core charge by setting rcore_flag=-1, not using the new core > > > correction at all (via pg), and even tried using the spin versions. > > > > > > > > > The SIESTA input file I've been using is given below. If anyone has > > > any ideas what I might be doing wrong, I'd really appreciate any help. > > > > > > Sincerely, > > > Sterling Paramore > > > > > > # > > > ----------------------------------------------------------------------- > > >---- -- # FDF for orthosilicic acid > > > # > > > # > > > ----------------------------------------------------------------------- > > >---- -- > > > > > > SystemName Orthosilicic Acid > > > SystemLabel SiO4H4 > > > > > > # > > > ### Species identification > > > # > > > > > > NumberOfSpecies 3 > > > NumberOfAtoms 9 > > > > > > %block ChemicalSpeciesLabel > > > 1 14 Si > > > 2 8 O > > > 3 1 H > > > %endblock ChemicalSpeciesLabel > > > > > > # > > > ### DFT Definitions > > > > > > XC.functional GGA > > > XC.authors RPBE > > > > > > # > > > ### Basis Definitions > > > # > > > > > > User-Basis .false. > > > #PAO.BasisSize SZ > > > PAO.BasisSize DZP > > > PAO.BasisType split > > > > > > PAO.EnergyShift 272 meV > > > > > > > > > # > > > #Boundary Conditions > > > # > > > > > > LatticeConstant 0. Ang > > > > > > %block LatticeVectors > > > 1.000 0.000 0.000 > > > 0.000 1.000 0.000 > > > 0.000 0.000 1.000 > > > %endblock LatticeVectors > > > > > > > > > AtomicCoordinatesFormat Ang > > > %block AtomicCoordinatesAndAtomicSpecies > > > -0.96890 0.87640 0.04850 1 > > > -1.83470 1.93920 1.08230 2 > > > 0.22500 1.75900 -0.81420 2 > > > -0.21160 -0.34630 0.98630 2 > > > -2.05870 0.15050 -1.06140 2 > > > -1.49050 -0.43130 -1.62660 3 > > > -0.28460 2.43880 -1.32300 3 > > > -1.14020 2.36020 1.64900 3 > > > 0.42990 0.14010 1.56310 3 > > > %endblock AtomicCoordinatesAndAtomicSpecies > > > > > > # > > > # MD Definitions > > > # > > > > > > WriteForces .true. > > > > > > MD.TypeOfRun CG > > > MD.NumCGsteps 10
