For every complex number $a$ one can define a function $Q_a$ on $\mathbb{C}$ ;

$Q_a(b)=$IP[ar, ai, br, bi],

$a,b\in \mathbb{C},ar=Re(a),ai=Im(a),br=Re(b),bi=Im(b)$,

IP[ar_, ai_, br_, 
  bi_] := (Norm[
     Exp[-(ar^2 + ai^2 + br^2 + bi^2)/2]*
      Exp[(br - I bi)*(ar + I ai)]]^2)/Pi

with IP defined as such.

For $a=1+i$, this is what the function would look like

ar = 1;(*Real Part*)
ai = 1;(*Imaginary part*)
 Plot3D[IP[ar, ai, br, bi], {br, -1, 3}, {bi, -1, 3}, 
  PlotRange -> Full](*Taken over all complex numbers br+bi I*)

enter image description here

Now, I have a different function $Q_{p.d}$, where

$Q_{p,d}(b)=$Norm[(br + bi I)]^(2 p) Exp[-Norm[(br + bi I) - d]^2]/(Pi p! LaguerreL[p, -Norm[d]^2])

$p \in \mathbb{R}, b,d\in \mathbb{C},br=Re(b),bi=Im(b)$

p = 10;
d=dr+di I;
Plot3D[Norm[(br + bi I)]^(2 p) Exp[-Norm[(br + bi I) - d]^2]/(Pi p! LaguerreL[p, -Norm[d]^2]), {br, 1, 5}, {bi, 1, 5}, PlotRange -> Full]

enter image description here

What I want

Approximate the function $Q_{p,d}$ for a fixed $p$ and $d$ using positive linear combinations of $Q_a$ i.e., Minimize the 3D volume $Q_{p,d}-\sum_{i=1}^\infty c_i Q_{a_i}$, where $c_i$ are real positive coeficients. I have to optimize for the complex numbers $a_i$ and the coefficients $c_i$ so that the difference is minimum.

What kind of techniques are used to tackle these types of problems?

  • 1
    $\begingroup$ What is the definition of Q[a]? What is the definition of Q[p,d]? In this post, there two complex variables a,b in the first function and there are three variables in the second function, p is real and b,d are complex. $\endgroup$
    – cvgmt
    Aug 26, 2022 at 10:23
  • $\begingroup$ Does my edit answer your queries? $\endgroup$
    – Dotman
    Aug 26, 2022 at 10:43

2 Answers 2


This seems a classic multidimensional expansion of a vector relative to a basis. I can not calculate an infinite expansion, but I need to truncate the expansion to a finite number of terms.

For practical reasons (I do not have the nerves to wait long) I give as an example a coarse expansion. You may later increase the number of basis function if needed.

We will expand your example:

p = 10;
d = 1 + I;
Qb[br_, bi_] = 
  Norm[(br + 
       bi I)]^(2 p) Exp[-Norm[(br + bi I) - d]^2]/(Pi p! LaguerreL[
       p, -Norm[d]^2]);

Plot3D[Qb[x, y], {x, 0, 5}, {y, 0, 5}]

enter image description here

We choose a basis a grid of your IP functions. To get a better expansion you must make the grid denser.

A single IP function looks like:

IP[ar_, ai_, br_, 
   bi_] := (Norm[
      Exp[-(ar^2 + ai^2 + br^2 + bi^2)/2]*
       Exp[(br - I bi)*(ar + I ai)]]^2)/Pi;
Plot3D[IP[1, 1, x, y], {x, -1, 3}, {y, -1, 3}]

enter image description here

The function is different from zero approx. in a circle of radius 2 around the center. Qb has a center of approx.: {3,3}, therefore we choose a grid from {1,2,3,4,5} in x and y direction:

bas[x_, y_] = 
  Flatten[Table[IP[bx, by, x, y], {bx, 1, 5}, {by, 1, 5}], 1];
centers = Flatten[Table[{bx, by}, {bx, 1, 5}, {by, 1, 5}], 1];
Plot3D[bas[x, y], {x, 0, 6}, {y, 0, 6}, PlotRange -> {0, 0.4}]

enter image description here

