ASKSAGE: Sage Q&A Forum - RSS feedhttps://ask.sagemath.org/questions/Q&A Forum for SageenCopyright Sage, 2010. Some rights reserved under creative commons license.Tue, 20 Apr 2021 20:55:34 +0200Mapping points on a cylinder to points on a cone.https://ask.sagemath.org/question/56745/mapping-points-on-a-cylinder-to-points-on-a-cone/My question has two parts.
1) How to map (not project) 3D coordinates from the surface of a cylinder to the surface of a sphere.
In this image the red squares are projections, the blue square is a (one of many possible) mapping.
![Contrasting projection with mapping](https://i.imgur.com/VLo7Nm6.png)
2) How to do this when the cylinder and cone are not concentric; and when their axes are not in the same plane.
By way of better explanation, this view represents a cutting tool (brown; middle left) rotating around its axis (black double-dashed; right) cutting the sides (orange; upper left) into the workpiece blank (green; left).
(The yellow shadow shows the resulting gear after 30 cuts and indexing.)
![Close up view showing the problem](https://i.imgur.com/OJFdQlV.jpg)
This one shows the same thing from the top and side to emphasise the non-orthogonal relationship between the cutting plane and the workpiece.
![Top and front views establishing context](https://i.imgur.com/FnlRsqo.jpg)
Given that I have a set of 3D points that represent the boundary of the cutting tool, I want to mathematically sweep those points around the axis through the workpiece and then map them into the coordinate space of the workpiece.
I've found a paper ([12MB pdf](https://core.ac.uk/display/10365043) that I believe develops the equations to do this, but my math education was long ago and never covered the Differential geometry used in that paper; and is therefore not strong enough to understand enough to be sure.
I would like help to verify that it does what I think it does; and if so, in working out which of the equations therein I need to program to apply the required transformation(s).
I realise that this problem is ill-defined mathematically as presented; and it is asking a lot for people to read a 16 page paper; but any help offered would be most gratefully received.
(I'm a retired engineer and a competent programmer and this project is a purely personal one with no commercial affiliations.)
Thanks Buk.BukTue, 20 Apr 2021 20:55:34 +0200https://ask.sagemath.org/question/56745/Cython / Sage: how to use Static Typing (cdef)?https://ask.sagemath.org/question/56519/cython-sage-how-to-use-static-typing-cdef/I'm looking at Sage developer guide for Cython:
https://doc.sagemath.org/html/en/developer/coding_in_cython.html
But I cannot find information about how to use static typing, which is what enables Cython (not necessarily in the context of Sage) to reach C-like performance. Since objects in Sage often have complicated types (expressions, polynomials), how would type declarations (cdef in Cython) work?
P.S. Does Cython have some type inference capability (like the "auto" keyword in C++) so the user can avoid typing long type names?greatpetTue, 06 Apr 2021 13:01:20 +0200https://ask.sagemath.org/question/56519/Error with matrix definition of a linear programhttps://ask.sagemath.org/question/55970/error-with-matrix-definition-of-a-linear-program/ According to the folloing tract `https://trac.sagemath.org/ticket/16714`
this command works
p = MixedIntegerLinearProgram()
x = p.new_variable()
p.add_constraint(A_matrix*x <= a_vector)
so why the following code returns an error ?
A=[[1,0,0,1,0,0,0,0,0,0,0,0],[0,1,0,0,0,0,1,0,0,0,0,0],[0,0,1,0,0,0,0,0,0,1,0,0],[0,0,0,0,1,0,0,1,0,0,0,0],[0,0,0,0,0,1,0,0,0,0,1,0],[0,0,0,0,0,0,0,0,1,0,0,1]]
U1=[1,1,1,1,1,1]
show(matrix(A))
show(vector(U1))
p = MixedIntegerLinearProgram(binary = true)
x = p.new_variable()
p.add_constraint(A*x <= U1)
p.show()
CyrilleMon, 01 Mar 2021 12:04:34 +0100https://ask.sagemath.org/question/55970/Integer programming bi-indexed variables.https://ask.sagemath.org/question/55872/integer-programming-bi-indexed-variables/Many integer programming programs as the following needed for Kemeny ranking
$$
\begin{array}{c}
\text{minimiser} \sum_{i,j \in \mathcal{A}}\omega_{i,j} x_{i,j}+\omega_{j,i} x_{j,i}\\\\
\text{sous les contraintes} \\\\
x_{i,j}+ x_{j,i} = 1, \forall i \not=j\\\\
x_{i,j}+ x_{j,k}+ x_{k,i}\geq 1, \forall i \not=j\not=k\not=i
\end{array}
$$
I want to construct a function wich takes the vector or the list $\omega$ and return the solution. But to construct such a function iot will be safer for me not to be obliged to assign on variable with a index to all two indexed variables.
Is there a way ? My problem is in the enumeration of the constraints for high indexes.CyrilleWed, 24 Feb 2021 10:41:01 +0100https://ask.sagemath.org/question/55872/Linear programming with a strange resulthttps://ask.sagemath.org/question/55430/linear-programming-with-a-strange-result/ I have encoutered a strange behavior in solving a linear program.
def sol_zero_sum_game(M=matrix([[1,-1], [-1,1]]),code=1) :
dim = M.nrows()
U=ones_matrix(dim,dim)
zsg=MixedIntegerLinearProgram(maximization=False, solver="GLPK")
x=zsg.new_variable(real=True,nonnegative=True, indices=[0..dim-1])
minM=min(min(M))
Id= identity_matrix(dim)
M1=(abs(minM)+1)*U+M
Bzsgl=M1*x
zsg.set_objective(sum(x[i] for i in range(dim)))
zsg.solve()
xx=zsg.get_values(x)
#show(xx)
for i in range(0,dim) :
zsg.add_constraint(Bzsgl[i]>=1)
if code==1 :
return zsg.show()
if code==2 :
return xx
The above code can display the program for the following `G` matrix :
G=matrix([[1,-1], [-1,1]])
sol_zero_sum_game(G,1)
The solution is given by :
G=matrix([[1,-1], [-1,1]])
sol_zero_sum_game(G,2)
which gives `{0: 0.0, 1: 0.0}` as a solution. But this is obviously false since after substitution of $(0,0)$ in the constrains it doesn't work. As this is the formalisation of a game I know that the solution is $(0.25, 0.25)$ as confirmed by an other software LIPS (see the screen capture below).
![image description](/upfiles/16114008054195529.png)
I have certainly made a mistake but I can't see where. Thanks for help.CyrilleSat, 23 Jan 2021 12:21:18 +0100https://ask.sagemath.org/question/55430/Linear Programminghttps://ask.sagemath.org/question/54672/linear-programming/I wonder if there is not a little display bug inside the SageMath
linear program since the display of the result here is given
before the display of the program (for some other programs
it is given inside).
# Nombre de contraintes et de variables.
nmd, mmd = 10, 6
# Coefficients pour chaque contrainte.
coeffs = [(250, 770, 360, 190, 230, -1),
(0, 0, 0, 0, 0, 1),
(31, 44, 14, 27, 3, 0),
(480, 1770, 800, 580, 160, 0),
(1, 0, 0, 0, 0, 0),
(0, 1, 0, 0, 0, 0),
(0, 0, 1, 0, 0, 0),
(0, 0, 0, 1, 0, 0),
(0, 0, 0, 0, 1, 0),
(0, 0, 0, 0, 0, 1)]
# Matrice des contraintes.
Amd = matrix(nmd, mmd, coeffs)
# Bornes inférieures pour les contraintes.
bmdmin = [0, 600, 30, 0, 0, 0, 0, 0, 0]
# Bornes supérieures pour les contraintes (oo = infini).
bmdmax = [0, 900, 1000, oo, oo, oo, oo, oo, oo, oo]
# Visualisons.
show(LatexExpr("A="), Amd)
show(LatexExpr("bmin="), bmdmin)
show(LatexExpr("bmax="), bmdmax)
# On crée le programme de maximisation.
md = MixedIntegerLinearProgram(maximization=False, solver="GLPK")
# Nouvelles variables: x_0 ... x_2; 'integer' dit si variables à valeurs entières
x = md.new_variable(integer=False, nonnegative=True, indices=[0 .. mmd - 1])
# Fonction linéaire pour les contraintes.
Bmd = Amd * x
# On fixe l'objectif.
md.set_objective(1.2*x[0] + 4.5*x[1] + 3.7*x[2] + 6.3*x[3] + 3.2*x[4])
# On construit les contraintes avec leurs bornes.
for i in range(0,4):
md.add_constraint(Bmd[i], min=bmdmin[i], max=bmdmax[i])
md.show()
md.solve()
md.get_values(x)
xmd = md.get_values(x)
# Le résultat.
show(xmd)
Output:
𝐴=(...)
𝑏𝑚𝑖𝑛=[0,600,30,0,0,0,0,0,0]
𝑏𝑚𝑎𝑥=[0,900,1000,+∞,+∞,+∞,+∞,+∞,+∞,+∞]
Minimization:
1.2 x_0 + 4.5 x_1 + 3.7 x_2 + 6.3 x_3 + 3.2 x_4
Constraints:
0.0 <= 250.0 x_0 + 770.0 x_1 + 360.0 x_2 + 190.0 x_3 + 230.0 x_4 - x_5 <= 0.0
600.0 <= x_5 <= 900.0
30.0 <= 31.0 x_0 + 44.0 x_1 + 14.0 x_2 + 27.0 x_3 + 3.0 x_4 <= 1000.0
0.0 <= 480.0 x_0 + 1770.0 x_1 + 800.0 x_2 + 580.0 x_3 + 160.0 x_4 <= inf
Variables:
x_0 is a continuous variable (min=0.0, max=+oo)
x_1 is a continuous variable (min=0.0, max=+oo)
x_2 is a continuous variable (min=0.0, max=+oo)
x_3 is a continuous variable (min=0.0, max=+oo)
x_4 is a continuous variable (min=0.0, max=+oo)
x_5 is a continuous variable (min=0.0, max=+oo)
{0:2.4000000000000004,1:0.0,2:0.0,3:0.0,4:0.0,5:600.0000000000001}
Of course I am aware that without the `show()` command there is no problem. But, for uniformity of presentation I need this command. As I could call it in an other cell the problem an be circumvented. But for economic reasons I would prefer to have the code and the result in the same cell. Nevertheless, this is just a remark.
**Edit to address a comment by slelievre.**
Inside a notebook this works perfectly. But inside a "Sage Cell" here is the result:
𝐴=(...)
𝑏𝑚𝑖𝑛=[0,600,30,0,0,0,0,0,0]
𝑏𝑚𝑎𝑥=[0,900,1000,+∞,+∞,+∞,+∞,+∞,+∞,+∞]
{0:2.4000000000000004,1:0.0,2:0.0,3:0.0,4:0.0,5:600.0000000000001}
Minimization:
1.2 x_0 + 4.5 x_1 + 3.7 x_2 + 6.3 x_3 + 3.2 x_4
Constraints:
0.0 <= 250.0 x_0 + 770.0 x_1 + 360.0 x_2 + 190.0 x_3 + 230.0 x_4 - x_5 <= 0.0
600.0 <= x_5 <= 900.0
30.0 <= 31.0 x_0 + 44.0 x_1 + 14.0 x_2 + 27.0 x_3 + 3.0 x_4 <= 1000.0
0.0 <= 480.0 x_0 + 1770.0 x_1 + 800.0 x_2 + 580.0 x_3 + 160.0 x_4 <= inf
Variables:
x_0 is a continuous variable (min=0.0, max=+oo)
x_1 is a continuous variable (min=0.0, max=+oo)
x_2 is a continuous variable (min=0.0, max=+oo)
x_3 is a continuous variable (min=0.0, max=+oo)
x_4 is a continuous variable (min=0.0, max=+oo)
x_5 is a continuous variable (min=0.0, max=+oo)
That is, the result, `{0:2.40...,1:0.0,2:0.0,3:0.0,4:0.0,5:600.00...}`,
is printed before the program's output, and sometimes inside it.CyrilleMon, 14 Dec 2020 18:13:04 +0100https://ask.sagemath.org/question/54672/Interactivelinearprogramming gives the wrong answerhttps://ask.sagemath.org/question/53982/interactivelinearprogramming-gives-the-wrong-answer/Here are two versions of the same linear programming problem.
The first one is simply the second one after eliminating the variable $z$.
- $\min$ { $-1.8 x - y + 600 \; | \; x \leq 30; \; y \leq 32 ; \; 2 x + 3 y \leq 90; \; 10 x + 6 y \leq 230$ }
- $\min$ { $4.2 x + 5 y + 6 z \; | \; x \leq 30; \; y \leq 32; \; 2 x + 3 y \leq 90; \; 10 x + 6 y \leq 230; \; x + y + z = 100$ }
Here is some SageMath code to solve the first one
%display latex
A = matrix([[1, 0], [0, 1], [2, 3], [10, 6]])
b = vector([30, 32, 90, 230])
c = vector([-1.8, -12])
# La forme du `problem_type` est obligatoire pour `InteractiveLPProblemStandardForm`
inc = InteractiveLPProblem(A, b, c, ["I_A", "I_B"],
problem_type="min",
objective_constant_term=600)
show(inc)
inc.plot()
inc = inc.standard_form()
inc = inc.run_simplex_method()
show(inc)
But the answer is obviously false. I have used the same type of code to solve the second one.
%display latex
A = matrix([[1, 0, 0], [0, 1, 0], [2, 3, 0], [10, 6, 0], [1, 1, 1]])
b = vector([30, 32, 90, 230, 100])
c = vector([4.2, 5, 6])
inc = InteractiveLPProblem(A, b, c, ["x", "y", "z"],
problem_type="min",
constraint_type=['<=', '<=', '<=', '<=', '=='],
objective_constant_term=600)
show(inc)
inc = inc.standard_form()
inc = inc.run_simplex_method()
show(inc)
Here the answer is infeasibility.
But with the following code I reach the good answer (verified by Lips <- Sourceforge,
which gives the same answer for both problems):
n = 5 # nombre de contraintes
m = 3 # nombre de variables
A = matrix(n, m, [1, 0, 0, 0, 1, 0, 2, 3, 0, 10, 6, 0, 1, 1, 1]) # les coefficients
bmin = [0, 0, 0, 0, 100] # bornes inférieures pour les contraintes
bmax = [30, 32, 90, 230, 100] # bornes supérieures pour les contraintes (oo = infini)
show(LatexExpr("A="), A)
show(LatexExpr("bmin="), bmin)
show(LatexExpr("bmax="), bmax)
sol = MixedIntegerLinearProgram(maximization=False, solver="GLPK") # on crée le programme
x = sol.new_variable(integer=False, indices=[0 .. m - 1]) # les nouvelles variables sont x_0 à x_11
B = A*x # la fonction linéaire pour les contraintes
sol.set_objective(4.2*x[0]+5*x[1]+6*x[2]+600) # fixe l'objectif
# on construit les contraintes avec leurs bornes
for i in range(n):
sol.add_constraint(B[i], min=bmin[i], max=bmax[i])
sol.show()
sol.solve()
xx = sol.get_values(x)
show(LatexExpr(r"I_A^\star ="),
xx[0], LatexExpr(r", I_B^\star ="),
xx[1], LatexExpr(r", E^\star ="),
xx[2], LatexExpr(r" \text{ et } C^\star ="),
sol.get_objective_value())
So I wonder if I made a mistake or if there is an error in the code.
I am perfectly aware that the `interactivelinearprogramming` package
is built for pedagogical reasons. But some things are not very obvious.
For instance it is sensitive to the order of the options.CyrilleTue, 20 Oct 2020 06:12:53 +0200https://ask.sagemath.org/question/53982/One step linear programming pivothttps://ask.sagemath.org/question/52769/one-step-linear-programming-pivot/The following question is not because I do not appreciate the numerous linear programming tools disponible from Sage or Python, but because I want to force my student to have a reflexion on how the Dantzig algorithm works by themselves
Here is some data
A=matrix([[1,2,3,4],[4,3,2,1],[2,4,3,1]])
b=vector([10,20,30])
c=vector([1,-2,3,-1])
With the help of this function
def tableau_(mat, cont, obj):
Id = identity_matrix(mat.nrows())
z = matrix(zero_vector(mat.ncols()-1))
return block_matrix([[A,II,matrix(b).transpose()],[matrix(c),matrix(z),zero_matrix(1,1)]],subdivide=False)
I can construct the following table
mm=matrix(tableau_(A, b, c))
show(mm)
now I want a simple step a the Dantzig pivot algorithm applied to mm throug a calling function like
`one_step_pivot(mm, mm[2][2])` where `mm[2][2]` is the pivot. I wonder if there is a simple tools in Sage or numpy or any other package that do the trick or if I would be oblige to write the code by myself nwhich will tazkes me a lot of time. By advance thanks for the help.CyrilleFri, 31 Jul 2020 08:29:40 +0200https://ask.sagemath.org/question/52769/Linear programming displaying Latex formated equationshttps://ask.sagemath.org/question/52724/linear-programming-displaying-latex-formated-equations/The following works fine:
A = ([1,0], [0, -1], [6,10], [-6,-10])
B = (18, -12,70, -70)
C = (-4.10, -8)
P = InteractiveLPProblem(A, B, C, ["x_1", "x_2"], variable_type=">=")
P
P.plot()
But when I run:
P.run_simplex_method()
SageMath says that the original problem is infeasible and should stop there.
But some LaTeX formatted equations appear after.
This could be troublesome for student so I would be glad to know
how to stop displaying anything after the diagnosis of infeasibility.
Is this possible? My students need to be confronted to infeasibility.CyrilleTue, 28 Jul 2020 16:55:19 +0200https://ask.sagemath.org/question/52724/Warning against a deprecated functionhttps://ask.sagemath.org/question/51576/warning-against-a-deprecated-function/ When programming the following code
def varelmat(M) :
return sum((vector(M)).apply_map((x-(sum(vector(M)))/len(vector(M)))^2))/len(vector(M))
and applying it
B= random_matrix(QQ, 7,7)
varelmat(B)
I end with the following warning
`/opt/sagemath-9.0/local/lib/python3.7/site-packages/sage/repl/ipython_kernel/__main__.py:3: DeprecationWarning: Substitution using function-call syntax and unnamed arguments is deprecated and will be removed from a future release of Sage; you can use named arguments instead, like EXPR(x=..., y=...)
See http://trac.sagemath.org/5930 for details.
IPKernelApp.launch_instance(kernel_class=SageKernel'
Could some one explain me how to keep my code save ?CyrilleWed, 27 May 2020 11:33:48 +0200https://ask.sagemath.org/question/51576/Some code gives error in sagemath 9 but OK in 8.9https://ask.sagemath.org/question/50204/some-code-gives-error-in-sagemath-9-but-ok-in-89/After I upgraded to sagemath 9.0 some of the code that used to work now gives error. The code was used to obtain an estimate of size of expression.
>sage --version
SageMath version 9.0, Release Date: 2020-01-01
>which sage
/bin/sage
It looks like sagemath 9 now uses python 3.0 while 8.9 used python 2 (since I had to change all my print statements to use () to make them work.
Here is an example of function that now gives an error. This is in file, say `bug_sage.py`
#!/usr/bin/env sage
from sage.all import *
def tree(expr):
if expr.operator() is None:
return expr
else:
return [expr.operator()]+map(tree, expr.operands())
var('x')
tree(x*e**((x*log(x) + 1)/log(x)))
Now from command line at Linux, I type `sage ./bug_sage.py` and it gives the error
Traceback (most recent call last):
File "./bug_sage.py", line 12, in <module>
print (tree(x*e**((x*log(x) + 1)/log(x))))
File "./bug_sage.py", line 9, in tree
return [expr.operator()]+map(tree, expr.operands())
TypeError: can only concatenate list (not "map") to list
In 8.9, no error is generated.
It looks like this is due to change in Python itself? This function is meant to generate list of all operands in expression in order to estimate the size of the expression. It is later used as follows
len(flatten(tree(anti)))
Any idea how to fix it? Or do you suggest better way to obtain size of expression (called leaf count in other CAS systems).
NasserMon, 09 Mar 2020 21:52:26 +0100https://ask.sagemath.org/question/50204/Should I avoid using unicode in names of variables?https://ask.sagemath.org/question/49370/should-i-avoid-using-unicode-in-names-of-variables/ Hello, Sage community!
Since Sage now supports Python 3, the possibility of using tildes and other unicode letters in names of variables is real. For example, I could define a constant `Gauß` instead of `Gauss` (with a "ß" instead of "ss"); or I could define a function called `Müller` (with two dots over the "u"). In Spanish and other languages, this is a more pressing matter. If I defined a function `design(howMany)`, in Spanish it would be `diseño(cuántos)`.
However, it seems quite strange to have tildes, betas, etc. in variable names. Is it a good programming behavior to avoid this kind of names?, or, should I use them sparely?
Thanks in advance!dsejasSun, 05 Jan 2020 02:23:51 +0100https://ask.sagemath.org/question/49370/Running Python code from Jupyter and transforming it in a Sage functionhttps://ask.sagemath.org/question/48995/running-python-code-from-jupyter-and-transforming-it-in-a-sage-function/ I would like to know if it is possible to run a python code directly within a notebook and transform it in a sagemath function. The problem is that I am nearly sure that is possible, but I do not know how to install it --- I am under windows an a poor programmer.
For instance if some one has to kindness to take as an exemple https://github.com/fordc5/PrisonersDilemmaCSWorkshop, it will be great.CyrilleWed, 11 Dec 2019 10:01:06 +0100https://ask.sagemath.org/question/48995/SEND+MORE=MONEY (one more time)https://ask.sagemath.org/question/48931/sendmoremoney-one-more-time/I have some heavy problems of formulation with Sagemath. It's certainly a formulation problem --- I will add the theoretrical aspect of the model at the end of the question :
I have corrected my code. May be some errors could stay. But even if I am wrong I need some information on how to write my program. You can find the source in the free paper https://pubsonline.informs.org/doi/pdf/10.1287/ited.2016.0163. There is also an Excel code on the internet. Here is my code
p = MixedIntegerLinearProgram(maximization=True,solver='PPL')
v = p.new_variable(nonnegative=True)# permet de définir v[i, j] aussi bien que v[i]
c = p.new_variable(nonnegative=True)# permet de définir v[i, j] aussi bien que v[i]
p.set_objective(v[0,0])
for i in range(7) : # nombre de lettres (8)
for j in range(9) : #nombre de chiffres
p.set_binary(v[i,j])
#S=sum(j*v[0,j] for j in range(9))
#E=sum(j*v[1,j] for j in range(9))
#N=sum(j*v[2,j] for j in range(9))
#D=sum(j*v[3,j] for j in range(9))
#M=sum(j*v[4,j] for j in range(9))
#O=sum(j*v[5,j] for j in range(9))
#R=sum(j*v[6,j] for j in range(9))
#Y=sum(j*v[7,j] for j in range(9))
for k in range(1,4) :
p.set_binary(c[k])
# Une lettre n'est représentée que par un chiffre
for j in range(8) :
p.add_constraint(0<= sum(v[i,j] for i in range(7)) <= 1)
# Un chiffre n'est associé qu'à une seule lettre
for i in range(7) :
p.add_constraint(sum(v[i,j] for j in range(9)) == 1)
# D + E = Y + 10 C_1
p.add_constraint(sum(j*v[3,j] for j in range(7))+sum(j*v[1,j] for j in range(9))==
sum(j*v[7,j] for j in range(7))+10*c[1])#
# C_1 + R + N = E + 10 C_2
p.add_constraint(c[1]+sum(j*v[6,j] for j in range(7))+sum(j*v[2,j] for j in range(9))==
sum(j*v[1,j] for j in range(9))+10*c[2])#
# C_2 + E + O = N + 10 C_3
p.add_constraint(c[2]+sum(j*v[5,j] for j in range(8))+sum(j*v[1,j] for j in range(9))==
sum(j*v[2,j] for j in range(9))+10*c[3])#
# C_3 + S + M = O + 10 C_4
p.add_constraint(c[3]+sum(j*v[0,j] for j in range(9)) + sum(j*v[4,j] for j in range(9))==
sum(j*v[5,j] for j in range(9))+10*c[4])#
# C_4 = M
p.add_constraint(c[4]==sum(j*v[4,j] for j in range(9)))
p.solve()
p.show()
xx=p.get_values(v)
show(xx)
for j in range(9) :
if xx[0,j]== 1 : print('S='+ `j`)
for j in range(9) :
if xx[1,j]== 1 : print('E='+ `j`)
for j in range(9) :
if xx[2,j]== 1 : print('N='+ `j`)
for j in range(9) :
if xx[3,j]== 1 : print('D='+ `j`)
for j in range(9) :
if xx[4,j]== 1 : print('M='+ `j`)
for j in range(9) :
if xx[5,j]== 1 : print('O='+ `j`)
for j in range(9) :
if xx[6,j]== 1 : print('R='+ `j`)
for j in range(9) :
if xx[7,j]== 1 : print('Y='+ `j`)
I have 8x10 = 80, v[i,j] variables v[i, j] and 4 c[j] variables. I have asked all to be binary. But when I look the printed model I see only 73 variables and from $x_{66}$ to $x_{73}$ they are not binary. More than that the (7,6) variables takes a value of 3/6. It's not particularly simple to associate the variables internal variables to the internal ones.
I need help! (Thanks by advance). The documentation is of no help for me.
MODELIZATION CLUES
In this so called "cryptarythm" one use 8 letters $S_{(1)}E_{(2)}N_{(3)}D_{(4)}M_{(5)}O_{(6)}R_{(7)}Y_{(8)}$ so we define the binary variables
$x_{ij} =
\begin{cases}
1 & \text{if letter } i (i= 1, \ldots 8) \text{ is associated to th digit } j (j = 0,\ldots, 9)\\\\
0 & \text{sinon}
\end{cases}$
We need also 4 variables $C_k$ for the retention.
From this one can write $S = \sum_{1}^9 j x_{1j }$, $E = \sum_{1}^9 j x_{2j }$ and so on. For the constraints one have 17 general constraints
--- each letter is associated to a unic digit $\sum_{j=0}^9 x_{ij} = 1, \text{ for }
i = 1, \ldots 8$
--- each digit is associated to a only one letter $\sum_{i=0}^8 x_{ij} \leq 1, \text{ for }
j = 1, \ldots 9$
Then one hase the pure addition constraints :
$D + E = Y + 10 C_1 \Longleftrightarrow \sum_{j=0}^9 j x_{2j} +\sum_{j=0}^9 j x_{4j}= \sum_{j=0}^9 j x_{8j} + 10 C_1$
$C_1 + R + N = Y + 10 C_2 \Longleftrightarrow C_1+\sum_{j=0}^9 j x_{3j} +\sum_{j=0}^9 j x_{4j}= \sum_{j=0}^9 j x_{2j} + 10 C_2$
$C_2 + E + O = N + 10 C_3 \Longleftrightarrow C_2+\sum_{j=0}^9 j x_{2j} +\sum_{j=0}^9 j x_{6j}= \sum_{j=0}^9 j x_{3j} + 10 C_3$
$C_3 + S + M = O + 10 C_4 \Longleftrightarrow C_3+\sum_{j=0}^9 j x_{1j} +\sum_{j=0}^9 j x_{5j}= \sum_{j=0}^9 j x_{6j} + 10 C_4$
$C_4 = M \Longleftrightarrow C_4 = \sum_{j=0}^9 j x_{5j}$
Here is the solution
$(S=9). (E=5). (N=6). (D=7).$
$+. (M=1). (O=0).(R=8). (E=5).$
$= (M=1). (O=0). (N=6). (E=5). (y=2)$CyrilleTue, 03 Dec 2019 16:01:23 +0100https://ask.sagemath.org/question/48931/The annulus problem in linear programminghttps://ask.sagemath.org/question/48924/the-annulus-problem-in-linear-programming/ I do not know if the following question is a question about Sagemath or a question about my incapacity to use `MixedIntegerLinearProgram`
The problem is to find the smallest annulus for a set of points. Here is a reference `https://www.ti.inf.ethz.ch/ew/courses/CG12/lecture/Chapter%2011%20and%2012.pdf`
# list of point
lis=[(4,1),(5,6),(2,2),(3,4),(6,3.2),(2.6,2.5),(5.8, 2.4),(2,3.4),(4.2,4),(5.6,1)]
from sage.plot.circle import Circle
# Evaluation of the norm of each point
d = [sqrt(lis[i][0]^2+lis[i][1]^2) for i in range(len(lis)-1)]
anmin = MixedIntegerLinearProgram(maximization=True, solver = "GLPK")# création du programme
nvar=len(lis)-1 # nombre de variables
x = anmin.new_variable(integer=False, indices=[0..nvar-1],nonnegative=True) # les nouvelles variables sont x[0] ... x[3]
anmin.set_objective(x[0]-x[1])# fixe l’objectif
# Constraints
for i in range(nvar) :
anmin.add_constraint(x[0]+ 2*(lis[i][0]*x[2]+lis[i][1]*x[3]), min=0, max=n(d[i]^2, digits = 2))
for i in range(nvar) :
anmin.add_constraint(x[1]+ 2*(lis[i][0]*x[2]+lis[i][1]*x[3]), min=n(d[i]^2, digits = 2), max=oo)
anmin.show()
anmin.solve()
xx=anmin.get_values(x)
show(xx)
The result is obviously false. It would be great if somebody helps me find the good answer.
CyrilleSun, 01 Dec 2019 16:58:04 +0100https://ask.sagemath.org/question/48924/SEND+MORE=MONEY (Milp programming)https://ask.sagemath.org/question/48795/sendmoremoney-milp-programming/ Her is the code (unhappily wrong) for the famous verbal arithmetic : SEND+MORE=MONEY coded by Simon (Puzzle—Verbal Arithmetic and Mastermind --- accessible in line)
p = MixedIntegerLinearProgram(maximization=True,solver='GLPK')
v = p.new_variable(nonnegative=True)
p.set_objective(v[0,0])
for i in range(9) :
for j in range(8) :
p.set_binary(v[i,j])
#S=sum(j*v[0,j] for i in range(0,7))
#E=sum(j*v[1,j] for i in range(0,7))
#N=sum(j*v[2,j] for i in range(0,7))
#D=sum(j*v[3,j] for i in range(0,7))
#M=sum(j*v[4,j] for i in range(0,7))
#O=sum(j*v[5,j] for i in range(0,7))
#R=sum(j*v[6,j] for i in range(0,7))
#Y=sum(j*v[7,j] for i in range(0,7))
# Une lettre n'est représentée que par un chiffre
for j in range(8) :
p.add_constraint(sum(v[i,j] for i in range(0,7)) <= 1)
# Un chiffre n'est associé qu'à une seule lettre
for i in range(7) :
p.add_constraint(sum(v[i,j] for j in range(0,8)) == 1)
# le modèle génère 80 variables, les retenues seront notées au delà
# D + E = Y + 10 C_1
p.add_constraint(sum(j*v[3,j] for i in range(0,7))+sum(j*v[1,j] for i in range(0,7))==sum(j*v[7,j] for i in range(0,7))+10*v[8,8])
# C_1 + R + N = E + 10 C_2
p.add_constraint(v[8,8]+sum(j*v[6,j] for i in range(0,7))+sum(j*v[2,j] for i in range(0,7))==sum(j*v[1,j] for i in range(0,7))+10*v[9,8])
# C_2 + E + O = N + 10 C_3
p.add_constraint(v[9,8]+sum(j*v[5,j] for i in range(0,7))+sum(j*v[1,j] for i in range(0,7))==sum(j*v[2,j] for i in range(0,7))+10*v[10,8])
# C_3 + S + M = O + 10 C_4
p.add_constraint(v[10,8]+sum(j*v[0,j] for i in range(0,7)) + sum(j*v[2,j] for i in range(0,7))==sum(j*v[7,j] for i in range(0,7))+10*v[11,8])
# C_4 = M
p.add_constraint(v[11,8]==sum(j*v[4,j] for i in range(0,7)))
p.show()
p.solve()
xx=p.get_values(v)
show(xx)
Obviously there is a mistake in my code that appears clearly with `p.show()`. The $j*V_{i,j}$ are evaluated outside of the context. Secondly, the code would be more readable if I could replace the sums by the letters. I think the evaluation must be postponed but I dont know how.
CyrilleMon, 18 Nov 2019 15:08:50 +0100https://ask.sagemath.org/question/48795/Mixed integer programming constraint definitionhttps://ask.sagemath.org/question/48228/mixed-integer-programming-constraint-definition/Sorry to ask so much questions but I am in hurry. I would like to know why this procedure doesn't work
> A=
> matrix(8,4,(1,1,1,-14,0,1,2,-8,-1,1,1,0,0,0,0,-1,1,0,0,0,0,1,0,0,0,0,0,0,0,1))
> b= matrix(8,1,(0,0,0,-1,0,0,0,0))
> sign=list('==' '==' '>=' '>=' '>='> '>=' '>=' '>=')
> x = vector(var('x',n=4, latex_name='x'))
>B=A*x
>p = MixedIntegerLinearProgram(maximization=False, solver = "GLPK")
> x = p.new_variable(integer=True)
> p.add_constraint(B[0] 'sign[0]' b[0])
for the last line I have tryed
> p.add_constraint(B[0] sign[0] b[0])
> p.add_constraint(B[0] == b[0])
without any success.CyrilleMon, 07 Oct 2019 21:23:21 +0200https://ask.sagemath.org/question/48228/On Linear Programming with a double sumhttps://ask.sagemath.org/question/48761/on-linear-programming-with-a-double-sum/With the help of **nbruin**, this is a perfect way to resolve a double sum problem.
c = matrix(3,3,(-1,2,-3,3,2,-1,2,1,3))
d = matrix(3,3,(1,2,3,3,2,1,2,1,3))
show(c)
show(d)
p = MixedIntegerLinearProgram(solver='GLPK')
v = p.new_variable(nonnegative=True)
for j in range(3):
p.set_objective(sum(c[i,j]*v[i,j] for i in range(0,3) for j in range(0,3)))
p.add_constraint(sum(d[i]*v[i,j] for i in range(0,3)) <= 1)
p.show()
p.solve()
vv=p.get_values(v)
show(vv)
w=v.items()
w
But 1) w is not sorted 2) I would like to construct a table with $x_{(a, b)}$ in one column and the optimal value of $x_c$ which corresponds to the $x_{(a, b)}$ . CyrilleFri, 15 Nov 2019 16:39:39 +0100https://ask.sagemath.org/question/48761/Linear programming sum and double sumhttps://ask.sagemath.org/question/48722/linear-programming-sum-and-double-sum/ I wonder why those syntax are invalid
p = MixedIntegerLinearProgram(solver='GLPK')
v = p.new_variable(nonnegative=True)
p.add_constraint(sum([v[i,1] for i in [0,3]) <= 1)
show(p)
or
p = MixedIntegerLinearProgram(solver='GLPK')
v = p.new_variable(nonnegative=True)
p.add_constraint(sum(sum([v[i,j] for i in [0,3]),for j in [1, 5]) <= 1)
show(p)CyrilleThu, 14 Nov 2019 07:24:32 +0100https://ask.sagemath.org/question/48722/Functions in GAPhttps://ask.sagemath.org/question/48711/functions-in-gap/ How can i write special numbers(Mersenne, Fibonacci etc.) or the average or the GCD's two different algorithm in GAP(Group, Algebra and Programming)?
hayyambeyTue, 12 Nov 2019 21:30:04 +0100https://ask.sagemath.org/question/48711/How do I write and implement a program in Sage for Windows?https://ask.sagemath.org/question/48422/how-do-i-write-and-implement-a-program-in-sage-for-windows/ I downloaded Sage 8.9 for Windows and want to use it to do some matrix computations. I do not have any prior experience in programming. So far I've been able to use the SageMath 8.9 Console to declare variables, define matrices, and multiply them together. But I want to write a program which does these computations successively and have SAGE implement it, rather than me typing in the commands one at a time.
I haven't been able to find how to do this in any tutorials. I don't understand where exactly I should be writing the code, how to save the file, how to get SAGE to read the program I wrote and implement it. D_SSat, 19 Oct 2019 05:29:30 +0200https://ask.sagemath.org/question/48422/Adapt the nauty_directg functionhttps://ask.sagemath.org/question/48275/adapt-the-nauty_directg-function/Dear reader,
For a research interest, I would like to generate the collection of non-isomorphic digraphs on a given graph, with a given subset of the edges fixed. Naturally, simply generating the full collection using digraphs.nauty_directg() and removing the undesirable members of the obtained collection would be functional, but it's quite the waste of computing time to take this approach. Hence, I was wondering whether I'd be able to get under the hood of the nauty_directg function to create my own adapted version, which simply skips over the orientation of the edges listed in some input argument. Would anyone be able to help me get started with this endeavor?
Best,
PepijnPepijnWissingFri, 11 Oct 2019 13:01:38 +0200https://ask.sagemath.org/question/48275/MixedIntegerLinearProgram strange behaviorhttps://ask.sagemath.org/question/48107/mixedintegerlinearprogram-strange-behavior/This is a reedition of my initial question
I find missleading the following comportment (I have add name="w") according the suggestion of dsejas
> %display latex p =
>
> MixedIntegerLinearProgram(maximization=False,
>
> solver = "GLPK")
>
> w = p.new_variable(integer=True,
>
> nonnegative=True,name="w")
>
> p.add_constraint(w[0] + w[1] + w[2] -
> 14*w[3] == 1)
>
> p.add_constraint(w[1] + 2*w[2] -
> 8*w[3] == 0)
>
> p.add_constraint(2*w[2] - 3*w[3] == 0)
>
> p.add_constraint(w[0] - w[1] - w[2] >=
> 0)
>
> p.add_constraint(w[3] >= 1)
>
> p.set_objective(w[2]+ 5*w[3])
>
> p.show()
but look at what is shown. The $x_i$'s are in the display. Then I have typpeset :
> b = p.get_backend()
>
> b.solver_parameter("simplex_or_intopt", "simplex_only")
>
> b.solver_parameter("verbosity_simplex", "GLP_MSG_ALL")
>
> ans = p.solve()
>
>ans
followed by
> val=p.get_values(w)
>
> val
and finally
> val0=p.get_values(w[2])
>
> val0
Look at w[2]'s value. It is not an integer value.CyrilleMon, 30 Sep 2019 22:13:18 +0200https://ask.sagemath.org/question/48107/Interactive linear programminghttps://ask.sagemath.org/question/48092/interactive-linear-programming/1) According to the documentation of "interactive simplex method" the following code seems legitimate
> %display latex A = ([1, 1, 3], [3, 1,
> 2], [2, 0, 2]) b = (1000, 1500, 1100)
> c = (10, 5, 2) P =
> InteractiveLPProblemStandardForm(A, b,
> c,["x_1", "x_2",
> "x_3"],slack_variables=["e_4", "e_5",
> "e_6"],
> problem_type="min",constraint_type=["<=",
> ">=","<="]) show(P)
seems correct. But I receive an error message
TypeError: __init__() got an unexpected keyword argument 'constraint_type'.
2) mai I write "="
3) If I add
> variable_type=[">=", ">="]
what must I enter if there is no restriction on a variable ?CyrilleSun, 29 Sep 2019 09:07:27 +0200https://ask.sagemath.org/question/48092/How to create a "sage program" with command line inputhttps://ask.sagemath.org/question/44708/how-to-create-a-sage-program-with-command-line-input/I am completely new to Sage and Pyhton, so this is probably a very trivial question, but I am not able to find any answer. I have a file script.sage
N=10
sage_command1
sage_command2
sage_command3
print(result)
I use this script via `sage script.sage` and it does what I want. How can I modify this so that I can pass the variable `N` in the command line? Something like `sage script.sage 10` and then it prints the result.
Thank you!RiccardoTue, 18 Dec 2018 23:35:35 +0100https://ask.sagemath.org/question/44708/Using multiple lines of pari/gp code in a Sage notebookhttps://ask.sagemath.org/question/41758/using-multiple-lines-of-parigp-code-in-a-sage-notebook/I am struggling to get multiple lines of pari/gp code working in a Jupyter Sage notebook.
When I enter:
%%gp
for(x=1,10,print(x))
it all works fine, however when I for instance write:
%%gp
for(x=1,10,{
print(x);
})
the system just 'hangs' and doesn't return any output. The same issue occurs in the Cocalc Sage cloud environment. I also tried 'pari/gp in your browser' and there it works fine.
Am I doing something wrong or isn't the multi-line pari/gp option supported in Sage?RuudHSun, 25 Mar 2018 22:58:41 +0200https://ask.sagemath.org/question/41758/Monkey-patching PARI callshttps://ask.sagemath.org/question/39807/monkey-patching-pari-calls/I'm currently trying to implement some number field computations on a very particular number field (specifically, $K = \mathbb{Q}(\zeta_{n-1},n^{1/(n-1)})$). The code I've written so far has some bottlenecks that profiling via %prun has identified, namely:
1. `{method '_nf_rnfeq' of 'sage.libs.cypari2.gen.gen' objects}`, which seems to compute the minimal polynomial of $K$ (currently I'm defining $K$ via `F.<u> = CyclotomicField(n-1)`, and then `K.<a> = F.extension(x^(n-1) - n)`, so $K$ is a relative extension. `_nf_rnfeq` seems to compute the absolute minimal polynomial of a given relative extension (via PARI).
2. `{method 'discriminant' of 'sage.rings.polynomial.polynomial_integer_dense_flint.Polynomial_integer_dense_flint' objects}`, which appears to just be the polynomial discriminant computation.
Assume I've done analysis by hand gives me a much faster computation of the absolute minimal polynomial for this particular number field. I'd then want a version of `_nf_rnfeq` that:
1. Checks if it's being called on my special case, where better algorithms exist. If so, run that.
2. Otherwise, run the traditional algorithm.
Of course, I'd prefer to make this modification without modifying the source code, as that seems to be an especially ugly solution.
One way to implement this is to change is by modifying the relevant method at run-time, via a technique known as "Monkey-Patching" (see [here](https://stackoverflow.com/questions/47503061/replacing-specific-function-calls/47503146#47503146)).
The issue I'm having now is that attempting to monkey-patch the PARI call gives me the error:
`TypeError: can't set attributes of built-in/extension type 'sage.libs.cypari2.gen.gen'`
The code I'm using (in the first cell of my IPython notebook) is below:
import sage.libs.cypari2.gen
orig_nf_rnfeq = sage.libs.cypari2.gen.gen._nf_rnfeq
def _nf_rnfeq(*args, **kwargs):
print("Rnfeq works!")
return orig_nf_rnfeq(*args, **kwargs)
sage.libs.cypari2.gen.gen._nf_rnfeq = _nf_rnfeqorangejakeMon, 27 Nov 2017 00:37:25 +0100https://ask.sagemath.org/question/39807/Can you help me with digits of pi?https://ask.sagemath.org/question/39517/can-you-help-me-with-digits-of-pi/ I want to put the first n digits of pi into a list. I have tried numerical approximation but I was not able to figure out if it could help me .Neither did I find a function for it.veritasMon, 13 Nov 2017 20:22:11 +0100https://ask.sagemath.org/question/39517/Can you help me with the digits of pi?https://ask.sagemath.org/question/39518/can-you-help-me-with-the-digits-of-pi/ I want to put the first n digits of pi into a list. I tried numerical approximation but it didn't give me solutions. Neither did I find function which shows the first n digits not the rounded version of piveritasMon, 13 Nov 2017 20:29:21 +0100https://ask.sagemath.org/question/39518/using sage library in C or C++?https://ask.sagemath.org/question/39196/using-sage-library-in-c-or-c/ Is there a way to run sage library in C/C++ or include the sage in C/C++? For example, I want to calculate minkowski sum, difference or decomposition in C or C++? studentboyTue, 17 Oct 2017 17:57:44 +0200https://ask.sagemath.org/question/39196/