Distinguish inexact rings

[wbhart]

Currently inexact and exact rings are not properly distinguished in our system. This makes it hard to write code that works for one kind of ring and not the other.

Also, we currently have Singular.Ring <: Ring and Singular.Field <: Field in Singular.jl, but not Singular.Field <: Singular.Ring.

These problems occur due to the lack of multiple inheritance in Julia. We formally require traits to work around this.

What I propose is that we parameterise the Ring and Field abstract types in AbstractAlgebra, Ring{M, E} and Field{M, E} where M is a symbol :generic, :nemo, :singular and where E is true if the ring/field is exact. Naturally we would still have Field{M, E} <: Ring{M, E}.

Of course RingElem and FieldElem would have to be similarly parameterised.

Each ring type would be required to belong to some specific Ring{M, E} or Field{M, E} and similarly for each ring/field element type.

Of course a polynomial ring over an inexact ring/field becomes itself an inexact ring which necessarily implies that all Ring/RingElem and Field/FieldElem concrete types in all our packages must similarly be parameterised, e.g.

mutable struct Poly{E, T <: RingElem{M, E} where M} <: RingElement{:generic, E} end

The only fly in the oitment currently is that types like Poly are parameterised on T <: RingElement which includes Julia Int and BigInt and so on, which don't have the M, E parameters. The change I am proposing would prevent us from using such types. There's no way to make Int belong to RingElem{M, E}

I don't have a workaround for this at present, but it's probably good that we see this problem now rather than later.

Of course making this change is rather major and would require a lot of effort over multiple packages. But it seems necessary in the limit.

[wbhart]

One possible way to handle the RingElement problem would be to define new types JuliaInt <: RingElement{:julia, true}, etc. and define all kinds of arithmetic for them with real Julia Int's so they for the most part behave just like system Ints.

Only generic objects like Poly and Mat etc. would actually use them and even then only internally. We could set up the polynomial module to return an actual system Int when one calls coeff on a Poly{<:JuliaInt}.

[wbhart]

The other issue that would need discussion is the extent to which Singular is expected to handle Singular rings as coefficient rings in future. If not then we probably need to distinguish not just RingElem{:singular} and RingElem{:nemo} etc, but also RingElem{:singular_coeff} and RingElem{:singular}, which could be problematic if we ever want one to belong to the other.

[thofma]

Do I understand correctly that this is a way to do traits without doing traits?

One question I have is: Does this scale? If we put exact/inexact in there, what other properties should be treated like this? Euclidean? Rings without coefficient explosion? I am not sure I understand the reason to push the exact/inexact into the type but not doing it for other properties. Just a thought a had after reading your proposal.

[wbhart]

Yeah this is just traits without explicit support for them. The only reason to put exact/inexact and module as traits and not anything else is that they are the main things holding us back as a project, I believe.

I don't really believe Euclidean can be a type level trait. I don't know the extent to which rings with/without coefficient explosion is holding us back. But if you have hit it as a severe limitation of the system many times, then we should absolute discuss including it. It's not easy to determine though is it? I mean, for example, is a residue ring over a ring without coefficient explosion a ring with or without coefficient explosion? I'm not entirely sure if it is or isn't really suitable as a type level trait.

At present, I don't believe there is any really neat way to handle arbitrary traits for all our types, so this is a proposal that does as little as it has to, in order to give us a substantial benefit. We could add an additional trait in future if it proves absolutely necessary, but after 6 years we should have a good handle on the ones that always seem to be killing us.

[thofma]

OK. (I was not aware that the lack of distinguishing exact and inexact is one of the two main things holding us back as a project.)

I guess this will work, but will make writing function signatures a bit painful. Maybe something like const ExactRingElem{T} = RingElem{T, :exact} would work?

What will be our definition of an exact/inexact field?

[wbhart]

I think the exact/inexact thing is holding us back at this point. When I said "main things holding us back", I meant "main things that can be addressed by this proposal", not main things unqualified, sorry. And even then that's only my opinion of course, which is why I made a proposal for discussion.

