I want to minimise S, which is given by,
$S = \sum_{i,j} \frac{1}{\sigma^2_{ij}}|V_{ij} - g_i g^*_j V^{model}_{ij}|^2$,
where $V_{ij}$ and $V^{model}_{ij}$ are the observed and model complex visibilities (this is a (M x 1) vector) both of which are known, $\sigma_{ij}$ is the error and $g_i$ and $g_j$ are the complex gains that I am minimising S for.
In my case I am only interested to solve for the phases of the complex gains assuming their amplitude is 1. I can therefore write the following:
$V^{model, corrupted}_{ij} = g_i g^*_j V^{model}_{ij} = e^{i (\phi_i - \phi_j)} V^{model}_{ij} = e^{i f \phi} V^{model}_{ij}$,
where f is a (M x N) matrix and $\phi$ is a (N x 1) vector of the gain phases (for N=4 antennas I have M=N(N-1)/2 visibilities), so that
$\begin{pmatrix} \phi_1 - \phi_2 \\ \phi_1 - \phi_3 \\ \phi_1 - \phi_4 \\ \phi_2 - \phi_3 \\ \phi_2 - \phi_4 \\ \phi_3 - \phi_4 \\ \end{pmatrix} = \begin{pmatrix} 1 & -1 & 0 & 0 \\ 1 & 0 & -1 & 0 \\ 1 & 0 & 0 & -1 \\ 0 & 1 & -1 & 0 \\ 0 & 1 & 0 & -1 \\ 0 & 0 & 1 & -1 \\ \end{pmatrix} \cdot \begin{pmatrix} \phi_1 & \phi_2 & \phi_3 & \phi_4 \end{pmatrix} = f\phi$
I can now minimise the expression S, using for example a gradient descent algorithm, for the vector $\phi$. However, as you can imagine as the number of antennas increases the problem becomes more and more computationally demanding.
I was wondering if there is a more efficient way to approach my problem, either analytic (doubtful) or with a more sophisticated algorithm that takes advantage of the fact that my solution has the form that I want (i.e. $g_j = e^{i \phi_j}$).
I assume that the $\sigma_{ij}$ are known.
Let $\Sigma$ denote the matrix with entries $\sigma_{ij}$. Let $U$ denote the matrix $$ U = \exp\operatorname{diag}(\phi_1,\phi_2,\phi_3) = \pmatrix{e^{i\phi_1 }\\ & e^{i \phi_2}\\ & & e^{i \phi_3}}. $$ Note that the matrix whose entries are $V_{ij}g_ig_j^*$ can be written as $UVU^*$. With that, your objective function can be expressed as $$ S = \sum_{i,j} \frac{1}{\sigma^2_{ij}}|V_{ij} - g_i g^*_j V^{model}_{ij}|^2 = \|\Sigma \odot (V - UV^{\text{model}}U^*)\|_F^2, $$ where $\odot$ denotes a Hadamard product. With that, a gradient with respect to the vector $\phi = (\phi_1,\phi_2,\phi_3)$ can be computed.
We have: $$ dS = d\|\Sigma\odot(V - UV^MU^*)\|_F^2 = 2 \operatorname{Re}\operatorname{tr}[(\Sigma\odot(V - UV^MU^*))^*d(\Sigma\odot(V - UV^MU^*))]. \\ d[\Sigma\odot(V - UV^MU^*)] = \Sigma \odot d(UV^MU^*).\\ d [UV^MU^*] = U(\phi + d\phi)V^MU^*(\phi + d\phi) - UV^MU^* \\ = U V^M[\operatorname{diag}(i e^{i\phi}d\phi)]^* + [\operatorname{diag}(i e^{i\phi}d\phi)]V^MU^*.\\ $$ To find the entries of the gradient, plug in $dx = e_j$ (the standard basis vectors). We have $dx = e_j \implies$ $$ d [UV^MU^*] = UV^M\operatorname{diag}(i e^{i\phi_j} e_j)^* + \operatorname{diag}(i e^{i\phi_j} e_j)V^MU^*\\ = -ie^{-i\phi_j} U (V^M e_j)e_j^T + ie^{i\phi_j}e_j(e_j^T V^M)V^MU^*.\\ d[\Sigma\odot(V - UV^MU^*)] = \Sigma \odot [-ie^{-i\phi_j} U (V^M e_j)e_j^T + ie^{i\phi_j}e_j(e_j^T V^M)V^MU^*].\\ dS =\\ 2 \operatorname{Re}\operatorname{tr}[(\Sigma\odot(V - UV^MU^*))^*(\Sigma \odot [-ie^{-i\phi_j} U (V^M e_j)e_j^T + ie^{i\phi_j}e_je_j^TV^MU^*])] \\ = 2 \operatorname{Re}\operatorname{tr}[\Sigma^{\odot2} \odot ((V - UV^MU^*)^*[-ie^{-i\phi_j} U (V^M e_j)e_j^T + ie^{i\phi_j}e_je_j^T V^MU^*])] \\ = 2 \operatorname{Re}[ -ie^{-i\phi_j} e_j^T( \Sigma^{\odot2} \odot [(V - UV^MU^*)^* U V^M ])e_j\\ \quad +ie^{-i\phi_j} e_j^T( \Sigma^{\odot2} \odot [(V - UV^MU^*)^* V^MU^*])e_j]. $$