Methods for polyhedron are ring dependant?

[jplab]

Hi everyone,

Currently, I am working on the ticket:

http://trac.sagemath.org/ticket/18442

which implements a barycentric subdivision of a polytope. For this I use the containment check. In the tests, it works fine for a regular pentagon over AA but not over RDF. The problem can be seen with the following code:

sage: P = polytopes.regular_polygon(5)
sage: a_vertex = P.vertices()[0]
sage: for facet in P.Hrepresentation(): print facet.contains(a_vertex), facet.interior_contains(a_vertex)

True False
True True
True False
True True
True True
sage: P = polytopes.regular_polygon(5, base_ring=RDF)
sage: a_vertex = P.vertices()[0]
sage: for facet in P.Hrepresentation(): print facet.contains(a_vertex), facet.interior_contains(a_vertex)

True True
True True
True True
True True
True True

(The first output is as expected: the vertex is not contained in the interior of exactly 2 facets. This is not the case in RDF)

Right now, I have the following question in mind: Does it make sense to have polyhedron over RDF? By which I mean: it can not even deal with containment of faces properly. Reasonably, the user knows that and will use a different ring... But what if?

I know that it is practical to have polyhedron defined over RDF for different reasons, but there is a fundamental change in the results (and pretty much any method doing computations) if the ring changes. I'm NOT asking to remove RDF as a possible ring. The reason I'm asking is that I have to use such containment check to implement the method barycentric subdivision for polyhedron, but then it breaks for polyhedron defined on RDF because of that. So there are methods that simply do not work when the ring is RDF.

Maybe this analogy could explain my thoughts:

If we compute the eigenvalues of a matrix defined over RDF, sage warns us that the results may go bad. It doesn't if the matrix is rational... Could such a warning be issued when working with RDF? And then put NotImplemented whenever it is known that RDF could cause trouble (and then fix it for that specific ring?).

I could simply set the barycentric subdivision for polyhedron over RDF not to be implemented, but I feel this would not answer the problem raised anyway.

Thanks!

[Jan Keitel]

In this particular case the problem can probably be resolved by being a bit stricter with the check in interior_contains:

Instead of checking

return self.polyhedron()._is_positive( self.eval(Vobj) )

one should have

return self.polyhedron()._is_positive( self.eval(Vobj) ) and not self.polyhedron()._is_zero( self.eval(Vobj) )

[Vincent Delecroix]

On 13/10/15 06:48, Jan Keitel wrote:
> In this particular case the problem can probably be resolved by being
> a bit stricter with the check in interior_contains:
>
> Instead of checking
>
> return self.polyhedron()._is_positive( self.eval(Vobj) )
>
> one should have
>
> return self.polyhedron()._is_positive( self.eval(Vobj) )
> and not
> self.polyhedron()._is_zero( self.eval(Vobj) )

This will not solve the problem! _is_positive(x) returns (x >= -1e-6)...
which is exactly the same thing as _is_nonneg(x)... You will also get
False positive and True negative, won't you? For example, you can be out
of the polytope but get from _is_positive that you are inside.

I would go for:
- a warning for non exact rings (using the method .is_exact())
- using _is_positive for non exact rings

Though, using intensively _is_positive or _is_nonneg might significantly
slows down the code (it is one more Python function call).

Vincent

[Jan Keitel]

Hi Vincent,

that's actually not true, have you checked it?

The code

elif Vobj.is_vertex():

            return (self.polyhedron()._is_positive( self.eval(Vobj) ) and

                    not self.polyhedron()._is_zero( self.eval(Vobj) ) )

inside interior_contains in representations.py gives

sage: P = polytopes.regular_polygon(5, base_ring=RDF)

sage: sage: a_vertex = P.vertices()[0]

sage: sage: for facet in P.Hrepresentation(): print facet.contains(a_vertex), facet.interior_contains(a_vertex) 

True True

True False

True True

True False

True True

Whether this is the correct way to go is an entirely different question - all I'm saying is that _is_positive checks (with fuzziness) whether that the argument is non-negative, but not whether it is non-zero. Maybe that is what should be changed.

[Vincent Delecroix]

No need to check, the code is explicitely as I said:


def _is_positive(self, x):
return x >= 1e-6

Here is an example:

sage: P = polytopes.regular_polygon(5, base_ring=RDF)

sage: v0,v1,v2,v3,v4 = [v.vector() for v in P.vertices()]

