d_Helmholtz_3D accuracy

[Ahmad_Falahatpisheh]

What is the accuracy of d_Helmholtz_3D? I would like to know what the residual error is after solving the equation.

Thanks.

[Alexander_K_Intel2 · Accepted solution]

Hi,

d_Helmholtz_3D is the direct solver of matrix correspond of 7-point grid Helmholtz equation. So its provide accuracy based on floating operations.

With best regards,

Alexander Kalinkin

[Ahmad_Falahatpisheh]

since it is double precision, does it mean the accuracy is 1E-16?

[Alexander_K_Intel2]

Not equal but about it.

With best regards,

Alexander Kalinkin

[Ahmad_Falahatpisheh]

Dear Alexander,

I checked the accuracy of d_Helmholtz_3D and it was much much larger than 1E-16. Since the solver uses a standard seven-point discretization, I verified the accuracy by the following code. (I have a uniform mesh for my problem.)
[cpp]for ( k=1; kGetTuple( i + j *NX + k *NX*NY); phi_im1jk = Phi->GetTuple( (i-1) + j *NX + k *NX*NY); phi_ip1jk = Phi->GetTuple( (i+1) + j *NX + k *NX*NY); phi_ijm1k = Phi->GetTuple( i + (j-1)*NX + k *NX*NY); phi_ijp1k = Phi->GetTuple( i + (j+1)*NX + k *NX*NY); phi_ijkm1 = Phi->GetTuple( i + j *NX + (k-1)*NX*NY); phi_ijkp1 = Phi->GetTuple( i + j *NX + (k+1)*NX*NY); X_i = X->GetTuple( i ); X_im1 = X->GetTuple( i-1 ); Y_j = Y->GetTuple( j ); Y_jm1 = Y->GetTuple( j-1 ); Z_k = Z->GetTuple( k ); Z_km1 = Z->GetTuple( k-1 ); rhs= f->GetTuple( i + j*NX+ k*NX*NY); res = ( phi_im1jk + phi_ip1jk - 2*phi_ijk)/pow( X_i - X_im1, 2 ) + ( phi_ijm1k + phi_ijp1k - 2*phi_ijk)/pow( Y_j - Y_jm1, 2 ) + ( phi_ijkm1 + phi_ijkp1 - 2*phi_ijk)/pow( Z_k - Z_km1, 2 ) + rhs; } } }[/cpp]
When I print res, the residual is about 1E-1. Is there something that I have to be careful when using the function? I need to have an accuracy about 1E-16. Please advise.

Thanks,
Ahmad

[Alexander_K_Intel2]

Hi Ahmad,

To verify it I need to have full example with rhs and boundary condition. Could you provide this example to me by e'mail or by private answer?

With best regards,

Alexander Kalinkin

[Ahmad_Falahatpisheh]

Alexander,

I found a bug in my code which misled me to the see large residuals. I fixed it and the error now is about 1E-15.

Thanks,
Ahmad

[Alexander_K_Intel2]

Hi Ahmad,

Nice to hear it, feel free to ask any question about PL in particular and MKL in general.

With best regards,

Alexander Kalinkin

[Ahmad_Falahatpisheh]

Hi Alexander,

I have another question. The solver is for uniform mesh. Does this mean that it has to have dx=dy=dz? Or we can have dx!=dy!=dz (constant dx, dy, dz everywhere in the domain)?

Thanks,
Ahmad

[Alexander_K_Intel2]

Hi Ahmad,

The uniform mesh mean that all mesh steps are equals in one direction, but mesh sizes for different dimension could be differ. For example hx=0.2, hy=0.5, hz=0.1.

With best regards,

Alexander Kalinkin

[Ahmad_Falahatpisheh]

I am writing a journal paper in which I have used d_Helmholtz. Regarding the 7-point grid Helmholtz equation, can I have the name of the method by which the system is solved?

Thanks,
Ahmad

[Alexander_K_Intel2]

Hi Ahmad,

The main information could be in paper prepared by us a several years ago so feel free to use it.

With best regards,

Alexander Kalinkin

[Ahmad_Falahatpisheh]

Thanks. You helped me a lot.
Best,
Ahmad

[Ahmad_Falahatpisheh]

Hi Alexander,

I didn't find the method by which the library solves the system. Is it gradient bi-conjugate, multigrid, overrelaxation, or Fourier? I would appreciate it.

Thanks,
Ahmad

[Alexander_K_Intel2]

Hi Ahmad,

Poisson library based on Fourier decomposition for elliptic problems with separable variables.

With best regards,

Alexander Kalinkin