I am looking to numerically solve the (complex) Time Domain Ginzburg Landau Equation. I wish to write a python simulator to observe the nucleation of fluxons over a square 2D superconductor domain (eventually 3D, cubic domain). I am using a fourth order Runge Kutta solver for this which I made in Python (the odeint solver does not handle complex variables). How do I implement the RK4 algorithm for complex variables? The equation I am looking to solve is.... ∂ψ/∂t = -p/i ∇2 - q/i |ψ|2 ψ + γ ψ Where p and q are chosen within the complex unit square (unit magnitude), and γ is real on [0,1]. Any help in implementing this would be appreciated. I'd also appreciate if you could point me to software or libraries that can advance my efforts.