Discussion :

[hengstepha]

Bonjour à tous, j'aimerais programmer une équation dynamique en C++ qui est:
f(0,x)=[K-exp(x)]+ pour x<u;
f(i,x)=0 pour x>=u;
f(i+1,x)=pn*f(i,x+sigma*sqrt(delta))+(1-pn)*f(i,x-sigma*sqrt(delta));
voila mon code:

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
  double partiepositive (double k,double e) // renvoie la partie positive de k-e
{
  if (k>e) return k-e;
  else 0;
}
 
double phi(double u, double K,double  T,double  sigma,double N,double r,double i,double x)//equiation dynamique
{
    double delta=T/N;
    double res;
    double mu=r-sigma*sigma/2;
    double pn=0.5+mu*sqrt(delta)/(2*sigma);
    res=partiepositive(K,exp(x));
     double w1=x+sigma*sqrt(delta);
   double w2=x-sigma*sqrt(delta);
 
    if((i==0)&&(x<u))
      return res;
    else
      {
    if((i==0)&&(x>=u))
      {return 0;}
    else
      {
        if(x>=u)
          return 0;
        else
          {   
        return pn*phi(u, K, T, sigma,N,r, i-1,w1)+(1-pn)*phi(u, K, T, sigma,N,r, i-1,w2);
          }
      }
      }
}

mon problème est que l'algorithme est très lent!!! à partir de i=27..
avez-vous des améliorations pour y remédier??
merci d'avance.

[Ulmo]

Tiens, un financier...

Pense à la balise code...

Et stocke plutôt les valeurs de f(.,.) dans un tableau. Actuellement tu calcules 2 fois f(i-2, x) au lieu d'une, par exemple (un fois pour x-sigma*sqrt(delta)+sigma*sqrt(delta) et une fois pour x+sigma*sqrt(delta)-sigma*sqrt(delta)).

[souviron34]

hou que c'est pas beau....

Bon c'est OK avec le langage, mais c'est pas esthéthique, pour du code de maths...

Déclare tes variables, et fait tes calculs.

Bon, ça c'est une chose (mais qui n'engage que moi)

Par contre, déjà tu calcules je sais combien de fois sqrt(delta)... Stockes le dans une variable.

Même chose pour sigma*sqrt(delta)

Au lieu de diviser par 2, mutiplies par 0.5

Ton K = exp, de même que toutes tes initialisations (et donc les calculs), tu les calcules à tous les coups, sans tester par rapport à u, alors que ton algo le dis

voilà pour l'instant

[hengstepha]

Tiens, un financier...

Pense à la balise code...

Et stocke plutôt les valeurs de f(.,.) dans un tableau. Actuellement tu calcules 2 fois f(i-2, x) au lieu d'une, par exemple (un fois pour x-sigma*sqrt(delta)+sigma*sqrt(delta) et une fois pour x+sigma*sqrt(delta)-sigma*sqrt(delta)).

tu veux que je stocke tous les valeurs de f(,) ou bien ce que les f(,) que je calcule plusieurs fois.

[hengstepha]

Tiens, un financier...

Pense à la balise code...

Et stocke plutôt les valeurs de f(.,.) dans un tableau. Actuellement tu calcules 2 fois f(i-2, x) au lieu d'une, par exemple (un fois pour x-sigma*sqrt(delta)+sigma*sqrt(delta) et une fois pour x+sigma*sqrt(delta)-sigma*sqrt(delta)).

je reformule...
je stocke tous les f(,) (je pense que c'est faisable mais faudrait un LONG tableau,ou des pointeurs de pointeurs.............. de pointeurs) ou seulement les f(,) que je calcule plusieurs fois( là c'est trop chaud pour un pauvre petit financier...)

ps: pas mal l'antithèse?

[Ulmo]

Pour calculer f(n, x), tu as besoin d'une valeur de f(n, .) (en x).
Donc de 2 valeurs de f(n-1, .).
Donc de 3 valeurs de f(n-2, .) puisque les incréments sont les mêmes vers le haut et vers le bas.
...
Donc de N+1 valeurs de f(0, .) que tu connais explicitement.