First we need to calculate the scalar product of all basis vectors: g=<bi|bj>, the so called metric tensor. To speed up the calculation we restrict the domain of integration. Toward this aim we pair each basis vector with its center, resulting in basc. With this we get g:

basc = Transpose[{bas[x, y], centers}];
g = Outer[
  NIntegrate[#1[[1]] #2, {x, #1[[2, 1]] - 3, #1[[2, 1]] + 
      3}, {y, #1[[2, 2]] - 3, #1[[2, 2]] + 3}] &, basc, bas[x, y], 1];

Next we need the scalar product between Qb and the basis vectors: rs

rs = Outer[
    NIntegrate[#1[[1]] #2[x, y], {x, #1[[2, 1]] - 3, #1[[2, 1]] + 
        3}, {y, #1[[2, 2]] - 3, #1[[2, 2]] + 3}] &, basc, {Qb}, 1] // 

We get the expansion coefficients by solving the linear problem: g. coefficients==rs:

coefs = LinearSolve[g, rs];

and the expansion:

exp = bas[x, y] . coefs;

We can plot the expansion together with the original Qb:

Plot3D[{exp, Qb[x, y]}, {x, 0, 5}, {y, 0, 5}]

enter image description here

  • $\begingroup$ This looks very nice!! Exactly what I was hoping for. Thank you so much $\endgroup$
    – Dotman
    Sep 6, 2022 at 12:07
  • $\begingroup$ Is there a way to restrict the final coefficients to be positive? $\endgroup$
    – Dotman
    Sep 6, 2022 at 12:07
  • $\begingroup$ What would the aim of this be? The sign of the coefficients depend on the sign of the basis function and these are arbitrary. $\endgroup$ Sep 6, 2022 at 12:30
  • $\begingroup$ I understand that they would be arbitrary. But negative solutions would be unphysical in my model. The coefficients are to be normalised to form a probability distribution. $\endgroup$
    – Dotman
    Sep 6, 2022 at 13:34
  • 1
    $\begingroup$ I do not see that the expansion coefficients have a physical meaning. They are only a tool to model a given function in an arbitrary basis. I do not know what you are doing, but probabilities often depend on the square of expansion coefficients. $\endgroup$ Sep 6, 2022 at 14:04

The Concept

You want a Generalized Fourier Transform.


Your use of IP and eyeballing $Q_a(b)$, tells me that you're using this as an inner product on $\mathbb{C}$. That's fine, but unrelated to the discussion of inner product below. Since you want to form $Q_{p,d}$ as a linear combination of $\mathcal{Q} := \{Q_{a_i}\}_{i\in \mathbb{N}}$, your vectors are the linear combinations of all such functions, in the vectorspace $V:=\mathop{span}\left(\mathcal{Q}\right)$. The inner product discussion below is for $V$.


For this to work, you need:

  1. $V\subseteq L^2(\mathbb{C})$, your functions need to be square integrable on $\mathbb{C}$,
  2. to find/define an inner product on $V$,
  3. a basis $P$ for $V$ which is orthonormal with respect to that inner product.

(1) Square Integrable Functions

For all $Q_a(b)$ (that is, for all $a\in\mathbb{C}$), you need $$\iint_{b\in\mathbb{C}}Q_a(b)\overline{Q_a(b)}\; db \le \infty.$$

Since $Q_a$ has real output, you don't need to conjugate and this is $$\iint_{b\in\mathbb{C}}Q_a(b)^2\; db \le \infty.$$

I tested $Q_0$ and got $$\iint_{b\in\mathbb{C}}Q_0(b)^2 \; db = \frac{1}{2\pi}.$$

Eyeballing $Q_a$, I assume they're all square integrable, so moving on.

(2) Inner Product on $V$

You need to choose a function $w:\mathbb{C}\to\mathbb{R}$ which will define an inner product on $V$ like this: $$\left\langle Q_{a_1}, Q_{a_2} \right\rangle_w := \iint_{b\in\mathbb{C}} w(b)Q_{a_1}(b)Q_{a_2}(b) \; db$$

I won't go into all the requirements needed for this to qualify as an inner product, but those conditions must be satisfied.

(3) Orthonormal Basis

Finally, you need to define a set of functions $\mathcal{P}=\{P_j\}_{j\in\mathbb{N}} \subseteq V$ (that is, each $P_j$ is a linear combination of your $Q_{a_i}$) such that

  1. $\mathop{span}{(\mathcal{P})} = V$, any linear combination of $\mathcal{Q}$ can be gotten as a linear combination of $\mathcal{P}$,
  2. $\left\langle P_j,P_j \right\rangle_w = 1$ for all $P_j \in \mathcal{P}$,
  3. $\left\langle P_j,P_k \right\rangle_w = 0$ for all $P_j,P_k \in \mathcal{P}$ where $P_j \ne P_k$.

Note that 1. is not saying that any $P_j$ can be gotten as a linear combination of $\mathcal{Q}$. That's already satisfied since $P_j \in V$, which is the set of all linear combinations of $\mathcal{Q}$. Instead, it's saying that every function in $V$ must be represented as a linear combination of $\mathcal{P}$. ($\mathcal{P}$ is large/complete enough that you can get to anything in $V$ using $\mathcal{P}$ as a basis.)

Conditions 2. says that each $P_j$ has unit norm and 3. says that each function in $\mathcal{P}$ is orthogonal to each other. Together, 1. and 3. ensure that $\mathcal{P}$ is a basis for $V$.

(4) The Transform

Finally, $$Q_{p,d} = \sum_{i\in\mathbb{N}} k_j P_j$$ where $$k_j = \frac{\left\langle Q_{p,d}, P_j \right\rangle_w}{\left\langle P_j, P_j \right\rangle_w}.$$ This is a linear combination of $\mathcal{Q}$ because every $P_j$ is a linear combination of $\mathcal{Q}$. If you need to represent $Q_{p,d}$ in terms of the $Q_{a_i} \in \mathcal{Q}$, you can expand each $P_j = \sum_{i\in\mathbb{N}} m_{i,j} Q_{a_i}$ and distribute each $k_j$ over the linear combination: $$k_j P_j = k_j \sum_{i\in\mathbb{N}} m_{i,j} Q_{a_i} = \sum_{i\in\mathbb{N}} k_j m_{i,j} Q_{a_i}$$ and then collect them in terms of each $i$ so that $$c_i := \sum_{j\in\mathbb{N}} k_j m_{i,j}.$$

Practical Considerations

In the best world, your $Q_a$ would initially be defined so that you can use $\mathcal{Q}$ as $\mathcal{P}$, and use the constant function $w(b)=1$ to define the inner product. You can't do that here though.

The norm associated with this inner product is $$\|Q_a\|_w := \sqrt{\left\langle Q_a, Q_a \right\rangle_w}.$$

Each $Q_a$ could be scaled so that it has unit norm for $w(b)=1$. However, these functions are not orthogonal to each other for this choice of $w$.

Consider the regular Fourier Transform. The orthonormal basis there is the set of all (scaled) integer frequency sine waves with phase parameter. That is, $$\mathcal{P} = \left\{\frac{1}{\sqrt{\pi}}\sin(ix+\phi_i)\right\}_{i\in\mathbb{N}}$$ is an orthonormal basis for their span, for any choice of $\{\phi_i\}_{i\in\mathbb{N}}$. Because of the relation between trig and the exponential function, the phase $\phi_i$ can be absorbed into the constant $c_i$ - how isn't important, only to note that your basis vectors/functions can be parameterized.

The more important thing is to note that these functions have different scale factors. That is, some have large wavelengths and some have tiny wavelengths. Your functions $Q_a$ are parameterized to move around via $a$, but it's not clear to me that they can be used to represent large and small scale features anywhere in the plane. I plotted several and visually they all seemed the same size as $a$ moved the plot around the plane, but $\|Q_a\|_w$ changed as they moved around (I used $w(b)=1$).

Notice also that the $\mathop{sin}$ functions take on positive and negative values and are thus able to cancel each other in order to be orthogonal.

I think you'll have the best chance of success if you build $\mathcal{P}$ as a wavelet basis, so that each $P_j$ takes on positive and negative values (via negative $m_{i,j}$) and they have different frequency content / scale factors independent of their translation across the plane. It's not clear to me that this can be done with your current definition of $Q_a$.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.