sage: v = (.5*v0 + .5*v1) - 1e-4*v2
sage: v in P
False

sage: v = (.5*v0 + .5*v1) - 1e-7*v2
sage: v in P
True

Everything is fine with base_ring=AA (and using exact vectors of course).

Vincent

[Jan Keitel]

Hi Vincent,

I don't know what you're checking, but the code in base_RDF.py says

return x>=-1e-6

So no, it is not checking true positivity.

[Nathann Cohen]

Hello,

Instead of checking

return self.polyhedron()._is_positive( self.eval(Vobj) )

one should have

return self.polyhedron()._is_positive( self.eval(Vobj) ) and not self.polyhedron()._is_zero( self.eval(Vobj) )

As Vincent said, this would just be a different kind of wrong. I personally don't mind if you want to turn the old wrong answers into new wrong answers, but perhaps what we should do here is protect the users against a false kind of security.

Indeed, anybody who actually plays with RDF matrices directly will soon notice that testing "==0" does not lead very far.

What would you think of adding a warning in _is_zero, only when it returns *true* ? If we then write the function your way -- by calling _is_zero -- then we would show a warning whenever some decision is taken based on the output of _is_zero ?

This way, all functions that call is_zero would show this warning (only once) on RDF, and ... well, perhaps add to it a suggestion to use an exact ring if results start looking weird ?

Nathann

[Vincent Delecroix]

Hi,

> I don't know what you're checking, but the code in base_RDF.py says

The point

(.5*v0 + .5*v1) - 1e-7*v2

is definitely not in the regular pentagon. But Sage answers that it is.

If you start taking barycenters, this is the kind of things you can run
into.

>> On 13/10/15 08:30, Jan Keitel wrote:
>>> that's actually not true, have you checked it?
>>> The code

What was actually not true in my initial statement?

Vincent

[Jan Keitel]

Hi Vincent,

you were missing the crucial minus in the code: Instead of checking that x>1e-6, it is checking whether x>-1e-6.

[Vincent Delecroix]

I did not miss it (but forgot it when I copied it in my e-mail). You
have cancellations in approximate rings:

sage: (1e60 - 1) - 1e60
0.000000000000000
sage: _.is_zero()
True

And you can check that this number is larger than -1e-6 if you wish. It
will still be wrong.

Vincent

[Volker Braun]

BTW the 1e-6 fuzzy zero is because that is used in cdd, the backend for RDF polyhedra.

You can't really distinguish positive from non-negative in the nonexact case, so I'm fine with just providing it "as is".   

[jplab]

Hi everyone,

Out of the discussion I could conclude the following:

-it would be a good idea to warn the user if some results may be wrong due to the used ring. Perhaps this deserves a ticket on its own.
-use is_positive for non exact ring

and I suggest:

-Perhaps containment should be done combinatorially in the case of an inexact ring (looking at the incidence matrix) and using the canonical ordering of vertices and facets?

Would that make sense?

I'm also fine with providing RDF polyhedra as they are, but for some methods I realized that it makes a difference in the way to implement them, which is good to know. In the case of barycentric subdivision it makes a huge difference in the computation time, the algorithmic method, and the truth of the results...

JP

[Nathann Cohen]

Hello,

-it would be a good idea to warn the user if some results may be wrong due to the used ring. Perhaps this deserves a ticket on its own.

+1

 

-use is_positive for non exact ring

I do not understand the aim of that.

-Perhaps containment should be done combinatorially in the case of an inexact ring (looking at the incidence matrix) and using the canonical ordering of vertices and facets?

It "may" be better in terms of speed. In terms of correctness, however, it is probably equivalent.

Nathann

[Volker Braun]

On Wednesday, October 14, 2015 at 8:27:44 AM UTC+2, jplab wrote:

-it would be a good idea to warn the user if some results may be wrong due to the used ring. Perhaps this deserves a ticket on its own.

I think thats a bit silly; really it is just a special case of "any result in RDF might me mathematically wrong due to fixed precision". That is true but such a fundamental fact that I don't expect it to be printed on the screen in every session where I use floating point numbers. 

[David Roe]

Agreed.  There are already some notes in the documentation for such rings warning of this kind of issue, but it might be worth expanding upon them.

David

--
You received this message because you are subscribed to the Google Groups "sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.
To post to this group, send email to sage-...@googlegroups.com.
Visit this group at http://groups.google.com/group/sage-devel.
For more options, visit https://groups.google.com/d/optout.