Travail de méthodes numériques

**Hyozanryu** · 03/01/2008, 15h50

Bonjour ^^
Je sais pas trop si je poste dans la bonne section donc veuillez m'excuser si je me suis trompé

Je suis étudiant et j'ai un travail de méthodes numériques à faire mais je ne suis pas très doué

Le travail à faire est de créer un programme réalisant les calculs de la diffraction de Fraunhofer. Le programme doit être en C et les graphes tracés avec gnuplot, le tout sous linux. Un sujet similaire au mien se trouve ici (seule la première partie "Rappels" m'a été distribuée).
Le programme que j'ai fait est très basique:une fonction transmitance et une fonction calculant l'intégrale en chaque point (c'est une simple méthode des rectangles).
Seulement, si j'augmente le pas d'intégration, l'ordinateur met beaucoup trop de temps pour faire les calculs

Je me retrouve limité à un pas de 1 et le résultat est totalement inexploitable...
En fait l'intégrale à calculer est complexe mais le résultat qui m'intéresse pour le tracé du graphe (l'intensité) est réel et correspond à AA* (A étant l'intégrale écrite sur le sujet). J'aimerai savoir s'il n'y a pas une méthode un peu plus précise que celle que j'utilises ^^' On m'a conseillé d'utiliser la transformée de Fourier seulement ça a pas l'air simple à programmer...

Voila, j'espères que quelqu'un pourra m'aider un peu ^^ Merci

Edit: Oops, j'oublie de poster le 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
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
 
#include<math.h>
#include<stdio.h>
#include<stdlib.h>
 
#define N 10
 
int transmitance(double lx, double ly, double x, double y);
 
double sommeReelle(double lx, double ly, double pas, double X, double Y, double constante);
 
double sommeImmaginaire(double lx, double ly, double pas, double X, double Y, double constante);
 
int main()
{ 
  double X[N], Y[N], I[N][N];
  double x, y;
  double lx, ly;
  int i, j;
  double pas =1;
  double f, lamda;
  double constante = (2*M_PI)/(lamda*f);
  FILE *fp;
 
  for(i=0;i<=N;i++)   //discrétisation des axes X et Y
  {
     X[i]=-N/2+i*pas;
     Y[i]=X[i];
  }
 
  fp = fopen("intensite.data","w");
 
  printf("Entrez lx et ly:");
  scanf("%lf %lf",&lx,&ly);
 
  printf("Entrez lamda et f:");
  scanf("%lf %lf",&lamda,&f);
 
  for(i=0;i<=N;i++)
  {
      for(j=0;j<=N;j++)
      {
I[i][j]=pow(sommeReelle(lx,ly,pas,X[i],Y[j],constante),2)+pow(sommeImmaginaire(lx,ly,pas,X[i],Y[j],constante),2);
fprintf(fp,"%lf %lf %lf \n",X[i],Y[j],I[i][j]);
      }
  }
 
  fclose(fp);
 
  return 0;
}
 
int transmitance(double lx, double ly, double x, double y)
{
  int valeurTransmitance = 0;
  if( (x<=lx/2)&&(x>=-lx/2)&&(y<=ly/2)&&(y>=-ly/2) )
    valeurTransmitance = 1;
 
return valeurTransmitance;
}
 
double sommeReelle(double lx, double ly, double pas, double X, double Y, double constante)
{
  int x, y;
  double somme = 0;
  for(x=-N/2;x<=N/2;x+=pas)
  {
    for(y=-N/2;y<=N/2;y+=pas);
    {
         somme=somme+pas*pas*transmitance(lx,ly,x,y)*cos(constante*(x*X+y*Y));
    }
  }
 
return somme;
}
 
double sommeImmaginaire(double lx, double ly, double pas, double X, double Y, double constante)
{
  int x, y;
  double somme = 0;
  for(x=-N/2;x<=N/2;x+=pas)
  {
     for(y=-N/2;y<=N/2;y+=pas)
     {
         somme=somme+pas*pas*transmitance(lx,ly,x,y)*sin(constante*(x*X+y*Y));
     }
  }
 
return somme;
}

**diogene** · 03/01/2008, 16h15

Il faut d'abord corriger qquelques erreurs :

Code :

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

1
2
3
4
5
6
int main()
{ 
  double X[N], Y[N], I[N][N];
...
  double f, lamda;
  double constante = (2*M_PI)/(lamda*f);

lamda et f ne sont pas initialisés et ont n'importe quoi comme valeur

Code :

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

1
2
3
4
5
  for(i=0;i<=N;i++)
  {
      for(j=0;j<=N;j++)
      {
I[i][j]=pow(sommeReelle(lx,ly,pas,X[i],Y[j],constante),2)+pow(sommeImmaginaire(lx,ly,pas,X[i],Y[j],constante),2);

i et j ne peuvent prendre la valeur N, sinon on sort des tableaux.
On verra ensuite pour accélérer le calcul

**Hyozanryu** · 03/01/2008, 16h21

Ok, je vais essayer en modifiant ça, merci ^^ Par contre juste une question, si je sors du tableau, pourquoi la compilation marche?

**Melem** · 03/01/2008, 16h24

Par contre juste une question, si je sors du tableau, pourquoi la compilation marche?

Parce que le code est syntaxiquement correct.

**Hyozanryu** · 03/01/2008, 16h49

Voilà j'ai essayé ^^ Le résultat est déja mieux, j'arrive à deviner la forme du sinus cardinal ^^ Merci!

Par contre il y a toujours le problème avec le pas d'intégration, j'ai essayé de le diminuer un peu (0.5) mais ça rallonge condérablement les calculs...

**diogene** · 03/01/2008, 16h55

Peux-tu donner les paramètres typiques du calcul lx, ly, lamda et f ?

**Hyozanryu** · 03/01/2008, 17h07

Damn ça marche plus

Quelque soit les valeurs que j'entre pour lamda et f le programme renvoie une intensité constante. Pourtant ça a marché une fois...
Logiquement il faudrait entrer des valeurs plausibles. En TP les dimensions de la fente sont de l'ordre du mm, f de la dizaine de cm et lamda de l'ordre du nm...
Dans le cas d'une fente fine de longueur lx (lx>>ly), on est censé obtenir en 2D un sinus cardinal qui s'annule tous les lamda*f/a.

Edit: Bon j'ai retrouvé le cas où ça avait marché: ne pas initialiser lamda et f, et calculer la constante juste après (donc le scanf qui suit ne sert à rien). Donc les valeurs sont prises n'importe comment je suppose ^^' je vais tenter de faire un printf pour les relever.

Re edit: Après plusieurs test ça fonctionne pour lx et ly de l'ordre de l'unité (5 et 1 par exemple) et une constante très petite (10^-21 et 1). C'est pas vraiment conforme à la réalité mais bon...

**diogene** · 03/01/2008, 18h08

Avant de songer à accélérer le calcul, il faut résoudre ce problème des valeurs des paramètres. Apparemment, tu ne sais pas dans quel système d'unités tu travailles. Je pense que la première étape est d'écrire ton intégrale en normalisant les grandeurs, Lx et Ly par rapport à lambda, x par rapport à Lx, y par rapport à Ly et X et Y par rapport à f. Sinon, il est très difficile de déterminer les plages et le pas de calcul de l'intégrale.

**Hyozanryu** · 03/01/2008, 18h34

Désolé mais je comprend pas très bien ce que tu veux dire par normaliser les grandeurs

Ca doit être fait dans les fonctions qui calculent l'intégrale? Il faut imposer des limites pour le choix des valeurs?
En tout cas merci de m'aider, c'est super gentil ^^

**diogene** · 03/01/2008, 18h53

Envoyé par Hyozanryu

Désolé mais je comprend pas très bien ce que tu veux dire par normaliser les grandeurs

Ca doit être fait dans les fonctions qui calculent l'intégrale? Il faut imposer des limites pour le choix des valeurs?
En tout cas merci de m'aider, c'est super gentil ^^

Ce que je veux dire, c'est exprimer l'intégrale en fonction de x/Lx, y/Ly , X/f et Y/f qui sont alors des grandeurs sans dimensions. Par exemple, la plage d'intégration pour ta fente devient -1/2 à 1/2 quelle que soit les dimensions physiques. Le pas d'intégration devient 1/10, 1/100,... selon la précision demandée, indépendamment des dimension réelles.

**Hyozanryu** · 03/01/2008, 20h04

J'ai fait ça mais je pense pas que ce soit juste, ça marche pas.

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 sommeReelle(double lx, double ly, double pas, double X, double Y, double constante, double f)
{
  int x, y;
  double a=x/lx, b=y/ly;
  double somme = 0;
  for(a=-1/2;a<=1/2;a+=pas)
  {
    for(b=-1/2;b<=1/2;b+=pas);
    {
         somme=somme+pas*pas*transmitance(lx,ly,x,y)*cos(constante*(a*(X/f)+b*(Y/f)));
    }
  }
 
return somme;
}

**diogene** · 04/01/2008, 10h48

Envoyé par Hyozanryu

J'ai fait ça mais je pense pas que ce soit juste, ça marche pas.

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 sommeReelle(double lx, double ly, double pas, double X, double Y, double constante, double f)
{
  int x, y;
  double a=x/lx, b=y/ly;
  double somme = 0;
  for(a=-1/2;a<=1/2;a+=pas)
  {
    for(b=-1/2;b<=1/2;b+=pas);
    {
         somme=somme+pas*pas*transmitance(lx,ly,x,y)*cos(constante*(a*(X/f)+b*(Y/f)));
    }
  }
 
return somme;
}

Ca ne peut pas marcher, x et y ne sont pas initialisés. D'ailleurs, ils sont inutiles puisque dans cette plage de a et b transmitance vaut 1.

**Hyozanryu** · 16/01/2008, 20h06

Bonsoir ^^
Désolé de pas avoir répondu plus tôt, j'avais des exams à réviser ^^' Bon j'ai pas trop eu le temps d'avancer dans ce fichu programme mais j'y ai un peu réfléchi et je vois pas comment faire...Je suis repartit de l'expression analytique de l'intégrale et j'ai appliqué le changement de variable. Seulement le résultat est pas bon du tout. Le problême c'est que j'arrive pas à avoir une expression dont toutes les valeurs sont sans simensions...Voilà ce que j'ai changé:

Code :

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

1
2
3
4
5
6
7
8
9
10
11
12
13
14
double sommeReelle(double lx, double ly, double pas, double X, double Y, double constante,double f)
{
  int a, b;
  double somme = 0;
 
   for(a=-N/(2*lx);a<=N/(2*lx);a+=pas)
    {
      for(b=-N/(2*ly);b<=N/(2*ly);b+=pas)
	{
	  somme = somme+pas*pas*transmitance(lx,ly,a*lx,b*ly)*cos(constante*(a*lx*(X/f)+b*ly*(Y/f)));
	}
    }
  return somme;
}

**souviron34** · 16/01/2008, 23h27

deja un petit probleme :

pas est un double, et tu l'ajoutes a un entier dans l'increment des boucles......

Travail de méthodes numériques

C

Discussions similaires

Partager

Partager