How to find the maximum of this oscillating function?

Question

Given $n \in\mathbb{N}$, $a_1, \dots, a_n \in\mathbb{C}$, $k \in \mathbb{R}$ and $x_1, \dots x_n \in \mathbb{R}^3$, let $\Phi : \mathbb{R}^3 \to \Bbb C$ be defined by

$$\Phi (y) := \sum_{i=1}^{n} {a_i \frac{e^{j k \| x_i - y \|_2}}{\| x_i - y \|_2}}$$

where $j := \sqrt{-1}$ is the imaginary unit. Note the singularities at $y \in \{ x_1, \dots x_n \}$.

Can somebody please point me to a more elegant method than brute forcing (just evaluating the expression at millions of points) to find $y$ for which $|\Phi|$ has the global maximum and is not near a singularity (distance > $k$)?

Answer 1

Follows a MATHEMATICA script which implements a maximization procedure to obtain the desired maximum. The nature of the problem, prone to many relative maxima, indicates to choose an evolutionary algorithm to perform the search. The charge points are given in X the coefficients $a_k$ are given as a[1],...,a[n], $k$ is represented by lambda. The objective function is handled as $\|\Phi(y)\|^2$ and the unknown $y$ as Y. Attached a plot showing in blue the X charge points and in red the solution for Y.

SeedRandom[2]
Clear[Phi]
n = 30;
r = 0.06;
parms = Join[Table[a[k] -> 1, {k, 1, n}], {lambda -> 5}];
X = RandomReal[2 {-1, 1}, {n, 3}];
Y = {x, y, z};
dots = Table[Sqrt[(X[[k]] - Y).(X[[k]] - Y)], {k, 1, n}];
Phi[j_] := Sum[a[k] Exp[j lambda dots[[k]]]/dots[[k]], {k, 1, n}];
objr = (Phi[I] + Phi[-I])/2 // ComplexExpand;
obji = (Phi[I] - Phi[-I])/(I 2) // ComplexExpand;
obj = objr^2 + obji^2;
restrs = Table[Sqrt[(X[[k]] - Y).(X[[k]] - Y)] - 2 Pi/lambda, {k, 1, n}];
restrs0 = Thread[restrs > 0];
case = Join[{obj}, restrs0] /. parms;
sol = Quiet@NMaximize[case, Y, Method -> "DifferentialEvolution"]

(******* To show the fitting to the restrictions *******)
restrs /. parms /. sol[[2]]
restrs0 = Thread[restrs > 0] /. parms /. sol[[2]]
(*******************************************************)

gr1 = Table[Graphics3D[{Blue, Sphere[X[[k]], r]}], {k, 1, n}];
gr2 = Graphics3D[{Red, Sphere[Y, r]} /. sol[[2]]];
Show[gr1, gr2]