This is a problem when I'm reading a paper. 


Equation:


$min\{\sum\_p(S\_p-I\_p)^2+\beta((\partial\_xS\_p-h\_p)^2+(\partial\_yS\_p-v\_p)^2) \} $


where $S,I,h,v$ are all $M\*N$ matrices and p stands for every element in the Matrix. $I,h,v$ are known.


The paper just mentioned "we diagonalize derivative operators after Fast Fourier Transform for speedup" and get the solution


$S=\mathscr{F}^{-1}\left(\frac{\mathscr{F}(I)+\beta(\mathscr{F}(\partial\_x)^\*\mathscr{F}(h)+\mathscr{F}(\partial\_y)^\*\mathscr{F}(v))}{\mathscr{F}(1)+\beta(\mathscr{F}(\partial\_x)^\*\mathscr{F}(\partial\_x)+\mathscr{F}(\partial\_y)^\*\mathscr{F}(\partial\_y)}\right)$


where $\mathscr{F}(1)$ stands for FFT of delta function. Plus, multiplication and division are all component-wise. "\*" means conjugation.


I've looked up some books but still can't get how this happen. I don't find any connection between DFT and minimization or quadratic function.


Here's the [paper](http://www.cse.cuhk.edu.hk/~leojia/papers/L0smooth_Siggraph_Asia2011.pdf), equation is on page 4.



How to use DFT to solve this minimization problem?
You need to minimize the objective function


$$f:S \mapsto \| S-I\|\_2^2+\beta(\| \partial\_xS-h\|^2\_2+\|\partial\_yS-v\|\_2^2) $$


where $\| \cdot \|\_2$ is the "entrywise" $\ell^2$ norm. This is done as usual: the global minimum $S$ must be a critical point and then the derivative must be zero at $S$. The differential of the objective function is given by


$$Df(S)(M) = 2\langle S-I,M\rangle + \beta\left( 2 \langle\partial\_x S -h,\partial\_x M\rangle + 2 \langle\partial\_y S - v,\partial\_y M \rangle \right)$$


where $\langle A,B \rangle = \sum\_p A\_p B\_p$ (the "entrywise" inner product). Now, that for any $M$ the last expression vanishes implies that


$$ S-I + \beta( \partial\_x^T \partial\_x S-\partial\_x^Th + \partial\_y^T \partial\_y S -\partial\_y^Tv) =0.$$


($\cdot^T$ is the adjoint operator). Now you can take the DFT on both sides,


$$ \mathscr{F}(S) - \mathscr{F}(I) + \beta \left( \mathscr{F}(\partial\_x^T\partial\_x S) - \mathscr{F}(\partial\_x^Th)+ \mathscr{F}(\partial\_y^T\partial\_y S) - \mathscr{F} (\partial\_y^Tv)\right)=0,$$
and then use the fact that the DFT diagonalizes the gradient $\mathscr{F}(\partial\_. S) = \mathscr{F}(\partial\_\cdot ) \mathscr{F} (S)$ and some algebra to find that
$$ (\mathscr{F}(1) + \beta (\mathscr{F}(\partial\_x)^\*\mathscr{F}(\partial\_x) + \mathscr{F}(\partial\_y)^\*\mathscr{F}(\partial\_y) ))\mathscr{F}(S)= \mathscr{F}(I) + \beta \left( \mathscr{F}(\partial\_x)^\*\mathscr{F}(h) + \mathscr{F}(\partial\_x)^\*\mathscr{F} (v)\right)$$
and solving for $S$ this gives the formula asked. (Remark: We could have written the objective function directly in terms of $\mathscr{F}(S)$ using Parseval and then optimize).


As for the speedup part, we would have to invert a very large ($NM \times NM$) matrix representing the operator 
$$ 1 + \beta( \partial\_x^T \partial\_x + \partial\_y^T \partial\_y )$$
if we do it directly in image space. 


