1!-*- mode: F90 -*-! 2!------------------------------------------------------------! 3! This file is distributed as part of the Wannier90 code and ! 4! under the terms of the GNU General Public License. See the ! 5! file `LICENSE' in the root directory of the Wannier90 ! 6! distribution, or http://www.gnu.org/copyleft/gpl.txt ! 7! ! 8! The webpage of the Wannier90 code is www.wannier.org ! 9! ! 10! The Wannier90 code is hosted on GitHub: ! 11! ! 12! https://github.com/wannier-developers/wannier90 ! 13!------------------------------------------------------------! 14 15program postw90 16 !! The postw90 program 17 use w90_constants, only: dp, eps6 18 use w90_parameters 19 use w90_io 20 21 use w90_kmesh 22 use w90_comms, only: on_root, num_nodes, comms_setup, comms_end, comms_bcast, comms_barrier 23 use w90_postw90_common, only: pw90common_wanint_setup, pw90common_wanint_get_kpoint_file, & 24 pw90common_wanint_param_dist, pw90common_wanint_data_dist 25 26 ! These modules deal with the interpolation of specific physical properties 27 ! 28 use w90_dos 29 use w90_berry 30 use w90_gyrotropic 31 use w90_spin 32 use w90_kpath 33 use w90_kslice 34 35 use w90_boltzwann 36 use w90_geninterp 37 38 implicit none 39 40 integer :: nkp, len_seedname 41 logical :: have_gamma 42 real(kind=dp) :: time0, time1, time2 43 character(len=9) :: stat, pos 44 logical :: wpout_found, werr_found, dryrun 45 character(len=50) :: prog 46 47 ! Put some descriptive comments here 48 ! 49 call comms_setup 50 51 library = .false. 52 ispostw90 = .true. 53 54 if (on_root) then 55 time0 = io_time() 56 prog = 'postw90' 57 call io_commandline(prog, dryrun) 58 len_seedname = len(seedname) 59 end if 60 call comms_bcast(len_seedname, 1) 61 call comms_bcast(seedname, len_seedname) 62 call comms_bcast(dryrun, 1) 63 64 if (on_root) then 65 ! If an error file (generated by postw90) exists, I delete it 66 ! Note: I do it only for the error file generated from the root node; 67 ! If error files generated by other nodes exist, I don't do anything for them 68 inquire (file=trim(seedname)//'.node_00000.werr', exist=werr_found) 69 if (werr_found) then 70 stdout = io_file_unit() 71 open (unit=stdout, file=trim(seedname)//'.node_00000.werr', status='old', position='append') 72 close (stdout, status='delete') 73 end if 74 75 inquire (file=trim(seedname)//'.wpout', exist=wpout_found) 76 if (wpout_found) then 77 stat = 'old' 78 else 79 stat = 'replace' 80 endif 81 pos = 'append' 82 83 stdout = io_file_unit() 84 open (unit=stdout, file=trim(seedname)//'.wpout', status=trim(stat), position=trim(pos)) 85 86 call param_write_header 87 if (num_nodes == 1) then 88#ifdef MPI 89 write (stdout, '(/,1x,a)') 'Running in serial (with parallel executable)' 90#else 91 write (stdout, '(/,1x,a)') 'Running in serial (with serial executable)' 92#endif 93 else 94 write (stdout, '(/,1x,a,i3,a/)') & 95 'Running in parallel on ', num_nodes, ' CPUs' 96 endif 97 end if 98 99 ! Read onto the root node all the input parameters from seendame.win, 100 ! as well as the energy eigenvalues on the ab-initio q-mesh from seedname.eig 101 ! 102 if (on_root) then 103 call param_read 104 call param_postw90_write 105 time1 = io_time() 106 write (stdout, '(1x,a25,f11.3,a)') & 107 'Time to read parameters ', time1 - time0, ' (sec)' 108 109 if (.not. effective_model) then 110 ! Check if the q-mesh includes the gamma point 111 ! 112 have_gamma = .false. 113 do nkp = 1, num_kpts 114 if (all(abs(kpt_latt(:, nkp)) < eps6)) have_gamma = .true. 115 end do 116 if (.not. have_gamma) write (stdout, '(1x,a)') & 117 'Ab-initio does not include Gamma. Interpolation may be incorrect!!!' 118 ! 119 ! Need nntot, wb, and bk to evaluate WF matrix elements of 120 ! the position operator in reciprocal space. Also need 121 ! nnlist to compute the additional matrix elements entering 122 ! the orbital magnetization 123 ! 124 call kmesh_get 125 time2 = io_time() 126 write (stdout, '(1x,a25,f11.3,a)') & 127 'Time to get kmesh ', time2 - time1, ' (sec)' 128 endif 129 130 ! GP, May 10, 2012: for the moment I leave this commented 131 ! since we need first to tune that routine so that it doesn't 132 ! print the memory information related to wannier90.x. 133 ! Note that the code for the memory estimation for the 134 ! Boltzwann routine is already there. 135 ! call param_memory_estimate 136 end if 137 138 if (dryrun) then 139 if (on_root) then 140 write (stdout, *) ' ' 141 write (stdout, *) ' ===============================' 142 write (stdout, *) ' DRYRUN ' 143 write (stdout, *) ' No problems found with win file' 144 write (stdout, *) ' ===============================' 145 endif 146 stop 147 endif 148 149 ! We now distribute a subset of the parameters to the other nodes 150 ! 151 call pw90common_wanint_param_dist 152 153 if (.not. effective_model) then 154 ! 155 ! Read files seedname.chk (overlap matrices, unitary matrices for 156 ! both disentanglement and maximal localization, etc.) 157 ! 158 if (on_root) call param_read_chkpt() 159 ! 160 ! Distribute the information in the um and chk files to the other nodes 161 ! 162 ! Ivo: For interpolation purposes we do not need u_matrix_opt and 163 ! u_matrix separately, only their product v_matrix, and this 164 ! is what is distributed now 165 ! 166 call pw90common_wanint_data_dist 167 ! 168 end if 169 170 ! Read list of k-points in irreducible BZ and their weights 171 ! 172 ! Should this be done on root node only? 173 ! 174 if (wanint_kpoint_file) call pw90common_wanint_get_kpoint_file 175 176 ! Setup a number of common variables for all interpolation tasks 177 ! 178 call pw90common_wanint_setup 179 180 if (on_root) then 181 time1 = io_time() 182 write (stdout, '(/1x,a25,f11.3,a)') & 183 'Time to read and process .chk ', time1 - time2, ' (sec)' 184 endif 185 ! 186 ! Now perform one or more of the following tasks 187 188 ! --------------------------------------------------------------- 189 ! Density of states calculated using a uniform interpolation mesh 190 ! --------------------------------------------------------------- 191 ! 192 if (dos .and. index(dos_task, 'dos_plot') > 0) call dos_main 193 194! find_fermi_level commented for the moment in dos.F90 195! if(dos .and. index(dos_task,'find_fermi_energy')>0) call find_fermi_level 196 197 ! -------------------------------------------------------------------- 198 ! Bands, Berry curvature, or orbital magnetization plot along a k-path 199 ! -------------------------------------------------------------------- 200 ! 201 if (kpath) call k_path 202 203 ! --------------------------------------------------------------------------- 204 ! Bands, Berry curvature, or orbital magnetization plot on a slice in k-space 205 ! --------------------------------------------------------------------------- 206 ! 207 if (kslice) call k_slice 208 209 ! -------------------- 210 ! Spin magnetic moment 211 ! -------------------- 212 ! 213 if (spin_moment) call spin_get_moment 214 215 ! ------------------------------------------------------------------- 216 ! dc Anomalous Hall conductivity and eventually (if 'mcd' string also 217 ! present in addition to 'ahe', e.g., 'ahe+mcd') dichroic optical 218 ! conductivity, both calculated on the same (adaptively-refined) mesh 219 ! ------------------------------------------------------------------- 220 ! 221 ! --------------------------------------------------------------- 222 ! Absorptive dichroic optical conductivity & JDOS on uniform mesh 223 ! --------------------------------------------------------------- 224 ! 225 ! ----------------------------------------------------------------- 226 ! Absorptive ordinary optical conductivity & JDOS on a uniform mesh 227 ! ----------------------------------------------------------------- 228 ! 229 ! ----------------------------------------------------------------- 230 ! Orbital magnetization 231 ! ----------------------------------------------------------------- 232 ! 233 if (berry) call berry_main 234 ! ----------------------------------------------------------------- 235 ! Boltzmann transport coefficients (BoltzWann module) 236 ! ----------------------------------------------------------------- 237 ! 238 if (on_root) then 239 time1 = io_time() 240 endif 241 242 if (geninterp) call geninterp_main 243 244 if (boltzwann) call boltzwann_main 245 246 if (gyrotropic) call gyrotropic_main 247 248 if (on_root .and. boltzwann) then 249 time2 = io_time() 250 write (stdout, '(/1x,a,f11.3,a)') & 251 'Time for BoltzWann (Boltzmann transport) ', time2 - time1, ' (sec)' 252 endif 253 254 ! I put a barrier here before calling the final time printing routines, 255 ! just to be sure that all processes have arrived here. 256 call comms_barrier 257 if (on_root) then 258 write (stdout, '(/,1x,a25,f11.3,a)') & 259 'Total Execution Time ', io_time(), ' (sec)' 260 if (timing_level > 0) call io_print_timings() 261 write (stdout, *) 262 write (stdout, '(/,1x,a)') 'All done: postw90 exiting' 263 close (stdout) 264 end if 265 266 call comms_end 267 268end program postw90 269 270