jga's profile - activity

Thank you, @Sébastien. I think what we're seeing is just a more contrived version of this: https://en.wikipedia.org/wiki/Wilkins... . Your latest experiments also point in this direction. I also did some experiments on my side, please see the original post for the code and a brief summary. I'm also accepting this answer because it provided a solid explanation as to why every method fails. Specially the last comment showing the instability.

Thanks for your reply! The polynomials are not irreducible and contain the 6 vertical lines x=19,20,21,22,23,24. I've added two lines in the original code to get rid of them. From what you wrote it seems like the issue goes all the way down to the "roots" algorithm. I gather this method relies on "Brent's method", which requires some kind of regularity from the function. Presumably these polynomials don't satisfy the required hypotheses. I will try to see if manually implementing some other root finding algorithms could help.

@Sébastien exactly. These three polynomials should define nice graphs passing through said points.This must be the case because for each fixed value of x there can be at most 31 values of y in the plot. These vary continuously as x does, so there's just no way it looks as the plot suggests. In my case this is very clear for x<10, where the plot looks well. The issue happens further away from the origin.

I added the code, @tmonteil

I have three huge degree 31 bivariate polynomials (20,000 characters long each) I want to plot, but I keep getting a lot of noise in my plot. I can't upload it, but the point is that in some regions I just get colorful noise. I've tried defining the polynomials over RealField(n) and increasing the number of plot_points, but neither of these approaches work. Any ideas on how to work around this? Thanks.

Edit: Tried using sympy's plot_implicit and it's so (SO!) slow. Then used numpy's contour_plot and it's fast, but has the same problem as sage.

Here's the code that produces the polynomials and plot. Be patient as it could be a bit slow (depending on your machine).

plot = Graphics()
m = 32-1
D = [(i,j) for i in range(0,m+1) for j in range(0,m+1) if i+j<m+1]

#Polygon
P = Polyhedron( vertices= [(0, 0), (32, 44), (23, 0), (10, 14), (2, 3)] )
points = P.integral_points()

plot_pts = point(points, rgbcolor=(0, 0, 0), size = 20).plot()
plot_np = P.plot(fill = False, point=False, line='black')

M = matrix(ZZ, len(points), len(D), 0)

for row_num, row in enumerate(points):
    for col_num, column in enumerate(D):
        i, j = row
        a, b= column
        #Matrix for interpolation:
        M[row_num, col_num] = (i^a)*(j^b) 

R = PolynomialRing(QQ, 2, 'xy') 
S = PolynomialRing(RealField(500), 2, 'uv') 
x, y = R.gens()       
u, v = S.gens() 
K = M.right_kernel() 
Kdim = K.dimension() 
print(Kdim) 
if Kdim > 0: 
     for l in range(Kdim): 
        K_basis = K.basis()[l] 
        #Writing the interpolating polynomial 
        f=0 
        for order, bidegree in enumerate(D): 
            a, b = bidegree 
            f += list(K_basis)[order] * u^a * v^b 
        F = f.factor()
        f = F[6][0]

        cols = ['red', 'blue', 'green'] 
        interpolation = implicit_plot(f, (u,-1,34), (v,-1,45), plot_points=100, color=cols[l]) 
        plot += interpolation 
plot += plot_pts + plot_np 
plot.show(figsize=10)

Edit 2: Using the mpmath library in Python and with the aid of Sébastien's code below I wrote a routine that allows us to control the root finding method and precision of our computations. I tried several methods, secant (default), newton, hailley, mnewton, etc. without success. Changing the precision and tolerance of the root finding function from mpmath didn't help either. I think this polynomial just behaves too wildly in the region of the plot.

Here's the code and relevant documentation for the "mpmath.findroot" function:

http://mpmath.org/doc/current/calculu... :

from sympy import *
import matplotlib.pyplot as plt
import mpmath as mp

mp.dps = 100
stepx = 0.1
stepy = 0.5
xrange = np.arange(7.5,12.5, stepx)
yrange = np.arange(0, 5, stepy)

def plot_roots_of_f(f):
    L = []
    for u in xrange:
        for v in yrange:
            Root = mp.findroot(f.eval({x:u}),v, method = 'mnewton', tol = E-60, verbose=False,verify=False)
            L.append([u,Root])
    return plt.plot(*zip(*L),linestyle='None', marker='.')

plot_roots_of_f(f)
plt.show()

I'm using Manjaro Linux and I don't remember exactly how I installed it. Probably using a package manager like Octopi or something. Perhaps I should ask this in the distro's forum? I thought it might be a sage problem because the counting command works well for small multiples of the polytope, the Ehrhart quasipolynomial command works and so on, but as you say it's probably a bad installation on my side.

I'm pretty sure it's installed. I actually worked around this issue by computing the Ehrhart quasipolynomial of the polytopes I am working with normaliz. Doesn't this require to have latte_int installed? The package that appears to be installed is "latte_integrale", which I assume is the same?

In any case, this query yields make: *** No rule to make target 'all-toolchain'. Stop. Any suggestion to deal with this issue? >.< I guess it's worth mentioning that this issue appeared recently. An OS update must have messed something up, I just don't know that.

Thanks for your answer!

Hello, I'm doing some polytope computations and the above problem appears as soon as the polytopes become too big. For example, I'm using:

from sage.interfaces.latte import *
p=Polyhedron(vertices=[(1/2,0),(0,1/3),(1,1)], backend='normaliz')
q1=10*p
q1.integral_points_count()
14
q2=100*p
q2.integral_points_count()
---------------------------------------------------------------------------
FeatureNotPresentError                    Traceback (most recent call last)
<ipython-input-6-ca4453575b9a> in <module>()
----> 1 q.integral_points_count()

/usr/lib/python2.7/site-packages/sage/geometry/polyhedron/base_QQ.pyc in integral_points_count(self, verbose, use_Hrepresentation, explicit_enumeration_threshold, preprocess, **kwds)
    216                 cdd=True,
    217                 verbose=verbose,
--> 218                 **kwds)

/usr/lib/python2.7/site-packages/sage/interfaces/latte.pyc in count(arg, ehrhart_polynomial, multivariate_generating_function, raw_output, verbose, **kwds)
    112     """
    113     # Check that LattE is present
--> 114     Latte().require()
    115 
    116     args = ['latte-count']

/usr/lib/python2.7/site-packages/sage/features/__init__.pyc in require(self)
    156         presence = self.is_present()
    157         if not presence:


  --> 158             raise FeatureNotPresentError(self, presence.reason, presence.resolution)
        159 
        160     def __repr__(self):
FeatureNotPresentError: LattE is not available.
Executable 'count' not found on PATH.
To install latte-count you can try to run 'sage -i latte_int'.
Further installation instructions might be available at [DELETED LINK].

Any help to resolve the issue would be greatly appreaciated!