Solving system of ODEs using Fourier transforms

I have a system of ODEs I'd like to solve. It seems in the literature that the best way of getting a complete solution is with Fourier transforms (and that under some conditions perturbative solutions do not work).
\begin{align} \frac{d P_{11}(t)}{dt} &= V(t) P_{21}- V(t)^* P_{12}\\ \frac{d P_{12}(t)}{dt} &= V(t)(P_{22}- P_{11})-c P_{12}\\ \frac{d P_{21}(t)}{dt} &= V^*(t)(P_{22}- P_{11})-c P_{21}\\ \frac{d P_{22}(t)}{dt} &= V^*(t) P_{12}- V(t) P_{21} \end{align} where ${P_{11},P_{12},P_{21},P_{22}}$ are functions of t that I'm solving for; c is a constant, V* is the complex conjugate of V; and $V = E_1 e^{-i \omega_1 t} + E_2 e^{-i \omega_2 t}$, with $E_1$ and $E_2$ being constants.

My plan is to Fourier transform this set of ODEs, solve the resultant linear algebra, and then inverse Fourier transform back to get my solution.

So first I apply some standard rules to Fourier transforms (convolution, $F(P(t) \cdot e^{i\omega_1 t}) = F(P)(\omega-\omega_1)$, and the rules for Fourier transformed derivatives $F(\dot{\rho}) = i \omega F(\rho)$).

But in doing so, I see that my system of equations (after transforming) has the form:

$$\begin{equation} i \omega \hat{\rho}_{11} + \hat{\rho}_{21}(\omega-\omega_1)-\hat{\rho}_{12}(\omega - \omega_2) = 0 \\ -(\hat{\rho}_{11}(\omega-\omega_1) -\hat{\rho}_{11}(\omega-\omega_2) + (i\omega + c)\hat{\rho}_{12}+ \hat{\rho}_{22}(\omega-\omega_1)+\hat{\rho}_{22}(\omega-\omega_2)) = 0 \end{equation} \\ ...$$ My four equations now have terms like $\hat{\rho}_{21}(\omega)$ and $\hat{\rho}_{21}(\omega-\omega_1)$. In the basic methods of Fourier transforming to solve different equations, we solve algebraically for the fouriered functions. But in this case it's not obvious how to solve this system algebraically. I've thought about manually solving this system by hand using substitution (of both the input arguments and the functions themselves), but haven't had any luck.

Any ideas on how I can move forward to solve the system?