The reason exact/inexact is so important is that different algorithms are needed for inexact rings/fields. The main thing I have in mind here is linear algebra and polynomial arithmetic over power series rings, both of which have come up in recent discussions with @fieker

At present we have no way to write linear algebra or polynomial routines that work differently for such rings. Typically one has to distinguish exact zero from approximate zero, do precision handling, manage "numerical" stability and a host of other things we can't currently do. The other issue is not being able to raise an exception when someone calls a function for a ring where a result is returned, but is mathematically wrong, due to the inexact nature of the ring.

I agree that some helpers like the ExactRingElem{T} you mentioned might ease the burden. It's a two edged sword. The more you try to hide from the developer, the more likely they get confused. But this particular suggestion seems ok to me.

As for definition of exact/inexact, that is a good question. I actually wonder if we will want to distinguish numerical inexact from other kinds of inexact. For example, one probably needs to handle BigFloats and arbs differently from PowerSeriesRing's.

But as a first step, I would define inexact as any ring in the system that stores an approximation to elements of a mathematical ring and maintains a precision. If we later decide that we want to write generics for multiple different kinds of inexact (I think we will in the limit, but not in the near future), then we can simply introduce additional values in the trait. Perhaps for this reason we should use :exact rather than true/false, or use a type called Exact or something like that. This way we can add more values later.

Of course we'd need to use Union's to write functions for multiple kinds of inexact if we ever did introduce additional values of that trait, but I don't see a problem with that.

I'm keen to not overengineer the proposal too. Unless there really are other really urgent problems that can actually be solved on the type level, it's probably best to just do the simplest thing to fix the problems we know about. I can't think of any other things that could be fixed by type level traits in this way that have repeatedly hit us over the years. But it is worth thinking about whether there are any others.

[wbhart]

By the way, I think very few actual functions will need to distinguish exact and inexact rings. The ones I can think of include those that do iszero checking, e.g. normalisation of polynomials, pivoting in linear algebra, printing of objects over series base rings, Euclidean division. There aren't that many. So although this proposal would be a massive imposition on our type system, I don't think it will interfere much with routine development.

Obviously anyone writing routines for numerical stable linear algebra and polynomial arithmetic will need this. But we don't have any of this right now, so it's not as if that is interfering with any of the things we have currently been working on.

[wbhart]

The proposal to define JuliaInt <: RingElement{:julia, :exact} and so on may actually have additional benefits. One is that we can define arithmetic for these however we want. So I think we get rid of the currently ugly hacks to handle exp, log, inv, sqrt, mod, etc internally differently than externally.

The second advantage, I think, is that we can get rid of the RingElement hack. We'd just need RingElem.

[fredrik-johansson]

To accommodate inexact rings, predicates (iszero, ==, etc.) really need to come in different, explicitly named versions.

• Representation-level: for example, equal intervals compare as equal, although they may represent different numbers. Mathematically meaningless, but useful, for example, if you want to cache inexact elements.
• One-directional: True implies mathematical truth but False is meaningless. This is how most Arb predicates work. It's the least invasive but requires some discipline to use correctly (all logical operations and if-statements must be designed around this), which in practice means generic code WILL be broken for inexact rings.
• Throwing: True or False imply certainty; uncertainty results in an exception. This has the advantage that code written assuming exact rings will be correct automatically for inexact rings. There will be lots of annoying exceptions though; generic code will need to be fixed to remove spurious exceptions.
• Triple-valued: returning True/False/Unknown (i.e. evaluating into an inexact boolean ring!). Also requires coding discipline, but may be more robust and lead to cleaner code in the long run.

I don't really see a way around explicitly and systematically distinguishing between these notions, which I'm doing in Calcium. Well, there is an alternative, and that is to have a Sage-like or Mathematica-like system where 0 == 1 and pigs fly.

[fieker]

To add more: you need, due to one-directionality possible iszero isnonzero In Magma there is also the even less well defined IsWeaklyZero Which is mainly useful, but not well defined...

On Thu, May 06, 2021 at 03:41:22AM -0700, Fredrik Johansson wrote:

