Electrical Conductivity of Molten Salts and Ionic Conduction in Electrolyte Solutions

Abstract

A microscopic description for the partial DC conductivities in molten salts has been discussed by using a Langevin equation for the constituent ions. The memory function γ(t) can be written as in the form of a decaying function with time. In order to solve the mutual relation between the combined-velocity correlation functions Zσ±(t) and the memory function γ(t) in a short time region, a new recursion method is proposed. Practical application is carried out for molten NaCl by using MD simulation. The fitted function is described by three kinds of Gaussian functions and their physical backgrounds are discussed. Also the electrical conductivity in aqueous solution of electrolyte has been obtained, based on a generalized Langevin equation for cation and anion in it. This treatment can connect and compare with the work of computer simulation. The obtained results for concentration dependence of electrical conductivity are given by a function of the square root of concentration. The electrophoretic effect and the relaxation one are also discussed.

1. Introduction

The phenomena of transport properties in ionic liquids are of great important in the industrial science and technology, as well as in physics and chemistry. In connection with these, a number of experimental and theoretical studies have been published until the present time [1, 2, 3]. Ionic liquids are mainly classified into two categories; one is a group of molten salts and the other is a large number of electrolytic solutions, in particular, aqueous solutions of electrolytes.

In the case of molten salts, Sundheim discovered that the ratio of the partial conductivities of cation and anion were always equal to their inverse mass ratio, namely, σ⁺(DC)/σ⁻(DC) = m⁻/m⁺ [4].

Later on, this golden rule or a unified rule was theoretically explained by our group [5, 6, 7, 8, 9]. Detailed procedure will be shown in what follows.

Paralleling to above discovery, a number of scientific studies in molten salts have been developed from 1960s by several researchers [10, 11].

In order to study the structural and transport properties in molten salts, experimental investigations and molecular dynamics simulations have also been carried out from mid-70s of the last century [12, 13, 14, 15, 16].

Following to these, we have been engaged in the study of transport properties in molten salts [6, 7, 8, 9, 17]. We have carried out a theoretical study on the electrical conductivity of molten salts, starting from the Langevin equation and the velocity correlation functions for the constituent ions. Subsequently this treatment was successful to obtain the golden rule σ⁺/σ⁻ = m⁻/m⁺ in a microscopic view point.

It remains, however, unclear how the adopted Langevin equation can be effectively solved within a short time region, under an appropriate memory function, because our former theory was only successful to get the partial conductivities.

We like to discuss more generally the correlation between the velocity correlation functions incorporated with the partial DC conductivities and some of useful memory functions which are closely related to the friction constants acting on cations and anions in molten salts.

Preceding the investigation for molten salts, on the other hand, there have been a number of studies for ionic solutions since the discovery of Faraday, in which a typical example is electrolytic solution. During such long-termed history of electrochemistry, it was well established by Kohlrausch that the experimental results on the ionic conductivities in dilute electrolytic solutions indicated the law of independent migration of ions, Λ_c = Λ₀ − kc^1/2, where Λ₀ being the conductivity in the dilute limit and c the concentration and k the constant specified by the electrolyte dissolved in water.

The beginning of the modern aspect, in particular, on the thermodynamic and transport properties in electrolytic solutions might be originated from Debye-Hückel theory [18].

In order to explain the ionic conductivity in electrolytic solution, successful works following to Debye-Hückel theory have been reported by Onsager [19], Prigogine [20], and Fuoss and his co-worker [21]. In these theories, Λ₀ is treated by the Stokes law and the concentration dependence is mainly explained by the electrophoretic effect and relaxation one. Therefore, these treatments are based on a kind of mixing of the microscopic and partially macroscopic view point.

Starting from the Liouville equation, statistical mechanics of irreversible process for the ionic conductivity in electrolytic solution have been developed by Davis and Résibois [22] and Friedman [23], although they did not derive any explicit expressions for the friction constant in terms of inter-particle interactions.

It has been required to investigate the static and dynamic properties of dissolved ions in aqueous solutions from the microscopic view point. Along this requirement, the technique of molecular dynamic simulation has been applied, using some qualified inter-particle potentials. Various theoretical attempts have been recently tried to establish the dynamical behaviors of dissolved ions in these solutions, which is able to discuss parallel with results obtained by MD simulation [24, 25, 26].

Chandra and Bagchi [27] have developed a new theoretical approach to study the ionic conduction in electrolytic solutions, based on the combination of the mode coupling theory and the generalized Langevin equation, and they were successful to obtain the Onsager equation. However, there still remains the task to obtain how to derive the theoretical formula for Λ₀ in terms of inter-particle potentials and corresponding pair distribution functions.

We will apply the linear response theory for the electrolytic solution and to obtain Λ₀ and the concentration dependence of the conductivity in terms of pair-wise potentials and pair distribution functions among ions and water molecules, which can compare parallel with dynamical properties of MD simulation [28].

In addition, we will also clarify how the electrophoretic and relaxation effects treated by many researchers are explained in a microscopic view point.

From these, we will see what is similar and what is different for the case of molten salts and that of electrolytic solutions.

2. Generalized Langevin equations for the cation and anion in a molten salt

Let us consider a molten salt composed of the density n⁺ = n⁻ = n₀ (= N/V₀), of the constituent ion’s masses m⁺ and m⁻, and of the charge z⁺ = − z⁻ = z = 1, where N being the total number of cation and/or anion in the volume V₀.

A golden rule, σ⁺(DC)/σ⁻(DC) = m⁻/m⁺, can be obtainable from a generalized Drude theory, as a law of motion under an electric field [5].

As an extension, the generalized Langevin equation for an arbitrary cation or anion in the system under an external field E is written as follows:

$${\mathrm{m}}^{\pm }\mathrm{d}{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right)/\mathrm{d}t=-{\mathrm{m}}^{\pm }{{\int }_{-\infty }}^t{\mathrm{\xi }}^{\pm }\left(t-t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right)\mathrm{d}t\prime +{{\mathbf{R}}_{\mathrm{i}}}^{\pm }\left(t\right)+{\mathrm{z}}^{\pm }\mathrm{e}\mathbf{E}$$E1

where ξ^±(t) and R_i^±(t) are the retarded friction function in relation to the friction force and the random fluctuating force, acting on the cation or anion i, respectively.

After taking the ensemble average, equations of time evolution based on Eq. (1) in respect to the partial ionic conductivities are then written as follows:

$${\mathrm{m}}^{\pm }\mathrm{d}<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right){{\mathbf{v}}_{\mathrm{j}}}^{\pm }\left(0\right)>/\mathrm{d}t=-{\mathrm{m}}^{\pm }{\int }_{-\infty }^t<{\mathrm{\xi }}^{\pm }\left(t-t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right){{\mathbf{v}}_{\mathrm{j}}}^{\pm }\left(0\right)>\mathrm{d}t\prime \left(\text{for}\mathrm{i}=\mathrm{j}\text{and}\mathrm{i}\ne \mathrm{j}\right)$$E2

and

$${\mathrm{m}}^{\pm }\mathrm{d}<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right){{\mathbf{v}}_{\mathrm{k}}}^{\mp }\left(0\right)>/\mathrm{d}t=-m^{\pm }{\int }_{-\infty }^t<{\mathrm{\xi }}^{\pm }\left(t-t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right){{\mathbf{v}}_{\mathrm{k}}}^{\mp }\left(0\right)>\left(\text{for}\mathrm{i}\ne \mathrm{k}\right)$$E3

And the equation of time evolution in relation to the diffusion constants of constituent ions is written as follows:

$${\mathrm{m}}^{\pm }\mathrm{d}<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(0\right)>/\mathrm{d}t=-{\mathrm{m}}^{\pm }{\int }_{-\infty }^t<{\mathrm{\xi }}^{\pm }\left(t-t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(0\right)>\mathrm{d}t\prime$$E4

As was previously illustrated [9], the retarded friction function ξ^±(t) cannot be independent for the averaging procedure and we have to define new memory functions as follows:

$$<{\mathrm{\xi }}^{\pm }\left(t-t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right){{\mathbf{v}}_{\mathrm{j}}}^{\pm }\left(0\right)>={{\mathrm{\gamma }}_{\mathrm{\sigma }}}^{\pm }\left(t\right)<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right){{\mathbf{v}}_{\mathrm{j}}}^{\pm }\left(0\right)>\left(\text{for}\mathrm{i}=\mathrm{j}\text{and}\mathrm{i}\ne \mathrm{j}\right)$$E5

and

$$<{\mathrm{\xi }}^{\pm }\left(t-t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right){{\mathbf{v}}_{\mathrm{k}}}^{\mp }\left(0\right)>={{\mathrm{\gamma }}_{\mathrm{\sigma }}}^{\pm }\left(t\right)<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right){{\mathbf{v}}_{\mathrm{k}}}^{\mp }\left(0\right)>\left(\text{for}\mathrm{i}\ne \mathrm{k}\right)$$E6

While, in the case of diffusion constants of constituent ions, that is, E = 0, we can define

$$<{\mathrm{\xi }}^{\pm }\left(t-t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(0\right)>={{\mathrm{\gamma }}_{\mathrm{D}}}^{\pm }\left(t\right)<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(0\right)>$$E7

It is emphasized that the memory functions γ_σ^±(t) is not equal to γ_D^±(t) as shown in previous paper [9]. In other words, the retarded friction function, ξ^±(t − t′), is a kind of vector function and is varied with the environment such as the existence of electric field E. Therefore, the memory function is varied in accordance with what sort of evolution is considered in the time-dependent correlation function [29].

Assuming that the ensemble average for the fluctuating force is zero and if we apply the following electric field,

$$\mathbf{E}\left(t\right)=Re{\mathbf{E}}_{0}exp\left(\mathrm{i\omega }t\right)$$E8

where Re means the real part and ω is the angular frequency, then the averaged ion’s velocity induced by this external filed is equal to

$$<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right)>=Re{\mathrm{\mu }}^{\pm }\left(\mathrm{\omega }\right){\mathrm{z}}^{\pm }\mathrm{e}\mathbf{E}\left(t\right)$$E9

where μ^±(ω) is the mobility of cation or anion.

Putting (9) into the equation of motion (1) after taking the ensemble average, we have

$${\mathrm{\mu }}^{\pm }\left(\mathrm{\omega }\right)=\left(1/m^{\pm }\right)\left[1/\left\{\mathrm{i\omega }+\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(\mathrm{\omega }\right)\right\}\right]$$E10

where

$$\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(\mathrm{\omega }\right)={{\int }_0}^{\infty }{\mathrm{\gamma }}^{\pm }\left(t\right)exp\left(-\mathrm{i\omega }t\right)\mathrm{d}t$$E11

Therefore, the current density is written as follows:

$${\mathrm{j}}^{\pm }\left(t\right)={\mathrm{nz}}^{\pm 2}{\mathrm{e}}^2<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right)>=Re{\mathrm{nz}}^{\pm 2}{\mathrm{e}}^2{\mathrm{\mu }}^{\pm }\left(\mathrm{\omega }\right)\mathbf{E}\left(t\right)$$E12

The partial conductivity is, then, equal to

$${\mathrm{\sigma }}^{\pm }\left(\mathrm{\omega }\right)={\mathrm{nz}}^{\pm 2}{\mathrm{e}}^2{\mathrm{\mu }}^{\pm }\left(\mathrm{\omega }\right)=\left({\mathrm{nz}}^{\pm 2}{\mathrm{e}}^2/{\mathrm{m}}^{\pm }\right)\left\{1/\left(\mathrm{i\omega }+\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(\mathrm{\omega }\right)\right)\right\}$$E13

and in the limit of ω = 0,

$${\mathrm{\sigma }}^{\pm }\left(\mathrm{DC}\right)={\mathrm{nz}}^{\pm 2}{\mathrm{e}}^2{\mathrm{\mu }}^{\pm }\left(0\right)=\left\{{\mathrm{nz}}^{\pm 2}{\mathrm{e}}^2/{\mathrm{m}}^{\pm }\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(0\right)\right\}$$E14

Therefore, $\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(0\right)$ is equal to the effective friction constant acting on each ion.

According to our previous studies [7, 8, 9], the following relation was recognized:

$$\stackrel{\sim }{\mathrm{\gamma }}{}^{+}\left(0\right)=\stackrel{\sim }{\mathrm{\gamma }}{}^{-}\left(0\right)\equiv \stackrel{\sim }{\mathrm{\gamma }}\left(0\right)$$E15

where $\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)$ is expressed as follows:

$$\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)={\left({\mathrm{\alpha }}^0/\text{3μ}\right)}^{1/2},\left(1/\mathrm{\mu }\right)=\left(1/{\mathrm{m}}^{+}\right)+\left(1/{\mathrm{m}}^{-}\right)$$E16

and

$${\mathrm{\alpha }}^0=n{{\int }_0}^{\infty }\left[{\partial }^2{\mathrm{\varphi }}^{+-}\left(r\right)/\partial r^2+\left(2/r\right)\left\{{\mathrm{\partial \varphi }}^{+-}\left(r\right)/\partial r\right\}\right]{\mathrm{ɡ}}^{+-}\left(r\right)·\text{4π}{r}^2\mathrm{d}r$$E17

ϕ⁺⁻(r) and ɡ⁺⁻(r) in this equation are the inter-ionic potential between cation and anion and the corresponding pair distribution function, respectively.

Therefore, we have a golden rule for the partial conductivities in a microscopic scale as follows:

$${\mathrm{\sigma }}^{+}\left(\mathrm{DC}\right)/{\mathrm{\sigma }}^{-}\left(\mathrm{DC}\right)={\mathrm{m}}^{-}/{\mathrm{m}}^{+}$$E18

In the following sections, as a numerical example, the MD simulation on molten NaCl at 1100 K is often utilized, for which the interionic potential functions suggested by Tosi and Fumi [30] for a study of solid alkali halides are applied. In order to make sure that the Tosi-Fumi potential for NaCl can be valid in the liquid state, we have estimated the partial pair distribution functions of molten NaCl liquid, ɡ_ij(r) (i,j = Na⁺, Cl⁻) as shown in Figure 1, which agree with those of experimental results obtained by Edwards et al. [31].

Figure 1.

Pair distribution functions, g_ij(r), for molten NaCl at 1148 K, obtained by MD simulation.

Using these ɡ_ij(r), we have also estimated the total neighboring numbers around arbitrary ions located at the distance r, which describe as n_ij = 4πʃ₀^rr²dr, as shown in Figure 2a–c.

Figure 2.

