8476080312

Committed 29 Mar 2024 01:46AM UTC coverage: 55.461% (-7.2%) from 62.648%

Build # 8476080312

Build Type

Pull #126

github

Committed by

d-cogswell

Commit Message

Adds a benchmark section to docs.

Pull Request Pull Request #126: v1.0.0

Run Details

2351 of 4239 relevant lines covered (55.46%)

2.22 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

66.89

/mpet/props_am.py

"""This module handles properties associated with the active materials.
Only helper functions are defined here.
Diffusion functions are defined in mpet.electrode.diffusion
Chemical potential functions are defined in mpet.electrode.materials"""
import types
import numpy as np

import mpet.geometry as geo
from mpet.config import constants
from mpet.utils import import_function


class muRfuncs():
    """ This class defines functions which describe the chemical
    potential of active materials.
    Each function describes a particular material. The
    chemical potential is directly related to the open circuit voltage
    (OCV) for solid solution materials.
    In each function, muR_ref is the offset (non-dimensional) chemical
    potential by which the returned value will be shifted (for
    numerical convenience).
    Each material function returns:
        muR -- chemical potential
        actR -- activity (if applicable, else None)
    """
    def __init__(self, config, trode, ind=None):
        """config is the full dictionary of
        parameters for the electrode particles, as made for the
        simulations. trode is the selected electrode.
        ind is optinally the selected particle, provided as (vInd, pInd)
        """
        self.config = config
        self.trode = trode
        self.ind = ind
        self.T = config['T']  # nondimensional
        # eokT and kToe are the reference values for scalings
        self.eokT = constants.e / (constants.k * constants.T_ref)
        self.kToe = 1. / self.eokT

        # If the user provided a filename with muRfuncs, try to load
        # the function from there, otherwise load it from this class
        filename = self.get_trode_param("muRfunc_filename")
        muRfunc_name = self.get_trode_param("muRfunc")
        if filename is None:
            # the function will be loaded from the materials folder
            muRfunc = import_function(None, muRfunc_name,
                                      mpet_module=f"mpet.electrode.materials.{muRfunc_name}")
        else:
            muRfunc = import_function(filename, muRfunc_name)

        # We have to make sure the function knows what 'self' is with
        # the types.MethodType function
        self.muRfunc = types.MethodType(muRfunc, self)

    def get_trode_param(self, item):
        """
        Shorthand to retrieve electrode-specific value
        """
        value = self.config[self.trode, item]
        # check if it is a particle-specific parameter
        if self.ind is not None and item in self.config.params_per_particle:
            value = value[self.ind]
        return value

    def get_muR_from_OCV(self, OCV, muR_ref):
        return -self.eokT*OCV + muR_ref

    ######
    # Helper functions
    ######

    def ideal_sln(self, y):
        """ Helper function: Should not be called directly from
        simulation. Call a specific material instead. """
        T = self.T
        muR = T*np.log(y/(1-y))
        return muR

    def reg_sln(self, y, Omga):
        """ Helper function """
        muR_IS = self.ideal_sln(y)
        enthalpyTerm = Omga*(1-2*y)
        muR = muR_IS + enthalpyTerm
        return muR

    def graphite_2param_homog(self, y, Omga, Omgb, Omgc, EvdW):
        """ Helper function """
        y1, y2 = y
        muR1 = self.reg_sln(y1, Omga)
        muR2 = self.reg_sln(y2, Omga)
        muR1 += Omgb*y2 + Omgc*y2*(1-y2)*(1-2*y1)
        muR2 += Omgb*y1 + Omgc*y1*(1-y1)*(1-2*y2)
        muR1 += EvdW * (30 * y1**2 * (1-y1)**2)
        muR2 += EvdW * (30 * y2**2 * (1-y2)**2)
        return (muR1, muR2)

    def graphite_1param_homog(self, y, Omga, Omgb):
        """ Helper function """
        width = 5e-2
        tailScl = 5e-2
        muLtail = -tailScl*1./(y**(0.85))
        muRtail = tailScl*1./((1-y)**(0.85))
        slpScl = 0.45
        muLlin = slpScl*Omga*4*(0.26-y)*step_down(y, 0.5, width)
        muRlin = (slpScl*Omga*4*(0.74-y) + Omgb)*step_up(y, 0.5, width)
        muR = muLtail + muRtail + muLlin + muRlin
        return muR

    def graphite_1param_homog_2(self, y, Omga, Omgb):
        """ Helper function """
        width = 5e-2
        tailScl = 5e-2
        slpScl = 0.45
        muLtail = -tailScl*1./(y**(0.85))
        muRtail = tailScl*1./((1-y)**(0.85))
        muLlin = (slpScl*Omga*12*(0.40-y)
                  * step_down(y, 0.49, 0.9*width)*step_up(y, 0.35, width))
        muRlin = (slpScl*Omga*4*(0.74-y) + Omgb)*step_up(y, 0.5, width)
        muLMod = (0.
                  + 40*(-np.exp(-y/0.015))
                  + 0.75*(np.tanh((y-0.17)/0.02) - 1)
                  + 1.0*(np.tanh((y-0.22)/0.040) - 1)
                  )*step_down(y, 0.35, width)
        muR = muLMod + muLtail + muRtail + muLlin + muRlin
        return muR

    def graphite_1param_homog_3(self, y, Omga, Omgb):
        """ Helper function with low hysteresis and soft tail """
        width = 5e-2
        tailScl = 5e-2
        muLtail = -tailScl*1./(y**(0.85))
        muRtail = tailScl*1./((1-y)**(0.85))
        muRtail = 1.0e1*step_up(y, 1.0, 0.045)
        muLlin = (0.15*Omga*12*(0.40-y**0.98)
                  * step_down(y, 0.49, 0.9*width)*step_up(y, 0.35, width))
        muRlin = (0.1*Omga*4*(0.74-y) + 0.90*Omgb)*step_up(y, 0.5, 0.4*width)
        muLMod = (0.
                  + 40*(-np.exp(-y/0.015))
                  + 0.75*(np.tanh((y-0.17)/0.02) - 1)
                  + 1.0*(np.tanh((y-0.22)/0.040) - 1)
                  )*step_down(y, 0.35, width)
        muR = 0.18 + muLMod + muLtail + muRtail + muLlin + muRlin
        return muR

    def non_homog_rect_fixed_csurf(self, y, ybar, B, kappa, ywet):
        """ Helper function """
        N = len(y)
        ytmp = np.empty(N+2, dtype=object)
        ytmp[1:-1] = y
        ytmp[0] = ywet
        ytmp[-1] = ywet
        dxs = 1./N
        curv = np.diff(ytmp, 2)/(dxs**2)
        muR_nh = -kappa*curv + B*(y - ybar)
        return muR_nh

    def non_homog_rect_variational(self, y, ybar, B, kappa):
        """ Helper function """
        # the taylor expansion at the edges is used
        N_2 = len(y)
        ytmp = np.empty(N_2+2, dtype=object)
        dxs = 1./N_2
        ytmp[1:-1] = y
        ytmp[0] = y[0] + np.diff(y)[0]*dxs + 0.5*np.diff(y,2)[0]*dxs**2
        ytmp[-1] = y[-1] + np.diff(y)[-1]*dxs + 0.5*np.diff(y,2)[-1]*dxs**2
        curv = np.diff(ytmp, 2)/(dxs**2)
        muR_nh = -kappa*curv + B*(y - ybar)
        return muR_nh

    def non_homog_round_wetting(self, y, ybar, B, kappa, beta_s, shape, r_vec):
        """ Helper function """
        dr = r_vec[1] - r_vec[0]
        Rs = 1.
        curv = geo.calc_curv(y, dr, r_vec, Rs, beta_s, shape)
        muR_nh = B*(y - ybar) - kappa*curv
        return muR_nh

    def general_non_homog(self, y, ybar):
        """ Helper function """
        ptype = self.get_trode_param("type")
        mod1var, mod2var = False, False
        if isinstance(y, np.ndarray):
            mod1var = True
            N = len(y)
        elif (isinstance(y, tuple) and len(y) == 2
                and isinstance(y[0], np.ndarray)):
            mod2var = True
            N = len(y[0])
        else:
            raise Exception("Unknown input type")
        if ("homog" not in ptype) and (N > 1):
            shape = self.get_trode_param("shape")
            if shape == "C3":
                if mod1var:
                    kappa = self.get_trode_param("kappa")
                    B = self.get_trode_param("B")
                    if self.get_trode_param("type") in ["ACR"]:
                        cwet = self.get_trode_param("cwet")
                        muR_nh = self.non_homog_rect_fixed_csurf(
                            y, ybar, B, kappa, cwet)
                    elif self.get_trode_param("type") in ["ACR_Diff"]:
                        muR_nh = self.non_homog_rect_variational(
                            y, ybar, B, kappa)
                elif mod2var:
                    kappa = self.get_trode_param("kappa")
                    B = self.get_trode_param("B")
                    cwet = self.get_trode_param("cwet")
                    muR_nh = self.non_homog_rect_fixed_csurf(
                        y, ybar, B, kappa, cwet)
            elif shape in ["cylinder", "sphere"]:
                kappa = self.get_trode_param("kappa")
                B = self.get_trode_param("B")
                beta_s = self.get_trode_param("beta_s")
                r_vec = geo.get_unit_solid_discr(shape, N)[0]
                if mod1var:
                    muR_nh = self.non_homog_round_wetting(
                        y, ybar, B, kappa, beta_s, shape, r_vec)
                elif mod2var:
                    muR1_nh = self.non_homog_round_wetting(
                        y[0], ybar[0], B, kappa, beta_s, shape, r_vec)
                    muR2_nh = self.non_homog_round_wetting(
                        y[1], ybar[1], B, kappa, beta_s, shape, r_vec)
                    muR_nh = (muR1_nh, muR2_nh)
        else:  # homogeneous particle
            if mod1var:
                muR_nh = 0*y
            elif mod2var:
                muR_nh = (0*y[0], 0*y[1])
        return muR_nh