Edit: The fourier transform of a derivative is known to be ${\widehat {f'\;}}(\xi )=2\pi i\xi {\hat {f}}(\xi )$. This is found by performing integration by parts and identifying that the function at the boundaries is zero. In my case, I'm not completely sure if $\lim_{x \to \infty} \hat{\rho}_{12}(\omega)$ (but, as discussed in the comments, I think that it's likely that there's some "energy" constraints limiting the maximum amount of frequencies allowed) And so, for now, I have omitted any potential contribution due to the boundaries.

• What are $E_1$ and $E_2$? What is $V^*$? – amsmath Sep 5 '18 at 5:49
• E1 and E2 are simply constants. V* is the complex conjugate of V – Steven Sagona Sep 5 '18 at 13:25
• @StevenSagona In your last equation I don't understand why there is $\hat{\rho}_{11}(0)$ and what is $\omega_{ab}$. – Delta-u Sep 5 '18 at 15:13
• I fixed $\omega_{ab}$. The original differential equation has a lot more constants in it, and I, when writing down my notes, did not properly convert everything to the correct constants. – Steven Sagona Sep 5 '18 at 15:18
• @StevenSagona Ok :-). I obtain similar equations (with $\omega+\omega_{1/2}$ instead of $-$ but I may have made an error somewhere) but I have no $\hat{\rho}(0)$. Where did it comes from ? – Delta-u Sep 5 '18 at 15:21

A partial answer, but hopefully somewhat useful:

First off, it would be instructive to know how you obtained this system, but in the absense of this context, we can still do a few things.

Notice that $$\frac{\text{d}}{\text{d} t}\left(P_{11}+P_{22}\right) = 0$$; hence, the system reduces to \begin{align} \frac{\text{d} x}{\text{d} t} &= 2 V\,z - 2 \,V^*y\\ \frac{\text{d} y}{\text{d} t} &= V\,x - c\,y\\ \frac{\text{d} z}{\text{d} t} &= V^*x - c\,z\\ \end{align} with $$x = P_{11}-P_{22}$$, $$y=P_{12}$$ and $$z=P_{21}$$. Now, when the frequencies $$\omega_1$$ and $$\omega_2$$ are commensurable, the matrix $$A(t) = \begin{pmatrix} 0 & -2 V^* & 2 V \\ V & -c & 0 \\ V^* & 0 & -c\end{pmatrix}$$ is periodic. In that case, you can try to attack the above system using Floquet theory, see e.g. here.

Alternatively, you can rewrite $$y = e^{-c t} Y$$, $$z = e^{-c t} Z$$ and $$V = e^{- c t} W$$ to obtain the even simpler system \begin{align} \frac{\text{d} x}{\text{d} t} &= 2 W\,Z - 2 W^*Y\\ \frac{\text{d} Y}{\text{d} t} &= W\,x \\ \frac{\text{d} Z}{\text{d} t} &= W^*x. \end{align} Unfortunately, $$W$$ is not periodic anymore; however, you might be able to make some headway with this formulation.

• Thanks. I'll add a note that explains how I get up to this point. I'll have to think about the commensurable requirement. The simulation is based on real frequencies of light in an experiement, so I'm not sure if you can say with certainty if the number is rational (but I'm sure it's fine to some approximation.) I'm aware that floquet theory is commonly used for time-dependent problems, but I'm not experienced enough with it to really make progress and am not aware of much. – Steven Sagona Oct 15 '18 at 22:07
• I totally agree: if the frequencies are (to be) experimentally obtained, then any rational relationships are questionable, unless there is some strong system symmetry. Also, I do have to admit that I have struggled to use Floquet theory to obtain explicit results, so its usefulness would maybe not be immediately clear -- maybe someone else would be able to help here. – Frits Veerman Oct 16 '18 at 8:52

Personally, I have never been a big fan of "Fourier transform" methods- or transform methods in general. Since these are linear equations with constant coefficients, I would write them in matrix form: $\frac{d}{dt}\begin{pmatrix}P_{11} \\ P_{12} \\ P_{21} \\ P_{22}\end{pmatrix}= \begin{pmatrix}0 & -V^* & V & 0 \\ -V & -c & 0 & V \\ -V^* & 0 & -c & V^* \\ 0 & V^* & V & 0 \end{pmatrix}\begin{pmatrix}P_{11} \\ P_{12} \\ P_{21} \\ P_{22}\end{pmatrix}$.

Now determine the eigenvalues and eigenvectors of that matrix.

• Maybe this isn't clear in the question, but V is a function of t, so I'm not sure if finding the eigenvectors or eigenvalues is sufficient to get a full solution (without some extra work maybe). I'll edit the question to make it more clear V is not a constant! – Steven Sagona Sep 5 '18 at 15:29
• If the coefficient matrix function $A(t)$ commutes with its integral for all $t$, then $\exp\left(\int_a^t A(s)\text ds\right)$ is a solution. I'm not sure if this is the case here, but it's worth checking if you're interested in a solution regardless of it's form. – AlexanderJ93 Sep 5 '18 at 20:10
• It's not clear to me what you mean here. What is a coefficient matrix function? Also, the answer should have exponential terms associated with it, such as $e^{i \omega_1 t}$ (because I know this is the form of the final answer) – Steven Sagona Sep 6 '18 at 15:51