Automatic determination of i2geo-coded geometric loci.
Use the i2geo File Format to code a geometric locus set. Sage will be used to (symbolically) compute its equation.
The allowed elements (currently) are: Free points, Midpoint(point-point), Point(on Circle, on line), Segment(point-point), Line(point-point, point-line -meaning a parallel), OrthogonalLine, Circle(center-radius, center-point, center-radius_as_segment), Intersect(object-object) and Locus.
Click evaluate. If you do not see evaluate, click %hide and shift+enter. Once the text field is generated, paste the the XML i2geo description in it
To avoid unexpected behaviour, RESTART THE WORKSHEET before every locus computation (Action -> Restart worksheet).
{{{id=6|
%hide
automatic_names(True)
@interact
def read_xml(XML_from_construction='Replace this text with the intergeo XML description of your construction'):
global texto
texto=XML_from_construction
///
}}}
Click evaluate. If you do not see evaluate, click %hide and shift+enter. This cell carries out all neccessary computations.
{{{id=8|
%hide
from xml.etree import ElementTree as ET
cons = ET.XML(texto)
"""next we process the i2geo locus construction"""
BoVar=[x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x13, x14, x15, x16, x17, x18, x19, x20, x21, x22, x23, x24, x25, x26, x27, x28, x29, x30, x31, x32, x33, x34, x35, x36, x37, x38, x39, x40, x41, x42, x43, x44, x45, x46, x47, x48, x49, x50, x51, x52, x53, x54, x55, x56, x57, x58, x59, x60, x61, x62, x63, x64, x65, x66, x67, x68, x69, x70, x71, x72, x73, x74, x75, x76, x77, x78, x79, x80, x81, x82, x83, x84, x85, x86, x87, x88, x89, x90, x91, x92, x93, x94, x95, x96, x97, x98, x99];BoVar.reverse()
Todo={}
def FreePoint(pnt,absc,orde):
"""Adds the point $pnt$ with coordinates $(absc,orde)$ to the geometric construction."""
Todo.update({pnt:{'coords':(absc,orde),'parents':Set([]),'type':'FreePoint','hist':['FreePoint',pnt,absc,orde],'eq':Set([])}})
def BoundedPoint(pnt, absc, orde):
Todo.update({pnt:{'coords':(absc,orde),'type':'BoundedPoint'}})
def is_point(p):
return ((p in Todo.keys()) and (Todo[p]['type']=='FreePoint' or Todo[p]['type']=='BoundedPoint'))
def x(obj):
if is_point(obj):
return Todo[obj]['coords'][0]
def y(obj):
if is_point(obj):
return Todo[obj]['coords'][1]
""" anc is a global variable used in the function ancestors. It must be set to nil before evaluating the function"""
anc=[]
def ancestors(n):
"""Returns a set consisting of $n$ and all its ancestors, if $n$ is a construction object. It returns the empty set in other case."""
global anc
res=anc
if n in Todo.keys():
res.append(n)
t=list(Todo[n]['parents'])
res.extend(t)
for i in t:
ancestors(i)
return Set(res)
def Line(n,p,q):
"""Constructs the line $n$ that goes through the points $p$ and $q$."""
vdir=(x(q) - x(p), y(q) - y(p))
eq=Set([(absc - x(p))*(y(q) - y(p)) - (orde - y(p))*(x(q) - x(p))])
pare=Set([p,q])
Todo.update({n:{'vdir':vdir,'eq':eq, 'parents':pare,'type':'Line','hist':['Line',n,p,q]}})
def Segment(n,p,q):
"""Constructs the line $n$ that goes through the points $p$ and $q$. Defined for GeoGebra compatibility."""
vdir=(x(q) - x(p), y(q) - y(p))
eq=Set([(absc - x(p))*(y(q) - y(p)) - (orde - y(p))*(x(q) - x(p))])
pare=Set([p,q])
Todo.update({n:{'vdir':vdir,'eq':eq, 'parents':pare,'type':'Line','hist':['Segment',n,p,q]}})
def Circle(n,c,p,*q):
"""Constructs the circle $n$ with center the point $c$ and radius (squared!) $p^2$ if $p$ is a number. If $p$ is a point, $n$is the circle
with center $c$ and passing through the point $p$. If there is a fourth argument, the point $q$, then so must be $p$ and the circle has center
$c$ and radius (squared) the distance between $p$ and $q$"""
pare=[c]
if is_point(p):
pare.append(p)
if len(q)==1:
radius2=(x(p)-x(q[0]))^2+(y(p)-y(q[0]))^2
pare.append(q[0])
else:
radius2=(x(p)-x(c))^2+(y(p)-y(c))^2
else:
radius2=p^2
eq=Set([(absc - x(c))^2 + (orde - y(c))^2 - radius2])
Todo.update({n:{'center':c,'radius2':radius2,'eq':eq, 'parents':Set(pare),'type':'Line','hist':['Circle',n,c,p,q]}})
def MidPoint(n,p,q):
"""Constructs the midpoint $n$ of points $p$ and $q$."""
if n in Todo.keys():
eq=Todo[n]['eq']
parents=Todo[n]['parents']
temp=Todo[n]['coords']
else:
eq=Set([])
parents=Set([])
temp=(BoVar.pop(),BoVar.pop())
eq=eq.union(Set([temp[0]-1/2*(x(p)+x(q)),temp[1]-1/2*(y(p)+y(q))]))
parents=parents.union(Set([p,q]))
Todo.update({n:{'coords':temp,'type':'BoundedPoint','eq':eq,'parents':parents,'hist':['MidPoint',n,p,q]}})
def PointOnObject(n,o):
"""Constructs a point $n$ lying on object $o$."""
if n in Todo.keys():
eq=Todo[n]['eq']
temp=Todo[n]['coords']
parents=Todo[n]['parents']
else:
eq=Set([])
parents=Set([])
temp=(BoVar.pop(),BoVar.pop())
flag_subs=False
for i in list(Todo[o]['eq']):
if (absc in list(i.variables()) or orde in list(i.variables())):
flag_subs=True
eq=eq.union(Set([i.subs(absc=temp[0],orde=temp[1])]))
parents=parents.union(Set([o]))
if flag_subs:
Todo.update({n:{'coords':temp,'parents':parents,'eq':eq,'type':'BoundedPoint','hist':['PointOnObject',n,o]}})
def ParallelLine(n,p,l):
"""Constructs the line $n$ parallel to the line $l$ passing through $p$."""
vdir=Todo[l]['vdir']
eq=Set([(absc - x(p))*vdir[1] - (orde - y(p))*vdir[0]])
pare=Set([p,l])
Todo.update({n:{'vdir':vdir,'eq':eq, 'parents':pare,'type':'Line','hist':['ParallelLine',n,p,l]}})
def PerpendicularLine(n,p,l):
"""Constructs the line $n$ perpendicular to the line $l$ passing through $p$."""
vdir=Todo[l]['vdir']
eq=Set([(absc - x(p))*vdir[0] + (orde - y(p))*vdir[1]])
vdir=(-vdir[1],vdir[0])
pare=Set([p,l])
Todo.update({n:{'vdir':vdir,'eq':eq, 'parents':pare,'type':'Line','hist':['PerpendicularLine',n,p,l]}})
def IntersectionObjectObject(n,o1,o2):
"""Constructs a point $n$ as intersection of objects $o1$ and $o2$."""
PointOnObject(n,o1)
PointOnObject(n,o2)
def Locus(n,t,*m):
"""Constructs the locus $n$ with tracer $t$. The third argument, if it exists, $m$, is the mover. If it does not exist, the locus is an
implicit one."""
######
locus_name = n
print "Locus name: " + locus_name
######
if len(m)==1:
mover=str(m[0])
implicit=False
else:
implicit=True
mover=''
if implicit:
parents=Set([t])
else:
parents=Set([t,mover])
vbl='';vblist=[]
anc=[];anc_de_t=ancestors(t)
anc=[];anc_de_m=ancestors_de_m=ancestors(mover)
anc_set=anc_de_t.union(anc_de_m)
for i in anc_set:
if Todo[i]['type']=='BoundedPoint':
vbl=vbl+','+str(x(i))
vblist.append(x(i))
vbl=vbl+','+str(y(i))
vblist.append(y(i))
vbl=vbl.replace(',','',1)
vblist.remove(x(t));vblist.remove(y(t))
Eqs=[]
for i in anc_set:
if Todo[i]['type']=='BoundedPoint':
for j in list(Todo[i]['eq']):
Eqs.append(j)
S = PolynomialRing(QQ, vbl)
S.inject_variables()
I=Eqs*S
H=I.elimination_ideal(eval(str(vblist)))
res=H.gens()
"""res=map(factor,H.gens())"""
eq=[]
a,b=str(x(t)),str(y(t))
for i in res:
i=str(i)
i=i.replace('^','**')
i=i.replace(a,'absc')
i=i.replace(b,'orde')
i=eval(i)
eq.append(i)
T=PolynomialRing(QQ,[absc,orde])
T.inject_variables()
J=eq*T
di=J.dimension()
print 'dimension=',di
if di==1:
type='Line'
eq=Set(eq)
elif di==0:
type='BoundedPoint'
temp=(BoVar.pop(),BoVar.pop())
eqq=[]
for i in eq:
eqq.append(i.subs({absc:temp[0],orde:temp[1]}))
eq=Set(eqq)
elif di==2:
type='Plane'
eq=Set([0])
else:
type='Nil'
eq=Set([1])
##########
print eq
##########
hist=['Locus',n,t]
if not implicit:
hist.append(mover)
diccio=[('implicit',implicit),('tracer',t),('mover',mover),('parents',parents),('eq',eq),('type',type),('hist',hist)]
if type=='BoundedPoint':
diccio.append(('coords',temp))
Todo.update({n:dict(diccio)})
def eval_cons(cons):
"""file format has two parts: elements (with coordinates) and constraints (with "algorithms")"""
elements = cons[0]
constraints = cons[1]
dict_segments_points = {'nameSegment1':'nameFirstEndPoint','nameSegment2':'nameSecondEndPoint'}
"""we read now all the elements from the constraints part"""
for i in constraints:
type = i.tag
if type == 'free_point':
name_free_point = i[0].text
for i in elements:
b = i.items()
type_e = i.tag
name_element = b[0][1]
if type_e == 'point' and name_element == name_free_point:
from fractions import Fraction
coord_x = sage_eval(str(Fraction(i[0][0].text)))
coord_y = sage_eval(str(Fraction(i[0][1].text)))
FreePoint(name_free_point, coord_x, coord_y)
elif type == 'circle_by_center_and_point':
name_circle = i[0].text
center = i[1].text
point_in_circle = i[2].text
Circle(name_circle, center, point_in_circle)
elif type == 'line_segment_by_points':
name_segment = i[0].text
point1 = i[1].text
point2 = i[2].text
Segment(name_segment, point1, point2)
dict_segments_points[name_segment+'1'] = point1
dict_segments_points[name_segment+'2'] = point2
elif type == 'point_on_circle':
name_point = i[0].text
name_circle = i[1].text
PointOnObject(name_point,name_circle)
elif type == 'point_on_line':
name_point = i[0].text
name_line = i[1].text
PointOnObject(name_point,name_line)
elif type == 'line_through_two_points':
name_line = i[0].text
point1 = i[1].text
point2 = i[2].text
Line(name_line,point1,point2)
elif type == 'circle_by_center_and_radius':
#radius given by a segment, not a number
name_circle = i[0].text
name_center = i[1].text
name_segment_radius = i[2].text
name_endpoint1 = dict_segments_points[name_segment_radius+'1']
name_endpoint2 = dict_segments_points[name_segment_radius+'2']
Circle(name_circle,name_center,name_endpoint1,name_endpoint2)
elif type == 'intersection_points_of_two_circles':
#original i2geo code includes two intersection points
name_point = i[0].text
name_circle1 = i[2].text
name_circle2 = i[3].text
IntersectionObjectObject(name_point,name_circle1,name_circle2)
elif type == 'intersection_points_of_circle_and_line':
name_point = i[0].text
name_circle = i[2].text
name_line = i[3].text
IntersectionObjectObject(name_point,name_circle,name_line)
elif type == 'locus_defined_by_point_on_circle':
name_locus = i[0].text
name_tracer = i[1].text
name_mover = i[2].text
Locus(name_locus,name_tracer,name_mover)
elif type == 'midpoint_of_two_points':
name_midpoint = i[0].text
point1 = i[1].text
point2 = i[2].text
MidPoint(name_midpoint,point1,point2)
elif type == 'line_parallel_to_line_through_point':
name_new_line = i[0].text
name_point = i[1].text
name_old_line = i[2].text
ParallelLine(name_new_line,name_point,name_old_line)
elif type == 'line_perpendicular_to_line_through_point':
name_new_line = i[0].text
name_point = i[1].text
name_old_line = i[2].text
PerpendicularLine(name_new_line,name_point,name_old_line)
elif type == 'point_intersection_of_two_lines':
name_point = i[0].text
name_line1 = i[1].text
name_line2 = i[2].text
IntersectionObjectObject(name_point, name_line1, name_line2)
elif type == 'locus_defined_by_point_on_line':
name_locus = i[0].text
name_tracer = i[1].text
name_mover = i[2].text
Locus(name_locus, name_tracer, name_mover)
elif type == 'locus_defined_by_point_on_line_segment':
name_locus = i[0].text
name_tracer = i[1].text
name_mover = i[2].text
Locus(name_locus, name_tracer, name_mover)
def get_eq_internal(n):
e=list(Todo[n]['eq'])
for i in e:
if Set(map(str,i.variables())).union(Set(['absc','orde']))!=Set(['absc','orde']):
return 'The object is parametrically defined'
if len(e)==1:
res=[]
for i in e[0].factor_list():
res.append(i[0])
return res
elif len(e)==0:
return 'The object has no equations'
else:
res=[]
for j in e:
restemp=[]
for i in j.factor_list():
restemp.append(i[0])
res.append(restemp)
return res
def get_eq(n):
e=list(Todo[n]['eq'])
for i in e:
if Set(map(str,i.variables())).union(Set(['absc','orde']))!=Set(['absc','orde']):
return 'The object is parametrically defined'
if len(e)==1:
res=[]
for i in e[0].factor_list():
res.append(i[0].subs({absc:X,orde:Y}))
return res
elif len(e)==0:
return 'The object has no equations'
else:
res=[]
for j in e:
restemp=[]
for i in j.factor_list():
restemp.append(i[0].subs({absc:X,orde:Y}))
res.append(restemp)
return res
def box_dim2d(e):
"""Returns adequate box dimensions for plotting the 2d-geometric object described by $e$, if $e$ is compact. Usa .roots(ring=RR) sobre los polinomios univariables para calcular las raices"""
e1,e2=diff(e,absc),diff(e,orde)
S = PolynomialRing(QQ, [absc,orde])
S.inject_variables()
I=[e,e2]*S
H=I.elimination_ideal([orde])
if H.is_principal() and not H.is_one() and not H.is_zero():
eq=gens(H)[0].subs({absc:X,orde:Y})
if eq==0:
m0,M1,px=0,0,0
else:
raices=eq.roots(ring=RR,multiplicities=False)
if raices==[]:
return False
M1,m0=max(raices),min(raices)
px=(M1-m0)/10
if px==0:
px=.5
else:
return False
I=[e,e1]*S
H=I.elimination_ideal([absc])
if H.is_principal() and not H.is_one() and not H.is_zero():
eq=gens(H)[0].subs({absc:X,orde:Y})
if eq==0:
m2,M3,py=0,0,0
else:
raices=eq.roots(ring=RR,multiplicities=False)
if raices==[]:
return False
M3,m2=max(raices),min(raices)
py=(M3-m2)/10
if py==0:
py=.5
else:
return False
"""px and py extend the drawing box"""
return [m0-px,M1+px,m2-py,M3+py]
def get_plot(l):
l=get_eq_internal(l)
if l==[]:
return 'There is nothing to plot'
else:
a=Graphics()
for i in l:
i=expand(i)
dim=box_dim2d(i)
if not (dim==False):
a=a+implicit_plot(i,(absc,dim[0],dim[1]),(orde,dim[2],dim[3]),aspect_ratio=1)
else:
a=a+implicit_plot(i,(absc,dims[0],dims[1]),(orde,dims[2],dims[3]),aspect_ratio=1)
return a
def get_plot_window(l):
l=get_eq_internal(l)
if l==[]:
return 'There is nothing to plot'
else:
a=Graphics()
for i in l:
i=expand(i)
dim=box_dim2d(i)
a=a+implicit_plot(i,(absc,dims[0],dims[1]),(orde,dims[2],dims[3]),aspect_ratio=1)
return a
eval_cons(cons)
///
}}}
Click evaluate. If you do not see evaluate, click %hide and shift+enter. Use the check boxes to show the equation and/or graph of the locus.
{{{id=13|
%hide
@interact
def ex(locus_name='Overwrite this text with your locus name, click in the menu area and check Show Equation',show_equation=checkbox(False,label='Show Equation'), plot_locus=checkbox(False,label='plot locus
(compact w/o singular points)'), window_locus=checkbox(False,label='plot locus
(in window)'), window_plot='[-5,5,-5,5]'):
global show_c
global ec
global dims
show_c=False
if show_equation:
ec=get_eq(locus_name)
print "Equation of locus: " + (str(ec).replace('[','')).replace(']','') + " = 0"
if plot_locus:
show(get_plot(locus_name))
if window_locus:
dims=eval(window_plot)
show(get_plot_window(locus_name))
///
}}}
{{{id=14|
///
}}}