본문으로 건너뛰기

Brillouin zone 샘플링과 smearing

앞 장에서 주기계의 고유상태에 k\mathbf k 라벨이 붙는 것을 보았다. 전자 밀도·총에너지 같은 물리량은 Brillouin zone(BZ) 전체에 대한 적분으로 정의되는데, 실제 계산은 유한 개의 k\mathbf k-점에서의 합으로 근사할 수밖에 없다. 이 장은 그 근사의 표준 도구 — Monkhorst–Pack 격자, 대칭 축약, 그리고 금속에서 필수가 되는 smearing — 를 정리하고, 차원(3D/2D/1D)과 물질(금속/절연체)에 따라 k-grid를 어떻게 설계하는지의 실무 규칙을 세운다. 이론 장이므로 입력 파일과 실행 섹션은 없다.

학습 목표

  • BZ 적분을 유한 k-grid 가중합으로 바꾸는 절차와 Monkhorst–Pack 격자의 구조를 이해한다.
  • 진공 방향은 k-점 1개로 충분한 이유와 이 튜토리얼의 1×1×641 \times 1 \times 64 grid의 근거를 설명한다.
  • 금속의 k-수렴이 절연체보다 느린 이유를 Fermi 면의 불연속 점유로 설명한다.
  • smearing 4종(Fermi–Dirac, Gaussian, Methfessel–Paxton, tetrahedron)의 특성과 용도를 구분한다.
  • smearing이 수렴 보조 수단이며 free energy와 E(σ0)E(\sigma \to 0)을 구분해야 함을 이해한다.

1. BZ 적분에서 유한 k-grid 합으로

주기계의 전자 밀도는 밴드와 BZ에 대한 이중 합·적분이다.

n(r)=nBZdkΩBZ  fnkψnk(r)2        nkwkfnkψnk(r)2n(\mathbf r) = \sum_n \int_{\mathrm{BZ}} \frac{d\mathbf k}{\Omega_{\mathrm{BZ}}}\; f_{n\mathbf k}\, |\psi_{n\mathbf k}(\mathbf r)|^2 \;\;\longrightarrow\;\; \sum_n \sum_{\mathbf k} w_{\mathbf k}\, f_{n\mathbf k}\, |\psi_{n\mathbf k}(\mathbf r)|^2

적분을 가중치 wkw_{\mathbf k}를 가진 유한 합으로 바꾸는 것이 k-샘플링이다. 표준은 Monkhorst–Pack(MP) 격자 — BZ를 각 역격자 방향으로 N1×N2×N3N_1 \times N_2 \times N_3 등분한 균일 격자다. 1차원으로 쓰면 NN-분할 격자점은 다음과 같다.

ki=2iN12N2πa,i=1,,Nk_i = \frac{2i - N - 1}{2N}\,\frac{2\pi}{a}, \qquad i = 1, \ldots, N

원래 정의의 격자점은 분할 수의 홀짝에 따라 위치가 달라진다. 홀수 분할이면 Γ점(k=0\mathbf k = 0)을 포함하고, 짝수 분할이면 격자가 반 칸 이동해 Γ를 포함하지 않는다. Γ-centered 격자는 분할 수와 무관하게 항상 Γ를 포함하도록 정렬한 변형이다. 대부분의 경우 둘의 수렴 성능은 비슷하지만, hexagonal 격자에서는 이동된 짝수 MP grid가 결정 대칭과 어긋난 k-집합을 만들 수 있어 Γ-centered를 쓰는 것이 안전한 표준이다.

균일 격자의 k-점들이 모두 독립은 아니다. 시간역전 대칭(εn(k)=εn(k)\varepsilon_n(\mathbf k) = \varepsilon_n(-\mathbf k))과 결정 점군 대칭으로 서로 등가인 점들을 묶으면, 대표점들의 집합인 **기약 BZ(irreducible BZ, IBZ)**와 각 대표점의 등가 개수에 비례하는 가중치 wkw_{\mathbf k}만 남는다. 코드가 자동으로 수행하며, 대칭이 높은 계에서는 실제 계산할 k-점 수가 명목 grid의 수 분의 일로 줄어든다. SIESTA에서는 %block kgrid_Monkhorst_Pack, VASP에서는 KPOINTS 파일이 이 grid를 지정한다.

