Posts for the month of February 2012

Electrical problem including metal diffusion

This post includes the effects of the counter electrode and the metal diffusion. The metal diffusion is included via a solution where the metal gradient doesn't change rapidly at the interface. This is not really very physical, but it will be useful when testing numerical solutions.

At the working electrode $W$:

$ i_T^W = -\kappa \vec{n}^W \cdot \nabla V_{IR}$

$ i_T^W = -c_{DL} \frac{\partial V_{DL}^W}{\partial t} + i_{F}^W $

$ i_F^W = \frac{c_{cu}^0}{ c_{cu}^{\infty}} i_0^W E^W$

$ D_{cu} \vec{n}^W \cdot \nabla c_{cu} = -\frac{i_F^W}{2 F}  $

At the counter electrode $C$:

$ i_T^C = -\kappa \vec{n}^C \cdot \nabla V_{IR}$

$ i_T^C = -c_{DL} \frac{\partial V_{DL}^C}{\partial t} + i_{F}^C $

$ i_F^C = i_0^C E^C$

where

$ E = \exp{\left(-\frac{\alpha F V_{DL}}{R T} \right)} -  \exp{\left(\frac{\left(2 -\alpha\right) F V_{DL}}{R T} \right)} $

The normals point out of the electrolyte.

In the bulk:

$ \nabla^2 V_{IR} = 0 $

$ \frac{\partial c_{cu}}{\partial t} = D_{cu} \nabla^2 c_{cu} $

and

$V_{CELL} = V_{DL}^W + V_{IR}^W - V_{DL}^C - V_{IR}^C$

$V_{APPLIED} = V_{DL}^W + V_{IR}^W$

$\delta_W V_{IR}^C + \delta_C V_{IR}^W = 0 $

$\delta_C + \delta_W = \delta$

where $\delta_C$ and $\delta_W$ are the distances from the counter and working electrodes to the reference electode, respectively. $\delta$ is the total distance between the electrodes.

The metal diffusion has an analytical solution

$ c_{cu} \left( x, t \right) = \frac{c_x^0}{k} \cos{\left(k x\right)} \exp{\left(-k D_{cu} t\right)} + c_x^0 \left(x - \delta \right) + c_{cu}^{\infty} $

where

$ c_x^0 = \frac{i_F^W}{2 F D_{cu}} $

$ k = \frac{\pi}{2 \delta} $

This solution is valid when $|c_{xt}^0| \ll 1$. This leads to two dimensionless quantities that must be small,

$\Gamma_0 = \frac{\delta i_0^W E^W}{c_{cu}^{\infty} 2 F D_{cu}} \ll 1$

and

$ \Gamma_1 = \Gamma_0 \frac{\alpha F V_{APPLIED}}{R T} \ll 1$

We can now derive an expression for $c_{cu}^0$.

$ c_{cu}^0 = \frac{i_F^W}{2 F D} \left[ \frac{ \exp{\left(-k D_{cu} t\right)} }{k} - \delta \right] + c_{cu}^{\infty} $

Using the above solution we can derive a new expression for $i_F^W$.

$ i_F^W = i_0^W \frac{E^W}{ 1 - \frac{1}{2 F D c_{cu}^{\infty}} \left( \frac{ \exp{\left(-k D_{cu} t\right)} }{k} - \delta \right) i_0^W E^W} $

If we do this the equations can be reduced to two coupled ODEs for $V_{DL}^C$ and $V_{DL}^W$:

$ -\frac{\kappa}{\delta_W} \left(V_{APPLIED} - V_{DL}^W\right) = -c_{DL} \frac{\partial V_{DL}^W}{\partial t} +  i_0^W \frac{E^W}{ 1 - \frac{1}{2 F D c_{cu}^{\infty}} \left( \frac{ \exp{\left(-k D_{cu} t\right)} }{k} - \delta \right) i_0^W E^W}$

and

$ \frac{\kappa}{\delta_W} \left(V_{APPLIED} - V_{DL}^W\right) = -c_{DL} \frac{\partial V_{DL}^C}{\partial t} +  i_0^C E^C$

