N.A. in output
MartinBeseda opened this issue · 7 comments
Hello,
I'm a completely new user of Duo and I'm puzzled by one thing.
I wanted to compute energy levels and spectroscopic constants of N2+ in SigmaU+ symmetry.
The problem is, I'm not getting N.A. symbol instead of the desired spectroscopic symols.
Equilibrium properties
imin = 6
r_imin / ang = 1.0595959596
V(r_imin) / cm-1 = -391.89938844
Has a single min? No
True r_e / ang = N.A.
Derivatives at true re
der0, cm-1 = N.A.
der1, cm-1/ang = N.A.
der2, cm-1/ang^2 = N.A.
der3, cm-1/ang^3 = N.A.
der4, cm-1/ang^4 = N.A.
Harmonic we, cm-1= N.A.
Rotat. B0, cm-1= N.A.
Anharm. const. xe= N.A.
Coriol. ae, cm-1= N.A.
Centr. De, cm-1= N.A.
Y00, cm-1= N.A.
Approximate J=0 vibrational energy levels (no couplings)
given by E(v, J=0) = V(re) + Y00 + we*(v+0.5) - we*xe*(v+0.5)^2
v
0 N.A.
1 N.A.
2 N.A.
3 N.A.
Spin = 0.5
|Lambda| = 1.0
Physical J_min = 0.5
Approximate J=J_min vibrational-rotational energy levels (no couplings)
given by E(v, J) = E(v, J=0) B0*J*(J+1)- ae*(v+0.5)*J*(J+1) - De*[J(J+1)]^2
v
0 N.A.
1 N.A.
2 N.A.
I suppose, that it's caused by a detection of multiple minima in my data. The problem is, there should be just one real minima, as can be seen in the plot.
Could you, please, tell me, what am I doing wrong in the input?
Complete output file
(Transcript of the input --->)
(N2+ 2SigmaU+)
masses 14 14
nstates 1
jrot -2 2
symmetry C2h
grid
npoints 100
range 0.85 5
type 0
end
vibrationalbasis
vmax 1
end
eigensolver
nroots 20
end
potential 1
name "2 SigmaU+"
symmetry u +
lambda 1
mult 2
spin 0.5
type grid
interpolationtype Cubicsplines
units angstrom eH
values
0.8 0.4663431999999972
0.85 0.2716072000000054
0.9 0.14446329999999818
0.95 0.06631509999999707
1 0.023365799999993442
1.05 0
1.1 -0.00121289999999874
1.15 0.008618799999993598
1.2 0.0252145999999982
1.25 0.04443120000000533
1.3 0.06253599999999437
1.35 0.0777470999999963
1.4 0.08978539999999668
1.45 0.099500699999993
1.5 0.1078013999999996
1.55 0.11524179999999262
1.6 0.12209389999999587
1.65 0.1284830000000028
1.7 0.13446980000000508
1.75 0.14007789999999432
1.8 0.14535119999999324
1.85 0.1503034999999926
1.9 0.15495169999999803
1.95 0.1593098000000026
2 0.16338980000000447
2.05 0.16720200000000318
2.1 0.17075570000000084
2.15 0.17405920000000208
2.2 0.1771202000000045
2.25 0.1799457999999987
2.3 0.18254299999999546
2.35 0.18491869999999722
2.4 0.18708049999999332
2.45 0.1890367000000026
2.5 0.19079659999999876
2.55 0.19237090000000023
2.6 0.1937710999999922
2.65 0.1950096000000059
2.7 0.1960983999999968
2.75 0.1970500000000044
2.8 0.19787619999999606
2.85 0.198588799999996
2.9 0.1991991000000013
2.95 0.19971769999999367
3 0.20015519999999754
3.05 0.20054960000000221
3.1 0.20084439999999404
3.15 0.2010855000000049
3.2 0.20128060000000403
3.25 0.20143659999999386
3.3 0.2015597999999983
3.35 0.20165559999999516
3.4 0.20172870000000387
3.45 0.20178300000000604
3.5 0.20182219999999518
3.55 0.20184899999999573
3.6 0.20186619999999778
3.65 0.20167109999999866
3.7 0.20167419999999936
3.75 0.20167290000000548
3.8 0.20166840000000263
3.85 0.2016617000000025
3.9 0.20165350000000615
3.95 0.20164450000000045
4 0.20163519999999835
4.05 0.20162589999999625
4.1 0.20161690000000476
4.15 0.20160850000000607
4.2 0.20160079999999425
4.25 0.20159400000000005
4.3 0.201588000000001
4.35 0.2015829999999994
4.4 0.20157890000000123
4.45 0.20157580000000053
4.55 0.20157229999999515
4.6 0.20157179999999641
4.65 0.20157209999999282
4.7 0.2015731000000045
4.75 0.20157480000000305
4.8 0.20157720000000268
4.85 0.20157999999999276
4.9 0.20158340000000408
4.95 0.2015872999999999
5 0.20159160000000043
30 0.20215919999999699
end
(<--- End of the input.)
____ ____ _ _ ___
| _ \ _ __ ___ __ _ _ __ __ _ _ __ ___ | _ \| | | |/ _ \
| |_) | '__/ _ \ / _` | '__/ _` | '_ ` _ \ | | | | | | | | | |
| __/| | | (_) | (_| | | | (_| | | | | | | | |_| | |_| | |_| |
|_| |_| \___/ \__, |_| \__,_|_| |_| |_| |____/ \___/ \___/
|___/
Please refer to:
Sergei N. Yurchenko, Lorenzo Lodi, Jonathan Tennyson and Andrey V. Stolyarov,
`DUO: a general program for calculating spectra of diatomic molecules',
Computer Physics Communication, (to be submitted), 2015.
Contacts: s.yurchenko@ucl.ac.uk; l.lodi@ucl.ac.uk; j.tennyson@ucl.ac.uk;
avstol@phys.chem.msu.ru
Values of physical constants used by DUO:
Planck constant h = 6.626069570000E-27 erg*second
Speed of light c = 29979245800.0000 cm/second
Bohr radius a0 = 0.5291772109200 angstroms
Hartree energy Eh = 219474.6313708000 cm^-1
Unified atomic mass unit u = 1.660538921000E-24 grams
Unified atomic mass unit u = 1822.8884861185961 me
Electron charge e = 1.602176565000E-19 Coulombs
Boltzmann constant kB = 1.380648800000E-16 erg/Kelvin
Avogadro constant = 6.022141290000E+23 mol^-1
Fine structure constant 1/alpha = 137.03599907466793
Electron mass me = 9.109382903261E-28 grams
a.u. of dipole e*a0 = 2.541746363812 debyes
1 eV = 8065.5442959967 cm^-1
Reduced mass is 7.00000000000000 atomic mass units (Daltons)
Reduced mass is 12760.2194028302 atomic units (electron masses)
Using a uniformly spaced grid with npoints = 100
rmin, rmax, step (in ang) = 0.85000000000000 5.00000000000000 0.04191919191919
rmin, rmax, step (in bohrs) = 1.60626720588030 9.44863062282531 0.07921579209035
Generate a grid representation for all Hamiltonian terms
...done!
Potential functions:
1 2 SigmaU+ type = GRID
r(Ang) 1
0.85000000 5.96108901E+04
0.89191919 3.53189271E+04
0.93383838 1.91364265E+04
0.97575758 8.99058349E+03
1.01767677 2.82452296E+03
1.05959596 -3.91899388E+02
1.10151515 -2.25374736E+02
1.14343434 1.51274429E+03
1.18535354 4.36117328E+03
1.22727273 7.82374126E+03
1.26919192 1.13348204E+04
1.31111111 1.45283006E+04
1.35303030 1.72424716E+04
1.39494949 1.94656822E+04
1.43686869 2.13149816E+04
1.47878788 2.29153448E+04
1.52070707 2.43544961E+04
1.56262626 2.56831723E+04
1.60454545 2.69279109E+04
1.64646465 2.81026012E+04
1.68838384 2.92150897E+04
1.73030303 3.02679469E+04
1.77222222 3.12667133E+04
1.81414141 3.22153268E+04
1.85606061 3.31149423E+04
1.89797980 3.39680011E+04
1.93989899 3.47762295E+04
1.98181818 3.55412089E+04
2.02373737 3.62643273E+04
2.06565657 3.69468352E+04
2.10757576 3.75898847E+04
2.14949495 3.81945240E+04
2.19141414 3.87617458E+04
2.23333333 3.92924672E+04
2.27525253 3.97875985E+04
2.31717172 4.02480264E+04
2.35909091 4.06746694E+04
2.40101010 4.10685374E+04
2.44292929 4.14307117E+04
2.48484848 4.17623801E+04
2.52676768 4.20648984E+04
2.56868687 4.23397119E+04
2.61060606 4.25883754E+04
2.65252525 4.28124940E+04
2.69444444 4.30136088E+04
2.73636364 4.31933546E+04
2.77828283 4.33532576E+04
2.82020202 4.34948664E+04
2.86212121 4.36196226E+04
2.90404040 4.37290901E+04
2.94595960 4.38245410E+04
2.98787879 4.39065934E+04
3.02979798 4.39826049E+04
3.07171717 4.40460431E+04
3.11363636 4.40955668E+04
3.15555556 4.41384255E+04
3.19747475 4.41740461E+04
3.23939394 4.42035995E+04
3.28131313 4.42279262E+04
3.32323232 4.42477068E+04
3.36515152 4.42636507E+04
3.40707071 4.42762234E+04
3.44898990 4.42860332E+04
3.49090909 4.42937717E+04
3.53282828 4.42976945E+04
3.57474747 4.43073206E+04
3.61666667 4.42916287E+04
3.65858586 4.42583558E+04
3.70050505 4.42624477E+04
3.74242424 4.42625506E+04
3.78434343 4.42612499E+04
3.82626263 4.42604591E+04
3.86818182 4.42589779E+04
3.91010101 4.42574437E+04
3.95202020 4.42557702E+04
3.99393939 4.42540594E+04
4.03585859 4.42523441E+04
4.07777778 4.42506594E+04
4.11969697 4.42490515E+04
4.16161616 4.42475431E+04
4.20353535 4.42461493E+04
4.24545455 4.42448973E+04
4.28737374 4.42437651E+04
4.32929293 4.42427836E+04
4.37121212 4.42419480E+04
4.41313131 4.42412548E+04
4.45505051 4.42407173E+04
4.49696970 4.42403199E+04
4.53888889 4.42400550E+04
4.58080808 4.42399173E+04
4.62272727 4.42399059E+04
4.66464646 4.42400113E+04
4.70656566 4.42402216E+04
4.74848485 4.42405411E+04
4.79040404 4.42409728E+04
4.83232323 4.42414664E+04
4.87424242 4.42420418E+04
4.91616162 4.42427079E+04
4.95808081 4.42434465E+04
5.00000000 4.42442421E+04
Equilibrium properties
imin = 6
r_imin / ang = 1.0595959596
V(r_imin) / cm-1 = -391.89938844
Has a single min? No
True r_e / ang = N.A.
Derivatives at true re
der0, cm-1 = N.A.
der1, cm-1/ang = N.A.
der2, cm-1/ang^2 = N.A.
der3, cm-1/ang^3 = N.A.
der4, cm-1/ang^4 = N.A.
Harmonic we, cm-1= N.A.
Rotat. B0, cm-1= N.A.
Anharm. const. xe= N.A.
Coriol. ae, cm-1= N.A.
Centr. De, cm-1= N.A.
Y00, cm-1= N.A.
Approximate J=0 vibrational energy levels (no couplings)
given by E(v, J=0) = V(re) + Y00 + we*(v+0.5) - we*xe*(v+0.5)^2
v
0 N.A.
1 N.A.
2 N.A.
3 N.A.
Spin = 0.5
|Lambda| = 1.0
Physical J_min = 0.5
Approximate J=J_min vibrational-rotational energy levels (no couplings)
given by E(v, J) = E(v, J=0) B0*J*(J+1)- ae*(v+0.5)*J*(J+1) - De*[J(J+1)]^2
v
0 N.A.
1 N.A.
2 N.A.
3 N.A.
Construct the J=0 matrix
Solving one-dimentional Schrodinger equations using : SINC
Vibrational (contracted) energies:
N Energy/cm State v
1 0.000000 [ 1 0 ] 2 SigmaU+
j = 2.0
Define the quanta book-keeping
The total number sigma/lambda states (size of the sigma-lambda submatrix) = 4
Sigma-Lambda basis set:
i jrot state spin sigma lambda omega
1 2.0 1 0.5 -0.5 1 0.5 2 SigmaU+
2 2.0 1 0.5 0.5 1 1.5 2 SigmaU+
3 2.0 1 0.5 -0.5 -1 -1.5 2 SigmaU+
4 2.0 1 0.5 0.5 -1 -0.5 2 SigmaU+
...done!
Contracted basis set:
i jrot ilevel ivib state v spin sigma lambda omega Name
1 2.0 1 1 1 0 0.5 -0.5 1 0.5 2 SigmaU+
2 2.0 2 1 1 0 0.5 0.5 1 1.5 2 SigmaU+
3 2.0 3 1 1 0 0.5 -0.5 -1 -1.5 2 SigmaU+
4 2.0 4 1 1 0 0.5 0.5 -1 -0.5 2 SigmaU+
Memory Report:
Active Arrays size (GiB)
GRID 0.000004
r-field 0.000003
contrfunc 0.000075
Total memory = 0.00008427 GiB
Maximal memory = 0.00020479 GiB ( 100000.0)
( 19 arrays contributing less than 1% are not shown)
Construct the hamiltonian matrix
...done!
Non-zero matrix elements of the coupled Sigma-Lambda matrix:
Threshold for printing is 2.22E-16 cm^-1
RV == Rotational-vibrational
SO == Spin-Orbit interaction
SS == Spin-Spin interaction
JS == J.S interaction (aka S-uncoupling)
LS == L.J interaction (aka L-uncoupling)
LS == L.S interaction (spin-electronic)
RV= 677.001; JS( 1 2)= -4.7036 ;
RV= 672.895;
RV= 672.895; JS( 3 4)= -4.7036 ;
RV= 677.001;
Eigenvalues for J = 2.0
J N Energy/cm State v lambda spin sigma omega parity
2.0 1 669.816022 1 0 1 0.5 0.5 1.5 + ||2 SigmaU+
2.0 2 680.080211 1 0 1 0.5 -0.5 0.5 + ||2 SigmaU+
J N Energy/cm State v lambda spin sigma omega parity
2.0 1 669.816022 1 0 -1 0.5 -0.5 -1.5 - ||2 SigmaU+
2.0 2 680.080211 1 0 -1 0.5 0.5 -0.5 - ||2 SigmaU+
Zero point energy (ZPE) = 0.000000
Memory Report:
Active Arrays size (GiB)
GRID 0.000004
r-field 0.000003
contrfunc 0.000075
Total memory = 0.00008338 GiB
Maximal memory = 0.00020479 GiB ( 100000.0)
( 25 arrays contributing less than 1% are not shown)
Timing data at 2019/04/23 15:38:11
Total time (seconds) Self time (seconds)
Timer Calls -------------------- -------------------
----- ----- Real CPU Real CPU
Duo 1. 0.0 0.0 0.0 0.0
Map on grid 1. 0.0 0.0 0.0 0.0
Solve vibrational part 1. 0.0 0.0 0.0 0.0
( 7 timers contributing less than 1% are not shown)
Dear Martin,
This feature is only working for analytical representations. It has not being implemented for the grid-type functions, sorry. It is on my todo list.
Cheers,
Sergey
Dear Sergey,
so, if I understand you well, you don't support piecewise-defined functions currently, do you? I was thinking about interpolating the dataset manually and providing the resulting function to DUO.
Cheers,
Martin
We do support piecewise-defined functions all right, just not for this feature of estimating the spectroscopic constants. The latter have never being a priority as the main purpose of Duo is a direct solution of the Schroedinger equation. Besides, spectroscopic constants, at least of the higher order, are heavily empirical. There is no much sense, from my experience, deriving them from the potential energy curves.
Cheers,
Sergey
Well, I see the sense of deriving them from PESs in two things: 1) Checking the PES comparing the constants with already known experients 2) Providing them to experimentalists, if they have no previous results.
And just off-topic question - which software did you use for the diatomic spectroscopic constants computations then?
Yes, I am aware of these applications. (1) Personally I prefer comparing energies not constants. (2) I do see sometimes (not very often) experimentalists trying to use ab initio spectroscopic constants, but not with much success.
RE: spectroscopic constants, I don't compute them at all. Spectroscopic constants only make sense as part of the perturbative/empirical methodology based on the effective rotational Hamiltonians for computing energies or line positions. We tend to use the variational methodology, based on the PECs (spin-orbit curves etc) directly. If I need to do a effective rotational Hamoltonian treatment I use PGOPHER.
Ok, I understand. Still, I'd like to compute them at least for this one case, but I understand, that it isn't your priority. Thank you for the prompt responses!
Martin
I will surely introduce this feature at some point. As a quick fix I can suggest to represent your potential energy curve as, e.g. a Morse potential, or even as a polynomial, or something like this
poten 1
name "X2Sigmau+"
symmetry u +
lambda 1
mult 2
type EMO
values
V0 0.00000000000000E+00
RE 1.08435981870765E+00
DE 4.44521643363054E+04
RREF -1.00000000000000E+00
PL 8.00000000000000E+00
PR 8.00000000000000E+00
NL 1.00000000000000E+00
NR 3.00000000000000E+00
B0 3.53434543259561E+00 fit
B1 3.41588567189263E-01 fit
B2 2.36501323107290E+00 fit
B3 -3.81206035003490E+00 fit
end
The spectroscopic constants only depend on the firs few derivatives, i.e. it is only the very bottom of the potential you need to fit accurately. However, looking at you potential, you do want to have more points around the equilibrium (below ~5000 cm-1) in order to have spectroscopic constants accurate. Otherwise the uncertainly will be too high. You can try using Duo to fit to ab initio data, see example below, but it will only work with a better coverage around equilibrium and accurate (smooth) ab initio data. A side note, you should use accurate masses or switch to "atoms" (see below).
Cheers,
Sergey
(N2+ 2SigmaU+)
atoms N N
nstates 1
jrot -2 2
symmetry C2v
grid
npoints 100
range 0.85 5
type 0
end
vibrationalbasis
vmax 1
end
eigensolver
nroots 50
end
poten 1
name "X2Sigmau+"
symmetry u +
lambda 1
mult 2
type EMO
values
V0 0.00000000000000E+00
RE 1.08435981870765E+00
DE 4.44521643363054E+04
RREF -1.00000000000000E+00
PL 8.00000000000000E+00
PR 8.00000000000000E+00
NL 1.00000000000000E+00
NR 3.00000000000000E+00
B0 3.53434543259561E+00 fit
B1 3.41588567189263E-01 fit
B2 2.36501323107290E+00 fit
B3 -3.81206035003490E+00 fit
end
DO_NOT_SHIFT_PECS
abinitio potential 1
name "2 SigmaU+"
symmetry u +
lambda 1
mult 2
spin 0.5
type grid
interpolationtype Cubicsplines
units angstrom eH
weighting ps1997 0.005 50000.0
values
0.8 0.4663431999999972
0.85 0.2716072000000054
0.9 0.14446329999999818
0.95 0.06631509999999707
1 0.023365799999993442
1.05 0
1.1 -0.00121289999999874
1.15 0.008618799999993598
1.2 0.0252145999999982
1.25 0.04443120000000533
1.3 0.06253599999999437
1.35 0.0777470999999963
1.4 0.08978539999999668
1.45 0.099500699999993
1.5 0.1078013999999996
1.55 0.11524179999999262
1.6 0.12209389999999587
1.65 0.1284830000000028
1.7 0.13446980000000508
1.75 0.14007789999999432
1.8 0.14535119999999324
1.85 0.1503034999999926
1.9 0.15495169999999803
1.95 0.1593098000000026
2 0.16338980000000447
2.05 0.16720200000000318
2.1 0.17075570000000084
2.15 0.17405920000000208
2.2 0.1771202000000045
2.25 0.1799457999999987
2.3 0.18254299999999546
2.35 0.18491869999999722
2.4 0.18708049999999332
2.45 0.1890367000000026
2.5 0.19079659999999876
2.55 0.19237090000000023
2.6 0.1937710999999922
2.65 0.1950096000000059
2.7 0.1960983999999968
2.75 0.1970500000000044
2.8 0.19787619999999606
2.85 0.198588799999996
2.9 0.1991991000000013
2.95 0.19971769999999367
3 0.20015519999999754
3.05 0.20054960000000221
3.1 0.20084439999999404
3.15 0.2010855000000049
3.2 0.20128060000000403
3.25 0.20143659999999386
3.3 0.2015597999999983
3.35 0.20165559999999516
3.4 0.20172870000000387
3.45 0.20178300000000604
3.5 0.20182219999999518
3.55 0.20184899999999573
3.6 0.20186619999999778
3.65 0.20167109999999866
3.7 0.20167419999999936
3.75 0.20167290000000548
3.8 0.20166840000000263
3.85 0.2016617000000025
3.9 0.20165350000000615
3.95 0.20164450000000045
4 0.20163519999999835
4.05 0.20162589999999625
4.1 0.20161690000000476
4.15 0.20160850000000607
4.2 0.20160079999999425
4.25 0.20159400000000005
4.3 0.201588000000001
4.35 0.2015829999999994
4.4 0.20157890000000123
4.45 0.20157580000000053
4.55 0.20157229999999515
4.6 0.20157179999999641
4.65 0.20157209999999282
4.7 0.2015731000000045
4.75 0.20157480000000305
4.8 0.20157720000000268
4.85 0.20157999999999276
4.9 0.20158340000000408
4.95 0.2015872999999999
5 0.20159160000000043
30 0.20215919999999699
end
FITTING
JLIST 0
itmax 20
fit_factor 1e-12
robust 0.0
output ESigma_Elander
energies (J parity NN energy ) (e-state v ilambda isigma omega weight)
0 + 1 0.00 1 0 1 0 0 1.0
end