Equa diff sous octave(matlab)

[inviteac2c37b6]

Bonjour, j'ai un soucis : je tente de résoudre des équations différentielles sous octave.
En fait je doit modéliser un choc parfaitement élastique de deux sphères grâce a une force de rappel d'un ressort. Donc si d(la distance entre les deux sphères : d=sqrt(x²+y²)) est inférieur à 2*r(le rayon des sphères) on à mx''=-Kx
et my''=-Ky. Et si la d>2r x''=0 et y''=0 (il faut que je fasse ça pour les deux boules donc ça vas donner un truck du genre xa,ya,x"a,y"a,xb,yb.....).

Donc la j'ai déjà un problème pour modéliser quand il n'y a pas de choc (d>2r), je pence que ces parce que je ne maitrise pas vraiment la résolution numérique (la formulation surtout) que ca marche pas.
Donc la pour l'instant je me suis simplifier la tache j'ai considérer que une sphère et une dimension : j'ai x''(t)=0, je définis x1(t)=x'(t) donc x1'(t)=0 . Je définis le vecteur y(t)=(x(t),(x1(t)) Donc dy=y'(t)=(x1(t),x1'(t))

Ainsi voici ce que donne ma fonction

function dy=fct(y,t)
dy=zeros(size(y))
dy(1)=y(2)
dy(2)=0

Et mon programme :

clear all
close all
clc
t=[0:1:10]
CI=[2,0]
Y=lsode(@fct,CI,t)
plot(t,Y,'*','r')

J'ai beau mettre ce que je veux dans CI (condition initiales), rien ne change, le plot n'affiche rien.
Donc voici mes problèmes :
-Pourquoi je n'obtiens pas le résultat attendus, une droite constante pour tout t.
-A quoi correspond les CI ? Cet a dire CI(1)= x(o) et CI(2)= x'(0) ou l'inverse ?
-Pour la suite du programme je ne vois pas comment je pourrais faire rentrer tout ca dans mon fichier fonction avec ce formalisme (boula a, boule b, x et y)

Merci d'avance pour votre aide.
-----

[inviteac2c37b6]

Bon la je viens de corriger un soucis avec le plot, j'ai enlever les options. Donc la ces bon j'obtiens ma trajectoire rectiligne, mais j'ai toujour ce problème pour modéliser le reste (dimension x et y, sphère a et b) dans la même fonction.

[inviteac2c37b6]

J'ai un peu avancer et la je suis vraiment confronter a ce problème de formalisme, je vous montre ce que j'ai écris

La fonction

function [dx,dy]=fct(x,y,t)
dx=zeros(size(x))
dy=zeros(size(y))
dx(1)=x(2)
dx(2)=0
dy(1)=y(2)
dy(2)=0

Le programme :
clear all
close all
clc
t=[0:1:10]
CI=[2,0,0,0]
[X,Y]=lsode(@fct,CI,t)
plot(X,Y,'o')

J’essaie ici de généraliser mon programme a deux dimensions, pour ensuite tracer y en fonction de x, dans le cas de mon problème je souhaite une trajectoire rectiligne et parallèle a l'axe des abscisses
Donc la dans mon esprit les conditions initiales correspondrais a cela CI(1) = x(0), CI(2) = x'(0), CI(3) =y(0) et CI(4) = y'(0), dans la réalité je ne sait pas, et je ne peux pas le vérifier car le programme ne marche pas.
Voici le message d’erreur :
error: invalid vector index = 2
error: evaluating assignment expression near line 6, column 6
error: called from `fct' in file `D:\3.0.5_gcc-4.3.0\programmes\fct.m'
error: lsode: evaluation of user-supplied function failed
error: lsode: inconsistent sizes for state and derivative vectors
error: near line 6 of file `D:\3.0.5_gcc-4.3.0\programmes\fctmain.m'
octave.exe:11:\3.0.5_gcc-4.3.0\programmes
>

Sachant que ligne 6 colone 6 correspond au ls ode.

Je conte généraliser ce programme aux deux sphères(a et b) et aux deux cas(collision pas collision a l'aide d'une boucle if dans la fonction), donc si j'arrive a formaliser pour les deux dimensions, je pourrais formaliser pour les deux dimension + les deux sphères, et je pence que la boucle if ne poseras pas de problèmes.

Merci de me dire ce qui ne vas pas, et comment corriger.

[inviteac2c37b6]

Alors pour palier a ce problème de formalisme (mettre les quatre équadif du second ordre(xa,xb,ya,yb) donc les 8 équadif de premier ordre dans le même fichier) j'ai séparer le tout en quatre fichier fonction + 1 fichier "main".
Un nouveau problème ce pose je vous explique :
Voici le fichier "main"


clear all
close all
clc

global d
CIxa=[5,1]
CIya=[5,1]
CIxb=[2,0]
CIyb=[2,0]
for t=0:0.1:100

Xa=lsode(@fctxa,CIxa,t)
Ya=lsode(@fctya,CIya,t)

Xb=lsode(@fctxb,CIxb,t)
Yb=lsode(@fctyb,CIyb,t)

d=sqrt(((Xa.^2)+(Ya.^2))-((Xb.^2)+(Yb.^2)))
end
plot(Xa,Ya,'o',Xb,Yb,'+')


Et voici l'un des quatre fichier fonction, ils sont similaires :


function dxa=fctxa(xa,t)
global d
dxa=zeros(size(xa))
if d>2

dxa(1)=xa(2)
dxa(2)=0
else
dxa(1)=xa(2)
dxa(2)=-xa(1)
end


Grâce a la boucle for il effectue à la suite les calculs ainsi : xa=, ya=, xb=, yb= ,d= et il recommence. Ces exactement ce que je recherche. Mais il y a un problème : je vous montre deux boucles vous allez comprendre

Xa =

5 1

Ya =

5 1

Xb =

2 0

Yb =

2 0

d =

6.4807 1.4142

Xa =

5 1

Ya =

5 1

Xb =

2 0

Yb =

2 0

d =

6.4807 1.4142

Et ces pareil pour tout t, il n'y a aucune variation de mes fonctions.
En fait je pence que les lsode ne fonctionne pas, il ne fait pas appel a mes fonctions dans la boucle for et ce contente de n'utiliser que les conditions initiales.
D’où viens le problème ? Sinon je pence toujours que si je pouvais tout compacter dans deux fichier : un fichier "main" et un fichier fonction je n'aurais pas de problème mais je ne sait pas comment faire.
Autre problème il me sort la variable d (distance entre les deux sphères) en deux dimensions, alors que cela devrais être en une dimension.

[A voir en vidéo sur Futura]