분할 수를 정하는 척도는 실공간 셀 길이와의 곱이다. 앞 장의 zone folding에서 보았듯 셀을 mm배 늘리면 BZ가 1/m1/m로 접히므로, Ni×aiN_i \times a_i(분할 수 × 셀 변 길이)를 일정하게 유지하면 k-공간 해상도가 유지된다. 단위셀에서 1×1×641 \times 1 \times 64가 필요했던 계라면, 같은 물질의 8배 supercell에는 1×1×81 \times 1 \times 8이면 같은 해상도다.

2. 차원별 원칙 — grid의 모양은 계의 주기성을 따른다

k-샘플링은 주기적인 방향에만 의미가 있다. 진공으로 분리된 방향은 이웃 상(image) 사이에 파동함수가 겹치지 않으므로 밴드가 그 방향으로 분산이 없고(평평하고), 어느 k\mathbf k에서 계산해도 같은 값이 나온다 — 여러 점을 찍는 것은 같은 계산의 반복일 뿐이다. 따라서 grid의 모양은 계의 차원을 그대로 반영한다.

grid 형태근거
3D bulkN×N×NN \times N \times N (셀 변 길이에 반비례하도록 조정)세 방향 모두 주기적
2D slabN×N×1N \times N \times 1표면 수직 방향은 진공 — 분산 없음
1D chain1×1×N1 \times 1 \times N사슬 축만 주기적, xyxy는 진공

이 튜토리얼의 1×1×641 \times 1 \times 64(챕터 01)가 정확히 셋째 행이다. xyxy 진공 방향은 Γ 1점, 주기 방향인 zz만 64점으로 조밀하게 샘플링한다. zz에 64점이 필요한 이유는 다음 절의 주제 — 이 사슬이 금속이기 때문이다.

금속과 절연체의 k-수렴은 질적으로 다르다. 절연체·반도체는 모든 밴드가 완전히 차 있거나 완전히 비어 있어 피적분 함수가 k\mathbf k에 대해 매끄럽고, 적은 k-점으로도 빠르게(사실상 지수적으로) 수렴한다. 금속은 Fermi 면이 BZ를 가로지르며 점유가 1에서 0으로 불연속하게 바뀌므로, 피적분 함수의 불연속을 유한 격자로 해상해야 해서 수렴이 느리다. 1D 금속 사슬에서는 Fermi "면"이 Fermi 으로 축소되어 그 근처 샘플링에 특히 민감하다 — 챕터 01 연습문제 1의 k-수렴 테스트에서 직접 확인하는 내용이다.

두 가지 구분을 함께 정리해 둔다. 첫째, 밴드 구조 계산의 k-경로는 SCF grid와 별개다 — SCF는 BZ 전체를 대표하는 균일 grid로 밀도를 수렴시키고, 밴드는 그 수렴된 밀도 위에서 고대칭 선분 경로를 따라 고유값만 계산한다(챕터 02). 둘째, transport 워크플로우에서 전극 계산의 k-grid는 device 계산과 일치시켜야 하는 제약이 따로 있고(챕터 05), TBtrans의 k-grid 지정에는 위치 규칙이라는 함정이 있다(챕터 07의 함정 2).

3. Smearing — 왜 필요한가

smearing된 점유 함수

그림 1. T=0T=0의 계단 점유 함수(검정)와 smearing된 점유(파랑·주황) — smearing 폭 σ\sigma가 클수록 Fermi면의 불연속이 부드러워져 k-적분이 빨리 수렴하지만, 물리량은 σ0\sigma\to0 외삽으로 판단해야 한다. (개념도)