Il suffit donc d'un tableau de taille n+1, mais il faut faire attention à l'ordre de remplissage pour n'écraser que les valeurs de l'étape précédente dont on n'a plus besoin.

Au passage, tu as une formule explicite avec la formule du binôme, mais ce n'est peut-être pas le but...

[hengstepha]

Pour calculer f(n, x), tu as besoin d'une valeur de f(n, .) (en x).
Donc de 2 valeurs de f(n-1, .).
Donc de 3 valeurs de f(n-2, .) puisque les incréments sont les mêmes vers le haut et vers le bas.
...
Donc de N+1 valeurs de f(0, .) que tu connais explicitement.

Il suffit donc d'un tableau de taille n+1, mais il faut faire attention à l'ordre de remplissage pour n'écraser que les valeurs de l'étape précédente dont on n'a plus besoin.

Au passage, tu as une formule explicite avec la formule du binôme, mais ce n'est peut-être pas le but...

a ok en faite on oublie la forme itérative ( on utilise pas l'équation dynamique). et on essai de trouver la forme générale?? si j'ai bien compris c'est bien cela?

[Ulmo]

Quand on peut c'est mieux, on gagne du temps. Ici il suffit de compter combien de fois il y aura f(0, x+k*h) dans f(n, x).

[salseropom]

Code :

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

1
2
3
4
5
  double partiepositive (double k,double e) // renvoie la partie positive de k-e
{
  if (k>e) return k-e;
  else 0;
}

Salut, ta fonction partiepositive marche, mais elle renvoie des warning. Je te propose la chose suivante

Code :

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

1
2
3
4
5
 
double partiepositive (double k,double e) // renvoie la partie positive de k-e
{
  return (k>e) ? k-e : 0.;
}

et peut-être que le mieux serait de faire une macro

Code :

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

1
2
 
#define MAX(a,b) ((a>b) ? (a) : (b))

et ensuite de faire partout, à la place de partiepositive(a,b) de mettre

Code :

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

1
2
 
MAX(k-e,0.);

je ne pense pas que ça ira bcp plus vite mais au moins tu auras enlevé qq warning lors de la compilation et je pense que le code sera plus lisible...

[hengstepha]

Quand on peut c'est mieux, on gagne du temps. Ici il suffit de compter combien de fois il y aura f(0, x+k*h) dans f(n, x).

ok
voila mais changement : j'ai sépare le cas pair et impair mais je ne trouve pas de bon résulats...


double fact(double n)
{
if (n<=1)
return 1;
else
return n*fact(n-1);
}

double binom(double n, double p)
{
return fact(n)/(fact(p)*(fact(n-p)));
}


double phiimpaire(double u, double K,double T,double sigma,int N,double r,double x)
{
int k=-1;
double * px=new double[N+1];
double delta=T/N;
double mu=r-sigma*sigma/2;
double pn=0.5+mu*sqrt(delta)/(2*sigma);
double j=sigma*sqrt(delta);
double res;

for (int i=0; i<floor(N/2)+1;i++)
{
px[i]=binom(N,i)*pow(pn,N-i)*pow(1-pn,i)*phiint(u,K,x+(N-2*i)*j);
}
for (int i=(int)(floor(N/2)+1);i<N+1;i++)
{
k++;
px[i]=binom(N,i)*pow(pn,N-i)*pow(1-pn,i)*phiint(u,K,x-(1+2*k)*j);
}
for(int i=0; i<N+1;i++)
{
res+=px[i];
}
delete [] px;
return res;

}

double phipaire(double u, double K,double T,double sigma,int N,double r,double x)
{
int k=0;
double *px=new double[N+1];
double delta=T/N;
double mu=r-sigma*sigma/2;
double pn=0.5+mu*sqrt(delta)/(2*sigma);
double j=sigma*sqrt(delta);
double res;

for (int i=0; i<N/2;i++)
{
px[i]=binom(N,i)*pow(pn,N-i)*pow(1-pn,i)*phiint(u,K,x+(N-2*i)*j);
}
px[N/2]=binom(N,N/2)*pow(pn,N-N/2)*pow(1-pn,N/2)*phiint(u,K,x);
for (int i=N/2+1;i<N+1;i++)
{
k++;
px[i]=binom(N,i)*pow(pn,N-i)*pow(1-pn,i)*phiint(u,K,x-(2*k)*j);
}
for(int i=0; i<N+1;i++)
{
res+=px[i];
}
delete [] px;
return res;

}
si tu peux reperer l'erreur sa serait cool

[hengstepha]

Quand on peut c'est mieux, on gagne du temps. Ici il suffit de compter combien de fois il y aura f(0, x+k*h) dans f(n, x).

j'ai remarqué que ma fonction binomiale marcche pour des petites valeurs maisne marchait pas pour les grandes par exemple n=400, c'est normal?

[diogene]

La première chose à faire est d'éliminer les calculs redondants (sqrt) ou à priori pas forcéments utiles (exp). Par exemple :

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
double calc(double u, double K,double  SigmaDel,double pn, double i,double x)//equiation dynamique
{
    if(x>=u) return 0;
    if(i==0)
    {
      double e = exp(x);
      return K>e? K-e:0;
    }
    return pn*calc(u, K,SigmaDel,pn,i-1,x+SigmaDel)+(1-pn)*calc(u, K,SigmaDel,pn, i-1,x-SigmaDel);
}
double phi(double u, double K,double  T,double  sigma,double N,double r,double i,double x)//equiation dynamique
{
    double sq = sqrt(T/N);
    double pn=0.5+(r-sigma*sigma/2)*sq/(2*sigma);
    return calc(u,K,sigma*sq,pn,i,x);
}

Ensuite, il faut abandonner la structure recursive pour une structure itérative

[hengstepha]

La première chose à faire est d'éliminer les calculs redondants (sqrt) ou à priori pas forcéments utiles (exp). Par exemple :

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
double calc(double u, double K,double  SigmaDel,double pn, double i,double x)//equiation dynamique
{
    if(x>=u) return 0;
    if(i==0)
    {
      double e = exp(x);
      return K>e? K-e:0;
    }
    return pn*calc(u, K,SigmaDel,pn,i-1,x+SigmaDel)+(1-pn)*calc(u, K,SigmaDel,pn, i-1,x-SigmaDel);
}
double phi(double u, double K,double  T,double  sigma,double N,double r,double i,double x)//equiation dynamique
{
    double sq = sqrt(T/N);
    double pn=0.5+(r-sigma*sigma/2)*sq/(2*sigma);
    return calc(u,K,sigma*sq,pn,i,x);
}

Ensuite, il faut abandonner la structure recursive pour une structure itérative

salut a toi
j'ai essayé ce que tu as donné mais pour un N >27 c'est très long....

[hengstepha]

La première chose à faire est d'éliminer les calculs redondants (sqrt) ou à priori pas forcéments utiles (exp). Par exemple :

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
double calc(double u, double K,double  SigmaDel,double pn, double i,double x)//equiation dynamique
{
    if(x>=u) return 0;
    if(i==0)
    {
      double e = exp(x);
      return K>e? K-e:0;
    }
    return pn*calc(u, K,SigmaDel,pn,i-1,x+SigmaDel)+(1-pn)*calc(u, K,SigmaDel,pn, i-1,x-SigmaDel);
}
double phi(double u, double K,double  T,double  sigma,double N,double r,double i,double x)//equiation dynamique
{
    double sq = sqrt(T/N);
    double pn=0.5+(r-sigma*sigma/2)*sq/(2*sigma);
    return calc(u,K,sigma*sq,pn,i,x);
}

Ensuite, il faut abandonner la structure recursive pour une structure itérative

merci du conseil mais le probleme pour passer à la forme itérative c'est que ici elle dépend de calc(u, K,SigmaDel,pn,i-1,x+SigmaDel) et de
calc(u, K,SigmaDel,pn,i-1,x-SigmaDel) pas très simple..

[Ulmo]

j'ai remarqué que ma fonction binomiale marcche pour des petites valeurs maisne marchait pas pour les grandes par exemple n=400, c'est normal?

Tu as des des erreurs de débordement : 400! ne tient pas dans un double. Tu ne peux aller que jusqu'à 170.
De toute façon, à partir de N = 20 tu peux remplacer efficacement l'arbre par une variable gaussienne. Tu as alors un calcul presque explicite (la formule de Black & Scholes).

[souviron34]

Tu as des des erreurs de débordement : 400! ne tient pas dans un double. Tu ne peux aller que jusqu'à 170.
De toute façon, à partir de N = 20 tu peux remplacer efficacement l'arbre par une variable gaussienne. Tu as alors un calcul presque explicite (la formule de Black & Scholes).

Complètement faux.. Il y a tout un thread sur le sujet :

http://www.developpez.net/forums/sho...d.php?t=278583

Ne confondez pas nombre de chiffres et valeur....

Sur des sytèmes 32 bits (double = 64 bits), on peut aller jusqu'à :

factorielle(1.97969 e+308)

(voir post bas de page 2..)

[Ulmo]

Complètement faux.. Il y a tout un thread sur le sujet :

http://www.developpez.net/forums/sho...d.php?t=278583

Ne confondez pas nombre de chiffres et valeur....

Sur des sytèmes 32 bits (double = 64 bits), on peut aller jusqu'à :

factorielle(1.97969 e+308)

(voir post bas de page 2..)

Plus précisemment jusqu'à :
factorielle(171) = 1.97969 e+308

Edit : D'ailleurs 171! ne tient pas dans un double, on s'arrête bien à 170!.
Et à 1754! avec des long double.

[souviron34]

milles excuses

j'avais affiché ceci en faisant une estimation rapide par dychotomie (pas en faisant le calcul complet)...

[diogene]

Si le sujet est toujours d'actualité, tu peux essayer ce code :

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
45
46
47
48
double valeur(double K, double x)
{
     double e = exp(x);
     return K>e? K-e:0;
}
// Les coefficients du binôme
long cnp(int n,int p)
{
 long C = 1;
 int i;
 if(p<0) return 0;
 for(i=1; i<=p; i++)C= C*(n-i+1)/i;
 return C ;
}
double calc2(double u, double K,double  SigmaDel,double pn, double n,double x)
{
  double del = (x-u)/SigmaDel/2;
  double expo = pow(pn,n);
  double fact = (1-pn)/pn;
  int q = n/2;
  int i;
  long coef;
  int m0;
  int i0;
  double result=0;
  if(x>=u) return 0;
  if(n==0)return valeur(K,x);
  i0 = ceil(n/2.0+del);
  i0 = i0<0 ? 0 : i0;
  for(i = i0; i <= n; i++)
  {
   m0 = ceil(n-i+2*del);
   m0 = m0<0? 0 : m0;
   coef =0;
   if(i<n-q  && m0<=i) coef = cnp(n,i)-cnp(n,m0-1);
   if(i>=n-q && m0<=n-i) coef = cnp(n,n-i)-cnp(n,m0-1);
   if(coef!=0)result += expo*valeur(K,x+(n-2*i)*SigmaDel)* coef;
   expo = expo*fact;
  }
  return result;
 
}
double phi(double u, double K,double  T,double  sigma,double N,double r,long i,double x)
{
    double sq = sqrt(T/N);
    double pn=0.5+(r-sigma*sigma/2)*sq/(2*sigma);
    return calc2(u,K,sigma*sq,pn,i,x);
}

D'après qq tests, il semble coller et être beaucoup plus rapide. On peut encore optimiser un peu ce code .