TABLE OF CONTENTS

::pwtk::pwi::supercell_fccHKL

SYNOPSIS

proc ::pwtk::pwi::supercell_fccHKL {args} {

USAGE

   ::pwtk::pwi::supercell_fccHKL  ?-alat1x1 ALAT?  ?-nk1 NK1?  ?-kshift "S1 S2 S2"?  ?-adsize SIZE_BOHR?  ?-vacuum VACUUM_SPECS?  hkl  supercell  Nlayers

PURPOSE

Generate CELL_PARAMETERS and K_POINTS for a supercell of the fcc(hkl) slab. Set also dipole correction if so requested by the -vacuum option (see below).

It is assumed that the slab is oriented along the XY plane and Z is the surface normal direction.

BEWARE

All lengths are in Bohr units!

OPTIONS

Further explanation of the -vacuum option:

* if "vacuum-specs" is a single number, it corresponds to the vacuum thickness, i.e.: 

       -vacuum 16.0

* if "vacuum-specs" is a list with 3 numbers, it corresponds to 

   "vac1st diplayer vac2nd" (see ::pwtk::pwi::dipoleCorr), i.e.:
       -vacuum {15.0 2.0 16.0}   
   In this case dipole correction is switched on. If any of vac1st,
   diplayer, or vac2nd is an empty string, the corresponding default value is used.

ARGUMENTS

In addition to the math::linearalgebra form, the supercell can be specified in the diagonal NxM or in the a1.b1_a2.b2 notation. The latter corresponds to the following expansion matrix: {{a1 b1} {a2 b2}}

SOURCE

    set narg 3
    set usage "?-alat1x1 value?  ?-nk1 nk1?  ?-kshift {s1 s2 s2}?  ?-adsize size?  ?-vacuum vacuum-specs?  hkl supercell Nlayers"
    set options {
        {alat1x1.arg  {}  "alat of the 2D (1x1) unitcell"}
        {nk1.arg      {}  "nk1 of (nk1 x nk2) k-mesh of the 1x1-fcc(hkl) unitcell"}
        {kshift.arg   {}  "a list of k-mesh offsets"}
        {adsize.arg   0.0 "size of the largest adsorbate to consider"}
        {vacuum.arg   {}  "either vacuum thickness or {vac1st diplayer vac2nd} list"}
    }
    ::pwtk::parseOpt_
    ::pwtk::checkOType_ -alat1x1 $opt(alat1x1) {number posreal} "positive real number"
    ::pwtk::checkOType_ -nk1     $opt(nk1)     {number posint}  "positive integer"
    ::pwtk::checkOType_ -kshift  $opt(kshift)  {numberlist binary} "list of binary numbers"
    ::pwtk::checkOType_ -adsize  $opt(adsize)  {number nonnegreal} "non-negative number"
    ::pwtk::checkOType_ -vacuum  $opt(vacuum)  {numberlist posreal} "list of positive numbers"

    lassign $args hkl supercell Nlayers
    
    # -alat1x1
    if { $opt(alat1x1) eq {} } {
        set opt(alat1x1) [::pwtk::pwi::alat]
        if { $opt(alat1x1) == 0.0 } {
            ::pwtk::error "cannot generate a supercell; alat1x1 is not set" 1
        }
    }
    puts "debug: alat = $opt(alat1x1)"
    
    set ku [::pwtk::input::cardGetFlags   K_POINTS trim]
    set kp [::pwtk::input::cardGetContent K_POINTS]

    # -nk1
    set nk1 $opt(nk1)
    if { $nk1 eq {} } {
        if { [string match *automatic* $ku] } {            
            set nk1 [lindex $kp 0]
        } else {
            set nk1 12
        }
    }
    
    # -kshift
    if { $opt(kshift) eq {} } {
        if { [string match *automatic* $ku] } {            
            set opt(kshift) [lrange $kp 3 end]
        } else {
            set opt(kshift) "1 1 1"
        }
    } elseif { [llength $opt(kshift)] != 3 } {
        ::pwtk::error "wrong value of the -kshift option; got \"$opt(kshift)\"
but expected a list of 3 binary numbers" 1
    }

    # -vacuum
    if { $opt(vacuum) == {} } {
        set vacuum 16.0
    } else {
        set vlen [llength $opt(vacuum)]
        if { $vlen == 1 } {
            set vacuum $opt(vacuum)
        } elseif { $vlen == 3 } {
            lassign $opt(vacuum) vac1st diplayer vac2nd
        } else {
            ::pwtk::error "wrong value of the -vacuum option; got \"$opt(vacuum)\"
but expected either a single number or a list of 3 numbers" 1
        }
    }

    set path_ [namespace path]
    namespace path ::math::linearalgebra

    # fcc bulk lattice parameter
    set a0_ [expr $opt(alat1x1)*sqrt(2)]

    # fcc(HKL) unitcell vectors
    switch -- [string trim [string map {{ } {}} $hkl] ()] {
        111 {
            set unitvec "1.0  0.0
                        -0.5  [expr 0.5*sqrt(3)]"
            set unitKmesh "$nk1 $nk1"
            set d12_bohr  [expr $a0_*sqrt(3)/3.0]
        }
        110 - 101 - 011 {
            set unitvec "1.0  0.0
                         0.0  [expr sqrt(2)]"
            set unitKmesh "$nk1 [expr $nk1/sqrt(2)]"
            set d12_bohr  [expr $a0_*sqrt(2)/4.0]
        }
        100 - 010 - 001 {
            set unitvec "1.0  0.00
                         0.0  1.00"
            set unitKmesh "$nk1 $nk1"
            set d12_bohr  [expr $a0_/2.0]
        }
        default {
            ::pwtk::error "unsupported Miller indices $hkl, must be 111, 110, or 100" 1
        }
    }

    # get the size of C vector in Bohr

    set fullslab_height [expr ($Nlayers-1)*$d12_bohr + $opt(adsize)]; # in Bohr
    
    ::pwtk::ifexist vacuum {
        # "vacuum" variable exist only for the case of NO dipole-correction
        # (see above the parsing of the -vacuum option)
        set Z_bohr [expr $fullslab_height + $vacuum]
    } else {
        # dipole-correction case
        set Z_bohr [dipoleCorr -vac1st $vac1st -diplayer $diplayer -vac2nd $vac2nd  $fullslab_height]
    }

    # generate a supercell
    lassign [::pwtk::supercell2D $unitvec $supercell $unitKmesh] supervec kp

    # get supercell alat (in units of alat1x1)
    lassign $supervec v1
    set alat [norm $v1]
    set alatBohr [expr $alat*$opt(alat1x1)]
    
    # rescale the supercell vectors so that norm(v1) = 1.0 and set Z 
    set supervec [::pwtk::latvec2Dto3D [scale [expr 1.0/$alat] $supervec] [expr $Z_bohr/$alatBohr]]
    lassign $supervec v1 v2 v3
    
    set nb    [norm $v2]
    set nc    [norm $v3]
    set kmesh [list {*}$kp {*}$opt(kshift)]
    
    SYSTEM "ibrav = 0, celldm(1) = [expr $alatBohr], A = "
    CELL_PARAMETERS (alat) [::pwtk::m2cell $supervec]
    K_POINTS (automatic) $kmesh

    # print the basic info to stdout

    set sv2D [::pwtk::latvec3Dto2D [::pwtk::matricize $supercell]]
    set sv3D [::pwtk::latvec2Dto3D $sv2D 1.0]    
    set area [expr abs([det $sv3D])]
    
    ::pwtk::print ++++ "Auto-generation of the fcc($hkl) SLAB supercell..."
    ::pwtk::print "
Number of expected layers in the slab:  $Nlayers
Unitcell k-mesh:  [lmap k [::pwtk::Kmeshize $unitKmesh 2] {expr round($k)}]
    
2D unitcell lattice vectors (in Bohr):

[show [scale $opt(alat1x1) [::pwtk::matricize $unitvec]] %15.10f]
2D supercell expansion matrix:  
                             
[show $sv2D %3.0f]
Supercell:
   alat = [format %.6f $alatBohr] Bohr = [format %.6f $alat] unitcell-alat
   
   Size of vector A  = [format %8.4f 1.0] alat
   Size of vector B  = [format %8.4f $nb] alat
   Size of vector C  = [format %8.4f $nc] alat
   
   2D area of supercell = [format %8.4f $area] of the (1x1)-unitcell area
   
   Expected slab thickness = [format %.4f $fullslab_height] Bohr
   Expected adsorbate size = $opt(adsize) Bohr"
    
    ::pwtk::ifexist vacuum {
        ::pwtk::print "   Vacuum thickness        = $vacuum Bohr\n"
    } else {
        ::pwtk::print "   Vacuum thickness        = [join [::pwtk::dipcorrDefault_ $vac1st $diplayer $vac2nd] { + }] Bohr\n"
    }
    ::pwtk::print "The following cell parameters and k-points cards were generated:\n"
    
    ::pwtk::input::cardPrint CELL_PARAMETERS
    ::pwtk::input::cardPrint K_POINTS

    namespace path $path_
}