QuantumLab-ZY / HamGNN

An E(3) equivariant Graph Neural Network for predicting electronic Hamiltonian matrix
GNU General Public License v3.0
54 stars 13 forks source link

Attempted to utilize the model from the Twist Bilayer MoS2 (TBM) Demo to predict angles for other configurations, but was unsuccessful. #7

Open newplay opened 8 months ago

newplay commented 8 months ago

Dear Yang Zhong: I'm tried use the HamGNN model to predict $21.79^\circ$ twist bilayer MoS2 but the result seems very different with the currect result. There's my input below: openmx input file (using poscar2openmx.py):

#
#      File Name      
#

System.CurrrentDirectory         ./    # default=./
System.Name                     MoS2_21_79_AA
DATA.PATH           /home/zjlin/openmx3.9/DFT_DATA19   # default=../DFT_DATA19
level.of.stdout                   1    # default=1 (1-3)
level.of.fileout                  1    # default=1 (0-2)
HS.fileout                   on       # on|off, default=off

#
# SCF or Electronic System
#

scf.XcType                  GGA-PBE    # LDA|LSDA-CA|LSDA-PW|GGA-PBE
scf.SpinPolarization        off        # On|Off|NC
scf.ElectronicTemperature   0      # default=300 (K)
scf.energycutoff           200.0       # default=150 (Ry)
scf.maxIter                 300         # default=40
scf.EigenvalueSolver        Band      # DC|GDC|Cluster|Band
scf.Kgrid                  6 6 6       # means 4x4x4
scf.Mixing.Type           rmm-diis     # Simple|Rmm-Diis|Gr-Pulay|Kerker|Rmm-Diisk
scf.Init.Mixing.Weight     0.10        # 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-7      # default=1.0e-6 (Hartree)

#
# MD or Geometry Optimization
#

MD.Type                      Nomd        # Nomd|Opt|NVE|NVT_VS|NVT_NH
                                       # Constraint_Opt|DIIS2|Constraint_DIIS2
MD.Opt.DIIS.History          4
MD.Opt.StartDIIS             5         # default=5
MD.maxIter                 100         # default=1
MD.TimeStep                1.0         # default=0.5 (fs)
MD.Opt.criterion          1.0e-4       # default=1.0e-4 (Hartree/bohr)

#
# MO output
#

MO.fileout                  off        # on|off, default=off
num.HOMOs                    2         # default=1
num.LUMOs                    2         # default=1

#
# 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
#
# Definition of Atomic Species
#
Species.Number       2
<Definition.of.Atomic.Species
S   S7.0-s2p2d1       S_PBE19
Mo   Mo7.0-s3p2d2       Mo_PBE19
Definition.of.Atomic.Species>

#
# Atoms
#
Atoms.Number          42
Atoms.SpeciesAndCoordinates.Unit   Ang # Ang|AU
<Atoms.SpeciesAndCoordinates           # Unit=Ang.
  1  Mo   0.0000000   0.0000000   7.4395018   7.00   7.00
  2  Mo   1.5951579   2.7628940   7.4395018   7.00   7.00
  3  Mo  -0.0000003   5.5257881   7.4395018   7.00   7.00
  4  Mo   3.1903157   5.5257881   7.4395018   7.00   7.00
  5  Mo   1.5951576   8.2886821   7.4395018   7.00   7.00
  6  Mo   4.7854736   8.2886821   7.4395018   7.00   7.00
  7  Mo   3.1903154  11.0515762   7.4395018   7.00   7.00
  8  Mo   0.0000000   0.0000000  13.7514776   7.00   7.00
  9  Mo   0.4557595   3.1575933  13.7514776   7.00   7.00
 10  Mo   3.4181958   4.3416910  13.7514776   7.00   7.00
 11  Mo   0.9115191   6.3151866  13.7514776   7.00   7.00
 12  Mo   3.8739553   7.4992843  13.7514776   7.00   7.00
 13  Mo   1.3672786   9.4727798  13.7514776   7.00   7.00
 14  Mo   4.3297149  10.6568776  13.7514776   7.00   7.00
 15  S   0.0000031   1.8419275   5.8746174   3.00   3.00
 16  S   1.5951610   4.6048216   5.8746174   3.00   3.00
 17  S   4.7854769   4.6048216   5.8746174   3.00   3.00
 18  S   0.0000028   7.3677156   5.8746174   3.00   3.00
 19  S   3.1903188   7.3677156   5.8746174   3.00   3.00
 20  S   1.5951607  10.1306097   5.8746174   3.00   3.00
 21  S   4.7854767  10.1306097   5.8746174   3.00   3.00
 22  S   0.0000031   1.8419275   9.0043862   3.00   3.00
 23  S   1.5951610   4.6048216   9.0043862   3.00   3.00
 24  S   4.7854769   4.6048216   9.0043862   3.00   3.00
 25  S   0.0000028   7.3677156   9.0043862   3.00   3.00
 26  S   3.1903188   7.3677156   9.0043862   3.00   3.00
 27  S   1.5951607  10.1306097   9.0043862   3.00   3.00
 28  S   4.7854767  10.1306097   9.0043862   3.00   3.00
 29  S   2.2788008   2.8944601  12.1865933   3.00   3.00
 30  S  -0.2278759   4.8679557  12.1865933   3.00   3.00
 31  S   2.7345604   6.0520534  12.1865933   3.00   3.00
 32  S   5.6969966   7.2361511  12.1865933   3.00   3.00
 33  S   0.2278837   8.0255490  12.1865933   3.00   3.00
 34  S   3.1903199   9.2096467  12.1865933   3.00   3.00
 35  S   3.6460794  12.3672400  12.1865933   3.00   3.00
 36  S   2.2788008   2.8944601  15.3163620   3.00   3.00
 37  S  -0.2278759   4.8679557  15.3163620   3.00   3.00
 38  S   2.7345604   6.0520534  15.3163620   3.00   3.00
 39  S   5.6969966   7.2361511  15.3163620   3.00   3.00
 40  S   0.2278837   8.0255490  15.3163620   3.00   3.00
 41  S   3.1903199   9.2096467  15.3163620   3.00   3.00
 42  S   3.6460794  12.3672400  15.3163620   3.00   3.00
