Un problème de Nguyen

Bonjour

sagemath via

_______________________

def norm(P):
    return P/(Linf*P)

def vecteur(P1,P2):
    return norm(P2)-norm(P1)

def perp(droite,point):
    Pinf=matM*droite
    return point.cross_product(Pinf)

def simple(vec):
    if len(vec)==3:
        num=gcd([numerator(vec[0]), numerator(vec[1]), numerator(vec[2])])
        den=lcm([denominator(vec[0]),denominator(vec[1]),denominator(vec[2])])
        facteur=num/den
        vec0=factor(vec[0]/facteur)
        vec1=factor(vec[1]/facteur)
        vec2=factor(vec[2]/facteur)
        return vector([vec0,vec1,vec2])
    
def distance2(P1,P2) :
    vec=vecteur(P1,P2)
    return vec*pyth*vec

def ver(P) :
    return vector([P[0]*(Linf*P),P[1]*(Linf*P),P[2]*(Linf*P),P*pyth*P])

def cercle(vec1,vec2,vec3):
    mat=transpose(matrix([ver(vec1),ver(vec2),ver(vec3)]))
    c1=det(mat[(1,2,3),:])
    c2=-det(mat[(0,2,3),:])
    c3=det(mat[(0,1,3),:])
    c4=-det(mat[(0,1,2),:])
    return vector([c1,c2,c3,c4])

def gram(cer1,cer2):
    return (cer1*matQ*cer1)*(cer2*matQ*cer2)-(cer1*matQ*cer2)^2

def centre(cer):  #TriangleCenter(P1,P2,P3,3)
    cen=matQ*cer
    return vector([cen[0], cen[1], cen[2]])

def solbar(vec):
    if len(vec)==3:
        vec0=vec[0].rhs()
        vec1=vec[1].rhs()
        vec2=vec[2].rhs()
        return vector([vec0,vec1,vec2])  
def cercle2(centre,point):
    return invQ*ver(norm(centre))- distance2(centre,point)* vector([Linf[0],Linf[1],Linf[2],0])    

#version coordonnées barycentriques
var('a b c' , domain='positive') #a=Distance(B,C) ; b=Distance(C,A) : c=Distance(A,B)
Sa=(b^2+c^2-a^2)/2 ; Sb=(a^2+c^2-b^2)/2 ; Sc=(a^2+b^2-c^2)/2 # notations de Conway
S=sqrt((a+b+c)*(b+c-a)*(c+a-b)*(a+b-c))/2  # S=2*aire(ABC) comme dans ETC
Linf=vector([1,1,1])
pyth=1/2*matrix([[0,-c^2,-b^2],[-c^2,0,-a^2],[-b^2,-a^2,0]]) #cercle circonscrit à ABC
matM=matrix([[a^2,-Sc,-Sb],[-Sc,b^2,-Sa],[-Sb,-Sa,c^2]]) #orthodir à multiplier par 1/S
matQ=-matrix([[a^2,-Sc,-Sb,-a^2*Sa],[-Sc,b^2,-Sa,-b^2*Sb],[-Sb,-Sa,c^2,-c^2*Sc],[-a^2*Sa,-b^2*Sb,-c^2*Sc,a^2*b^2*c^2]])/(2*S^2) #quartique de l’espace des cycles
invQ=matrix([[0,c^2,b^2,1],[c^2,0,a^2,1],[b^2,a^2,0,1],[1,1,1,0]])

I=vector([a,b,c])
Ia=vector([-a,b,c])
A=vector([1,0,0]);B=vector([0,1,0]);C=vector([0,0,1])

AB=A.cross_product(B);BC=B.cross_product(C);CA=C.cross_product(A)
Ia=simple(BC.cross_product(perp(BC,I)));Ib=simple(CA.cross_product(perp(CA,I)));Ic=simple(AB.cross_product(perp(AB,I)))
inscrit=cercle(Ia,Ib,Ic)

Ap=norm(B)+norm(C);Bp=norm(A)+norm(C);Cp=norm(A)+norm(B)
euler=cercle(Ap,Bp,Cp)
J=simple(centre(euler))/8

print('test de tangence =>',factor(gram(inscrit,euler)))
var('x,y,z')
M=vector([x,y,z])
sol=solve([inscrit*ver(M)==0,euler*ver(M)==0,Linf*M==1],x,y,z)
Fe=simple(solbar(sol[0]))
print('I =',I)
print('J =',J)
print('Fe =',Fe)
AI=A.cross_product(I)
var('x y z')
M=vector([x,y,z])
sol=solve([cercle(A,B,C)*ver(M)==0,AI*M==0,Linf*M==1],x,y,z)
D=simple(solbar(sol[0]))
print (D)
K=norm(C)+norm(B)
L=norm(I)+norm(D)
M=norm(I)+norm(A)
KL=K.cross_product(L)
KLinf=Linf.cross_product(KL)
t=M.cross_product(KLinf)
print(factor(t*Fe))

________________________________

La droite (t) passe par Fe.

$\square$

Cordialement