[mpmath] r1241 committed - Fixed test_exp_pade, added precision in exp_pade....

2 views
Skip to first unread message

mpm...@googlecode.com

unread,
Oct 30, 2011, 7:24:15 AM10/30/11
to mpmath-...@googlegroups.com
Revision: 1241
Author: mario.pernici
Date: Sun Oct 30 04:23:48 2011
Log: Fixed test_exp_pade, added precision in exp_pade.

test_exp_pade did not test at all exp_pade; put a condition on the
matrix for the similarity transformation to avoid using high precision.

Added precision in exp_pade depending on the order of the matrix,
because tests show that without doing this the precision of the result
decreases with the order of the matrix.


http://code.google.com/p/mpmath/source/detail?r=1241

Modified:
/trunk/mpmath/matrices/calculus.py
/trunk/mpmath/tests/test_linalg.py

=======================================
--- /trunk/mpmath/matrices/calculus.py Thu Feb 3 05:08:12 2011
+++ /trunk/mpmath/matrices/calculus.py Sun Oct 30 04:23:48 2011
@@ -111,9 +111,13 @@
42.0927851137247

"""
- A = ctx.matrix(A)
if method == 'pade':
- return ctx._exp_pade(A)
+ prec = ctx.prec
+ ctx.prec += 2*A.rows
+ A = ctx.matrix(A)
+ res = ctx._exp_pade(A)
+ return res
+ A = ctx.matrix(A)
prec = ctx.prec
j = int(max(1, ctx.mag(ctx.mnorm(A,'inf'))))
j += int(0.5*prec**0.5)
=======================================
--- /trunk/mpmath/tests/test_linalg.py Sat Jan 8 15:52:03 2011
+++ /trunk/mpmath/tests/test_linalg.py Sun Oct 30 04:23:48 2011
@@ -221,20 +221,26 @@
def test_exp_pade():
for i in range(3):
dps = 15
- extra = 5
+ extra = 15
mp.dps = dps + extra
dm = 0
- while not dm:
- m = randmatrix(3)
+ N = 3
+ dg = range(1,N+1)
+ a = diag(dg)
+ expa = diag([exp(x) for x in dg])
+ # choose a random matrix not close to be singular
+ # to avoid adding too much extra precision in computing
+ # m**-1 * M * m
+ while abs(dm) < 0.01:
+ m = randmatrix(N)
dm = det(m)
m = m/dm
- a = diag([1,2,3])
a1 = m**-1 * a * m
+ e2 = m**-1 * expa * m
mp.dps = dps
e1 = expm(a1, method='pade')
mp.dps = dps + extra
- e2 = m * a1 * m**-1
- d = e2 - a
+ d = e2 - e1
#print d
mp.dps = dps
assert norm(d, inf).ae(0)

Reply all
Reply to author
Forward
0 new messages