Discussion :

[cpalperou]

Bonjour,
je dois trouver numériquement les solutions de l'équation différentielle suivante:

d(alpha)/dz = (1/L(z))*(1-E(z)*exp(alpha*fac(z))

J'ai bien sûr les vecteurs L(z), E(z), fac(z) et z lui même
Je n'ai pas l'ombre du spectre du fantôme d'une idée! A part utiliser la fonction ode45

Merci de votre aide

[cpalperou]

Voici ce que j'ai essayé:
dans un 1er fichier, intitulé allol.m :

Code :

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

1
2
3
 
% extrait de ce code: 
[T,alp]=ode45(@(alpha) myfunode(alpha,Ldif,zdif,fac_exp),zdif,1e-10);

et dans un 2ème intitulé myfunode.m:

Code :

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

1
2
3
4
 
function yfun = myfunode(alpha,L,Ec,fac)
 
    yfun=(1./L)*(1-Ec.*exp(alpha,fac));

et, quand je lance dans matlab le prog allol.m, voici l'erreur qu'il me retourne:

Code :

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

1
2
3
4
5
6
7
8
9
10
11
12
 
 Error using @(alpha)myfunode(alpha,Ldif,zdif,fac_exp)
Too many input arguments.
 
Error in odearguments (line 88)
f0 = feval(ode,t0,y0,args{:});   % ODE15I sets args{1} to yp0.
 
Error in ode45 (line 114)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
 
Error in allol (line 113)
    [T,alp]=ode45(@(alpha) myfunode(alpha,Ldif,zdif,fac_exp),zdif,1e-10);

[Invité]

Bonjour,

La fonction odefun dont le handle est passé en premier argument de la fonction ode doit prendre deux paramètres, or toi tu ne lui en fournis qu'un.

Code :

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

1
2
% extrait de ce code: 
[T,alp]=ode45(@(z,alpha) myfunode(alpha,Ldif,zdif,fac_exp),zdif,1e-10);

[cpalperou]

Merci pour votre réponse winjerome.
J'ai effectué la correction que vous apportez. J'ai aussi trouvé une petite bétise dans ce que j'avais écrit: dans l'appel à ode45 je n'utilisais pas les bons vecteurs (zdif au lieu de Ec):

Code :

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

1
2
% extrait de ce code: 
[T,alp]=ode45(@(z,alpha) myfunode(alpha,Ldif,E,fac_exp),zdif,1e-10);

et dans la fonction myfunode, j'ai dû changer quelques * par .*:

Code :

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

1
2
function yfun = myfunode(alpha,L,Ec,fac)
    yfun=(1./L).*(1-Ec.*exp(alpha.*fac));

et j'ai encore une erreur:

Code :

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

1
2
3
4
5
6
7
8
Error using odearguments (line 91)
@(Z,ALPHA)MYFUNODE(ALPHA,LDIF,E,FAC_EXP) must return a column vector.
 
Error in ode45 (line 114)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
 
Error in allol (line 114)
    [T,alp]=ode45(@(z,alpha)myfunode(alpha,Ldif,E,fac_exp),zdif,1e-10);

[Invité]

Comme dit dans le message d'erreur : must return a column vector. : ta fonction odefun doit retourner un vecteur colonne.

Par contre en voyant ton premier message, L, E et fac dépendent de z, il faudrait alors je pense utiliser l'astuce de l'exemple 3 de le documentation.

[cpalperou]

Re,
après quelques messages d'erreurs m'indiquant que les inputs de la fonction myfunode doivent être de même dimensions, je me suis assuré de ce point.
Un message m'a également indiqué que la condition initiale devait être un vecteur de la dimension de mes inputs (E, Ldif,...). Ce qu j'ai fait.
Voici mes codes:
allol.m:

Code :

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

1
2
3
% extrait de ce code: 
I0 = zeros(1,numel(E))+1e-10;
[T,alp]=ode45(@(z,alpha) myfunode(alpha,Ldif,E,fac_exp),zdif,I0);

et myfunode.m:

Code :

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

1
2
3
4
5
function yfun = myfunode(alpha,L,Ec,fac)
L=L';
Ec=Ec';
fac=fac';
yfun=(1./L).*(1-Ec.*exp(alpha.*fac));

J'ai lancé le programme il y a maintenant 30 minutes et il est toujours en train d'itérer.
Je pense l'interrompre et chercher une autre solution.
Winjerome, oui, en effet L,E et fac dépendent de z. Mon équa diff se rapproche de celle de l'exemple 3 de l'aide. La mienne est du genre:
y'+f(t)*exp[g(t)*y]=h(t)

Je ne vois pas comment m'inspirer de cet exemple.

[cpalperou]

Lorsque j'indique une condition initiale unique:
I0=1e-10
voici l'erreur obtenue:

Code :

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

1
2
3
4
5
6
7
8
9
10
Error using odearguments (line 93)
@(Z,ALPHA)MYFUNODE(ALPHA,LDIF,E,FAC_EXP) returns a vector of length 45000, but the length of initial
conditions vector is 1. The vector returned by @(Z,ALPHA)MYFUNODE(ALPHA,LDIF,E,FAC_EXP) and the initial
conditions vector must have the same number of elements.
 
Error in ode45 (line 114)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
 
Error in RefAngle4 (line 114)
    [T,alp]=ode45(@(z,alpha)myfunode(alpha,Ldif,E,fac_exp),zdif,I0);

[Invité]

Tes variables L,E et fac sont des tableaux de valeurs qui correspondent à certaines valeurs de z.

Par exemple (non représentatif de tes valeurs à toi, c'est juste un exemple)

• pour la valeur z = 1, tu prendras les valeurs d'indice 10 : L(10), E(10), fac(10)
• pour la valeur z = 3, tu prendras les valeurs d'indice 15 : L(15), E(15), fac(15)

Maintenant la fonction ode45, comment fonctionne t-elle ?
Itérativement, elle va donner des valeurs (z, alpha) à ta fonction myfunode afin d'obtenir en retour la valeur yfun qui correspond à ces deux entrées selon ta formule de départ : (1/L(z))*(1-E(z)*exp(alpha*fac(z)).
Mais, L, E et fac n'étant que des tableaux et non des fonctions, il est fort probable que la fonction ode45 te fournisse des valeurs de z pour lesquelles tu n'aies pas de correspondance exacte, mais des valeurs approchantes. C'est là que la fonction interp1 de l'exemple intervient : elle va effectuer le calcul des valeurs L(z), E(z) et fac(z).

Et tu n'as qu'un valeur de condition initiale à donner ici : alpha(z₀)

[cpalperou]

Merci WinJérome,
je n'y arrive toujours pas (j'ai beaucoup de difficultés avec les fonctions handle qu'utilisent ode45 ou fmincon par exemple).
Si je comprends bien:
j'ai mes vecteurs L, E et fac qui sont défini pour certaines valeurs de z. Ces valeurs de z sont regroupés dans un vecteur que j'ai nommé zdif. C'est sur cet intervalle que je veux le résultat de la résolution de mon équa dif.
Ainsi, à zdif(1) correspond L(1), E(1), fac(1) ....
je dois trouver la solution de l'ode d(alpha)/dz = (1/L(z))*(1-E(z)*exp(alpha*fac(z)) et le résultat, alpha(z) doit être donné sur l'intervalle zdif.
C'est à dire que je veux qu'à zdif(1) corresponde alpha(1), à zdif(2) corresponde alpha(2), ...
J'ai modifié mes codes selon winbjérome ainsi:
allol.m:

Code :

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

1
2
I0=1e-10;
    [T,alp]=ode45(@(z,alpha) myfunode(alpha,Ldif,E,fac_exp,zdif),zdif,I0);

myfunode.m:

Code :

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

1
2
3
4
5
6
function yfun = myfunode(alpha,L,Ec,fac,zdif)
zi = [zdif(1):-.001:zdif(end)]';
L = interp1(zdif,L,zi);
Ec = interp1(zdif,Ec,zi);
fac = interp1(zdif,Ec,zi);
yfun=(1./L).*(1-Ec.*exp(alpha.*fac));

et j'obtiens l'erreur suivante:

Code :

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

1
2
3
4
5
6
7
8
9
10
Error using odearguments (line 93)
@(Z,ALPHA)MYFUNODE(ALPHA,LDIF,E,FAC_EXP) returns a vector of length 45000, but the length of initial
conditions vector is 1. The vector returned by @(Z,ALPHA)MYFUNODE(ALPHA,LDIF,E,FAC_EXP) and the initial
conditions vector must have the same number of elements.
 
Error in ode45 (line 114)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
 
Error in RefAngle4 (line 114)
    [T,alp]=ode45(@(z,alpha)myfunode(alpha,Ldif,E,fac_exp),zdif,I0);

comprends pas!!

[cpalperou]

Bien, j'ai modifié mon programme myfunode.m ainsi :

Code :

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

1
2
3
4
5
function yfun = myfunode(z,alpha,L,Ec,fac,zdif)
L = interp1(zdif,L,z);
Ec = interp1(zdif,Ec,z);
fac = interp1(zdif,Ec,z);
yfun=(1./L).*(1-Ec.*exp(alpha.*fac));

et :

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
Index exceeds matrix dimensions.
 
Error in interp1 (line 124)
        V = V(idx,:);
 
Error in myfunode (line 7)
fac = interp1(zdif,Ec,z);
 
Error in @(z,alpha)myfunode(z,alpha,Ldif,E,fac_exp,zdif)
 
 
Error in odearguments (line 88)
f0 = feval(ode,t0,y0,args{:});   % ODE15I sets args{1} to yp0.
 
Error in ode45 (line 114)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
 
Error in RefAngle4 (line 115)
    [T,alp]=ode45(@(z,alpha) myfunode(z,alpha,Ldif,E,fac_exp,zdif),Tspan,I0);

?????

[cpalperou]

Bien, j'ai trouvé une erreu dans mon programme myfunode.m:

Code :

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

1
2
3
4
5
function yfun = myfunode(z,alpha,L,Ec,fac,zdif)
L = interp1(zdif,L,z);
Ec = interp1(zdif,Ec,z);
fac = interp1(zdif,Ec,z);
yfun=(1./L).*(1-Ec.*exp(alpha.*fac));

En gras, c'est pas Ec que je dois indiquer mais bien sûr fac!
et le programme ne renvoie plus d'erreur: cool .....

mais les alpha obtenus sont négatifs! ils devraient être positifs!!

[Invité]

Que retourne

Code :

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

whos L Ec fac zdif

?