We numerically solve for the electric potential on a box which obeys Ohm's law: \begin{equation}\nabla\cdot(\overset{\scriptscriptstyle\leftrightarrow}{\sigma} \nabla \Phi) = 0, \label{Eqn:OhmsLaw}\end{equation} where $\overset{\scriptscriptstyle\leftrightarrow}{\sigma}$ is the 2D conductivity tensor, and $\Phi$ are the electric potential and 2D current density respectively. The diagonal and off-diagonal elements of the conductivity tensor are given by the functions defined in the Mathematica code. They look like this:
The first plot is $\sigma_{xx}$ as a function of $r$ and second plot is $\sigma_{xy}$ and everything is radially symmetric. When $\sigma_{xx}$ is zero in that tiny interval the electric potential $\hat{z}\cdot((\nabla\sigma_{xy})\times(\nabla\Phi))$, contours of electric potential are parallel to contour of constant $\sigma_{xy}$ and since it is radially symmetric they should be circles.
Now here is the Mathematica code :
omegactimestau[B_] := 4.8 10^-10 B/(1.7 10^-29 3 10^10) 3 10^-4/10^8;
x1[r_] := -12.000000000000016` + 20.000000000000025` r;
(*Define conductivity tensor*)
sigmaxy2[x_, y_, B_] :=Piecewise[{{-1., Inequality[0.45, LessEqual, r, Less, 0.55]}, {0., r < 0},
{1., Inequality[0.6499999999999999, LessEqual, r, Less, 0.75]},
{-12.000000000000016 + 20.000000000000025*r, Inequality[0.55, LessEqual, r, Less,
0.6499999999999999]}, {18.849565896451544 - 32.1787625422147/Sqrt[1 + 4.*r^2], r >= 0.75}},
22.918285831105766 - 32.17876254221469/Sqrt[1 + 4.*r^2]]
sigmaxx2[x_, y_, B_] := 1/omegactimestau[B] (10^-3 + Piecewise[{
{Sqrt[x1[r]] (Sqrt[1 - x1[r]]), 0.6 <= r <= 0.65},
{-Sqrt[x1[r]] (Sqrt[-1 - x1[r]]), 0.55 <= r < 0.6},
{0, 0.45 <= r < 0.55},
{(37.69913179290309` - 64.3575250844294`/Sqrt[1 + 4.` r^2] - 2),
r > 0.75},
{Abs[
45.83657166221153` - 64.35752508442938`/Sqrt[1 + 4.` r^2]] -
2, 0 <= r < 0.45}
}] )/. {r -> Sqrt[(x - Lx/2)^2 + (y - Ly/2)^2]};
(*Solve Ohms law*)
gauss = 1 ;
tesla = 10^4 gauss ;
Lx = 15; Ly = 15;
rad = 0.297986 ;
rectangle = ImplicitRegion[0 <= x <= Lx && 0 <= y <= Ly, {x, y}];
frefine = Function[{vertices, area}, Block[{x1, y1}, {x1, y1} = Mean[vertices];
If[(x1 - Lx/2)^2 + (y1 - Ly/2)^2 < 100 rad^2 && area > 0.001, True,
False]]];
mesh = ToElementMesh[rectangle, MaxBoundaryCellMeasure -> 10^-2,
MeshRefinementFunction -> frefine];
<< NDSolve`FEM`
(*eqn[B_]=Inactive[Div][{{sigmaxx[x,y,B],sigmaxy[x,y,B]},{-sigmaxy[x,y,B],\
sigmaxx[x,y,B]}}.Inactive[Grad][phi[x,y],{x,y}],{x,y}];*)
reg = mesh;
brange = Range[0.5, 2, 0.5]*tesla;
sol1 = NDSolveValue[{Inactive[
Div][{{sigmaxx2[x, y, #], sigmaxy2[x, y, #]}, {-sigmaxy2[x, y, #],
sigmaxx2[x, y, #]}} . Inactive[Grad][phi[x, y], {x, y}], {x, y}] ==
0, phi[Lx, y] == 5, phi[0, y] == -5}, phi, {x, y} \[Element] reg] & /@
brange;
(*Plots*)
ContourPlot[
sol1[x + Lx/2, y + Ly/2] /. sol1 -> sol1[[2]], {x, -1, 1}, {y, -1, 1},
AspectRatio -> Automatic, ColorFunction -> ColorData["SunsetColors"],
ImageSize -> 650, PlotRangePadding -> None, FrameStyle -> Black,
BaseStyle -> FontSize -> 22, PerformanceGoal -> "Quality", Contours -> 100,
PlotRange -> All, PlotPoints -> 100]
Having obtained a solution for phi[x, y] ($\equiv \Phi(x,y)$) I want to calculate $\bar{\Phi}$ defined as
\begin{equation}
\bar{\Phi} ({r})= \int\,\frac{d^2 q}{2\pi^2}\tilde{\Phi}({q}) J_{0}\left(R_{c}|{q}|\right)e^{i {q}\cdot{r}},
\label{convphi}
\end{equation}
where $\tilde{\Phi}(q)$ is the Fourier transform of $\Phi(x,y)$, $J_{0}$ is the zero order Bessel function of first kind, $R_{c}$ is some constant parameter (for the code above it is of the order 0.1). Essentially the steps are first take inverse Fourier transform of solution obtained from NDSolve and multiply by the Bessel function and take Fourier transform of result.
In Mathematica there is no 2d Fourier transform equivalent for NFourierTransform.
Is there some alternative for the 2d case?
Edit : This question is edited after being deleted by me because I realized that my previous question had no good answer or that my question was bad. I apologize to the stack exchange and the answer below for bad practice.
ToElementMeshbefore you loadNDSolve`FEM`. Fixing that problem leads to another error "The PDE coefficient (big mess containingPiecewisefunctions) does not evaluate to a numeric matrix", which I'm not sure how to fix. $\endgroup$NIntegrateto evaluate$$\frac{1}{2\pi} \iint \Phi(x,y) e^{i(q_x x + q_y y)}\, dx\, dy?$$ $\endgroup$