Atoms.SpeciesAndCoordinates>
Atoms.UnitVectors.Unit             Ang #  Ang|AU
<Atoms.UnitVectors                     # unit=Ang.
       6.3806000   5.5258000   0.0000000
      -1.5952000   8.2887000   0.0000000
       0.0000000   0.0000000  50.0000000
Atoms.UnitVectors>

pbsnode shell script:

#! /bin/bash
#PBS -N MoS2_21_79
#PBS -e MoS2_21_79.err
#PBS -o MoS2_21_79.log
#PBS -V
#PBS -l nodes=1:ppn=16,walltime=9999:00:00

cd $PBS_O_WORKDIR
mpiexec -np 16 /home/zjlin/HamGNN/openmx_postprocess/openmx_postprocess MoS2_21_79_AA_1.dat > openmx.std

graph_data_gen.yaml:

ao_max: 19
graph_data_save_path: '/home/zjlin/ML_work/HamGNN/Bilayer_MoS2_Demo/database/work_dir/dataset/graph/MoS2_21_79/AA/'
read_openmx_path: '/home/zjlin/HamGNN/openmx_postprocess/read_openmx'
max_SCF_skip: 200
scfout_paths: '/home/zjlin/ML_work/HamGNN/Bilayer_MoS2_Demo/database/openmx/MoS2_21_79/AA' # The directory of the .scfout file calculated by openmx/openmx_postprocess, or a wildcard directory name to match multiple directories
dat_file_name: 'MoS2_21_79_AA_1.dat'
std_file_name: null  # None if no openmx computation is performed
scfout_file_name: 'overlap.scfout' # If the openmx self-consistent Hamiltonian is not required as the target, "overlap.scfout" can be used instead.
soc_switch: False # generate graph_data.npz for SOC (True) or Non-SOC (False) Hamiltonian

band_cal.yaml:

nao_max: 19
graph_data_path: '/home/zjlin/ML_work/HamGNN/Bilayer_MoS2_Demo/database/work_dir/dataset/graph/MoS2_21_79/AA/graph_data.npz'
hamiltonian_path: '/home/zjlin/ML_work/HamGNN/Bilayer_MoS2_Demo/database/bilayer_MoS2/Examples/Moire_twisted_bilayer_MoS2/version_0/target_hamiltonian.npy'
nk: 60          # the number of k points
save_dir: '/home/zjlin/ML_work/HamGNN/Bilayer_MoS2_Demo/database/work_dir/result/MoS2_21_79/AA/' # The directory to save the results
strcture_name: 'tw_MoS2_21_79_AA'  # The name of each cif file saved is strcture_name_idx.cif after band calculation
soc_switch: False
auto_mode: False # If the auto_mode is used, users can omit providing k_path and label, as the program will automatically generate them based on the crystal symmetry.
k_path: [[0.,0.,0.],[0.5,0.,0.0],[0.33333333,0.33333333,0.],[0.,0.,0.]]
label: ['$G$','$M$','$K$','$G$'] # The lable for each k points in K_path