금속의 k-수렴이 느린 근본 원인은 T=0T = 0 점유 함수가 계단이라는 것이다: fnk=θ(EFεnk)f_{n\mathbf k} = \theta(E_F - \varepsilon_{n\mathbf k}). 유한 grid에서 k-점을 하나 추가하거나 옮길 때마다 어떤 상태가 Fermi 준위를 넘나들며 점유가 0과 1 사이를 불연속하게 점프하고, 총에너지가 grid에 따라 요동한다. 처방은 점유를 폭 σ\sigma의 매끄러운 함수로 바꾸는 것이다.

fnk=f ⁣(εnkμσ)f_{n\mathbf k} = f\!\left(\frac{\varepsilon_{n\mathbf k} - \mu}{\sigma}\right)

점유가 매끄러워지면 피적분 함수도 매끄러워져 성긴 grid로도 안정적인 값이 나온다. 대가는 결과가 σ\sigma에 의존하게 된다는 것이다.

또 하나 유의할 것은 변분 일관성이다. 매끄러운 점유와 변분적으로 일관된 양은 내부 에너지 EE가 아니라 일반화 엔트로피 항을 포함한 free energy다.

F=EσS[f]F = E - \sigma\, S[f]

힘과 응력도 FF의 미분으로 계산해야 원자를 움직였을 때의 에너지 변화와 일치한다. 코드 출력에 free energy가 따로 찍히는 이유가 이것이다 — SIESTA stdout의 FreeEng 열(챕터 01), VASP OSZICAR의 F(챕터 03). FFEE의 차이는 σ\sigma가 커질수록 벌어지므로, 그 크기가 곧 smearing 오염의 척도다.

4. Smearing의 종류와 특성

방법점유 함수특징VASPSIESTA
Fermi–Dirac[1+ex]1[1 + e^{x}]^{-1}물리적 전자 온도 T=σ/kBT = \sigma/k_B의 열적 점유 그대로 — 유일하게 σ\sigma에 물리적 해석 존재. 꼬리가 길어 같은 σ\sigma 대비 오염 큼ISMEAR = -1OccupationFunction FD + ElectronicTemperature
Gaussianerf 기반무난한 범용. 금속/반도체 판별이 안 된 시스템, relax에 안전한 기본값ISMEAR = 0
Methfessel–PaxtonHermite 다항 보정EEσ\sigma 의존을 고차로 억제 — 금속 relax의 표준. 점유가 음수가 되거나 1을 넘을 수 있어(비물리적) 점유를 해석하는 분석에는 주의ISMEAR = 1, 2OccupationFunction MP
tetrahedron (+Blöchl 보정)점유 함수 아님 — BZ를 사면체로 분할해 선형 보간 적분smearing 파라미터 자체가 불필요. 정밀 총에너지·DOS에 최적. 힘·응력에는 변분적 일관성이 없어 relax에 부적합, k-점 최소 개수 필요ISMEAR = -5미지원

두 코드의 파라미터 대응은 단위만 다르다. VASP은 SIGMA를 eV로 직접 지정하고, SIESTA는 ElectronicTemperature를 온도(또는 에너지 단위)로 지정하며 σ=kBT\sigma = k_B T 관계로 환산된다. 300 K25.9 meV300\ \mathrm{K} \approx 25.9\ \mathrm{meV}이고, 챕터 03의 SIGMA = 0.05 eV는 약 580 K580\ \mathrm{K}에 해당한다. SIESTA는 점유 함수를 지정하지 않으면 Fermi–Dirac을 쓴다.

SIESTA의 디폴트 smearing

SIESTA의 디폴트는 OccupationFunction FD + ElectronicTemperature 300 K, 즉 σ25.9 meV\sigma \approx 25.9\ \mathrm{meV}의 Fermi–Dirac smearing이다. 아무것도 지정하지 않은 챕터 01 계산도 이 설정으로 돌고 있으며, stdout의 FreeEng 열이 총에너지와 미세하게 다른 값을 가지는 이유가 이 엔트로피 항이다.

Methfessel–Paxton의 음의 점유

