def calc_omega(cp):
cp.insert
a=[]
for i in range(len(cp)):
ptmp = []
tmp = 0
for j in range(len(cp)):
if j != i:
row = []
row.insert(0,1/(cp[i]-cp[j]))
row.insert(1,-cp[j]/(cp[i]-cp[j]))
ptmp.insert(tmp,row)
tmp += 1
p=[1]
for j in range(len(cp)-1):
p = conv(p,ptmp[j])
pint = numpy.polyint(p)
arow = []
for j in range(len(cp)):
arow.append(numpy.polyval(pint,cp[j]))
a.append(arow)
return a
python类polyint()的实例源码
def test_polyint_type(self):
# Ticket #944
msg = "Wrong type, should be complex"
x = np.ones(3, dtype=np.complex)
assert_(np.polyint(x).dtype == np.complex, msg)
msg = "Wrong type, should be float"
x = np.ones(3, dtype=np.int)
assert_(np.polyint(x).dtype == np.float, msg)
def test_polyint_type(self):
# Ticket #944
msg = "Wrong type, should be complex"
x = np.ones(3, dtype=np.complex)
assert_(np.polyint(x).dtype == np.complex, msg)
msg = "Wrong type, should be float"
x = np.ones(3, dtype=np.int)
assert_(np.polyint(x).dtype == np.float, msg)
def oscillation_indicator():
""" Generate the oscillation indicator matrix
"""
? = zeros([N+1, N+1])
for a in range(1,N+1):
?Dera = ?Der[a]
for p, m in product(range(N+1), range(N+1)):
antiderivative = polyint(?Dera[p] * ?Dera[m])
?[p,m] += antiderivative(1) - antiderivative(0)
return ?
def basis_polys():
""" Returns the basis polynomials and their derivatives and antiderivatives
"""
nodes, _, _ = quad()
? = [lagrange(nodes,eye(N+1)[i]) for i in range(N+1)]
?Der = [[polyder(?_p, m=a) for ?_p in ?] for a in range(N+1)]
?Int = [polyint(?_p) for ?_p in ?]
return ?, ?Der, ?Int
def test_polyint_type(self):
# Ticket #944
msg = "Wrong type, should be complex"
x = np.ones(3, dtype=np.complex)
assert_(np.polyint(x).dtype == np.complex, msg)
msg = "Wrong type, should be float"
x = np.ones(3, dtype=np.int)
assert_(np.polyint(x).dtype == np.float, msg)
def test_polyint_type(self):
# Ticket #944
msg = "Wrong type, should be complex"
x = np.ones(3, dtype=np.complex)
assert_(np.polyint(x).dtype == np.complex, msg)
msg = "Wrong type, should be float"
x = np.ones(3, dtype=np.int)
assert_(np.polyint(x).dtype == np.float, msg)
def test_polyint_type(self):
# Ticket #944
msg = "Wrong type, should be complex"
x = np.ones(3, dtype=np.complex)
assert_(np.polyint(x).dtype == np.complex, msg)
msg = "Wrong type, should be float"
x = np.ones(3, dtype=np.int)
assert_(np.polyint(x).dtype == np.float, msg)
def test_polyint_type(self):
# Ticket #944
msg = "Wrong type, should be complex"
x = np.ones(3, dtype=np.complex)
assert_(np.polyint(x).dtype == np.complex, msg)
msg = "Wrong type, should be float"
x = np.ones(3, dtype=np.int)
assert_(np.polyint(x).dtype == np.float, msg)
