multiple calls to poisson solver produces odd results

[ringlenscl_msu_edu]

Greetings,

I'm have difficulties with the 3D Poisson solver. I'm trying to implement it in a particle-in-cell code that requires Poisson's equation to be solved repeatedly on a constant grid with changes to the boundary conditions (possibly) and charge densities on the internal mesh points.

According to the MKL guide, in the Sequence of Invoking PL Routines section, I should be able to initialize in the beginning, then loop the commit and solve routines, and clean up after the I'm done.

If I only perform one iteration, the correct solution is produced. However, if I try to loop, without changing any parameters, and look at the solution afterwards it is different if nx=ny=nz is not satisfied. If nx=ny=nz IS satisfied, then a given number of iterations produces the same answer (as expected). The code is given below. If anyone notices any problems, please let me know.

ThreeDCartPotSolver solver(64,64,64,1,1,1); //initialize solver class

solver.setY0BoundPotsConst(10); //set boundary points constant (all dirichlet)

solver.setY1BoundPotsConst(10);

solver.setX0BoundPotsConst(-10);

solver.setX1BoundPotsConst(-10);

for(int s=0; s<100; s++) //loop the solve routine (without altering boundary conditions).

{

solver.zeroSolutionVec(); //writes zeroes solution vector (f from MKL manual)

solver.commitSolver(); //commit

solver.calcPotential(); //solve

}

solver.writePotData("cartpot.dat"); //writes out f to data file

[Alexander_K_Intel2]

Hi,

As for me you have some problems with changed array ipar but I can't check it in your interfaces. Could you provide testcase that I could compiler and run on my side?

With best regards,

Alexander Kalinkin

[ringlenscl_msu_edu]

Hi Alexander,

I suspected that it might be something like that. I'll attach the files if you want to take a look at it. Thanks for your help.

Ryan

[Alexander_K_Intel2]

Also check dpar[0]-dpar[2]. These values have different meaning before and after commit step, before commit step they store length of correspondent domain edge but after they store mesh size. So, if you call in loop {commit,solve} check that these parameters set as length of domain edge before each commit step.

With best regards,

Alexander Kalinkin

[Alexander_K_Intel2 · Accepted solution]

Well, I've reproduced your issue on MKL example by calling commit and solve step in loop. This issue could be resolved by 2 ways:

First variant:

loop

{

set f = f_old;

dpar[0]=(bx-ax);

dpar[1]=(by-ay);

dpar[2]=(bz-az);

commit();

solve();

}

Second variant:

loop

{

set f = f_old;

commit();

solve();

ipar[0]=9; // more detailed inIntel Math Kernel Library ReferenceManual

}

If these ways doesn't resolve your issue please wrote it here and I will try to find solution of problem on your example

With best regards,

Alexander Kalinkin

[ringlenscl_msu_edu]

Hi Alexander,

Thanks for your help. I chose the first variant, and it worked perfectly. That probably would have taken me quite some time to find on my own. Thanks, again!

Ryan

[Alexander_K_Intel2]

Hi Ryan,

You are welcome, feel free to ask any question about Poisson library or other parts of MKL.

With best regards,

Alexander Kalinkin

[ringlenscl_msu_edu]

Hi Alexander,

I just wanted to update this thread with another issue I've found by looping over commit and solve. I noticed when I ran this loop thousands of times a memory leak developed. It disappeared when I changed the cycle to initialize->commit->solve->free. I wanted to avoid initializing at every iteration as I didn't know how much overhead was involved in this process.

The "Typical Order of Invoking PL Routines" figure in the MKL manual seems to indicate that initializing at the beginning of every iteration isn't required, so maybe this memory leak is a bug?

Thanks,

Ryan

[Alexander_K_Intel2]

Hi Rian,

You are right, if you choose first variant (dpar[0]=bx-ax and so on) the memory leak appeared because it is not recommended way of use Poisson library in loop, it is just a trick. Please choose second variant (set ipar[0]=9) then there is not memory leaks and result of PL is correct.

With best regards,

Alexander Kalinkin