本文整理汇总了Python中chemreac.ReactionDiffusion.calc_efield方法的典型用法代码示例。如果您正苦于以下问题:Python ReactionDiffusion.calc_efield方法的具体用法?Python ReactionDiffusion.calc_efield怎么用?Python ReactionDiffusion.calc_efield使用的例子?那么, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在类chemreac.ReactionDiffusion
的用法示例。
在下文中一共展示了ReactionDiffusion.calc_efield方法的1个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的Python代码示例。
示例1: integrate_rd
# 需要导入模块: from chemreac import ReactionDiffusion [as 别名]
# 或者: from chemreac.ReactionDiffusion import calc_efield [as 别名]
def integrate_rd(D=-3e-1, t0=0.0, tend=7., x0=0.1, xend=1.0, N=1024,
base=0.5, offset=0.25, nt=25, geom='f',
logt=False, logy=False, logx=False, random=False,
nstencil=3, lrefl=False, rrefl=False,
num_jacobian=False, method='bdf', plot=False,
savefig='None', atol=1e-6, rtol=1e-6, random_seed=42,
surf_chg=(0.0, 0.0), sigma_q=101, sigma_skew=0.5,
verbose=False, eps_rel=80.10, use_log2=False):
"""
A negative D (diffusion coefficent) denotes:
mobility := -D
D := 0
A positive D calculates mobility from Einstein-Smoluchowski relation
"""
assert 0 <= base and base <= 1
assert 0 <= offset and offset <= 1
if random_seed:
np.random.seed(random_seed)
n = 2
if D < 0:
mobility = -D
D = 0
else:
mobility = electrical_mobility_from_D(D, 1, 298.15)
print(D, mobility)
# Setup the grid
logb = (lambda arg: log(arg)/log(2)) if use_log2 else log
_x0 = logb(x0) if logx else x0
_xend = logb(xend) if logx else xend
x = np.linspace(_x0, _xend, N+1)
if random:
x += (np.random.random(N+1)-0.5)*(_xend-_x0)/(N+2)
# Setup the system
stoich_active = []
stoich_prod = []
k = []
rd = ReactionDiffusion(
n, stoich_active, stoich_prod, k, N,
D=[D, D],
z_chg=[1, -1],
mobility=[mobility, -mobility],
x=x,
geom=geom,
logy=logy,
logt=logt,
logx=logx,
nstencil=nstencil,
lrefl=lrefl,
rrefl=rrefl,
auto_efield=True,
surf_chg=surf_chg,
eps_rel=eps_rel, # water at 20 deg C
faraday_const=1,
vacuum_permittivity=1,
use_log2=use_log2
)
# Initial conditions
sigma = (xend-x0)/sigma_q
sigma = [(1-sigma_skew)*sigma, sigma_skew*sigma]
y0 = np.vstack(pair_of_gaussians(
rd.xcenters, [base+offset, base-offset], sigma, logy, logx, geom, use_log2)).transpose()
if logy:
y0 = sigm(y0)
if plot:
# Plot initial E-field
import matplotlib.pyplot as plt
plt.figure(figsize=(6, 10))
rd.calc_efield((rd.expb(y0) if logy else y0).flatten())
plt.subplot(4, 1, 3)
plt.plot(rd.xcenters, rd.efield, label="E at t=t0")
plt.plot(rd.xcenters, rd.xcenters*0, label="0")
# Run the integration
tout = np.linspace(t0, tend, nt)
integr = run(rd, y0, tout,
atol=atol, rtol=rtol, sigm_damp=True,
C0_is_log=logy,
with_jacobian=(not num_jacobian), method=method)
Cout = integr.Cout
if verbose:
print(integr.info)
# Plot results
if plot:
def _plot(y, ttl=None, **kwargs):
plt.plot(rd.xcenters, y, **kwargs)
plt.xlabel((('log_%s({})' % ('2' if use_log2 else 'e')) if logx else '{}').format('x / m'))
plt.ylabel('C / M')
if ttl:
plt.title(ttl)
for i in range(nt):
#.........这里部分代码省略.........