(a) gNa-Cl(r) and nNa-Cl(r) for molten NaCl at 1148 K, obtained by MD simulation. (b) gNa-Na(r) and nNa-Na(r) for molten NaCl at 1148 K, obtained by MD simulation. (c) gCl-Cl(r) and nCl-Cl(r) for molten NaCl at 1148 K, obtained by MD simulation.

The nearest neighbor number is defined as n_ij(r₁), where r₁ is the position of the first minimum of ɡ_ij(r).

Then, the nearest neighbors around a Na⁺ are nearly equal to 5.0, since the distance r₁ is taken at the minimum position of ɡ_Na-Cl(r) as shown in Figure 2a.

The application of Tosi-Fumi potentials in the MD simulations for viscosity and electrical conductivity is also valid to reproduce their experimental results [5].

Therefore, the following MD simulations for molten NaCl must be reliable to see their microscopic view.

3. Linear response theory for the partial conductivities

On the other hand, according to our previous investigations [6, 7, 8, 9, 17, 29], the partial DC conductivities σ⁺(DC) and σ⁺(DC) are expressed as follows,

$${\mathrm{\sigma }}^{+}\left(\mathrm{DC}\right)={\mathrm{\sigma }}^{++}+{\mathrm{\sigma }}^{+-}=\left(1/{\text{3k}}_{\mathrm{B}}\mathrm{T}\right){{\int }_0}^{\infty }<{\mathbf{j}}^{+}\left(t\right)\mathbf{j}\left(0\right)>\mathrm{d}t$$E19

$${\mathrm{\sigma }}^{-}\left(\mathrm{DC}\right)={\mathrm{\sigma }}^{--}+{\mathrm{\sigma }}^{+-}=\left(1/{\text{3k}}_{\mathrm{B}}\mathrm{T}\right){{\int }_0}^{\infty }<{\mathbf{j}}^{-}\left(t\right)\mathbf{j}\left(0\right)>\mathrm{d}t$$E20

where

$${\mathrm{\sigma }}^{\pm \pm }=\left(1/{\text{3k}}_{\mathrm{B}}\mathrm{T}\right){{\int }_0}^{\infty }<{\mathbf{j}}^{\pm }\left(t\right){\mathbf{j}}^{\pm }\left(0\right)>\mathrm{d}t$$E21

$${\mathrm{\sigma }}^{+-}=\left(1/{\text{3k}}_{\mathrm{B}}\mathrm{T}\right){{\int }_0}^{\infty }<{\mathbf{j}}^{+}\left(t\right){\mathbf{j}}^{-}\left(0\right)>\mathrm{d}t$$E22

and

$$\mathbf{j}\left(t\right)={\mathbf{j}}^{+}\left(t\right)+{\mathbf{j}}^{-}\left(t\right)$$E23

where

