Discussion :

[Rhun sans coeur]

Bonjour
Je suis en peine avec la library fftw et fft de Matlab.
Je doit reprendre dans une modélisation l'équation w = laplacien(phi). Je connais w, et dois remonter à phi. Je souhaite travailler en spectral pour avoir des conditions de bords périodiques.
Le problème est en 2D.
L'idée est la suivante :
FFT(w)=FFT(laplacien(phi))=FFT( -( (kx)²+(ky)² )*(phi) )
Et donc :
phi =iFFT ( - FFT (w) / ( (kx)²+(ky)² ) ).
Dans le principe, ca marche, avec kx et ky ls vecteurs d'ondes associés aux deux dimentions x et y.

Sauf que. Je suis complètement incapable de réaliser la division de FFT(w) par les bons vecteurs d'ondes, car je n'arrive pas à réaliser la surface de kx et ky correctement. Et ceci à cause de l'ordre dans lequel la fonction fft place les coefficients d'amplitudes et à cause de la parité de ma grille.

Ma routine de test ressemble à cela pour l'instant : elle est en 1D car avant de faire gros, il faut faire petit.

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
 
%test fft1D
lx=1000;
LX=2^nextpow2(lx);
dx=1;
%calcul d'une perturbation initiale m,n en A.
i=1:lx;
A(i)=sin(1*i*pi*2/lx);
subplot(1,6,1)
plot(A);
 
%calcul du laplacien de la perturbation A.
DELA=diff(A,2);
subplot(1,6,2)
plot(DELA);
 
%calcul de la fft du laplacien de la perturbation A.
FTDELA=fft(DELA,LX)/lx;
subplot(1,6,3)
plot(abs(FTDELA),'+');
 
%calcul de la fft de A a partir de la fft du laplacien de A.
k = 1/(dx*2)*linspace(0,1,LX/2+1);
m= 1/(dx*2)*linspace(-1,0,LX/2+1);
 
for i=2:LX/2+1;
FTA(i)=-FTDELA(i)/(k(i)*k(i)*(4*pi*pi));
end;
for  i=LX/2+2:LX
	FTA(i)= -FTDELA(i)/(m(i-LX/2)*m(i-LX/2)*(4*pi*pi));
end
subplot(1,6,4)
plot(abs(FTA),'+');
 
%vérification de la tete de ce que l'on devrait avoir en fft de A
FTA2=fft(A,LX)/lx;
subplot(1,6,5)
plot(abs(FTA2),'+');
 
%calcul de la perturbation A à nouveau.
B=ifft(FTA,LX)*lx;
subplot(1,6,6)
plot(abs(B));

C'est réalisé sous QTOctave, je suis débutant avec fourier et apparement : la fonction finale ne ressemble pas à la fonction initiale.

[phryte]

Bonjour.
Je mettrais :

Code :

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

plot(real(B));axis([0 lx -1 1]);grid

[FR119492]

Salut!

Je doit reprendre dans une modélisation l'équation w = laplacien(phi). Je connais w, et dois remonter à phi. Je souhaite travailler en spectral pour avoir des conditions de bords périodiques.

Pauvre Joseph Fourier! Il doit se retournerdans sa tombe. Pourquoi utiliser la FFT qui est justement faite pour les cas où les conditions ne sont pas périodiques. Tu as le choix entre trois excellentes méthodes:

• les différences finies;
• les éléments finis;
• les séries de Fourier (sous leur forme originale).

Jean-Marc Blanc

[Rhun sans coeur]

@phryte
En fait : ici les plots sont juste pour voir où j'en suis. La routine n'a pas pour but de tracer quoique ce soit, mais de recalculer un champ.

@Jean-Marc Blanc
Alors , il y a peut être des choses que je n'ai pas comprises. La première méthode que j'ai utilisé pour résoudre cela est basée sur une récurrence. Cela doit être la méthode des différences finies.
w(i,j)= -4*phi(i,j)+phi(i+1,j)+phi(i-1,j)+phi(i,j+1)+phi(i,j-1)
Cela correspond à une dérivée forward et backward pour faire le laplacien. Si on connait phi au bord de la grille, on peut recalculer phi partout. Mais connaître phi tout autour veut dire connaitre les conditions de bord. Une manière de les connaître, c'est de les fixer. Et j'arrive a des conditions de bords fixes.
Or mon système est périodique, mais ne peut pas supporter ces conditions de bords fixes.
C'est pour ça que j'ai pense à passer en spectral. Et donc à la fft. Où me plante-je donc ?