!-*- mode: F90 -*-! !------------------------------------------------------------! ! This file is distributed as part of the Wannier90 code and ! ! under the terms of the GNU General Public License. See the ! ! file `LICENSE' in the root directory of the Wannier90 ! ! distribution, or http://www.gnu.org/copyleft/gpl.txt ! ! ! ! The webpage of the Wannier90 code is www.wannier.org ! ! ! ! The Wannier90 code is hosted on GitHub: ! ! ! ! https://github.com/wannier-developers/wannier90 ! !------------------------------------------------------------! program postw90 !! The postw90 program use w90_constants, only: dp, eps6 use w90_parameters use w90_io use w90_kmesh use w90_comms, only: on_root, num_nodes, comms_setup, comms_end, comms_bcast, comms_barrier use w90_postw90_common, only: pw90common_wanint_setup, pw90common_wanint_get_kpoint_file, & pw90common_wanint_param_dist, pw90common_wanint_data_dist ! These modules deal with the interpolation of specific physical properties ! use w90_dos use w90_berry use w90_gyrotropic use w90_spin use w90_kpath use w90_kslice use w90_boltzwann use w90_geninterp implicit none integer :: nkp, len_seedname logical :: have_gamma real(kind=dp) :: time0, time1, time2 character(len=9) :: stat, pos logical :: wpout_found, werr_found, dryrun character(len=50) :: prog ! Put some descriptive comments here ! call comms_setup library = .false. ispostw90 = .true. if (on_root) then time0 = io_time() prog = 'postw90' call io_commandline(prog, dryrun) len_seedname = len(seedname) end if call comms_bcast(len_seedname, 1) call comms_bcast(seedname, len_seedname) call comms_bcast(dryrun, 1) if (on_root) then ! If an error file (generated by postw90) exists, I delete it ! Note: I do it only for the error file generated from the root node; ! If error files generated by other nodes exist, I don't do anything for them inquire (file=trim(seedname)//'.node_00000.werr', exist=werr_found) if (werr_found) then stdout = io_file_unit() open (unit=stdout, file=trim(seedname)//'.node_00000.werr', status='old', position='append') close (stdout, status='delete') end if inquire (file=trim(seedname)//'.wpout', exist=wpout_found) if (wpout_found) then stat = 'old' else stat = 'replace' endif pos = 'append' stdout = io_file_unit() open (unit=stdout, file=trim(seedname)//'.wpout', status=trim(stat), position=trim(pos)) call param_write_header if (num_nodes == 1) then #ifdef MPI write (stdout, '(/,1x,a)') 'Running in serial (with parallel executable)' #else write (stdout, '(/,1x,a)') 'Running in serial (with serial executable)' #endif else write (stdout, '(/,1x,a,i3,a/)') & 'Running in parallel on ', num_nodes, ' CPUs' endif end if ! Read onto the root node all the input parameters from seendame.win, ! as well as the energy eigenvalues on the ab-initio q-mesh from seedname.eig ! if (on_root) then call param_read call param_postw90_write time1 = io_time() write (stdout, '(1x,a25,f11.3,a)') & 'Time to read parameters ', time1 - time0, ' (sec)' if (.not. effective_model) then ! Check if the q-mesh includes the gamma point ! have_gamma = .false. do nkp = 1, num_kpts if (all(abs(kpt_latt(:, nkp)) < eps6)) have_gamma = .true. end do if (.not. have_gamma) write (stdout, '(1x,a)') & 'Ab-initio does not include Gamma. Interpolation may be incorrect!!!' ! ! Need nntot, wb, and bk to evaluate WF matrix elements of ! the position operator in reciprocal space. Also need ! nnlist to compute the additional matrix elements entering ! the orbital magnetization ! call kmesh_get time2 = io_time() write (stdout, '(1x,a25,f11.3,a)') & 'Time to get kmesh ', time2 - time1, ' (sec)' endif ! GP, May 10, 2012: for the moment I leave this commented ! since we need first to tune that routine so that it doesn't ! print the memory information related to wannier90.x. ! Note that the code for the memory estimation for the ! Boltzwann routine is already there. ! call param_memory_estimate end if if (dryrun) then if (on_root) then write (stdout, *) ' ' write (stdout, *) ' ===============================' write (stdout, *) ' DRYRUN ' write (stdout, *) ' No problems found with win file' write (stdout, *) ' ===============================' endif stop endif ! We now distribute a subset of the parameters to the other nodes ! call pw90common_wanint_param_dist if (.not. effective_model) then ! ! Read files seedname.chk (overlap matrices, unitary matrices for ! both disentanglement and maximal localization, etc.) ! if (on_root) call param_read_chkpt() ! ! Distribute the information in the um and chk files to the other nodes ! ! Ivo: For interpolation purposes we do not need u_matrix_opt and ! u_matrix separately, only their product v_matrix, and this ! is what is distributed now ! call pw90common_wanint_data_dist ! end if ! Read list of k-points in irreducible BZ and their weights ! ! Should this be done on root node only? ! if (wanint_kpoint_file) call pw90common_wanint_get_kpoint_file ! Setup a number of common variables for all interpolation tasks ! call pw90common_wanint_setup if (on_root) then time1 = io_time() write (stdout, '(/1x,a25,f11.3,a)') & 'Time to read and process .chk ', time1 - time2, ' (sec)' endif ! ! Now perform one or more of the following tasks ! --------------------------------------------------------------- ! Density of states calculated using a uniform interpolation mesh ! --------------------------------------------------------------- ! if (dos .and. index(dos_task, 'dos_plot') > 0) call dos_main ! find_fermi_level commented for the moment in dos.F90 ! if(dos .and. index(dos_task,'find_fermi_energy')>0) call find_fermi_level ! -------------------------------------------------------------------- ! Bands, Berry curvature, or orbital magnetization plot along a k-path ! -------------------------------------------------------------------- ! if (kpath) call k_path ! --------------------------------------------------------------------------- ! Bands, Berry curvature, or orbital magnetization plot on a slice in k-space ! --------------------------------------------------------------------------- ! if (kslice) call k_slice ! -------------------- ! Spin magnetic moment ! -------------------- ! if (spin_moment) call spin_get_moment ! ------------------------------------------------------------------- ! dc Anomalous Hall conductivity and eventually (if 'mcd' string also ! present in addition to 'ahe', e.g., 'ahe+mcd') dichroic optical ! conductivity, both calculated on the same (adaptively-refined) mesh ! ------------------------------------------------------------------- ! ! --------------------------------------------------------------- ! Absorptive dichroic optical conductivity & JDOS on uniform mesh ! --------------------------------------------------------------- ! ! ----------------------------------------------------------------- ! Absorptive ordinary optical conductivity & JDOS on a uniform mesh ! ----------------------------------------------------------------- ! ! ----------------------------------------------------------------- ! Orbital magnetization ! ----------------------------------------------------------------- ! if (berry) call berry_main ! ----------------------------------------------------------------- ! Boltzmann transport coefficients (BoltzWann module) ! ----------------------------------------------------------------- ! if (on_root) then time1 = io_time() endif if (geninterp) call geninterp_main if (boltzwann) call boltzwann_main if (gyrotropic) call gyrotropic_main if (on_root .and. boltzwann) then time2 = io_time() write (stdout, '(/1x,a,f11.3,a)') & 'Time for BoltzWann (Boltzmann transport) ', time2 - time1, ' (sec)' endif ! I put a barrier here before calling the final time printing routines, ! just to be sure that all processes have arrived here. call comms_barrier if (on_root) then write (stdout, '(/,1x,a25,f11.3,a)') & 'Total Execution Time ', io_time(), ' (sec)' if (timing_level > 0) call io_print_timings() write (stdout, *) write (stdout, '(/,1x,a)') 'All done: postw90 exiting' close (stdout) end if call comms_end end program postw90