The material parameters are:

$i_0^W = i_0^C = 40\;\text{A / m}^2$

$\alpha = 0.4$

$F = 9.6485 \times 10^{-4}\;\text{J / V / mol}$

$V_{APPLIED} = -0.275\;\text{V}$

$R = 8.314\;\text{J / K / mol}$

$T = 298.0\;\text{K}$

$\kappa = 15.26\;\text{A / V / m}$

$\delta = 100 \times10^{-9}\;\text{m}$

$\delta_W = \delta / 2 \;\text{m}$

$C_{DL} = 0.3\;\text{A s / V / m}^2$

$c_{cu}^{\infty} = 250.0\;\text{mol / m}^3$

$\Gamma_1 = 0.012024 $

$\Gamma_2 = -0.0515 $

$D_{cu} = 5 \times 10^{-10} \;\text{m}^2\text{ / s}$

Integrating this gives:

$V_{DL}^W = -0.27499$

$V_{DL}^C = 0.0686$

$i_T^C = 2865.65$

$V_{IR}^W =  -9.38943 \times 10^{-6}$

$V_{IR}^C =  9.389434 \times 10^{-6}$

$V_{CELL} = -0.3436$

The code for this is source:trunk/moffat/electrical/misc/1D_metal.py.

The evolution of the $V_{DL}$:

source:trunk/moffat/electrical/misc/voltageVersusTimeCopper.png

and the cupric concentrations at various times:

source:trunk/moffat/electrical/misc/copper.png

  • Posted: 2012-02-21 16:59 (Updated: 2012-03-22 13:26)
  • Author: wd15
  • Categories: (none)
  • Comments (0)

Symmetric Electrical Problem

Extending blog:MoffatElectrical1 by making the system totally symmetric

The changes are:

14c14
< alpha = 0.5
---
> alpha = 0.4
30c30
<     return i0 * (numpy.exp(-alpha * Fbar * V) - numpy.exp((1 - alpha) * Fbar * V))
---
>     return i0 * (numpy.exp(-alpha * Fbar * V) - numpy.exp((2 - alpha) * Fbar * V))
34c34
<                       - (1 - alpha) * Fbar * numpy.exp((1 - alpha) * Fbar * V))
---
>                       - (2 - alpha) * Fbar * numpy.exp((2 - alpha) * Fbar * V))

gives:

$V_{DL}^W = -0.255$

$V_{DL}^C = 0.255$

$i_T^C = 5833$

$V_{IR}^W = -0.0191 $

$V_{IR}^C =  0.0191$

$V_{CELL} = -0.55$

The code for this is source:trunk/moffat/electrical/misc/1D_symmetric.py.

The evolution of the $V_{DL}$:

source:trunk/moffat/electrical/misc/voltageVersusTimeSymmetry.png

  • Posted: 2012-02-14 18:22 (Updated: 2012-03-22 13:28)
  • Author: wd15
  • Categories: (none)
  • Comments (0)

Electrical Problem including counter electrode

Including the effects of the counter electrode.

At the working electrode $W$:

$ i_T^W = -\kappa \vec{n}^W \cdot \nabla V_{IR}$

$ i_T^W = -c_{DL} \frac{\partial V_{DL}^W}{\partial t} + i_{F}^W $

$ i_F^W = i_0^W \left[ \exp{\left(-\frac{\alpha F V_{DL}^W}{R T} \right)} -  \exp{\left(\frac{\left(2 -\alpha\right) F V_{DL}^W}{R T} \right)}  \right]$

At the counter electrode $C$:

$ i_T^C = -\kappa \vec{n}^C \cdot \nabla V_{IR}$

$ i_T^C = -c_{DL} \frac{\partial V_{DL}^C}{\partial t} + i_{F}^C $

$ i_F^C = i_0^C \left[ \exp{\left(-\frac{\alpha F V_{DL}^C}{R T} \right)} -  \exp{\left(\frac{\left(2 -\alpha\right) F V_{DL}^C}{R T} \right)}  \right]$

