Python math 模块，erf() 实例源码
我们从Python开源项目中，提取了以下39个代码示例，用于说明如何使用math.erf()。
def get_winrates(user_id_1: int, user_id_2: int, ndchar: NDChar, amplified: bool) -> tuple or None:
stats_1 = await get_character_stats(user_id_1, ndchar, amplified)
stats_2 = await get_character_stats(user_id_2, ndchar, amplified)
if not stats_1.has_wins or not stats_2.has_wins:
return None
m2_minus_m1 = stats_2.mean - stats_1.mean
sum_var = stats_1.var + stats_2.var
erf_arg = m2_minus_m1 / math.sqrt(2*sum_var)
if m2_minus_m1 > 0:
winrate_of_1_if_both_finish = (1.0 + math.erf(erf_arg))/2.0
else:
winrate_of_1_if_both_finish = (1.0 - math.erf(-erf_arg))/2.0
both_finish_prob = stats_1.winrate * stats_2.winrate
neither_finish_prob = (1-stats_1.winrate)*(1-stats_2.winrate)
winrate_of_1 = winrate_of_1_if_both_finish*both_finish_prob + (stats_1.winrate - both_finish_prob)
winrate_of_2 = (1.0-winrate_of_1_if_both_finish)*both_finish_prob + (stats_2.winrate - both_finish_prob)
return winrate_of_1, winrate_of_2, neither_finish_prob
def blackcox_pd(equity, extasset, sigma):
"""Compute the probability of default for external assets following a
Geometric Brownian Motion and the Black and Cox model.
Parameters:
equity (float): equity
extasset (float): external assets
sigma (float): volatility of the Geometric Browninan Motion
Returns:
probability of default
"""
if equity <= 0.0:
return 1.0
if equity >= extasset:
return 0.0
else:
#return 1 + (- 1/2 * (1 + math.erf((-math.log(1 - equity/extasset) - sigma**2/2) /
# (math.sqrt(2) * sigma)) )
# + (extasset/equity)/2 * (1 + math.erf((math.log(1 - equity/extasset) - sigma**2/2) /
# (math.sqrt(2) * sigma)) ) )
return (1/2 * (1 + math.erf((math.log(1 - equity/extasset) + sigma**2/2) /
(math.sqrt(2) * sigma)) ) +
(extasset/(extasset - equity))/2 * (1 + math.erf((math.log(1 - equity/extasset) - sigma**2/2) /
(math.sqrt(2) * sigma)) ) )
def black_scholes_numba_opt(price, strike, t, mr, sig_sig_two, vol, call, put):
P = float( price [0] )
S = strike [0]
T = t [0]
a = log(P / S)
b = T * mr[0]
z = T * sig_sig_two[0]
c = 0.25 * z
y = 1./sqrt(z)
w1 = (a - b + c) * y
w2 = (a - b - c) * y
d1 = 0.5 + 0.5 * erf(w1)
d2 = 0.5 + 0.5 * erf(w2)
Se = exp(b) * S
res = P * d1 - Se * d2
call [0] = res
put [0] = res - P + Se
def black_scholes_numba_opt(price, strike, t, mr, sig_sig_two):
P = price
S = strike
T = t
a = log(P / S)
b = T * mr
z = T * sig_sig_two
c = 0.25 * z
y = 1./sqrt(z)
w1 = (a - b + c) * y
w2 = (a - b - c) * y
d1 = 0.5 + 0.5 * erf(w1)
d2 = 0.5 + 0.5 * erf(w2)
Se = exp(b) * S
r = P * d1 - Se * d2
return complex(r, r - P + Se)
def black_scholes_numba_opt(price, strike, t, mr, sig_sig_two):
P = price
S = strike
T = t
a = log(P / S)
b = T * mr
z = T * sig_sig_two
c = 0.25 * z
y = 1./sqrt(z)
w1 = (a - b + c) * y
w2 = (a - b - c) * y
d1 = 0.5 + 0.5 * erf(w1)
d2 = 0.5 + 0.5 * erf(w2)
Se = exp(b) * S
r = P * d1 - Se * d2
return complex(r, r - P + Se)
def __get_snd_score(self, x):
return (1.0 + math.erf(x / math.sqrt(2.0))) / 2.0
def lognormal_pd(equity, extasset, sigma):
"""Compute the probability of default for external assets following a
Geometric Brownian Motion.
Such probability of default is the correct probability of default to use
for the NEVA interbank valuation function with external assets following a
Geometric Brownian Motion, implemented in `exante_en_merton_gbm`. See
Eq. (16a) in [1].
Parameters:
equity (float): equity
extasset (float): external assets
sigma (float): volatility of the Geometric Browninan Motion
Returns:
probability of default
References:
[1] P. Barucca, M. Bardoscia, F. Caccioli, M. D'Errico, G. Visentin,
S. Battiston, G. Caldarelli. Network Valuation in Financial Systems,
https://arxiv.org/abs/1606.05164
"""
if equity >= extasset:
#print('wow1')
return 0.0
else:
#print('wow2', (sigma**2 / 2 + math.log(1.0 - equity/extasset)) / (math.sqrt(2) * sigma))
#print('wow2', sigma**2 / 2, math.log(1.0 - equity/extasset), (math.sqrt(2) * sigma))
return 1/2 * (1 + math.erf((sigma**2 / 2 + math.log(1 - equity/extasset))
/ (math.sqrt(2) * sigma)))
def lognormal_cav_aext(equity, extasset, ibliabtot, sigma):
"""Compute the conditional expected endogenous recovery for external
assets following a Geometric Brownian Motion.
Such conditional expected endogenous recovery is the correct conditional
expected endogenous recovery to use for the NEVA interbank valuation
function with external assets following a Geometric Brownian Motion,
implemented in `end_lognormal_dr`. See Eq. (16b) in [1].
Parameters:
equity (float): equity
extasset (float): external assets
ibliabtot (float): total interbank liabilities
sigma (float): volatility of the Geometric Browninan Motion
Returns:
conditional expected endogenous recovery
References:
[1] P. Barucca, M. Bardoscia, F. Caccioli, M. D'Errico, G. Visentin,
S. Battiston, G. Caldarelli. Network Valuation in Financial Systems,
https://arxiv.org/abs/1606.05164
"""
out = 0.0
if extasset > equity:
tmp_sigma_1 = sigma**2 / 2
tmp_sigma_2 = math.sqrt(2) * sigma
out += 1/2 * (1 + math.erf((math.log(1 - equity/extasset) - tmp_sigma_1)
/ tmp_sigma_2))
if extasset > equity + ibliabtot:
out -= 1/2 * (1 + math.erf((math.log(1 - (equity + ibliabtot)/extasset) -
tmp_sigma_1)
/ tmp_sigma_2))
return extasset * out
def compute_relative_viscosity(self, state):
phi = state.get_crystal_fraction()
if self._strain_rate == 1.0:
# needle-like B particles from Cimarelli et al., 2011
# self.phi_max_2 = 0.44
delta_1 = 4.45
gama_1 = 8.55
phi_star_1 = 0.28
epsilon_1 = 0.001
f = (1. - epsilon_1) * math.erf(min(25., (
(math.sqrt(math.pi) / (2. * (1. - epsilon_1))) * (phi / phi_star_1) * (
1. + (math.pow((phi / phi_star_1), gama_1))))))
relative_viscosity = (1. + math.pow((phi / phi_star_1), delta_1)) / (
math.pow((1. - f), (2.5 * phi_star_1)))
return relative_viscosity
if self._strain_rate == 0.0001:
# needle-like, B particles from Cimarelli et al., 2011
# self.phi_max = 0.36
delta_1 = 7.5
gama_1 = 5.5
phi_star_1 = 0.26
epsilon_1 = 0.0002
f = (1. - epsilon_1) * math.erf(min(25., (
(math.sqrt(math.pi) / (2. * (1. - epsilon_1))) * (phi / phi_star_1) * (
1. + (math.pow((phi / phi_star_1), gama_1))))))
relative_viscosity = (1. + math.pow((phi / phi_star_1), delta_1)) / (
math.pow((1. - f), (2.5 * phi_star_1)))
return relative_viscosity
def compute_relative_viscosity(self, state):
phi = state.get_crystal_fraction()
if self._strain_rate == 1.0:
# for spheres, A particles from Cimarelli et al., 2011
# self.phi_max = 0.61,
delta_1 = 11.4
gama_1 = 1.6
phi_star_1 = 0.67
epsilon_1 = 0.01
f = (1. - epsilon_1) * math.erf(min(25., (
(math.sqrt(math.pi) / (2. * (1. - epsilon_1))) * (phi / phi_star_1) * (
1. + (math.pow((phi / phi_star_1), gama_1))))))
relative_viscosity = (1. + math.pow((phi / phi_star_1), delta_1)) / (
math.pow((1. - f), (2.5 * phi_star_1)))
return relative_viscosity
if self._strain_rate == 0.0001:
# spheres A particles from Cimarelli et al., 2011
# self.phi_max_1 = 0.54,
delta_1 = 11.48
gama_1 = 1.52
phi_star_1 = 0.62
epsilon_1 = 0.005
f = (1. - epsilon_1) * math.erf(min(25., (
(math.sqrt(math.pi) / (2. * (1. - epsilon_1))) * (phi / phi_star_1) * (
1. + (math.pow((phi / phi_star_1), gama_1))))))
relative_viscosity = (1. + math.pow((phi / phi_star_1), delta_1)) / (
math.pow((1. - f), (2.5 * phi_star_1)))
return relative_viscosity
def normal_cdf(x, mu=0,sigma=1):
return (1 + math.erf((x - mu) / math.sqrt(2) / sigma)) / 2
def _step(self, action):
# the first element of action is the actual current action
current_action = action[0]
observation, reward, done, info = self.cartpole._step(current_action)
if not done:
# We add the newly predicted observations to the list before checking predictions
# in order to give the agent a chance to predict the observations that they
# are going to get _this_ round.
self.predicted_observations.append(action[1:])
if self.iteration > TIME_BEFORE_BONUS_ALLOWED:
for i in xrange(min(NUM_PREDICTED_OBSERVATIONS, len(self.predicted_observations))):
l2dist = np.sqrt(np.sum(np.square(np.subtract(
self.predicted_observations[-(i + 1)][i],
observation
))))
bonus = CORRECT_PREDICTION_BONUS * (1 - math.erf(l2dist))
reward += bonus
self.iteration += 1
return observation, reward, done, info
def erf(x):
#Python 2.7 defines math.erf(), but we need to cater for older versions.
a = 8*(math.pi-3)/(3*math.pi*(4-math.pi))
x2 = x*x
y = -x2 * (4/math.pi + a*x2) / (1 + a*x2)
return math.copysign(math.sqrt(1 - math.exp(y)), x)
def erf_inv(x):
# Above erf formula inverted analytically
a = 8*(math.pi-3)/(3*math.pi*(4-math.pi))
y = math.log(1-x*x)
z = 2/(math.pi*a) + y/2
return math.copysign(math.sqrt(math.sqrt(z*z - y/a) - z), x)
def phi(q):
# Cumlative distribution function for the standard Gaussian law: quantile -> probability
return 0.5*(1+erf(q/math.sqrt(2)))
def lognorm_opt(med=0, spread=0, step=.01, just_norm=False, shaking=False):
p_norm = (math.erf((shaking-med)/(math.sqrt(2) * spread)) + 1)/2
return p_norm * 100
def phi(x):
"""Cumulative distribution function for the standard normal distribution
Taken from python math module documentation"""
return (1.0 + math.erf(x / math.sqrt(2.0))) / 2.0
def _erf(x):
''' Currently not needed, but might be useful for a CUDA implementation '''
ax = _np.abs(x)
if ax < 0.5:
t = x*x
top = ((((.771058495001320e-04*t-.133733772997339e-02)*t+.323076579225834e-01)*t+.479137145607681e-01)*t+.128379167095513e+00) + 1.0
bot = ((.301048631703895e-02*t+.538971687740286e-01)*t+.375795757275549e+00)*t + 1.0
return x * (top / bot)
if ax < 4.0:
top = ((((((-1.36864857382717e-07*ax+5.64195517478974e-01)*ax+7.21175825088309e+00)*ax+4.31622272220567e+01)*ax+1.52989285046940e+02)*ax+3.39320816734344e+02)*ax+4.51918953711873e+02)*ax + 3.00459261020162e+02
bot = ((((((1.0*ax+1.27827273196294e+01)*ax+7.70001529352295e+01)*ax+2.77585444743988e+02)*ax+6.38980264465631e+02)*ax+9.31354094850610e+02)*ax+7.90950925327898e+02)*ax + 3.00459260956983e+02
erf = 0.5 + (0.5 - _np.exp(-x * x) * top / bot)
if x < 0.0:
erf = -erf
return erf
if ax < 5.8:
x2 = x*x
t = 1.0 / x2
top = (((2.10144126479064e+00*t+2.62370141675169e+01)*t+2.13688200555087e+01)*t+4.65807828718470e+00)*t + 2.82094791773523e-01
bot = (((9.41537750555460e+01*t+1.87114811799590e+02)*t+9.90191814623914e+01)*t+1.80124575948747e+01)*t + 1.0
erf = (.564189583547756e0 - top / (x2 * bot)) / ax
erf = 0.5 + (0.5 - _np.exp(-x2) * erf)
if x < 0.0:
erf = -erf
return erf
return _np.sign(x)
def _gaussian_integral(x, mu, sigma):
sq_norm = 0.70710678118654757 / sigma # sq_norm = sqrt(0.5/sigma**2)
d = x - mu
return 0.5 * (_math.erf((d + 0.5) * sq_norm) - _math.erf((d - 0.5) * sq_norm))
def stdnormdist(x):
#Cumulative standard normal distribution function
#Lifted from https://docs.python.org/3.2/library/math.html
return (1.0 + math.erf(x / math.sqrt(2.0))) / 2.0
def _step(self, action):
# the first element of action is the actual current action
current_action = action[0]
observation, reward, done, info = self.cartpole._step(current_action)
if not done:
# We add the newly predicted observations to the list before checking predictions
# in order to give the agent a chance to predict the observations that they
# are going to get _this_ round.
self.predicted_observations.append(action[1:])
if self.iteration > TIME_BEFORE_BONUS_ALLOWED:
for i in xrange(min(NUM_PREDICTED_OBSERVATIONS, len(self.predicted_observations))):
l2dist = np.sqrt(np.sum(np.square(np.subtract(
self.predicted_observations[-(i + 1)][i],
observation
))))
bonus = CORRECT_PREDICTION_BONUS * (1 - math.erf(l2dist))
reward += bonus
self.iteration += 1
return observation, reward, done, info
def normal(x):
y = 0.5 + 0.5*math.erf(x / 1.414213)
return y
# how much of the distribution is < x?
def real_normal(x):
y = 0.5 + 0.5*math.erf(x / 1.414213)
return y
def phi(n):
"""
Computes the cumulative distribution function of the standard,
normal distribution for values <= n.
>>> round(phi(0), 3)
0.5
>>> round(phi(-1), 3)
0.159
>>> round(phi(+1), 3)
0.841
"""
return (1+erf(n/sqrt(2)))/2
def percentage_of_area_under_std_normal_curve_from_zcore(z_score):
return .5 * (math.erf(z_score / 2 ** .5) + 1)
def StandardNormalCdf(x):
"""Evaluates the CDF of the standard Normal distribution.
See http://en.wikipedia.org/wiki/Normal_distribution
#Cumulative_distribution_function
Args:
x: float
Returns:
float
"""
return (math.erf(x / ROOT2) + 1) / 2
def get_winrate(rating_1: Rating, rating_2: Rating):
delta_mu = rating_1.mu - rating_2.mu
if delta_mu >= 0:
beta = trueskill.global_env().beta
denom = sqrt(2 * (2 * beta * beta + rating_1.sigma * rating_1.sigma + rating_2.sigma * rating_2.sigma))
return (erf(delta_mu/denom) + 1.0)/2.0
else:
return 1.0 - get_winrate(rating_2, rating_1)
def StandardNormalCdf(x):
"""Evaluates the CDF of the standard Normal distribution.
See http://en.wikipedia.org/wiki/Normal_distribution
#Cumulative_distribution_function
Args:
x: float
Returns:
float
"""
return (math.erf(x / ROOT2) + 1) / 2
def StandardNormalCdf(x):
"""Evaluates the CDF of the standard Normal distribution.
See http://en.wikipedia.org/wiki/Normal_distribution
#Cumulative_distribution_function
Args:
x: float
Returns:
float
"""
return (math.erf(x / ROOT2) + 1) / 2
def norm_cdf(x, mu=0, sigma=1):
"""
Parameters
----------
x :
mu : float
distribution's mean
sigma : float
distribution's standard deviation
Returns
-------
float
pdf or cdf value, depending on input flag ``f``
Notes
-----
http://stackoverflow.com/questions/809362/how-to-calculate-cumulative-normal-distribution-in-python
Examples
--------
Compares total absolute error for 100 values
>>> from scipy.stats import norm
>>> sum( [abs(Util.norm_cdf(x) - norm.cdf(x)) for x in range(100)])
3.3306690738754696e-16
"""
y = 0.5 * (1 - math.erf(-(x - mu)/(sigma * math.sqrt(2.0))))
if y > 1: y = 1
return y
def normalCDF(x, mu = 0, sigma = 1.0):
return (1 + math.erf((x - mu) / math.sqrt(2) / sigma)) / 2
def _step(self, action):
# the first element of action is the actual current action
current_action = action[0]
observation, reward, done, info = self.cartpole._step(current_action)
if not done:
# We add the newly predicted observations to the list before checking predictions
# in order to give the agent a chance to predict the observations that they
# are going to get _this_ round.
self.predicted_observations.append(action[1:])
if self.iteration > TIME_BEFORE_BONUS_ALLOWED:
for i in xrange(min(NUM_PREDICTED_OBSERVATIONS, len(self.predicted_observations))):
l2dist = np.sqrt(np.sum(np.square(np.subtract(
self.predicted_observations[-(i + 1)][i],
observation
))))
bonus = CORRECT_PREDICTION_BONUS * (1 - math.erf(l2dist))
reward += bonus
self.iteration += 1
return observation, reward, done, info
def black_scholes ( nopt, price, strike, t, rate, vol, call, put ):
mr = -rate
sig_sig_two = vol * vol * 2
for i in range(nopt):
P = float( price [i] )
S = strike [i]
T = t [i]
a = log(P / S)
b = T * mr
z = T * sig_sig_two
c = 0.25 * z
y = invsqrt(z)
w1 = (a - b + c) * y
w2 = (a - b - c) * y
d1 = 0.5 + 0.5 * erf(w1)
d2 = 0.5 + 0.5 * erf(w2)
Se = exp(b) * S
call [i] = P * d1 - Se * d2
put [i] = call [i] - P + Se
def black_scholes( nopt, price, strike, t, rate, vol, call, put):
mr = -rate
sig_sig_two = vol * vol * 2
for i in range(nopt):
P = float( price[i] )
S = strike [i]
T = t [i]
a = log(P / S)
b = T * mr
z = T * sig_sig_two
c = 0.25 * z
y = 1./sqrt(z)
w1 = (a - b + c) * y
w2 = (a - b - c) * y
d1 = 0.5 + 0.5 * erf(w1)
d2 = 0.5 + 0.5 * erf(w2)
Se = exp(b) * S
call [i] = P * d1 - Se * d2
put [i] = call [i] - P + Se
def F(x, u, std_dev):
return (1.0 + erf((x-u) / (std_dev * sqrt(2.0)))) / 2.0
def F(x, u, std_dev):
return (1.0 + erf((x-u) / (std_dev * sqrt(2.0)))) / 2.0
def F(x, u, std_dev):
return 0.5*(1 + erf((x-u)/(std_dev*(2**0.5))))
def __call__(self, temp, strain_rate, material, **overrides):
"""Solves for the yield stress and flow stress at the given condition
Args:
temp(float): Temperature, in degrees Kelvin
strain_rate(float): Strain rate, in sec**-1
material(str): Key for material type
Keyword Args:
**overrides(dict): Passed as a chain of keyword arguments. These
arguments override any material property
Return:
flow_stress(float): Flow stress in ??Pa??
yield_stress(float): Yield stress in ??Pa??
"""
g_modu, t_norm, psi_norm = self.__pre__(temp, strain_rate, material,
**overrides)
s_o = self.get_option('s_o')
s_inf = self.get_option('s_inf')
kappa = self.get_option('kappa')
beta = self.get_option('beta')
y_o = self.get_option('y_o')
y_inf = self.get_option('y_inf')
y_1 = self.get_option('y_1')
y_2 = self.get_option('y_2')
if psi_norm == 0.0:
erf_psi_norm = 0.0
else:
erf_psi_norm = erf(kappa * t_norm * log(psi_norm**(-1)))
# end
glide_flow_stress = s_o - (s_o - s_inf) * erf_psi_norm
shock_flow_stress = s_o * psi_norm**beta
glide_yield_stress = y_o - (y_o - y_inf) * erf_psi_norm
shock_yield_stress = y_1 * psi_norm**y_2
flow_stress = max((glide_flow_stress, shock_flow_stress))
yield_stress = max((glide_yield_stress,
min((shock_yield_stress, shock_flow_stress))))
flow_stress *= g_modu
yield_stress *= g_modu
return flow_stress, yield_stress