Source code for bmlite.P2D.domains

"""
Contains classes to construct the battery for P2D simulations. Each class reads
in keyword arguments that define parameters relevant to its specific domain.
For example, the area and temperature are `Battery` level parameters because
they are the same everywhere, but the discretizations `Nx` and `Nr` may be
different for the anode, separator, and cathode domains.

"""

import numpy as np


[docs] class Battery: def __init__(self, **kwargs) -> None: """ A class for battery-level attributes. Parameters ---------- **kwargs : dict, required Keyword arguments to set the battery-level attributes. The required keys and descriptions are given below: ====== ================================================= Key Description [units] (type) ====== ================================================= cap nominal battery capacity [A*h] (*float*) temp temperature [K] (*float*) area area normal to current collectors [m2] (*float*) ====== ================================================= """ self.cap = kwargs.get('cap') self.temp = kwargs.get('temp') self.area = kwargs.get('area') self.update()
[docs] def update(self): """ Updates any secondary/dependent parameters. At present, this method does not do anything for the `Battery` class. """ pass
[docs] class Electrolyte: def __init__(self, **kwargs) -> None: """ A class for the electrolyte attributes and methods. Parameters ---------- **kwargs : dict, required Keyword arguments to set the electrolyte attributes. The required keys and descriptions are given below: ========== ================================================ Key Description [units] (type) ========== ================================================ li_0 initial Li+ concentration [kmol/m3] (*float*) D_deg `D` degradation factor [-] (*float*) t0_deg `t0` degradation factor [-] (*float*) kappa_deg `kappa` degradation factor [-] (*float*) gamma_deg `gamma` degradation factor [-] (*float*) material class name from `bmlite.materials` [-] (*str*) ========== ================================================ """ self.Li_0 = kwargs.get('Li_0') self.D_deg = kwargs.get('D_deg') self.t0_deg = kwargs.get('t0_deg') self.kappa_deg = kwargs.get('kappa_deg') self.gamma_deg = kwargs.get('gamma_deg') self.material = kwargs.get('material')
[docs] def update(self) -> None: """ Updates any secondary/dependent parameters. For the `Electrolyte` class, this only initializes the material class. """ from .. import materials ElyteMaterial = getattr(materials, self.material) self._material = ElyteMaterial()
[docs] def get_D(self, C_Li: float | np.ndarray, T: float) -> float | np.ndarray: """ Calculate the lithium ion diffusivity in the electrolyte solution at concentration `C_Li` and temperature `T`. Parameters ---------- C_Li : float | 1D array Lithium ion concentration in the electrolyte [kmol/m3]. T : float Battery temperature [K]. Returns ------- D : float | 1D array Lithium ion diffusivity in the electrolyte [m2/s]. """ return self.D_deg * self._material.get_D(C_Li, T)
[docs] def get_t0(self, C_Li: float | np.ndarray, T: float) -> float | np.ndarray: """ Calculate the lithium ion transference number at concentration `C_Li` and temperature `T`. Parameters ---------- C_Li : float | 1D array Lithium ion concentration in the electrolyte [kmol/m3]. T : float Battery temperature [K]. Returns ------- t0 : float | 1D array Lithium ion transference number [-]. """ return self.t0_deg * self._material.get_t0(C_Li, T)
[docs] def get_kappa(self, C_Li: float | np.ndarray, T: float) -> float | np.ndarray: """ Calculate the electrolyte conductivity at concentration `C_Li` and temperature `T`. Parameters ---------- C_Li : float | 1D array Lithium ion concentration in the electrolyte [kmol/m3]. T : float Battery temperature [K]. Returns ------- kappa : float | 1D array Electrolyte conductivity [S/m]. """ return self.kappa_deg * self._material.get_kappa(C_Li, T)
[docs] def get_gamma(self, C_Li: float | np.ndarray, T: float) -> float | np.ndarray: """ Calculate the electrolyte thermodynamic factor at concentration `C_Li` and temperature `T`. Parameters ---------- C_Li : float | 1D array Lithium ion concentration in the electrolyte [kmol/m3]. T : float Battery temperature [K]. Returns ------- gamma : float | 1D array Thermodynamic factor [-]. """ return self.gamma_deg * self._material.get_gamma(C_Li, T)
[docs] class Electrode: def __init__(self, name: str, **kwargs) -> None: """ A class for the electrode-specific attributes and methods. This class is used to build both the anode and cathode for P2D simulations. Parameters ---------- name : str Name of the electrode, must be either 'anode' or 'cathode'. **kwargs : dict, required Keyword arguments to set the electrode attributes. The required keys and descriptions are given below: ========== ======================================================= Key Description [units] (type) ========== ======================================================= Nx number of `x` discretizations [-] (*int*) Nr number of `r` discretizations [-] (*int*) thick electrode thickness [m] (*float*) R_s represenatative particle radius [m] (*float*) eps_el electrolyte volume fraction [-] (*float*) eps_CBD carbon binder domain volume fraction [-] (*float*) eps_void void volume fraction [-] (*float*) p_sol solid Bruggeman factor, `eps_s**p_sol` [-] (*float*) p_liq liquid Bruggeman factor, `eps_el**p_liq` [-] (*float*) alpha_a Butler-Volmer anodic symmetry factor [-] (*float*) alpha_c Butler-Volmer cathodic symmetry factor [-] (*float*) Li_max max solid-phase Li concentration [kmol/m3] (*float*) x_0 initial solid Li intercalation fraction [-] (*float*) i0_deg `i0` degradation factor [-] (*float*) Ds_deg `Ds` degradation factor [-] (*float*) material class name from `bmlite.materials` [-] (*str*) csvfile path to CSV file with OCV data (optional) (*str*) submodels `submodels` classes to include (*dict[dict]*) ========== ======================================================= """ from . import submodels if name not in ['anode', 'cathode']: raise ValueError("'name' must be either 'anode' or 'cathode'.") self._name = name self.Nx = kwargs.get('Nx') self.Nr = kwargs.get('Nr') self.thick = kwargs.get('thick') self.R_s = kwargs.get('R_s') self.eps_s = kwargs.get('eps_s') self.eps_el = kwargs.get('eps_el') self.eps_CBD = kwargs.get('eps_CBD') self.p_sol = kwargs.get('p_sol') self.p_liq = kwargs.get('p_liq') self.alpha_a = kwargs.get('alpha_a') self.alpha_c = kwargs.get('alpha_c') self.Li_max = kwargs.get('Li_max') self.x_0 = kwargs.get('x_0') self.i0_deg = kwargs.get('i0_deg') self.Ds_deg = kwargs.get('Ds_deg') self.material = kwargs.get('material') self.csvfile = kwargs.get('csvfile', None) self.update() self._submodels = {} all_submodels = kwargs.get('submodels', {}) if 'Hysteresis' in all_submodels: Hysteresis, opt = submodels.Hysteresis, all_submodels['Hysteresis'] self._submodels['Hysteresis'] = Hysteresis(self, **opt)
[docs] def update(self) -> None: """ Updates any secondary/dependent parameters. For the `Electrode` class, this initializes the material class, and sets the following: * Void-phase volume fraction [-]: `eps_void = 1 - eps_s - eps_el` * Activate material volume fraction [-]: `eps_AM = eps_s - eps_CBD` * Solid-phase conductivity [S/m]: `sigma_s = 10 * eps_s` * Specific particle surface area [m2/m3]: `A_s = 3 * eps_AM / R_s` """ import inspect from .. import materials self.eps_void = 1. - self.eps_s - self.eps_el self.eps_AM = self.eps_s - self.eps_CBD self.sigma_s = 10. * self.eps_s self.A_s = 3. * self.eps_AM / self.R_s if self.eps_void < -np.finfo(float).eps: raise ValueError('eps_s + eps_el > 1.0') self.eps_void = max(self.eps_void, 0.) if self.eps_AM <= 0.: raise ValueError('eps_s <= eps_CBD.') Material = getattr(materials, self.material) if 'csvfile' in inspect.signature(Material).parameters: self._material = Material(self.alpha_a, self.alpha_c, self.Li_max, csvfile=self.csvfile) else: self._material = Material(self.alpha_a, self.alpha_c, self.Li_max)
[docs] def get_Ds(self, x: float | np.ndarray, T: float, fluxdir: float | np.ndarray) -> float | np.ndarray: """ Calculate the lithium diffusivity in the solid phase given the local intercalation fraction `x` and temperature `T`. Parameters ---------- x : float | 1D array Lithium intercalation fraction [-]. T : float Battery temperature [K]. fluxdir : float | 1D array Lithiation direction: +1 for lithiation, -1 for delithiation, 0 for rest. Used for direction-dependent parameters. Ensure the zero case is handled explicitly or via a default (lithiating or delithiating). Returns ------- Ds : float | 1D array Lithium diffusivity in the solid phase [m2/s]. """ return self.Ds_deg * self._material.get_Ds(x, T, fluxdir)
[docs] def get_i0(self, x: float | np.ndarray, C_Li: float | np.ndarray, T: float, fluxdir: float | np.ndarray) -> float | np.ndarray: """ Calculate the exchange current density given the surface intercalation fraction `x` at the particle surface, the local lithium ion concentration `C_Li`, and temperature `T`. The input types for `x` and `C_Li` should both be the same (i.e., both float or both 1D arrays). Parameters ---------- x : float | 1D array Lithium intercalation fraction at `r = R_s` [-]. C_Li : float | 1D array Lithium ion concentration in the local electrolyte [kmol/m3]. T : float Battery temperature [K]. fluxdir : float | 1D array Lithiation direction: +1 for lithiation, -1 for delithiation, 0 for rest. Used for direction-dependent parameters. Ensure the zero case is handled explicitly or via a default (lithiating or delithiating). Returns ------- i0 : float | 1D array Exchange current density [A/m2]. """ return self.i0_deg * self._material.get_i0(x, C_Li, T, fluxdir)
[docs] def get_Eeq(self, x: float | np.ndarray) -> float | np.ndarray: """ Calculate the equilibrium potential given the surface intercalation fraction `x` at the particle surface. Parameters ---------- x : float | 1D array Lithium intercalation fraction at `r = R_s` [-]. Returns ------- Eeq : float | 1D array Equilibrium potential [V]. """ return self._material.get_Eeq(x)
[docs] def get_Mhyst(self, x: float | np.ndarray) -> float | np.ndarray: """ Calculate the hysteresis magnitude given the surface intercalation fraction `x` at the particle surface. Parameters ---------- x : float | 1D array Lithium intercalation fraction at `r = R_s` [-]. Returns ------- M_hyst : float | 1D array Hysteresis magnitude [V]. """ return self._material.get_Mhyst(x)
[docs] def make_mesh(self, xshift: float = 0., pshift: int = 0) -> None: """ Determines/sets the `x` locations for all of the "minus" interfaces `xm`, "plus" interfaces `xp`, and control volume centers `x` based on the electrode thickness and `Nx` discretization. At present, only a uniform mesh is supported. Also determines/sets the `r` locations for all of the "minus" interfaces `rm`, "plus" interfaces `rp`, and control volume centers `r` based on the representative particle radius and `Nr` discretization. At present, only a uniform mesh is supported. Parameters ---------- xshift : float, optional pshift : int, optional See Also -------- batmods.mesh.x_ptr, batmods.mesh.xr_ptr, batmods.mesh.uniform_mesh """ from ..mesh import x_ptr, xr_ptr, uniform_mesh # Mesh locations self.xm, self.xp, self.x = uniform_mesh(self.thick, self.Nx, xshift) self.rm, self.rp, self.r = uniform_mesh(self.R_s, self.Nr) # Pointers # [[ptr_an], [ptr_sep], [ptr_ca]] # ptr_an and ptr_ca -> [[Li_ed(0->R_s)], phi_ed, Li_el, phi_el, ...] self.ptr = {} self.ptr['xs'] = 0 + pshift self.ptr['r_off'] = 1 self.ptr['phis'] = self.ptr['xs'] + self.Nr self.ptr['ce'] = self.ptr['phis'] + 1 self.ptr['phie'] = self.ptr['ce'] + 1 xvars = ['phis', 'ce', 'phie'] last_xvar = 'phie' # Submodels only support new x variables (like hysteresis... not xr) for model in self._submodels.values(): new_xvar = model.make_mesh(last_xvar, pshift) xvars.append(new_xvar) last_xvar = new_xvar submodel_count = len(self._submodels) self.ptr['x_off'] = self.Nr + 3 + submodel_count self.ptr['start'] = pshift self.ptr['size'] = self.ptr['x_off']*self.Nx self.ptr['shift'] = self.ptr['size'] x_ptr(self, xvars) xr_ptr(self, ['xs'])
def sv0(self, el: object) -> np.ndarray: start = self.ptr['start'] size = self.ptr['size'] sv0 = np.zeros(size) sv0[self.xr_ptr['xs'].flatten() - start] = self.x_0 sv0[self.x_ptr['phis'] - start] = self.phi_0 sv0[self.x_ptr['ce'] - start] = el.Li_0 sv0[self.x_ptr['phie'] - start] = el.phi_0 for model in self._submodels.values(): model.sv0(sv0) return sv0 def algidx(self) -> np.ndarray: algidx = np.hstack([self.x_ptr['phis'], self.x_ptr['phie']]) for model in self._submodels.values(): model.algidx(algidx) return np.sort(algidx) def to_dict(self, soln: object) -> dict: xs = np.zeros([soln.t.size, self.Nx, self.Nr]) for i in range(soln.t.size): xs_tmp = soln.y[i, self.xr_ptr['xs'].flatten()] xs[i, :, :] = xs_tmp.reshape([self.Nx, self.Nr]) ed_soln = { 'x': self.x, 'r': self.r, 'xs': xs, 'cs': xs*self.Li_max, 'phis': soln.y[:, self.x_ptr['phis']], 'ce': soln.y[:, self.x_ptr['ce']], 'phie': soln.y[:, self.x_ptr['phie']], } for model in self._submodels.values(): outputs = model.to_dict(soln) ed_soln.update(outputs) return ed_soln def _boundary_voltage(self, soln) -> np.ndarray: """ Calculate and return the boundary voltage at all solution times. Parameters ---------- soln : Solution A P2D solution instance, step or cycle. Returns ------- voltage_V : np.ndarray Boundary voltage, in volts. """ if self._name == 'anode': idx = 0 elif self._name == 'cathode': idx = -1 else: raise ValueError("Electrode name not in {'anode', 'cathode'}.") V_ptr = self.x_ptr['phis'][idx] return soln.y[:, V_ptr] def _boundary_current(self, soln) -> np.ndarray: """ Calculate and return the boundary current at all solution times. Parameters ---------- soln : Solution A P2D solution instance, step or cycle. Returns ------- current_A : np.ndarray Boundary current, in amps. """ sim = soln._sim bat, an, ca = sim.bat, sim.an, sim.ca c, T = sim.c, bat.temp # calculate the boundary current using solid-phase cons. of charge ed_soln = self.to_dict(soln) phis = ed_soln['phis'] xs_R = ed_soln['xs'][:, :, -1] phie = ed_soln['phie'] ce = ed_soln['ce'] if 'Hysteresis' in self._submodels: hyst = ed_soln['hyst'] Hyst = self.get_Mhyst(xs_R)*hyst else: Hyst = 0. eta = phis - phie - (self.get_Eeq(xs_R) + Hyst) fluxdir = -np.sign(eta) i0 = self.get_i0(xs_R, ce, T, fluxdir) sdot = i0 / c.F * ( np.exp( self.alpha_a*c.F*eta / c.R / T) - np.exp(-self.alpha_c*c.F*eta / c.R / T) ) if self._name == 'anode': i_ext = sdot[:, 0]*an.A_s*c.F*(an.xp[0] - an.xm[0]) \ - an.sigma_s*an.eps_s**an.p_sol \ * (phis[:, 1] - phis[:, 0]) / (an.x[1] - an.x[0]) elif self._name == 'cathode': i_ext = -sdot[:, -1]*ca.A_s*c.F*(ca.xp[-1] - ca.xm[-1]) \ - ca.sigma_s*ca.eps_s**ca.p_sol \ * (phis[:, -1] - phis[:, -2]) / (ca.x[-1] - ca.x[-2]) current_A = i_ext*bat.area return current_A
[docs] class Separator: def __init__(self, **kwargs) -> None: """ A class for the separator attributes and methods. This class is used to build both the separator for P2D simulations. Parameters ---------- **kwargs : dict, required Keyword arguments to set the separator attributes. The required keys and descriptions are given below: ======== ======================================================== Key Description [units] (type) ======== ======================================================== Nx number of `x` discretizations [-] (*int*) thick electrode thickness [m] (*float*) eps_el electrolyte volume fraction [-] (*float*) p_liq liquid Bruggeman factor, `eps_el**p_liq` [-] (*float*) ======== ======================================================== """ self.Nx = kwargs.get('Nx') self.thick = kwargs.get('thick') self.eps_el = kwargs.get('eps_el') self.p_liq = kwargs.get('p_liq') self.update()
[docs] def update(self) -> None: """ Updates any secondary/dependent parameters. For the `Separator` class, this sets the following: * Solid-phase volume fraction [-]: `eps_s = 1 - eps_el` """ self.eps_s = 1 - self.eps_el
[docs] def make_mesh(self, xshift: float = 0., pshift: int = 0) -> None: """ Determines/sets the `x` locations for all of the "minus" interfaces `xm`, "plus" interfaces `xp`, and control volume centers `x` based on the electrode thickness and `Nx` discretization. At present, only a uniform mesh is supported. Parameters ---------- xshift : float, optional pshift : int, optional See Also -------- batmods.mesh.x_ptr, batmods.mesh.xr_ptr, batmods.mesh.uniform_mesh """ from ..mesh import x_ptr, uniform_mesh # Mesh locations self.xm, self.xp, self.x = uniform_mesh(self.thick, self.Nx, xshift) # Pointers # [[ptr_an], Li_el, phi_el, ..., [ptr_ca]] self.ptr = {} self.ptr['ce'] = 0 + pshift self.ptr['phie'] = self.ptr['ce'] + 1 self.ptr['x_off'] = 2 self.ptr['shift'] = self.Nx * self.ptr['x_off'] x_ptr(self, ['ce', 'phie'])
def sv0(self, el: object) -> np.ndarray: import numpy as np return np.tile([el.Li_0, el.phi_0], self.Nx) def algidx(self) -> np.ndarray: return self.x_ptr['phie'] def to_dict(self, soln: object) -> dict: sep_soln = { 'x': self.x, 'ce': soln.y[:, self.x_ptr['ce']], 'phie': soln.y[:, self.x_ptr['phie']], } return sep_soln