The normals point out of the electrolyte.

In the bulk:

$ \nabla^2 V_{IR} = 0 $

and

$V_{CELL} = V_{DL}^W + V_{IR}^W - V_{DL}^C - V_{IR}^C$

$V_{APPLIED} = V_{DL}^W + V_{IR}^W$

$\delta_W V_{IR}^C + \delta_C V_{IR}^W = 0 $

$\delta_C + \delta_W = \delta$

where $\delta_C$ and $\delta_W$ are the distances from the counter and working electrodes to the reference electode, respectively. $\delta$ is the total distance between the electrodes.

If we do this the equations can be reduced to two coupled ODEs for $V_{DL}^C$ and $V_{DL}^W$:

$ -\frac{\kappa}{\delta_W} \left(V_{APPLIED} - V_{DL}^W\right) = -c_{DL} \frac{\partial V_{DL}^W}{\partial t} +  i_0^W \left[ \exp{\left(-\frac{\alpha F V_{DL}^W}{R T} \right)} -  \exp{\left(\frac{\left(2 -\alpha\right) F V_{DL}^W}{R T} \right)}  \right]$

and

$ \frac{\kappa}{\delta_W} \left(V_{APPLIED} - V_{DL}^W\right) = -c_{DL} \frac{\partial V_{DL}^C}{\partial t} +  i_0^C \left[ \exp{\left(-\frac{\alpha F V_{DL}^C}{R T} \right)} -  \exp{\left(\frac{\left(2 -\alpha\right) F V_{DL}^C}{R T} \right)}  \right]$

The material parameters are:

$i_0^W = i_0^C = 40\;\text{A / m}^2$

$\alpha = 0.4$

$F = 9.6485 \times 10^{-4}\;\text{J / V / mol}$

$V_{APPLIED} = -0.275\;\text{V}$

$R = 8.314\;\text{J / K / mol}$

$T = 298.0\;\text{K}$

$\kappa = 15.26\;\text{A / V / m}$

$\delta = 100 \times10^{-6}\;\text{m}$

$\delta_W = 50 \times10^{-6}\;\text{m}$

$C_{DL} = 0.3\;\text{A s / V / m}^2$

Integrating this gives:

$V_{DL}^W = -0.2666$

$V_{DL}^C = 0.0667$

$i_T^C = 2546$

$V_{IR}^W =  -0.008345$

$V_{IR}^C =  0.008345$

$V_{CELL} = -0.350$

The code for this is source:trunk/moffat/electrical/misc/1D_counter.py.

The evolution of the $V_{DL}$:

source:trunk/moffat/electrical/misc/voltageVersusTimeCounter.png

  • Posted: 2012-02-14 11:43 (Updated: 2012-03-22 13:29)
  • Author: wd15
  • Categories: (none)
  • Comments (0)

A few changes to Bill's code

Bill wants me to make four changes to his code, see source:trunk/kirkendall/A-V_phase-field.py@1457.

He needs to make four changes.

  • Update $V_m$ to be an expression involving $\sigma_{xx}$. This is done in r1458.
  • Save the data from the simulations for later use and read the saved files and plot the results independently from the simulations. See diff:@1458:1460. This now writes a data file directory of the form data/A-V_phase-field_17_36_53Thu02Feb2012. The data files are stored as numpy zip files .npz and numbered by their time step value. plotdata.py reads and plots from the latest data directory created.

The name of an alternative path can be specified as an argument

$ python plotdata.py data/A-V_phase-field_17_33_42Thu02Feb2012

This will use the specified directory, while

$ python plotdata.py

will use the latest created data directory. plotdata.py plots the data for the four main fields using the stored data.

  • Add an improved method for extracting the 0.5 level set. See r1461. This is done using scipy's interp1d. After plotting the four main fields for every stored timestep, plotdata.py will throw up a plot of the interface position.
  • Posted: 2012-02-02 11:45 (Updated: 2012-02-02 18:11)
  • Author: wd15
  • Categories: (none)
  • Comments (0)