To accommodate inexact rings, predicates (iszero, ==, etc.) really need to come in different, explicitly named versions.

• Representation-level: for example, equal intervals compare as equal, although they may represent different numbers. Mathematically meaningless, but useful, for example, if you want to cache inexact elements.
• One-directional: True implies mathematical truth but False is meaningless. This is how most Arb predicates work. It's the least invasive but requires some discipline to use correctly (all logical operations and if-statements must be designed around this), which in practice means generic code WILL be broken for inexact rings.
• Throwing: True or False imply certainty; uncertainty results in an exception. This has the advantage that code written assuming exact rings will be correct automatically for inexact rings. There will be lots of annoying exceptions though; generic code will need to be fixed to remove spurious exceptions.
• Triple-valued: returning True/False/Unknown (i.e. evaluating into an inexact boolean ring!). Also requires coding discipline, but may be more robust and lead to cleaner code in the long run.

I don't really see a way around explicitly and systematically distinguishing between these notions, which I'm doing in Calcium. Well, there is an alternative, and that is to have a Sage-like or Mathematica-like system where 0 == 1 and pigs fly.

-- You are receiving this because you were mentioned. Reply to this email directly or view it on GitHub: https://github.com/Nemocas/AbstractAlgebra.jl/issues/860#issuecomment-833423263

[wbhart]

Yes, I absolutely believe you are right @fredrik-johansson

I actually have no idea how to handle @fieker 's proposal for linear algebra over an inexact Euclidean domain (K[[x]] for a field K, taken to some undefined precision, to be precise; excuse the pun). Perhaps you have some ideas.

Obviously pivoting depends on detecting zero elements. But there's no way to know if you have a mathematical zero in K[[x]] at any precision. The big problem I see is that exact division in any power series ring results in lower precision. Eventually you are asking whether O(x^0) is zero. For that reason you may need to increase the precision to get a meaningful answer. But then you hit the precision cap.

Similarly, @fieker wants to work in a ring where he does not know a precision cap in advance. He wants to be able to lift to some undefined precision, which effectively means removing the precision cap on our AbsSeries rings R. But this then leaves values like zero(R) and one(R) that were previously used, now at the wrong precision (they are taken at the precision cap). And without a precision cap to be used as a default precision, one cannot currently even define zero(R) or one(R) or R(2) for that matter.

We could add exact elements, I guess, but then as you know, the can is just kicked down the road to divexact of exact elements.

The main motivation for this particular proposal here is simply so that we can define an inexact Euclidean ring K[[x]] and then annotate generic functions, which definitely don't work with such a ring, to only accept exact inputs. This prevents the ugly situation of users getting an answer, believing it to be correct and then publishing papers based on a faulty assumption about what our system can actually do.

If you have any suggestions, please do chip in @fredrik-johansson

I'm really only happy to do it at all if it is done properly. But currently, I don't see what that even is. This proposal is hopefully one small step in that direction.

[fieker]

On Thu, May 06, 2021 at 03:59:54AM -0700, wbhart wrote:

Yes, I absolutely believe you are right @fredrik-johansson

I actually have no idea how to handle @fieker 's proposal for linear algebra over an inexact Euclidean domain K[[x]] to be precise. Perhaps you have some ideas.

Obviously pivoting depends on detecting zero elements. But there's no Nope, not really. We're oding euclidean operations (think HNF). For stability you'd select elements of minimal valuation. Furthermore, HNF implies only unimodular operations, so this is even mostly stable, in case of row and column operations (SNF) the result is best possible with no unneccessary loss of anything. Avi (and Tristan and Caruso(?) did this over the p-adics, but the theory is the same.

This is totally different fromt he real case.

way to know if you have a mathematical zero in K[[x]] at any precision. The big problem I see is that exact division in any power series ring results in lower precision. Eventually you are asking whether O(x^0) is zero. For that reason you may need to increase the precision to get a meaningful answer. But then you hit the precision cap.

