$Cast AIv, $Cast VLev
Считаю задачу методом ParticleInCell, сейчас использую следующую неявную разностную схему для скорости:
v_{1}-v_{0}=\frac{e\delta t}{2m}\left(E_{0}+E_{1}+\frac{1}{c}\left(\left[v_{0}\times H_{0}\right]+\left[v_{1}\times H_{1}\right]\right)\right)
и для координаты:
r_{1}=r_{0}+0.5\delta t\left(v_{0}+v_{1}\right),
где \detla t — шаг по времени, e --заряд, E,H — электрическое и магнитное поле, c - -скор света. 0 — значение в текущий момент времени, 1 — в последующий.
в коде в немного оптимизированном виде она выглядит вот так:
//высчитываю некоторый предварительный результат избавившись от величин с 0м индексом
u0(1) = v0(1)+H0(3)*v0(2)-H0(2)*v0(3)+E0(1)
u0(2) = -H0(3)*v0(1)+v0(2)+H0(1)*v0(3)+E0(2)
u0(3) = H0(2)*v0(1)-H0(1)*v0(2)+v0(3)+E0(3)
xi0=x0+0.5*dt*v0(1)
//и приступаю к основному счёту
C0=1./(1.+H1(1)**2+H1(2)**2+H1(3)**2)
v1(1)=(1+H1(1)**2)*(u0(1)+E1(1))*C0&
+(H1(3)+H1(2)*H1(1))*(u0(2)+E1(2))*C0&
+(H1(3)*H1(1)-H1(2))*(u0(3)+E1(3))*C0
v1(2)=-(H1(3)-H1(2)*H1(1))*(u0(1)+E1(1))*C0&
+(1+H1(2)**2)*(u0(2)+E1(2))*C0&
+(H1(1)+H1(3)*H1(2))*(u0(3)+E1(3))*C0
v1(3)=(H1(3)*H1(1)+H1(2))*(u0(1)+E1(1))*C0&
+(-H1(1)+H1(3)*H1(2))*(u0(2)+E1(2))*C0&
+(1+H1(3)**2)*(u0(3)+E1(3))*C0
x1=xi0 + 0.5*dt*v1(1)
последний фрагмент итерирую несколько раз и уже за 2-3 итерации он сходится. Но эта штука работает для чисто декатовой (X*Y) геометрии. Возникла необходимость посчитать всё это в полноценной цилиндрической (R*Z) геометрии а значит необходимо учитывать центробежную силу, действующую на частицы, двигающиеся в азимутальном направлении. Например, для метода LeapFrog есть такая формула для скорости, которая применяется после обычной ф-лы для нахождения новой скорости:
v1(1)=v0(1)
v1(2)=v0(2)+v0(3)**2/R*dt-1.5*v0(3)**2*v0(2)/R**2*dt**2
v1(3)=v0(3)-v0(2)*v0(3)/R*dt+v0(3)/R**2*(v0(2)**2-0.5*v0(3)**2)*dt**2