The calcutation of overlay matrix clearly complete.

There're result: image by HamGNN

image by OpenMX

Is there any documentation or action of mine that contains errors?

Best regards, TzuChing

QuantumLab-ZY commented 8 months ago

Can you show me the parameters in the config.yaml file that are used for predicting the Hamiltonian of this structure? Did you use the .ckpt model weights file I provided in the repository to predict this band structure, or did you train your own .ckpt file?

QuantumLab-ZY commented 8 months ago

You can adjust the displayed range of the y-axis for the predicted bands to [-4, 4] and then compare them with the target bands. We only compare energy bands near the Fermi surface as predictions for very high empty states are unreliable; fortunately, we usually do not use these high energy bands.

newplay commented 8 months ago

Did you use the .ckpt model weights file I provided in the repository to predict this band structure, or did you train your own .ckpt file?

I used the .ckpt model from your dataset in the MoS2 demo, so I think the configuration is the same as yours. and if change the ylim to [-4,4]: image

QuantumLab-ZY commented 8 months ago

I noticed that your hamiltonian_path parameter in band_cal.yaml is '/home/zjlin/ML_work/HamGNN/Bilayer_MoS2_Demo/database/bilayer_MoS2/Examples/Moire_twisted_bilayer_MoS2/version_0/target_hamiltonian.npy',The target Hamiltonian here is fake because you haven't actually performed the SCF calculation, it is from the openmx_postprocess. Therefore, you should use "prediction_hamiltonian.npy" to calculate the energy bands.

newplay commented 8 months ago

I noticed that your hamiltonian_path parameter in band_cal.yaml is '/home/zjlin/ML_work/HamGNN/Bilayer_MoS2_Demo/database/bilayer_MoS2/Examples/Moire_twisted_bilayer_MoS2/version_0/target_hamiltonian.npy',The target Hamiltonian here is fake because you haven't actually performed the SCF calculation, it is from the openmx_postprocess. Therefore, you should use "prediction_hamiltonian.npy" to calculate the energy bands.

Thank you. I trust that's the essence of this issue. I will attempt it later.

Best regards, TzuChing

newplay commented 8 months ago

I noticed that your hamiltonian_path parameter in band_cal.yaml is '/home/zjlin/ML_work/HamGNN/Bilayer_MoS2_Demo/database/bilayer_MoS2/Examples/Moire_twisted_bilayer_MoS2/version_0/target_hamiltonian.npy',The target Hamiltonian here is fake because you haven't actually performed the SCF calculation, it is from the openmx_postprocess. Therefore, you should use "prediction_hamiltonian.npy" to calculate the energy bands.

Thank you. I trust that's the essence of this issue. I will attempt it later.

Best regards, TzuChing

I apologize for the unfortunate news: I changed the target_Hamiltonian.npy to prediction_Hamiltonian.npy, but the result remains the same as depicted in the figure above.

QuantumLab-ZY commented 8 months ago

Hi, TzuChing. I'm not sure if there is anything wrong with the graph_data.npz you got, but you can try the graph_data.npz I gave here for the bilayer MoS2 with 21.79 degrees. graph_data.npz.zip

newplay commented 8 months ago

Hi, TzuChing. I'm not sure if there is anything wrong with the graph_data.npz you got, but you can try the graph_data.npz I gave here for the bilayer MoS2 with 21.79 degrees. graph_data.npz.zip

Thank you , will try it immediately.

newplay commented 8 months ago

Dear Yang Zhong: May I kindly inquire if you have had the opportunity to execute the 'graph_data' file? I am still experiencing consistent results on my end. If convenient, could you please share the outcomes of your execution? Best regards, TzuChing

newplay commented 8 months ago

Dear Yang Zhong: Apologies for my oversight; I mistakenly executed 'band_cal' without utilizing .ckpt for prediction. I regret any time and effort that you may have been wasted. Moving forward, I will ensure to carefully review the requirements and inspect the source code before seeking assistance. Thank you for your support, and I apologize once again. Best regards, TzuChing

QuantumLab-ZY commented 8 months ago

Dear TzuChing,

No problem at all! These things happen, and I understand that troubleshooting can sometimes be a complex process. If you encounter any other issues, feel free to contact me.

