← Back to branch summary

~siesta-pseudos-bases/siesta/trunk-psml

« back to all changes in this revision

Viewing changes to Util/TS/TBtrans/m_tbt_save.F90

  • Committer: Alberto Garcia
  • Date: 2019-09-02 14:09:43 UTC
  • mfrom: (427.6.323 trunk)
  • Revision ID: albertog@icmab.es-20190902140943-mzmbe1jacgefpgxw
Sync to trunk-776 (notably nc/soc wavefunction support)


Show diffs side-by-side

added added

removed removed

Lines of Context:
Util/TS/TBtrans/m_tbt_save.F90
311
311
  end subroutine cdf_precision_cmplx
312
312
 
313
313
  subroutine init_cdf_save(fname,TSHS,r,btd,ispin, &
314
 
      N_Elec, Elecs, raEl, rElpd, btd_El, &
 
314
      N_Elec, Elecs, raEl, roElpd, btd_El, &
315
315
      nkpt, kpt, wkpt, NE, Eta, &
316
316
      a_Dev, a_Buf, sp_dev_sc, &
317
317
      save_DATA )
346
346
    integer, intent(in) :: ispin
347
347
    integer, intent(in) :: N_Elec
348
348
    type(Elec), intent(in) :: Elecs(N_Elec)
349
 
    type(tRgn), intent(in) :: raEl(N_Elec), rElpd(N_Elec), btd_El(N_Elec)
 
349
    type(tRgn), intent(in) :: raEl(N_Elec), roElpd(N_Elec), btd_El(N_Elec)
350
350
    integer, intent(in) :: nkpt
351
351
    real(dp), intent(in), target :: kpt(3,nkpt), wkpt(nkpt)
352
352
    integer, intent(in) :: NE
363
363
    type(hNCDF) :: ncdf, grp
364
364
    type(dict) :: dic
365
365
    logical :: exist, sme, isGamma
366
 
    integer :: iEl, jEl, i, nnzs_dev, N_eigen
 
366
    integer :: iEl, jEl, i, nnzs_dev, N_eigen, no_e
367
367
    integer :: prec_DOS, prec_T, prec_Teig, prec_J, prec_COOP, prec_DM
368
368
    type(OrbitalDistribution) :: fdit
369
369
    real(dp) :: mem
372
372
    type(tRgn) :: a_Dev_sort, r_tmp
373
373
#ifdef TBT_PHONON
374
374
    character(len=*), parameter :: T_unit = 'g0'
 
375
    character(len=*), parameter :: COHP_unit = 'Ry'
375
376
#else
376
377
    character(len=*), parameter :: T_unit = 'G0'
 
378
    character(len=*), parameter :: COHP_unit = 'Ry/Ry'
377
379
#endif
378
380
#ifdef MPI
379
381
    integer :: MPIerror
603
605
       mem = mem + calc_mem(NF90_INT, a_Buf%n)
604
606
    end if
605
607
 
606
 
    if ( 'DOS-Gf' .in. save_DATA ) then
607
 
       dic = dic // ('info'.kv.'Density of states')//('unit'.kv.'1/Ry')
608
 
       call ncdf_def_var(ncdf,'DOS',prec_DOS,(/'no_d','ne  ','nkpt'/), &
609
 
            atts = dic , chunks = (/r%n,1,1/) , compress_lvl = cmp_lvl )
610
 
       mem = mem + calc_mem(prec_DOS, r%n, NE, nkpt)
611
 
    end if
612
 
 
613
608
    dic = dic // ('info'.kv.'k point')//('unit'.kv.'b')
614
609
    call ncdf_def_var(ncdf,'kpt',NF90_DOUBLE,(/'xyz ','nkpt'/), &
615
610
         atts = dic)
633
628
#endif
634
629
    call ncdf_def_var(ncdf,'eta',NF90_DOUBLE,(/'one'/), atts = dic)
635
630
 
 
631
    if ( 'DOS-Gf' .in. save_DATA ) then
 
