Я пытаюсь решить уравнение Пуассона в прямоугольной области с помощью быстрого преобразования Фурье-косинуса с библиотекой FFTW3.
Граничное условие четырех сторон является нулевым граничным условием Неймана.
d^u/dx^2+d^2u/dy^2=-f
f=cos(x)+cos(y)
, домен [-pi pi,-pi pi]
u=f
,f
удовлетворяют граничному условию Неймана.U=F/(lamda(i)+lamda(j))
lamda(k)=2*(1-cos(i/(n0-1)))
U
инвертируется, чтобы стать маленьким u
,Тем не менее, я все еще получаю неправильный результат по сравнению с точным решением.
Не могли бы вы мне помочь? Огромное спасибо.
Вот мой код
#include <stdio.h>
#include <math.h>
#include <fftw3.h>
#include <iostream>
using namespace std;
int main() {
int n0=64;
int n1=64;
int N=n0*n1;
double pi = 3.14159265359;
double L = 2.*pi;
double dx = L/(n0);
double *in1=new double [N];
double *in2=new double [N];
double *out1=new double [N];
double *out2=new double [N];
double *X=new double [N];
double *Y=new double [N];fftw_plan p, q;
int i,j;
p = fftw_plan_r2r_2d(n0,n1, in1, out1, FFTW_REDFT00, FFTW_REDFT00, FFTW_ESTIMATE);
q = fftw_plan_r2r_2d(n0,n1, in2, out2, FFTW_REDFT00, FFTW_REDFT00, FFTW_ESTIMATE);
for(i = 0;i <n0;i++){
X[i] =-pi+i*dx ;
for(j = 0;j<n1;j++){
Y[j] = -pi+j*dx ;
in1[i*n0 + j]= cos(X[i]) + cos(Y[j]) ; // row major ordering
}
}
fftw_execute(p);
double *lamda=new double [n0];
for (i=0;i<n0;i++){
lamda[i]=cos(pi*i/(n0-1));
}
for ( i= 0; i< n0; i++){ // f = g / ( lamda(ii)+lamda(jj) )
for( j = 0; j < n1; j++){
double fact=0;
in2[i*n0 + j]=0;
fact=2*(2-lamda[i]-lamda[j]);
if(fact!=0){
in2[i*n0 + j] = out1[i*n0 + j]/fact;
}else{
in2[i*n0 + j] =0;
}
}
}
fftw_execute(q);
double erl1 = 0.;
for ( i = 0; i < n0; i++) {
for( j = 0; j < n1; j++){
cout<< i <<" "<< j<<" "<< cos(X[i])+cos(Y[j])<<" "<< dx*dx*out2[i*n0+j]/(2.*(double)(n0-1))/(2.0*(double)(n1-1)) <<" "<< endl;
}
}
fftw_destroy_plan(p); fftw_destroy_plan(q); fftw_cleanup();
delete [] in1;
delete [] in2;
delete [] out1;
delete [] out2;
delete [] X;
delete [] Y;
return 0;
}
Задача ещё не решена.