有効遮蔽媒質法の応用と拡張：より現実的な電気化学系の記述に向けて大谷実

エネルギー材料シミュレーショングループ
ナノシステム研究部門
産業技術総合研究所

常行班：第一原理分子動力学法、熱物理とダイナミクス

2011/10/07 東大工学部

Outline
Introduction

Examples of application of electrochemistry
Why do we need the effective screening medium method?

MethodWhat is the ESM method?
How to apply bias potential to metal/water interface

Application

Model of the electrochemical system
Electrochemical reactions (a few steps in the electrolysis of water)

Extension

Limitation of the current ESM method
Toward more realistic modeling

SummaryApplied electrochemistryBattery
Manganese dry cell
Lead battery
NiCd, NiH secondary battery
Fuel cell
Lithium secondary battery
primary: only dischargeablesecondary: chargable
Capacitor
Electrolytic condenser
Double layer condenser
Photovoltaic cellc-Si, a-Si solar cell
Dye sensitized solar cellSensorpH meterion selective concentration
meter
glucose, etc. (using enzyme)gas (oxygen, etc.)
Separator
DNA, protein (usingelectromigration)
Electroplating
Cathodic protection
Electrolysis
Aluminum, Copper, etc.
Water, salt, etc.
Organic chemicals
tetraethyl lead

Electrolysis vs. Fuel cell

Cathode:Anode:$%&$2H+ + 2e− → H2
2H2 O → O2 + 4H+ + 4e−!"#-.+,!'()*+,!

Ｈ%#&%#$＋ !"###"$!＋ "#%&$!＋
Ｈ"$'()*%

!

+,-.)*%!"#!"#

%&%&$

$%&$

Hydrogen fuel cell$＋Electrolysis of water

O2 + 4H+ + 4e− → 2H2 O
H2 → 2H+ + 2e−

Goal of our study• Develop a new calculation method tosimulate electrochemical systems.

• Reproduce electrochemical reactions on anelectrode/electrolyte interface• Calculate the free energy difference andactivation energy of electrochemical
reactions.

• Design or predict high efficiency, durable andcheap catalysts for the fuel cell or electrolysisof water.

Things to be considered:!

Ion distribution!

Screening effect of water!

Interaction btw water & metal!

Electronic structure!

Bias potential!

electrode

Electrode/electrolyte interface










V
εf

Electric double layer

Effective Screening Medium (ESM) method
M.O. and O. Sugino, PRB 73, 115407 (2006)µV∞z
Electric double layer

Effective screening medium method
M.O. and O. Sugino, PRB 73, 115407 (2006)

Slab

waterelectrodeelectrodeSlab model

VεfA slab model for an electrode! Put water molecules! Introduce a counter electrodewith an imaginary battery!

µV∞z
Electric double layer

ESM method

M.O. and O. Sugino, PRB 73, 115407 (2006)

Slabwaterelectrode

Electrochemical cell

ESM method

M.O. and O. Sugino, PRB 73, 115407 (2006)

Electrochemical cellESM�r = 1conventional formalism:
Solved PE with PBC
�(r) = 1�ρtot (r � )V (r) = dr|r − r � |SlabESMwater�r = ∞ESM formalism:
Model dependent BCs
�(r) : model dependent
V (r) =

�drG(r, r � )ρtot (r � )ESM method

M.O. and O. Sugino, PRB 73, 115407 (2006)

Electrochemical cell!
+!
+Slab ! water + ESM!
+ESMx

�r = 1zbare coulomb part∇[�(r)∇]V (r) = −4πρtot (r)�V (g� , z1 ) = 0∂z V (g� , z → −∞) = 0g� = (gx , gy )4π −g� |z−z� |G(g� , z, z ) =e2g�
4π −g� (2z1 −z−z� )Image charge − e2g�part�
V (g� , z) = dzG(g� , z, z � )ρtot (g� , z � )�Advantages of the ESM method• Various boundary conditions areapplicable.

• We do not need any correction or

modification of the coulomb interaction.

• Easy to implement without additional
calculation costs.

• SIESTA, PWscf, OpenMX are available.

Advantages of the ESM method• Various boundary conditions areapplicable.

Batteryany correction or• We do not need!"#$electrodemodification of the coulomb interaction.

