sagemath escribir ecuaciones descargar como cientifico python math numpy linear-algebra arbitrary-precision

python - escribir - numpy álgebra lineal de precisión arbitraria



python cientifico pdf (5)

Entonces necesitas 350 dígitos de precisión. No lo conseguirás con los números flotantes IEEE (qué numpy está usando). Puede obtenerlo con el programa bc:

$ bc -l bc 1.06 Copyright 1991-1994, 1997, 1998, 2000 Free Software Foundation, Inc. This is free software with ABSOLUTELY NO WARRANTY. For details type `warranty''. scale=350 e(-800) .<hundreds of zeros>00366

Tengo una matriz 2d numpy [tamaño mediano / grande - digamos 500x500]. Quiero encontrar los valores propios del exponente a nivel de elemento de él. El problema es que algunos de los valores son bastante negativos (-800, -1000, etc.) y sus exponentes tienen un flujo inferior (lo que significa que están muy cerca de cero, por lo que numpy los trata como cero). ¿Hay alguna forma de usar precisión arbitraria en numpy?

La forma en que lo sueño:

import numpy as np np.set_precision(''arbitrary'') # <--- Missing part a = np.array([[-800.21,-600.00],[-600.00,-1000.48]]) ex = np.exp(a) ## Currently warns about underflow eigvals, eigvecs = np.linalg.eig(ex)

He buscado una solución con gmpy y mpmath en vano. Cualquier idea será bienvenida.


No tengo experiencia con Numpy en particular, pero agregar un punto decimal con una cantidad establecida de ceros puede ayudar. Por ejemplo, use 1.0000 en lugar de 1. En las secuencias de comandos python normales donde he tenido este problema, esto ha ayudado, así que a menos que su problema sea causado por una rareza en numpy y no tenga nada que ver con python, esto debería ayudar.

¡Buena suerte!


Por lo que yo sé, Numpy no admite una precisión superior a doble (float64), que es el valor predeterminado si no se especifica.

Intenta usar esto: http://code.google.com/p/mpmath/

Lista de características (entre otras)

Aritmética:

  • Números reales y complejos con precisión arbitraria
  • Ilimitado exponente de tamaños / magnitudes

SymPy puede calcular la precisión arbitraria:

from sympy import exp, N, S from sympy.matrices import Matrix data = [[S("-800.21"),S("-600.00")],[S("-600.00"),S("-1000.48")]] m = Matrix(data) ex = m.applyfunc(exp).applyfunc(lambda x:N(x, 100)) vecs = ex.eigenvects() print vecs[0][0] # eigen value print vecs[1][0] # eigen value print vecs[0][2] # eigen vect print vecs[1][2] # eigen vect

salida:

-2.650396553004310816338679447269582701529092549943247237903254759946483528035516341807463648841185335e-261 2.650396553004310816338679447269582701529092549943247237903254759946483528035516341807466621962539464e-261 [[-0.9999999999999999999999999999999999999999999999999999999999999999999999999999999999999994391176386872] [ 1]] [[1.000000000000000000000000000000000000000000000000000000000000000000000000000000000000000560882361313] [ 1]]

puedes cambiar 100 en N (x, 100) a otra precisión, pero, como probé 1000, el cálculo de eigen vect falló.


En los sistemas de 64 bits, hay un numpy.float128 numpy.float128. (También creo que hay un float96 float96 en los sistemas de 32 bits). Mientras que numpy.linalg.eig no admite flotantes de 128 bits, scipy.linalg.eig (más o menos) sí.

Sin embargo, nada de esto va a importar , a la larga. Cualquier solucionador general para un problema de valor propio va a ser iterativo, en lugar de exacto, ¡por lo que no obtendrás nada manteniendo la precisión extra! np.linalg.eig funciona para cualquier forma, pero nunca devuelve una solución exacta .

Si siempre estás resolviendo matrices de 2x2, es trivial escribir tu propio solucionador que debería ser más exacto. Voy a mostrar un ejemplo de esto al final ...

Independientemente, avanzando en contenedores de memoria sin sentido precisos:

import numpy as np import scipy as sp import scipy.linalg a = np.array([[-800.21,-600.00],[-600.00,-1000.48]], dtype=np.float128) ex = np.exp(a) print ex eigvals, eigvecs = sp.linalg.eig(ex) # And to test... check1 = ex.dot(eigvecs[:,0]) check2 = eigvals[0] * eigvecs[:,0] print ''Checking accuracy..'' print check1, check2 print (check1 - check2).dot(check1 - check2), ''<-- Should be zero''

Sin embargo, notarás que lo que obtienes es idéntico a simplemente hacer np.linalg.eig(ex.astype(np.float64) . De hecho, estoy bastante seguro de que eso es lo que hace scipy , mientras que numpy plantea un error que hacerlo en silencio. Sin embargo, podría estar bastante equivocado ...

Si no quiere usar scipy, una solución alternativa es reescalar elementos después de la exponenciación pero antes de resolver los valores propios, convertirlos en flotantes "normales", resolver los valores propios y luego refundir cosas como float128 después y reescalar.

P.ej

import numpy as np a = np.array([[-800.21,-600.00],[-600.00,-1000.48]], dtype=np.float128) ex = np.exp(a) factor = 1e300 ex_rescaled = (ex * factor).astype(np.float64) eigvals, eigvecs = np.linalg.eig(ex_rescaled) eigvals = eigvals.astype(np.float128) / factor # And to test... check1 = ex.dot(eigvecs[:,0]) check2 = eigvals[0] * eigvecs[:,0] print ''Checking accuracy..'' print check1, check2 print (check1 - check2).dot(check1 - check2), ''<-- Should be zero''

Finalmente, si solo resuelve matrices de 2x2 o 3x3, puede escribir su propio solucionador, que devolverá un valor exacto para esas formas de matrices.

import numpy as np def quadratic(a,b,c): sqrt_part = np.lib.scimath.sqrt(b**2 - 4*a*c) root1 = (-b + sqrt_part) / (2 * a) root2 = (-b - sqrt_part) / (2 * a) return root1, root2 def eigvals(matrix_2x2): vals = np.zeros(2, dtype=matrix_2x2.dtype) a,b,c,d = matrix_2x2.flatten() vals[:] = quadratic(1.0, -(a+d), (a*d-b*c)) return vals def eigvecs(matrix_2x2, vals): a,b,c,d = matrix_2x2.flatten() vecs = np.zeros_like(matrix_2x2) if (b == 0.0) and (c == 0.0): vecs[0,0], vecs[1,1] = 1.0, 1.0 elif c != 0.0: vecs[0,:] = vals - d vecs[1,:] = c elif b != 0: vecs[0,:] = b vecs[1,:] = vals - a return vecs def eig_2x2(matrix_2x2): vals = eigvals(matrix_2x2) vecs = eigvecs(matrix_2x2, vals) return vals, vecs a = np.array([[-800.21,-600.00],[-600.00,-1000.48]], dtype=np.float128) ex = np.exp(a) eigvals, eigvecs = eig_2x2(ex) # And to test... check1 = ex.dot(eigvecs[:,0]) check2 = eigvals[0] * eigvecs[:,0] print ''Checking accuracy..'' print check1, check2 print (check1 - check2).dot(check1 - check2), ''<-- Should be zero''

Esta devuelve una solución verdaderamente exacta, pero solo funcionará para matrices de 2x2. ¡Sin embargo, es la única solución que realmente se beneficia con la precisión adicional!