A multi-force-field stack design and validation framework for asphalt molecular dynamics: Part I. Theoretical background and general-purpose force-field comparison)

Taeyoung Yun

doi:10.22702/jkai.2026.16.1.6

Preview

Journal of the Korean Asphalt Institute. 30 June 2026. 65-78
https://doi.org/10.22702/jkai.2026.16.1.6

A multi-force-field stack design and validation framework for asphalt molecular dynamics: Part I. Theoretical background and general-purpose force-field comparison)

아스팔트 분자동역학 시뮬레이션을 위한 다중 역장 스택 설계 및 검증 프레임워크: Part I. 이론적 배경과 범용 역장 비교

Taeyoung Yun¹²^*

윤 태영¹²^*

¹Research Fellow, Highway and Transportation Research Institute, Korea Institute of Civil Engineering and Building Technology (Ph.D.)

²Professor, University of Science and Technology (Ph.D.)

¹한국건설기술연구원 도로연구실 연구위원(박사)

²과학기술연합대학원대학교 정교수(박사)

^{*Corresponding Author}

License (open-access, https://creativecommons.org/licenses/by-nc/4.0/):

This is an Open-Access article distributed under the terms of the Creative Commons Attribution Non-Commercial License (https://creativecommons.org/licenses/by-nc/4.0) which permits unrestricted non-commercial use, distribution, and reproduction in any medium, provided the original work is properly cited.

ABSTRACT

This paper is Part I of a two-part series that develops a force-field (FF) selection and validation framework for atomistic molecular dynamics (MD) simulations of asphalt binder systems that involve SARA hydrocarbons, polymeric modifiers, mineral aggregate surfaces, and aqueous ionic environments under oxidative aging. Because no single FF rigorously covers this heterogeneity, FF requirements are derived at the intersection of two axes: a usage-environment axis (binder multi-component nature, mineral diversity, water and ions, oxidative aging) and a theory-and-methodology axis (bonded and non-bonded energy terms, mixing rules, 1-4 scaling compatibility, and long-range electrostatic management). Part I establishes the theory-and-methodology axis by systematizing the decomposition of bonded and non-bonded MD energy, the partial-charge models and reusable FF artifacts (AM1-BCC, RESP, CHELPG, CM1A/CM5, NPA), the non-bonded summation methods (cutoff, Ewald/PPPM, slab correction, and QEq), and the mixing-rule and 1-4 scaling conditions that govern multi-FF stacking. On this basis, Part I provides a comparative evaluation of the major candidate force fields (GAFF/GAFF2 + AM1-BCC, OPLS-AA/OPLS4, COMPASS/COMPASS III, AMBER ff14SB, CHARMM CGenFF, INTERFACE FF, and ReaxFF) in terms of functional form, charge scheme, mixing-rule compatibility, and primary application. The companion Part II derives the asphalt-specific FF requirements, operationalizes the resulting 3+1 lane partitioning (organic binder, mineral interface, water and ions, and reactive oxidation) into a scenario-based operating matrix, and validates the Lane 1 baseline through parity-plot comparisons of GAFF2 against commercial COMPASS III.

Keywords

Asphalt molecular dynamics

Multi-force-field stack

Force-field comparison

Partial-charge model and FF artifact

Mixing-rule and 1-4 scaling compatibility

MAIN

1. 서 론
2. MD 에너지 구성과 FF 아티팩트의 이론적 배경
2.1 결합 및 비결합 에너지의 분해
2.2 부분 전하 모델과 FF 아티팩트
2.3 비결합 에너지 계산 방법론
2.4 조합 소재 대상 Multi-FF 필요성
2.5 혼합 규칙, 1-4 scaling과 FF 스택
3. FF 비교
3.1 GAFF / GAFF2 + AM1-BCC
3.2 OPLS-AA / OPLS4
3.3 COMPASS / COMPASS III
3.4 AMBER ff14SB
3.5 CHARMM CGenFF
3.6 INTERFACE Force Field (IFF)
3.7 ReaxFF
4. 결 론

1. 서 론

아스팔트 바인더는 SARA 4성분(saturate, aromatic, resin, asphaltene)의 복잡한 탄화수소 혼합물(Li and Greenfield, 2014)에 고분자 화학 첨가제, 광물 골재 표면, 그리고 도로 환경의 수분･이온까지 포함하는 다성분･다상(multi-component, multi-phase) 시스템이다. 아스팔트 바인더의 물성을 예측하기 위한 분자동역학(molecular dynamics, MD) 활용 시뮬레이션 연구는 Zhang and Greenfield(2007) 이후로 활발히 수행되어 왔지만, 특정 첨가제와 광물을 고려한 아스팔트 바인더 시스템에 대한 MD 결과의 신뢰성은 힘장(force field, FF) 선택에 결정적으로 의존하기 때문에, 단일 FF로 이러한 화학적･물리적 이질성 전체를 엄밀하게 다루는 것은 원리적으로 불가능하다.

이와 같은 한계를 극복하기 위하여 화학이나 소재과학 분야에서는 유기물(organic)･광물(mineral)･물(water)/이온(ions)･반응성(reactive) 등 서로 다른 방법 또는 Lane에 서로 다른 FF를 매칭한 다중 FF 스택(multi-FF stack)이 사실상 표준이 되고 있지만, 그 설계 원리와 검증 절차가 학술적으로 일반화된 프레임워크로 정리된 사례는 흔하지 않다. 특히 다양한 종류의 복잡한 소재를 사용하는 아스팔트 혼합물에 대해서는 이러한 프레임워크가 체계적으로 제안된 바가 없다.

본 연구에서는 다양한 종류의 첨가물과 같이 사용되면서도 골재와 같은 광물에 부착되어 수분이나 나트륨과 같은 이온에 노출되는 아스팔트 혼합물의 물성을 MD로 예측하는 데에 활용될 수 있는, 기존의 엄밀한 이론과 최근 연구 논문에 근거한 다중 FF 스택의 설계 및 검증 프레임워크를 제안하고자 하였다. 본 프레임워크는 아스팔트 바인더에 대한 일반 MD 연구 전반에도 적용할 수 있지만, 이른바 소재 역설계와 같은 닫힌-루프 자율 시뮬레이션(closed-loop autonomous simulation) 파이프라인 구축을 통하여 일관성 있는 결과를 얻기 위해서도 필요성이 높다.

본 연구는 내용의 분량과 논리적 완결성을 고려하여 두 편(Part I, Part II)으로 나누어 제시한다. 제I편에 해당하는 본 논문은 다중 FF 스택을 구성하기 위한 이론적 기반과 범용 역장의 비교에 초점을 둔다. 즉, 분자동역학 포텐셜 에너지의 결합･비결합 항 분해, 부분 전하 모델과 재사용 가능한 FF 아티팩트, 비결합 에너지의 수치 계산 방법론, 그리고 서로 다른 FF를 동일 시뮬레이션에서 결합하기 위한 혼합 규칙과 1-4 scaling 호환성 등 이론･방법론 축을 체계화하고(제2장), 이를 토대로 아스팔트･탄화수소 시뮬레이션에 활용 가능한 주요 역장들을 함수 형태･전하 모델･혼합 규칙 호환성･주요 적용 분야의 관점에서 비교한다(제3장). 이어지는 제II편에서는 본 편에서 정리한 이론･방법론 축을 아스팔트 혼합물의 사용 환경 축과 교차하여 FF 요구사항을 도출하고, 그 결과로 제안되는 3+1 Lane 분담 구조를 시나리오별 운영 매트릭스로 구체화하며, Lane 1 베이스라인을 상용 COMPASS III와의 parity-plot 비교로 검증한다. 본 논문(Part I)의 구성은 다음과 같다. 제2장에서는 MD 에너지 구성과 FF 아티팩트의 이론적 배경을, 제3장에서는 주요 FF의 비교를 제시하고, 제4장에서 본 편의 결론을 도출한다.

2. MD 에너지 구성과 FF 아티팩트의 이론적 배경

우선 제2장에서는 여러 연구에서 단편적으로 활용되는 FF 이론 요소를 단일 프레임워크로 정리하였다. 본 연구에서 아스팔트 혼합물을 위하여 특정 FF 조합 아티팩트(artifact)의 생성･재사용을 제안하는 이유, 다수의 비결합 방법론이 동시에 필요한 물리적 근거, 그리고 유기물･광물･이온 트랙이 단일 FF로 통합될 수 없는 이론적 근거를 제시하였다.

2.1 결합 및 비결합 에너지의 분해

MD 시뮬레이션에서 시스템의 전체 포텐셜 에너지(E_potential)는 결합에너지(E_bonded)와 비결합에너지(E_non-bonded)로 분해될 수 있다. 결합 에너지는 공유결합 토폴로지에 의해 정의되어 결합으로 묶인 원자 쌍･삼중･사중에만 적용되며, 결합 신장(bond stretching), 결합각 변형(angle bending), 이면각 회전(dihedral torsion), 평면성(improper torsion) 등의 항으로 구성되어 각각 식 (1), (2), (3), (4)와 같이 표현된다.

(1)

E_{b o n d} = \sum k_{b} (r - r_{0})^{2}

(2)

E_{a n g l e} = \sum k_{θ} (θ - θ_{0})^{2}

(3)

E_{d i h e d r a l} = \sum V_{n} [1 + c o s (n φ - δ)]

(4)

E_{i m p r o p e r} = \sum k_{χ} (χ - χ_{0})^{2}

여기에서 파라미터들(k_b, r₀, k_θ, θ₀, V_n)은 원자들의 종류 및 조합특성에 따라 결정된다. 비결합 에너지는 모든 원자 쌍에 적용되어 결합 토폴로지와 무관하며, 두 가지 다른 물리적 기원을 갖는다. 비결합 에너지 중 van der Waals(Lennard-Jones 12-6) 항은 식 (5)와 같이 표현된다.

(5)

E_{L J} = 4 ε [(\frac{σ}{r})^{12} - (\frac{σ}{r})^{6}]

여기에서 (σ/r)¹² 반발항은 전자 구름 중첩 시 Pauli 배타 원리에 의한 양자역학적 반발력을 나타내고, (σ/r)⁶ 인력항은 London 분산력으로 순간 쌍극자-유도 쌍극자 상호작용을 나타낸다. ε는 포텐셜 우물 깊이 또는 상호작용 세기, σ는 충돌 직경 또는 원자 크기를 의미하며, r^-6 감쇠로 인해 통상 12 Å의 cutoff로 충분히 처리할 수 있는 단거리력이다. 또다른 비결합 에너지 중 Coulomb 항은 식 (6)과 같이 표현된다.

(6)

E_{c o u l} = \sum_{i < j} \frac{q_{i} q_{j}}{4 π ε_{0} r_{i j}}

이는 부분 전하 q_i의 곱에 비례하며 r^-1로 감쇠하는 장거리력이다. Cutoff만으로는 충분히 수렴하지 않는 대표적인 장거리력인 r^-1 감쇠율은 다음 장에서 구체적으로 설명된다.

2.2 부분 전하 모델과 FF 아티팩트

비결합 에너지인 장거리력은 아스팔트 바인더와 골재사이의 부착 거동을 지배하는 핵심인자이며, 이러한 장거리 전자기적 상호작용은 연속적 전자밀도 함수 $ρ (r)$ 에 의하여 결정되는 전하 분포에 의하여 결정된다. 전통적인 MD에서는 이를 원자핵 위치에 존재하는 점전하(point charge) $q_{i}$ 로 근사하며. 이 근사가 부분 전하이다. 그러나 부분 전하는 관측 가능한 물리량이라기보다 전자밀도 분포를 근사적으로 표현하기 위한 모델 의존적 파라미터에 가깝다. 따라서 부분 전하를 얻기 위해서는 상대적으로 비용이 큰 양자화학 계산 또는 이에 준하는 전하 할당(charge assignment) 절차가 필요하다. 한편 이러한 계산을 반복하는 것은 비효율적이므로 계산 결과를 저장해놓고 사용가능한 행태의 FF 아티팩트의 관리하는 것이 합리적이다.

장거리력을 기술하는 대표적인 부분 전하 모델로는 AM1-BCC 전하 모델(Jakalian et al., 2002)이 널리 사용된다. 이 모델을 활용하는 방법은 크게 두 단계로 구성되는데, 첫 번째 단계는 AM1(Austin Model 1; Dewar et al., 1985)이라는 NDDO 근사 기반 반경험적 분자 궤도 계산이다. 이 단계에서는 Hartree-Fock의 경험적 근사방법에 근거하여 자체 일관장(self-consistent field, SCF)으로 전자 구조를 풀어 멀리켄 점유(Mulliken population) 전하를 얻는다. 이것이 FF 파라미터를 설정하는 소프트웨어로 알려진 Antechamber가 내부적으로 반경험 양자역학(semiempirical quantum mechanics, SQM)을 호출하는 이유이다. 일반적으로 SCF는 분자당 수십~수백 초가 소요되고 수렴 실패 가능성이 존재하므로 수렴기준을 완화하고 초기 추측 변경이 가능한 강건한(robust) 프로파일을 활용하기도 한다. 두 번째 단계는 결합 전하 보정(bond charge correction, BCC)인데, 앞서 결정한 AM1 멀리켄 점유 전하가 갖는 오차를 결합 종류별 경험 보정 파라미터로 수정하는 과정이다. 여기에서 FF 아티팩트에 저장되는 핵심 데이터는 (i) 양자화학 계산 결과로 시뮬레이션 중 재계산 불가한 부분 전하 $q_{i}$ , (ii) Antechamber의 GAFF2의 화학 환경 인지(chemical-environment-aware) 원자 타이핑 결과인 원자 타입 $f f_{t y p e}$ , (iii) 원자 타입에 의해 결정되는 LJ 파라미터 $ε \cdot σ$ , (iv) 원자 타입 쌍･삼중･사중에 의해 결정되는 결합 파라미터(k, $r_{0}$ 등) 이다. Antechamber + SQM 실행에는 분자당 수백 초가 소요되므로, 한 번 계산하여 JSON 등 구조화된 아티팩트로 저장하고 시뮬레이션마다 재사용하는 것이 계산 시간 뿐만 아니라 계산의 연속성을 확보할 수 있는 방법이라 할 수 있다. ABCG2는 GAFF2의 BCC 파라미터를 개선한 후속 모델(He et al., 2020)이다. 이외에도 부분 전하를 할당하는 방법으로 분자 표면의 정전기 퍼텐셜(electrostatic potential, ESP)에 직접 피팅하여 형태 이성질체 앙상블(conformer ensemble)을 샘플링함으로써 더 높은 정확도를 확보하는 것으로 알려져 있으며, 결과의 비결정론적 성격으로 자동화 친화도가 낮은 구속 정전기 퍼텐셜(Restrained Electrostatic Potential, RESP; Bayly et al., 1993), 분자 표면 격자점에서의 ESP 피팅을 단순화한 변형인 CHELPG(CHarges from Electrostatic Potentials using a Grid-based method; Breneman and Wiberg, 1990), OPLS 계열에서 채택되는 전하 모델(charge model)로 CM1A･CM5(Charge Model 1A･5; Storer et al., 1995; Marenich et al., 2012), 자연 결합 궤도(natural bond orbital, NBO) 기반의 대안 분석법으로 자연 점유 분석(Natural Population Analysis, NPA; Reed et al., 1985)이 있다.

2.3 비결합 에너지 계산 방법론

비결합 에너지 계산에 다양한 방법론이 동시에 필요한 근본 원인은 2.1절에서 언급된 LJ( $r^{- 6}$ )와 Coulomb( $r^{- 1}$ )의 감쇠 비대칭 때문이다. LJ는 빠른 $r^{- 6}$ 감쇠로 cutoff 12 Å에서 충분히 수렴하지만, 3차원 주기적 시스템에서 Coulomb 합 $\sum \frac{q_{i} q_{j}}{r_{i j}}$ 는 조건부로 수렴한다. 이는 합산 순서에 따라 결과가 달라지는 것을 의미하므로 cutoff만으로는 정확히 계산할 수 없고, 시뮬레이션 셀의 기하학(non-periodic, 3D periodic, 2D periodic slab, reactive)에 따라 서로 다른 수학적 처리가 필요하다. 앞선 방법으로 부분전하가 결정되면 시뮬레이션 도중 부분전하를 동적으로 합산하는데 사용될 수 있는 방법은 다음과 같다.

(가) Cutoff 단순 절단(coul/cut): r > $r_{c u t}$ 일 때 E_coul = 0으로 처리한다. 계산이 빠르지만 Coulomb의 장거리 기여를 무시하므로 벌크 물성 평가에는 부정확하다. 주기적 이미지가 없는 가스상 단일 분자에 한해 합리적이므로, 아스팔트 바인더를 구성하는 개별 분자에 대한 계산을 위한 단일 분자 진공 모드(single-molecule vacuum mode)에만 적용해야 한다.

(나) Ewald 분할 → PPPM(coul/long + pppm). Ewald summation의 핵심 아이디어는 장거리 합을 두 개의 빠르게 수렴하는 합으로 분할하는 방법이다: E_coul = E_short-range(실공간) + E_long-range(역공간). 실공간 항은 차폐된(screened) Coulomb $\frac{e r f c (α r)}{r}$ 로 단거리에서 빠르게 감쇠하고, 역공간 항은 부드러운 장거리 부분으로 Fourier 변환을 통해 역공간에서 효율적으로 계산된다. PPPM(Particle-Particle Particle-Mesh)은 역공간 항을 격자에서 FFT로 계산하는 효율적 변형으로, 순수 Ewald의 $O (N^{3} / 2)$ 대비 개선되는 $O (N l o g N)$ 스케일링을 제공한다. 모든 3D 주기적 벌크･계면 시스템에서 표준적으로 적용할 수 있으므로, 아스팔트 바인더에 우선 적용을 고려할 수 있다.

(다) Slab 보정(kspace_modify slab). x･y는 주기적이지만 z는 비주기적인 2D 주기적 시스템(유기물-광물 계면 등)에서, 표준 3D PPPM은 z 방향 이미지 전하를 비물리적으로 포함한다. Yeh and Berkowitz(1999)는 이 문제를 해결하기 위해 z 방향 진공 갭을 약 3배 정도로 충분히 확장하고 표준 3D Ewald 식에 쌍극자(dipole) 보정 항을 추가하는 방법을 제안했다. 두 점전하 사이 쌍힘(pair force) 비교에서 엄밀한 2D Ewald 결과와 매우 잘 일치하면서도 계산 시간 단축이 가능함을 보였는데, 아스팔트 바인더와 골재사이에서 생성되는 유기물-광물 계면 lane에서 이 slab 보정이 필수적이라고 할 수 있다.

(라) QEq 동적 전하(qeq/reax). ReaxFF는 고정 전하를 사용하지 않고, 매 시간 간격(timestep)마다 전기음성도 균등화(electronegativity equalization)으로 $q_{i}$ 를 재계산한다. 이는 결합 형성･절단 시 전하 재분배를 허용하기 위한 필수 절차이지만, 매 스텝의 추가 반복 계산 비용이 크게 증가하므로 ReaxFF를 검증 lane으로 한정하거나 제한된 반응에만 적용하는 것이 바람직하다. 결정구조의 고정(freeze) 가정도 IFF의 설계 범위와 직결된 실용적인 선택으로 고려될 수 있다. IFF는 광물-유기물 계면 에너지에 피팅되었지만, 광물 벌크 포논(phonon)을 정확히 재현하도록 설계되지 않았다. 따라서 LAYER_BULKFF 시뮬레이션에서 광물 원자를 fix spring/self + rigid freeze로 고정하여 정확도를 확보할 수 있다. 이로써 유기물-광물 계면 접착 에너지･접촉각･흡착 구조가 정확히 재현될 수 있으며, 포논-흡착 결합(phonon-adsorbate coupling), 표면 재구성, 열전도는 정확하게 재현되지 않을 수 있다. 광물 내부 동역학까지 필요한 응용에서는 CLAYFF 등 광물 전용 FF로 광물도 동적으로 다루어야 하나, 아스팔트 바인더 물성 계산에서는 표면 접착이 지배적이므로 이 방법이 바람직하다.

2.4 조합 소재 대상 Multi-FF 필요성

유기 분자는 명확한 공유결합 토폴로지를 가지며 결합이 생성･절단되지 않는다는 MD 가정 하에서 AM1-BCC가 충분히 정확한 부분 전하를 제공하므로 GAFF2 + AM1-BCC가 아스팔트 바인더와 같은 범용 유기물 FF로 적합하다고 할 수 있다. 반면 골재의 주요 성분이 되는 광물 표면( $S i O_{2}$ , $A l_{2} O_{3}$ 등)은 이온성 결합과 공유 결합이 혼합된 시스템이며, 계면에서의 흡착･접착이 핵심 물성인 반면 범용 FF(UFF; Rappé et al., 1992)는 벌크 구조에 피팅되어 계면 에너지가 부정확하다. INTERFACE FF(IFF)(Heinz et al., 2013)는 광물의 LJ 파라미터를 표면 에너지･접촉각 등 거시 규모(macro-scale) 측정값에 열역학적으로 일관성있게 피팅한 점에 있는데, 실리카 표면의 유기물 흡착이 정전기에 의해 지배되도록 설계된다. 이온( $N a^{+}$ , $K^{+}$ , $C a^{2 +}$ , $C l^{-}$ 등)은 정수 형식 전하(formal charge)를 가지므로 AM1-BCC가 아닌 형식 전하를 그대로 사용하지만, 이온의 비결합 에너지에서 Coulomb이 지배적이더라도 LJ 파라미터가 수화 자유에너지(hydration free energy)와 격자 상수(lattice constant)를 동시에 결정한다. Joung and Cheatham(2008)은 TIP3P/SPC/E/TIP4P-Ew 물 모델과 호환되도록 9개 단원자 이온의 LJ $ε \cdot σ$ 를 Lennard-Jones 공간 전체에서 스캔하면서 수화 자유에너지(hydration free energy) + 격자에너지/상수(lattice energy/constant) + 제1 수화각 반경(first-shell radius)을 동시에 균형 있게 재현하도록 피팅한 표준 이온 파라미터 세트를 제시했는데, FF(UFF, GAFF2)의 이온 파라미터로 수화 구조를 정확히 재현하기 어려운 경우 활용될 수 있다. 이온은 도로에서 결빙 방지나 차량 운행 중 발생되는 화학물질들을 고려할 때 필요하다. Fig. 1은 반응성까지 포함하여 호환성을 우선으로 고려하는 본 프레임워크의 4-Lane 구조도를 나타낸다.

https://cdn.apub.kr/journalsite/sites/jkai/2026-016-01/N0850160106/images/jkai_2026_161_65_F1.jpg

Fig. 1.

4-Lane structure of the multi-force-field stack

2.5 혼합 규칙, 1-4 scaling과 FF 스택

서로 다른 FF의 원자가 만나는, 예를 들어 유기물 C ↔ 광물 Si의 이종 FF(cross-FF) 인터페이스에서 LJ 파라미터는 혼합 규칙(mixing rule)으로 조합된다. Lorentz-Berthelot(L-B) 규칙은 식 (7)의 산술 평균과 식 (8)의 기하 평균으로 정의된다.

(7)

σ_{i j} = \frac{σ_{i} + σ_{j}}{2}

(8)

ε_{i j} = \sqrt{ε_{i} \cdot ε_{j}}

IFF는 L-B mixing을 가정하고 피팅되었으므로 IFF host로 AMBER･CHARMM 변종을 채택하는 유기물-광물 계면 lane에서는 반드시 산술평균 혼합(arithmetic mixing)을 강제해야 한다. 반면에 OPLS-AA가 기하평균 혼합(geometric mixing)을 채택한다는 점이 GAFF･AMBER와의 cross-FF 사용을 어렵게 만드는 이유이며, 이 호환성 문제가 닫힌 계(closed system) 적용 원리에 반영되어야 한다.

결합으로 가까이 연결된 원자 쌍의 비결합 상호작용은 스케일링된다. 1-2 결합(직접 결합) 및 1-3(각도) 거리는 이미 결합 항에 포함되어 LJ･Coulomb 모두 0으로 제외되며, 1-4(이면각) 거리는 AMBER(Cornell et al., 1995) 계열 표준으로 LJ를 0.5, Coulomb을 5/6 ≈ 0.8333로 스케일링하며, 1-5 이상은 전체(1.0)를 적용한다. 이 1-4 scaling 값은 AMBER 계열 FF의 관례이며, GAFF2의 이면각(dihedral) 파라미터가 이 값을 가정하고 피팅되었기 때문에 GAFF2 사용 시 0.5/0.8333 외 값을 임의로 적용하지 말아야 한다. OPLS는 LJ･Coulomb 모두 0.5를 사용하므로 이면각 파라미터의 매칭이 다른 점도 고려되어야 한다.

LAMMPS thermo 출력으로 추적되는 에너지 항은 E_bond + E_angle + E_edihed + E_imp로써 이들은 전체 결합에너지 E_mol를 구성하며, E_vdwl + E_coul + E_long은 비결합 에너지 E_pair를 구성한다. 본 논문에서 제안된 프레임워크에서는 이 모든 항들이 아티팩트 메타데이터에 기록되어 활용된다.

3. FF 비교

3.1 GAFF / GAFF2 + AM1-BCC

GAFF는 Wang et al.(2004)이 일반 유기 소분자(약물･리간드･작은 유기물)를 다루기 위해 개발한 AMBER 호환 FF이다. 원자 타이핑(atom typing)은 혼성화(hybridization)･방향족성(aromaticity)･EWG 인접성 기반으로 정의되며, 결합항은 조화함수, 이면각은 Fourier 형태이다. GAFF2의 BCC 파라미터는 He et al.(2020)에서 442개 유기 용질(organic solute)의 수용액 SFE에 피팅하여 최적화되었으며, 이 새로운 전하모형은 ABCG2(AM1-BCC-GAFF2)로 정의되었다. GAFF의 한계는 Caleman et al.(2012)의 액상 유기물(organic liquid)에 대한 벤치마크에서 확인된 바 있다. 146개 액상 유기물과 1,200개 이상의 실험 측정값에 대해 GAFF･OPLS-AA는 밀도를 일관적으로 약간 과소평가(slope a = 0.96, $R^{2}$ = 97%)하며 $Δ H_{v a p}$ 을 과대평가(slope a = 1.07)한 것으로 확인되었다. 또한 고정 점전하(fixed-point-charge) 모델은 분극(polarization) 효과 부재로 다양한 용매에서 액체물성과 용매화 자유에너지를 동시에 재현하기 어렵다는 한계가 있기도 하다.

AM1-BCC 전하 모델은 Jakalian et al.(2000, 2002)이 제시한 방법으로, 2,755 유기분자의 훈련 데이터셋을 통해 309개의 unique BCC 파라미터를 도출했고, AMBER + AM1-BCC 조합은 수소결합 이량체(hydrogen-bond dimer) 에너지를 ab initio 대비 0.95 kcal/mol RMS, monofunctional isostere의 상대 SFE를 실험 대비 0.69 kcal/mol 이내로 재현했으며, 검증 전반에서 상관계수는 0.96 이상을 유지한 것으로 알려져 있다.

3.2 OPLS-AA / OPLS4

OPLS-AA는 Jorgensen et al.(1996)가 개발한 액상 유기체 시뮬레이션용 전 원자(all-atom) FF이다. 결합･결합각은 AMBER에서 거의 그대로 가져왔고, 이면각은 50개 이상의 organic 분자･이온의 RHF/6-31G*//RHF/6-31G* 회전 에너지 프로파일에 피팅하여 평균 형태 에너지(conformational energy) 오차 0.2 kcal/mol 이하로 매우 우수한 것으로 알려져 있다. 비결합 파라미터는 Monte Carlo 시뮬레이션과 결합하여 알칸(alkane), 알케인(alkene), 알코올(alcohol), 에테르(ether), 아세탈(acetal), 티올(thiol), 설파이드(sulfide), 아민(amine) 등 34개 유기 액체의 열역학적･구조적 특성에 맞춰 개발되었다. OPLS의 혼합 규칙은 기하평균이기 때문에 산술평균을 사용하는 Lorentz-Berthelot인 GAFF/AMBER과 호환되지 않으므로 이종 FF 사용 시 주의가 필요하다. 부분 전하와 LJ 파라미터에 대해서 비교적 정확한 값을 제시하는 LigParGen(Dodda et al., 2017)은 SMILES/MOL/PDB 입력에 대해 OPLS-AA/1.14*CM1A(-LBCC) 파라미터를 자동 생성하기 때문에 활용이 편리한 장점이 있다. 황 함유 작용기에 대한 정확도를 개선한 상용 버전으로 OPLS4(Lu et al., 2021)가 있다.

3.3 COMPASS / COMPASS III

COMPASS(Sun, 1998)는 ab initio 기반으로 응축상(condensed-phase) 응용을 최적화한 class II FF이다. 28개 분자류･격리 분자･액상분자･결정구조에 대한 검증에서 구조･형태･진동수･열역학적 특성을 동시에 정확히 예측한다고 보고되었다. COMPASS의 교차항(cross-term)과 9-6 LJ form은 GAFF/OPLS의 Class I form과 호환되지 않으므로 이종 FF 사용이 어렵다. COMPASS/PCFF는 BIOVIA(Materials Studio) 의존성이 강하기 때문에 활용에 제약이 있다.

3.4 AMBER ff14SB

AMBER ff14SB는 단백질의 주사슬(backbone)･곁사슬(side chain) 파라미터를 NMR과 QM 데이터에 더 정밀하게 피팅한 표준 단백질 FF이다. 아스팔트 응용에서는 작은 organic 분자가 시뮬레이션 대상의 대부분이므로 GAFF/GAFF2가 소분자(small molecule)를 처리하고 ff14SB는 직접 사용되지 않지만, 두 FF는 동일한 함수 형태를 갖으므로 호환이 가능하다.

3.5 CHARMM CGenFF

CGenFF(Vanommeslaeghe et al., 2010)는 CHARMM36 호환 범용(generic) FF로, 원자 타이핑･유사 기반 파라미터 할당･벌점 점수를 통해 신뢰성을 긍정적으로 평가한 바 있다.

3.6 INTERFACE Force Field (IFF)

Heinz et al.(2013) 가 제시한 IFF는 무기-유기 계면 시뮬레이션을 위한 열역학적으로 일관성이 확보된 FF이다. IFF는 원자 척도에서의 물리화학적 이해를 바탕으로 모든 파라미터(특히 원자 전하(atomic charge)와 반데르발스 상수(vdW constant))를 할당하며, 각 화합물의 거시 척도 물성을 측정값과 비교하여 검증한다. 이 접근법은 다른 파라미터화 프로토콜(protocol)에서 발생하는 벌크･표면 물성(bulk･surface property)의 최대 두 자리수 오차를 제거하고, 열역학적으로 일관성을 도입함으로써 파라미터 전이성(transferability)을 향상시킨다.

IFF의 핵심적 강점은 호스트 FF(host FF; PCFF, COMPASS, CHARMM, AMBER, GROMACS, OPLS-AA)와 함수 형태･혼합 규칙을 공유하여 organic-inorganic 시뮬레이션을 단일 프레임워크에서 수행할 수 있다는 점이다. Mishra et al.(2021)은 무수석고(anhydrite), 반수석고(bassanite), 석고(gypsum)를 포함하는 황산칼슘(calcium sulfate) 에 대해 1~5% 정확도로 격자(lattice)･표면･수화･기계적･열적물성을 성공적으로 재현한 바 있다. 나아가 Kanhaiya et al.(2023)은 IFF를 산화물(oxide)･수산화물(hydroxide) 및 유기 하이브리드 소재로 확장하여 마이크로미터 규모까지 구조･계면 물성을 일관되게 재현할 수 있음을 보고하였다. Fig. 2는 5종 유기물 FF와 5종 IFF 사이의 혼합 규칙 호환성을 나타내는데, 예를 들어 그림 아래 설명된 바와 같이 Closed L-B match는 적정한 조합이 될 수 있는 조합을 나타낸다.

https://cdn.apub.kr/journalsite/sites/jkai/2026-016-01/N0850160106/images/jkai_2026_161_65_F2.jpg

Fig. 2.

Mixing-rule compatibility between 5 organic FFs and 5 IFF host variants

3.7 ReaxFF

van Duin et al.(2001)이 탄화수소(hydrocarbon)용으로 발표한 ReaxFF는 결합 차수(bond order)에 의존하는 에너지 함수와 전하 평형화(charge equilibration, QEq)을 결합하여 결합 형성･파괴를 시뮬레이션 도중 자연스럽게 다룰 수 있게 한다. Chenoweth et al.(2008)은 탄화수소 산화(hydrocarbon oxidation)을 위해 C/H/O 파라미터를 추가했고, 이후 N/S/금속 등으로 확장되었다. Senftle et al.(2016)은 ReaxFF의 개발과 적용을 종합적으로 검토했는데, ReaxFF는 (i) 매우 짧은 시간 간격에 의한 큰 계산 비용, (ii) 파라미터 셋 사이의 전이성 부족이 제약이라 할 수 있다.

Table 1은 아스팔트/탄화수소 시뮬레이션에서 고려할 만한 주요 FF의 특성을 나타낸다. 함수 형태 컬럼에는 비결합 형태, 혼합 규칙이 표기되어 있다.

Table 1.

Comparison summary of the major force fields

FF	Functional form	Charge / Param	Primary application
GAFF / GAFF2	Harmonic + Fourier; 12-6 LJ; L-B mixing; 1.0 fs	AM1-BCC	Organic, hydrocarbon, SARA
OPLS-AA / OPLS4	AMBER form, 4-term Fourier; geometric mixing	CM1A / OPLS4	Organic liquid, drug
COMPASS III	Class II, 4-term anharmonic, 9-6 LJ; 6th-power	ab initio	Polymer, hydrocarbon
PCFF	Class II (isomorphic with COMPASS)	ab initio	Polymer, some inorganics
CVFF	Class I+ (with Morse option)	Empirical	Peptide, IFF host
AMBER ff14SB	Isomorphic with GAFF	RESP	Proteins & nucleic acids
CHARMM CGenFF	Harmonic + Urey-Bradley + 12-6 LJ	Analogy + penalty	Drug, lipid
MMFF94	Buffered 14-7 vdW	QM-fit	Drug docking
UFF	Harmonic + 12-6 LJ; generic	Generic	Cross-domain estimate only
DREIDING	Harmonic, generic	Empirical	Polymer (deprecated)
INTERFACE FF	Host FF compatible (PCFF/AMBER/CHARMM)	Born-extended	Minerals, metals, oxides, cement
ReaxFF	Bond-order + QEq;	Dynamic charge	Reactive (oxidation, combustion)

4. 결 론

제I편은 아스팔트 바인더･골재 계면의 분자동역학(MD) 시뮬레이션을 위한 다중 역장(FF) 스택을 설계하기에 앞서, 이를 뒷받침하는 이론･방법론 축을 체계화하고 주요 범용 역장을 비교하는 것을 목적으로 수행되었다. 이를 위해 포텐셜 에너지를 결합 항(결합 신축･각도 굽힘･이면각･부적절 이면각)과 비결합 항(반데르발스･정전기)으로 분해하고, 부분 전하가 관측 가능한 물리량이 아닌 모델 의존적 파라미터임을 근거로 양자화학 계산 결과를 재사용 가능한 FF 아티팩트로 관리하는 전략을 정리하였다. 또한 LJ(r^-6)와 Coulomb(r^-1)의 감쇠 비대칭에서 비롯되는 비결합 에너지 계산 방법론(cutoff, Ewald/PPPM, slab 보정, QEq)을 시뮬레이션 셀의 기하학과 연계하여 구분하였다.

주요 결과는 다음과 같다. 첫째, 서로 다른 FF를 동일 시뮬레이션에서 결합하기 위한 호환 조건으로 혼합 규칙(Lorentz-Berthelot 대 기하평균)과 1-4 scaling 체계의 일치가 핵심임을 정리하고, GAFF2/AMBER 계열과 INTERFACE-AMBER가 함수 형태･혼합 규칙･1-4 scaling을 공유하여 닫힌 계에서 결합 가능한 호환군을 형성함을 확인하였다. 둘째, 아스팔트･탄화수소 시뮬레이션에 활용 가능한 주요 역장(GAFF/GAFF2 + AM1-BCC, OPLS-AA/OPLS4, COMPASS/COMPASS III, AMBER ff14SB, CHARMM CGenFF, INTERFACE FF, ReaxFF)을 함수 형태･전하 모델･혼합 규칙 호환성･주요 적용 분야의 관점에서 비교하여(Table 1), 유기상･무기 계면･반응성 영역에 각각 적합한 후보군과 그 호환성･제약을 명확히 하였다.

본 편은 개별 FF의 장단점을 단순 나열하는 기존 리뷰와 달리, 다중 FF가 동일 시뮬레이션 내에서 호환･공존하기 위한 이론적 조건과 호환군을 명시함으로써 다중 FF 스택 설계의 이론적 토대를 제공한다. 이를 바탕으로 제II편에서는 본 편에서 정리한 이론･방법론 축을 아스팔트 혼합물의 사용 환경 축(바인더 다성분성･골재 광물 다양성･수분･이온･산화 노화)과 교차하여 FF 요구사항을 도출하고, 그 결과로 제안되는 3+1 Lane 분담 구조를 시나리오별 운영 매트릭스로 구체화하며, Lane 1 베이스라인(GAFF2 + AM1-BCC)을 상용 COMPASS III와의 parity-plot 비교로 검증한다.

Acknowledgements

본 연구는 한국건설기술연구원의 2026년 주요사업(과제번호: 20260120-001 분자동역학과 소재정보학을 활용한 도로 건설 재료 배합기술 개발 연구)의 재원으로 수행된 연구 결과입니다.

References

Bayly, C.I., Cieplak, P., Cornell, W.D. and Kollman, P.A. (1993). “A well-behaved electrostatic potential based method using charge restraints for deriving atomic charges: The RESP model”, J. Phys. Chem., 97, pp. 10269-10280.

10.1021/j100142a004

Breneman, C.M. and Wiberg, K.B. (1990). “Determining atom-centered monopoles from molecular electrostatic potentials. The need for high sampling density in formamide conformational analysis”, J. Comput. Chem., 11, pp. 361-373.

10.1002/jcc.540110311

Caleman, C., van Maaren, P.J., Hong, M., Hub, J.S., Costa, L.T. and van der Spoel, D. (2012). “Force field benchmark of organic liquids”, J. Chem. Theory Comput., 8, pp. 61-74.

10.1021/ct200731v22241968PMC3254193

Chenoweth, K., van Duin, A.C.T. and Goddard, W.A. (2008). “ReaxFF reactive force field for molecular dynamics simulations of hydrocarbon oxidation”, J. Phys. Chem. A, 112, pp. 1040-1053.

10.1021/jp709896w

Cornell, W.D., Cieplak, P., Bayly, C.I., Gould, I.R., Merz, K.M., Ferguson, D.M., Spellmeyer, D.C., Fox, T., Caldwell, J.W. and Kollman, P.A. (1995). “A second generation force field for the simulation of proteins, nucleic acids, and organic molecules”, J. Am. Chem. Soc., 117, pp. 5179-5197.

10.1021/ja00124a002

Dewar, M.J.S., Zoebisch, E.G., Healy, E.F. and Stewart, J.J.P. (1985). “Development and use of quantum mechanical molecular models. 76. AM1: A new general purpose quantum mechanical molecular model”, J. Am. Chem. Soc., 107, pp. 3902-3909.

10.1021/ja00299a024

Dodda, L.S., Cabeza de Vaca, I., Tirado-Rives, J. and Jorgensen, W.L. (2017). “LigParGen web server: an automatic OPLS-AA parameter generator for organic ligands”, Nucleic Acids Res., 45, pp. W331-W336.

10.1093/nar/gkx31228444340PMC5793816

He, X., Man, V.H., Yang, W., Lee, T.-S. and Wang, J. (2020). “A fast and high-quality charge model for the next generation general AMBER force field”, J. Chem. Phys., 153, 114502.

10.1063/5.001905632962378PMC7728379

Heinz, H., Lin, T.J., Mishra, R.K. and Emami, F.S. (2013). “Thermodynamically consistent force fields for the assembly of inorganic, organic, and biological nanostructures: The INTERFACE force field”, Langmuir, 29, pp. 1754-1765.

10.1021/la3038846

Jakalian, A., Bush, B.L., Jack, D.B. and Bayly, C.I. (2000). “Fast, efficient generation of high-quality atomic charges. AM1-BCC model: I. Method”, J. Comput. Chem., 21, pp. 132-146.

10.1002/(SICI)1096-987X(20000130)21:2<132::AID-JCC5>3.0.CO;2-P

Jakalian, A., Jack, D.B. and Bayly, C.I. (2002). “Fast, efficient generation of high-quality atomic charges. AM1-BCC model: II. Parameterization and validation”, J. Comput. Chem., 23, pp. 1623-1641.

10.1002/jcc.10128

Jorgensen, W.L., Maxwell, D.S. and Tirado-Rives, J. (1996). “Development and testing of the OPLS all-atom force field on conformational energetics and properties of organic liquids”, J. Am. Chem. Soc., 118, pp. 11225-11236.

10.1021/ja9621760

Joung, I.S. and Cheatham, T.E. III (2008). “Determination of alkali and halide monovalent ion parameters for use in explicitly solvated biomolecular simulations”, J. Phys. Chem. B, 112, pp. 9020-9041.

10.1021/jp800161418593145PMC2652252

Kanhaiya, K., Nathanson, M., in ‘t Veld, P.J., Zhu, C., Nikiforov, I., Tadmor, E.B., Choi, Y.K., Im, W., Mishra, R.K. and Heinz, H. (2023). “Accurate force fields for atomistic simulations of oxides, hydroxides, and organic hybrid materials up to the micrometer scale”, J. Chem. Theory Comput., 19, pp. 8293-8322.

10.1021/acs.jctc.3c00750

Li, D.D. and Greenfield, M.L. (2014). “Viscosity, relaxation time, and dynamics within a model asphalt of larger molecules”, J. Chem. Phys., 140, 034507.

10.1063/1.4848736

Lu, C., Wu, C., Ghoreishi, D., Chen, W., Wang, L., Damm, W., Ross, G.A., Dahlgren, M.K., Russell, E., Von Bargen, C.D., Abel, R., Friesner, R.A. and Harder, E.D. (2021). “OPLS4: Improving force field accuracy on challenging regimes of chemical space”, J. Chem. Theory Comput., 17, pp. 4291-4300.

10.1021/acs.jctc.1c00302

Marenich, A.V., Jerome, S.V., Cramer, C.J. and Truhlar, D.G. (2012). “Charge Model 5: An extension of Hirshfeld population analysis for the accurate description of molecular interactions in gaseous and condensed phases”, J. Chem. Theory Comput., 8, pp. 527-541.

10.1021/ct200866d

Mishra, R.K., Kanhaiya, K., Winetrout, J.J., Flatt, R.J. and Heinz, H. (2021). “Force field for calcium sulfate minerals to predict structural, hydration, and interfacial properties”, Cem. Concr. Res., 139, 106262.

10.1016/j.cemconres.2020.106262

Rappé, A.K., Casewit, C.J., Colwell, K.S., Goddard, W.A. and Skiff, W.M. (1992). “UFF, a full periodic table force field for molecular mechanics and molecular dynamics simulations”, J. Am. Chem. Soc., 114, pp. 10024-10035.

10.1021/ja00051a040

Reed, A.E., Weinstock, R.B. and Weinhold, F. (1985). “Natural population analysis”, J. Chem. Phys., 83, pp. 735-746.

10.1063/1.449486

Senftle, T.P., Hong, S., Islam, M.M., Kylasa, S.B., Zheng, Y., Shin, Y.K., Junkermeier, C., Engel-Herbert, R., Janik, M.J., Aktulga, H.M., Verstraelen, T., Grama, A. and van Duin, A.C.T. (2016). “The ReaxFF reactive force-field: development, applications and future directions”, npj Comput. Mater., 2, 15011.

10.1038/npjcompumats.2015.11

Storer, J.W., Giesen, D.J., Cramer, C.J. and Truhlar, D.G. (1995). “Class IV charge models: A new semiempirical approach in quantum chemistry”, J. Comput.-Aided Mol. Des., 9, pp. 87-110.

10.1007/BF00117280

Sun, H. (1998). “COMPASS: An ab initio force-field optimized for condensed-phase applications”, J. Phys. Chem. B, 102, pp. 7338-7364.

10.1021/jp980939v

van Duin, A.C.T., Dasgupta, S., Lorant, F. and Goddard, W.A. (2001). “ReaxFF: A reactive force field for hydrocarbons”, J. Phys. Chem. A, 105, pp. 9396-9409.

10.1021/jp004368u

Vanommeslaeghe, K., Hatcher, E., Acharya, C., Kundu, S., Zhong, S., Shim, J., Darian, E., Guvench, O., Lopes, P., Vorobyov, I. and MacKerell, A.D. (2010). “CHARMM general force field”, J. Comput. Chem., 31, pp. 671-690.

10.1002/jcc.2136719575467PMC2888302

Wang, J., Wolf, R.M., Caldwell, J.W., Kollman, P.A. and Case, D.A. (2004), “Development and testing of a general amber force field”, J. Comput. Chem., 25, pp. 1157-1174.

10.1002/jcc.20035

Yeh, I.-C. and Berkowitz, M.L. (1999). “Ewald summation for systems with slab geometry”, J. Chem. Phys., 111, pp. 3155-3162.

10.1063/1.479595

Zhang, L. and Greenfield, M.L. (2007). “Analyzing properties of model asphalts using molecular simulation”, Energy Fuels, 21, pp. 1712-1716.

10.1021/ef060658j

Journal of the Korean Asphalt Institute ISSN:2234-0785(Print) 2635-9553(Online) 한국아스팔트학회지

Preview

A multi-force-field stack design and validation framework for asphalt molecular dynamics: Part I. Theoretical background and general-purpose force-field comparison)

ABSTRACT

MAIN

(1)

(2)

(3)

(4)

(5)

(6)

Fig. 1.

4-Lane structure of the multi-force-field stack

(7)

(8)

Fig. 2.

Mixing-rule compatibility between 5 organic FFs and 5 IFF host variants

Table 1.

Comparison summary of the major force fields

Acknowledgements

References