Re: open shell systems ( No.1 ) |
- Date: 2005/12/04 02:57
- Name: T.Ozaki
- Hi,
> Question: is this set up for a charged molecule calculation really correct (I seem to > add one extra electron to the system, at least in the Atoms.SpeciesAndCoordinates > section of the .dat file)?
Your setting is correct. The total number of electons should be neutral in the specification of "Atoms.SpeciesAndCoordinates", where each number is used for only the generation of the initial charge density. You can control the actual total number electons by the keyword, "scf.system.charge". Also you can easily check the total number of electrons by summing the Mulliken charge up after the calculation to confirm it.
> If the settings are correct, together with set.system.charge = 1, I seem to get wrong > dipole momentum of the ammonium ion: it is of the same order as the dipole momentum of > the ammonia molecule (nh3). As far as I can grasp, nh4+ is a symmetrical molecule and > can not have any sizable dipole momentum (other than from numerical error). > The situation with so4(2-) ion seems to be even worse: the dipole momentum is very > large, almost 1e3, which is completely unbelievable. Most probably I my setting are > too wrong here.
I have also tried the two systems and obtained completely the zero dipole moment for both the systems as we expect. Possible source of errors you met can be listed as folllows:
(1) Could you obtain the convergent results ? (2) scf.lapack.dste=dstevx is better for a symmetric system. The lapack routine "dstedc" and "dstegr" tends to fail to correctly calculate symmetric systems. (3) Remind that OpenMX uses the regular real space grid for the numerical integration. So, the total energy and the dipole moment depend on the rotation of system if you use a lower cutoff energy (scf.energycutoff), while, of course, the dependency on the rotation of system is reduced, as the cutoff energy increases. To avoid this dependency, I constructed the structures by modifying "Methane.dat" file in the "work" directory, since the relative structural orientation to the regular real space grid is symmetrical. Also, you will see that OpenMX gives a reasonable dipole moment for many molecules, and how the dipole moment depends on the choice of basis functions at http://netserver.aip.org/cgi-bin/epaps?ID=E-JCPSA6-121-301439
> Question 2: if the settings are right, does it mean that the "default" dft method in > nh3.dat file example can handle only closed electronic shells and can not be used > for "open shell" calculations, such as molecular ions?
You can try both "closed electronic shells" (spin-unpolarized) and "open shell" (spin-polarized) calculations.
Best regards,
T.Ozaki
|
Re: open shell systems ( No.2 ) |
- Date: 2005/12/05 16:35
- Name: peter
- Dear Prof. Ozaki,
below I am attaching the .dat file for my nh4 calculation. I tried to modify the settings according to your reccomendations, but failed to obtain zero dipole momentum. In fact my result for the dipole momentum (8.12) is very stable with respect to the settings modifications. Please have a look at the .dat file: ----------------------------------------------------------- # # SCF calculation of a methane molecule by the LDA # and the cluster method #
# # File Name # System.CurrrentDirectory ./ System.Name NH4 level.of.stdout 1 # default=1 (1-3) level.of.fileout 1 # default=1 (0-2)
# # Definition of Atomic Species #
Species.Number 2 <Definition.of.Atomic.Species H H4.0-s1 H_TM N N4.5-s1p1 N_TM_PCC Definition.of.Atomic.Species>
# # Atoms # Atoms.Number 5 Atoms.SpeciesAndCoordinates.Unit Ang # Ang|AU <Atoms.SpeciesAndCoordinates 1 N -0.976418 -0.975766 -0.97872 2.5 2.5 0.0 0.0 1 2 H -0.965599 0.057111 -0.990461 0.5 0.5 0.0 0.0 1 3 H -0.00646 -1.32998 -0.950363 0.5 0.5 0.0 0.0 1 4 H -1.48788 -1.30529 -0.143907 0.5 0.5 0.0 0.0 1 5 H -1.44574 -1.3249 -1.83015 0.5 0.5 0.0 0.0 1 Atoms.SpeciesAndCoordinates> Atoms.UnitVectors.Unit Ang # Ang|AU <Atoms.UnitVectors 10.0 0.0 0.0 0.0 10.0 0.0 0.0 0.0 10.0 Atoms.UnitVectors>
# # SCF or Electronic System #
scf.XcType LDA # LDA|LSDA-CA|LSDA-PW|GGA-PBE scf.SpinPolarization off # On|Off|NC scf.ElectronicTemperature 300.0 # default=300 (K) scf.energycutoff 200.0 # default=150 (Ry) scf.maxIter 100 # default=40 scf.EigenvalueSolver cluster # DC|GDC|Cluster|Band scf.lapack.dste dstevx # dstegr|dstedc|dstevx, default=dstegr scf.Kgrid 1 1 1 # means n1 x n2 x n3 scf.Mixing.Type rmm-diis # Simple|Rmm-Diis|Gr-Pulay|Kerker|Rmm-Diisk scf.Init.Mixing.Weight 0.30 # default=0.30 scf.Min.Mixing.Weight 0.001 # default=0.001 scf.Max.Mixing.Weight 0.400 # default=0.40 scf.Mixing.History 7 # default=5 scf.Mixing.StartPulay 5 # default=6 scf.criterion 1.0e-8 # default=1.0e-6 (Hartree)
scf.system.charge 1
# # 1D FFT #
1DFFT.NumGridK 900 # default=900 1DFFT.NumGridR 900 # default=900 1DFFT.EnergyCutoff 2500.0 # default=3DFFT.EnergyCutoff*3.0 (Ry)
# # Orbital Optimization #
orbitalOpt.Method off # Off|Unrestricted|Restricted orbitalOpt.InitCoes Symmetrical # Symmetrical|Free orbitalOpt.initPrefactor 0.1 # default=0.1 orbitalOpt.scf.maxIter 15 # default=12 orbitalOpt.MD.maxIter 7 # default=5 orbitalOpt.per.MDIter 20 # default=1000000 orbitalOpt.criterion 1.0e-4 # default=1.0e-4 (Hartree/borh)^2
# # output of contracted orbitals #
CntOrb.fileout off # on|off, default=off Num.CntOrb.Atoms 1 # default=1 <Atoms.Cont.Orbitals 1 Atoms.Cont.Orbitals> # # SCF Order-N #
orderN.HoppingRanges 6.0 # default=5.0 (Ang) orderN.NumHoppings 2 # default=2
# # MD or Geometry Optimization #
MD.Type Opt # Nomd|Opt|NVE|NVT_VS|NVT_NH # Constraint_Opt|DIIS MD.maxIter 1 # default=1 MD.TimeStep 1.0 # default=0.5 (fs) MD.Opt.criterion 1.0e-4 # default=1.0e-4 (Hartree/bohr)
# # restart using a restart file, *.rst #
scf.restart off # on|off,default=off
# # MO output #
MO.fileout off # on|off num.HOMOs 1 # default=1 num.LUMOs 1 # default=1 MO.Nkpoint 1 # default=1 <MO.kpoint 0.0 0.0 0.0 MO.kpoint>
# # DOS and PDOS #
Dos.fileout off # on|off, default=off Dos.Erange -10.0 10.0 # default = -20 20 Dos.Kgrid 1 1 1 # default = Kgrid1 Kgrid2 Kgrid3
# # output Hamiltonian and overlap #
HS.fileout off # on|off, default=off
|
Re: open shell systems ( No.3 ) |
- Date: 2005/12/05 16:40
- Name: peter
- Please have a look at the end of the openmx output for the .dat file above:
---------------------- ******************* MD= 1 SCF=11 ******************* <Poisson> Poisson's equation using FFT... <Set_Hamiltonian> Hamiltonian matrix for VNA+dVH+Vxc... <Cluster> Eigenvalue problem... Atom 1 MulP 2.9671 2.9671 sum 5.9343 Atom 2 MulP 0.2582 0.2582 sum 0.5163 Atom 3 MulP 0.2582 0.2582 sum 0.5164 Atom 4 MulP 0.2582 0.2582 sum 0.5165 Atom 5 MulP 0.2582 0.2582 sum 0.5165 <DFT> Total spin S = 0.000000000000 <DFT> Mixing_weight= 0.348032419010 <DFT> Uele = -6.086293721744 dUele = 0.000000002256 <DFT> NormRD = 0.000000000000 Criterion = 0.000000010000 <MD= 1> Force calculation Force calculation #1 Force calculation #2 Force calculation #3 Force calculation #5 <MD= 1> Double-counting Energy Force calculation #6 Force calculation #7
******************************************************* Dipole moment (Debye) *******************************************************
Absolute D 8.12725258
Dx Dy Dz Total -4.68953265 -4.68622176 -4.70104709 Core -42.20947180 -42.18122895 -42.30895582 Electron 37.51993916 37.49500719 37.60790873
******************************************************* Kohn-Sham's energies (Hartree) at MD = 1 *******************************************************
Uele = -6.086293721744 Uxc0 = 0.400301653350 Uxc1 = 0.400301653350 UH0 = -20.588765127127 UH1 = 1.578292054581 Ucore = 12.127643937873 Utot = -12.168519549718
******************************************************* Computational times (s) at MD = 1 *******************************************************
Set_OLP_Kin = 0.17895 Set_Nonlocal = 0.40981 Set_Hamiltonian = 0.71516 Poisson = 2.76630 diagonalization = 0.00572 Mixing_DM = 0.00060 Force = 0.87290 Correction_Energy = 0.65550
******************************************************* MD or geometry opt. at MD = 1 *******************************************************
<Steepest_Descent> SD_scaling= 0.930797663119 <Steepest_Descent> |Maximum force| (Hartree/Bohr) = 0.088589990501 <Steepest_Descent> Criterion (Hartree/Bohr) = 0.000100000000
atom= 1, XYZ(ang) Fxyz(a.u.)= -0.9764 -0.9758 -0.9787 0.0000 0.0000 -0.0000 atom= 2, XYZ(ang) Fxyz(a.u.)= -0.9656 0.0607 -0.9905 -0.0009 -0.0886 0.0010 atom= 3, XYZ(ang) Fxyz(a.u.)= -0.0031 -1.3312 -0.9503 -0.0830 0.0303 -0.0024 atom= 4, XYZ(ang) Fxyz(a.u.)= -1.4897 -1.3064 -0.1410 0.0438 0.0282 -0.0715 atom= 5, XYZ(ang) Fxyz(a.u.)= -1.4474 -1.3261 -1.8331 0.0402 0.0299 0.0729
outputting data on grids to files... ---------------------------------------------
As you can see, the molecule has a very large dipole momentum.
Another experiment I did is the following. I took your Methane.dat file from the work directory of your distribution. After scf.system.charge=1 setting, openmx finds the dipole moment 5.4e-5 (vs. 1.5e-7 for the neutral molecule). The results seems to be fine, does it mean that a small change in the orientation of the molecule (Methane.dat vs. our nh4.dat) can lead to such a drastic error in the scf calculation?
Could you please provide me with your .dat file for nh4+ calculation? Are there other reccomendations we could follow?
Thank you very much for your promt answers, Peter
|
Re: open shell systems ( No.4 ) |
- Date: 2005/12/06 01:18
- Name: peter
- note added:
when nitrogen atom is placed instead of carbon right in the methane.dat file, the dipole momentum of the nh4 molecule appears to be zero after openmx calculation. If a molecular dynamics is used to do geometry optimization, the bond length is fixed, but a small dipole momentum appears.
I tried to calculate hydroxile ion oh-. The .dat file for the calculation i obtained from your water.dat by removing hydrogen atom. The dipole momentum of the ion appears 0.68-0.77 depending on the basis set and a few settings. The "experimental value" can be found on the net and are all larger than one, e.g. 1.58D (http://www.colby.edu/chemistry/webmo/OH-.html).
An attempt to calculate ch3o- ion also leads to a giant dipole momentum, while the dipole momentum for the neutral molecule (ch3oh) is fine.
What is the source for this instability? Does it mean that all ions calculations are unstable with respect to the relative orientation of a molecular ion vs. the FFT "cavity"?
regards, Peter
|
Re: open shell systems ( No.5 ) |
- Date: 2005/12/06 02:35
- Name: T.Ozaki
- Hi,
Thanks for your detailed reports.
I noticed from your calculations a thing I did not notice for calculation of the dipole moment of charged systems.
For a charged system, the dipole moment is not well defined if the background contribution is not taken into account.
For instance, imagine a charged linear molecule "A1-A2-A3" along x-axis, whose x-coordinate is X-d, X, and X+d for A1, A2, and A3, respectively. And the charge is Z, 0, and Z for A1, A2, and A3, respectively. Then, the dipole moment is calculated as
Z(X+d) + 0(X) + Z(X-d) = 2ZX
We cleary see that the dipole moment of the charged molecule depends on the coordinate X, and in fact this is the source of the problem you met, while the background contribution was zero (?un)fortunately in my structures.
The problem is resolved by taking into account the contribution of the back ground charge. Then, the total system becomes neutral, and the coordinate dependency of the dipole moment disappears.
The following is the result including the contribution of the background charge calculated by your input file.
******************************************************* Dipole moment (Debye) *******************************************************
Absolute D 0.00071255
Dx Dy Dz Total 0.00041165 0.00057932 -0.00005157 Core -42.20947180 -42.18122895 -42.30895582 Electron 37.51993916 37.49500719 37.60790873 Back ground 4.68994430 4.68680108 4.70099552
We see that the back ground contribution compensates other contributions. Please download the following two modified files, in which the contribution of the background charge is taken into account, and replace the corresponding files in OpenMX2.73 (note that the two files are valid for only OpenMX2.73). http://staff.aist.go.jp/t-ozaki/bugfixed/05Dec06/Correction_Energy.c http://staff.aist.go.jp/t-ozaki/bugfixed/05Dec06/openmx_common.h
Regards,
T.Ozaki
|
Re: open shell systems ( No.6 ) |
- Date: 2005/12/08 00:56
- Name: peter
- Dear Prof. Ozaki,
thank you for your explanations. I played with the settings in the . dat file and managed to find very good dipole momentum for small organic molecules. The next test I did is the molecular polarization. I applied the electric field in the orthogonal directions and derived the polarizability tensor for a few simple molecules (up to 6 atoms). The common observation is that the polarizability seems to be 30-50% lower than I expected from experimental data. This holds even for molecules with good numerical values of dipole momentum. To be specific: I diagonalized the polarizability tensor to compare with the triples of polarizability values given in the literature. If the tensor components where not found, I used the mean value.
Polarizability must be sensitive to inclusion of higher AO in the basis set. Have you performed polarizability calculations with openmx? Is there a place I could read about it and figure out the proper settings?
thank you for your answers, Peter
|
Re: openmx runtime questions ( No.7 ) |
- Date: 2005/12/08 20:01
- Name: T.Ozaki
- Hi,
The calculation of polarizability has never been reported using OpenMX. I guess that a relatively high quality basis set will be required for the convergent result.
Regards,
T.Ozaki
|
|