This thread is locked.Only browsing is available.
Top Page > Browsing
open shell systems
Date: 2005/12/02 18:06
Name: peter

Dear Prof. Ozaki

thank you very much for your explanation. Let me ask another question. In order to try openmx for ion calculations, I got your example .dat file for ammonia molecule, added there another hydrogen (using the same basis set), set the system charge to 1 and performed the calculation. The inital occupations for the spin states of the added hydrogen atom (Atoms.SpeciesAndCoordinates of the .dat file) were set to 0.5 0.5.

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)?

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.

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?

Sorry if the answer to the question is covered somewhere in the manual, I tried my best to find it there.

many thanks to you again for your answers,
Peter
メンテ
Page: [1]

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
メンテ

Page: [1]