OrderN / CONQUEST-release

Full public release of large scale and linear scaling DFT code CONQUEST
http://www.order-n.org/
MIT License
96 stars 25 forks source link

Merge post-processing utility #92

Closed davidbowler closed 2 years ago

davidbowler commented 2 years ago

New post-processing utility for output: structure, charge density, band density, STM simulation, DOS.

tsuyoshi38 commented 2 years ago

Let me ask you a question. I thought this branch was made from the latest "develop" branch. If so, I cannot find XC_LibXC_v5_module.f90. Is it okay?

davidbowler commented 2 years ago

This branch started some time ago, so is behind the latest develop, which explains the lack of LibXC5. However, when it is merged into develop it will only add to the code, and not remove the new changes.

tsuyoshi38 commented 2 years ago

Thank you for your reply. I am relieved.

tsuyoshi38 commented 2 years ago

A few comments:

  1. When reading the manual for the post-processing part, I thought we should make an input file (or standard input) specific for the utility PostProcessCQ. It may be better to explain that the input file for this utility is also Conquest_input.
  2. When I tried to make a cube file for total charge density, I noticed that I have no file "chden.***" because IO.Iprint was 1 to run CONQUEST. I think it is better to use IO.DumpChargeDensity instead of increasing IO.Iprint. But...
  3. I am sorry, but it is possible that I did not consider carefully the timing when we should print out chden. files. I thought it was not so important to print out these files, if we have matrix files (Kmatrix2. etc.). I have noticed that we use "if (flag_DumpChargeDensity .or. iprint_SC > 2)" in H_matrix_module, while in SelfCon_module, we use "if (flag_DumpChargeDensity .and. iprint_SC > 1)" . I am not sure whether these are okay or not.
  4. I think we should have a clear idea when we print out chden.*** files.

I just wonder for 3 or 4, we should change them now or later.

tsuyoshi38 commented 2 years ago

I have found that Conquest will not produce Process0000001WF.dat etc. if I use the following Conquest_input.

If the lines

IO.min_wf_E -1.0

IO.max_wf_E 1.0

are activated, I could get the files.

Is it a correct behavior?

---- Conquest_input ------- IO.outputWF T IO.maxnoWF 6

IO.min_wf_E -1.0

IO.max_wf_E 1.0

%block WaveFunctionsOut 62 63 64 65 66 67 %endblock ---- Conquest_input -------

davidbowler commented 2 years ago

No, this is not correct behaviour - thank you for finding it! I can see where the problem is (in initial_read_module.f90) and will fix it.

tsuyoshi38 commented 2 years ago

Thank you very much for the update. I have confirmed that the problems I pointed out above have been solved.

I have some questions, comments or suggestions for DOS.

  1. I think users sometimes want to set only maximum value (IO.max_DOS_E) with the default setting for the minimum value (the lowest eigenvalue). But, it does not work. In such cases, IO.min_DOS_E is set to 0. (initial default value)
  2. I wonder the default value of minimum value should be a little lower (by 4xsigma, for example) than the lowest eigenvalue so that we can see the edge of the DOS from the lowest eigen value.
  3. I think the default value of sigma should be at least 3-4 times bigger than the grid-spacing for the energy. It is especially important when the levels are discrete. I think it is also important to check whether the integrated DOS correctly shows the number of electrons or not.
  4. Related to 3, I wonder it is better to show the integrated DOS simultaneously. (I mean, for E< IO.min_DOS_E, we simply count the number of states, and then sum up the electrons for each bin. )

I am sorry, but I have not checked the STM image. Others look perfect to me.

davidbowler commented 2 years ago

I have now updated DOS as discussed: improving default parameters (1001 bins, min/max eigenvalues, sigma); added integrated DOS (including importing smearing from Conquest for Fermi or Methfessel-Paxton).

I have also refactored the storage to account for bands output by Conquest and bands to be processed (which do not have to be the same set, though the second has to be a subset of the first).

davidbowler commented 2 years ago

Also the output of chden.nnn files from Conquest is now turned on by IO.DumpChargeDensity and happens at SCF, and optionally at each SCF step if IO.Iprint_SC is larger than 2.

tsuyoshi38 commented 2 years ago

Also the output of chden.nnn files from Conquest is now turned on by IO.DumpChargeDensity and happens at SCF, and optionally at each SCF step if IO.Iprint_SC is larger than 2.

This sounds very reasonable to me.

tsuyoshi38 commented 2 years ago

Thank you very much for the update of DOS procedure. It looks wonderful, but I have found a small problem and suggestion.

  1. In the subroutine process_dos, E_DOS_min or E_DOS_max are shifted by 2xsigma_DOS. But, this should be done AFTER the setting of sigma_DOS, I think.
  2. I think integrated DOS should not consider the occupation. I am sorry if this is because my comment was unclear. I thought, when the lower boundary is higher than the lowest eigenvalues, it is better to add the number of states below the lower boundary (E_DOS_min) to the integrated DOS. (It may be a little difficult to consider this when E_DOS_min is within the dispersion of some band.)

I will send an e-mail showing the corresponding part.

tsuyoshi38 commented 2 years ago

Thank you very much for the update of DOS related tools. Let me ask you one more. You have changed the keywords for DOS, from IO.min_DOS_E etc. to Process.min_DOS_E. I agree that it should be Process.***, but I have noticed that the explanation in manual has not been changed. Please change them.

As far as I have checked, others are fine and let me approve this pull request. Thank you very much for your efforts.

davidbowler commented 2 years ago

I have made the changes to the documentation of the parameter names (Process.XXX) and added an option to allow integrated DOS to be local (i.e. just within the energy window) or total (i.e. from lowest eigenvalue).

I've also realised that the band density and STM are not (yet) compatible with MSSF, so I've added a check to quit if post-processing is attempted on these areas with MSSF.