Bill Allombert on Wed, 13 May 2009 14:55:42 +0200 |
[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]
Re: Another problem with matrix inversion |
On Wed, May 13, 2009 at 08:50:53AM +0200, Lorenz Minder wrote: > Hi, > > > PARI only know perform polynomial and linear algebra over a field, so > > it assumes that moduli are prime. > > Yes, of course. My point is that it would be better if it worked also > if the base-ring is not a field. Since even in that case, this is a > meaningful question, a useful answer would be preferable, in my opinion. There is the major issue that module over a ring do not generally have a basis. > I did a sketchy implementation of what I outlined in the other mail, and > a patch is attached. This patch only works for small integers right now, > as I have only modified Flm_gauss(), but not FpM_gauss(). It's not > yet production code, only a prototype for discussion. Please do not change Flm and FpM function: they are meant to be restricted to prime modulus. You will need to create new functions ZnM_xxx to handle matrices over Z/nZ. > It seems to work fine for what I need it, and so I'd like to > ask if a more complete and tested patch among the same lines would be > considered for inclusion in PARI. I trust that you can write a PARI function ZnM_gauss that would work correctly. The issue is how do you plan to integrate that with GP and keep some consistency ? i.e. it is always annoying when a feature (support for rings) is only supported for a particular ring (Z/nZ) and for a particular function. And furthermore what application do you have in mind ? > Baring any mistakes I may have made, the code should work for any > modulus, and claim that a matrix is singular only if it is indeed > singular. > > As is, the code in the patch reduces to the original algorithm if p is > a prime. For composite moduli, O(log(m)) copies of the matrix are kept > on the stack in the worst case. > > Incidentially, this patch seems to fix a nasty bug with Flm_gauss() when > the right hand side is a t_COL instead of a t_MAT; but I could not find > any caller using it that way in any case. (And I _might_ be > misunderstanding the old code, but I think that is unlikely.) A function named FpM_gauss is not supposed to handle arguments that are not t_MAT. > I have now also noticed that there is matsolvemod(), which seems to work > fine, too. So alternatively it would be possible to implement an > algorithm on top of modsolvemod(); would that be better or does it have > other drawbacks such as slower performance? (Is it possible to give a > bunch of right hand sides to matsolvemod() at once?) I would not know. Cheers, Bill.