additional• Easy to implement without Capacitor
calculationIsolatedslab costs.

electrodeelectrode!"#$!"#$PWscf, OpenMX are available.
• SIESTA,

Implementation of the ESM
M.O. and O. Sugino, PRB 73, 115407 (2006)

Only the electrostatic parts are modified:• Ewald energy & force (each MD step)• Local (long-rang) part (each MD step)• Hartree potential (each SCF step)

Conventional g-component

gy

g-component in ESM methodgygzzgxgxElectrolysis of water
2H2 O → 2H2 + O2HEROER

H3 O+ + e− → H2 O + Had

H2 O → (H2 O)adHad + Had → H2(H2 O)ad → (OH)ad + H3 O+ + e−or

H3 O + Had + e → H2 O + H2+−Direct simulation(OH)ad + (OH)ad → Oad + H2 O
Oad + Oad → O2

Constrained dynamics (Blue-moon ensemble)

Simulation condition
32 water molecules
36 Pt atoms (3-layers)
GGA-PBE
Plane wave - Ultrasoft pp
Temperature: 80℃

r= 
metal
vacuum with ze
barrierpotentials zw

H₂OwithH₃O+Bias potential ( -/+ )Pt(111) slab
vacuumH3 O+OH−r = 1

Simulation modelsystem
A3A2A1B
C1C2aC2bC3HERQESM(e/cell) component on surface-0.7
H65O32
H3O+0.35
H65O32
H3O+0.0
H65O32
H3O+0.0
H64O320.0
H63O32
OH0.35
H63O32
OH0.35
H63O322OH
0.7
H63O322OH

OER

in water

HER
PZC
H3O+OER

H3O+O. Sugino, et al., Surf. Sci. 601, 5237 (2007)
T. Ikeshoji, et al., Phys. Chem. Chem. Phys. in press
M. Otani et al., J. Phys. Soc. Jpn. 77, 024802 (2008)
M. Otani, et al., Phys. Chem. Chem. Phys. 10, 1609 (2008)

MD simulations with bias potential
ESMGrotthuss mechanism

Charge redistribution∆ρ = ρ(charged) − ρ(neutral)⎫⎜⎜ polarization of water⎬⎜⎜⎜⎭ Large charge redistribution occurs
at the interface

Accumulated charge on
metal surface

Electrolysis of water
2H2 O → 2H2 + O2HEROER

H3 O+ + e− → H2 O + Had

H2 O → (H2 O)adHad + Had → H2(H2 O)ad → (OH)ad + H3 O+ + e−or

H3 O + Had + e → H2 O + H2+−Direct simulation(OH)ad + (OH)ad → Oad + H2 OOad + Oad → O2

Constrained dynamics (Blue-moon ensemble)

Hydrogen evolution reaction
2H + 2e → H2−+Hydrogen evolution+−HO+e→ H2 O + Had
3Volmer step:e-Heyrovsky step:Electrode

ElectrodeH3 O+ + Had + e− → H2 O + H2

e-ElectrodeH2Tafel step: Had + Had → H2H2ElectrodeElectrochemical reaction (cathode)Hydronium ion diffuses vis
the Grotthus mechanism
H atom is adsorbed on a top
site of Pt surface.

ePtPt

Q=0.95 (e/cell)Charge population [e]Charge population [e]Charge transfer-0.05

Charge population of top
most Pt layer averaged by
the number of Pt atoms

-0.06
-0.07
-0.08

Pt-0.09
-0.10
0.40

H0.30

Charge population of theadsorbed hydrogen atom.

0.20
0.10
0.00
01

2

3

4

567

Simulation time [ps]• Charge population of the Pt surface gradually decreases until thereaction and suddenly increases at the reaction.

• Charge population of the hydrogen atom slightly decreases before
the reaction and suddenly decreases at adsorptionElectron transfer from Pt to wateree

Apply bias potentialeee e HOMO

PtPt
LUMO

1412Energy (eV)108642¡f0-2-4

5.27

5.28

5.29
5.3
5.31
Simulation time (ps)

5.32

0.2
0.18
0.16
0.14
0.12
0.1
0.08
0.06
0.04
0.02
0eee
DOS weighted by GPOP (eV<1)

LUMO

LUMO
HOMO

