Copy from here https://mathematica.stackexchange.com/questions/284809/why-does-the-minimum-eigenvalue-change-dramatically-when-one-basis-function-is-a
I have a basis set which describes with high accuracy the ground state of this system:$$H=-\frac{1}{2}\Delta-\frac{1}{r}+\frac{25}{8}\rho^2-5/2$$where $r=(\rho,z,\phi)$ is a coordinate in the cylindrical system.
Ground state from the literature $E_{min}=-1.3803975$
The anisotropic gaussian basis set which I use to solve this task:$$\psi_j=e^{-b_{0j}\;z^2}e^{-a_{0j}\;\rho^2}$$where $a_{0j}$ and $b_{0j}$ are parameters
$a_{0j}$={1.25, 1.25003, 1.25006, 1.251393, 1.252726, 1.2587955, 1.264865, 1.260499, 1.303186, 1.546232, 2.337185, 4.507307, 9.426298, 20.270036, 43.977048, 95.675215, 208.32642, 453.736619, 988.32236, 2152.803005, 4689.357171, 10214.647096}
$b_{0j}$={0.026284, 0.0417685, 0.057253, 0.09098300000000001, 0.124713, 0.19818550000000001, 0.271658,0.417292, 0.61588, 1.008925, 1.949878, 4.24735, 9.251849, 20.15297, 43.898488,95.622497, 208.291042, 453.712878`, 988.306428, 2152.792314,4689.349997, 10214.642282}
This basis gives the energy of the ground state $E_{min}=-1.3803971$, which agrees perfectly with the value from the literature.
Now if we add only one element to the list $a_{0j}$ and list $b_{0j}$ (1.2626819999999999 and 0.344475 respectively), then the absolute value of minimum energy becomes huge ($E_{min}=-7.1824042$).
Why is this happening?
It is known from theory that an increase in the number of basis functions should either increase the accuracy or not change the values of the eigenvalues.
The code had written in Wolfram Mathematica. In the code I have renamed $\rho ≡ r$
ClearAll["Global`*"]
b00 = {0.026284`, 0.0417685`, 0.057253`, 0.09098300000000001`,
0.124713`, 0.19818550000000001`, 0.271658`,(*0.344475`,*)0.417292`,
0.61588`, 1.008925`, 1.949878`, 4.24735`, 9.251849`, 20.15297`,
43.898488`, 95.622497`, 208.291042`, 453.712878`, 988.306428`,
2152.792314`, 4689.349997`, 10214.642282`};
a00 = {1.25`, 1.25003`, 1.25006`, 1.251393`, 1.252726`, 1.2587955`,
1.264865`,(*1.2626819999999999`,*)1.260499`, 1.303186`, 1.546232`,
2.337185`, 4.507307`, 9.426298`, 20.270036`, 43.977048`,
95.675215`, 208.32642`, 453.736619`, 988.32236`, 2152.803005`,
4689.357171`, 10214.647096`};
nmax = Dimensions[b00][[1]];
a0 [j_] := a00[[j]];
b0 [j_] := b00[[j]];
Psi[r_, z_, j_] := Exp[-b0[j]*z^2]*Exp[-a0[j]*r^2];
Kk[r_, z_, j1_, j2_] :=
FullSimplify[
Psi[r, z, j2]*
Laplacian[Psi[r, z, j1], {r, \[Theta], z}, "Cylindrical"]*r*2*Pi];
Kx[j1_, j2_] := -(1/2)*
NIntegrate[
Kk[r, z, j1, j2], {r, 0, Infinity}, {z, -Infinity, Infinity}];
Kxx = Table[Kx[j1, j2], {j1, 1, nmax}, {j2, 1, nmax}];
Ka[r_, z_, j1_, j2_] :=
FullSimplify[Psi[r, z, j2]*25/8*r^2*Psi[r, z, j1]*r*2*Pi];
KA[j1_, j2_] :=
KA[j1, j2] =
KA[j2, j1] =
NIntegrate[
Ka[r, z, j1, j2], {r, 0, Infinity}, {z, -Infinity, Infinity}];
KAx = Table[KA[j1, j2], {j1, 1, nmax}, {j2, 1, nmax}];
Ks[r_, z_, j1_, j2_] :=
FullSimplify[Psi[r, z, j2]*(-5/2)*Psi[r, z, j1]*r*2*Pi];
KS[j1_, j2_] :=
KS[j1, j2] =
KS[j2, j1] =
NIntegrate[
Ks[r, z, j1, j2], {r, 0, Infinity}, {z, -Infinity, Infinity}];
KSx = Table[KS[j1, j2], {j1, 1, nmax}, {j2, 1, nmax}];
VP1[r_, z_] := -(1/Sqrt[r^2 + z^2]);
Px1[j1_, j2_] :=
Px1[j1, j2] =
Px1[j2, j1] =
NIntegrate[
Psi[r, z, j2]*VP1[r, z]*Psi[r, z, j1]*r*2*Pi, {r, 0,
Infinity}, {z, -Infinity, Infinity}];
Pxx1 = Table[Px1[j1, j2], {j1, 1, nmax}, {j2, 1, nmax}];
Bb[j1_, j2_] :=
Bb[j1, j2] =
Bb[j2, j1] =
NIntegrate[
Psi[r, z, j2]*Psi[r, z, j1]*r*2*Pi, {r, 0,
Infinity}, {z, -Infinity, Infinity}];
Bxx = Table[Bb[j1, j2], {j1, 1, nmax}, {j2, 1, nmax}];
EE1 = Eigenvalues[{Kxx + Pxx1 + KAx + KSx, Bxx}];
Min[EE1]
Out[48]= -1.3804
