TABLE OF CONTENTS
ABINIT/conducti [ 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-2024 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)
SOURCE
25 #if defined HAVE_CONFIG_H 26 #include "config.h" 27 #endif 28 29 #include "abi_common.h" 30 31 program conducti 32 33 use defs_basis 34 use m_xmpi 35 use m_errors 36 use m_abicore 37 #if defined HAVE_MPI2 38 use mpi 39 #endif 40 use m_conducti 41 42 use m_io_tools, only : open_file 43 use m_time, only : timein 44 use m_fstrings, only : sjoin, itoa 45 use m_nctk, only : nctk_test_mpiio 46 use m_paw_optics,only : linear_optics_paw 47 48 implicit none 49 50 !Local variables------------------------------- 51 !scalars 52 integer,parameter :: master=0 53 integer :: incpaw,nproc,comm,inunt,my_rank,mpierr 54 real(dp) :: tcpu,tcpui,twall,twalli 55 character(len=fnlen) :: filnam,filnam_out 56 character(len=500) :: msg 57 !arrays 58 real(dp) :: tsec(2) 59 60 ! ********************************************************************************* 61 62 !Change communicator for I/O (mandatory!) 63 call abi_io_redirect(new_io_comm=xmpi_world) 64 65 !Start, MPI init, ... 66 call xmpi_init() 67 call timein(tcpui,twalli) 68 69 !Initialize memory profiling if it is activated 70 !if a full abimem.mocc report is desired, set the argument of abimem_init to "2" instead of "0" 71 !note that abimem.mocc files can easily be multiple GB in size so don't use this option normally 72 #ifdef HAVE_MEM_PROFILING 73 call abimem_init(0) 74 #endif 75 76 !Test if the netcdf library supports MPI-IO 77 call nctk_test_mpiio(print_warning=.false.) 78 79 !Extract MPI parallelization data 80 comm = xmpi_world 81 nproc = xmpi_comm_size(comm) 82 my_rank = xmpi_comm_rank(comm) 83 84 85 !Read some input data 86 if (my_rank==master) then 87 88 ! Read data file name 89 write(std_out,'(a)')' Please, give the name of the data file ...' 90 read(std_in, '(a)') filnam 91 write(std_out,'(2a)')' The name of the data file is: ',trim(filnam) 92 93 ! Read type of calculation 94 if (open_file(filnam,msg,newunit=inunt,form='formatted')==0) then 95 rewind(inunt) 96 read(inunt,*) incpaw 97 close(inunt) 98 else 99 ABI_ERROR('Error opening data file!') 100 end if 101 102 ! Read output file name 103 write(std_out,'(a)')' Give the name of the output file ...' 104 read(std_in, '(a)') filnam_out 105 write(std_out,'(2a)')' The name of the output file is: ',filnam_out 106 107 end if 108 109 !Broadcast input data 110 call xmpi_bcast(incpaw,master,comm,mpierr) 111 call xmpi_bcast(filnam,master,comm,mpierr) 112 call xmpi_bcast(filnam_out,master,comm,mpierr) 113 114 !Call main routine 115 if (incpaw==1) then 116 if (my_rank==master) then 117 call conducti_nc(filnam,filnam_out) 118 end if 119 elseif (incpaw==2) then 120 call conducti_paw(filnam,filnam_out) 121 elseif (incpaw==3) then 122 if (my_rank==master) then 123 call linear_optics_paw(filnam,filnam_out) 124 end if 125 elseif (incpaw==4) then 126 call conducti_paw(filnam,filnam_out) 127 call conducti_paw_core(filnam,filnam_out,with_absorption=.true.,with_emissivity=.true.) 128 elseif (incpaw==5) then 129 call conducti_paw_core(filnam,filnam_out,with_absorption=.true.) 130 elseif (incpaw==6) then 131 call conducti_paw_core(filnam,filnam_out,with_absorption=.true.,with_emissivity=.true.) 132 elseif (incpaw==42) then 133 call conducti_paw(filnam,filnam_out,varocc=1) 134 else 135 ABI_ERROR(sjoin("Wrong incpaw:", itoa(incpaw))) 136 end if 137 138 !End, memory cleaning 139 call timein(tcpu,twall) 140 tsec(1)=tcpu-tcpui 141 tsec(2)=twall-twalli 142 if (my_rank==0) then 143 write(std_out, '(a,a,a,f13.1,a,f13.1)' ) & 144 & '-',ch10,'- Proc. 0 individual time (sec): cpu=',tsec(1),' wall=',tsec(2) 145 end if 146 147 call abinit_doctor("__conducti") 148 149 call xmpi_end() 150 151 end program conducti