MP smearing은 에너지 정확도를 위해 점유 함수에 음수 영역을 허용한다. 총에너지·힘에는 문제가 없지만, 점유수 자체를 물리량으로 읽는 분석(부분 전하, 자화 분해 등)에서는 인위적인 값이 섞일 수 있다. gap이 있는 시스템에 큰 σ\sigma의 MP를 쓰면 gap 가장자리 상태에 음의 점유가 걸리는 전형적 함정이 있다 — 반도체가 될 수 있는 시스템에는 Gaussian이 안전하다.

5. 실무 규칙

용도별 선택.

  • 확실한 금속의 relax — Methfessel–Paxton + 작은 σ\sigma (0.05–0.2 eV). 매끄러운 점유가 주는 안정적인 힘이 관건이다.
  • 금속/반도체 판별 전, 또는 relax 일반 — Gaussian + 작은 σ\sigma. 챕터 03의 ISMEAR = 0, SIGMA = 0.05가 이 규칙의 적용이다.
  • 고정 구조의 정밀 총에너지·DOS — 코드가 지원하면 tetrahedron(+Blöchl 보정)이 최선이다. smearing 오염 자체가 없다. 단 relax에는 쓰지 않는다.
  • 유한 온도 물성, NEGF 전극의 점유 — Fermi–Dirac. σ\sigma가 물리 파라미터가 되는 유일한 경우다(아래 참조).

DOS처럼 스펙트럼의 모양을 보는 계산은 총에너지보다 조밀한 k-grid를 요구한다는 점도 함께 기억한다 — 총에너지는 적분이라 오차가 상쇄되지만, DOS는 에너지별로 해상해야 하기 때문이다.

smearing은 수렴 보조 수단이지 물리 파라미터가 아니다. 우리가 원하는 답은 σ0\sigma \to 0, k-grid 무한 조밀 극한의 T=0T = 0 에너지다. σ\sigma는 그 극한에 유한한 grid로 접근하기 위한 수치 장치일 뿐이며, σ\sigma를 키우면 k-수렴은 빨라지지만 답 자체가 오염된다 — 두 파라미터는 함께 검사해야 한다. 예외는 Fermi–Dirac 하나다. 실제 전자 온도의 열적 점유를 의도적으로 넣는 경우(예: 유한 온도 물성, NEGF 전극의 Fermi 함수)에만 σ\sigma가 물리 파라미터가 된다.

오염의 크기는 3절에서 본 FFEE의 차이로 가늠하고, 필요하면 σ0\sigma \to 0 외삽으로 제거한다. Gaussian smearing에서는 E(σ0)(F+E)/2E(\sigma \to 0) \approx (F + E)/2라는 외삽 관계가 성립하며, VASP OSZICAR의 E0이 바로 이 외삽값이다 — 챕터 03에서 "구조 간 에너지 비교는 E0으로 한다"고 한 근거가 이것이다. FFEE의 차이가 무시할 수 없이 크면 σ\sigma가 너무 크다는 신호다.

k-수렴 테스트 절차. 새 시스템에서는 항상 한 번 수행한다.

  1. 기저·cutoff 등 나머지 파라미터를 고정한다.
  2. 주기 방향의 분할 수만 증가 시퀀스(예: 1D에서 kzk_z = 16, 32, 64, 128)로 바꿔가며 총에너지를 계산한다.
  3. 원자당 총에너지를 표로 만들고, 인접한 두 설정의 차이가 목표 임계(관례적으로 원자당 1 meV) 이하로 내려가는 최소 grid를 채택한다.
  4. 금속이면 σ\sigma 두 값(예: 표준값과 그 절반)에서 같은 표를 반복해, 채택한 grid가 σ\sigma 선택에 둔감한지 확인한다.

정리하는 표의 틀은 다음과 같다 (값은 각자의 계산으로 채운다).

kzk_zEE (eV/atom)직전 대비 차이 (meV/atom)판정
16
32
641 미만채택
128확인용

에너지 차이(결합 에너지, 상대 안정성)를 비교하는 계산이라면 비교 대상들이 같은 grid·같은 smearing을 쓰는 것이 절대값 수렴보다 먼저다 — 오차 상쇄를 살리는 배치다. 이 절차의 실습이 챕터 01 연습문제 1이다.