Electrochemical reaction (anode)

Water dissociation:
2H2 O + Pt → Pt-OH− + H3 O+

e-++++++++Pt

Pt

QESM=-0.35 (e/cell)Electrochemical reaction (anode)

Water formation:
Pt-OH− + H3 O+ → 2H2 O + Pte-++++

++++

Pt

Pt

Q=-0.35 (e/cell)Reversible reaction of water splitting!"#$%&'()

Q=0 (e/cell)

Q=-0.35 (e/cell)

Q=-0.35 (e/cell)

Q=-0.70 (e/cell)ac(d)
muheishehehis, itsa
hisge
Pt
acof
en
hissoits
his
be-tials in the solution except the electrode and vacuum interface
regions after giving σ = 7 µC/cm2 was ca. 4.42 V against theFermi energy of Pt. The same potential for pure water system,
which gives PZC, is ca. 4.84 V. This means electrode potential is ca. 0.42 V vs. PZC and corresponds to 0.6 ∼ 0.7V vs.
SHE, taking into account the experimentally determined PZC
31∗ The estimated potential, 0.6 ∼ 0.7of Pt (0.23profile
V vs. SHE).
DensityofV vs. SHE, is close to the experimentally observed
potential
Activation
barrier of waterOH
distance
15,16(0.7 ∼ 0.8 V).

Activation barrier hight & potential
splitting∆E = RT ∆ log(ρ) = 0.02 eV"!"#$!$#"$)$,!)$)$*+$(a)(b)=>9?@4A$@>@A?@$B7A'4CC7D84C

"!"#$)$&'$!$&'#"$)$,!)$)$*+$

3#,$456789:7;

Current@E88*;'

onns,
all!"#$!$#!%&'($)$!)$)$*+$

-./0

06.$3!< -.2$0

(c)

1.10

1."0

Fig. 2 Two-dimensional density profiles of two O–H distances in
the hydrogen transfer reaction. Different definition of OH distances,r1 and r2 , are indicated in each graph. (a) hydrogen transfer between
H2 O and OH on the Pt surface with σ = 0, (b) water electrolysis toform H3 O+ and OH on the Pt surface with σ = 0.35 µC/cm2 , (c)Voltage+electrode potential= 0.6~0.7 V(exp. = 0.7~0.8V)

What we have achieved are...
direct simulation ofelectrochemical reactions
potential dependence ofinterface structure andvibrational frequencies
calculation of activationenergies (using Blue-moon
ensemble method)

Toward more realistic modeling...

!= vacuum with ze
barrierpotentials zw

H₂O
with
H₃O+Pt(111) slab
vacuum! = !

What we have achieved are...
direct simulation ofelectrochemical reactions
potential dependence ofinterface structure andvibrational frequencies
calculation of activationenergies (using Blue-moon
ensemble method)

Toward more realistic modeling...

!= vacuum with zebarrier
potentials zw!= H₂Owith
H₃O+H₂O
with
H₃O+

Pt(111) slab

Pt(111) slab

vacuumvacuum

! = !

! = !

Toward more realistic modeling...
�(z) = ea(z−z1 )�r = ∞

�r = 1

�r = 1�r → ∞

z1

z1

ρ(r)

ρ(r)

Discontinuity

Smooth ESM

No discontinuity

