wann_write_xyz Subroutine

private subroutine wann_write_xyz()

Uses

  • proc~~wann_write_xyz~~UsesGraph proc~wann_write_xyz wann_write_xyz module~w90_utility w90_utility proc~wann_write_xyz->module~w90_utility module~w90_parameters w90_parameters proc~wann_write_xyz->module~w90_parameters module~w90_io w90_io proc~wann_write_xyz->module~w90_io module~w90_constants w90_constants module~w90_utility->module~w90_constants module~w90_parameters->module~w90_io module~w90_parameters->module~w90_constants module~w90_io->module~w90_constants

Arguments

None

Calls

proc~~wann_write_xyz~~CallsGraph proc~wann_write_xyz wann_write_xyz proc~io_date io_date proc~wann_write_xyz->proc~io_date proc~io_file_unit io_file_unit proc~wann_write_xyz->proc~io_file_unit

Contents

Source Code


Source Code

  subroutine wann_write_xyz()
    !=====================================!
    !                                     !
    ! Write xyz file with Wannier centres !
    !                                     !
    !=====================================!

    use w90_io, only: seedname, io_file_unit, io_date, stdout
    use w90_parameters, only: translate_home_cell, num_wann, wannier_centres, &
      lenconfac, real_lattice, recip_lattice, iprint, &
      num_atoms, atoms_symbol, atoms_pos_cart, &
      num_species, atoms_species_num
    use w90_utility, only: utility_translate_home

    implicit none

    integer          :: iw, ind, xyz_unit, nsp, nat
    character(len=9) :: cdate, ctime
    real(kind=dp)    :: wc(3, num_wann)

    wc = wannier_centres

    if (translate_home_cell) then
      do iw = 1, num_wann
        call utility_translate_home(wc(:, iw), real_lattice, recip_lattice)
      enddo
    endif

    if (iprint > 2) then
      write (stdout, '(1x,a)') 'Final centres (translated to home cell for writing xyz file)'
      do iw = 1, num_wann
        write (stdout, 888) iw, (wc(ind, iw)*lenconfac, ind=1, 3)
      end do
      write (stdout, '(1x,a78)') repeat('-', 78)
      write (stdout, *)
    endif

    xyz_unit = io_file_unit()
    open (xyz_unit, file=trim(seedname)//'_centres.xyz', form='formatted')
    write (xyz_unit, '(i6)') num_wann + num_atoms
    call io_date(cdate, ctime)
    write (xyz_unit, *) 'Wannier centres, written by Wannier90 on'//cdate//' at '//ctime
    do iw = 1, num_wann
      write (xyz_unit, '("X",6x,3(f14.8,3x))') (wc(ind, iw), ind=1, 3)
    end do
    do nsp = 1, num_species
      do nat = 1, atoms_species_num(nsp)
        write (xyz_unit, '(a2,5x,3(f14.8,3x))') atoms_symbol(nsp), atoms_pos_cart(:, nat, nsp)
      end do
    end do
    close (xyz_unit)

    write (stdout, '(/a)') ' Wannier centres written to file '//trim(seedname)//'_centres.xyz'

    return

888 format(2x, 'WF centre', i5, 2x, '(', f10.6, ',', f10.6, ',', f10.6, ' )')

  end subroutine wann_write_xyz