def step_down(x, xc, delta):
    return 0.5*(-np.tanh((x - xc)/delta) + 1)


def step_up(x, xc, delta):
    return 0.5*(np.tanh((x - xc)/delta) + 1)

1	"""This module handles properties associated with the active materials.
2	Only helper functions are defined here.
3	Diffusion functions are defined in mpet.electrode.diffusion
4	Chemical potential functions are defined in mpet.electrode.materials"""
5	import types	4✔
6	import numpy as np	4✔
7
8	import mpet.geometry as geo	4✔
9	from mpet.config import constants	4✔
10	from mpet.utils import import_function	4✔
11
12
13	class muRfuncs():	4✔
14	""" This class defines functions which describe the chemical
15	potential of active materials.
16	Each function describes a particular material. The
17	chemical potential is directly related to the open circuit voltage
18	(OCV) for solid solution materials.
19	In each function, muR_ref is the offset (non-dimensional) chemical
20	potential by which the returned value will be shifted (for
21	numerical convenience).
22	Each material function returns:
23	muR -- chemical potential
24	actR -- activity (if applicable, else None)
25	"""
26	def __init__(self, config, trode, ind=None):	4✔
27	"""config is the full dictionary of
28	parameters for the electrode particles, as made for the
29	simulations. trode is the selected electrode.
30	ind is optinally the selected particle, provided as (vInd, pInd)
31	"""
32	self.config = config	4✔
33	self.trode = trode	4✔
34	self.ind = ind	4✔
35	self.T = config['T'] # nondimensional	4✔
36	# eokT and kToe are the reference values for scalings
37	self.eokT = constants.e / (constants.k * constants.T_ref)	4✔
38	self.kToe = 1. / self.eokT	4✔
39
40	# If the user provided a filename with muRfuncs, try to load
41	# the function from there, otherwise load it from this class
42	filename = self.get_trode_param("muRfunc_filename")	4✔
43	muRfunc_name = self.get_trode_param("muRfunc")	4✔
44	if filename is None:	4✔
45	# the function will be loaded from the materials folder
46	muRfunc = import_function(None, muRfunc_name,	4✔
47	mpet_module=f"mpet.electrode.materials.{muRfunc_name}")
48	else:
49	muRfunc = import_function(filename, muRfunc_name)	×
50
51	# We have to make sure the function knows what 'self' is with
52	# the types.MethodType function
53	self.muRfunc = types.MethodType(muRfunc, self)	4✔
54
55	def get_trode_param(self, item):	4✔
56	"""
57	Shorthand to retrieve electrode-specific value
58	"""
59	value = self.config[self.trode, item]	4✔
60	# check if it is a particle-specific parameter
61	if self.ind is not None and item in self.config.params_per_particle:	4✔
62	value = value[self.ind]	4✔
63	return value	4✔
64
65	def get_muR_from_OCV(self, OCV, muR_ref):	4✔
66	return -self.eokT*OCV + muR_ref	4✔
67
68	######
69	# Helper functions
70	######
71
72	def ideal_sln(self, y):	4✔
73	""" Helper function: Should not be called directly from
74	simulation. Call a specific material instead. """
75	T = self.T	4✔
76	muR = T*np.log(y/(1-y))	4✔
77	return muR	4✔
78
79	def reg_sln(self, y, Omga):	4✔
80	""" Helper function """
81	muR_IS = self.ideal_sln(y)	4✔
82	enthalpyTerm = Omga(1-2y)	4✔
83	muR = muR_IS + enthalpyTerm	4✔
84	return muR	4✔
85
86	def graphite_2param_homog(self, y, Omga, Omgb, Omgc, EvdW):	4✔
87	""" Helper function """
88	y1, y2 = y	4✔
89	muR1 = self.reg_sln(y1, Omga)	4✔
90	muR2 = self.reg_sln(y2, Omga)	4✔
91	muR1 += Omgby2 + Omgcy2(1-y2)(1-2*y1)	4✔
92	muR2 += Omgby1 + Omgcy1(1-y1)(1-2*y2)	4✔
93	muR1 += EvdW * (30 * y1*2 (1-y1)**2)	4✔
94	muR2 += EvdW * (30 * y2*2 (1-y2)**2)	4✔
95	return (muR1, muR2)	4✔
96
97	def graphite_1param_homog(self, y, Omga, Omgb):	4✔
98	""" Helper function """
99	width = 5e-2	×
100	tailScl = 5e-2	×
101	muLtail = -tailScl1./(y*(0.85))	×
102	muRtail = tailScl1./((1-y)*(0.85))	×
103	slpScl = 0.45	×
104	muLlin = slpSclOmga4(0.26-y)step_down(y, 0.5, width)	×
105	muRlin = (slpSclOmga4(0.74-y) + Omgb)step_up(y, 0.5, width)	×
106	muR = muLtail + muRtail + muLlin + muRlin	×
107	return muR	×
108
109	def graphite_1param_homog_2(self, y, Omga, Omgb):	4✔
110	""" Helper function """
111	width = 5e-2	×
112	tailScl = 5e-2	×
113	slpScl = 0.45	×
114	muLtail = -tailScl1./(y*(0.85))	×
115	muRtail = tailScl1./((1-y)*(0.85))	×
116	muLlin = (slpSclOmga12*(0.40-y)	×
117	* step_down(y, 0.49, 0.9width)step_up(y, 0.35, width))
118	muRlin = (slpSclOmga4(0.74-y) + Omgb)step_up(y, 0.5, width)	×
119	muLMod = (0.	×
120	+ 40*(-np.exp(-y/0.015))
121	+ 0.75*(np.tanh((y-0.17)/0.02) - 1)
122	+ 1.0*(np.tanh((y-0.22)/0.040) - 1)
123	)*step_down(y, 0.35, width)
124	muR = muLMod + muLtail + muRtail + muLlin + muRlin	×
125	return muR	×
126
127	def graphite_1param_homog_3(self, y, Omga, Omgb):	4✔
128	""" Helper function with low hysteresis and soft tail """
129	width = 5e-2	×
130	tailScl = 5e-2	×
131	muLtail = -tailScl1./(y*(0.85))	×
132	muRtail = tailScl1./((1-y)*(0.85))	×
133	muRtail = 1.0e1*step_up(y, 1.0, 0.045)	×
134	muLlin = (0.15Omga12(0.40-y*0.98)	×
135	* step_down(y, 0.49, 0.9width)step_up(y, 0.35, width))
136	muRlin = (0.1Omga4(0.74-y) + 0.90Omgb)step_up(y, 0.5, 0.4width)	×
137	muLMod = (0.	×
138	+ 40*(-np.exp(-y/0.015))
139	+ 0.75*(np.tanh((y-0.17)/0.02) - 1)
140	+ 1.0*(np.tanh((y-0.22)/0.040) - 1)
141	)*step_down(y, 0.35, width)
142	muR = 0.18 + muLMod + muLtail + muRtail + muLlin + muRlin	×
143	return muR	×
144
145	def non_homog_rect_fixed_csurf(self, y, ybar, B, kappa, ywet):	4✔
146	""" Helper function """
147	N = len(y)	4✔
148	ytmp = np.empty(N+2, dtype=object)	4✔
149	ytmp[1:-1] = y	4✔
150	ytmp[0] = ywet	4✔
151	ytmp[-1] = ywet	4✔
152	dxs = 1./N	4✔
153	curv = np.diff(ytmp, 2)/(dxs**2)	4✔
154	muR_nh = -kappacurv + B(y - ybar)	4✔
155	return muR_nh	4✔
156
157	def non_homog_rect_variational(self, y, ybar, B, kappa):	4✔
158	""" Helper function """
159	# the taylor expansion at the edges is used
160	N_2 = len(y)	×
161	ytmp = np.empty(N_2+2, dtype=object)	×
162	dxs = 1./N_2	×
163	ytmp[1:-1] = y	×
164	ytmp[0] = y[0] + np.diff(y)[0]dxs + 0.5np.diff(y,2)[0]dxs*2	×
165	ytmp[-1] = y[-1] + np.diff(y)[-1]dxs + 0.5np.diff(y,2)[-1]dxs*2	×
166	curv = np.diff(ytmp, 2)/(dxs**2)	×
167	muR_nh = -kappacurv + B(y - ybar)	×
168	return muR_nh	×
169
170	def non_homog_round_wetting(self, y, ybar, B, kappa, beta_s, shape, r_vec):	4✔
171	""" Helper function """
172	dr = r_vec[1] - r_vec[0]	4✔
173	Rs = 1.	4✔
174	curv = geo.calc_curv(y, dr, r_vec, Rs, beta_s, shape)	4✔
175	muR_nh = B(y - ybar) - kappacurv	4✔
176	return muR_nh	4✔
177
178	def general_non_homog(self, y, ybar):	4✔
179	""" Helper function """
180	ptype = self.get_trode_param("type")	4✔
181	mod1var, mod2var = False, False	4✔
182	if isinstance(y, np.ndarray):	4✔
183	mod1var = True	4✔
184	N = len(y)	4✔
185	elif (isinstance(y, tuple) and len(y) == 2	4✔
186	and isinstance(y[0], np.ndarray)):
187	mod2var = True	4✔
188	N = len(y[0])	4✔
189	else:
190	raise Exception("Unknown input type")	×
191	if ("homog" not in ptype) and (N > 1):	4✔
192	shape = self.get_trode_param("shape")	4✔
193	if shape == "C3":	4✔
194	if mod1var:	4✔
195	kappa = self.get_trode_param("kappa")	4✔
196	B = self.get_trode_param("B")	4✔
197	if self.get_trode_param("type") in ["ACR"]:	4✔
198	cwet = self.get_trode_param("cwet")	4✔
199	muR_nh = self.non_homog_rect_fixed_csurf(	4✔
200	y, ybar, B, kappa, cwet)
201	elif self.get_trode_param("type") in ["ACR_Diff"]:	×
202	muR_nh = self.non_homog_rect_variational(	×
203	y, ybar, B, kappa)
204	elif mod2var:	×
205	kappa = self.get_trode_param("kappa")	×
206	B = self.get_trode_param("B")	×
207	cwet = self.get_trode_param("cwet")	×
208	muR_nh = self.non_homog_rect_fixed_csurf(	×
209	y, ybar, B, kappa, cwet)
210	elif shape in ["cylinder", "sphere"]:	4✔
211	kappa = self.get_trode_param("kappa")	4✔
212	B = self.get_trode_param("B")	4✔
213	beta_s = self.get_trode_param("beta_s")	4✔
214	r_vec = geo.get_unit_solid_discr(shape, N)[0]	4✔
215	if mod1var:	4✔
216	muR_nh = self.non_homog_round_wetting(	4✔
217	y, ybar, B, kappa, beta_s, shape, r_vec)
218	elif mod2var:	4✔
219	muR1_nh = self.non_homog_round_wetting(	4✔
220	y[0], ybar[0], B, kappa, beta_s, shape, r_vec)
221	muR2_nh = self.non_homog_round_wetting(	4✔
222	y[1], ybar[1], B, kappa, beta_s, shape, r_vec)
223	muR_nh = (muR1_nh, muR2_nh)	4✔
224	else: # homogeneous particle
225	if mod1var:	4✔
226	muR_nh = 0*y	4✔
227	elif mod2var:	4✔
228	muR_nh = (0y[0], 0y[1])	4✔
229	return muR_nh	4✔
230
231
232	def step_down(x, xc, delta):	4✔
233	return 0.5*(-np.tanh((x - xc)/delta) + 1)	×
234
235
236	def step_up(x, xc, delta):	4✔
237	return 0.5*(np.tanh((x - xc)/delta) + 1)	×

TRI-AMDD / mpet / 8476080312

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous