Discussion :

[chardonnette]

Bonjour à tous!
Je sollicite une aide bienveillante pour un problème de compréhension (ou syntaxe!) avec ODE45

Voici la fonction ....

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
function dydt=fxncubic(t,y)
dydt=zeros(size(y));
global x0 b0 x1 b1 T 
%paramètre temporel
E=0.1;
 
dydt(1)=(-y(1).^3-T.*y(1)+y(2))./E
dydt(2)=y(1)-x1
dydt=[dydt(1);dydt(2)]

....que j'utilise comme ceci :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

[t,y]=ode45('fxncubic',[0 10],[x0,b0]);

xo et bo sont définis dans mon script. Mon problème est le suivant :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
???  In an assignment  A(I) = B, the number of elements in B and
 I must be the same.
 
Error in ==> fxncubic at 8
dydt(1)=(-y(1).^3-T.*y(1)+y(2))./E

+ les autres messages d'erreurs liés à ODE
J'ai eu beau retourner le problème das tout les sens, je ne vois pas en quoi tout ces vecteurs n'ont pas la même taille!

Je vous remercie à l'avance pour votre aide (ou lecture ^^)

[Jerome Briot]

Quelles sont les dimensions de la variable T ?
Pourquoi avoir autant de variables globales ?

Et plus généralement, nous expliquer clairement l'objet de ton étude...
Quelle est l'équation que tu essaies de résoudre ?
A quel phénomène physique est-elle liée ?

[chardonnette]

Alors, l'objet d'étude est la modélisation d'un cycle cardique (diastole+systole) et de ses dysfonctionnements, un projet vraiment passionnant! (pour ceux que ça intéresse je peux toujours fournir les documents du sujet).Je n'en suis qu'au début et débute plus ou moins avec matlab, mais ça s'annonce assez simpa!

Ensuite, T est un scalaire, et pour les variables globales...a vrai dire je ne sais même pas pourquoi j'en ais mis autant également ^^, je met toutes les constantes la dedans pour ne pas être embêtée.

Sinon , le plus simple c'est que je mette le script en entier. Je ne voulais pas vous embêter avec des lignes inutile!

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
%choix du paramètre T la tension de la fibre musculaire
T=input('valeur de T : ');
x=-2:0.1:2;
b=-x.^3+T.*x;
plot(b,x); grid on; 
 
%hold on
 
 
%paramètre temporel
E=0.1;
 
% durée du cycle
t1=0
tf=10
 
%condition initiale sur la taille de la fibre musculaire systole
%x1=input('valeur de x1 : ')
x0=-1.3
x1=-1.3
%condition initiale sur le potentiel électrochimique diastole
b0=-x0^3+T*x0
 
%condition initiale sur le potentiel électrochimique systole
b1=-x1^3+T*x1
 
%conditions initiales sur leurs dérivées respectives
xp1=-(x1.^3-T*x1+b1)/E;
bp1=b1-x1;
 
%résolution avec ODE 45
[t,y]=ode45('fxncubic',[t1 tf],[x0,xp1,b0,bp1]);
 
%graphes
pause
plot(y(:,1),y(:,2))

[Jerome Briot]

Je ne suis pas un spécialiste de ce genre de problème mais essaye ceci :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
 
function test
 
global T x1 E
 
%choix du paramètre T la tension de la fibre musculaire
T = input('valeur de T : ');
x = -2:0.1:2;
b = -x.^3+T.*x;
 
figure
subplot(211)
plot(b,x); grid on; 
 
%hold on
 
 
%paramètre temporel
E = 0.1;
 
% durée du cycle
t1 = 0;
tf = 10;
 
%condition initiale sur la taille de la fibre musculaire systole
%x1=input('valeur de x1 : ')
x0 = -1.3;
x1 = -1.3;
%condition initiale sur le potentiel électrochimique diastole
b0 = -x0^3+T*x0;
 
%condition initiale sur le potentiel électrochimique systole
b1 = -x1^3+T*x1;
 
%conditions initiales sur leurs dérivées respectives
xp1 = -(x1.^3-T*x1+b1)/E;
bp1 = b1-x1;
 
%résolution avec ODE 45
[t,y]=ode45('fxncubic',[t1 tf],[x0,xp1,b0,bp1]);
 
%graphes
subplot(212)
plot(y(:,1),y(:,2))

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
function dydt=fxncubic(t,y)
global T x1 E
 
dydt=zeros(size(y));
% global x0 b0 x1 b1 T 
%paramètre temporel
% E=0.1;
 
dydt(1) = (-y(1).^3-T.*y(1)+y(2))./E;
dydt(2) = y(1)-x1;
% dydt=[dydt(1);dydt(2)]

ou beaucoup plus proprement sans les variables globales :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
function test
 
% global T x1 E
 
%choix du paramètre T la tension de la fibre musculaire
T = input('valeur de T : ');
x = -2:0.1:2;
b = -x.^3+T.*x;
 
figure
subplot(211)
plot(b,x); grid on; 
 
%hold on
 
 
%paramètre temporel
E = 0.1;
 
% durée du cycle
t1 = 0;
tf = 10;
 
%condition initiale sur la taille de la fibre musculaire systole
%x1=input('valeur de x1 : ')
x0 = -1.3;
x1 = -1.3;
%condition initiale sur le potentiel électrochimique diastole
b0 = -x0^3+T*x0;
 
%condition initiale sur le potentiel électrochimique systole
b1 = -x1^3+T*x1;
 
%conditions initiales sur leurs dérivées respectives
xp1 = -(x1.^3-T*x1+b1)/E;
bp1 = b1-x1;
 
%résolution avec ODE 45
[t,y]=ode45(@(t,y) fxncubic(t,y,T,x1,E),[t1 tf],[x0,xp1,b0,bp1]);
 
%graphes
subplot(212)
plot(y(:,1),y(:,2))

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
function dydt=fxncubic(t,y,T,x1,E)
% global T x1 E
 
dydt=zeros(size(y));
% global x0 b0 x1 b1 T 
%paramètre temporel
% E=0.1;
 
dydt(1) = (-y(1).^3-T.*y(1)+y(2))./E;
dydt(2) = y(1)-x1;
% dydt=[dydt(1);dydt(2)]