I am writing to inquire if it is possible somehow to extract directly the finite differences from the Intel MKL Helmholtz solver. I am using the solver for the Poisson equation, in order to solve for the stream function Psi(x,y,z) of a fluid. In this case I am not interested in the direct value of Psi itself, but rather the gradients of Psi. Seeing as the capability exists to use Neumann boundary conditions, and it is stated in the solver description that a five-point finite difference solver is used, than I am assuming this data may be somehow stored already within the solver during execution time.
I could possibly also use the same solver, however now solving for the gradient variables -df/dx instead of -f, but this will probably end up being exactly as expensive as calculating the finite difference gradients which I am currently doing, possibly even more so due to the more expensive calculation of the gradients on the boundary. Hence why I wonder if there is a way simply to extract directly the gradient. Perhaps I am also overseeing a simple way to use the Poisson solver directly?
Any help would be appreciated.
I don't think this is possible and definitely it cannot be done through the current interface. Based on the information you gave, I'd rather consider other existing ways of approximating your problem so that the gradient of the stream function will be a proper (maybe, the main) unknown of the problem rather than a side-effect of computing the stream function. Of course, you can employ a post-processing scheme to compute the gradients from the stream function but these tend to lack robustness/are unstable.
Yes I already gave some thought to formulating the problem such that the gradient is the unknown of the problem, this is in fact what I meant in the second paragraph of my post. This may however end up being more expensive than the problem itself was originally. Thank you for your quick response anyway.
Alternatively, I have read elsewhere that the curl of a given field can also be directly extracted by performing an operation on the field in fourier space such as described here:
I have seen some rough examples for Matlab, so I thought I might check here if somebody has any experience carrying this out with the MKL Fast Poisson solver. Of course it would require preparation of my input arrays, however I cannot find any solid examples online. Any help would be appreciated!