arbitrary precision linear algebra
Hello. I have a written Python program which currently uses numpy to perform linear algebra operations. Specifically, I do matrix*matrix, matrix*vector, numpy.linalg.inv(matrix), and linalg.eig(matrix) operations. Now I am interested in allowing arbitrary precision. I have tried gmpy, bigfloat, mpmath, and decimal but I have been unable to easily implement any with my current program. I suspect I have to change some commands but I am unsure what. My question is which of the arbitrary precision implementations will most easily handle linear algebra? I don't care about speed, just ease of use. Online tutorials for arbitrary precision linear algebra operations would be useful. For example, it looks like mpmath can handle matrix operations http://fredrik-j.blogspot.com/search?q=matrix but I was unable to find a clear tutorial. The tutorials for most of the arbitrary precision implementations demonstrate simple scalar examples. Thanks in advance -- http://mail.python.org/mailman/listinfo/python-list
Re: arbitrary precision linear algebra
On Mar 2, 8:42 am, Ben123 wrote: > Hello. I have a written Python program which currently uses numpy to > perform linear algebra operations. Specifically, I do matrix*matrix, > matrix*vector, numpy.linalg.inv(matrix), and linalg.eig(matrix) > operations. Now I am interested in allowing arbitrary precision. I > have tried gmpy, bigfloat, mpmath, and decimal but I have been unable > to easily implement any with my current program. I suspect I have to > change some commands but I am unsure what. > > My question is which of the arbitrary precision implementations will > most easily handle linear algebra? I don't care about speed, just ease > of use. Online tutorials for arbitrary precision linear algebra > operations would be useful. > > For example, it looks like mpmath can handle matrix > operationshttp://fredrik-j.blogspot.com/search?q=matrix > but I was unable to find a clear tutorial. The tutorials for most of > the arbitrary precision implementations demonstrate simple scalar > examples. > > Thanks in advance Hello again. I forgot to mention I'm using Python 2.6.4 mpmath 0.17 bigfloat 0.2.1 gmp 5.01 gmpy2 2.0.0a1 mpfr 3.0.0 all on Ubuntu x64 -- http://mail.python.org/mailman/listinfo/python-list
Re: arbitrary precision linear algebra
On Mar 2, 9:04 am, Arthur Mc Coy <[email protected]> wrote: > What do you mean by "arbitrary precision" ? Each method of calculating > of something has its own precision... If you are unfamiliar with arbitrary precision, I'm referring to http://en.wikipedia.org/wiki/Arbitrary-precision_arithmetic Suppose I find the eigenvalues of a matrix and the eigenvalues range from 1 to 0.0001. This can be handled by numpy in Python because the smallest eigenvalue is larger than then numerical precision of 1E-19. However, if the range of eigenvalues is 1 to 1E-40, then I will need to increase the precision of all calculations leading up to finding the eigenvalues. I am working with complex valued matrices and I expect to get real eigenvalues back (based on the physics of the system). The precision of numpy is apparent from the imaginary component of the eigenvalues I find, currently 1E-19 or 1E-20. I need better precision for small eigenvalues. In case you are curious, the complex-valued matrices are 20x20. Thanks -- http://mail.python.org/mailman/listinfo/python-list
Re: arbitrary precision linear algebra
On Mar 2, 10:21 am, Arthur Mc Coy <[email protected]> wrote: > On Mar 2, 5:26 pm, Ben123 wrote: > > > > > > > > > > > On Mar 2, 9:04 am, Arthur Mc Coy <[email protected]> wrote: > > > > What do you mean by "arbitrary precision" ? Each method of calculating > > > of something has its own precision... > > > If you are unfamiliar with arbitrary precision, I'm referring > > tohttp://en.wikipedia.org/wiki/Arbitrary-precision_arithmetic > > > Suppose I find the eigenvalues of a matrix and the eigenvalues range > > from 1 to 0.0001. This can be handled by numpy in Python because the > > smallest eigenvalue is larger than then numerical precision of 1E-19. > > However, if the range of eigenvalues is 1 to 1E-40, then I will need > > to increase the precision of all calculations leading up to finding > > the eigenvalues. > > > I am working with complex valued matrices and I expect to get real > > eigenvalues back (based on the physics of the system). The precision > > of numpy is apparent from the imaginary component of the eigenvalues I > > find, currently 1E-19 or 1E-20. I need better precision for small > > eigenvalues. > > > In case you are curious, the complex-valued matrices are 20x20. > > > Thanks > > You probably have to change the method of finding eigenvalues. > Which one do you use? Power or algebraic ? I'm not sure what you mean by this. As I mentioned, in Python I am using linalg.eig() from numpy on complex matrices. I have not investigated how this is implemented. > Do you use Gaussian method to simplify matrices ? No > > Languages can't support infinitely large or small numbers, so try to > multiply the inner variables by 10^n to increase their values if this > will not involve on the method. For example, I did this when was > calculating geometric means of computer benchmarks. Currently I have values between 1 and 1E-300 (not infinitely small). I don't see how scaling by powers of 10 will increase precision. > In such way you will be storing the number of zeros as n. Are you saying python cares whether I express a number as 0.001 or scaled by 10^5 to read 100? If this is the case, I'm still stuck. I need the full range of eigenvalues from 1 to 1E-300, so the entire range could be scaled by 1E300 but I would still need better precision than 1E19 > > Yes, interesting what are you calculating. -- http://mail.python.org/mailman/listinfo/python-list
Re: arbitrary precision linear algebra
On Mar 2, 11:28 am, Robin Becker wrote: > On 02/03/2011 16:39, Ben123 wrote: > ... > > > > > > > > >> Languages can't support infinitely large or small numbers, so try to > >> multiply the inner variables by 10^n to increase their values if this > >> will not involve on the method. For example, I did this when was > >> calculating geometric means of computer benchmarks. > > > Currently I have values between 1 and 1E-300 (not infinitely small). I > > don't see how scaling by powers of 10 will increase precision. > > >> In such way you will be storing the number of zeros as n. > > > Are you saying python cares whether I express a number as 0.001 or > > scaled by 10^5 to read 100? If this is the case, I'm still stuck. I > > need the full range of eigenvalues from 1 to 1E-300, so the entire > > range could be scaled by 1E300 but I would still need better precision > > than 1E19 > > ... > > If you enter a number as 1e-19 then python will treat as a float; by default I > think that float is IEEE double precision so you're getting a 48 bit mantissa > (others may have better details). That means you've already lost any idea of > arbitrary precision. > > When you say you have numbers like 1E-300 are those actually numerically zero > or > have you some valid inputs that vary over a huge range. It should be possible > to > compute determinants/inverses etc to arbitrary precision as those are known to > be polynomial functions of the input elements. However, eigenvalues are not. I meant that I expect the eigenvalues to be in a range from 1 to 1E-300. The values in the matrix are from sine and cosine values and have range [-10,10] for the real and imaginary parts. Thus, I could specify the matrix elements to arbitrary precision prior to diagonalization. However, I'm looking for the resulting eigenvalues to have precision to 1E-300 > > eg > > [0 2] > [1 0] > > has eigenvalues +/- sqrt(2) so even though you can represent the matrix in > finite precision the eigenvalues require infinite precision. Here's an example of a matrix I am diagonalizing: from numpy import * exmpl=array([[ 9.39582700e-04 +1.21760986e-21j, 3.32958159e-04 -6.03359588e-04j, \ -3.71551527e-04 -1.77143102e-04j, 2.87113322e-04 -9.84374562e-04j], \ [ 3.32958159e-04 +6.03359588e-04j, 5.05441331e-04 +6.80604208e-21j, \ -1.79123354e-05 -3.01368276e-04j, 7.33866807e-04 -1.64459142e-04j], \ [ -3.71551527e-04 +1.77143102e-04j, -1.79123354e-05 +3.01368276e-04j, \ 1.80324968e-04 -2.63622461e-21j, 7.20508934e-05 +4.43394732e-04j], \ [ 2.87113322e-04 +9.84374562e-04j, 7.33866807e-04 +1.64459142e-04j, \ 7.20508934e-05 -4.43394732e-04j, 1.11903651e-03 -6.61744490e-21j]]) The eigenvalues I get back using linalg.eig(exmpl)[0] are 2.74438550e-03 +7.67536143e-20j 6.54324852e-12 +2.03438929e-20j 1.39471366e-13 +4.68335525e-20j 1.76591707e-12 +2.20243218e-20j]) Here I would want to have better precision in the eigenvalues, a smaller imaginary component. To use your example, I'd like to increase the number of digits shown in the eigenvalue: from numpy import * test=array([[0, 2],[1, 0]]) linalg.eig(test)[0] but using import mpmath from mpmath import * mp.dps = 50 from numpy import * test=array([[0, 2],[1, 0]]) linalg.eig(test)[0] does not help > > Eigenvalues are roots of a polynomial in the elements and root solving may > require an infinite number of steps so it will be difficult with arbitrary > matrices to keep arbitrary precision. > -- > Robin Becker -- http://mail.python.org/mailman/listinfo/python-list
Re: arbitrary precision linear algebra
On Mar 2, 12:22 pm, geremy condra wrote: > On Wed, Mar 2, 2011 at 10:21 AM, geremy condra wrote: > > On Wed, Mar 2, 2011 at 6:42 AM, Ben123 wrote: > >> Hello. I have a written Python program which currently uses numpy to > >> perform linear algebra operations. Specifically, I do matrix*matrix, > >> matrix*vector, numpy.linalg.inv(matrix), and linalg.eig(matrix) > >> operations. Now I am interested in allowing arbitrary precision. I > >> have tried gmpy, bigfloat, mpmath, and decimal but I have been unable > >> to easily implement any with my current program. I suspect I have to > >> change some commands but I am unsure what. > > >> My question is which of the arbitrary precision implementations will > >> most easily handle linear algebra? I don't care about speed, just ease > >> of use. Online tutorials for arbitrary precision linear algebra > >> operations would be useful. > > >> For example, it looks like mpmath can handle matrix operations > >>http://fredrik-j.blogspot.com/search?q=matrix > >> but I was unable to find a clear tutorial. The tutorials for most of > >> the arbitrary precision implementations demonstrate simple scalar > >> examples. > > >> Thanks in advance > > > Have you looked at Sage[0]? I don't know for a fact, but you should be > > able to define a matrix over RealField(precision_in_bits) and then > > take the eigenvalue of it. I don't know if it will actually produce > > the precision you need though. > > > Geremy Condra > > Apologies, forgot the links: > > http://www.sagemath.org/doc/constructions/linear_algebra.htmlhttp://www.sagemath.org/doc/reference/sage/rings/complex_field.html > > Geremy Condra I've already implemented the algorithm in Mathematica for exactly this reason - increasing precision there is trivial. I was interested in Python because it is free and more portable. I knew of Sage but wasn't sure if it was more appropriate than Python. -- http://mail.python.org/mailman/listinfo/python-list
Re: arbitrary precision linear algebra
On Mar 2, 12:22 pm, geremy condra wrote: > On Wed, Mar 2, 2011 at 10:21 AM, geremy condra wrote: > > On Wed, Mar 2, 2011 at 6:42 AM, Ben123 wrote: > >> Hello. I have a written Python program which currently uses numpy to > >> perform linear algebra operations. Specifically, I do matrix*matrix, > >> matrix*vector, numpy.linalg.inv(matrix), and linalg.eig(matrix) > >> operations. Now I am interested in allowing arbitrary precision. I > >> have tried gmpy, bigfloat, mpmath, and decimal but I have been unable > >> to easily implement any with my current program. I suspect I have to > >> change some commands but I am unsure what. > > >> My question is which of the arbitrary precision implementations will > >> most easily handle linear algebra? I don't care about speed, just ease > >> of use. Online tutorials for arbitrary precision linear algebra > >> operations would be useful. > > >> For example, it looks like mpmath can handle matrix operations > >>http://fredrik-j.blogspot.com/search?q=matrix > >> but I was unable to find a clear tutorial. The tutorials for most of > >> the arbitrary precision implementations demonstrate simple scalar > >> examples. > > >> Thanks in advance > > > Have you looked at Sage[0]? I don't know for a fact, but you should be > > able to define a matrix over RealField(precision_in_bits) and then > > take the eigenvalue of it. I don't know if it will actually produce > > the precision you need though. > > > Geremy Condra > > Apologies, forgot the links: > > http://www.sagemath.org/doc/constructions/linear_algebra.htmlhttp://www.sagemath.org/doc/reference/sage/rings/complex_field.html > > Geremy Condra I'm not sufficiently familiar with Sage, but from http://www.sagemath.org/doc/constructions/linear_algebra.html "currently Sage does not implement multiprecision numerical eigenvalues/eigenvectors" I'll ask on the Sage forums about this. In the mean time, I'm still trying to get arbitrary precision linear algebra in Python -- http://mail.python.org/mailman/listinfo/python-list
Re: arbitrary precision linear algebra
On Mar 2, 4:48 pm, Nobody wrote: > On Wed, 02 Mar 2011 06:42:22 -0800, Ben123 wrote: > > Hello. I have a written Python program which currently uses numpy to > > perform linear algebra operations. Specifically, I do matrix*matrix, > > matrix*vector, numpy.linalg.inv(matrix), and linalg.eig(matrix) > > operations. Now I am interested in allowing arbitrary precision. I > > have tried gmpy, bigfloat, mpmath, and decimal but I have been unable > > to easily implement any with my current program. I suspect I have to > > change some commands but I am unsure what. > > numpy is implemented in C, and is limited to the language's primitive > types. It provides a "float96" type which I would guess is actually a > "long double" (80 bits on x86). You can't extend numpy's precision by > using a separate arbitary-precision library. Hello. To clarify, I don't need numpy explicitly. I'm looking for an arbitrary precision library which can also do linear algebra operations: matrix*matrix, matrix*vector, matrix inverse, and eigenvalues of matrix. -- http://mail.python.org/mailman/listinfo/python-list
Re: arbitrary precision linear algebra
On Mar 2, 1:34 pm, geremy condra wrote: > On Wed, Mar 2, 2011 at 10:47 AM, Ben123 wrote: > > On Mar 2, 12:22 pm, geremy condra wrote: > >> On Wed, Mar 2, 2011 at 10:21 AM, geremy condra wrote: > >> > On Wed, Mar 2, 2011 at 6:42 AM, Ben123 wrote: > >> >> Hello. I have a written Python program which currently uses numpy to > >> >> perform linear algebra operations. Specifically, I do matrix*matrix, > >> >> matrix*vector, numpy.linalg.inv(matrix), and linalg.eig(matrix) > >> >> operations. Now I am interested in allowing arbitrary precision. I > >> >> have tried gmpy, bigfloat, mpmath, and decimal but I have been unable > >> >> to easily implement any with my current program. I suspect I have to > >> >> change some commands but I am unsure what. > > >> >> My question is which of the arbitrary precision implementations will > >> >> most easily handle linear algebra? I don't care about speed, just ease > >> >> of use. Online tutorials for arbitrary precision linear algebra > >> >> operations would be useful. > > >> >> For example, it looks like mpmath can handle matrix operations > >> >>http://fredrik-j.blogspot.com/search?q=matrix > >> >> but I was unable to find a clear tutorial. The tutorials for most of > >> >> the arbitrary precision implementations demonstrate simple scalar > >> >> examples. > > >> >> Thanks in advance > > >> > Have you looked at Sage[0]? I don't know for a fact, but you should be > >> > able to define a matrix over RealField(precision_in_bits) and then > >> > take the eigenvalue of it. I don't know if it will actually produce > >> > the precision you need though. > > >> > Geremy Condra > > >> Apologies, forgot the links: > > >>http://www.sagemath.org/doc/constructions/linear_algebra.htmlhttp://w... > > >> Geremy Condra > > > I'm not sufficiently familiar with Sage, but from > >http://www.sagemath.org/doc/constructions/linear_algebra.html > > > "currently Sage does not implement multiprecision numerical > > eigenvalues/eigenvectors" > > > I'll ask on the Sage forums about this. In the mean time, I'm still > > trying to get arbitrary precision linear algebra in Python > > I'd suggest you read that slightly more carefully. It's speaking > specifically of RR and CC, ie, double-width reals and complex values. > Using RealField and ComplexField- the arbitrary precision real and > complex fields- seems to be working. Using the earlier example: > > sage: M1 = Matrix(RealField(1000), [[0, 2], [1, 0]]) > sage: M2 = Matrix(RR, [[0, 2], [1, 0]]) > sage: M1.eigenvalues() > [1.414213562373095048801688724209698078569671875376948073176679737990732478 > 462107038850387534327641572735013846230912297024924836055850737212644121497 > 099935831413222665927505592755799950501152782060571470109559971605970274534 > 596862014728517418640889198609552329230484308714321450839762603627995251407 > 99, > -1.414213562373095048801688724209698078569671875376948073176679737990732478 > 462107038850387534327641572735013846230912297024924836055850737212644121497 > 099935831413222665927505592755799950501152782060571470109559971605970274534 > 596862014728517418640889198609552329230484308714321450839762603627995251407 > 99] > sage: M2.eigenvalues() > [1.41421356237310, -1.41421356237310] > > Converting the first of the latter values to an element of > RealField(1000) yields much what I would expect from higher precision > arithmetic: > > R = RealField(1000) > sage: x = M1.eigenvalues()[0] > sage: y = R(M2.eigenvalues()[0]) > sage: x > 1.4142135623730950488016887242096980785696718753769480731766797379907324784 > 621070388503875343276415727350138462309122970249248360558507372126441214970 > 999358314132226659275055927557999505011527820605714701095599716059702745345 > 968620147285174186408891986095523292304843087143214508397626036279952514079 9 > sage: y > 1.41421356237309514547462185873882845044136047363281250 > 000 > 000 > 000 0 > > So, while I don't know for a fact that it's using the precision you > need, it certainly does seem to be using high precision arithmetic > here. Furthermore, repeating it for various precisions seems to > increase the difference, as would be expected from better > approximations, and the number of digits in the result is consistent > with the interpretation that it is using the specified precision. > > All of this to say that it seems to be doing what you want it to do. > > Geremy Condra Hello. I was able to install Sage 4.6.1 and get your example to work. I will update this thread once I get my program to work with Sage. -- http://mail.python.org/mailman/listinfo/python-list
