-
Notifications
You must be signed in to change notification settings - Fork 16
/
Copy pathpost_process.f90
105 lines (91 loc) · 5.6 KB
/
post_process.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
#include "alias.inc"
subroutine post_process(PINPT, PPRAM, PPRAM_FIT, PKPTS, EDFT, PWGHT, PGEOM, NN_TABLE, PINPT_BERRY, PINPT_DOS, PRPLT)
use parameters, only: hopping, incar, energy, spmat, &
params, poscar, kpoints, weight, &
dos, berry, gainp, replot, unfold
use set_default, only: init_params
use mpi_setup
use reorder_band
use total_energy
#ifdef MKL_SPARSE
use sparse_tool, only: feast_initialize
#endif
implicit none
type(incar) :: PINPT
type(params ) :: PPRAM_FIT
type(params ), dimension(PINPT%nsystem) :: PPRAM
type(kpoints), dimension(PINPT%nsystem) :: PKPTS
type(energy ), dimension(PINPT%nsystem) :: EDFT
type(energy ), dimension(PINPT%nsystem) :: ETBA
type(weight ), dimension(PINPT%nsystem) :: PWGHT
type(poscar ), dimension(PINPT%nsystem) :: PGEOM
type(hopping), dimension(PINPT%nsystem) :: NN_TABLE
type(berry ), dimension(PINPT%nsystem) :: PINPT_BERRY
type(dos ), dimension(PINPT%nsystem) :: PINPT_DOS
type(replot ), dimension(PINPT%nsystem) :: PRPLT
type(unfold ), dimension(PINPT%nsystem) :: PUFLD
type(gainp ) :: PKAIA ! temp
integer*4 i, ik, my_ik
integer*4 mpierr
integer*8 buffer, buffer_limit, buffer2
integer*4, allocatable :: feast_ne(:,:), ne_found(:)
integer*4, allocatable :: ourjob(:), ourjob_disp(:)
real*8, allocatable :: E(:,:)
if(PINPT%flag_tbfit_finish) then
! this print_mode applies to the read_input routine, and make not to print input settings as it reads
! INCAR-TB again. This constraint make the output report redundunt
print_mode = 99
endif
buffer_limit = int8(2)**int8(31) - int8(1)
do i = 1, PINPT%nsystem
write(message,'( A)')' ' ; write_msg
write(message,'( A)')' #=======================================================' ; write_msg
if(PINPT%flag_tbfit_finish) then
write(message,'(2A)')' START POST-PROCESSING PROCEDURE: ',trim(PGEOM(i)%gfilenm) ; write_msg
else
write(message,'(1A)')' START POST-PROCESSING PROCEDURE: ' ; write_msg
endif
write(message,'( A)')' #=======================================================' ; write_msg
call read_input(PINPT,PPRAM(i),PKPTS(i), PGEOM(i), PWGHT(i), EDFT(i), NN_TABLE(i), &
PINPT_DOS(i), PINPT_BERRY(i), PKAIA, PRPLT(i), PUFLD(i), i)
if(PINPT%flag_tbfit_finish) then
call init_params(PPRAM(i), PINPT)
PPRAM(i) = PPRAM_FIT ! just use fitted parameter if flag_tbfit_finish
endif
if(PRPLT(i)%flag_replot) cycle
if(PINPT%flag_get_band .or. PINPT%flag_get_berry_curvature) then
if(.not. PINPT%flag_print_energy_singlek) then
call allocate_ETBA(PGEOM(i), PINPT, PKPTS(i), ETBA(i), .FALSE., 0)
call get_eig(NN_TABLE(i), PKPTS(i)%kpoint, PKPTS(i)%nkpoint, PINPT, PPRAM(i), ETBA(i)%E, ETBA(i)%V, ETBA(i)%SV, PGEOM(i)%neig, &
PGEOM(i)%init_erange, PGEOM(i)%nband, PINPT%flag_get_orbital, PINPT%flag_sparse, .true., PINPT%flag_phase)
elseif(PINPT%flag_print_energy_singlek) then ! called by PRTSEPK tag
call allocate_ETBA(PGEOM(i), PINPT, PKPTS(i), ETBA(i), .FALSE., 0)
call get_eig_sepk(NN_TABLE(i), PKPTS(i)%kpoint, PKPTS(i)%nkpoint, PINPT, PPRAM(i), ETBA(i)%E, PGEOM(i)%neig, &
PGEOM(i)%init_erange, PGEOM(i)%nband, PINPT%flag_get_orbital, PINPT%flag_sparse, .true., PINPT%flag_phase)
endif
call get_total_energy(ETBA(i), PINPT, PKPTS(i), PGEOM(i), .true., .true.)
if(PINPT%flag_get_band_order .or. PINPT%flag_get_band_order_print_only) then
call get_ordered_band(ETBA(i), PKPTS(i), PGEOM(i), PWGHT(i), PINPT, .false., PPRAM(i)%flag_use_overlap)
endif
call print_band_structure(PKPTS(i), ETBA(i), EDFT(i), PGEOM(i), PINPT, PPRAM(i), PWGHT(i))
! POST PROCESSING
if(PINPT%flag_get_berry_curvature) call get_berry_curvature(NN_TABLE(i), PINPT, PINPT_BERRY(i), PGEOM(i), PKPTS(i), ETBA(i))
if(PINPT%flag_get_circ_dichroism) call get_circular_dichroism(NN_TABLE(i), PINPT, PGEOM(i), PKPTS(i), ETBA(i))
if(PINPT%flag_plot_stm_image) call plot_stm_image(PINPT,PGEOM(i),PKPTS(i), ETBA(i), PPRAM(i)%flag_use_overlap)
if(PINPT%flag_plot_eigen_state) call plot_eigen_state(PINPT,PGEOM(i),PKPTS(i), ETBA(i), PPRAM(i)%flag_use_overlap)
if(PINPT%flag_get_effective_ham) call get_eig_downfold(PINPT,PPRAM(i),PKPTS(i),PGEOM(i),NN_TABLE(i)) ! NOTE: THIS IS EXPERIMENTAL, BUT WORKS ANYWAY.
endif
#ifdef MPI
call MPI_Barrier(mpi_comm_earth,mpierr)
#endif
! In the following routines, they call "get_eig" to get eigenvalues & eigenvectors in their routine
if(PINPT%flag_get_dos) call get_dos(NN_TABLE(i), PINPT, PPRAM(i), PINPT_DOS(i), PGEOM(i), PKPTS(i))
if(PINPT%flag_get_zak_phase) call get_zak_phase(NN_TABLE(i), PINPT, PPRAM(i), PINPT_BERRY(i), PGEOM(i), PKPTS(i))
if(PINPT%flag_get_wcc) call get_wcc(NN_TABLE(i), PINPT, PPRAM(i), PINPT_BERRY(i), PGEOM(i), PKPTS(i))
if(PINPT%flag_get_z2) call get_z2(NN_TABLE(i), PINPT, PPRAM(i), PINPT_BERRY(i), PGEOM(i), PKPTS(i))
if(PINPT%flag_get_parity) call get_parity(NN_TABLE(i), PINPT, PPRAM(i), PINPT_BERRY(i), PGEOM(i), PKPTS(i))
if(PINPT%flag_get_symmetry) call get_symmetry_eig(NN_TABLE(i), PINPT, PPRAM(i), PINPT_BERRY(i), PGEOM(i), PKPTS(i))
if(PINPT%flag_get_unfold) call get_unfold(PINPT, PPRAM(i), PUFLD(i))
enddo
return
endsubroutine