study_conics



This a study of a conic (real) from coefficients
of equation  ax^2+bxy+cy^2+dx+ey+f
It founds: type of conic, center, half-axis, vertexes, foci, eccentricity, axes equation.
Some foci and vertexes  are not real (for simplification of block)
and use of mnewton (if solve not works).

conic(a,b,c,d,e,f):=block([m,m1,detm,detm1,c1,c2,centro,s,k,sem,assi,eq,vert],
m:matrix([a,b/2,d/2],[b/2,c,e/2],[d/2,e/2,f]),
m1:matrix([a,b/2],[b/2,c]),
detm:determinant(m),detm1:determinant(m1),
if detm#0  then print("Real Conic"),
detm2:b^2-4*a*c,
if detm2<0 then print("Ellipse"),
if detm2<0 and a=c and b=0 then print("Circle"),
if detm2>0 and a#c then print("Hyperbola"),
if detm2>0 and  a=c then print("Equilater Hyperbola"),
if detm1#0 then 
(c1:a*x+b/2*y+d/2,
c2:b/2*x+c*y+e/2, 
(if detm1 #0 then centro:first(solve([c1,c2],[x,y])),
print("center",centro),
if b=0 and a =c and detm2<0 then print("radius",1/2*sqrt(d^2+e^2-4*f)),
load(eigen),
s:first(eivals(m1)),
k:-detm/detm1,
sem:[k/first(s),k/second(s)],
k1:float(sqrt(abs(sem))),print("half-axis",sqrt(abs(sem)),"->", k1),
if detm2<0 and second(sem)<first(sem) then ecc:sqrt(1-abs(second(sem)/first(sem))),
if detm2<0 and second(sem)>first(sem) then ecc:sqrt(1-abs(first(sem)/second(sem))),
if detm2>0 and second(sem)<first(sem) then ecc:sqrt(1+abs(second(sem)/first(sem))),
if detm2>0 and second(sem)>first(sem) then ecc:sqrt(1+abs(first(sem)/second(sem))),
assi:1/2*b*(x-first(centro))^2-(a-c)*(x-first(centro))*(y-second(centro))-
b/2*(y-second(centro))^2,
print("axes-equations",rhs(assi)),
eq:a*x^2+b*x*y+c*y^2+d*x+e*y+f,
vert:(solve([eq,assi],[x,y])),
print("eccentricity",ecc,"->", float(ecc)),
print("vertexes",vert),
load(mnewton),
v1:mnewton([assi,eq],[x,y],[first(sem),second(sem)]),
v2:mnewton([assi,eq],[x,y],[first(sem),-second(sem)]),
v3:mnewton([assi,eq],[x,y],[-first(sem),second(sem)]),
v4:mnewton([assi,eq],[x,y],[-first(sem),second(sem)]),
print(v1,v2,v3,v4),
k2:abs(first(sem)-second(sem)),
feq:k2=rhs((x-first(centro))^2+(y-second(centro))^2),
fuochi:solve([assi,feq],[x,y]),
print("foci",fuochi),
f1:mnewton([assi,feq],[x,y],[first(sem),second(sem)]),
f2:mnewton([assi,feq],[x,y],[first(sem),-second(sem)]),
f3:mnewton([assi,feq],[x,y],[-first(sem),second(sem)]),
f4:mnewton([assi,feq],[x,y],[-first(sem),-second(sem)]),
print(f1,f2,f3,f4))),
if detm1=0 then (print("Parabola"),
eq:a*x^2+b*x*y+c*y^2+d*x+e*y+f,
asse:a*diff(eq,x)+b/2*diff(eq,y),
print("axis-equation",asse),
ve:first(solve([eq,asse],[x,y])),
print("vertex",ve),
p:sqrt(-detm/(a+c)^3),
print("Parameter p:",p,"->",float(p)),
fp:rhs(expand((x-first(ve))^2+(y-second(ve))^2)),
fp1:fp=(p/2)^2,
fu:solve([fp1,asse],[x,y]),
print("focus",fu)),
return(done))$

Example:
conic(1,2,1,-13,-11,32);
Real Conic
Parabola
axis-equation4*y+4*x-24
vertex[x=1,y=5]
Parameter p:1/2^(3/2)->0.35355339059327
focus[[x=7/8,y=41/8],[x=9/8,y=39/8]]

conic(9,0,4,-18,16,-11);
Real Conic
Ellipse
center[x=1,y=-2]
half-axis[3,2]->[3.0,2.0]
axes-equations-5*(x-1)*(y+2)
eccentricitysqrt(5)/3->0.74535599249993
vertexes[[x=1,y=1],[x=1,y=-5],[x=-1,y=-2],[x=3,y=-2]]
foci[[x=1,y=-sqrt(5)-2],[x=1,y=sqrt(5)-2],[x=1-sqrt(5),y=-2],[x=sqrt(5)+1,y=-2]]