Modification of the Verlet algorithms for the pendulum problem

Question

I'm trying to write a program to integrate the motion equations of the pendulum in the damped and forced case, that is, following this equation: $$ \frac{d^2\theta}{dt^2}=-\frac{g}{L}\sin(\theta)-\mu\frac{d\theta}{dt}+F\sin(\Omega t) $$ where $\mu$ is the friction coefficient, $F$ a multiplicative constant for the external force, which in turn is a periodic function with frequency $\Omega$.

I can define the usual system of ODEs without any difficult and use an integrator of the Runge Kutta type, but I encounter some difficulties while trying to use a symplectic integrator, say the Verlet integrators.

The position Verlet I used in different programs is: $$ x^{n+\frac{1}{2}}=x^n+\frac{h}{2}v^n\\ v^{n+1}=v^n+ha(x^{n+\frac{1}{2}})\\ x^{n+1}=x^{n+\frac{1}{2}}+\frac{h}{2}v^{n+1} $$ where $a(...)%$ is the acceleration, defined in the first equation of my question.

Just for the sake of clarity, here is the velocity Verlet: $$ v^{n+\frac{1}{2}}=v^n+\frac{h}{2}a(x^n)\\ x^{n+1}=x^n+hv^{n+\frac{1}{2}}\\ v^{n+1}=v^{n+\frac{1}{2}}+\frac{h}{2}a(x^{n+1}) $$

Here is my problem: in the expression of $a(...)$ the velocity appears in the damping term. But in both the Verlet algorithms I have no clue of how I can pass the velocity as an argument of the function.

In the position Verlet, I suppose I need to know $v^{n+\frac{1}{2}}$, while in the velocity Verlet i can use the initial conditions to kick the algorithm in first place, but in the third equation I think I need to use $v^{n+1}$, that is actually the result of the right hand side.

As I stated, I have no clue of how to resolve this issue, and actually, I'm not even sure if a solution is possible at all. I would be glad for any kind of answer, even a hint would be very kind.

Answer 1

I think that the answer and comments (including references to previous related questions), although correct and also interesting (I have discovered a kind of predictor-corrector variant of the velocity Verlet I was not aware of), are missing a point which could be relevant for your problem.

If the force depends linearly on velocity, the implicit step, which, for example, in the case of velocity Verlet would be $$ v_{n+1}=v_n + \frac{h}{2}\left(~~a_{n+1}(x_{n+1},v_{n+1},t_{n+1})+a_n(x_{n},v_{n},t_n)~~\right), $$ can be easily transformed into an explicit step, by solving a linear equation. Indeed, if $a_n = \tilde a_n(x_n,t_n) + w v_n$, the previous implicit equation for $v_{n+1}$ becomes: $$ v_{n+1}=v_n + \frac{h}{2}\left(~~a_n(x_{n},v_{n})+\tilde a_{n+1}(x_{n+1},t_{n+1}) + w v_{n+1}~~\right) $$ giving $$ v_{n+1} = \frac{v_n+\frac{h}{2}(a_n+\tilde a_{n+1})}{1-\frac{hw}{2} } $$

The resulting algorithm keeps most of the qualities of the usual velocity Verlet.