632
       dic = dic // ('info'.kv.'Density of states')//('unit'.kv.'1/Ry')
 
633
       call ncdf_def_var(ncdf,'DOS',prec_DOS,(/'no_d','ne  ','nkpt'/), &
 
634
            atts = dic , chunks = (/r%n,1,1/) , compress_lvl = cmp_lvl )
 
635
       mem = mem + calc_mem(prec_DOS, r%n, NE, nkpt)
 
636
    end if
 
637
 
636
638
    ! Clean-up dictionary
637
639
    call delete(dic)
638
640
 
705
707
 
706
708
    end if
707
709
 
 
710
    dic = dic // ('unit'.kv.'1/Ry')
708
711
    if ( 'DM-Gf' .in. save_DATA ) then
709
 
       dic = dic // ('info'.kv.'Green function density matrix')//('unit'.kv.'1/Ry')
 
712
       dic = dic // ('info'.kv.'Green function density matrix')
710
713
       call ncdf_def_var(ncdf,'DM',prec_DM,(/'nnzs','ne  ','nkpt'/), &
711
714
           atts = dic , chunks = (/nnzs_dev/) , compress_lvl=cmp_lvl)
712
715
       mem = mem + calc_mem(prec_DM, nnzs_dev, NE, nkpt)
713
716
    end if
714
717
    if ( 'COOP-Gf' .in. save_DATA ) then
715
 
       dic = dic // ('info'.kv.'Crystal orbital overlap population')//('unit'.kv.'1/Ry')
 
718
       dic = dic // ('info'.kv.'Crystal orbital overlap population')
716
719
       call ncdf_def_var(ncdf,'COOP',prec_COOP,(/'nnzs','ne  ','nkpt'/), &
717
720
            atts = dic , chunks = (/nnzs_dev/) , compress_lvl=cmp_lvl)
718
721
       mem = mem + calc_mem(prec_COOP, nnzs_dev, NE, nkpt)
719
722
    end if
720
723
    if ( 'COHP-Gf' .in. save_DATA ) then
721
 
       dic = dic // ('info'.kv.'Crystal orbital Hamilton population')//('unit'.kv.'Ry/Ry')
 
724
       dic = dic // ('info'.kv.'Crystal orbital Hamilton population')//('unit'.kv.COHP_unit)
722
725
       call ncdf_def_var(ncdf,'COHP',prec_COOP,(/'nnzs','ne  ','nkpt'/), &
723
726
            atts = dic , chunks = (/nnzs_dev/) , compress_lvl=cmp_lvl)
724
727
       mem = mem + calc_mem(prec_COOP, nnzs_dev, NE, nkpt)
748
751
       call ncdf_put_var(grp,'a_down',raEl(iEl)%r)
749
752
       mem = mem + calc_mem(NF90_INT, raEl(iEl)%n)
750
753
 
751
 
       call ncdf_def_dim(grp,'n_btd',btd_El(iEl)%n)
752
 
       call ncdf_def_dim(grp,'no_down',rElpd(iEl)%n)
753
 
 
754
754
       ! Save generic information about electrode
755
755
       dic = dic//('info'.kv.'Bloch expansion')
756
756
       call ncdf_def_var(grp,'bloch',NF90_INT,(/'xyz'/), atts = dic)
757
 
       call ncdf_put_var(grp,'bloch',Elecs(iEl)%Bloch)
 
757
       call ncdf_put_var(grp,'bloch',Elecs(iEl)%Bloch%B)
 
758
       mem = mem + calc_mem(NF90_INT, 3)
 
759
 
 
760
       call ncdf_def_dim(grp,'no_down',roElpd(iEl)%n)
758
761
 
759
762
       dic = dic//('info'.kv.'Downfolding region orbital pivot table')
760
 
       call ncdf_def_var(grp,'pivot',NF90_INT,(/'no_down'/), atts = dic)
761
 
       call ncdf_put_var(grp,'pivot',rElpd(iEl)%r)
762
 
       mem = mem + calc_mem(NF90_INT, rElpd(iEl)%n)
 
