TABLE OF CONTENTS


ABINIT/conducti [ Programs ]

[ Top ] [ Programs ]

NAME

 conducti

FUNCTION

 This program computes the elements of the optical frequency dependent
 conductivity tensor and the conductivity along the three principal axes
 from the Kubo-Greenwood formula.

COPYRIGHT

 Copyright (C) 2006-2018 ABINIT group (FJ,SMazevet)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

INPUTS

  (main routine)

OUTPUT

  (main routine)

PARENTS

CHILDREN

      abi_io_redirect,abimem_init,abinit_doctor,conducti_nc,conducti_paw
      conducti_paw_core,destroy_mpi_enreg,emispec_paw,linear_optics_paw
      timein,xmpi_end,xmpi_init

SOURCE

 33 #if defined HAVE_CONFIG_H
 34 #include "config.h"
 35 #endif
 36 
 37 #include "abi_common.h"
 38 
 39 program conducti
 40 
 41  use defs_basis
 42  use defs_abitypes
 43  use m_xmpi
 44  use m_errors
 45  use m_abicore
 46 #if defined HAVE_MPI2
 47  use mpi
 48 #endif
 49  use m_conducti
 50 
 51  use m_io_tools,  only : open_file
 52  use m_time,      only : timein
 53  use m_fstrings,  only : sjoin, itoa
 54  use m_mpinfo,    only : destroy_mpi_enreg
 55  use m_paw_optics,only : linear_optics_paw
 56 
 57 !This section has been created automatically by the script Abilint (TD).
 58 !Do not modify the following lines by hand.
 59 #undef ABI_FUNC
 60 #define ABI_FUNC 'conducti'
 61 !End of the abilint section
 62 
 63  implicit none
 64 
 65 !Arguments -----------------------------------
 66 
 67 !Local variables-------------------------------
 68  integer :: incpaw,nproc,comm,inunt,my_rank
 69  real(dp) :: tcpu,tcpui,twall,twalli
 70  real(dp) :: tsec(2)
 71  character(len=fnlen) :: filnam,filnam_out
 72  character(len=500) :: msg
 73  type(MPI_type) :: mpi_enreg_seq ! needed for calling rwwf
 74 
 75 ! *********************************************************************************
 76 
 77 !Change communicator for I/O (mandatory!)
 78  call abi_io_redirect(new_io_comm=xmpi_world)
 79 
 80  call xmpi_init()
 81 
 82 !Initialize memory profiling if it is activated
 83 !if a full abimem.mocc report is desired, set the argument of abimem_init to "2" instead of "0"
 84 !note that abimem.mocc files can easily be multiple GB in size so don't use this option normally
 85 #ifdef HAVE_MEM_PROFILING
 86  call abimem_init(0)
 87 #endif
 88 
 89  call timein(tcpui,twalli)
 90 
 91  comm = xmpi_world
 92  nproc = xmpi_comm_size(comm); my_rank = xmpi_comm_rank(comm)
 93 
 94  if ( nproc > 1) then
 95    MSG_ERROR("conducti is not parallelized: run with one processor.")
 96  end if
 97 
 98 !Read data file name
 99  write(std_out,'(a)')' Please, give the name of the data file ...'
100  read(std_in, '(a)')filnam
101  write(std_out,'(a)')' The name of the data file is :',trim(filnam)
102 
103  if (open_file(filnam,msg,newunit=inunt,form='formatted') /= 0) then
104    MSG_ERROR(msg)
105  end if
106  rewind(inunt)
107  read(inunt,*) incpaw
108  close(inunt)
109 
110  write(std_out,'(a)')' Give the name of the output file ...'
111  read(std_in, '(a)') filnam_out
112  write(std_out,'(a)')' The name of the output file is :',filnam_out
113 
114  if (incpaw==1) then
115    call conducti_nc(filnam,filnam_out,mpi_enreg_seq)
116  elseif (incpaw==2) then
117    call conducti_paw(filnam,filnam_out,mpi_enreg_seq)
118  elseif (incpaw==3) then
119    call linear_optics_paw(filnam,filnam_out,mpi_enreg_seq)
120  elseif (incpaw==4) then
121    call conducti_paw(filnam,filnam_out,mpi_enreg_seq)
122    call conducti_paw_core(filnam,filnam_out,mpi_enreg_seq)
123    call emispec_paw(filnam,filnam_out,mpi_enreg_seq)
124  elseif (incpaw==5) then
125    call conducti_paw_core(filnam,filnam_out,mpi_enreg_seq)
126  elseif (incpaw==6) then
127    call emispec_paw(filnam,filnam_out,mpi_enreg_seq)
128  else
129    MSG_ERROR(sjoin("Wrong incpaw:", itoa(incpaw)))
130  end if
131 
132  call destroy_mpi_enreg(mpi_enreg_seq)
133 
134  call timein(tcpu,twall)
135 
136  tsec(1)=tcpu-tcpui
137  tsec(2)=twall-twalli
138  write(std_out, '(a,a,a,f13.1,a,f13.1)' ) &
139 & '-',ch10,'- Proc.   0 individual time (sec): cpu=',tsec(1),'  wall=',tsec(2)
140 
141  call abinit_doctor("__conducti")
142 
143  call xmpi_end()
144 
145  end program conducti