Similarly, @fieker wants to work in a ring where he does not know a precision cap in advance. He wants to be able to lift to some undefined precision, which effectively means removing the precision cap on our AbsSeries rings R. But this then leaves values like zero(R) and one(R) that were previously used at the wrong precision (they are taken at the precision cap). I basically want to sanction the behaviour of the underlying implementation: the precision cap is rarely ever used, binary operations look at the elements, not the ring. In applications like the Galois groups you cannot give a (meaningful) upper bound on the precision. And even if you can, you only want to work with the minimal precision for this step. This is exactly the same as anu Newton/Hensel for reals: do as little as possible and increase...

We could add exact elements, I guess, but then as you know, the can is just kicked down the road to divexact(exact_zero(R), f) for some inexact element f. I am not in favour of exact elements, they produce unpredictable runtime.

The main motivation for this particular proposal here is simply so that we can define an inexact Euclidean ring K[[x]] and then annotate generic functions which definitely don't work with such a ring to only accept exact inputs. This prevents the ugly situation of users getting an answer, believing it to be correct and the publishing papers based on a faulty assumption about what our system can actually do. ... like in numerics? But as I said: here, using HNF techniques, the rref is loss-less. So unless the input precision is too low (in which case you get false zeros), the result is as correct as the input.

If you have any suggestions, please do chip in @fredrik-johansson

I'm really only happy to do it at all if it is done properly. But currently, I don't see what that even is. This proposal is hopefully one small step in that direction.

-- You are receiving this because you were mentioned. Reply to this email directly or view it on GitHub: https://github.com/Nemocas/AbstractAlgebra.jl/issues/860#issuecomment-833433548

[wbhart]

On Thu, May 06, 2021 at 03:59:54AM -0700, wbhart wrote: Yes, I absolutely believe you are right @fredrik-johansson I actually have no idea how to handle @fieker 's proposal for linear algebra over an inexact Euclidean domain K[[x]] to be precise. Perhaps you have some ideas. Obviously pivoting depends on detecting zero elements. But there's no Nope, not really. We're oding euclidean operations (think HNF). For stability you'd select elements of minimal valuation.

Right, but this is just another way of saying, "least zero" element. And my original objection still applies. O(x^0) is always the element with least valuation, but mathematically it may not have actually been the element with least valuation in K[[x]].