763
       call ncdf_def_var(grp,'pivot_down',NF90_INT,(/'no_down'/), atts = dic)
 
764
       call ncdf_put_var(grp,'pivot_down',roElpd(iEl)%r)
 
765
       mem = mem + calc_mem(NF90_INT, roElpd(iEl)%n)
 
766
 
 
767
       call ncdf_def_dim(grp,'n_btd',btd_El(iEl)%n)
763
768
 
764
769
       dic = dic//('info'.kv.'Blocks in BTD downfolding for the pivot table')
765
770
       call ncdf_def_var(grp,'btd',NF90_INT,(/'n_btd'/), atts = dic)
766
771
       call ncdf_put_var(grp,'btd',btd_El(iEl)%r)
767
772
       mem = mem + calc_mem(NF90_INT, btd_El(iEl)%n)
768
773
 
 
774
       no_e = Elecs(iEl)%o_inD%n
 
775
       call ncdf_def_dim(grp,'no_e',no_e)
 
776
 
 
777
       dic = dic//('info'.kv.'Orbital pivot table for self-energy')
 
778
       call ncdf_def_var(grp,'pivot',NF90_INT,(/'no_e'/), atts = dic)
 
779
       call ncdf_put_var(grp,'pivot',Elecs(iEl)%o_inD%r)
 
780
       mem = mem + calc_mem(NF90_INT, no_e)
 
781
 
769
782
       dic = dic//('info'.kv.'Chemical potential')//('unit'.kv.'Ry')
770
783
       call ncdf_def_var(grp,'mu',NF90_DOUBLE,(/'one'/), atts = dic)
771
784
       call ncdf_put_var(grp,'mu',Elecs(iEl)%mu%mu)
785
798
       call ncdf_def_var(grp,'eta',NF90_DOUBLE,(/'one'/), atts = dic)
786
799
       call ncdf_put_var(grp,'eta',Elecs(iEl)%Eta)
787
800
 
 
801
       dic = dic//('info'.kv.'Accuracy of the self-energy')//('unit'.kv.'Ry')
 
802
      call ncdf_def_var(grp,'Accuracy',NF90_DOUBLE,(/'one'/), atts = dic)
 
803
      call ncdf_put_var(grp,'Accuracy',Elecs(iEl)%accu)
 
804
      call delete(dic)
 
805
 
 
806
       mem = mem + calc_mem(NF90_DOUBLE, 4) ! mu, kT, eta, Accuracy
 
807
 
788
808
       call delete(dic)
789
809
 
790
810
       if ( ('DOS-Elecs' .in. save_DATA) .and. .not. Elecs(iEl)%out_of_core ) then
805
825
          dic = dic//('info'.kv.'Unit cell')//('unit'.kv.'Bohr')
806
826
          call ncdf_def_var(grp,'cell',NF90_DOUBLE,(/'xyz','xyz'/), &
807
827
               atts = dic)
 
828
          mem = mem + calc_mem(NF90_DOUBLE, 3, 3)
808
829
          
809
830
          dic = dic//('info'.kv.'Atomic coordinates')