$${\mathbf{j}}^{+}\left(t\right)=\sum _{\mathrm{i}=1}^{\mathrm{n}}{\mathrm{z}}^{+}\mathrm{e}{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(t\right),{\mathbf{j}}^{-}\left(t\right)=\sum _{\mathrm{k}=1}^{\mathrm{n}}{\mathrm{z}}^{-}\mathrm{e}{{\mathbf{v}}_{\mathrm{k}}}^{-}\left(t\right)$$E24

Considering the ensemble averages of (19) and (20), it is convenient to define the velocity correlation functions Z_σ⁺(t) and Z_σ⁻(t) as follows:

$${{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(t\right)\equiv <{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(t\right){{\mathbf{v}}_{\mathrm{j}}}^{+}\left(0\right)>-<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(t\right){{\mathbf{v}}_{\mathrm{k}}}^{-}\left(0\right)>$$E25

and

$${{\mathrm{Z}}_{\mathrm{\sigma }}}^{-}\left(t\right)\equiv <{{\mathbf{v}}_{\mathrm{k}}}^{-}\left(t\right){{\mathbf{v}}_{\mathrm{l}}}^{-}\left(0\right)>-<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(t\right){{\mathbf{v}}_{\mathrm{k}}}^{-}\left(0\right)>$$E26

where < > means the ensemble average.

Using (25) and (26), the partial DC conductivities (19) and (20) are written, respectively, as follows:

$${\mathrm{\sigma }}^{+}\left(\mathrm{DC}\right)=\left({\mathrm{nz}}^2{\mathrm{e}}^2/{\text{3k}}_{\mathrm{B}}\mathrm{T}\right){{\int }_0}^{\infty }{{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(t\right)\mathrm{d}t$$E27

$${\mathrm{\sigma }}^{-}\left(\mathrm{DC}\right)=\left({\mathrm{nz}}^2{\mathrm{e}}^2/{\text{3k}}_{\mathrm{B}}\mathrm{T}\right){{\int }_0}^{\infty }{{\mathrm{Z}}_{\mathrm{\sigma }}}^{-}\left(t\right)\mathrm{d}t$$E28

On the other hand, combining Eqs. (25) or (26) and (1), we have

$$\partial \left\{{{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(t\right)\right\}/\partial t=-{{\int }_0}^t{\mathrm{\gamma }}^{+}\left(t-s\right){{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(s\right)\mathrm{d}s$$E29

and/or

$$\partial \left\{{{\mathrm{Z}}_{\mathrm{\sigma }}}^{-}\left(t\right)\right\}/\partial t=-{{\int }_0}^t{\mathrm{\gamma }}^{-}\left(t-s\right){{\mathrm{Z}}_{\mathrm{\sigma }}}^{-}\left(s\right)\mathrm{d}s$$E30

Taking the Laplace transformation of ∂{Z_σ⁺(t)}/ ∂t in (29) as follows,

$$\mathcal{L}\left[\partial \left\{{{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(t\right)\right\}/\partial t\right]\equiv {{\int }_0}^{\infty }exp\left(-\mathrm{i\omega }t\right)\left[\partial \left\{{{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(t\right)/\partial t\right\}\right]\mathrm{d}t$$E31

$$={{\left[exp\left(-\mathrm{i\omega }t\right){{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(t\right)\right]}_0}^{\infty }+\mathrm{i\omega }{{\int }_0}^{\infty }exp\left(-\mathrm{i\omega }t\right){{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(t\right)\mathrm{d}t$$E32

$$=-{{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(0\right)+\mathrm{i\omega }{\stackrel{\sim }{\mathrm{Z}}}_{\mathrm{\sigma }}{}^{+}\left(\mathrm{\omega }\right)$$E33

Here, we have used an evident condition Z_σ⁺(t = ∞) = 0.

On the other hand, the right hand side of (29) is given by the following expressions:

$$\begin{array}{ll}\mathcal{L}\left\{-{{\int }_0}^{\infty }{\mathrm{\gamma }}^{+}\left(s\right){{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(t-s\right)\mathrm{d}s\right\}& =-{{\int }_0}^{\infty }exp\left\{-\mathrm{i\omega }\left(t-s\right)\right\}{\mathrm{\gamma }}^{+}\left(t-s\right)\mathrm{d}\left(t-s\right)\\ & {{\int }_0}^{\infty }exp\left(-\mathrm{i\omega }s\right){{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(s\right)\mathrm{d}s=-\stackrel{\sim }{\mathrm{\gamma }}{}^{+}\left(\mathrm{\omega }\right){\stackrel{\sim }{\mathrm{Z}}}_{\mathrm{\sigma }}{}^{+}\left(\mathrm{\omega }\right)\end{array}$$E34

Therefore we have,

$$-{{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(0\right)+\mathrm{i\omega }{\stackrel{\sim }{\mathrm{Z}}}_{\mathrm{\sigma }}{}^{+}\left(\mathrm{\omega }\right)=-\stackrel{\sim }{\mathrm{\gamma }}{}^{+}\left(\mathrm{\omega }\right){\stackrel{\sim }{\mathrm{Z}}}_{\mathrm{\sigma }}{}^{+}\left(\mathrm{\omega }\right)$$E35

$${\stackrel{\sim }{\mathrm{Z}}}_{\mathrm{\sigma }}{}^{+}\left(\mathrm{\omega }\right)={{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(0\right)/\left\{\mathrm{i\omega }+\stackrel{\sim }{\mathrm{\gamma }}{}^{+}\left(\mathrm{\omega }\right)\right\}=\left({\text{3k}}_{\mathrm{B}}\mathrm{T}/{\mathrm{m}}^{+}\right)/\left\{\mathrm{i\omega }+\stackrel{\sim }{\mathrm{\gamma }}{}^{+}\left(\mathrm{\omega }\right)\right\}$$E36

In a similar way, we have,

$$\stackrel{\sim }{\mathrm{Z}}{}_{\mathrm{\sigma }}^{-}\left(\mathrm{\omega }\right)={{\mathrm{Z}}_{\mathrm{\sigma }}}^{-}\left(0\right)/\left\{\mathrm{i\omega }+\stackrel{\sim }{\mathrm{\gamma }}{}^{-}\left(\mathrm{\omega }\right)\right\}=\left({\text{3k}}_{\mathrm{B}}\mathrm{T}/{\mathrm{m}}^{-}\right)/\left\{\mathrm{i\omega }+\stackrel{\sim }{\mathrm{\gamma }}{}^{-}\left(\mathrm{\omega }\right)\right\}$$E37

If an appropriate memory function γ(t), which is valid for both cation and anion in the system, is considered and its Laplace transformation is inserted into either (36) or (37), then we can get the partial AC conductivities.

4. Microscopic representation for the Z_σ⁺(t) and Z_σ⁻(t) in a molten salt

We have already shown the microscopic expressions for Z_σ⁺(t) and Z_σ⁻(t) as Taylor expansion forms in a molten salt in which the inter-ionic potential between cation and anion and the corresponding pair distribution function are concerned by Koishi et al. [7]. In these combined velocity correlation functions, it can be shown that the odd power terms of the time t have vanishing coefficients which, it turns out, is related to the fact that any positions and their differentiations with time are uncorrelated in an ensemble average. In facts, the velocity auto-correlation function can be expressed in terms of even powers of the time t [32, 33].

The short-time expansion forms of Z_σ⁺(t) and Z_σ⁻(t) are actually shown in the following forms:

$${{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(t\right)=\left({\text{3k}}_{\mathrm{B}}\mathrm{T}/{\mathrm{m}}^{+}\right)\left[1-\left(t^2/2\right)\left({\mathrm{\alpha }}^0/\text{3μ}\right)+\left(\text{over}{t}^4\right)\right]$$E38

and

$${{\mathrm{Z}}_{\mathrm{\sigma }}}^{-}\left(t\right)=\left({\text{3k}}_{\mathrm{B}}\mathrm{T}/{\mathrm{m}}^{-}\right)\left[1-\left(t^2/2\right)\left({\mathrm{\alpha }}^0/\text{3μ}\right)+\left(\text{over}{t}^{\mathit{4}}\right)\right]$$E39

Thus, the partial conductivities for cation and anion in a molten salt are written as in the following Kubo-formulae:

$${\mathrm{\sigma }}^{+}\left(\mathrm{DC}\right)=\left({\mathrm{n}}_0{\mathrm{e}}^2/{\mathrm{m}}^{+}\right){{\int }_0}^{\infty }\left\{1-\left(t^2/2\right)\left({\mathrm{\alpha }}^0/\text{3μ}\right)+\left(\text{over}{t}^4\right)\right\}\mathrm{d}t$$E40

and

$${\mathrm{\sigma }}^{-}\left(\mathrm{DC}\right)=\left({\mathrm{n}}_0{\mathrm{e}}^2/{\mathrm{m}}^{-}\right){{\int }_0}^{\infty }\left\{1-\left(t^2/2\right)\left({\mathrm{\alpha }}^0/\text{3μ}\right)+\left(\text{over}{t}^4\right)\right\}\mathrm{d}t$$E41

Using (14), (16), (40) and (41), we have a very interesting relation written in the following form:

$$1/\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)={{\int }_0}^{\infty }\left\{1-\left(t^2/2\right){\left(\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)\right)}^2+\left(\text{over}{t}^4\right)\right\}\mathrm{d}t$$E42

However, it is generally difficult to obtain Z_σ^±(t) from appropriate memory functions, by using the well-known method in statistical mechanics [33].

Under these circumstances, we explore a new method to solve Langevin Eqs. (29) and (30), in order to clarify a detailed correlation between γ(t) and Z_σ^±(t) within the short time region, which will be shown in later section.

5. Method of continued-fraction based on Mori formulae

Many years ago, Mori [34, 35] had generalized the Langevin equation starting from the Hamilton’s canonical equation of motion in a system of a monatomic liquid with the component’s mass as m. Along his theory, Copley and Lovesey [36] have concluded that the memory function in the generalized Langevin equation could be expressed as follows:

$$\partial \left\{{\mathrm{\gamma }}_{\mathrm{n}}\left(t\right)\right\}/\partial t=-{{\int }_0}^t{\mathrm{\gamma }}_{\mathrm{n}+1}\left(t-s\right){\mathrm{\gamma }}_{\mathrm{n}}\left(s\right)\mathrm{d}s\mathrm{n}=1,2,3,\dots$$E43

where γ_n(t) is the n-th stage memory function and the first stage memory function is equal to γ(t) in Eqs. (29) and (30). The Fourier-Laplace transform of the above equation provides the following continued-fraction representation,

$${\stackrel{\sim }{\mathrm{\gamma }}}_{\mathrm{n}}\left(\mathrm{\omega }\right)=-{\mathrm{\delta }}_{\mathrm{n}}/\left[\mathrm{\omega }+\stackrel{\sim }{\mathrm{\gamma }}{}_{\mathrm{n}+1}\left(\mathrm{\omega }\right)\right]$$E44

where the Mori coefficient δ_n is equal to γ_n(0).

The method of Copley and Lovesey [36] was able to express the short time expansion for the velocity correlation function Z(t) (= < v_i(t) v_j(0)>) described as in the following form:

$$\mathrm{Z}\left(t\right)={\mathrm{Z}}_0\left\{1-\left(t^2/2!\right){\mathrm{Z}}_2+\left(t^4/4!\right){\mathrm{Z}}_4-\left(t^6/6!\right){\mathrm{Z}}_6+\dots \right\}$$E45

Thus, they provided the following relations if several δ_n’s are known:

$${\mathrm{Z}}_0=\left({\text{3k}}_{\mathrm{B}}T/\mathrm{m}\right),{\mathrm{Z}}_2={\mathrm{Z}}_0{\mathrm{\delta }}_1,{\mathrm{Z}}_4={\mathrm{Z}}_0{\mathrm{\delta }}_1\left({\mathrm{\delta }}_1+{\mathrm{\delta }}_2\right),{\mathrm{Z}}_6={\mathrm{Z}}_0\left\{\left({\mathrm{\delta }}_1+{\mathrm{\delta }}_2\right)+{\mathrm{\delta }}_2{\mathrm{\delta }}_3\right\},\dots$$E46

Therefore, the problem is ascribed to the derivation of δ_n’s. Because of a hard task in such repeating calculations, it is difficult to obtain a number of δ_n’s. However, several applications along these procedures have been carried out [37, 38].

Instead of the method of continued-fraction described in the above, we will provide a simple but new method to obtain the mutual relation between the combined velocity correlation function Z_σ^±(t) and γ(t) in a short time region, in the following section.

6. Recursion formulae for Z_σ^±(t) and γ(t)

Here, we provide a new and useful method to solve the Langevin equation based on recursion process [29]. Its detail is shown below.

Let us consider a Langevin equation for an evolution function being equivalent to (29) and (30), as follows:

$$\mathrm{dy}\left(t\right)/\mathrm{d}t={{\int }_0}^t\mathrm{q}\left(t-s\right)\mathrm{y}\left(s\right)\mathrm{d}s$$E47

The power expansion for q(t) is defined as follows:

$$\mathrm{q}\left(t\right)=\sum _{\mathrm{n}=0}^{\infty }\left({\mathrm{q}}_{\mathrm{n}}/\mathrm{n}!\right)t^{\mathrm{n}}\left({\mathrm{q}}_{\mathrm{n}}={\mathrm{q}}^{\left(\mathrm{n}\right)}\left(0\right)\right)$$E48

and the corresponding expansion formula for y(t) is written as follows:

$$\mathrm{y}\left(t\right)=\sum _{\mathrm{m}=0}^{\infty }\left({\mathrm{y}}_{\mathrm{m}}/\mathrm{m}!\right)t^{\mathrm{m}}\left({\mathrm{y}}_{\mathrm{m}}={\mathrm{y}}^{\left(\mathrm{m}\right)}\left(0\right)\right)$$E49

Putting (48) and (49) into the right hand side of Eq. (47), we have

$$\begin{array}{ll}{{\int }_0}^t\mathrm{q}\left(t-s\right)\mathrm{y}\left(s\right)\mathrm{d}s& =\sum _{\mathrm{n},\mathrm{m}=0}^{\infty }\left({\mathrm{q}}_{\mathrm{n}}/\mathrm{n}!\right)\left({\mathrm{y}}_{\mathrm{m}}/\mathrm{m}!\right){{\int }_0}^t{\left(t-s\right)}^{\mathrm{n}}{s}^{\mathrm{m}}\mathrm{d}s\\ & =\sum _{\mathrm{n},\mathrm{m}=0}^{\infty }\left({\mathrm{q}}_{\mathrm{n}}/\mathrm{n}!\right)\left({\mathrm{y}}_{\mathrm{m}}/\mathrm{m}!\right)t^{\left(\mathrm{n}+\mathrm{m}+1\right)}{{\int }_0}^1{\left(1-p\right)}^{\mathrm{n}}{p}^{\mathrm{m}}\mathrm{d}p\\ & =\sum _{\mathrm{n},\mathrm{m}=0}^{\infty }\left({\mathrm{q}}_{\mathrm{n}}/n!\right)\left({\mathrm{y}}_{\mathrm{m}}/m!\right)t^{\left(\mathrm{n}+\mathrm{m}+1\right)}\mathrm{B}\left(\mathrm{n}+1,\mathrm{m}+1\right)\\ & =\sum _{\mathrm{n},\mathrm{m}=0}^{\infty }\left({\mathrm{q}}_{\mathrm{n}}/n!\right)\left({\mathrm{y}}_{\mathrm{m}}/m!\right)t^{\left(\mathrm{n}+\mathrm{m}+1\right)}\left\{\mathrm{\Gamma }\left(\mathrm{n}+1\right)\mathrm{\Gamma }\left(\mathrm{m}+1\right)\right\}/\left\{\mathrm{\Gamma }\left(n+\mathrm{m}+2\right)\right\}\\ & =\sum _{\mathrm{n},\mathrm{m}=0}^{\infty }\left({\mathrm{q}}_{\mathrm{n}}/\mathrm{n}!\right)\left({\mathrm{y}}_{\mathrm{m}}/\mathrm{m}!\right)t^{\left(\mathrm{n}+\mathrm{m}+1\right)}\left(\mathrm{n}!\mathrm{m}!\right)/\left(\mathrm{n}+\mathrm{m}+1\right)!\\ & =\sum _{\mathrm{n},\mathrm{m}=0}^{\infty }(\left\{\left({\mathrm{q}}_{\mathrm{n}}{\mathrm{y}}_{\mathrm{m}}\right)t^{\left(\mathrm{n}+\mathrm{m}+1\right)}/\left(\mathrm{n}+\mathrm{m}+1\right)!\right\}\\ & =\sum _{\mathrm{k}=1}^{\infty }\left(z_{\mathrm{k}}/\mathrm{k}!\right)t^{\mathrm{k}}\end{array}$$E50

where B(n + 1, m + 1) and Γ(n + 1) mean the beta-function and the gamma-function, respectively, and

$$z_{\mathrm{k}}=\sum _{\mathrm{k}=\mathrm{n}+\mathrm{m}+1}{\mathrm{q}}_{\mathrm{n}}{\mathrm{y}}_{\mathrm{m}}$$E51

On the other hand, the left hand side of Eq. (47) is equal to the following formulae:

$$\mathrm{y}\prime \left(t\right)=\left\{\sum _{\mathrm{k}=0}^{\infty }\left({\mathrm{y}}_{\mathrm{k}}/\mathrm{k}!\right)t^{\mathrm{k}}\right\}\prime =\sum _{\mathrm{k}=0}^{\infty }\left({\mathrm{y}}_{\mathrm{k}+1}/\mathrm{k}!\right)t^{\mathrm{k}}$$E52

Compare both expressions (50) and (52), we can get the recursion formulae as follows,

$${\mathrm{y}}_1=0;{\mathrm{y}}_{\mathrm{k}+1}=\left\{\sum _{\mathrm{m}=0}^{\mathrm{k}-1}\left({\mathrm{q}}_{\mathrm{k}-\mathrm{m}-1}{\mathrm{y}}_{\mathrm{m}}\right)\right\}\left(\mathrm{k}=1,2,\dots \right)$$E53

Therefore, Eq. (49) is practically expressed in the following series:

$$\begin{array}{cc}{\mathrm{y}}_1& =0;{\mathrm{y}}_2={\mathrm{y}}_0{\mathrm{q}}_0;{\mathrm{y}}_3={\mathrm{y}}_0{\mathrm{q}}_1+{\mathrm{y}}_1{\mathrm{q}}_0={\mathrm{y}}_0{\mathrm{q}}_1;\\ {\mathrm{y}}_4& ={\mathrm{q}}_2{\mathrm{y}}_0+{\mathrm{q}}_1{\mathrm{y}}_1+{\mathrm{q}}_0{\mathrm{y}}_2={\mathrm{y}}_0\left({{\mathrm{q}}_0}^2+{\mathrm{q}}_2\right);\\ {\mathrm{y}}_5& ={\mathrm{y}}_0{\mathrm{q}}_3+{\mathrm{y}}_1{\mathrm{q}}_2+{\mathrm{y}}_2{\mathrm{q}}_1+{\mathrm{y}}_3{\mathrm{q}}_0={\mathrm{y}}_0\left({\text{2q}}_0{\mathrm{q}}_1+{\mathrm{q}}_3\right)\end{array}$$E54

and so on.

And vice versa, q_n’s are expressed as follows:

$$\begin{array}{cc}{\mathrm{q}}_1& =\left(1/{\mathrm{y}}_0\right){\mathrm{y}}_3;{\mathrm{q}}_2=\left({\mathrm{y}}_4/{\mathrm{y}}_0\right)-{\mathrm{q}}_0\left({\mathrm{y}}_2/{\mathrm{y}}_0\right)=\left(1/{\mathrm{y}}_0\right)\left\{{\mathrm{y}}_4-{\mathrm{q}}_0{\mathrm{y}}_2\right)\};\\ {\mathrm{q}}_3& =\left(1/{\mathrm{y}}_0\right)\left\{{\mathrm{y}}_5-\left({\mathrm{y}}_2{\mathrm{y}}_3/y_0\right)-{\mathrm{q}}_0{\mathrm{y}}_3\right\}\end{array}$$E55

and so on.

This method can be immediately applicable in the following way, comparing with Eqs. (38) and (39).

$${\mathrm{q}}_0=-\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^2$$E56

$$\mathrm{y}\left(t\right)={\mathrm{y}}_0\left[1-\left(t^2/2!\right)\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^2+\dots \right]$$E57

where

$${\mathrm{y}}_0=\left({\text{3k}}_{\mathrm{B}}\mathrm{T}/{\mathrm{m}}^{\pm }\right)\left(={{\mathrm{Z}}_{\mathrm{\sigma }}}^{\pm }\left(0\right)\equiv {{\mathrm{Z}}_0}^{\pm }\right)$$E58

7. Fluctuation dissipation theorem on the Laplace transformation of γ(t)

Considering Eqs. (56) and (57), the memory function γ(t) can be taken as the following form:

$$\mathrm{\gamma }\left(t\right)=\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^2\mathrm{f}\left(t\right)$$E59

where f(0) = 1.

The Laplace transformation of (59) in the long wavelength limit is then written as follows:

$$\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)=\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^2{{\int }_0}^{\infty }\mathrm{f}\left(t\right)\mathrm{d}t$$E60

Therefore, we have immediately,

$${{\int }_0}^{\infty }\mathrm{f}\left(t\right)\mathrm{d}t=1/\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)$$E61

On the other hand, the memory function and its Laplace transformation are described as in the following forms, by using the fluctuation dissipation theorem [6, 7, 8, 9],

$$\mathrm{\gamma }\left(t\right)=\left(1/{\text{3μk}}_{\mathrm{B}}\mathrm{T}\right)<{\mathbf{R}}_{\mathrm{i}}\left(t\right){\mathbf{R}}_{\mathrm{j}}\left(0\right)>$$E62

and

$$\stackrel{\sim }{\mathrm{\gamma }}\left(\mathrm{\omega }\right)=\left(1/{\text{3μk}}_{\mathrm{B}}\mathrm{T}\right){{\int }_0}^{\infty }exp\left(-\mathrm{i\omega }t\right)<{\mathbf{R}}_{\mathrm{i}}\left(t\right){\mathbf{R}}_{\mathrm{j}}\left(0\right)>\mathrm{dt}$$E63

The most simplest expression for < R_i(t) R_j(0) > can be taken as in the following form:

$$<{\mathbf{R}}_{\mathrm{i}}\left(t\right){\mathbf{R}}_{\mathrm{j}}\left(0\right)>=<{{\mathbf{R}}_{\mathrm{ij}}}^2>\mathrm{h}\left(t\right)$$E64

where < R_ij² > = < R_i(0) R_j(0)>.

Putting (64) into (62) and using (59), we have

$$\mathrm{\gamma }\left(t\right)=\left(1/{\text{3μk}}_{\mathrm{B}}\mathrm{T}\right)<{{\mathbf{R}}_{\mathrm{ij}}}^2>\mathrm{h}\left(t\right)=\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^2\mathrm{f}\left(t\right)$$E65

This equation gives h(t) ∝ f(t), and if we take both functions are identical, then

$$<{{\mathbf{R}}_{\mathrm{ij}}}^2>=\left(1/{\text{3μk}}_{\mathrm{B}}T\right)\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^2$$E66

Putting this relation into (62), we obtain again the relation (59), which indicates that the assumption, h(t) = f(t), is exactly justified.

Therefore, the general form for the memory function γ(t) is always written in the form of Eq. (59).

8. Former theories of velocity correlation functions in molten salts

Various analytic forms for memory functions were proposed [7, 8, 12, 39, 40, 41, 42, 43] and all these functions are qualitatively useful to obtain the combined velocity correlation functions, although some of these theories cannot predict the result obtained by MD simulation.

For example, if we use an approximate form for the memory function as

$$\begin{array}{ll}{{\mathrm{Z}}_{\mathrm{\sigma }}}^{\pm }\left(\mathrm{t}\right)=& \left({\text{3k}}_{\mathrm{B}}\mathrm{T}/{\mathrm{m}}^{\pm }\right)exp\left\{-\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)t/2\right\}\left[cos\left(\surd \text{3}\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)t/2\right)\right.\\ & \left.+\left\{\left(\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)/2\right)/\left(\surd \text{3}\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)/2\right)\right\}sin\left(\surd \text{3}\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)t/2\right)\right]\\ & =\left({\text{3k}}_{\mathrm{B}}\mathrm{T}/{\mathrm{m}}^{\pm }\right)\left[1-\left(t^2/2!\right)\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^2-\left(t^3/3!\right)\text{3}\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^3/8+\left(\text{over}{t}^4\right)\right]\end{array}$$E67

As shown in our previous results [29], the calculated Z_σ⁺(t) for cation by using Eq. (67) agrees with that of MD simulation [7] qualitatively and semi-quantitatively.

However, the time expansion forms of Z_σ^±(t) are essentially equal to the even powers expansion forms, which contradicts to the expression of (67). It is, therefore, necessary to seek an appropriate memory function which can be expanded as the even powers of the time t, even though the obtained result is numerically very close to the expression of γ(t) = $\stackrel{\sim }{\mathrm{\gamma }}$(0)² exp{−$\stackrel{\sim }{\mathrm{\gamma }}$(0)t}.

9. Application of recursion method for the derivation of γ(t) from Z_σ^±(t)

So far, we are successful to obtain the mutual relation between γ(t) and Z_σ^±(t) within a short time region to satisfy the Langevin equations in molten salts.

There are several works to obtain the auto-velocity correlation functions in monatomic liquids from appropriate memory functions γ(t) [39, 41, 42].

However, it is not known what sorts of model functions are suitable for the combined velocity correlation function Z_σ^±(t) until the present time. In order to elucidate this question, we will try to calculate the coefficients y_m’s of simulated Z_σ^±(t) of molten NaCl in a short time region, and from these the corresponding γ(t) will be obtained.

Previously we have already carried out the MD simulation for the combined velocity correlation functions Z_σ^±(t) [7].

We try two types of power expansion forms in order to fit the combined correlation functions Z_σ^±(t) by MD simulation. One is an arbitrary expansion form given by the even power series of the time t, which is theoretically exact for the combined correlation function. Another one is the series of even and odd powers for higher order terms over t² one. Practical reason for the use of latter case will be given below.

In the case of the utilization of only even powers, it was quite difficult to get to the simulated Z_σ^±(t) even if the power’s number is taken up to 36th order of time t.

On the other hand, we can get an agreement if we use even and odd serial powers over t² up to t⁹. This fact encourages us that the combined velocity correlation functions Z_σ^±(t) in molten systems must be practically analyzed in terms of even and odd powers of the time over t².

Therefore, the method utilizing the odd and even power series has a more rapid convergence for obtaining Z_σ^±(t), in comparison with the method utilizing only even power series.

The fitting parameters, which are equal to y_m’s, are obtained by the non-linear least mean square method as so-called Levenberg-Marquart method [44].

The primary value in this non-linear least mean square method is inferred by utilization of simplex method.

It is inevitable that the coefficients of y_m’s (m = 3, 4, …) are slightly variable because of the termination effect in the expansion form. But we have no difficulty to elucidate γ(t) in an appropriate short time range.

By using these obtained y_m’s, it is immediately possible to obtain q_n’s. And thereafter we can get a fitted curve indicating the curve of γ(t) within a short time region. In this figure, the fitting curve of γ(t) is obtained for the time range of 0 < t < 5.0 × 10⁻¹⁴ seconds, from the expansion form of Z_σ^±(t) up to t¹⁵.

It is therefore emphasized that the utilization of odd terms within the short time region is necessary for the derivation of q_n’s from the y_m’s obtained by MD simulation.

For references, several analytic functional forms describing γ(t) can also be given. The following two-types of functional forms are known as model functions being suitable for the auto-velocity correlation functions in liquids.

$$\left(\mathbf{a}\mathbf{-}\mathbf{1}\right)\mathrm{\gamma }\left(t\right)=\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^{2}sech\left\{\left(\mathrm{\pi }/2\right)\stackrel{\sim }{\mathrm{\gamma }}\left(0\right)t\right\}$$E68

$$\left(\mathbf{a}\mathbf{-}\mathbf{2}\right)\mathrm{\gamma }\left(t\right)=\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^{2}exp\left\{-\left(\mathrm{\pi }/4\right)\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^2{t}^2\right)\}$$E69

The γ(t) is expressed by the form of $\stackrel{\sim }{\mathrm{\gamma }}$(0)²exp{−$\stackrel{\sim }{\mathrm{\gamma }}$(0)t} agrees, at least within the short time region, with that of MD simulation.

However, an inevitable fact is that the theoretical memory function must be an expansion form of only even powers of the time, even though it is numerically close to the exponentially decaying function which includes the odd powers.

Is it possible to get a model function to fit the obtained curve of γ(t) by MD simulation? To answer this question, we have carried out the fitting procedure by using a combination of poly-Gaussian functions [29]. Practically, the following form composed of three kinds of Gaussian functions is good enough to reproduce the obtained curve of γ(t) under the condition of Eq. (61) for molten NaCl at 1100 K,

$$\mathrm{\gamma }\left(t\right)=\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^2\sum _{\mathrm{i}=1}^3{\mathrm{a}}_{\mathrm{i}}exp\left\{-\left(\mathrm{\pi }/4\right){\mathrm{b}}_{\mathrm{i}}\stackrel{\sim }{\mathrm{\gamma }}{\left(0\right)}^2{t}^2\right)\}$$E70

where

$$\sum _{\mathrm{i}=1}^3{\mathrm{a}}_{\mathrm{i}}=1,\text{and}\left\{{\left({\mathrm{b}}_2{\mathrm{b}}_3\right)}^{1/2}{\mathrm{a}}_1+{\left({\mathrm{b}}_3{\mathrm{b}}_1\right)}^{1/2}{\mathrm{a}}_2+{\left({\mathrm{b}}_1{\mathrm{b}}_2\right)}^{1/2}{\mathrm{a}}_3\right\}/{\left({\mathrm{b}}_1{\mathrm{b}}_2{\mathrm{b}}_3\right)}^{1/2}=1$$E71

Using (70) and (71), we could reproduce the obtained curve of γ(t) by MD simulation in molten NaCl at 1100 K. And these are approximated to as {a₁ = 0.2, a₂ = 0.3 and a₃ = 0.5}, which values correspond to the existing fractions of each short range configuration i = 1, i = 2, and i = 3, respectively. And values of {b₁ = 97.50, b₂ = 6.52, and b₃ = 0.38} correspond to their structural decaying speeds, respectively.

According to Figure 2a, the averaged nearest neighbor’s number around the Na⁺ ion is equal to 5.0. Any local coordination numbers around a Na⁺ are possible to be 4, 5, and 6 under the condition of density fluctuation in the liquid state.

It is possible to consider that stable short range configurations seem to be two types. One is the case of cubic structure-type configuration having with the coordination of 6 chlorine ions around the centered sodium ion as shown in Figure 3a, which is similar to the solid type configuration with a sort of lengthen fluctuation of the interionic distance.

Figure 3.

(a) A stable short range configuration of 6 Cl⁻ ions around a Na⁺ ion. (b) Another stable short range configuration of 4 Cl⁻ ions around a Na⁺ ion.

The other is close to a tetrahedral coordination of chlorine ions around the centered sodium ion as shown in Figure 3b.

For simplicity, here we assume that the decaying or releasing of these two types of rather stable short range configurations is nearly the same, then the combined configurational decaying is given by i = 3 and b₃.

On the other hand, there exist two types of rather unstable short range configurations as shown in Figure 4a and b, respectively, in which the surrounded Cl⁻ ions around a Na⁺ ion are spatially asymmetric.

Figure 4.

(a) A rather unstable short range configuration of 5 Cl⁻ ions around a Na⁺ ion. (b) Another unstable short range configuration of 4 Cl⁻ ions around a Na⁺ ion.

Totally, the local configuration types of Cl⁻ ions around a centered Na⁺ ion are listed in Table 1.

Degree of stability	Configuration type
	Coordination of 4 Cl^− ions	Coordination of 5 Cl^− ions	Coordination of 6 Cl^− ions	Existing probability, ai
i = 1	0.2			0.2
i = 2		0.3		0.3
i = 3	0.15		0.35	0.5

Table 1.

Local configuration types of Cl⁻ ions around a centered Na⁺ ion.

10. Discussion and conclusions in the case of molten salts

As shown in the previous section, the combined velocity correlation functions Z_σ^±(t) can be analyzed in terms of odd and even powers over t² in their expansion forms and the corresponding memory function includes the terms of odd and even powers in its expansion form.

In addition, it is emphasized that the γ(t) obtained from the simulated Z_σ⁺(t) agrees completely with that from Z_σ⁻(t). This fact means that the memory functions for cation and anion are identical and Eq. (15) is automatically justified by the present new type of experiment such as computer simulation.

In conclusion, we have newly obtained the mutual relation between the memory function γ(t) and the combined velocity correlation function Z_σ^±(t), by using a recursion method to solve the Langevin equation and it may be applicable for finding a suitable memory function in all liquid matters.

11. Generalized Langevin equation in electrolytic solution

Hereafter, we will consider the strong electrolytic solution composed of N⁺ cations, N⁻ anions and X water molecules in a volume V_M. For simplicity, we take that N⁺ = N⁻ = N and ions charges are equal to z⁺ = − z⁻ = z. Then the number densities of ions and water molecules are equal to n⁺ = n⁻ = n = N/V_M and x = X/V_M, respectively. And furthermore we assume that the dissociation of electrolyte is complete under the condition of N ≪ X.

In the present system, a generalized Langevin equation for the cation (or anion) i under an external field E is written as follows:

$${\mathrm{m}}^{\pm }\mathrm{d}{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right)/\mathrm{d}t=-{\mathrm{m}}^{\pm }{{\int }_0}^t{\mathrm{\gamma }}^{\pm }\left(t-t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(\mathrm{t}’\right)\mathrm{d}t\prime +{\mathrm{z}}^{\pm }\mathrm{e}\left(\mathbf{E}+\mathbf{F}\right)+{\mathrm{z}}^{\pm }\mathrm{e}{\mathbf{\epsilon }}_{\mathrm{i}}\left(t\right)$$E72

where γ^±(t) is the memory function incorporating with the friction force acting on its cation (or anion). F is the induced internal field yielded by the change of ion’s distribution which is resulted from the applying external field E, and ε_i(t) is the random fluctuating force acting on the ion i.

According to Berne and Rice [16], the internal field F induced by the asymmetric ion’s distribution in an ionic melt is expressed as follows:

$$\mathbf{F}=-\mathrm{\delta }·\mathbf{E}=-\text{4πn}/{\text{3k}}_{\mathrm{B}}\mathrm{T}{{\int }_{\mathrm{d}}}^{\infty }\left\{{\mathrm{d\varphi }}^{+-}\left(r\right)/\mathrm{d}r\right\}{\mathrm{ɡ}}^{+-}\left(r\right)r^3\mathrm{d}r·\mathbf{E}$$E73

where ɡ⁺⁻(r) is the pair distribution function between cation and anion, and d is the hard-core contact distance between cation and anion. Hereafter, we will use this result.

If we take E = E₀ e^−iωt, then the ensemble average for v_i^±(t) is written in the following form:

$$<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right)>=Re{\mathrm{\mu }}^{\pm }\left(\mathrm{\omega }\right){\mathrm{z}}^{\pm }\mathrm{e}{\mathbf{E}}_0{\mathrm{e}}^{-\mathrm{i\omega }t}$$E74

Inserting (74) into (72) and taking ensemble average under the assumption of <ε_i(t) > = 0, we have

$${\mathrm{m}}^{\pm }<\mathrm{d}{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right)>/\mathrm{d}t=-{\mathrm{m}}^{\pm }{{\int }_0}^{\mathrm{t}}{\mathrm{\gamma }}^{\pm }\left(t-t\prime \right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\prime \right)\mathrm{d}t\prime +{\mathrm{z}}^{\pm }\mathrm{e}\left(1–\mathrm{\delta }\right)\mathbf{E}$$E75

Therefore,

$${\mathrm{\mu }}^{\pm }\left(\mathrm{\omega }\right)=\left(1–\mathrm{\delta }\right)/{\mathrm{m}}^{\pm }\left\{-\mathrm{i\omega }+\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(\mathrm{\omega }\right)\right\}$$E76

where

$$\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(\mathrm{\omega }\right)={{\int }_0}^{\infty }{\mathrm{\gamma }}^{\pm }\left(t\right){\mathrm{e}}^{\mathrm{i\omega }t}\mathrm{d}t$$E77

The dc current density j^± is then written as follows:

$${\mathbf{j}}^{\pm }={\mathrm{nz}}^{\pm }\mathrm{e}<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right)>{}_{\mathrm{\omega }=0}={\mathrm{nz}}^2{\mathrm{e}}^2{\mathrm{\mu }}^{\pm }\left(0\right){\mathbf{E}}_0={\mathrm{nz}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right){\mathbf{E}}_0/{\mathrm{m}}^{\pm }\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(0\right)$$E78

On the other hand, j^± is expressed as j^± = σ^±E₀, where σ^± being equal to the partial conductivity for cation or anion. Therefore, σ^± is written as follows:

$${\mathrm{\sigma }}^{\pm }={\mathrm{nz}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right)/{\mathrm{m}}^{\pm }\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(0\right)$$E79

The Laplace transformation of the memory function in the long wavelength limit $\stackrel{\sim }{\mathrm{\gamma }}$^±(0) in Eq. (79) will be obtained in later section.

In the next section, we will discuss velocity correlation functions.

12. Linear response theory for electrolytic solutions

Eq. (79) is also obtainable from the following simplified Langevin equation:

$${\mathrm{m}}^{\pm }\mathrm{d}{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right)/\mathrm{d}t=-{\mathrm{m}}^{\pm }\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(0\right){{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right)+{\mathrm{z}}^{\pm }\mathrm{e}\left(\mathbf{E}+\mathbf{F}\right)+{\mathrm{z}}^{\pm }\mathrm{e}{\mathbf{\epsilon }}_{\mathrm{i}}\left(t\right)$$E80

Its derivation can be easily seen in a standard book of statistical physics.

Starting from Eq. (80) with an infinitesimal external field E, it is also easily obtainable the following Kubo-Green formulae for the partial conductivities σ⁺ and σ⁻ [6, 7, 28]:

$${\mathrm{\sigma }}^{+}=\left(1/{\text{3k}}_{\mathrm{B}}\mathrm{T}\right){{\int }_0}^{\infty }<{\mathbf{j}}^{+}\left(t\right)·\mathbf{j}\left(0\right)>\mathrm{d}t$$E81

and

$${\mathrm{\sigma }}^{-}=\left(1/{\text{3k}}_{\mathrm{B}}\mathrm{T}\right){{\int }_0}^{\infty }<{\mathbf{j}}^{-}\left(t\right)·\mathbf{j}\left(0\right)>\mathrm{d}t$$E82

where the current densities j^±(t) and j(t) are defined by the following expressions:

$${\mathbf{j}}^{\pm }\left(t\right)={\mathrm{z}}^{\pm }\mathrm{e}\sum _{\mathrm{n}}{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(t\right)\text{and}\mathbf{j}\left(t\right)={\mathbf{j}}^{+}\left(t\right)+{\mathbf{j}}^{-}\left(t\right)$$E83

In order to obtain the partial conductivities based on Eqs. (81) and (82), it is necessary to study the velocity correlation functions, < v_i⁺(t) v_j⁺(0)>, < v_k⁻(t) v_l⁻(0) > and < v_i⁺(t) v_k⁻(0)>.

In the next section, we will discuss velocity correlation functions described in terms of inter-molecular (or ionic) potentials and pair distribution functions in order to obtain the $\stackrel{\sim }{\mathrm{\gamma }}$^±(0).

13. Short time expansion of velocity correlation functions in electrolytic solutions

The short time expansion of velocity correlation function, < v_i⁺(t) v_j⁺(0) > for cation is written as

$$<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(t\right){{\mathbf{v}}_{\mathrm{j}}}^{+}\left(0\right)>=<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right){{\mathbf{v}}_{\mathrm{j}}}^{+}\left(0\right)>+t^2/2!<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right){{{\mathbf{v}}^{¨}}_{\mathrm{j}}}^{+}\left(0\right)>+\left(\text{higher order over}{t}^4\right)$$E84

In the present aqueous solution of electrolyte, the total Hamiltonian of the system is written as follows:

$$\mathrm{H}=\sum _{\mathrm{i}=1}^{\mathrm{N}+}{{\mathbf{p}}_{\mathrm{i}}}^{+2}/{\text{2m}}^{+}+\sum _{\mathrm{k}=1}^{\mathrm{N}-}{{\mathbf{p}}_{\mathrm{i}}}^{-2}/2m^{-}+\sum _{\mathrm{q}=1}^{\mathrm{X}}{{\mathbf{p}}_{\mathrm{q}}}^{\text{w2}}/{\text{2m}}^{\mathrm{w}}+\mathrm{V}$$E85

where

$$\begin{array}{ll}\mathrm{V}=& \sum _{\mathrm{i}\ne \mathrm{j}}^{\mathrm{N}+}{\mathrm{\varphi }}^{++}\left(|{{\mathbf{r}}_{\mathrm{i}}}^{+}-{{\mathbf{r}}_{\mathrm{j}}}^{+}|\right)+\sum _{\mathrm{k}\ne \mathrm{l}}^{\mathrm{N}-}{\mathrm{\varphi }}^{--}\left(|{{\mathbf{r}}_{\mathrm{k}}}^{—}-{{\mathbf{r}}_{\mathrm{l}}}^{-}|\right)+\sum _{\mathrm{i},\mathrm{k}}^{\mathrm{N}+,\mathrm{N}-}{\mathrm{\varphi }}^{+-}\left(|{{\mathbf{r}}_{\mathrm{i}}}^{+}-{{\mathbf{r}}_{\mathrm{k}}}^{-}|\right)\\ & +\sum _{\mathrm{i},\mathrm{q}}^{\mathrm{N}+,\mathrm{X}}{\mathrm{\varphi }}^{+\mathrm{w}}\left(|{{\mathbf{r}}_{\mathrm{i}}}^{+}-{{\mathbf{r}}_{\mathrm{q}}}^{\mathrm{w}}|\right)+\sum _{\mathrm{k},\mathrm{q}}^{\mathrm{N}-,\mathrm{X}}{\mathrm{\varphi }}^{-\mathrm{w}}\left(|{{\mathbf{r}}_{\mathrm{k}}}^{-}-{{\mathbf{r}}_{\mathrm{q}}}^{\mathrm{w}}|\right)+\sum _{\mathrm{q},\mathrm{s}}^{\mathrm{X}}{\mathrm{\varphi }}^{\mathrm{w}\mathrm{w}}\left(|{{\mathbf{r}}_{\mathrm{q}}}^{\mathrm{w}}-{{\mathbf{r}}_{\mathrm{s}}}^{\mathrm{w}}|\right)\end{array}$$E86

Since the water molecule is not spherical in its molecular configuration, it is difficult to define the position of r_q^w. However, we tentatively assume that its position is located at the center of oxygen atom in the H₂O molecule.

From the Poisson’s equation of motion,

$$\begin{array}{ll}{\mathbf{p}}_i^{+}\ddot{\mathbf{p}}{}_i^{+}=& -\sum _{\mathrm{j}=1}^{\mathrm{N}+}\left({{\mathbf{p}}_{\mathrm{i}}}^{+}{{\mathbf{p}}_{\mathrm{j}}}^{+}/{\mathrm{m}}^{+}\right){\partial }^2\mathrm{V}/\partial {{\mathbf{r}}_{\mathrm{i}}}^{+}\partial {{\mathbf{r}}_{\mathrm{j}}}^{+}-\sum _{\mathrm{k}=1}^{\mathrm{N}-}\left({{\mathbf{p}}_{\mathrm{i}}}^{+}{{\mathbf{p}}_{\mathrm{k}}}^{-}/{\mathrm{m}}^{-}\right){\partial }^2\mathrm{V}/\partial {{\mathbf{r}}_{\mathrm{k}}}^{-}\partial {{\mathbf{r}}_{\mathrm{i}}}^{+}\\ & -\sum _{\mathrm{q}=1}^{\mathrm{X}}\left({{\mathbf{p}}_{\mathrm{i}}}^{+}{{\mathbf{p}}_{\mathrm{q}}}^{\mathrm{w}}/{\mathrm{m}}^{\mathrm{w}}\right){\partial }^2\mathrm{V}/\partial {{\mathbf{r}}_{\mathrm{q}}}^{\mathrm{w}}\partial {{\mathbf{r}}_{\mathrm{i}}}^{+}\end{array}$$E87

and

$$\begin{array}{ll}{\mathbf{p}}_i^{+}\ddot{\mathbf{p}}{}_{\mathrm{i}\ne \mathrm{i}’}^{+}=& -\sum _{\mathrm{j}}^{\mathrm{N}+}\left({{\mathbf{p}}_{\mathrm{i}’}}^{+}{{\mathbf{p}}_{\mathrm{j}}}^{+}/{\mathrm{m}}^{+}\right){\partial }^2\mathrm{V}/\partial {{\mathbf{r}}_{\mathrm{i}}}^{+}\partial {{\mathbf{r}}_{\mathrm{j}}}^{+}-\sum _{\mathrm{k}=1}^{\mathrm{N}-}\left({{\mathbf{p}}_{\mathrm{i}}}^{+}{{\mathbf{p}}_{\mathrm{k}}}^{-}/{\mathrm{m}}^{-}\right){\partial }^2\mathrm{V}/\partial {{\mathbf{r}}_{\mathrm{k}}}^{-}\partial {{\mathbf{r}}_{\mathrm{i}}}^{+}\\ & -\sum _{\mathrm{q}=1}^{\mathrm{X}}\left({{\mathbf{p}}_{\mathrm{i}}}^{+}{{\mathbf{p}}_{\mathrm{q}}}^{\mathrm{w}}/{\mathrm{m}}^{\mathrm{w}}\right){\partial }^2\mathrm{V}/\partial {{\mathbf{r}}_{\mathrm{q}}}^{\mathrm{w}}\partial {{\mathbf{r}}_{\mathrm{i}}}^{+}\end{array}$$E88

Since the second derivative of the potential term V is independent for the product of momenta, all other terms other than i = j in (87) must vanish on averaging. And in a similar way, the meaningful terms of (88) for averaging must be also equal to the case i ≠ i’ = j. Therefore, taking the ensemble averages for (87) and (88), we have

$$<{\mathbf{p}}_i^{+}\ddot{\mathbf{p}}{}_j^{+}>=<{\mathbf{p}}_i^{+}\ddot{\mathbf{p}}{}_i^{+}>+<{\mathbf{p}}_i^{+}\ddot{\mathbf{p}}{}_{\mathrm{j}\ne \mathrm{i}}^{+}>=-{\mathrm{k}}_{\mathrm{B}}\mathrm{T}\left\{\mathrm{n}<{\mathrm{\varphi }}^{+-}>+\mathrm{x}<{\mathrm{\varphi }}^{+\mathrm{w}}>\right\}$$E89

where

$$<{\mathrm{\varphi }}^{+-}>={{\int }_0}^{\infty }\left\{{\partial }^2{\mathrm{\varphi }}^{+-}\left(r\right)\partial r^2+\left(2/r\right)\partial {\mathrm{\varphi }}^{+-}\left(r\right)\partial r\right\}{\mathrm{ɡ}}^{+-}\left(\mathrm{r}\right)\text{4π}{r}^2\mathrm{d}r$$E90

and

$$<{\mathrm{\varphi }}^{+\mathrm{w}}>={{\int }_0}^{\infty }\left\{{\partial }^2{\mathrm{\varphi }}^{+\mathrm{w}}\left(r\right)\partial r^2+\left(2/\mathrm{r}\right)\partial {\mathrm{\varphi }}^{+\mathrm{w}}\left(r\right)\partial r\right\}{\mathrm{ɡ}}^{+\mathrm{w}}\left(r\right)\text{4π}{r}^2\mathrm{d}r$$E91

In this equation, ɡ^+w(r) is the pair distribution function between cation and water molecule.

It is emphasized that there is no contribution from ϕ⁺⁺(r) in Eq. (89) because of the cancelation by the terms of i = j and i ≠ j in < ${\mathbf{p}}_i^{+}\ddot{\mathbf{p}}{}_j^{+}$> [7].

Insertion of (89) into (84) gives us the following form:

$$\begin{array}{ll}{\mathrm{Z}}^{+}\left(t\right)& \equiv <{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(t\right){{\mathbf{v}}_{\mathrm{j}}}^{+}\left(0\right)>\\ & =<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right){{\mathbf{v}}_{\mathrm{j}}}^{+}\left(0\right)>-\left(t^2/2!\right){\mathrm{k}}_{\mathrm{B}}\mathrm{T}\left\{\left(\mathrm{n}<{\mathrm{\varphi }}^{+-}>+\mathrm{x}<{\mathrm{\varphi }}^{+\mathrm{w}}>\right)/{\mathrm{m}}^{+}{\mathrm{m}}^{+}\right\}\\ & +\left(\text{higher order over}{t}^4\right)\end{array}$$E92

In a similar way, the term < ${\mathbf{p}}_i^{+}\ddot{\mathbf{p}}{}_k^{-}$> can be described as follows:

$$<{\mathbf{p}}_i^{+}\ddot{\mathbf{p}}{}_k^{-}>={\mathrm{k}}_{\mathrm{B}}\mathrm{T}\mathrm{n}<{\mathrm{\varphi }}^{+-}>$$E93

Using this relation, the distinct velocity correlation function is written as follows:

$$\begin{array}{ll}<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(t\right){{\mathbf{v}}_{\mathrm{k}}}^{-}\left(0\right)>& =<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right){{\mathbf{v}}_{\mathrm{k}}}^{-}\left(0\right)>+t^2/2!<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right){{{\mathbf{v}}^{¨}}_{\mathrm{k}}}^{-}\left(0\right)>\\ & +\left(\text{higher order over}{t}^4\right)\\ & =<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right){{\mathbf{v}}_{\mathrm{k}}}^{-}\left(0\right)>+\left(t^2/2!\right){\mathrm{k}}_{\mathrm{B}}\mathrm{T}\mathrm{n}<{\mathrm{\varphi }}^{+-}>/{\mathrm{m}}^{+}{\mathrm{m}}^{-}\\ & +\left(\text{higher order over}{t}^4\right)\end{array}$$E94

Using (92) and (94), the combined velocity correlation function Z_σ⁺(t)(= < j⁺(t)j(0)>/n²z²e²) incorporation with the partial conductivity σ⁺ is therefore expressed as follows:

$$\begin{array}{ll}{{\mathrm{Z}}_{\mathrm{\sigma }}}^{+}\left(t\right)& \equiv <{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(t\right){{\mathbf{v}}_{\mathrm{j}}}^{+}\left(0\right)>-<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(t\right){{\mathbf{v}}_{\mathrm{k}}}^{-}\left(0\right)>\\ & =\left({\text{3k}}_{\mathrm{B}}\mathrm{T}/{\mathrm{m}}^{+}\right)\left\{1–\left(t^2/2!\right)\left(\mathrm{n}<{\mathrm{\varphi }}^{+-}>/\text{3μ}+\mathrm{x}<{\mathrm{\varphi }}^{+\mathrm{w}}>/{\text{3m}}^{+}\right)\right.\\ & \left.+\left(\text{higher order over}{t}^4\right)\right\}\end{array}$$E95

where μ is equal to the reduced mass of m⁺ and m⁻. In deriving (95), we have assumed the initial conditions as follows:

$$<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right){{\mathbf{v}}_{\mathrm{j}}}^{+}\left(0\right)>=<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right){{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right)>=\left({\text{3k}}_{\mathrm{B}}\mathrm{T}/{\mathrm{m}}^{+}\right)\text{and}<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right){{\mathbf{v}}_{\mathrm{k}}}^{-}\left(0\right)>=0$$E96

These initial conditions are confirmed by our own molecular dynamic simulation, which will be shown in the later section. In a similar way, we have

$$\begin{array}{ll}{{\mathrm{Z}}_{\mathrm{\sigma }}}^{-}\left(t\right)& \equiv <{{\mathbf{v}}_{\mathrm{i}}}^{-}\left(t\right){{\mathbf{v}}_{\mathrm{j}}}^{-}\left(0\right)>-<{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right){{\mathbf{v}}_{\mathrm{j}}}^{-}\left(0\right)>\\ & =\left({\text{3k}}_{\mathrm{B}}\mathrm{T}/{\mathrm{m}}^{-}\right)\left\{1–\left(t^2/2!\right)\left(\mathrm{n}<{\mathrm{\varphi }}^{+-}>/\text{3μ}+\mathrm{x}<{\mathrm{\varphi }}^{-\mathrm{w}}>/{\text{3m}}^{-}\right)\right.\\ & \left.+\left(\text{higher order over}{t}^4\right)\right\}\end{array}$$E97

where

$$<{\mathrm{\varphi }}^{-\mathrm{w}}>={{\int }_0}^{\infty }\left\{{\partial }^2{\mathrm{\varphi }}^{-\mathrm{w}}\left(r\right)\partial r^2+\left(2/r\right)\partial {\mathrm{\varphi }}^{-\mathrm{w}}\left(r\right)\partial r\right\}{\mathrm{ɡ}}^{-\mathrm{w}}\left(r\right)\text{4π}{r}^2\mathrm{d}r$$E98

ɡ^−w(r) of this equation means the pair distribution function between anion and water molecule. And it is also emphasized that the contribution from ϕ⁻⁻(r) to < v_i⁻(t) v_j⁻(0) > is also vanished to be zero.

It is impossible to obtain the partial conductivities by the insertion of (95) and (97) into (81) and (82), because we knew only the terms up to t² in their explicit forms. However, these equations can be utilized for the derivation of $\stackrel{\sim }{\mathrm{\gamma }}$^±(0) as shown in the next section.

14. Derivation of $\stackrel{\sim }{\mathrm{\gamma }}$^±(0) in electrolytic solutions

According to the fluctuation dissipation theorem applied for the present system with the condition of no external field or of infinitesimal external field, the Laplace transformation of the memory function γ^±(t) and that of the ensemble average of time correlation function for the fluctuating random force <ε_i(t)ε_i(0) > have the following relation [25, 28]:

$$\begin{array}{ll}{\mathrm{m}}^{\pm 2}<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(0\right)·{{\mathbf{v}}_{\mathrm{j}}}^{\pm }\left(0\right)>\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(\mathrm{\omega }\right)& ={\mathrm{m}}^{\pm 2}<{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(0\right)·{{\mathbf{v}}_{\mathrm{i}}}^{\pm }\left(0\right)>\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(\mathrm{\omega }\right)\\ & ={\mathrm{z}}^2{\mathrm{e}}^2{{\int }_0}^{\infty }<{\mathbf{\epsilon }}_{\mathrm{i}}\left(t\right){\mathbf{\epsilon }}_{\mathrm{i}}\left(0\right)>{\mathrm{e}}^{\mathrm{i\omega }t}\mathrm{d}t\end{array}$$E99

The fluctuation dissipation theorem tells us the following relation:

$$\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(\mathrm{\omega }\right)=\left(1/3{\mathrm{m}}^{\pm }{\mathrm{k}}_{\mathrm{B}}\mathrm{T}\right)\left({\mathrm{z}}^2{\mathrm{e}}^2\right){{\int }_0}^{\infty }<{\mathbf{\epsilon }}_{\mathrm{i}}\left(t\right){\mathbf{\epsilon }}_{\mathrm{i}}\left(0\right)>{\mathrm{e}}^{\mathrm{i\omega }t}\mathrm{d}t$$E100

In the long wavelength limit, this relation is expressed by

$$\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(0\right)=\left(1/3{\mathrm{m}}^{\pm }{\mathrm{k}}_{\mathrm{B}}\mathrm{T}\right)\left({\mathrm{z}}^2{\mathrm{e}}^2\right){{\int }_0}^{\infty }<{\mathbf{\epsilon }}_{\mathrm{i}}\left(t\right){\mathbf{\epsilon }}_{\mathrm{i}}\left(0\right)>\mathrm{d}t$$E101

Let us go back to the memory function γ^±(t) and assume a combined exponential decay functions for it as follows, although this assumption is not exactly consistent with Eq. (84), but technically acceptable one as discussed in the case of molten salt [29],

$${\mathrm{\gamma }}^{\pm }\left(t\right)={{\mathrm{\gamma }}_0}^{\pm }\left(t\right)(\left\{={{\mathrm{\gamma }}_{00}}^{\pm }exp\left(-{{\mathrm{\beta }}_0}^{\pm }t\right)\right\}+{{\mathrm{\gamma }}_1}^{\pm }\left(t\right)\left\{={{\mathrm{\gamma }}_{01}}^{\pm }exp\left(-{{\mathrm{\beta }}_1}^{\pm }t\right)\right\}$$E102

In this equation, the pre-exponential factor γ₀₀^± is subject to the interactions between the central ion and surrounding water molecules. The decaying constants are related to the time dependence of its configuration change. The pre-exponential factor, γ₀₁^±, is equal to the magnitude of memory function at t = 0 in respect to the friction force acting on the central cation or anion caused by interactions between its central ion and neighboring ions. In other words, the first term on the right hand side of this equation means the case of dilute limit of electrolytic solution and the second one is equal to the effective friction caused by the addition of a moderate amount of electrolyte. Therefore, the first term is related to either <ϕ^+w > or < ϕ^−w >, while the second one is related to the term <ϕ⁺⁻ > .

Using (94) and (96), γ₀₀^± and γ₀₁^± are expressed as follows:

$${{\mathrm{\gamma }}_{00}}^{\pm }+{{\mathrm{\gamma }}_{01}}^{\pm }=-\ddot{\mathrm{Z}}{}_{\mathrm{\sigma }}^{\pm }\left(0\right)/\mathrm{Z}{}_{\mathrm{\sigma }}^{\pm }\left(0\right)=\mathrm{x}<{\mathrm{\varphi }}^{\pm \mathrm{w}}>/{\text{3m}}^{\pm }+\mathrm{n}<{\mathrm{\varphi }}^{+-}>/\text{3μ}$$E103

In the dilute limit of n ≪ x, we have

$$\left({{\mathrm{\gamma }}_{00}}^{\pm }\right)=\mathrm{x}<{\mathrm{\varphi }}^{\pm \mathrm{w}}>/{\text{3m}}^{\pm }$$E104

And then we have

$$\left({{\mathrm{\gamma }}_{01}}^{\pm }\right)=\mathrm{n}<{\mathrm{\varphi }}^{+-}>/\text{3μ}$$E105

At the dilution limit of electrolyte where the contribution of γ₁^±(t) can be neglected, the Laplace transformation of γ₀^±(t) in the long wavelength limit is then described as follows:

$$\stackrel{\sim }{\mathrm{\gamma }}{}_0^{\pm }\left(0\right)={{\mathrm{\gamma }}_{00}}^{\pm }/{{\mathrm{\beta }}_0}^{\pm }=\mathrm{x}<{\mathrm{\varphi }}^{\pm \mathrm{w}}>/{\text{3m}}^{\pm }{{\mathrm{\beta }}_0}^{\pm }=\left(1/3{\mathrm{m}}^{\pm }{\mathrm{k}}_{\mathrm{B}}\mathrm{T}\right)\left({\mathrm{z}}^2{\mathrm{e}}^2\right){{\int }_0}^{\infty }<{\mathbf{\epsilon }}_{\text{i0}}\left(t\right){\mathbf{\epsilon }}_{\text{i0}}\left(0\right)>\mathrm{d}t$$E106

where the auto-correlation function of random fluctuating force <ε_i0(t)ε_i0(0) > is only related to either <ϕ^+w > or < ϕ^−w >.

As seen in Eq. (79), the Laplace transformation of memory function in the long wavelength limit, $\stackrel{\sim }{\mathrm{\gamma }}{}_0^{\pm }$(0), corresponds to effective friction constants for cation and anion, which means that the auto-correlation function of random fluctuating force.

<ε_i0(t) ε_i0(0) > may be represented by an exponential decaying function with the time constant of $\stackrel{\sim }{\mathrm{\gamma }}{}_0^{\pm }$(0) as follows:

$$<{\mathbf{\epsilon }}_{\text{i0}}\left(t\right){\mathbf{\epsilon }}_{\text{i0}}\left(0\right)>=<{\mathbf{\epsilon }}_{\text{i0}}\left(0\right){\mathbf{\epsilon }}_{\text{i0}}\left(0\right)>exp\left\{-\stackrel{\sim }{\mathrm{\gamma }}{}_0^{\pm }\left(0\right)t\right\}$$E107

Insertion of (107) into (106) gives us

$$<{\mathbf{\epsilon }}_{\text{i0}}\left(0\right){\mathbf{\epsilon }}_{\text{i0}}\left(0\right)>=3{\mathrm{m}}^{\pm }{\mathrm{k}}_{\mathrm{B}}\mathrm{T}{\left(\stackrel{\sim }{\mathrm{\gamma }}{}_0^{\pm }\left(0\right)\right)}^2$$E108

Therefore, we obtain

$$\mathrm{\gamma }{}_0^{\pm }\left(t\right)={\left(\stackrel{\sim }{\mathrm{\gamma }},{}_0^{\pm },\left(0\right)\right)}^{2}exp\left(-\stackrel{\sim }{\mathrm{\gamma }}{}_0^{\pm }\left(0\right)t\right)$$E109

Compare (106) and (109) we have

$${{\mathrm{\beta }}_0}^{\pm }=\stackrel{\sim }{\mathrm{\gamma }}{}_0^{\pm }\left(0\right)={\left({{\mathrm{\gamma }}_{00}}^{\pm }\right)}^{1/2}={\left(\mathrm{x}<{\mathrm{\varphi }}^{\pm \mathrm{w}}>/{\text{3m}}^{\pm }\right)}^{1/2}$$E110

By the analogy with this relation, we can infer the following relation:

$${{\mathrm{\beta }}_1}^{\pm }=\stackrel{\sim }{\mathrm{\gamma }}{}_1^{\pm }\left(0\right)={\left({{\mathrm{\gamma }}_{01}}^{\pm }\right)}^{1/2}={\left(\mathrm{n}<{\mathrm{\varphi }}^{\pm \mathrm{w}}>/\text{3μ}\right)}^{1/2}$$E111

Therefore, Eq. (102) is explicitly written as follows:

$$\begin{array}{ll}{\mathrm{\gamma }}^{\pm }\left(t\right)& ={{\mathrm{\gamma }}_0}^{\pm }\left(t\right)+{{\mathrm{\gamma }}_1}^{\pm }\left(t\right)=\left(\mathrm{x}<{\mathrm{\varphi }}^{\pm \mathrm{w}}>/{\text{3m}}^{\pm }\right)exp\left[-{\left(\mathrm{x}<{\mathrm{\varphi }}^{\pm \mathrm{w}}>/{\text{3m}}^{\pm }\right)}^{1/2}t\right]\\ & +\left(\mathrm{n}<{\mathrm{\varphi }}^{+-}>/\text{3μ}\right)exp\left[-{\left(\mathrm{n}<{\mathrm{\varphi }}^{+-}>/\text{3μ}\right)}^{1/2}t\right]\end{array}$$E112

And the Laplace transformation of this equation in the long wavelength limit is equal to

$$\stackrel{\sim }{\mathrm{\gamma }}{}^{\pm }\left(0\right)=\left[{\left(\mathrm{x}<{\mathrm{\varphi }}^{\pm \mathrm{w}}>/{\text{3m}}^{\pm }\right)}^{1/2}\right]+\left[{\left(\mathrm{n}<{\mathrm{\varphi }}^{+-}>/\text{3μ}\right)}^{1/2}\right]$$E113

15. Partial conductivities σ⁺ and σ⁻

Putting Eq. (113) into (79), we obtain the formulae of the partial conductivities, σ⁺ and σ⁻, which are expressed in terms of the pair distribution functions and pair potentials as follows [28],

$${\mathrm{\sigma }}^{+}={\mathrm{nz}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right)/{\mathrm{m}}^{+}\stackrel{\sim }{\mathrm{\gamma }}{}^{+}\left(0\right)={\mathrm{nz}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right)/{\mathrm{m}}^{+}\left[{\left(\mathrm{x}<{\mathrm{\varphi }}^{+\mathrm{w}}>/{\text{3m}}^{+}\right)}^{1/2}+{\left(\mathrm{n}<{\mathrm{\varphi }}^{+-}>/\text{3μ}\right)}^{1/2}\right]$$E114

and

$${\mathrm{\sigma }}^{-}={\mathrm{nz}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right)/{\mathrm{m}}^{-}\stackrel{\sim }{\mathrm{\gamma }}{}^{-}\left(0\right)={\mathrm{nz}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right)/{\mathrm{m}}^{-}\left[{\left(\mathrm{x}<{\mathrm{\varphi }}^{-\mathrm{w}}>/{\text{3m}}^{-}\right)}^{1/2}+{\left(\mathrm{n}<{\mathrm{\varphi }}^{+-}>/\text{3μ}\right)}^{1/2}\right]$$E115

If the concentration c is defined as the number of moles of electrolyte in the unit volume (actually taken as 1 cc), then the number density n is equal to cN_A, where N_A being the Avogadro’s number. Then, the partial conductivities, σ⁺ and σ⁻, are written as follows:

$${\mathrm{\sigma }}^{+}={\mathrm{ne\mu }}^{+}=\mathrm{c}{\mathrm{N}}_{\mathrm{A}}{\mathrm{z}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right)/{\mathrm{m}}^{+}\left[{\left(\mathrm{x}<{\mathrm{\varphi }}^{+\mathrm{w}}>/{\text{3m}}^{+}\right)}^{1/2}+{\left({\mathrm{cN}}_{\mathrm{A}}<{\mathrm{\varphi }}^{+-}>/\text{3μ}\right)}^{1/2}\right]$$E116

and

$${\mathrm{\sigma }}^{-}={\mathrm{ne\mu }}^{-}=\mathrm{c}{\mathrm{N}}_{\mathrm{A}}{\mathrm{z}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right)/{\mathrm{m}}^{-}\left[{\left(\mathrm{x}<\mathrm{\varphi }{}^{-\mathrm{w}}>/{\text{3m}}^{-}\right)}^{1/2}+{\left({\mathrm{cN}}_{\mathrm{A}}<{\mathrm{\varphi }}^{+-}>/\text{3μ}\right)}^{1/2}\right]$$E117

In these equations, μ⁺ and μ⁻ are called as the mobility of cation and anion.

The partial molar conductance Λ⁺ and Λ⁻ are defined as Λ^± = σ^±/c. Then the total conductance Λ_c is written as follows:

$$\begin{array}{ll}{\mathrm{\Lambda }}_{\mathrm{c}}=& {\mathrm{\Lambda }}^{+}+{\mathrm{\Lambda }}^{-}={\mathrm{N}}_{\mathrm{A}}{\mathrm{z}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right)\times \left[1/{\mathrm{m}}^{+}\left\{{\left(\mathrm{x}<{\mathrm{\varphi }}^{+\mathrm{w}}>/{\text{3m}}^{+}\right)}^{1/2}+{\left({\mathrm{cN}}_{\mathrm{A}}<{\mathrm{\varphi }}^{+-}>/\text{3μ}\right)}^{1/2}\right\}\right.\\ & \left.+1/{\mathrm{m}}^{-}\left\{{\left(\mathrm{x}<\mathrm{\varphi }{}^{-\mathrm{w}}>/{\text{3m}}^{-}\right)}^{1/2}+{\left({\mathrm{cN}}_{\mathrm{A}}<{\mathrm{\varphi }}^{+-}>/\text{3μ}\right)}^{1/2}\right\}\right]\end{array}$$E118

Under the condition of n(=cN_A) ≪ x, they are approximated to as follows:

$${\mathrm{\Lambda }}^{+}={\mathrm{N}}_{\mathrm{A}}{\mathrm{z}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right){\left(3/{\mathrm{xm}}^{+}<{\mathrm{\varphi }}^{+\mathrm{w}}>\right)}^{1/2}\left[1–{\left({\mathrm{cN}}_{\mathrm{A}}{\mathrm{m}}^{+}<{\mathrm{\varphi }}^{+-}>/\mathrm{\mu x}<{\mathrm{\varphi }}^{+\mathrm{w}}>\right)}^{1/2}\right]$$E119

and

$${\mathrm{\Lambda }}^{-}={\mathrm{N}}_{\mathrm{A}}{\mathrm{z}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right){\left(3/{\mathrm{xm}}^{-}<\mathrm{\varphi }{}^{-\mathrm{w}}>\right)}^{1/2}\left[1–{\left({\mathrm{cN}}_{\mathrm{A}}{\mathrm{m}}^{-}<{\mathrm{\varphi }}^{+-}>/\mathrm{\mu x}<\mathrm{\varphi }{}^{-\mathrm{w}}>\right)}^{1/2}\right]$$E120

From Eqs. (119) and (120), we have a form of Λ_c = (Λ⁺ + Λ⁻) ≃ Λ₀ + Λ₁– kc^1/2. Λ₀ and k are written as follows:

$${\mathrm{\Lambda }}_0={\mathrm{N}}_{\mathrm{A}}{\mathrm{z}}^2{\mathrm{e}}^2\left\{{\left(3/{\mathrm{xm}}^{+}<{\mathrm{\varphi }}^{+\mathrm{w}}>\right)}^{1/2}+{\left(3/{\mathrm{xm}}^{-}<\mathrm{\varphi }{}^{-\mathrm{w}}>\right)}^{1/2}\right\}$$E121

$${\mathrm{\Lambda }}_1={\mathrm{N}}_{\mathrm{A}}{\mathrm{\delta z}}^2{\mathrm{e}}^2\left\{{\left(3/{\mathrm{xm}}^{+}<{\mathrm{\varphi }}^{+\mathrm{w}}>\right)}^{1/2}+{\left(3/{\mathrm{xm}}^{-}<\mathrm{\varphi }{}^{-\mathrm{w}}>\right)}^{1/2}\right\}$$E122

and

$$\mathrm{k}={\mathrm{N}}_{\mathrm{A}}{\mathrm{z}}^2{\mathrm{e}}^2\left(1–\mathrm{\delta }\right){\left({\text{3N}}_{\mathrm{A}}<{\mathrm{\varphi }}^{+-}>/\mathrm{\mu }\right)}^{1/2}\left\{\right(1/\mathrm{x}<{\mathrm{\varphi }}^{+\mathrm{w}}>+\left(1/\mathrm{x}<{\mathrm{\varphi }}^{-\mathrm{w}}>\right)\}$$E123

As seen in these expressions, Λ₀ means the conductance in the dilution limit of electrolyte and Λ₁ is the correction term appeared by the so-called relaxation effect. The last term kc^1/2 is composed of the so-called electrophoretic effect and the combined term of both effects.

In the case of a moderate concentration of electrolyte, in particular, of relatively weak electrolyte, we have to take account of the degree of association between the opposite ions into the expression for the partial conductivities.

16. Pair potentials in electrolytic solution

A number of research works to obtain the model potentials in electrolytic solutions have been presented since the Debye-Hückel theory [18]. In particular, various qualified model potentials, which satisfy the experimental data such as the hydration free energy and the enthalpies in condensed and gas phases, have recently been proposed in order to carry out the molecular dynamic simulation. It is not our intension to compare or evaluate these potentials and therefore we like to refer only some of these for our interests [24, 25, 26, 27, 45]. It may be possible to estimate these potentials by using wave mechanical approach. In fact the ion-water molecule interactions were obtained by such an elaborating method [46, 47, 48].

The essential point for these model potentials in electrolytic solutions is that the dielectric character should be concerned. According to Sack [49], the water-molecules around the constituent ion are strongly oriented and the ion’s orientating ability to neighboring water-molecules decreases with increasing of the distance between the ion and those water-molecules. Oka [50] also estimated the change of effective dielectric constant as a function of distance between the ion and water-molecule.

We propose the following model function to satisfy these results as follows:

$$\mathrm{\epsilon }\left(r\right)=1+\left({\mathrm{\epsilon }}_0-1\right)\left[1-exp\left\{-\left(r-{\mathrm{r}}_0\right)/\mathrm{\kappa }\right\}\right]$$E124

where ε₀ (=78.35) is the dielectric constant of water. Other parameters are numerically equal to r₀ = 5 A and κ = 3.44 A⁻¹, respectively.

The insertion of this dielectric function ε(r) for the long range part of inter-particle potential is not necessary in molecular dynamic simulations. The dynamics produces automatically the configuration of constituents to satisfy the dielectric behavior.

On the analogy of the inter-ionic potentials in molten salts, ϕ⁺⁻(r) in aqueous solution, where the dipole-dipole and dipole-quadrupole dispersion forces are neglected, may be given as follows:

$${\mathrm{\varphi }}^{+-}\left(r\right)=\left({\mathrm{z}}^{+}{\mathrm{z}}^{-}{\mathrm{e}}^2\right)/\left\{r\mathrm{\epsilon }\left(r\right)\right\}+{\mathrm{A}}^{+-}exp\left[{\mathrm{B}}^{+-}\left({{\mathrm{d}}_{\mathrm{i}}}^{+}+{{\mathrm{d}}_{\mathrm{j}}}^{-}-r\right)\right]$$E125

where A⁺⁻ is a constant in relation to the magnitude of repulsive force between cation i and anion j. B⁺⁻ the softness parameter and (d_i⁺ + d_j⁻) is the hard core contact between cation i and anion j. A⁺⁻ and B⁺⁻ are also given in the literature [27]. The difference between this expression and that of ionic crystal or of molten salt is only ascribed to whether the introduction of the dielectric function ε(r) is necessary or not.

For simplicity, the pair potentials ϕ^{+ w}(r) and ϕ^{− w}(r) are assumed to be a combined form of the repulsive potential and the charge-dipole potential. Then the pair potential between cation and water molecule centered at the oxygen atom, ϕ^{+ w}(r), is written as follows:

$${\mathrm{\varphi }}^{+\mathrm{w}}\left(r\right)={{\mathrm{\varphi }}_{\mathrm{rep}}}^{+\mathrm{w}}\left(r\right)-{\mathrm{z}}^{+}\text{e}\text{μ}cos\text{θ}\left(1–{\text{3l}}^2/8r^2\right)/\left\{r^2\mathrm{\epsilon }\left(r\right)\right\}$$E126

where ϕ_rep^+w(r) is repulsive potential and its explicit form will be given be later. μ is the dipole of water molecule and l its length. θ is the dipole’s angle in the direction of cation.

It is well-known that the above expression is converted to the following form according to Boltzmann law,

$${\mathrm{\varphi }}^{+\mathrm{w}}\left(r\right)={{\mathrm{\varphi }}_{\mathrm{rep}}}^{+\mathrm{w}}\left(r\right)-{\mathrm{z}}^{+2}{\mathrm{e}}^2{\mathrm{\mu }}^2\left(1–{\text{3l}}^2/8r^2\right)/\left\{{\text{3k}}_{\mathrm{B}}\mathrm{T}{r}^2\mathrm{\epsilon }\left(r\right)\right\}$$E127

On the other hand, a modified Lennard-Jones potential for water molecule, ϕ^ww(r), is useful and is written as follows [45]:

$${\mathrm{\varphi }}^{\mathrm{ww}}\left(r\right)=\text{4C}\left[{\left({\mathrm{d}}^{\mathrm{w}}/r\right)}^{12}-{\left({\mathrm{d}}^{\mathrm{w}}/\mathrm{r}\right)}^6\right]-{\text{2μ}}^2/r^3$$E128

In this equation, the term 4C(d^w/r)¹² is equal to the repulsive part and the parameters C and d^w for water molecule in the gaseous state are equal to 230.9 k_B and 2.824 Å, respectively.

The repulsive part of inter-ionic potential for ϕ⁺⁺(r) may be approximately described as the form of A⁺⁺exp[B⁺⁺(d⁺ + d⁺ − r)] similar to the repulsive one in Eq. (125), since its interaction occurs around the distance of close contact where the dielectric behavior of neighboring water molecules must be neglected.

Now let us assume that the repulsive potential ϕ_rep^+w(r) is represented by the root mean square of 4C(d^w/r)¹² and A⁺⁺exp[B⁺⁺ (d⁺ + d⁺- r)] as follows:

$${{\mathrm{\varphi }}_{\mathrm{rep}}}^{+\mathrm{w}}\left(\mathrm{r}\right)={\left\{\text{4C}{\left({\mathrm{d}}^{\mathrm{w}}/r\right)}^{12}{\mathrm{A}}^{++}exp\left[{\mathrm{B}}^{++}\left({\mathrm{d}}^{+}+{\mathrm{d}}^{+}-r\right)\right]\right\}}^{1/2}$$E129

Insertion of (129) into Eq. (127) gives us the following expression,

$${\mathrm{\varphi }}^{+\mathrm{w}}\left(r\right)={\left\{\text{4C}{\left({\mathrm{d}}^{\mathrm{w}}/r\right)}^{12}{\mathrm{A}}^{++}exp\left[{\mathrm{B}}^{++}\left({\mathrm{d}}^{+}+{\mathrm{d}}^{+}-r\right)\right]\right\}}^{1/2}-{\mathrm{z}}^2{\mathrm{e}}^2{\mathrm{\mu }}^2\left(1–{\text{3l}}^2/{\text{8r}}^2\right)/\left\{{\text{3k}}_{\mathrm{B}}\mathrm{T\epsilon }\left(r\right)r^2\right\}$$E130

In a similar way, the inter-particle potential between anion and water molecule is expressed as follows:

$${\mathrm{\varphi }}^{-\mathrm{w}}\left(r\right)={\left\{\text{4C}{\left({\mathrm{d}}^{\mathrm{w}}/r\right)}^{12}{\mathrm{A}}^{--}exp\left[{\mathrm{B}}^{--}\left({\mathrm{d}}^{-}+{\mathrm{d}}^{-}-r\right)\right]\right\}}^{1/2}-{\mathrm{z}}^2{\mathrm{e}}^2{\mathrm{\mu }}^2\left(1–{\text{3l}}^2/8r^2\right)/\left\{{\text{3k}}_{\mathrm{B}}\mathrm{T\epsilon }\left(r^2\right)\right\}$$E131

The dipole moment of water molecule is known to be μ = 0.38 (in the unit of e times 1 Å = 1.6 × 10⁻²⁹ C·m) and l ≒0.5 Å. Therefore, all parameters in (130) and (131) are known. According to Bopp et al. [51], the repulsive parts in (130) and (131) are converted to the exponential decaying functions similar to the repulsive part in (125) [46, 47].

Under these circumstances, it is possible to use either our empirical expressions (130) and (131), or to apply the inter-particle potentials derived by Bopp et al. [51]. It is also possible to estimate the repulsion terms in (130) and (131) by using wave mechanical approach. In fact, the ion-water molecule interactions were obtained by such an elaborating method [33, 52]. However, we will use the above empirical potentials for numerical application, for simplicity.

17. Momentum conservation and the tag of water molecules by ion’s movement

We will investigate the tag of water molecules by ion’s moving in the electrolytic solutions from the view point of equation of motion under an applied field E [28].

Under this situation, the second law of motion for the cation i can be written as follows:

$${\mathrm{m}}^{+}\mathrm{d}{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(t\right)/\mathrm{d}t=\sum _{\mathrm{j}\ne \mathrm{i}}^{\mathrm{N}+\text{and}\mathrm{N}-}{\mathbf{f}}_{\mathrm{ij}}+\sum _{\mathrm{q}\ne \mathrm{i}}^{\mathrm{X}}{\mathbf{F}}_{\mathrm{iq}}+{\mathrm{z}}^{+}\mathrm{e}\mathbf{E}$$E132

f_ij is the force acting on the ion i from the ion j and F_iq is that from the water molecule q.

At the time of steady state, τ, after applying the external field E, we have

$${\mathrm{m}}^{+}{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(\mathrm{\tau }\right)={{\int }_0}^{\mathrm{\tau }}\sum _{\mathrm{j}\ne \mathrm{i}}^{\mathrm{N}+\text{and}\mathrm{N}-}{\mathbf{f}}_{\mathrm{ij}}\mathrm{d}t+{{\int }_0}^{\mathrm{\tau }}\sum _{\mathrm{q}\ne \mathrm{i}}^{\mathrm{X}}{\mathbf{F}}_{\mathrm{iq}}\mathrm{d}t+{{\int }_0}^{\mathrm{\tau }}{\mathrm{z}}^{+}\mathrm{e}\mathbf{E}\mathrm{d}t+{\mathrm{m}}^{+}{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right)$$E133

In a similar way, we have

$${\mathrm{m}}^{-}\underset{\mathrm{k}}{\overset{-}{\mathbf{v}}}\left(\mathrm{\tau }\right)={{\int }_0}^{\mathrm{\tau }}\sum _{\mathrm{l}\ne \mathrm{k}}^{\mathrm{N}-\text{and}\mathrm{N}+}{\mathbf{f}}_{\mathrm{kl}}\mathrm{d}t+{{\int }_0}^{\mathrm{\tau }}\sum _{\mathrm{q}\ne \mathrm{k}}^{\mathrm{X}}{\mathbf{F}}_{\mathrm{kq}}\mathrm{d}t+{{\int }_0}^{\mathrm{\tau }}{\mathrm{z}}^{-}\mathrm{e}\mathbf{E}\mathrm{d}t+{\mathrm{m}}^{-}{{\mathbf{v}}_{\mathrm{k}}}^{-}\left(0\right)$$E134

and

$${\mathrm{m}}^{\mathrm{w}}{{\mathbf{v}}_{\mathrm{q}}}^{\mathrm{w}}\left(\mathrm{\tau }\right)={{\int }_0}^{\mathrm{\tau }}\sum _{\mathrm{q}\ne \left(\mathrm{k}\text{and}\mathrm{i}\right)}^{\mathrm{X}}\left({\mathbf{F}}_{\mathrm{iq}}+{\mathbf{F}}_{\mathrm{kq}}\right)\mathrm{d}t+{\mathrm{m}}^{\mathrm{w}}{{\mathbf{v}}_{\mathrm{q}}}^{\mathrm{w}}\left(0\right)$$E135

In a unit volume, the total summation of the ensemble averages of these momenta is written as follows:

$$\begin{array}{l}\mathrm{n}<{\mathrm{m}}^{+}{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(\mathrm{\tau }\right)>+\mathrm{n}<{\mathrm{m}}^{-}{{\mathbf{v}}_{\mathrm{k}}}^{-}\left(\mathrm{\tau }\right)>+{\mathrm{n}}_{\mathrm{w}}<{\mathrm{m}}^{\mathrm{w}}{{\mathbf{v}}_{\mathrm{q}}}^{\mathrm{w}}\left(\mathrm{\tau }\right)>\\ =1/{\mathrm{V}}_{\mathrm{M}}{{\int }_0}^{\mathrm{\tau }}\sum \left({\mathbf{f}}_{\mathrm{ij}}+{\mathbf{f}}_{\mathrm{ji}}\right)\mathrm{d}t+1/{\mathrm{V}}_{\mathrm{M}}{{\int }_0}^{\mathrm{\tau }}\sum \left({\mathbf{F}}_{\mathrm{iq}}+{\mathbf{F}}_{\mathrm{kq}}+{\mathbf{F}}_{\mathrm{qi}}+{\mathbf{F}}_{\mathrm{qk}}\right)\mathrm{d}t\\ +1/{\mathrm{V}}_{\mathrm{M}}{{\int }_0}^{\mathrm{\tau }}\left({\mathrm{z}}^{+}+{\mathrm{z}}^{-}\right)\mathrm{e}\mathbf{E}\mathrm{d}t+\mathrm{n}<{\mathrm{m}}^{+}{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(0\right)>+\mathrm{n}<{\mathrm{m}}^{-}{{\mathbf{v}}_{\mathrm{k}}}^{-}\left(0\right)>\\ +{\mathrm{n}}_{\mathrm{w}}<{\mathrm{m}}^{\mathrm{w}}{{\mathbf{v}}_{\mathrm{q}}}^{\mathrm{w}}\left(0\right)>\end{array}$$E136

where n_w is the number density of water molecules.

The summation of last three terms on the right hand side of this equation is equal to zero, because there is no external force at t = 0. All other terms on the right hand side of this equation are equal to zero by considering the law of action and reaction and charge neutrality condition.

Therefore, we have

$$\mathrm{n}<{\mathrm{m}}^{+}{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(\mathrm{\tau }\right)>+\mathrm{n}<{\mathrm{m}}^{-}{{\mathbf{v}}_{\mathrm{k}}}^{-}\left(\mathrm{\tau }\right)>+{\mathrm{n}}_{\mathrm{w}}<{\mathrm{m}}^{\mathrm{w}}{{\mathbf{v}}_{\mathrm{q}}}^{\mathrm{w}}\left(\mathrm{\tau }\right)>=0$$E137

This equation indicates that the partial conductivity ratio < v_i⁺(τ)>/< v_k⁻(τ) > is not equal to the inverse mass ratio m⁻/m⁺, which is essentially different from the case of molten salts.

Some of water molecules may be simultaneously pulled by the dissolved ions under an external field E. Here, we neglect the relative time-relaxation for velocities of particles undergoing the co-operative motion. Taking the numbers of pulled water-molecules by each cation and anion, as x⁺ and x⁻, we have

$${\mathrm{nx}}^{+}+{\mathrm{nx}}^{-}+{\mathrm{x}}^{\mathrm{r}}={\mathrm{n}}_{\mathrm{w}}$$E138

Here, x^r is equal to the number density of un-pulled water molecules.

Since, the movements of remainder water molecules under the external field must be isotropic, we have x^r < m^wv^w(τ) > = 0. Then n_w < m^wv_q^w(τ) > is expressed as follows:

$${\mathrm{n}}_{\mathrm{w}}<{\mathrm{m}}^{\mathrm{w}}{{\mathbf{v}}_{\mathrm{q}}}^{\mathrm{w}}\left(\mathrm{\tau }\right)>={\mathrm{x}}^{+}<{\mathrm{m}}^{\mathrm{w}}{{\mathbf{v}}_{\mathrm{i}}}^{+}\left(\mathrm{\tau }\right)>+{\mathrm{x}}^{-}<{\mathrm{m}}^{\mathrm{w}}{{\mathbf{v}}_{\mathrm{k}}}^{-}\left(\mathrm{\tau }\right)>$$E139

Insertion of this equation into (66) gives us the following relation:

$$<\left({\mathrm{m}}^{+}+{\mathrm{x}}^{+}{\mathrm{m}}^{\mathrm{w}}\right){\mathbf{v}}^{+}\left(\mathrm{\tau }\right)>+<\left({\mathrm{m}}^{-}+{\mathrm{x}}^{-}{\mathrm{m}}^{\mathrm{w}}\right){\mathbf{v}}^{-}\left(\mathrm{\tau }\right)>=0$$E140

Hereafter, we omit the suffix of ion i or k.

Therefore, we have

$$\mid <{\mathbf{v}}^{+}\left(\mathrm{\tau }\right)>\mid /\mid <{\mathbf{v}}^{-}\left(\mathrm{\tau }\right)>\mid =\left({\mathrm{m}}^{-}+{\mathrm{x}}^{-}{\mathrm{m}}^{\mathrm{w}}\right)/\left({\mathrm{m}}^{+}+{\mathrm{x}}^{+}{\mathrm{m}}^{\mathrm{w}}\right)$$E141

We cannot apply the above treatment for H⁺ and OH⁻ ions, because their conduction mechanisms differ from that of all other dissolved ions. Their mechanisms are known as the Grotthus-type conduction which is a kind of hopping conduction of electrons or holes [3].

It is, however, straightforward to obtain the following relation for all dissolved ions in their dilute limits except for H⁺ and OH⁻ ones,

$$\mid <{\mathbf{v}}^{+}\left(\mathrm{\tau }\right)>\mid \left({\mathrm{m}}^{+}+{\mathrm{x}}^{+}{\mathrm{m}}^{\mathrm{w}}\right)=\mid <{\mathbf{v}}^{-}\left(\mathrm{\tau }\right)>\mid \left({\mathrm{m}}^{-}+{\mathrm{x}}^{-}{\mathrm{m}}^{\mathrm{w}}\right)$$E142

This relation seems to be valid for all aqueous solutions of equivalent electrolytes in the dilution limit.

Using Eqs. (114) and (115), Eq. (142) for the dilution limit of electrolytic solution is expressed as follows:

$$\begin{array}{ll}{\mathrm{\sigma }}^{+}/{\mathrm{\sigma }}^{-}& =\left(\text{mass of an anion plus masses of water molecules pulled}\mathrm{by}\mathrm{its}\text{anion}\right)/\\ & \left(\text{mass of}\mathrm{a}\text{cation plus masses of water molecules pulled}\mathrm{by}\mathrm{its}\text{cation}\right)\\ & =\left({\mathrm{m}}^{-}+{\mathrm{x}}^{-}{\mathrm{m}}^{\mathrm{w}}\right)/\left({\mathrm{m}}^{+}+{\mathrm{x}}^{+}{\mathrm{m}}^{\mathrm{w}}\right)={\left({\mathrm{m}}^{-}<{\mathrm{\varphi }}^{-\mathrm{w}}>\right)}^{1/2}/{\left({\mathrm{m}}^{+}<{\mathrm{\varphi }}^{+\mathrm{w}}>\right)}^{1/2}\end{array}$$E143

This equation may correspond to the inverse mass ratio for the partial conductivities of molten salt [6].

18. Numerical results in electrolytic solutions

According to the theoretical results we have discussed so far, the pair distribution functions appear in the essential equations [28]. Therefore, how to obtain the pair distribution functions is one of the matters of vital importance.

There are several standard theoretical methods to obtain the pair distribution functions in molecular liquids from the knowledge of inter-particle potentials [33]. In the calculation of site-site distribution function for such a molecular liquid, the reference interaction-site model (RISM) approximation proposed by Chandler and Anderson [52] seems to be useful. Until the present time, the extension of RISM approximation, in order to obtain the potentials of mean force and also the site-site pair distribution functions ɡ^μν(r)‘s in electrolytic solutions, has been carried out by several authors [53, 54, 55]. These attempts cover the insufficient experimental knowledge for pair distribution functions ɡ⁺⁻(r), ɡ^{+ w}(r) and ɡ^{− w} (r).

However, we will use the ɡ^μν(r)‘s in aqueous solution of sodium chloride obtained by our own MD simulation. The essential numerical procedure of MD simulation in this study is same as our previous works of molten salts [56]. The procedure of MD simulation of electrolyte aqueous solution will be briefly described as follows for reader’s benefit. In MD for the electrolyte aqueous solution, the rigid body models (TIP4P) [57] are used for water molecules. The interactions between constituent TIP4P water molecules are expressed as the charged L-J type potentials, as,

$${}_{\mathit{ij}}\left(r\right)=\frac{z_i{z}_j{e}^2}{r}+\frac{A}{r^{12}}-\frac{B}{r^6}$$E144

The interactions between alkali metal cation and halide anion, TIP4P- alkali metal anion, and TIP4P – halide anion are expressed as [58]:

$${}_{\mathit{ij}}\left(r\right)=\frac{z_i{z}_j{e}^2}{r}+\frac{C}{r^9}-\frac{D}{r^6}$$E145

In (144) and (145), i and j stand for the constituent atoms; e is the elementary charge. The used charges for the constituent species z_i and the interaction parameters are taken from the literature; TIP4P – TIP4P [57]; TIP4P – alkali metal cation, TIP4P – halide anon, between alkali metal cation, between halide anion, and between alkali metal cation and halide anion [58]. The Ewald method is used for the calculation of the Coulomb interaction. For the structure calculation, MD is performed in NTP constant condition [59, 60, 61] under the pressure of 1 atm at 283 K. MD is performed for 50,000 steps with 0.1 fs one time step in 1.1% NaCl aqueous solution. MD cell contains about 10,000 molecules (i.e. 30,000 atoms) for the calculation of the structure and the velocity autocorrelation function for alkali halide aqueous solution. The numbers of the constituent ions in the MD cell are listed in Table 2.

Solute	Water (TIP4P)	Cation	Anion
Li^+ Cl^−	10,000	112	112
Na^+ Cl^−	10,000	112	112
K^+ Cl^−	10,000	112	112

Table 2.

The numbers of ions in MD cell.

The main part of MD is performed using SIGRESS ME package (Fujitsu) at the supercomputing facilities in Kyushu University.

The obtained figures of ɡ_ij(r) are shown in Figures 5–8. And using these data, we have estimated the numbers of water molecules involved within a sphere of radius r from the centered ion, n_ij(r) (i = Li⁺, Na⁺, K⁺ and Cl⁻; j = oxygen of water molecule) = 4πʃ₀^r ɡ_ij(r)r² dr, which are also figured in them.

Figure 5.

Pair distribution function of water molecules around a Li + ion, gLi-w(r) in the electrolyte solution of LiCl. And numbers of locating water molecules around a Li + ion within the sphere of the length r centered at its Li + ion, nLi-w(r), in its solution obtained by MD simulation.

Figure 6.

Pair distribution function of water molecules around a Na + ion, gNa-w(r) in the electrolyte solution of NaCl. And numbers of locating water molecules around a Na + ion within the sphere of the length r centered at its Na + ion, nNa-w(r), in its solution obtained by MD simulation.

Figure 7.

Pair distribution function of water molecules around a K+ ion, gK-w(r) in the electrolyte solution of KCl. And numbers of locating water molecules around a K+ ion within the sphere of the length r centered at its K+ ion, nK-w(r), in its solution obtained by MD simulation.

Figure 8.

Pair distribution function of water molecules around a Cl⁻ ion, gCl-w(r) in the electrolyte solution of MCl (M = Li, Na and K). And numbers of locating water molecules around a Cl⁻ ion within the sphere of the length r centered at its Cl⁻ ion, nCl-w(r), in its solution obtained by MD simulation.

Using Eq. (143), that is, σ⁺/σ⁻ = (m⁻ + x⁻m^w)/(m⁺ + x⁺m^w), and taking an assumption that the pulling water molecules for Na⁺ ion is equal to 6.0 although its plausible justification seems to be difficult, then we obtain the pulling water molecules for other ions as shown in Table 3, in which the hydration numbers seen in a text book [62] and our results obtained by MD simulation, for reference.

Ions	Pulling water molecules, x^+ or x^−	Hydration numbers in the text book [36]	Hydration numbers obtained from MD simulations
Li^+	7.6	4.3 ± 0.6	4.1
Na^+	6.0*	5.6 ± 1.7	5.7
K^+	2.8	5.5 ± 1.3	6.4
Cl^−	2.8	6.0 ± 0.7	6.5

Table 3.

Numbers of pulling water molecules, x⁺ or x⁻ and hydration numbers.

^*

Assumption that the pulling number x⁺ of Na⁺ ion is equal to be 6.0 and also that the pulling numbers of water molecules for Cl⁻ are not changed even for that the pairing positive ions are different.

Using these pulling numbers for the constituent ions, we can estimate the term, (m⁻ + x⁻m^w)/(m⁺ + x⁺m^w) as shown in Table 4. As seen in this table, agreements for both terms are satisfactory, which is a kind of proof for the assumption x⁺ is equal to 6.0.

Electrolyte	σ^+/σ^−	(m^− + x^−m^w)/(m^++ x^+m^w)
Li^+ Cl^−	0.595	0.598
Na^+ Cl^−	0.659	0.655
K^+ Cl^−	0.963	0.960

Table 4.

The ratio of ionic conductivity and the calculation results by using Table 3.

It is emphasized that the pulling number of water molecules by moving ion has no relation to the hydration number of water molecules as seen in Table 3. The hydration of water molecules around electrolytic ions is originated essentially by the thermodynamic stability which is related not only to the interaction energies among ions and water molecules but also to the configuration entropy terms. This is because that the pulling number is not always related to the hydration one.

19. Discussion on the electrical conductivities in electrolytic solutions

The present theory seems essentially comparable to the treatments developed by Onsager [19], Fuoss et al. [21], Prigogine [20], Friedman [23], Chandra and Bagchi [27], and Matsunaga and Tamaki [28].

Friedman [23] used a technique of diagram expansion starting from Kubo-Green formula for the conductivity of electrolytic solution and the obtained expression was also written in the form of Λ_c = (Λ⁺ + Λ⁻) = Λ₀ + Λ₁ – kc^1/2. However, his theory is very much sophisticated and too mathematical to understand with a physical insight.

Recent theoretical work carried out by Chandra and Bagchi [27] is basically started from a Kubo-Green type theory, that is, the partial conductivities are derived from velocity correlation functions. Their treatment seems to be a modernized and beautiful and therefore it is very much appreciable. However, the friction force of their theory involves various terms which make it difficult to calculate practically the partial conductivities. In fact, there still remains the task to represent a microscopic formula for Λ₀.

The present treatment is easily to understand in view of physical insight and is successful for deriving the formula of Λ₀.

The short-time expansion forms for < v_i⁺(t) v_i⁺(0)>, < v_k⁻(t) v_k⁻(0) > and < v_i⁺(t) v_k⁻(0) > are expressed in terms inter-particle potentials and corresponding pair distribution functions as seen in (95) and (97). In the case of molten salts, all these velocity correlation functions yield some physical quantities in relation to a part of partial conductivities [6]. In the present case, however, Z_σ⁺(t) and Z_σ⁻(t) play its role. Such an essential difference between the case of molten salt and the electrolytic solution may be ascribed to the difference in the momentum conservation of the system.