Surface de Clebsch
(droites, points d'intersection, pts d'Eckardt)
Les calculs suivants sont les résultats de demandes successives de mon collègue LG Vidiani qui effectuait des recherches sur la surface de Clebsch.
Des animations de la surface sont visibles à cette page
Des représentations Povray de la surface, d'un double-six sont visibles à cette page
>
restart:
#================= Cubique 1 (Clebsch) =====================
S:=81*(x^3+y^3+z^3)-189*(x^2*y+x^2*z+y^2*x+y^2*z+z^2*x+z^2*y)+54*(x*y*z)+126*(x*y+x*z+y*z)-9*(x^2+y^2+z^2)-9*(x+y+z)+1 ;
resz:=subs(x=a*z+p,y=b*z+q,S):# droites à z non cst
resx:=subs(y=a*x+p,z=b*x+q,S):# droites à x non cst
resx0:=subs(x=p,z=q,S):# droites à x cst et z cst
Hx:=collect(resx,x):
Hx0:=collect(resx0,y):
> Hz:=collect(resz,z):
>
Cz:=[seq(coeff(Hz,z,i),i=0..3)]:
Cx:=[seq(coeff(Hx,x,i),i=0..3)]:
Cx0:=[seq(coeff(Hx0,y,i),i=0..3)]:
#-------------- solz----------------
solutionsz :=map(allvalues,{solve({seq(Cz[i]=0,i=1..4)},{a,b,p,q})}):
solz:=map(u -> collect(subs(op(u),[a*z+p,b*z+q,z]),[x,y,z]),solutionsz):
print(seq(op([cat(`droite n = `,i,` `),solz[i]]),i=1..nops(solz)));
printf("%s %a\n","nbre droites non horizontales = ",nops(solz));
#-------------- solx0----------------
solutionsx0 :=map(allvalues,{solve({seq(Cx0[i]=0,i=1..4)},{p,q})}):
solx0:=map(u -> radnormal(subs(op(u),[p,y,q])),solutionsx0):
print(seq(op([cat(`droite n = `,nops(solz)+i,` `),solx0[i]]),i=1..nops(solx0)));
printf("%s %a\n","nbre droites horizontales avec x cste = ",nops(solx0));
#-------------- solx----------------
solutionsx :=map(allvalues,{solve({seq(Cx[i]=0,i=1..4)},{a,b,p,q})}):
# ==== j'élimine les solutions qui sont à z non constant =========
solutionsxbis:=NULL:
for i to nops(solutionsx) do
if member(b=0 , solutionsx[i]) then solutionsxbis:=solutionsxbis,solutionsx[i] fi; od:
solutionsxbis:={solutionsxbis}:
solx:=map(u -> radnormal(subs(op(u),[x,a*x+p,b*x+q])),solutionsxbis):
print(seq(op([cat(`droite n = `,nops(solz)+nops(solx0)+i,` `),solx[i]]),i=1..nops(solx)));
printf("%s %a\n", "nbre droites horizontales avec x non cste = ",nops(solx));
printf("%s %a\n","nombre de droites contenues dans la surface = ",nops(solx)+nops(solz)+nops(solx0));
nbre droites non horizontales = 22
nbre droites horizontales avec x cste = 0
nbre droites horizontales avec x non cste = 5
nombre de droites contenues dans la surface = 27
II. Représentation graphique
>
with(plots):
S:=81*(x^3+y^3+z^3)-189*(x^2*y+x^2*z+y^2*x+y^2*z+z^2*x+z^2*y)+54*(x*y*z)+126*(x*y+x*z+y*z)-9*(x^2+y^2+z^2)-9*(x+y+z)+1 :
deb:=-1: fin := 1.5:debz:=-1:finz:=1.5:
r:=0.005*(finz-debz):
g_droites:= tubeplot({op(solz)},radius=r,z=debz..finz,numpoints=2, color=magenta,style=patchnogrid),
tubeplot({op(solx)},radius=r,x=deb..fin,numpoints=2, color=black,style=patchnogrid),
tubeplot({op(solx0)},radius=r,y=deb..fin,numpoints=2, color=black,style=patchnogrid):
surf:=implicitplot3d(S,x=deb..fin,y=deb..fin,z=debz..finz,style=patchnogrid,grid=[50,50,50]):
>
plotsetup(inline); #window
display([surf,g_droites],scaling=constrained,view=[deb..fin,deb..fin,debz..finz],
orientation=[60,125],lightmodel=light2);
# light=[30,40,0.9,0.9,0],light=[120,20,0.8,0.6,0]2
> restart:
>
#------------ Utilitaires à compiler au préalable ----------------
#### recherche du numéro d'ordre dans l'ensemble "droites" #####
numero:=proc(dr)
local j:
for j from 1 while (j<= nops(droites)) and evalb(droites[j]<>dr) do od:
j;
end:
#### Numéros d'ordre d'une partie de l'ensemble "droites" ######
# ------- 1 pour faire afficher 0 sinon --------------
ordre_partie:= proc(partie,affichage::{0,1})
local ordre,i:
ordre:=[seq(0,i=1..nops(partie))]:
for i from 1 to 6 do
ordre[i]:=numero(partie[i]);
od:
if affichage = 1 then
printf("Numéros des droites %a",ordre);
>
for i from 1 to nops(ordre) do print(partie[i],"___",droites[ordre[i]]) od;
fi;
ordre;
end:
III. Droites et double-six
Deux exécutions de la recherche des droites ne donnent pas la liste des droites dans le même ordre, ordre qui n'a d'ailleurs aucune logique apparente
a. Ensemble des droites avec ordre fixé pour la suite
>
droites:= [ [x, -x, -1/3], [x, -x+1/3, 0], [x, -x+2/3, 1/3],
b. Equations des droites du double-six de "base" et leurs numéros
>
############ double-six initial
Numéros d'ordre dans l'ensemble "droites" des droites [a1,.., a6] et [b1,.., b6]
>
numa:=ordre_partie(a,1);# 1 pour afficher, 0 pour résultat seul
IV. Points d'intersection d'un ensemble de droites ( en particulier, les 27 droites )
>
#### 0 pas d'affichage, 1 affichage des points
# ====== faire effectuer l'instruction suivante désirée ============
>
#### si on travaille sur la solution calculée
>
#### si on travaille sur l'ordre pris dans fichier de données
>
#### pour un double-six de la liste,
>
#### pour le double-six calculé
[(-3+sqrt(5))*z+7/6-1/2*sqrt(5), sqrt(5)*z+1/2-1/6*sqrt(5), z],
[0, 3*z-1/3, z],
[-1/5*sqrt(5)*z+1/6+1/10*sqrt(5), (1+3/5*sqrt(5))*z+1/6+1/30*sqrt(5), z],
[-sqrt(5)*z+1/2+1/6*sqrt(5), (-3-sqrt(5))*z+7/6+1/2*sqrt(5), z],
[x, 3*x-1/3, 0],
[(-3/4-1/4*sqrt(5))*z+1/4-1/12*sqrt(5),(-5/4-3/4*sqrt(5))*z+1/12+1/12*sqrt(5), z],
[1/3, -z+2/3, z], [0, -z+1/3, z], [-1/3, -z, z],
[(-5/4+3/4*sqrt(5))*z+1/12-1/12*sqrt(5), (-3/4+1/4*sqrt(5))*z+1/4+1/12*sqrt(5), z],
[1/5*sqrt(5)*z+1/6-1/10*sqrt(5), (1-3/5*sqrt(5))*z+1/6-1/30*sqrt(5), z],
[1/3*z+1/9, 0, z], [0, 1/3*z+1/9, z],
[(1-3/5*sqrt(5))*z+1/6-1/30*sqrt(5), 1/5*sqrt(5)*z+1/6-1/10*sqrt(5), z],
[(-3/4+1/4*sqrt(5))*z+1/4+1/12*sqrt(5), (-5/4+3/4*sqrt(5))*z+1/12-1/12*sqrt(5), z],
[-z, -1/3, z], [-z+1/3, 0, z], [-z+2/3, 1/3, z],
[(-5/4-3/4*sqrt(5))*z+1/12+1/12*sqrt(5), (-3/4-1/4*sqrt(5))*z+1/4-1/12*sqrt(5), z],
[x, 1/3*x+1/9, 0],
[(-3-sqrt(5))*z+7/6+1/2*sqrt(5), -sqrt(5)*z+1/2+1/6*sqrt(5), z],
[(1+3/5*sqrt(5))*z+1/6+1/30*sqrt(5), -1/5*sqrt(5)*z+1/6+1/10*sqrt(5), z],
[3*z-1/3, 0, z], [sqrt(5)*z+1/2-1/6*sqrt(5), (-3+sqrt(5))*z+7/6-1/2*sqrt(5), z]]:
############""" affichage éventuel ##############
#seq(op([printf("Droite n° %d",i),print(droites[i])]),i=1..27):
a:=[[(-3/4-1/4*sqrt(5))*z+1/4-1/12*sqrt(5), (-5/4-3/4*sqrt(5))*z+1/12+1/12*sqrt(5), z], [x, 1/3*x+1/9, 0], [(1-3/5*sqrt(5))*z+1/6-1/30*sqrt(5), 1/5*sqrt(5)*z+1/6-1/10*sqrt(5), z], [1/3*z+1/9, 0, z], [0, 3*z-1/3, z], [sqrt(5)*z+1/2-1/6*sqrt(5), (-3+sqrt(5))*z+7/6-1/2*sqrt(5), z]]:
b:=[[x, 3*x-1/3, 0], [(-5/4-3/4*sqrt(5))*z+1/12+1/12*sqrt(5), (-3/4-1/4*sqrt(5))*z+1/4-1/12*sqrt(5), z], [0, 1/3*z+1/9, z], [1/5*sqrt(5)*z+1/6-1/10*sqrt(5), (1-3/5*sqrt(5))*z+1/6-1/30*sqrt(5), z], [(-3+sqrt(5))*z+7/6-1/2*sqrt(5), sqrt(5)*z+1/2-1/6*sqrt(5), z], [3*z-1/3, 0, z]]:
numb:=ordre_partie(b,0);
ptscom:=proc(tot,affichage::{0,1})
local i,j,n,sol,pts :
n:=nops(tot):pts:=NULL:
for i from 1 to n-1 do
for j from i+1 to n do
sol:=[solve({tot[i,1]-tot[j,1], tot[i,2]-tot[j,2],
tot[i,3]-tot[j,3]}, {x,y,z})]:
if sol<>[] then
sol:=op(sol);;
if member('x=x',sol) then sol:=subs('x=x'=NULL,sol)
elif member('y=y',sol) then sol:=subs('y=y'=NULL,sol)
elif member('z=z',sol) then sol:=subs('z=z'=NULL,sol)
fi;
pts:=pts,[radnormal(subs(sol,tot[i]),'rationalized'),{i,j}]
fi;
od;
od:
pts:=[pts]:
if affichage =1 then
printf(" Nombre de droites = %d\n",nops(tot));
printf("Nombre de points d'intersection = %a %s",nops(pts),
"--- [ point,[n° dr1 dans partie, n° dr2 dans partie]]");
print(pts);
fi;
pts;
end:
ptsdroites:=ptscom([op(solz),op(solx0),op(solx)],0):
ptsdroites:=ptscom(droites,0):
#[11,26,1,25,9,21], [8,18,23,15,5,10] est le premier élément de la liste
ptsSix:=ptscom(map(u ->droites[u],[11,26,1,25,9,21, 8,18,23,15,5,10] ),1);
ptssix:=ptscom([op(a),op(b)],1):
Début  
V. Nomemclature Schläfli des 15 droites complémentaires du double-six
>
c:='c':
for i from 1 to 6 do
dr1:=droites[numa[i]]:
for j from i+1 to 6 do
dr2:=droites[numb[j]]:
trouve:=false:
for k from 1 to nops(ptsdroites) while not(trouve) do
tamp:=ptsdroites[k,2]:
if member(numa[i],tamp) and
evalb({numa[i],numb[j]}<>tamp) then
tamp1:=op(subs(numa[i]=NULL,tamp)):
for t from 1 while t<= nops(ptsdroites) and
evalb({numb[j],tamp1}<>ptsdroites[t,2]) do od;
if t<nops(ptsdroites)+1 then trouve:=true;
c[10*i+j]:=op(subs(numa[i]=NULL,
numb[j]=NULL,ptsdroites[t,2])):
fi;
fi;
od:
od;
od:
>
#--- affichage des c_ij, n° et équations -----
for i from 1 to 6 do
for j from i+1 to 6 do
printf("c%a%a = droites[%a] ; ",i,j,c[10*i+j])
od:od:
for i from 1 to 6 do
for j from i+1 to 6 do
printf("c%a%a = %a ; ",i,j,droites[c[10*i+j]])
od:od:
c12 = droites[2] ; c13 = droites[7] ;
c14 = droites[10] ; c15 = droites[19] ;
c16 = droites[25] ; c23 = droites[21] ;
c24 = droites[24] ; c25 = droites[6] ;
c26 = droites[12] ; c34 = droites[1] ;
c35 = droites[11] ; c36 = droites[18] ;
c45 = droites[13] ; c46 = droites[20] ;
c56 = droites[3] ;
c12 = [x, -x+1/3, 0] ;
c13 = [-5^(1/2)*z+1/2+1/6*5^(1/2), (-3-5^(1/2))*z+7/6+1/2*5^(1/2), z] ;
c14 = [1/3, -z+2/3, z] ; c15 = [-z, -1/3, z] ;
c16 = [(1+3/5*5^(1/2))*z+1/6+1/30*5^(1/2), -1/5*5^(1/2)*z+1/6+1/10*5^(1/2), z] ;
c23 = [-z+2/3, 1/3, z] ;
c24 = [(-3-5^(1/2))*z+7/6+1/2*5^(1/2), -5^(1/2)*z+1/2+1/6*5^(1/2), z] ;
c25 = [-1/5*5^(1/2)*z+1/6+1/10*5^(1/2), (1+3/5*5^(1/2))*z+1/6+1/30*5^(1/2), z] ;
c26 = [-1/3, -z, z] ; c34 = [x, -x, -1/3] ; c35 = [0, -z+1/3, z] ;
c36 = [(-3/4+1/4*5^(1/2))*z+1/4+1/12*5^(1/2), (-5/4+3/4*5^(1/2))*z+1/12-1/12*5^(1/2), z] ;
c45 = [(-5/4+3/4*5^(1/2))*z+1/12-1/12*5^(1/2), (-3/4+1/4*5^(1/2))*z+1/4+1/12*5^(1/2), z] ;
c46 = [-z+1/3, 0, z] ; c56 = [x, -x+2/3, 1/3] ;
>
# les indices ij des droites c_ij
les_c_indices:=map( x -> op(x),[indices(c)]);
# --------- les n° des droites c_ij ---------
les15num:=map( x -> op(x),[entries(c)]);
verif:={op(les15num),op(numa),op(numb)};nops(verif);
VI. Les autres doubles-six
( faire compiler "les_c_indices" et "les15num" ci-dessus)
a. recherche
Il y a 15 Nij et 20 Nijk
avec Nij=[a_i, b_i ,c_{j,k}, c_{j,l}, c_{j,m}, c_{j,n}/
a_j, b_j, c_{i,k}, c_{i,l}, c_{i,m,} c_{i,n}]
et Nijk= [a_i , a_j, a_k, c_{m,n}, c_{l,n}, c_{l,m}/
c_{j,k}, c_{i,k}, c_{i,j}, b_l, b_m, b_n]
#### Nijk =[ a = [a_i , a_j, a_k, c_{m,n}, c_{l,n}, c_{l,m}] ; b = [ c_{j,k}, c_{i,k}, c_{i,j}, b_l, b_m, b_n]
>
D_S:=[]:
########### les Nijk ###########
for i from 1 to 6 do
for j from i+1 to 6 do
for k from j+1 to 6 do
b_ind:=subs(i=NULL,j=NULL,k=NULL,[seq(t,t=1..6)]);
Nijk[i,j,k]:=[ [numa[i],numa[j],numa[k],
c[10*b_ind[2]+b_ind[3]],
c[10*b_ind[1]+b_ind[3]],
c[10*b_ind[1]+b_ind[2]]],
[c[10*j+k], c[10*i+k], c[10*i+j],
numb[b_ind[1]],numb[b_ind[2]],numb[b_ind[3]]] ]:
D_S:=[op(D_S),Nijk[i,j,k]]:
od;
od;
od:
nds:=nops(D_S):
# --------- affichage des N_ijk ? --------------
> #printf("Nombre de N_ijk = %d",nds);D_S;
### Nij = [ a = [a_i, b_i ,c_{j,k}, c_{j,l}, c_{j,m}, c_{j,n}] ; b= [ a_j, b_j, c_{i,k}, c_{i,l}, c_{i,m,} c_{i,n}]
>
########### les Nij ###########
Nij:='Nij':
for i from 1 to 6 do
for j from i+1 to 6 do
reste:=subs(i=NULL,j=NULL,[seq(t,t=1..6)]);
Nij[i,j]:=[ [numa[i],numb[i],
c[10*min(j,reste[1])+max(j,reste[1])],
c[10*min(j,reste[2])+max(j,reste[2])],
c[10*min(j,reste[3])+max(j,reste[3])],
c[10*min(j,reste[4])+max(j,reste[4])]],
[numa[j],numb[j],
c[10*min(i,reste[1])+max(i,reste[1])],
c[10*min(i,reste[2])+max(i,reste[2])],
c[10*min(i,reste[3])+max(i,reste[3])],
c[10*min(i,reste[4])+max(i,reste[4])]] ]:
D_S:=[op(D_S),Nij[i,j]]:
od:od:
# --------- affichage des N_ij ? --------------
> #printf("Nombre de N_ij = %d",nops(D_S)-nds);D_S[nds+1..nops(D_S)];
>
#-------- affichage final des doubles-sis -----------
printf("Nombre de doubles-six = %d : [a1, ..a6],[b1,..b6]]",nops(D_S)); D_S;
Nombre de doubles-six = 35 : [[a1, ..a6],[b1,..b6]]
b. Vérification
>
#--- le premier est le double-six de base
six_ini:=[[ordre_partie(a,0),ordre_partie(b,0)],op(D_S)] :
touind := {{1, 8}, {1, 9}, {1, 10}, {1, 11}, {1, 12}, {2, 7}, {2, 9}, {2, 10}, {2, 11}, {2, 12}, {3, 7}, {3, 8}, {3, 10}, {3, 11}, {3, 12}, {4, 7}, {4, 8}, {4, 9}, {4, 11}, {4, 12}, {5, 7}, {5, 8}, {5, 9}, {5, 10}, {5, 12}, {6, 7}, {6, 8}, {6, 9}, {6, 10}, {6, 11}}:
>
# verif des six_ini
for i from 1 to nops(six_ini) do
printf("\nSIX_INI[%a] = %a",i,six_ini[i]);
tamp:=ptscom(map(u -> droites[u],map(op,six_ini[i])),0):
# élimination de l'équation
ind:=convert(map( u ->op(subsop(1=NULL,u)),tamp),set);
reste:= touind minus ind:
if nops(tamp)<30 then
printf("\n%s : %a","Couple de droites non sécantes = ",reste);
print(map( w -> [ droites[six_ini[i,1,w[1]]],
droites[six_ini[i,2,w[2]-6]]] ,reste)):
else printf("\n %s\n","Toute droite ai coupe toute droite bj, i<>j"):
fi;
od:
VII. Points de Eckardt à distance finie
>
# Eckardt
pts:=ptscom(droites,0):
p:=nops(pts):
Ec:=NULL:
for i from 1 to p-1 do
tamp:=pts[i,1]: dr:={pts[i,2]}; app:=1;
for j from i+1 to p do
if pts[j,1]=tamp then
app:=app+1:dr:=dr union {pts[j,2]};
fi;
od:
if app=3 then Ec:=Ec,[pts[i,1],dr] fi;
od:
Eckardt:=[Ec]:
printf("%s %a","Nombre de points de Eckardt à distance finie =",nops(Eckardt));
print(Eckardt);
Nombre de points de Eckardt à distance finie = 7
VIII. Détermination des triplets de droites parallèles ( points de Eckardt à l'infini)
>
#tot:=[op(solz),op(solx0),op(solx)]:
tot:= droites: #### si on travaille sur l'ordre pris dans fichier de données
pv:=proc(a,b) [a[2]*b[3]-a[3]*b[2],a[3]*b[1]-a[1]*b[3],a[1]*b[2]-a[2]*b[1]] end:
dir:=NULL:
for i from 1 to nops(tot) do
dir:=dir,expand(subs(x=0,y=0,z=0,tot[i])-subs(x=1,y=1,z=1,tot[i]))
od:
dir:=[dir]:
dr_para:=NULL:
n:=nops(dir):
for i from 1 to n-1 do
cpt:=1: tamp:=i :
for j from i+1 to n do
if pv(dir[i],dir[j])=[0,0,0] then cpt:=cpt+1: tamp:=tamp,j: fi:
od:
if cpt=3 then dr_para:=dr_para,[dir[i],{tamp}] fi;
od:
dr_para:=[dr_para]:
printf("%a %s",nops(dr_para),"triplets de droites parallèles ou points de Eckardt à l'infini");
seq(print(dr_para[i]),i=1..nops(dr_para));
3 triplets de droites parallèles ou points de Eckardt à l'infini