핵심 요약

개념내용
BZ 적분BZdkkwk\int_{\mathrm{BZ}} d\mathbf k \to \sum_{\mathbf k} w_{\mathbf k} — MP 균일 격자 + IBZ 대칭 축약
Γ-centered항상 Γ 포함. hexagonal 셀에서는 필수적 선택
grid 밀도 척도Ni×aiN_i \times a_i 일정 — 셀이 커지면 분할 수는 반비례로 줄인다
grid 모양주기 방향만 샘플링 — 3D는 N3N^3, slab은 N×N×1N \times N \times 1, 1D chain은 1×1×N1 \times 1 \times N
금속 vs 절연체절연체는 빠른 수렴, 금속은 Fermi 면 불연속으로 느림 — 조밀 grid + smearing
smearing계단 점유를 폭 σ\sigma로 완화. 변분 일관량은 F=EσSF = E - \sigma S
종류FD(물리적 온도), Gaussian(범용), MP(금속 relax, 음의 점유 주의), tetrahedron(정밀 에너지·DOS, relax 부적합)
코드 대응VASP ISMEAR/SIGMA ↔ SIESTA OccupationFunction/ElectronicTemperature (σ=kBT\sigma = k_B T)
원칙smearing은 수렴 보조 수단 — 목표는 σ0\sigma \to 0 극한. E0 = Gaussian 외삽값

연습문제

  1. 1D toy 모델의 k-수렴. 1D tight-binding 밴드 ε(k)=2tcos(ka)\varepsilon(k) = -2t\cos(ka), half-filling(EF=0E_F = 0)에서 밴드 에너지 aπ0π/aε(k)f(ε(k))dk\frac{a}{\pi}\int_0^{\pi/a} \varepsilon(k)\, f(\varepsilon(k))\, dk를 균일 k-grid 합으로 계산하는 코드를 짜라. σ=0\sigma = 0(계단 점유)과 Gaussian smearing(σ=0.1t\sigma = 0.1t)에서 grid 수를 늘려가며 수렴 속도를 비교하라.
  2. Γ 포함 여부. 1절의 MP 격자점 공식에서 NN이 홀수일 때만 k=0k = 0이 포함됨을 확인하고, 짝수 NN의 예로 Γ-centered 격자와 점 배치가 어떻게 달라지는지 그림으로 그려 보라.
  3. σ\sigma 외삽. 어떤 금속 계산에서 σ\sigma = 0.05, 0.1, 0.2 eV에 대해 FFEE를 얻었다고 하자. 3절의 F=EσSF = E - \sigma S 구조를 바탕으로, FFEEσ0\sigma \to 0에서 같은 값으로 수렴해야 하는 이유와 (F+E)/2(F+E)/2가 좋은 외삽인 이유를 설명하라.
  4. 단위 환산. SIGMA = 0.2 eV는 몇 K의 전자 온도에 해당하는가? 반대로 실온(300 K)의 kBTk_B T는 몇 meV인가? 이 비교로부터, relax에 쓰는 smearing 폭이 대부분 실제 실험 온도보다 훨씬 높은 가상의 온도에 해당함을 확인하라.
  5. grid 설계. 한 변 4 Å의 정육면체 셀을 가진 3D 금속 bulk에 12312^3 grid가 적절했다면, 같은 물질의 slab(면내 셀은 동일, 수직 방향 진공 20 Å)에는 어떤 grid를 제안하겠는가? 진공 방향에 k-점을 여러 개 쓰면 낭비 외에 결과가 달라지는 것이 있는지도 논하라. 여기까지가 Fundamentals의 이론 준비다. 나노 소자에서 전류가 어떻게 흐르는지의 개념적 뼈대 — 챕터 00 — 양자 수송 개관 — 부터 본편이 시작된다.

Ref: Monkhorst & Pack, Phys. Rev. B 13, 5188 (1976); Methfessel & Paxton, Phys. Rev. B 40, 3616 (1989); Blöchl, Jepsen & Andersen, Phys. Rev. B 49, 16223 (1994).