Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Refactor into solve_dirac_eigenproblem() #14

Merged
merged 13 commits into from
Oct 12, 2023
Merged

Conversation

certik
Copy link
Contributor

@certik certik commented Sep 21, 2023

Test using:

fpm test  test_dft_dirac

TODO:

  • Use solve_dirac_eigenproblem in schroed_dirac_solver.f90 in total_energy, test using fpm test test_coulomb_dirac.

@certik
Copy link
Contributor Author

certik commented Sep 26, 2023

TODO: finish the integration in schroed_dirac_solver.f90, more work is needed to convert or create some variables that are needed.

@certik
Copy link
Contributor Author

certik commented Sep 26, 2023

Here is the current status:

$ fpm test  test_coulomb_dirac
schroed_dirac_solver.f90               failed.
[  4%] Compiling...
././src/schroed_dirac_solver.f90:117:79:

  117 |         call solve_dirac_eigenproblem(Nb, Nq, Lmin2, Lmax, alpha, alpha_j, xe, xiq1, &
      |                                                                               1
Error: Rank mismatch in argument 'xiq_gj' at (1) (rank-2 and rank-1)
<ERROR> Compilation failed for object " src_schroed_dirac_solver.f90.o "
<ERROR> stopping due to failed compilation
STOP 1

We need to get all the arguments to fit well, but it's close.

@certik
Copy link
Contributor Author

certik commented Sep 28, 2023

Now the error is:

$ fpm test  test_coulomb_dirac
Project is up to date
  Z  rmax   Ne       a  p Nq DOFs
At line 208 of file ././src/feutils.f90
Fortran runtime error: Array bound mismatch for dimension 1 of array 'xq' (53/8)

Error termination. Backtrace:
#0  0x105473c07
#1  0x1054747a7
#2  0x105474be3
#3  0x104e94537
#4  0x104e724eb
#5  0x104ee0c4f
#6  0x104d96187
#7  0x104d973b7
<ERROR> Execution for object " test_coulomb_dirac " returned exit code  2
<ERROR> *cmd_run*:stopping due to failed executions
STOP 2

The fpm test test_dft_dirac works.

@certik
Copy link
Contributor Author

certik commented Oct 11, 2023

Both now run locally for me:

$ fpm test  test_coulomb_dirac
Project is up to date
  Z  rmax   Ne       a  p Nq DOFs
 92  50.0    7   100.0 23 53  322   -16817.941478426688

 Comparison of calculated and reference energies

 Total energy:
                   E               E_ref     error
 -16817.941478426688 -16817.941478425906  7.82E-10
 Eigenvalues:
   n                   E               E_ref     error
   1  -4861.198023119627  -4861.198023119372  2.55E-10
   2  -1257.395890258162  -1257.395890257889  2.73E-10
   3  -1089.611420919940  -1089.611420919875  6.55E-11
   4  -1257.395890257914  -1257.395890257889  2.55E-11
   5   -539.093341793847   -539.093341793890  4.37E-11
   6   -489.037087678211   -489.037087678200  1.09E-11
   7   -539.093341793920   -539.093341793890  2.91E-11
   8   -476.261595161170   -476.261595161155  1.46E-11
   9   -489.037087678265   -489.037087678200  6.55E-11
  10   -295.257844100335   -295.257844100397  6.18E-11
  11   -274.407758840138   -274.407758840065  7.28E-11
  12   -295.257844100372   -295.257844100397  2.55E-11
  13   -268.965877827133   -268.965877827130  3.64E-12
  14   -274.407758840051   -274.407758840065  1.46E-11
  15   -266.389447187808   -266.389447187816  7.28E-12
  16   -268.965877827151   -268.965877827130  2.18E-11
  17   -185.485191678570   -185.485191678552  1.82E-11
  18   -174.944613583484   -174.944613583462  2.18E-11
  19   -185.485191678556   -185.485191678552  3.64E-12
  20   -172.155252323741   -172.155252323737  3.64E-12
  21   -174.944613583466   -174.944613583462  3.64E-12
  22   -170.828937049904   -170.828937049879  2.55E-11
  23   -172.155252323744   -172.155252323737  7.28E-12
  24   -170.049934288581   -170.049934288552  2.91E-11
  25   -170.828937049868   -170.828937049879  1.09E-11
  26   -127.093638842613   -127.093638842631  1.82E-11
  27   -121.057538029556   -121.057538029549  7.28E-12
  28   -127.093638842634   -127.093638842631  3.64E-12
  29   -119.445271987152   -119.445271987141  1.09E-11
  30   -121.057538029556   -121.057538029549  7.28E-12
  31   -118.676410324359   -118.676410324351  7.28E-12
  32   -119.445271987148   -119.445271987141  7.28E-12
  33   -118.224144624928   -118.224144624903  2.55E-11
  34   -118.676410324362   -118.676410324351  1.09E-11
  35   -117.925825597205   -117.925825597293  8.73E-11
  36   -118.224144624917   -118.224144624903  1.46E-11
  37    -92.440787600950    -92.440787600943  7.28E-12
  38    -88.671749052017    -88.671749052017  0.00E+00
  39    -92.440787600939    -92.440787600943  3.64E-12
  40    -87.658287631901    -87.658287631893  7.28E-12
  41    -88.671749052020    -88.671749052017  3.64E-12
  42    -87.173966671959    -87.173966671948  1.09E-11
  43    -87.658287631890    -87.658287631893  3.64E-12
  44    -86.888766390963    -86.888766390941  2.18E-11
  45    -87.173966671951    -87.173966671948  3.64E-12
  46    -86.700519572765    -86.700519572809  4.37E-11
  47    -86.888766390948    -86.888766390941  7.28E-12
 Eigenfunctions saved in data_coulomb_dirac.txt
$ fpm test  test_dft_dirac    
Project is up to date
 SCF iteration:           1
 SCF iteration:           2
 SCF iteration:           3
 SCF iteration:           4
 SCF iteration:           5
 SCF iteration:           6
 SCF iteration:           7
 SCF convergence error:   2.4217968974735413     
 SCF iteration:           8
 SCF convergence error:   6.0885791899636388E-003
 SCF iteration:           9
 SCF convergence error:   1.8715260521275923E-003
 SCF iteration:          10
 SCF convergence error:   9.0274908870924264E-005
 SCF iteration:          11
 SCF convergence error:   2.0544301150948741E-005
 SCF iteration:          12
 SCF convergence error:   1.0121515515493229E-005
 SCF iteration:          13
 SCF convergence error:   1.9830640667350963E-006
 SCF iteration:          14
 SCF convergence error:   1.0529329301789403E-006
 SCF iteration:          15
 SCF convergence error:   2.0726292859762907E-007
 SCF iteration:          16
 SCF convergence error:   1.1235169949941337E-007
 SCF iteration:          17
 SCF convergence error:   2.7750502340495586E-008
 SCF iteration:          18
 SCF convergence error:   9.4005372375249863E-009
 SCF iteration:          19
 SCF convergence error:   1.3582393876276910E-008
 SCF iteration:          20
 SCF convergence error:   4.0745362639427185E-009
 SCF iteration:          21
 SCF convergence error:   5.2386894822120667E-010
 SCF iteration:          22
 SCF convergence error:   2.8194335754960775E-009
 SCF iteration:          23
 SCF convergence error:   5.1295501179993153E-010
 SCF iteration:          24
 SCF convergence error:   8.9676177594810724E-010
 Comparison of calculated and reference energies

 Total energy:
               E           E_ref     error
 -28001.13232549 -28001.13232549  3.50E-09

 Eigenvalues:
   n               E           E_ref     error
   1  -4223.41902046  -4223.41902046  7.39E-10
   2   -789.48978233   -789.48978233  3.64E-10
   3   -761.37447597   -761.37447597  1.16E-09
   4   -622.84809456   -622.84809456  9.65E-10
   5   -199.42980564   -199.42980564  3.32E-10
   6   -186.66371312   -186.66371312  5.05E-10
   7   -154.70102667   -154.70102667  5.29E-10
   8   -134.54118029   -134.54118029  5.30E-10
   9   -128.01665738   -128.01665738  3.25E-10
  10    -50.78894806    -50.78894806  4.39E-10
  11    -45.03717129    -45.03717129  4.76E-10
  12    -36.68861049    -36.68861049  5.20E-10
  13    -27.52930624    -27.52930624  5.32E-10
  14    -25.98542891    -25.98542891  3.98E-10
  15    -13.88951423    -13.88951423  4.42E-10
  16    -13.48546969    -13.48546969  4.65E-10
  17    -11.29558710    -11.29558710  5.25E-10
  18     -9.05796425     -9.05796425  5.15E-10
  19     -7.06929563     -7.06929563  5.03E-10
  20     -3.79741623     -3.79741623  5.41E-10
  21     -3.50121718     -3.50121718  4.89E-10
  22     -0.14678838     -0.14678838  4.82E-10
  23     -0.11604716     -0.11604717  4.93E-10
  24     -1.74803995     -1.74803995  6.13E-10
  25     -1.10111900     -1.10111900  5.50E-10
  26     -0.77578418     -0.77578418  5.36E-10
  27     -0.10304081     -0.10304082  4.89E-10
  28     -0.08480202     -0.08480202  4.94E-10
  29     -0.16094728     -0.16094728  2.96E-10

@certik
Copy link
Contributor Author

certik commented Oct 11, 2023

But this doesn't work yet:

$ fpm test  test_harmonic_dirac
Project is up to date
  Z  rmax   Ne       a  p Nq DOFs
 D(           1 ,           1 ) is exactly zero.
 The factorization has been completed, but the block diagonal
 matrix D is exactly singular, so the solution could not be
 computed.
ERROR STOP DSYSV ERROR.

Error termination. Backtrace:
#0  0x105203c07
#1  0x1052047a7
#2  0x105205a9b
#3  0x104d797cb
#4  0x104d2d82f
#5  0x104d14eff
#6  0x104d845ab
#7  0x104c381a3
#8  0x104c3949f
<ERROR> Execution for object " test_harmonic_dirac " returned exit code  1
<ERROR> *cmd_run*:stopping due to failed executions
STOP 1

The potential has to be handled correctly there, right now it assumes Coulomb.

@certik
Copy link
Contributor Author

certik commented Oct 11, 2023

Everything works for me locally now.

TODO for subsequent PRs

  • do the same refactoring for the Schroedinger solver
  • Lmin, Lmax refactoring (currently they are set in two different places)
  • further cleanup of the total_energy function

@certik certik marked this pull request as ready for review October 12, 2023 09:34
real(dp) :: E_dirac_shift
rmin = 0
allocate(xe(Ne+1), xiq(Nq), xiq1(Nq), wtq(Nq), wtq1(Nq), V(Nq, Ne), xiq_lob(p+1), wtq_lob(p+1))
allocate(phihq(Nq, p+1), dphihq(Nq, p+1), in(p+1, Ne), ib(p+1, Ne), xin(p+1), xq1(Nq, Ne))
allocate(phihq(Nq, p+1), dphihq(Nq, p+1), in(p+1, Ne), ib(p+1, Ne), xin(p+1), xq1(Nq,1))
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

1 is correct here, consistent with DFT. Inside solve_dirac_eigenproblem only the first element is used. This should eventually be converted to just a 1D array.

@certik certik merged commit e7747e6 into atomic-solvers:main Oct 12, 2023
4 checks passed
@certik certik deleted the refactor branch October 12, 2023 09:38
@certik certik mentioned this pull request Oct 12, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants