Intersection de deux droites
Bonjour
Dans le cadre de mon travail je suis confronté à un problème. Je déroule une méthode que je décris ci-dessous qui est très coûteuse en temps de calcul et j'aurais voulu la simplifier.
Mon but : j'ai deux droites définies ci-dessous.
a+A*t=x
b+B*t=y
c+C*t=z
et
a'+A'*t'=x
b'+B'*t'=y
c'+C'*t'=z
Je résous pour obtenir t et t' :
a+A*t=a'+A'*t'
b+B*t=b'+B'*t'
t =-(a*B'-b*A'-a'*B'+b'*A')/(A*B'-A'*B)
t'=-(a*B-b*A-a'*B+b'*A)/(A*B'-A'*B)
Je vérifie que si c+C*t=c'+C'*t' alors intersection des droites.
Maintenant c'est là que ça se complique, comme je sais que les droites ne s'intersectent pas je fais tourner l'une d'elle suivant un vecteur défini d'un certain angle jusqu'à intersection.
Et c'est cet angle que je cherche ! Je le trouve en faisant une boucle mais je cherche à le calculer directement pour gagner en temps de calcul.
Voici comment je procède :
Ma droite "fixe":
a+A*t=x
b+B*t=y
c+C*t=z
La "mobile"
a'+A'*t'=x
b'+B'*t'=y
c'+C'*t'=z
Mon vecteur suivant lequel la droite tourne {vx, vy, vz} et d'angle k
Ma matrice de rotation :
[vx^2*(-cos(k)+1)+cos(k), -vz*sin(k)+vx*vy*(-cos(k)+1), vy*sin(k)+vx*vz*(-cos(k)+1)]
[vz*sin(k)+vx*vy*(-cos(k)+1), vy^2*(-cos(k)+1)+cos(k), -vx*sin(k)+vy*vz*(-cos(k)+1)]
[-vy*sin(k)+vx*vz*(-cos(k)+1), vx*sin(k)+vy*vz*(-cos(k)+1), vz^2*(-cos(k)+1)+cos(k)]
La nouvelle position de la droite "mobile"
a'=a'*(vx^2*(-cos(k)+1)+cos(k))+b'*(-vz*sin(k)+vx*vy*(-cos(k)+1))+c'*(vy*sin(k)+vx*vz*(-cos(k)+1))
b'=a'*(vz*sin(k)+vx*vy*(-cos(k)+1))+b'*(vy^2*(-cos(k)+1)+cos(k))+c'*(-vx*sin(k)+vy*vz*(-cos(k)+1))
c'=a'*(-vy*sin(k)+vx*vz*(-cos(k)+1))+b'*(vx*sin(k)+vy*vz*(-cos(k)+1))+c'*(vz^2*(-cos(k)+1)+cos(k))
A'=A'*(vx^2*(-cos(k)+1)+cos(k))+B'*(-vz*sin(k)+vx*vy*(-cos(k)+1))+C'*(vy*sin(k)+vx*vz*(-cos(k)+1))
B'=A'*(vz*sin(k)+vx*vy*(-cos(k)+1))+B'*(vy^2*(-cos(k)+1)+cos(k))+C'*(-vx*sin(k)+vy*vz*(-cos(k)+1))
C'=A'*(-vy*sin(k)+vx*vz*(-cos(k)+1))+B'*(vx*sin(k)+vy*vz*(-cos(k)+1))+C'*(vz^2*(-cos(k)+1)+cos(k))
t =-(a*B'-b*A'-a'*B'+b'*A')/(A*B'-A'*B)
t'=-(a*B-b*A-a'*B+b'*A)/(A*B'-A'*B)
Dans ma boucle je fais évoluer k pour vérifier : c+C*t=c'+C'*t'
Mais ce que je voudrais c'est connaître k par calcul direct ! Pour cela je cherche à résoudre : c+C*t=c'+C'*t' pour en sortir k
J'ai lancé le calcul avec différent logiciel de calcul formel mais je n'ai toujours pas de réponse après 96h de traitement.
Quand j'essaye de résoudre à la main j'arrive très vite à des longueurs de formule hallucinante !
Est-ce que je m'y prends mal ? Y a-t-il une méthode plus simple ?
Merci pour votre aide !
Dans le cadre de mon travail je suis confronté à un problème. Je déroule une méthode que je décris ci-dessous qui est très coûteuse en temps de calcul et j'aurais voulu la simplifier.
Mon but : j'ai deux droites définies ci-dessous.
a+A*t=x
b+B*t=y
c+C*t=z
et
a'+A'*t'=x
b'+B'*t'=y
c'+C'*t'=z
Je résous pour obtenir t et t' :
a+A*t=a'+A'*t'
b+B*t=b'+B'*t'
t =-(a*B'-b*A'-a'*B'+b'*A')/(A*B'-A'*B)
t'=-(a*B-b*A-a'*B+b'*A)/(A*B'-A'*B)
Je vérifie que si c+C*t=c'+C'*t' alors intersection des droites.
Maintenant c'est là que ça se complique, comme je sais que les droites ne s'intersectent pas je fais tourner l'une d'elle suivant un vecteur défini d'un certain angle jusqu'à intersection.
Et c'est cet angle que je cherche ! Je le trouve en faisant une boucle mais je cherche à le calculer directement pour gagner en temps de calcul.
Voici comment je procède :
Ma droite "fixe":
a+A*t=x
b+B*t=y
c+C*t=z
La "mobile"
a'+A'*t'=x
b'+B'*t'=y
c'+C'*t'=z
Mon vecteur suivant lequel la droite tourne {vx, vy, vz} et d'angle k
Ma matrice de rotation :
[vx^2*(-cos(k)+1)+cos(k), -vz*sin(k)+vx*vy*(-cos(k)+1), vy*sin(k)+vx*vz*(-cos(k)+1)]
[vz*sin(k)+vx*vy*(-cos(k)+1), vy^2*(-cos(k)+1)+cos(k), -vx*sin(k)+vy*vz*(-cos(k)+1)]
[-vy*sin(k)+vx*vz*(-cos(k)+1), vx*sin(k)+vy*vz*(-cos(k)+1), vz^2*(-cos(k)+1)+cos(k)]
La nouvelle position de la droite "mobile"
a'=a'*(vx^2*(-cos(k)+1)+cos(k))+b'*(-vz*sin(k)+vx*vy*(-cos(k)+1))+c'*(vy*sin(k)+vx*vz*(-cos(k)+1))
b'=a'*(vz*sin(k)+vx*vy*(-cos(k)+1))+b'*(vy^2*(-cos(k)+1)+cos(k))+c'*(-vx*sin(k)+vy*vz*(-cos(k)+1))
c'=a'*(-vy*sin(k)+vx*vz*(-cos(k)+1))+b'*(vx*sin(k)+vy*vz*(-cos(k)+1))+c'*(vz^2*(-cos(k)+1)+cos(k))
A'=A'*(vx^2*(-cos(k)+1)+cos(k))+B'*(-vz*sin(k)+vx*vy*(-cos(k)+1))+C'*(vy*sin(k)+vx*vz*(-cos(k)+1))
B'=A'*(vz*sin(k)+vx*vy*(-cos(k)+1))+B'*(vy^2*(-cos(k)+1)+cos(k))+C'*(-vx*sin(k)+vy*vz*(-cos(k)+1))
C'=A'*(-vy*sin(k)+vx*vz*(-cos(k)+1))+B'*(vx*sin(k)+vy*vz*(-cos(k)+1))+C'*(vz^2*(-cos(k)+1)+cos(k))
t =-(a*B'-b*A'-a'*B'+b'*A')/(A*B'-A'*B)
t'=-(a*B-b*A-a'*B+b'*A)/(A*B'-A'*B)
Dans ma boucle je fais évoluer k pour vérifier : c+C*t=c'+C'*t'
Mais ce que je voudrais c'est connaître k par calcul direct ! Pour cela je cherche à résoudre : c+C*t=c'+C'*t' pour en sortir k
J'ai lancé le calcul avec différent logiciel de calcul formel mais je n'ai toujours pas de réponse après 96h de traitement.
Quand j'essaye de résoudre à la main j'arrive très vite à des longueurs de formule hallucinante !
Est-ce que je m'y prends mal ? Y a-t-il une méthode plus simple ?
Merci pour votre aide !
Connectez-vous ou Inscrivez-vous pour répondre.
Réponses
Il y a une infinité de rotations qui font s'intersecter les deux droites. Pour chacune, le point d'intersection est différent. Donc il te suffit de choisir un point sur D(la droite fixe), et une rotation (*) qui amène un point de (D') sur celui-ci.
Avant de calculer, penser géométriquement des questions géométriques abrège le travail.
Cordialement.
(*) je n'ai pas trop compris quelles rotations tu envisageais, ton exposé est trop "calculs sans explication". Et la matrice donne une rotation vectorielle, pas euclidienne.
Schématiquement ce que je voudrais faire.
J'impose la rotation d'une droite suivant un axe quelconque mais défini, passant par l'origine pour qu'elle intersecte avec une autre droite. Et ce que je voudrais connaître c'est l'angle parcouru autour de cette axe depuis son point de départ jusqu’à l'intersection. Et du coup il n'y a plus une infinité de solutions mais une seule, à n demi-tours près.
Cordialement.
> il n'y a plus une infinité de solutions mais une seule, à n demi-tours près.
Non. La droite mobile en tournant autour de l'axe fixe engendre un hyperboloïde de révolution à une nappe. L'autre droite coupe en général cet hyperboloïde en deux points, ou ne le coupe pas (c'est ce qui se passe pour une droite et une quadrique). Et il n'y a aucune raison pour que les deux solutions diffèrent par un demi-tour.
Cordialement.
Il me semble pourtant avoir compris l'explication que tu donnes du problème, et ça me semble coller avec les calculs que tu fais.
Une remarque sur la matrice de rotation que tu écris. Je devine à partir de cette écriture que le vecteur $v$ est de norme 1 : $v_x^2+v_y^2+v_z^2=1$. C'est bien ça ? Ensuite je vois au moins une erreur dans l'écriture de ta matrice : des signes faux sur la diagonale. En effet la somme des coefficients diagonaux d'une matrice de rotation (sa trace) est égale à $1+2\cos(k)$, si $k$ est l'angle de rotation.
Pour mener le calcul, je ne m'y prendrais pas comme ça.
Tu as une paramétrisation $M(t)$ de la droite fixe et une paramétrisation $M'(t')$ de la droite que tu fais tourner autour de l'axe passant par l'origine $O$ dirigé par le vecteur $\vec v$. Soit $P(t)$ la projection ortogonale de $M(t)$ sur cet axe, $P'(t')$ celle de $M'(t')$ (ça se calcule sans difficulté). Tu peux alors trouver $t'=\varphi(t)$ tel que $P(t)=P'(\varphi(t))$. Il suffit d'exprimer que le produit scalaire de $\overrightarrow{OM}(t)$ avec $\vec v$ est égal à celui de $\overrightarrow{OM'}(t')$ avec $\vec v$, et on en tire $t'$ comme fonction affine $\varphi$ de $t$.
Il reste ensuite à exprimer que la distance de $P(t)$ à $M(t)$ est égale à la distance de $P(t)$ à $M'(\varphi(t))$, ce qui fait une équation de degré deux en $t$, qui quand les dieux sont favorables admet deux racines réelles $t_1,t_2$. Reste ensuite à déterminer les angles $\widehat{M'(\varphi(t_i))P(t_i)M(t_i)}$, ce qui se fait sans problème .
Tout à fait, le vecteur v est de norme 1. Je viens de vérifier la matrice et je ne trouve pas l’erreur de signe. Sinon j'ai suivi votre raisonnement jusqu'au point d'exprimer la distance de P(t) à M(t). Je n'ai pas saisi de quelle manière l'exprimer, pouvez-vous m'éclairer, merci beaucoup.
Cordialement.
Au temps pour moi, ta formule est correcte.
Je bloque sur la résolution de l'équation de degré deux en t.
Voici mes étapes de calcul:
$M(t)$
$a+A.t$
$b+B.t$
$c+C.t$
$M'(t')$
$a'+A'.t'$
$b'+B'.t'$
$c'+C'.t'$
Axe de rotation[size=large]
$vx$
$vy$
$vz$
[/size]
Voici les projections orthogonales de $M(t)$ et $M'(t')$ sur l'axe de rotation, je les ai exprimés sous la même forme que $M(t)$ et $M'(t')$
par contre les vecteurs ne sont pas normalisés, est-ce un problème?
$P(t)$ ,la projection orthogonale de $M(t)$ sur l'axe de rotation
[size=large]
$vx.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})+(a-vx.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})).t$
$vy.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})+(b-vy.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})).t$
$vz.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})+(c-vz.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})).t$
[/size]
$P'(t')$ ,la projection orthogonale de $M'(t')$ sur l'axe de rotation
[size=large]
$vx.(\frac{(vx.a'+vy.b'+vz.c')}{\sqrt{vx^2+vy^2+vz^2}})+(a'-vx.(\frac{(vx.a'+vy.b'+vz.c')}{\sqrt{vx^2+vy^2+vz^2}})).t'$
$vy.(\frac{(vx.a'+vy.b'+vz.c')}{\sqrt{vx^2+vy^2+vz^2}})+(b'-vy.(\frac{(vx.a'+vy.b'+vz.c')}{\sqrt{vx^2+vy^2+vz^2}})).t'$
$vz.(\frac{(vx.a'+vy.b'+vz.c')}{\sqrt{vx^2+vy^2+vz^2}})+(c'-vz.(\frac{(vx.a'+vy.b'+vz.c')}{\sqrt{vx^2+vy^2+vz^2}})).t'$
[/size]
Vous demandiez d'exprimer $t'=\varphi(t)$
$(a+A.t).vx+(b+B.t).vy+(c+C.t).vz=(a'+A'.t').vx+(b'+B'.t').vy+(c'+C'.t').vz$
[size=x-large]
$t'=\frac{((A.vx+B.vy+C.vz).t+a.vx+b.vy+c.vz-a'.vx-b'.vy-c'.vz)}{(A'.vx+B'.vy+C'.vz)}$
[/size]
Distance $P(t)$ à $M(t)$
[size=large]
$((vx.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})+(a-vx.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})).t)-(a+A.t))^2$
$+((vy.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})+(b-vy.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})).t)-(b+B.t))^2$
$+((vz.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})+(c-vz.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})).t)-(c+C.t))^2$
[/size]
Distance $P(t)$ à $M'(\varphi(t))$
[size=petite]
$((vx.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})+(a-vx.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})).t)-(a'+A'.(\frac{((A.vx+B.vy+C.vz).t+a.vx+b.vy+c.vz-a'.vx-b'.vy-c'.vz)}{(A'.vx+B'.vy+C'.vz)})))^2$
$+((vy.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})+(b-vy.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})).t)-(b'+B'.(\frac{((A.vx+B.vy+C.vz).t+a.vx+b.vy+c.vz-a'.vx-b'.vy-c'.vz)}{(A'.vx+B'.vy+C'.vz)})))^2$
$+((vz.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})+(c-vz.(\frac{(vx.a+vy.b+vz.c)}{\sqrt{vx^2+vy^2+vz^2}})).t)-(c'+C'.(\frac{((A.vx+B.vy+C.vz).t+a.vx+b.vy+c.vz-a'.vx-b'.vy-c'.vz)}{(A'.vx+B'.vy+C'.vz)})))^2$
[/size]
C'est ici mon problème.
je me retrouve avec une égalité à résoudre trop longue. Suis-je sur la bonne voie?
En effet, j'ai la tête dans le guidon depuis quelques jours, je n'ai pas vu cette évidence ! Je suis peut-être un peu vieux jeu mais j'aime prendre la main sur ce que je fais. Je me sers du système de calcul formel pour simplifier si possible les expressions. Mais malheureusement j'arrive encore à des longueurs de formule stratosphérique pour $t1$ et $t2$, est-ce normal ?
Cordialement.
Il vaut mieux comprendre ce qu'on fait, et faire faire les calculs pénibles par un logiciel de ton choix.
PS. Pourquoi cherches-tu des formules, qui seront forcément horribles ? Au vu de tes précédentes interventions, j'ai l'impression que tu serais plutôt intéressé par une procédure qui prend en entrée les paramétrisation des deux droites et le vecteur v et qui recrache les angles de rotation qui amènent la deuxième droite à rencontrer la première - ou indique que ce n'est pas possible.
Au contraire, je cherche a simplifier au maximum les étapes de calcul. Je me suis servi d'un support 3D pour comprendre votre démarche en plus d'une feuille de calcul Excel pour valider les calculs. Et je me trompe peut-être mais je n'obtiens pas le résultat théorique avec votre méthode. Voici mes calculs:
$M(t)$
$100-0.173648177667.t$
$50$
$80-0.984807753012.t$
$M'(t')$
$80$
$50+t'$
$70$
L'axe de rotation
$-0.117755093199$
$-0.076023486696$
$0.990128359101$
Je ressors des angles de:
$16.4651472754$
$70.7452175857$
Pour un théorique suivant ma méthode de calcul confirmé en 3D de:
$15.362279395489$
J'ai fait un code à partir de ce que j'ai expliqué. Je pourrais le commenter si besoin est. C'est du SageMath :
Et voila l'éxécution sur ton exemple :
donne en sortie
angle : 71.1583011665874 degrés
angle : -15.3622793954808 degrés
Un autre exemple :
angle : 17.8931910329837 degrés
angle : 163.804833385816 degrés
CPU times: user 19.1 ms, sys: 279 µs, total: 19.4 ms
Wall time: 18.4 ms
Votre code change-t-il profondément si l'axe de rotation ne passe pas par l'origine mais est défini par un point et une direction ?