Furthermore, HNF implies only unimodular operations, so this is even mostly stable, in case of row and column operations (SNF) the result is best possible with no unneccessary loss of anything. Avi (and Tristan and Caruso(?) did this over the p-adics, but the theory is the same. This is totally different fromt he real case.

Here is my view on this. Suppose someone asks solve(A, b) for a matrix A and vector b over your ring. Suppose that every element of A and b are given to the same precision. What can you prove about the result? Ideally you'd like that the result is given to the same precision, if that is even possible.

So what have your students actually proved? What was their model of power series?

I totally get that you can look at any specific mathematical binary operation (other than division) and say that if both elements are given to precision n then the result can be given to precision n. But how is division handled? How are zero and one handled? What happens if you need to raise the precision at some point in the computation?

As a former mathematician, I am trained to be very skeptical of people's claims about their computer algebra systems. I am particularly skeptical of claims such as "so and so is doing blah in system blah".

In addition, I don't know how general this can be made. Should I only accept the answers they gave over the one specific ring they worked with, or is this a generally applicable technique?

Anyhow, your statement that this is different to numerical calculation is at least suggesting that maybe we need :exact, :adic and :numeric instead of just bistate :exact and :inexact in the current proposal.

way to know if you have a mathematical zero in K[[x]] at any precision. The big problem I see is that exact division in any power series ring results in lower precision. Eventually you are asking whether O(x^0) is zero. For that reason you may need to increase the precision to get a meaningful answer. But then you hit the precision cap. Similarly, @fieker wants to work in a ring where he does not know a precision cap in advance. He wants to be able to lift to some undefined precision, which effectively means removing the precision cap on our AbsSeries rings R. But this then leaves values like zero(R) and one(R) that were previously used at the wrong precision (they are taken at the precision cap). I basically want to sanction the behaviour of the underlying implementation: the precision cap is rarely ever used, binary operations look at the elements, not the ring. In applications like the Galois groups you cannot give a (meaningful) upper bound on the precision. And even if you can, you only want to work with the minimal precision for this step. This is exactly the same as anu Newton/Hensel for reals: do as little as possible and increase...

Not a big deal, but I don't understand the word "sanction" in this context, so I am probably missing your point.

I agree the precision cap is rarely used. But our generic code would currently, definitely use it. At the very least we agree that our current generic code will not Just Work TM.

Regarding raising the precision: for linear algebra for example, who would raise the precision here? Would it be the person who called the linear algebra step, or could it happen in the linear algebra step itself?

In other words, who is responsible for determining if the output of HNF over an inexact ring is mathematically meaningful/to sufficient precision? Is it the person who called the routine, or the responsibility of the routine itself?

What exactly would the routine have to guarantee mathematically?

We could add exact elements, I guess, but then as you know, the can is just kicked down the road to divexact(exact_zero(R), f) for some inexact element f. I am not in favour of exact elements, they produce unpredictable runtime.

I didn't like them much, but @fredrik-johansson uses them in various contexts I believe. He still likes them I think.

The main motivation for this particular proposal here is simply so that we can define an inexact Euclidean ring K[[x]] and then annotate generic functions which definitely don't work with such a ring to only accept exact inputs. This prevents the ugly situation of users getting an answer, believing it to be correct and the publishing papers based on a faulty assumption about what our system can actually do. ... like in numerics?

I wish this were confined to numerics only. Plenty of existing systems place the burden on the user for checking if the output of calls to their functions are meaningful. That is still the case with us, even in AbstractAlgebra, as perfect and without bugs, error or design inconsistencies as it is.

But as I said: here, using HNF techniques, the rref is loss-less. So unless the input precision is too low (in which case you get false zeros), the result is as correct as the input.

Because it only uses ring operations? There are plenty of divexact calls in our current implementation, which makes me a little skeptical.

[wbhart]

What I am really asking at the point I talk about "mathematical proofs" above is: what identity is true and to which precision, after calling rref, or hnf. And how is equality defined in this identity?

Regarding raising the precision cap, the only place I can see this as being useful is where you run an rref with entries at some precision, decide you don't like the answer, then completely recompute the input from scratch at a higher precision (from elements of infinite precision coming from your mathematical problem) and re-run the computation. It seems that creating a new ring at the new precision is not really problematic in this case.

If there are examples where raising the precision cap for some other reason is actually desirable, it might be useful to have those.

[fieker]

On Thu, May 06, 2021 at 04:55:38AM -0700, wbhart wrote:

On Thu, May 06, 2021 at 03:59:54AM -0700, wbhart wrote: Yes, I absolutely believe you are right @fredrik-johansson I actually have no idea how to handle @fieker 's proposal for linear algebra over an inexact Euclidean domain K[[x]] to be precise. Perhaps you have some ideas. Obviously pivoting depends on detecting zero elements. But there's no Nope, not really. We're oding euclidean operations (think HNF). For stability you'd select elements of minimal valuation.

Right, but this is just another way of saying, "least zero" element. And my original objection still applies. O(x^0) is always the element with least valuation, but mathematically it may not have actually been the element with least valuation in K[[x]].

Furthermore, HNF implies only unimodular operations, so this is even mostly stable, in case of row and column operations (SNF) the result is best possible with no unneccessary loss of anything. Avi (and Tristan and Caruso(?) did this over the p-adics, but the theory is the same. This is totally different fromt he real case.

Here is my view on this. Suppose someone asks solve(A, b) for a matrix A and vector b over your ring. Suppose that every element of A and b are given to the same precision. What can you prove about the result? Ideally you'd like that the result is given to the same precision, if that is even possible.

So what have your students actually proved? What was their model of power series?

Not my students, that is Xavier Caruso (Sage) and Tristan Vaccon(?) The claim is there are unimodular transformations U and V s.th. UAV = diag(a_i), val(a_1) \le ... \le val(a_n) U and V unimodular means, in particular, that UAV has exactly the same precision as A and that inv(U) and inv(V) also have the same precision and are unimodular.

Thus the unavoidable precision loss in the solution is precisely val(a_n).

We had a talk about this in Lautern while Avi was here.

I totally get that you can look at any specific mathematical binary operation (other than division) and say that if both elements are given to precision n then the result can be given to precision n. But how is division handled? How are zero and one handled? What happens if you need to raise the precision at some point in the computation? If I change the precision in the middle of a single rref call, I deserve what I get. However, if a perform a Newton lifting this makes sense: I increase the precision used to map exact elements in and mathematically know that the next iteration will fill the additional power series coefficients with the correct values - regardless of the previous ones. (meaning if I do add non-zero additional digits I still get the same result)

As a former mathematician, I am trained to be very skeptical of people's claims about their computer algebra systems. I am particularly skeptical of claims such as "so and so is doing blah in system blah". I'd rather to not ahve this pointless remarks in a public discussion. It is really not my fault that imprecise rings exist and are useful and offend many computer algebra people.

In addition, I don't know how general this can be made. Should I only accept the answers they gave over the one specific ring they worked with, or is this a generally applicable technique? The theory as far as my memory works is valid for all rings with a non-archimedian valuation, since they are euclidean with the property that RING operations (so no generic division) do not increase errors/ loose information. The definition should include power series over an (exact) field. I am confident that a lot also works for power series over a non-achimedia ring, but I am not claiming that.

Anyhow, your statement that this is different to numerical calculation is at least suggesting that maybe we need :exact, :adic and :numeric instead of just bistate :exact and :inexact in the current proposal.

If we want to branch on type: yes.

I agree the precision cap is rarely used. But our generic code would currently, definitely use it. At the very least we agree that our current generic code will not Just Work TM.

Our current generic code does not change the precision so it will continue to work. What would happen is that M some matrix (fixed with some precision) A = M^10 change precision on the ring B = M^10 might have A ne B I do not know if this is a problem

Regarding raising the precision: for linear algebra for example, who would raise the precision here? Would it be the person who called the linear algebra step, or could it happen in the linear algebra step itself? As usual: if you change the input, same as for floats, you have to reset it. Of course, Dixon lifting (or similar) might use this internally - or not.

In other words, who is responsible for determining if the output of HNF over an inexact ring is mathematically meaningful/to sufficient precision? Is it the person who called the routine, or the responsibility of the routine itself? I cannot change mathematics. If you hand me 2 elements, one with valuation, then the division will kill absolute precision. This is unavoidable. However as the author of the solve function I can guarantee to have the precision loss at the matheamtically unavoidable minimum.

If this is enough for the user: cannot say, don't know their intend. Howver, what I did was best possible, Why is this different from the current situation wrt. element arithmetic? Binary operations to the best they can - to the mathematical limit. If this is enough: depends on your application.

What exactly would the routine have to guarantee mathematically? I cannot speak for all possible functions, however I can design and write a linear algebra solver that, for a given input, will produce the best possible output, ie a solution with the largest precision determined by the given input. I can also get weakter versions of this that are faster (e.g assuming the input has uniform precision or so)

We could add exact elements, I guess, but then as you know, the can is just kicked down the road to divexact(exact_zero(R), f) for some inexact element f. I am not in favour of exact elements, they produce unpredictable runtime.

I didn't like them much, but @fredrik-johansson uses them in various contexts I believe. He still likes them I think.

They have their uses... and I've used them before, fell into traps, did not use them, fell into other traps...

The main motivation for this particular proposal here is simply so that we can define an inexact Euclidean ring K[[x]] and then annotate generic functions which definitely don't work with such a ring to only accept exact inputs. This prevents the ugly situation of users getting an answer, believing it to be correct and the publishing papers based on a faulty assumption about what our system can actually do. ... like in numerics?

I wish this were confined to numerics only. Plenty of existing systems place the burden on the user for checking if the output of calls to their functions are meaningful. That is still the case with us, even in AbstractAlgebra, as perfect and without bugs, error or design inconsistencies as it is.

But as I said: here, using HNF techniques, the rref is loss-less. So unless the input precision is too low (in which case you get false zeros), the result is as correct as the input.

Because it only uses ring operations? There are plenty of divexact calls in our current implementation, which makes me a little skeptical.

-- You are receiving this because you were mentioned. Reply to this email directly or view it on GitHub: https://github.com/Nemocas/AbstractAlgebra.jl/issues/860#issuecomment-833463859

[fieker]

On Thu, May 06, 2021 at 05:15:55AM -0700, wbhart wrote:

What I am really asking at the point I talk about "mathematical proofs" above is: what identity is true and to which precision, after calling rref, or hnf. And how is equality defined in this identity?

Regarding raising the precision cap, the only place I can see this as being useful is where you run an rref with entries at some precision, decide you don't like the answer, then completely recompute the input from scratch at a higher precision (from elements of infinite precision coming from your mathematical problem) and re-run the computation. It seems that creating a new ring at the new precision is not really problematic in this case.

If there are examples where raising the precision cap for some other reason is actually desirable, it might be useful to have those.

Newton lifting: given an exact polynomial f in Q(x)[y] I want the power-series expansion (in x) of the roots of f. I do not know in advance how many coefficients I need. (I cannot know this, not my lazyness)

In the initial steps of the lifting I'd expect the precision to be smaller than the degree of the coefficients of f, so truncation will occur to make things fast Later, the coefficients will probably be exact. Unless they are rational....

-- You are receiving this because you were mentioned. Reply to this email directly or view it on GitHub: https://github.com/Nemocas/AbstractAlgebra.jl/issues/860#issuecomment-833475306

[wbhart]

Do Caruso and Vaccon give an algorithm? Or is it just an existence claim?

And if they have an algorithm that works, under what model can they guarantee this?

As I said, I am a former mathematician who is naturally very skeptical about people's claims about what they are doing. In pure mathematical circles this skepticism is seen as a virtue!

To my current level of understanding, a function that "does the best it can" is only heuristic. One has to not only state the precision to which the result is given, but what is given to that precision, i.e. what identity. One cannot then transform that identity algebraically and expect that it still holds to the same precision. This approach has to be carried through the entire computation, otherwise the entire result is heuristic.

That has its place, but lots of people are allergic to it, as you already noted.

[wbhart]

Also, regarding your question about what is different compared with other kinds of rings: I would raise a not-implemented error for implementations over those rings unless and until we can prove something specific there as well. And we should of course preferentially be using the stuff @fredrik-johansson has implemented for numerics.

As things currently are in our system, users get the idea that they get meaningful answers when they don't, precisely because we have had no way to tell them that they don't. This is exactly what my proposal on this GitHub issue is all about!

We also have similar issues with non-integral domains in the system, e.g. the isunit issue for polynomials over Z/nZ. I think we also have irreducibility tests and factorisation routines that give answers over these rings which are completely meaningless. It's unfortunate that the non-integral domain issue can't be solved at the type level.

[rfourquet]

Maybe I'm missing something obvious, but probably worth being explicit about it:

What I propose is that we parameterise the Ring and Field abstract types in AbstractAlgebra, Ring{M, E} and Field{M, E} where M is a symbol :generic, :nemo, :singular and where E is true if the ring/field is exact. Naturally we would still have Field{M, E} <: Ring{M, E}.

So although this proposal would be a massive imposition on our type system, I don't think it will interfere much with routine development.

So why not using usual function ("holy") traits? Instead of imposing a type parameter, you impose an isexact(x) trait, which in many cases could be resolved at compile time. To be less intrusive, there could also be an isexact(x::Whatever) = true default, and why not an isexact(x) = let b = base_ring(x); b === Union{} ? true : isexact(b) or something like that. And in routines needing this information, instead of dispatching on the type with an E parameter, you can either have an if: routine(x) = isexact(x) ? foo(x) : bar(x) or if you prefer routine(x) = routine(x, isexact(x)) and then define two two-arguments methods. Also, a (small) problem with requiring E as a type parameter as the interface is that for Ring{M,E} you know the parameter is at the second position, but for a general subtype, you don't really know what are its type parameter and in what order they are (unless you impose all subtypes to have E as there second paremeter, which is clearly annoying, in particular for types which don't even have one type parameter in the first place). So anyway there wouldn't seem to be a way around having to define a trait like isexact(x::Ring{M,E}) where {M,E} = E and have to use the isexact trait while defining routine...

[wbhart]

Because I can't inherit from them. Consider the Singular system, where we want SingularRing <: Ring, SingularField <: Field and SingularField <: SingularRing.

The same thing happens with isexact in the end.

My point I guess is that it is much more complicated in the Singular.jl setting, and since we should implement this solution there, why not use it for the other situation that has been a major problem since the beginning that can definitely be decided at compiler time.

[wbhart]

Here's a more specific rationale.

Suppose you wish to write a generic function in Singular.jl (or Oscar.jl) which only works for two Singular RingElem's. How could you distinguish the Singular version from the AbstractAlgebra version of the function if they both accepted RingElem's and only the Singular one checked issingular on their types subsequently? This would result in an ambiguity.

I can't think of an easy example of this sort of thing in the case of isexact, admittedly. And maybe it does impose enough of a burden on the type system to do that functionally.

[wbhart]

I'm actually not sure we will want the M and E type parameters in the first two positions. I think we probably will want them always as the last two parameters. That way if you don't care about M and E you can just write Poly{T} or RingElem{T} in your function signature, which means the least imposition on the average developer.

But you are right that when people use these type parameters, they will need to get the order correct. I think tricks like @thofma's will make this easier though. And I really don't think we want lots of type level traits. I'd stick with just the two if possible.

[wbhart]

Also (and I am not 100% sure I'm right about this), the only types that will need the extra parameters are generic ones and Singular.jl ones. All other types will just declare themselves to belong to RingElem{M, E} for some appropriate values of M and E.

[wbhart]

I guess we should just try to solve the problems we have in the simplest possible way with the least disruption. So I guess we should just implement an isexact method for all types and use "if" statements to raise not implemented errors for functions that don't support inexact arithmetic.

So I now propose:

isexact(R) - return true if the ring is exact isnumeric(R) - return true if the ring uses floating point (e.g. BigFloat), though not sure if Arb would be considered numeric

Eventually we should introduce formal interfaces for rings like Arbs, padic rings, etc.

I don't know how to solve the Singular issue without an extra type parameter. I guess we just put up with things how they are for now.

[a-kulkarn]

Mind if I jump in on the discussion?

Here is my view on this. Suppose someone asks solve(A, b) for a matrix A and vector b over your ring. Suppose that every element of A and b are given to the same precision. What can you prove about the result? Ideally you'd like that the result is given to the same precision, if that is even possible.

The short answer to this is, you know the output up to the condition number of the matrix A. If your ring admits a notion of "singular value", the condition number equals largest singular value/smallest singular value.

The theory as far as my memory works is valid for all rings with a non-archimedian valuation, since they are euclidean with the property that RING operations (so no generic division) do not increase errors/ loose information. The definition should include power series over an (exact) field. I am confident that a lot also works for power series over a non-achimedean ring, but I am not claiming that.

Certainly for anything with a uniformizer. If one replaces "divisible by a power of the uniformizer" with "lives in a power of the maximal ideal", one might run into problems like how to compute an SVD of a matrix like

[p, 0]
[t, x] 

over the ring Zp[[t,x]], with maximal ideal <p, t, x>.

Why is having an "inexact" trait a better solution than writing specialized methods for the specific cases of interest (Power series, Laurent series, extensions of Qp, and flavors of the Real and Complex fields)? Most of the time it seems as though the algorithm in each case is usually different enough to merit a separate implementation. (Not to mention the semantics of what the results of the function actually mean.) Is there a specific example of where the trait architecture solves a problem specialization can't?