To limit the electronic structure computational cost, it is ad#14.6#14.5COOH
are bonded antiOCH
together.
Theparameters
are#12.3chosen so that
vantageous
to limit as intermuch as possible the size of the discretizaCH CONH
#13.8
#13.6
kinetic
energy,tion grid. It is thus important for the dielectric function not totowavesolvation
freeelectrostatic
energies
inref. 32.good agreement with exp
DFT/GGA
energies froma finer grid than thepseudopotentials andfunctions.
and therequire
exchangeandThis implicitly introduces an upper limit for the parameter ".
values.
The possibility of defining an interface solute–soFor theapplications presented
in this article,
we use the Hactrostatic
energyE,esexperimental cohesion energy (9.9 kcal/mol ). However, thismann pseudopotential and PBE exchange correlation funcvalue could be based
computed in uponprinciple byinitio supercell density isosu
molecular
cavity—anan abelectronic
In pseudopotentialbecalculations, the valence
electronicenergytional.
functional
calculation of liquid water. Among the many possible different
density # can be very low in the core of the atoms. To avoid
choices, by# " 0.0004," " 1.3 or # " et0.0006," " in1.7 satisfy
this isodensity peffects on ! , we add a finite density to thefunction proposed# usedbeen
Foresmanal.
their%["] asspurious
thesolutionofcriterion. For the first set of parameters, Figure 1 illustrates theto define ! at the atoms centers. The exact form of this additional
22 methanol molecule comshapeof the dielectric
function around aIt
density is not important because it is only needed
toensure
thecontinuummodel(IPCM).
isjustifiedby the ideaa(z−z1)
pared to the cavities used in PCM.
saturation of ! [#] to 1 inside the atomic cores. The ionic forces are
�(z) = ecomputed using the Hellmann–Feynman theorem.
No additionalis located in regions where the electronic denssolvent
forces are introduced by the dielectric function because it does not
Numerical
Results
explicitly depend on the atomic positions.
solute is small.
However,
to be able to represent accurat
In the determination of (2)the parameters, it is natural to choose
Table 1 shows the solvent-interaction energies—the differences! " 78.36 for water (experimental dielectric constant at low
regulargrid,
it isthe importantforto medium—have
a smoothfunbetweencalculations
withand1
withoutus
dielectric=∞�→∞frequency
and ambientconditions), the�only
solvent
consideredin�
=1
�
=r
rr
obtained for a series rof small neutral organic molecules. All thethis article. The parameters # and " control the magnitude of theThesolute–solvent
interface
being
anyway much mvalueshavebeenobtainedby fully optimizing
the molecularsolvation effect. By increasing # and/or decreasing ", one iny of thegeometries. The introduction of the dielectric medium lowers thecreases computation
the solvation effect. We define the parameters # , " usingsimple
surfacemadeinterlocking
sphetotalaenergies.
This effectis more important
for theof
mostpolar
the following criterion: we require the decreaseplicatedin energy due tothan
SmoothESMzzmolecules.
Table
1
also
shows
thevalues
obtainedfortheelectro1 molecule to be equal to the1
the dielectric medium for a single water
static solvationenergy in thePCM approach.
The larger
is noto limit
ourmodelto values
a fordiscontinuous
cohesion energy of liquid water at room temperature.
We usereason
thethe solvent effects obtained with our model can be explained byvariation needs to be
the choiceof models
our criterion to define# and ". Indeed,
we implicitly As a dielectric
function as inthedescribedabove.
include nonelectrostatic effects into the potential $, while in PCM,
ergy functional. Using
anadditional repulsion–
dispersion term comeslocal
from othermodel
models.
!["] we introduce
thefollowing
A good agreement
between the two approaches could be obtained
Cavitymodel
ostatic energy can be
if desired by tuning the parameters # and ". Table 1 also shows
332aSmooth ESM

3130

60!

0!

!

!

00000"#2)!
"#r$/"
1(#( (10$Organic Compounds (in [kcal/mol] in favor of second! !fortautomer/conformer).
#Some" #r$$#1$
1$2) .
2
1 $ #"#r$/"0 $Table 2. Change in Relative Stability Induced by Dielectric Solvent

2(3)

dr.

2-Hydroxypyridine/2-pyridone
3-Hydroxypyrazole/4-pyrazolin-3-one
Acetic acid syn-/anticonformer

This workOthers

5.4
8.4
3.3

4.6a
8.7bc
3.2, 3.4,d 4.8eThis smooth function
depends &onF.threeparameters:!(, "0J-L FattebertGygi,JCC(2001)
Experimental
value fromthe asymptotic
value
ofref.!33.! for " 3 0. "0 denotes the critic
Figure 1. Contour plot of the dielectric function ! with # " 0.0004
DFT/PCM from ref. 34.
" " 1.3
CH OH in atheplane containing COO. The contour theespectandto" forgives
from ref.interface
35.
in the middleDFT/PCMof thesolute–solvent, while ) cana03bc

Summary• Develop a new calculation method tosimulate electrochemical systems.

• Reproduce electrochemical reactions on anelectrode/electrolyte interface• Calculate the free energy difference andactivation energy of electrochemical
reactions.

• Some new extensions are applied to the ESM
method