810
831
          call ncdf_def_var(grp,'xa',NF90_DOUBLE,(/'xyz ','na_u'/), &
827
848
       
828
849
       ! Now we will only add information that is calculated
829
850
       if ( iEl == N_Elec ) then
830
 
          ! check if all are calculated
831
 
          if ( ('DOS-A-all' .nin. save_DATA) .and. &
832
 
               ('T-all'.nin. save_DATA) ) cycle
 
851
         ! check if all are calculated
 
852
         if ( ('DOS-A-all' .nin. save_DATA) .and. &
 
853
             ('T-all'.nin. save_DATA) ) cycle
833
854
       end if
834
855
 
 
856
       ! All have this unit
 
857
       dic = ('unit'.kv.'1/Ry')
 
858
 
835
859
       if ( 'DM-A' .in. save_DATA ) then
836
 
          dic = dic//('info'.kv.'Spectral function density matrix')//('unit'.kv.'1/Ry')
 
860
          dic = dic//('info'.kv.'Spectral function density matrix')
837
861
          call ncdf_def_var(grp,'DM',prec_DM,(/'nnzs','ne  ','nkpt'/), &
838
862
              atts = dic , chunks = (/nnzs_dev/) , compress_lvl=cmp_lvl)
839
863
          mem = mem + calc_mem(prec_DM, nnzs_dev, NE, nkpt)
840
864
       end if
841
865
 
842
866
       if ( 'DOS-A' .in. save_DATA ) then
843
 
          dic = dic//('info'.kv.'Spectral function density of states')// &
844
 
               ('unit'.kv.'1/Ry')
 
867
          dic = dic//('info'.kv.'Spectral function density of states')
845
868
          call ncdf_def_var(grp,'ADOS',prec_DOS,(/'no_d','ne  ','nkpt'/), &
846
869
               atts = dic, chunks = (/r%n,1,1/) , compress_lvl=cmp_lvl)
847
870
          mem = mem + calc_mem(prec_DOS, r%n, NE, nkpt)
848
871
       end if
849
872
 
850
873
       if ( 'COOP-A' .in. save_DATA ) then
851
 
          dic = dic//('info'.kv.'Crystal orbital overlap population')//('unit'.kv.'1/Ry')
 
874
          dic = dic//('info'.kv.'Crystal orbital overlap population')
852
875
          call ncdf_def_var(grp,'COOP',prec_COOP,(/'nnzs','ne  ','nkpt'/), &
853
876
               atts = dic , chunks = (/nnzs_dev/) , compress_lvl=cmp_lvl)
854
877
          mem = mem + calc_mem(prec_COOP, nnzs_dev, NE, nkpt)
855
878
       end if
856
879
 
857
880
       if ( 'COHP-A' .in. save_DATA ) then
858
 
          dic = dic//('info'.kv.'Crystal orbital Hamilton population')//('unit'.kv.'Ry/Ry')
 
881
          dic = dic//('info'.kv.'Crystal orbital Hamilton population')//('unit'.kv.COHP_unit)
859
882
          call ncdf_def_var(grp,'COHP',prec_COOP,(/'nnzs','ne  ','nkpt'/), &
860
883
               atts = dic , chunks = (/nnzs_dev/) , compress_lvl=cmp_lvl)
861
884
          mem = mem + calc_mem(prec_COOP, nnzs_dev, NE, nkpt)
862
885
       end if
863
886
 
 
887
       ! All quantities here are transmissions.
 
888
       dic = dic//('unit'.kv.T_unit)
 
889
 
864
890
       if ( 'orb-current' .in. save_DATA ) then
865
 
          dic = dic//('info'.kv.'Orbital transmission')//('unit'.kv.T_unit)
 
891
          dic = dic//('info'.kv.'Orbital transmission')
866
892
          call ncdf_def_var(grp,'J',prec_J,(/'nnzs','ne  ','nkpt'/), &
867
893
               atts = dic , chunks = (/nnzs_dev/) , compress_lvl=cmp_lvl)
868
894
          mem = mem + calc_mem(prec_J, nnzs_dev, NE, nkpt)
869
895
       end if
870
896
 
871
 
       ! All quantities here are transmissions.
872
 
       dic = dic//('unit'.kv.T_unit)
873
 
 
874
897
       tmp = trim(Elecs(iEl)%name)
875
898
       do jEl = 1 , N_Elec
876
899
          if ( ('T-all' .nin. save_DATA ) .and. &
1990
2013
      call ncdf_get_var(grp, var, r3)
1991
2014
 
1992
2015
      ! Immediately sum all orbitals and convert to 1/eV
1993
 
!$OMP parallel do default(shared), collapse(2), &
1994
 
!$OMP&  private(ik,ie,io,DOS)
 
2016
!$OMP parallel do default(shared), private(ik,ie,io,DOS)
1995
2017
      do ik = 1, nkpt
1996
2018
        do ie = 1, NE
1997
2019
          DOS = 0._dp