왜 진동수를 계산하는가
최적화 뒤에 진동수를 계산하는 이유는 크게 셋이다.
- 임계점의 성격 확인: 모든 진동수가 양수이면 안정한 최소, 허수 진동수가 하나 있으면 전이 상태다.
- 열역학 보정량 산출: 영점 진동에너지(ZPE), 열용량 보정, 엔트로피 → 깁스 자유에너지 G.
- IR·Raman 스펙트럼 예측: 실험과 직접 비교 가능한 진동 스펙트럼.
기본 사용법: AnFreq와 NumFreq
방식은 두 가지다.
| 키워드 | 방식 | 특징 |
|---|---|---|
Freq 또는 AnFreq | 해석적 미분 | 빠르고 정확. HF, DFT, MP2, CASSCF에서 지원. |
NumFreq | 수치 미분 | 모든 방법에서 사용 가능. 6 × Natom배의 그래디언트 계산. |
# 해석적 진동수 (가능한 경우 권장)
! B3LYP D4 def2-TZVP Freq
# 수치 진동수 (CCSD(T) 등 해석적 미분이 없는 경우)
! CCSD(T) cc-pVDZ NumFreq
진동수 계산은 그래디언트가 0인 점(stationary point)에서만 의미가 있다. 최적화 안 한 구조에서 돌리면 회전·병진 모드가 0이 아닌 값으로 새어 나오고, 그게 모든 진동수에 오차로 섞인다. 반드시 같은 메소드/기저집합으로 먼저 최적화해 둬야 한다.
최적화와 동시에
가장 흔한 패턴은 최적화와 진동수를 한 줄에 묶는 것이다.
! B3LYP D4 def2-TZVP RIJCOSX def2/J TightOpt Freq
* xyzfile 0 1 mol.xyz
최적화가 끝나면 ORCA가 같은 레벨로 진동수 계산을 자동으로 이어 간다. 허수 진동수가 나오면 그 모드를 따라 구조를 살짝 비틀어 다시 최적화하면 대개 풀린다.
출력 해석
진동수 계산이 끝나면 출력에 이런 표가 뜬다.
-----------------------
VIBRATIONAL FREQUENCIES
-----------------------
Scaling factor for frequencies = 1.000000000 (already applied!)
0: 0.00 cm**-1
1: 0.00 cm**-1
2: 0.00 cm**-1
3: 0.00 cm**-1
4: 0.00 cm**-1
5: 0.00 cm**-1
6: 1659.34 cm**-1
7: 3805.21 cm**-1
8: 3914.47 cm**-1
처음 6개(선형 분자는 5개)는 병진·회전 모드라 늘 0에 가까운 값이 찍힌다. 그 다음부터가 진짜 진동 모드다. 위 예에서 물 분자는 굽힘(1659), 대칭 신축(3805), 반대칭 신축(3914) 세 모드를 가진다.
정규 좌표 시각화
각 모드의 변위 벡터는 “NORMAL MODES” 섹션에 찍힌다.
------------
NORMAL MODES
------------
0 1 2 3 4
0 0.000000 0.000000 -0.000000 0.000000 -0.000000
1 0.000000 0.000000 0.000000 -0.000000 0.000000
...
정규 좌표를 눈으로 보고 싶으면 ORCA의 orca_pltvib 유틸리티로 진동 trajectory를 만들면 된다.
# 6번 모드의 trajectory 추출 (XYZ로 저장)
orca_pltvib mol.hess 6
허수 진동수 처리
출력에서 허수 진동수는 음수로 표기된다.
6: -456.78 cm**-1 ***imaginary mode***
7: 1234.56 cm**-1
8: ...
전이 상태도 아닌데 허수 진동수가 나왔다면, 이렇게 푼다.
- 허수 모드의 정규 좌표 방향으로 모든 원자를 약간(예: ± 0.1 Å) 밀어낸 새 좌표를 만든다.
orca_pltvib로 첫 진폭의 구조를 뽑으면 편하다. - 이 구조에서 다시 최적화한다. 임계값을
VeryTightOpt로 올려도 좋다. - 다시 진동수를 계산해 전부 양수인지 확인한다.
±10 cm⁻¹ 정도의 아주 작은 음의 진동수는 적분 격자나 SCF 임계값에서 오는 수치 잡음인 경우가 많다.
DefGrid3 VeryTightSCF로 격자와 SCF 정밀도를 올려 다시 돌리면 대개 사라진다.
열역학량과 깁스 자유에너지
진동수 계산이 끝나면 ORCA가 열역학량을 자동으로 계산해 출력한다. 기본 조건은 298.15 K, 1 atm이다.
-------------------------
THERMOCHEMISTRY AT 298.15K
-------------------------
Temperature ... 298.15 K
Pressure ... 1.00 atm
Total Mass ... 18.01 AMU
Zero point energy ... 0.02143521 Eh 13.45 kcal/mol
Thermal vibrational correction ... 0.00010256 Eh 0.06 kcal/mol
Thermal rotational correction ... 0.00141572 Eh 0.89 kcal/mol
Thermal translational correction ... 0.00141572 Eh 0.89 kcal/mol
-------------------------------------------------
Total thermal energy -76.39541321 Eh
Total enthalpy -76.39447892 Eh
Final Gibbs free energy -76.42157345 Eh
-------------------------------------------------
G-E(el) ... 0.00219978 Eh 1.38 kcal/mol
주요 항목의 뜻은 이렇다.
| 항목 | 의미 |
|---|---|
| Zero point energy | 영점 진동에너지(ZPE). 절대 영도에서도 남는 진동. |
| Total thermal energy | U = Eel + ZPE + thermal corrections. |
| Total enthalpy | H = U + pV. |
| Final Gibbs free energy | G = H − TS. 반응 자유에너지 계산의 핵심 값. |
온도·압력 변경
%freq
Temp 373.15 # 100 °C
Pressure 1.0 # atm
end
기본 RRHO(rigid-rotor harmonic oscillator) 모델은 저주파수(< 100 cm⁻¹) 모드에서 엔트로피를 과대평가한다.
Grimme의 quasi-RRHO 보정을 켜면 자유에너지가 한결 정확해진다.
%freq QuasiRRHO true CutOffFreq 35 end으로 켠다.
특히 컨포머 탐색이나 회전 자유도가 큰 분자에서 차이가 크다.
IR · Raman 스펙트럼 그리기
진동수 계산은 IR 강도를 자동으로 같이 출력한다. 라만 강도는 !NumFreq나 polarizability 계산이 같이 켜져 있어야 나온다.
스펙트럼을 그리려면 orca_mapspc 유틸리티를 쓴다.
# IR 스펙트럼 (.dat 파일 생성)
orca_mapspc mol.out ir -w25
# Raman 스펙트럼
orca_mapspc mol.out raman -w50
# 출력 파일 mol.out.ir.dat 또는 mol.out.raman.dat을 gnuplot, matplotlib 등으로 시각화
-w 옵션은 로렌츠 / 가우스 라인폭(cm⁻¹)을 정한다.
orca_mapspc를 인수 없이 실행하면 전체 옵션 목록이 나온다.
자주 쓰는 옵션
%freq
CentralDiff true # NumFreq에서 양측 차분 (기본 단측). 정확도 향상, 시간 2배.
Increment 0.005 # 변위 크기 (bohr). 기본 0.005.
Restart true # 중단된 NumFreq 재시작
ProjectTR true # 병진/회전 모드 투영 제거 (기본 true)
Mass2016 true # IUPAC 2016 원자 질량 사용
QuasiRRHO true # quasi-RRHO 엔트로피 보정
CutOffFreq 35.0 # quasi-RRHO 컷오프 (cm^-1), 기본 35
end
DFT가 내놓는 조화 진동수는 실험치보다 보통 5~10% 높게 나온다. 문헌과 비교할 때는 함수별 권장 스케일링 인자를 곱해 보정하는 게 좋다. 예를 들어 B3LYP/6-31G(d)는 약 0.961, B3LYP/def2-TZVP는 약 0.965다. NIST CCCBDB(cccbdb.nist.gov)에 종합 표가 있다.