Best wishes, Yang Zhong

Dear Yang Zhong: Apologies for my oversight; I mistakenly executed 'band_cal' without utilizing .ckpt for prediction. I regret any time and effort that you may have been wasted. Moving forward, I will ensure to carefully review the requirements and inspect the source code before seeking assistance. Thank you for your support, and I apologize once again. Best regards, TzuChing

newplay commented 8 months ago

Dear Yang Zhong,

I trust this message finds you well. Following the execution of 'HamGNN --config config.yaml,' the 'prediction_hamiltonian.npy' was successfully generated. However, I encountered an issue during the subsequent processing:

Traceback (most recent call last):
  File "/home/zjlin/HamGNN/utils_openmx/band_cal.py", line 438, in <module>
    main()
  File "/home/zjlin/HamGNN/utils_openmx/band_cal.py", line 365, in main
    w, v = eigh(a=HK[ik], b=SK[ik])
  File "/home/zjlin/anaconda3/envs/HamGNN/lib/python3.9/site-packages/scipy/linalg/_decomp.py", line 593, in eigh
    raise LinAlgError('The leading minor of order {} of B is not '
numpy.linalg.LinAlgError: The leading minor of order 557 of B is not positive definite. The factorization of B could not be completed, and no eigenvalues or eigenvectors were computed.

Upon investigation, I used the code snippet below to check the positive definiteness of matrices A and B:

A = np.linalg.eigvals(HK[ik])
B = np.linalg.eigvals(SK[ik])
if A.all() > 0 :
    print(f"A is a positive definite matrix")
elif B.all() > 0:
    print(f"B is a positive definite matrix")
else:
    print(f"A and B are not positive definite")

Surprisingly, the output consistently indicates 'B is a positive definite matrix.' I have also confirmed that the Hamiltonian matrix HK for the k-point is not Hermitian, as indicated by:

is_hermitian = np.allclose(HK[ik], HK[ik].conj().T)
if is_hermitian:
    print(f"For k-point {ik}: HK is a Hermitian matrix.")
else:
    print(f"For k-point {ik}: HK is NOT a Hermitian matrix.")

The output for this check is "For k-point 0: HK is NOT a Hermitian matrix."

I am seeking your guidance on resolving this issue. Additionally, I would appreciate clarification on whether the Hamiltonian matrix should be Hermitian in HamGNN. In quantum mechanics, it is generally understood that the Hamiltonian should be Hermitian.

Best regards, TzuChing

QuantumLab-ZY commented 8 months ago

Hi, TzuChing The first thing that comes to my mind is whether the graph_data.npz you used matches the target_hamiltonian.npz or prediction_hamiltonian.npz. If they don't match, it will result in non-Hermitian HK and SK constructions. I'm not sure if all the data you used is from the calculations for the same structure. You might be using another graph_data.npz file by mistake. Matrix 'B' here refers to the overlap matrix S, which we do not predict but directly obtain from openmx calculations. In extremely rare cases, the overlap matrix computed by openmx for certain structures may not be positive definite; I have encountered such situations before. Theoretically, HamGNN predicts Hamiltonian matrices that strictly satisfy Hermiticity. Before outputting the Hamiltonian matrix, HamGNN will perform symmetrization for Hamiltonian matrix. You can check if target_hamiltonian.npz has the Hermiticity problem. I am unsure of what minimum error np.allclose function allows when comparing two matrices; you can directly compare the differences between these two matrices, HK[ik]-HK[ik].conj().T, to see how large the error of Hermiticity is. Best wishes, Yang Zhong

Dear Yang Zhong,

I trust this message finds you well. Following the execution of 'HamGNN --config config.yaml,' the 'prediction_hamiltonian.npy' was successfully generated. However, I encountered an issue during the subsequent processing:

Traceback (most recent call last):
  File "/home/zjlin/HamGNN/utils_openmx/band_cal.py", line 438, in <module>
    main()
  File "/home/zjlin/HamGNN/utils_openmx/band_cal.py", line 365, in main
    w, v = eigh(a=HK[ik], b=SK[ik])
  File "/home/zjlin/anaconda3/envs/HamGNN/lib/python3.9/site-packages/scipy/linalg/_decomp.py", line 593, in eigh
    raise LinAlgError('The leading minor of order {} of B is not '
numpy.linalg.LinAlgError: The leading minor of order 557 of B is not positive definite. The factorization of B could not be completed, and no eigenvalues or eigenvectors were computed.

Upon investigation, I used the code snippet below to check the positive definiteness of matrices A and B:

A = np.linalg.eigvals(HK[ik])
B = np.linalg.eigvals(SK[ik])
if A.all() > 0 :
    print(f"A is a positive definite matrix")
elif B.all() > 0:
    print(f"B is a positive definite matrix")
else:
    print(f"A and B are not positive definite")

Surprisingly, the output consistently indicates 'B is a positive definite matrix.' I have also confirmed that the Hamiltonian matrix HK for the k-point is not Hermitian, as indicated by:

is_hermitian = np.allclose(HK[ik], HK[ik].conj().T)
if is_hermitian:
    print(f"For k-point {ik}: HK is a Hermitian matrix.")
else:
    print(f"For k-point {ik}: HK is NOT a Hermitian matrix.")

The output for this check is "For k-point 0: HK is NOT a Hermitian matrix."

I am seeking your guidance on resolving this issue. Additionally, I would appreciate clarification on whether the Hamiltonian matrix should be Hermitian in HamGNN. In quantum mechanics, it is generally understood that the Hamiltonian should be Hermitian.

Best regards, TzuChing

newplay commented 8 months ago

Hi, TzuChing The first thing that comes to my mind is whether the graph_data.npz you used matches the target_hamiltonian.npz or prediction_hamiltonian.npz. If they don't match, it will result in non-Hermitian HK and SK constructions. I'm not sure if all the data you used is from the calculations for the same structure. You might be using another graph_data.npz file by mistake. Matrix 'B' here refers to the overlap matrix S, which we do not predict but directly obtain from openmx calculations. In extremely rare cases, the overlap matrix computed by openmx for certain structures may not be positive definite; I have encountered such situations before. Theoretically, HamGNN predicts Hamiltonian matrices that strictly satisfy Hermiticity. Before outputting the Hamiltonian matrix, HamGNN will perform symmetrization for Hamiltonian matrix. You can check if target_hamiltonian.npz has the Hermiticity problem. I am unsure of what minimum error np.allclose function allows when comparing two matrices; you can directly compare the differences between these two matrices, HK[ik]-HK[ik].conj().T, to see how large the error of Hermiticity is. Best wishes, Yang Zhong

Dear Yang Zhong, I trust this message finds you well. Following the execution of 'HamGNN --config config.yaml,' the 'prediction_hamiltonian.npy' was successfully generated. However, I encountered an issue during the subsequent processing:

Traceback (most recent call last):
  File "/home/zjlin/HamGNN/utils_openmx/band_cal.py", line 438, in <module>
    main()
  File "/home/zjlin/HamGNN/utils_openmx/band_cal.py", line 365, in main
    w, v = eigh(a=HK[ik], b=SK[ik])
  File "/home/zjlin/anaconda3/envs/HamGNN/lib/python3.9/site-packages/scipy/linalg/_decomp.py", line 593, in eigh
    raise LinAlgError('The leading minor of order {} of B is not '
numpy.linalg.LinAlgError: The leading minor of order 557 of B is not positive definite. The factorization of B could not be completed, and no eigenvalues or eigenvectors were computed.

Upon investigation, I used the code snippet below to check the positive definiteness of matrices A and B:

A = np.linalg.eigvals(HK[ik])
B = np.linalg.eigvals(SK[ik])
if A.all() > 0 :
    print(f"A is a positive definite matrix")
elif B.all() > 0:
    print(f"B is a positive definite matrix")
else:
    print(f"A and B are not positive definite")

Surprisingly, the output consistently indicates 'B is a positive definite matrix.' I have also confirmed that the Hamiltonian matrix HK for the k-point is not Hermitian, as indicated by:

is_hermitian = np.allclose(HK[ik], HK[ik].conj().T)
if is_hermitian:
    print(f"For k-point {ik}: HK is a Hermitian matrix.")
else:
    print(f"For k-point {ik}: HK is NOT a Hermitian matrix.")

The output for this check is "For k-point 0: HK is NOT a Hermitian matrix." I am seeking your guidance on resolving this issue. Additionally, I would appreciate clarification on whether the Hamiltonian matrix should be Hermitian in HamGNN. In quantum mechanics, it is generally understood that the Hamiltonian should be Hermitian. Best regards, TzuChing

Indeed, you are correct! I employed the 'graph_data.npz' from your resources, but for the prediction process, I utilized my own approach. It is now functioning correctly! I sincerely appreciate your invaluable assistance and guidance throughout this process. Thank you very much!