TABLE OF CONTENTS
- ABINIT/abitk
- abitk/abitk_show_help
- abitk/get_path_cryst
- abitk/get_path_ebands
- abitk/get_path_ebands_cryst
- abitk/parse_skw_params
ABINIT/abitk [ Programs ]
NAME
abitk
FUNCTION
Command line interface to perform very basic post-processing of output files (mainly netcdf files). Use `abitk --help` to get list of possible commands.
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (MG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
(main program)
OUTPUT
(main routine)
SOURCE
23 #if defined HAVE_CONFIG_H 24 #include "config.h" 25 #endif 26 27 #include "abi_common.h" 28 29 program abitk 30 31 use defs_basis 32 use m_abicore 33 use m_xmpi 34 use m_errors 35 use m_hdr 36 use m_ebands 37 use m_crystal 38 use m_kpts 39 #ifdef HAVE_NETCDF 40 use netcdf 41 #endif 42 use m_nctk 43 44 use defs_datatypes, only : ebands_t 45 use m_build_info, only : abinit_version 46 use m_fstrings, only : sjoin, strcat, basename, itoa 47 use m_io_tools, only : open_file, enforce_fortran_io 48 use m_specialmsg, only : herald 49 use m_symtk, only : matr3inv 50 use m_numeric_tools, only : arth 51 use m_bz_mesh, only : kpath_new, kpath_t 52 use m_unittests, only : tetra_unittests, kptrank_unittests, tetra_zinv_convergence 53 use m_argparse, only : get_arg, get_arg_list, parse_kargs 54 use m_common, only : ebands_from_file, crystal_from_file 55 use m_parser, only : geo_t, geo_from_poscar_path 56 use m_phgamma, only : find_ewin 57 58 implicit none 59 60 !Arguments ---------------------------- 61 !Local variables----------------------- 62 !scalars 63 integer,parameter :: master = 0 64 integer :: ii, nargs, comm, my_rank, nprocs, prtvol, fform, rdwr, prtebands, spin 65 integer :: kptopt, nshiftk, new_nshiftk, chksymbreak, nkibz, nkbz, intmeth, lenr !occopt, 66 integer :: ndivsm, abimem_level, ierr, ntemp, ios, itemp, use_symmetries, ltetra 67 real(dp) :: spinmagntarget, extrael, doping, step, broad, abimem_limit_mb, fs_ewin !, tolsym, tsmear 68 logical :: is_metal 69 character(len=500) :: command, arg, msg, ptgroup 70 character(len=fnlen) :: path, other_path, out_path !, prefix 71 type(hdr_type) :: hdr 72 type(ebands_t) :: ebands, ebands_kpath, other_ebands 73 type(edos_t) :: edos 74 type(gaps_t) :: gaps 75 type(jdos_t) :: jdos 76 type(crystal_t) :: cryst, other_cryst 77 type(geo_t) :: geo 78 type(kpath_t) :: kpath 79 !arrays 80 integer :: kptrlatt(3,3), new_kptrlatt(3,3), ngqpt(3) 81 real(dp) :: skw_params(4), tmesh(3) 82 real(dp),allocatable :: bounds(:,:), kTmesh(:), mu_e(:) 83 real(dp),allocatable :: shiftk(:,:), new_shiftk(:,:), wtk(:), kibz(:,:), kbz(:,:) 84 real(dp),allocatable :: n_ehst(:,:,:) 85 86 !******************************************************* 87 88 ! Change communicator for I/O (mandatory!) 89 call abi_io_redirect(new_io_comm=xmpi_world) 90 91 ! Initialize MPI 92 call xmpi_init() 93 comm = xmpi_world; my_rank = xmpi_comm_rank(comm); nprocs = xmpi_comm_size(comm) 94 95 if (my_rank == master) then 96 #ifdef FC_NAG 97 open(unit=ab_out,file="abitk.out",form='formatted',status='unknown', action="write", recl=ABI_RECL, iomsg=msg, iostat=ios) 98 #else 99 open(unit=ab_out,file="abitk.out",form='formatted',status='unknown', action="write", iomsg=msg, iostat=ios) 100 #endif 101 ABI_CHECK(ios == 0, msg) 102 rewind (unit=ab_out) 103 else 104 close(std_out) 105 if (open_file(NULL_FILE, msg, unit=std_out, action="write") /= 0) then 106 ABI_ERROR(msg) 107 end if 108 end if 109 110 ! Initialize memory profiling if it is activated 111 ! if a full abimem.mocc report is desired, set the argument of abimem_init to "2" instead of "0" 112 ! note that abimem.mocc files can easily be multiple GB in size so don't use this option normally 113 ABI_CHECK(get_arg("abimem-level", abimem_level, msg, default=0) == 0, msg) 114 ABI_CHECK(get_arg("abimem-limit-mb", abimem_limit_mb, msg, default=20.0_dp) == 0, msg) 115 #ifdef HAVE_MEM_PROFILING 116 call abimem_init(abimem_level, limit_mb=abimem_limit_mb) 117 #endif 118 119 nargs = command_argument_count() 120 ABI_CHECK(get_arg("prtvol", prtvol, msg, default=0) == 0, msg) 121 122 ! Command line options. 123 do ii=1,command_argument_count() 124 call get_command_argument(ii, arg) 125 if (arg == "--version") then 126 write(std_out,"(a)") trim(abinit_version); goto 100 127 128 else if (arg == "--enforce-fortran-io") then 129 call enforce_fortran_io(.True.) 130 131 else if (arg == "-h" .or. arg == "--help") then 132 ! Document the options. 133 call abitk_show_help() 134 goto 100 135 end if 136 end do 137 138 call get_command_argument(1, command) 139 140 select case (command) 141 case ("from_poscar") 142 call get_command_argument(2, path) 143 geo = geo_from_poscar_path(path, comm) 144 call geo%print_abivars(std_out) 145 call geo%free() 146 147 !case ("to_poscar") 148 ! call get_command_argument(2, path) 149 ! call get_path_cryst(path, cryst, comm) 150 ! call prtposcar(fcart, fnameradix, natom, ntypat, rprimd, typat, ucvol, xred, znucl) 151 152 !case ("to_cif") 153 ! call get_command_argument(2, path) 154 ! call get_path_cryst(path, cryst, comm) 155 ! call prt_cif(brvltt, ciffname, natom, nsym, ntypat, rprimd, spgaxor, spgroup, spgorig, symrel, tnon, typat, xred, znucl) 156 157 case ("hdr_print") 158 ABI_CHECK(nargs > 1, "FILE argument is required.") 159 call get_command_argument(2, path) 160 call hdr_read_from_fname(hdr, path, fform, comm) 161 ABI_CHECK(fform /= 0, "fform == 0") 162 rdwr = 3; if (prtvol > 0) rdwr = 4 163 call hdr%echo(fform, rdwr, unit=std_out) 164 165 case ("ibz") 166 ! Print list of kpoints in the IBZ with the corresponding weights 167 call get_path_cryst(path, cryst, comm) 168 169 call parse_kargs(kptopt, kptrlatt, nshiftk, shiftk, chksymbreak) 170 ABI_CHECK(any(kptrlatt /= 0), "kptrlatt or ngkpt must be specified") 171 172 call kpts_ibz_from_kptrlatt(cryst, kptrlatt, kptopt, nshiftk, shiftk, nkibz, kibz, wtk, nkbz, kbz, & 173 new_kptrlatt=new_kptrlatt, new_shiftk=new_shiftk) 174 new_nshiftk = size(new_shiftk, dim=2) 175 176 write(std_out, "(/, a)")" Input_kptrlatt | New_kptrlatt" 177 do ii=1,3 178 write(std_out, "(2(a, 3(i0,1x)), a)")" [", kptrlatt(ii, :), "], [ ", new_kptrlatt(ii, :), "]" 179 end do 180 write(std_out, "(/, a)")" Input_shiftk | New_shiftk" 181 do ii=1,max(nshiftk, new_nshiftk) 182 if (ii <= nshiftk .and. ii <= new_nshiftk) then 183 write(std_out, "(2(a, 3(f3.1, 1x)), a)")" [", shiftk(:, ii), "], [", new_shiftk(:, ii), "]" 184 else if (ii <= nshiftk .and. ii > new_nshiftk) then 185 write(std_out, "((a, 3(f3.1, 1x)), a)")" [", shiftk(:, ii), "] ......" 186 end if 187 end do 188 write(std_out,"(/,a)")"# List of kpoints in the IBZ with the corresponding weights:" 189 write(std_out,"(a,i0)")" kptopt ", kptopt 190 write(std_out,"(a,i0,/,a)")" nkpt ", nkibz, " kpt" 191 do ii=1,nkibz 192 write(std_out, "(3(es11.4),a,es11.4)") kibz(:, ii), " # wtk ", wtk(ii) 193 end do 194 195 !case ("testkgrid") 196 !call get_path_cryst(path, cryst, comm) 197 !call parse_kargs(kptopt, kptrlatt, nshiftk, shiftk, chksymbreak) 198 !call testkgrid(bravais, iout, kptrlatt, kptrlen, msym, nshiftk, nsym, prtkpt, rprimd, shiftk, symafm, symrel, vacuum) 199 200 case ("crystal_print") 201 call get_path_cryst(path, cryst, comm) 202 call cryst%print(unit=std_out, prtvol=prtvol) 203 204 case ("crystal_abivars") 205 call get_path_cryst(path, cryst, comm) 206 call cryst%print_abivars(std_out) 207 208 case ("ebands_print", "ebands_xmgrace", "ebands_gnuplot") 209 call get_path_ebands(path, ebands, comm) 210 if (command == "ebands_print") then 211 call ebands_print(ebands, unit=std_out, prtvol=prtvol) 212 else 213 prtebands = 1; if (command == "ebands_gnuplot") prtebands = 2 214 call ebands_write(ebands, prtebands, basename(path)) 215 end if 216 217 case ("ebands_gaps") 218 call get_path_ebands_cryst(path, ebands, cryst, comm) 219 call ebands_print_gaps(ebands, std_out) 220 221 case ("ebands_edos", "ebands_jdos") 222 call get_path_ebands_cryst(path, ebands, cryst, comm) 223 ABI_CHECK(get_arg("intmeth", intmeth, msg, default=2) == 0, msg) 224 ABI_CHECK(get_arg("step", step, msg, default=0.02 * eV_Ha) == 0, msg) 225 ABI_CHECK(get_arg("broad", broad, msg, default=0.06 * eV_Ha) == 0, msg) 226 227 if (command == "ebands_edos") then 228 edos = ebands_get_edos(ebands, cryst, intmeth, step, broad, comm) 229 call edos%print(std_out, header="Electron DOS") 230 out_path = strcat(basename(path), "_EDOS") 231 call wrtout(std_out, sjoin("Writing electron DOS to file:", out_path)) 232 call edos%write(out_path) 233 234 else if (command == "ebands_jdos") then 235 NOT_IMPLEMENTED_ERROR() 236 jdos = ebands_get_jdos(ebands, cryst, intmeth, step, broad, comm, ierr) 237 !call jdos%write(strcat(basename(path), "_EJDOS")) 238 end if 239 240 case ("ebands_bxsf") 241 call get_path_ebands_cryst(path, ebands, cryst, comm) 242 if (ebands_write_bxsf(ebands, cryst, strcat(basename(path), "_BXSF")) /= 0) then 243 ABI_ERROR("Cannot produce file for Fermi surface in BXSF format. Check log file for info.") 244 end if 245 246 !case ("ebands_nesting") 247 !call get_path_ebands_cryst(path, ebands, cryst, comm) 248 !if (ebands_write_nesting(ebands, cryst, filepath, prtnest, tsmear, fermie_nest, qpath_vertices, errmsg) /= 0) then 249 ! ABI_ERROR("Cannot produce file for nesting factor. Check log file for info.") 250 !end if 251 252 case ("skw_kpath") 253 ! Get energies on the IBZ from path file 254 call get_path_ebands_cryst(path, ebands, cryst, comm) 255 256 ! Generate k-path 257 ABI_CHECK(get_arg("ndivsm", ndivsm, msg, default=20) == 0, msg) 258 !ABI_COMMENT("Using hard-coded k-path because nkpath not present in input file.") 259 !nbounds = 5 260 ABI_MALLOC(bounds, (3, 5)) 261 bounds = reshape([zero, zero, zero, half, zero, zero, zero, half, zero, zero, zero, zero, zero, zero, half], [3,5]) 262 !ABI_CHECK(get_args_from_string("bounds", bounds, msg, default="[0,0,0,1/2,0,0,0,1/2]") == 0, msg) 263 kpath = kpath_new(bounds, cryst%gprimd, ndivsm) 264 ABI_FREE(bounds) 265 266 call kpath%print(header="Interpolating energies on k-path", unit=std_out) 267 268 ! Interpolate band energies with star-functions 269 call parse_skw_params(skw_params) 270 !ABI_CHECK(get_arg_list("einterp", ivec9, lenr, msg, exclude="ngkpt", want_len=9) == 0, msg) 271 !if (nint(dtset%einterp(1)) == 1) skw_params = dtset%einterp 272 ebands_kpath = ebands_interp_kpath(ebands, cryst, kpath, skw_params, [1, ebands%mband], comm) 273 274 !call wrtout(std_out, sjoin(" Writing interpolated bands to:", path) 275 ABI_CHECK(get_arg("prtebands", prtebands, msg, default=2) == 0, msg) 276 call ebands_write(ebands_kpath, prtebands, path) 277 278 !ebands_kmesh = ebands_interp_kmesh(ebands, cryst, skw_params, intp_kptrlatt, intp_nshiftk, intp_shiftk, & 279 ! band_block, comm, out_prefix) 280 !call ebands_free(ebands_kmesh) 281 282 case ("skw_compare") 283 ! Get energies on the IBZ from path 284 call get_path_ebands_cryst(path, ebands, cryst, comm) 285 286 ! Get ab-initio energies for the second file (assume k-path!) 287 call get_path_ebands_cryst(other_path, other_ebands, other_cryst, comm, argpos=3) 288 289 ! Interpolate band energies on the path with star-functions 290 kpath = kpath_new(other_ebands%kptns, other_cryst%gprimd, -1) 291 292 call parse_skw_params(skw_params) 293 ebands_kpath = ebands_interp_kpath(ebands, cryst, kpath, skw_params, [1, ebands%mband], comm) 294 295 ! Compare gaps 296 ABI_CHECK(get_arg("is-metal", is_metal, msg, default=.False.) == 0, msg) 297 if (.not. is_metal) then 298 write(std_out, "(2a)")" Will try to compare gaps. Use --is-metal option to skip this check.",ch10 299 call ebands_print_gaps(other_ebands, std_out, header="Ab-initio gaps") 300 call ebands_print_gaps(ebands_kpath, std_out, header="SKW interpolated gaps") 301 end if 302 303 !ABI_CHECK(get_arg("prtebands", prtebands, msg, default=2) == 0, msg) 304 !call ebands_write(ebands_kpath, prtebands, path) 305 306 ! Write EBANDS file 307 #ifdef HAVE_NETCDF 308 NCF_CHECK(ebands_ncwrite_path(other_ebands, cryst, "abinitio_EBANDS.nc")) 309 NCF_CHECK(ebands_ncwrite_path(ebands_kpath, other_cryst, "skw_EBANDS.nc")) 310 #endif 311 312 call wrtout(std_out, & 313 ch10//" Use `abicomp.py ebands abinitio_EBANDS.nc skw_EBANDS.nc -p combiplot` to compare the bands with AbiPy.", & 314 newlines=2) 315 316 case ("ebands_mu_T") 317 ! Get energies on the IBZ from filepath 318 call get_path_ebands_cryst(path, ebands, cryst, comm) 319 320 ABI_CHECK(get_arg("extrael", extrael, msg, default=zero) == 0, msg) 321 if (extrael == zero) then 322 ABI_CHECK(get_arg("doping", doping, msg, default=zero) == 0, msg) 323 ! Units of eph_doping is e_charge / cm^3 324 extrael = - doping * cryst%ucvol * (Bohr_meter * 100) ** 3 325 write(std_out, *)"Adding doping", doping 326 end if 327 328 !ABI_CHECK(get_arg("occopt", occopt, msg, default=3) == 0, msg) 329 !ABI_CHECK(get_arg("tsmear", tsmear, msg, default=tol2) == 0, msg) 330 ABI_CHECK(get_arg("spinmagntarget", spinmagntarget, msg, default=-99.99_dp) == 0, msg) 331 332 ebands%nelect = ebands%nelect + extrael 333 gaps = ebands_get_gaps(ebands, ierr) 334 335 ABI_CHECK(get_arg_list("tmesh", tmesh, lenr, msg, default_list=[5._dp, 59._dp, 6._dp] ) == 0, msg) 336 ntemp = nint(tmesh(3)) 337 ABI_CHECK_IGEQ(ntemp, 1, "ntemp <= 0") 338 ABI_MALLOC(kTmesh, (ntemp)) 339 kTmesh = arth(tmesh(1), tmesh(2), ntemp) * kb_HaK 340 341 ABI_MALLOC(mu_e, (ntemp)) 342 call ebands_get_muT_with_fd(ebands, ntemp, kTmesh, spinmagntarget, prtvol, mu_e, comm) 343 !mu_e = 6.715 * eV_Ha 344 345 call gaps%print(unit=std_out, header="KS gaps", kTmesh=kTmesh, mu_e=mu_e) 346 !stop 347 348 ABI_MALLOC(n_ehst, (2, ebands%nsppol, ntemp)) 349 call ebands_get_carriers(ebands, ntemp, kTmesh, mu_e, n_ehst) 350 351 !write(msg, "(a16,a32,a32)") 'Temperature [K]', 'e/h density [cm^-3]', 'e/h mobility [cm^2/Vs]' 352 do spin=1,ebands%nsppol 353 do itemp=1,ntemp 354 write(std_out, "(a, 2f16.2, 2e16.2)")& 355 " T (K), mu_e (eV), nh, ne", kTmesh(itemp) / kb_HaK, mu_e(itemp) * Ha_eV, & 356 n_ehst(2,spin,itemp) / cryst%ucvol / Bohr_cm**3, & 357 n_ehst(1,spin,itemp) / cryst%ucvol / Bohr_cm**3 358 end do 359 end do 360 361 ABI_CHECK(get_arg("intmeth", intmeth, msg, default=2) == 0, msg) 362 ABI_CHECK(get_arg("step", step, msg, default=0.02 * eV_Ha) == 0, msg) 363 ABI_CHECK(get_arg("broad", broad, msg, default=0.06 * eV_Ha) == 0, msg) 364 365 edos = ebands_get_edos(ebands, cryst, intmeth, step, broad, comm) 366 call edos%print(std_out, header="Electron DOS") 367 call edos%get_carriers(ntemp, kTmesh, mu_e, n_ehst) 368 369 !write(msg, "(a16,a32,a32)") 'Temperature [K]', 'e/h density [cm^-3]', 'e/h mobility [cm^2/Vs]' 370 do spin=1,ebands%nsppol 371 do itemp=1,ntemp 372 write(std_out, "(a, 2f16.2, 2e16.2)")& 373 " T (K), mu_e (eV), nh, ne", kTmesh(itemp) / kb_HaK, mu_e(itemp) * Ha_eV, & 374 n_ehst(2, itemp, spin) / cryst%ucvol / Bohr_cm**3, & 375 n_ehst(1, itemp, spin) / cryst%ucvol / Bohr_cm**3 376 end do 377 end do 378 379 ABI_FREE(kTmesh) 380 ABI_FREE(mu_e) 381 ABI_FREE(n_ehst) 382 383 !case ("ebands_dope") 384 385 ! ==================== 386 ! Tools for developers 387 ! ==================== 388 !case ("denpot_f2nc") 389 ! call get_command_argument(2, path) 390 ! call get_command_argument(3, other_path) 391 ! call denpot_f2nc(path, other_path) 392 393 ! Get important dimensions from the first header and rewind the file. 394 !call hdr_fort_read(new%hdr_ref, unt, fform) 395 !if (dvdb_check_fform(fform, "read_dvdb", msg) /= 0) then 396 ! ABI_ERROR(sjoin("While reading:", path, ch10, msg)) 397 !end if 398 !! Fortran IO 399 !do ispden=1,hdr1%nspden 400 ! read(units(jj), err=10, iomsg=msg) (v1(ifft), ifft=1,cplex*nfft) 401 ! write(ount, err=10, iomsg=msg) (v1(ifft), ifft=1,cplex*nfft) 402 !end do 403 !if (new%debug) call new%hdr_ref%echo(fform, 4, unit=std_out) 404 !call hdr1%free() 405 406 !! first order potentials are always written because the eph code requires them 407 !! the files are small (much much smaller that 1WFK, actually we should avoid writing 1WFK) 408 !rdwrpaw=0 409 !call appdig(pertcase,dtfil%fnameabo_pot,fi1o) 410 !! TODO: should we write pawrhoij1 or pawrhoij. Note that ioarr writes hdr%pawrhoij 411 !call fftdatar_write_from_hdr("first_order_potential",fi1o,dtset%iomode,hdr,& 412 !ngfftf,cplex,nfftf,dtset%nspden,vtrial1,mpi_enreg) 413 414 case ("find_ewin") 415 ltetra = 2 416 ABI_CHECK(get_arg("ltetra", ltetra, msg, default=2) == 0, msg) 417 call get_path_ebands_cryst(path, ebands, cryst, comm) 418 ! Use Q-mesh == K-mesh for double delta. 419 call find_ewin(ebands%nkpt, ebands%kptns, cryst, ebands, ltetra, fs_ewin, comm) 420 421 ! =========== 422 ! Unit tests 423 ! =========== 424 425 case ("tetra_unit_tests") 426 ABI_CHECK(get_arg("ptgroup", ptgroup, msg, default="m-3m") == 0, msg) 427 ABI_CHECK(get_arg_list("ngqpt", ngqpt, lenr, msg, default_list=[20, 20, 20]) == 0, msg) 428 ABI_CHECK(get_arg("use-symmetries", use_symmetries, msg, default=1) == 0, msg) 429 call tetra_unittests(ptgroup, ngqpt, use_symmetries, prtvol, comm) 430 431 case ("tetra_zinv_convergence") 432 ABI_CHECK(get_arg("ptgroup", ptgroup, msg, default="m-3m") == 0, msg) 433 !ABI_CHECK(get_arg_list("ngqpt", ngqpt, lenr, msg, default_list=[20, 20, 20]) == 0, msg) 434 ABI_CHECK(get_arg("use-symmetries", use_symmetries, msg, default=1) == 0, msg) 435 call tetra_zinv_convergence(ptgroup, use_symmetries, comm) 436 437 case ("kptrank_unit_tests") 438 ABI_CHECK(get_arg("ptgroup", ptgroup, msg, default="m-3m") == 0, msg) 439 ABI_CHECK(get_arg_list("ngqpt", ngqpt, lenr, msg, default_list=[20, 20, 20]) == 0, msg) 440 ABI_CHECK(get_arg("use-symmetries", use_symmetries, msg, default=1) == 0, msg) 441 call kptrank_unittests(ptgroup, ngqpt, use_symmetries, comm) 442 443 case default 444 call abitk_show_help() 445 ABI_ERROR(sjoin("Invalid command:", command)) 446 end select 447 448 ! Deallocate memory to make memcheck happy. 449 call hdr%free() 450 call cryst%free() 451 call other_cryst%free() 452 call kpath%free() 453 call ebands_free(ebands) 454 call ebands_free(ebands_kpath) 455 call ebands_free(other_ebands) 456 call edos%free() 457 call jdos%free() 458 call gaps%free() 459 460 ABI_SFREE(kibz) 461 ABI_SFREE(wtk) 462 ABI_SFREE(kbz) 463 ABI_SFREE(new_shiftk) 464 465 call abinit_doctor("__abitk") 466 467 100 call xmpi_end() 468 469 contains
abitk/abitk_show_help [ Functions ]
[ Top ] [ abitk ] [ Functions ]
NAME
abitk_show_help
FUNCTION
Show command line help
INPUTS
OUTPUT
SOURCE
485 subroutine abitk_show_help() 486 487 write(std_out,"(a)")" --version Show version number and exit." 488 write(std_out,"(a)")" -h, --help Show this help and exit." 489 write(std_out,"(a)")" -v Increase verbosity level" 490 491 write(std_out,"(2a)")ch10,"=== HEADER ===" 492 write(std_out,"(a)")"hdr_print FILE [--prtvol 0] Print ABINIT header." 493 494 write(std_out,"(2a)")ch10,"=== KPOINTS ===" 495 write(std_out,"(a)")"ibz FILE --ngkpt 2 2 2 or --kptrlatt [--kptopt 1] [--shiftk 0.5 0.5 0.5] [--chksymbreak 1]" 496 497 write(std_out,"(2a)")ch10,"=== CRYSTAL ===" 498 write(std_out,"(a)")"crystal_print FILE Print info on crystalline structure." 499 write(std_out,"(a)")"from_poscar POSCAR_FILE Read POSCAR file, print abinit variables." 500 501 write(std_out,"(2a)")ch10,"=== ELECTRONS ===" 502 write(std_out,"(a)")"ebands_print FILE Print info on electron band structure." 503 write(std_out,"(a)")"ebands_xmgrace FILE Produce XMGRACE file with electron bands." 504 write(std_out,"(a)")"ebands_gnuplot FILE Produce GNUPLOT file with electron bands." 505 write(std_out,"(a)")"ebands_edos FILE --intmeth, --step, --broad Compute electron DOS." 506 write(std_out,"(a)")"ebands_bxsf FILE Produce BXSF file for Xcrysden." 507 !write(std_out,"(a)")"ebands_mu_t FILE --occopt --tsmear --extrael Change number of electron, compute new Fermi level." 508 write(std_out,"(a)")"ebands_gaps FILE Print info on gaps" 509 !write(std_out,"(a)")"ebands_jdos FILE --intmeth, --step, --broad Compute electron DOS." 510 write(std_out,"(a)")"skw_kpath FILE Interpolate band structure along a (hardcoded) k-path." 511 write(std_out,"(a)")"skw_compare IBZ_WFK KPATH_WFK Use e_nk from IBZ_WFK to interpolate on the k-path in KPATH_WFK." 512 513 write(std_out,"(2a)")ch10,"=== DEVELOPERS ===" 514 write(std_out,"(a)")"tetra_unit_tests Run unit tests for tetrahedron routines." 515 write(std_out,"(a)")"kptrank_unit_tests Run unit tests for kptrank routines." 516 517 end subroutine abitk_show_help
abitk/get_path_cryst [ Functions ]
[ Top ] [ abitk ] [ Functions ]
NAME
get_path_cryst
FUNCTION
INPUTS
OUTPUT
SOURCE
600 subroutine get_path_cryst(path, cryst, comm) 601 602 !Arguments ------------------------------------ 603 !scalars 604 character(len=*),intent(out) :: path 605 type(crystal_t),intent(out) :: cryst 606 integer,intent(in) :: comm 607 608 !Arguments ---------------------------- 609 !Local variables----------------------- 610 integer :: nargs 611 612 nargs = command_argument_count() 613 ABI_CHECK(nargs > 1, "FILE argument is required.") 614 call get_command_argument(2, path) 615 cryst = crystal_from_file(path, comm) 616 617 end subroutine get_path_cryst
abitk/get_path_ebands [ Functions ]
[ Top ] [ abitk ] [ Functions ]
NAME
get_path_ebands
FUNCTION
INPUTS
OUTPUT
SOURCE
568 subroutine get_path_ebands(path, ebands, comm) 569 570 !Arguments ------------------------------------ 571 !scalars 572 character(len=*),intent(out) :: path 573 type(ebands_t),intent(out) :: ebands 574 integer,intent(in) :: comm 575 576 !Arguments ---------------------------- 577 !Local variables----------------------- 578 integer :: nargs 579 580 nargs = command_argument_count() 581 ABI_CHECK(nargs > 1, "FILE argument is required.") 582 call get_command_argument(2, path) 583 ebands = ebands_from_file(path, comm) 584 585 end subroutine get_path_ebands
abitk/get_path_ebands_cryst [ Functions ]
[ Top ] [ abitk ] [ Functions ]
NAME
get_path_ebands_cryst
FUNCTION
INPUTS
OUTPUT
SOURCE
532 subroutine get_path_ebands_cryst(path, ebands, cryst, comm, argpos) 533 534 !Arguments ------------------------------------ 535 !scalars 536 character(len=*),intent(out) :: path 537 type(ebands_t),intent(out) :: ebands 538 type(crystal_t),intent(out) :: cryst 539 integer,intent(in) :: comm 540 integer,intent(in),optional :: argpos 541 542 !Arguments ---------------------------- 543 !Local variables----------------------- 544 integer :: nargs, apos 545 546 apos = 2; if (present(argpos)) apos = argpos 547 nargs = command_argument_count() 548 ABI_CHECK(nargs >= apos, "FILE argument is required.") 549 call get_command_argument(apos, path) 550 cryst = crystal_from_file(path, comm) 551 ebands = ebands_from_file(path, comm) 552 553 end subroutine get_path_ebands_cryst
abitk/parse_skw_params [ Functions ]
[ Top ] [ abitk ] [ Functions ]
NAME
parse_skw_params
FUNCTION
INPUTS
OUTPUT
SOURCE
632 subroutine parse_skw_params(params) 633 634 !Arguments ---------------------------- 635 !Local variables----------------------- 636 real(dp),intent(out) :: params(4) 637 638 integer :: lpratio 639 real(dp) :: rcut, rsigma 640 641 !******************************************************* 642 643 ABI_CHECK(get_arg("lpratio", lpratio, msg, default=5) == 0, msg) 644 ABI_CHECK(get_arg("rcut", rcut, msg, default=zero) == 0, msg) 645 ABI_CHECK(get_arg("rsigma", rsigma, msg, default=zero) == 0, msg) 646 params = [one, dble(lpratio), rcut, rsigma] 647 648 end subroutine parse_skw_params