swh:1:snp:122bde0cb0e54f3d002c308e151c63f07e45e6be
Tip revision: fff29182352fedf9c5323dedecb66a7cb8aa1f6f authored by Hanno Rein on 04 March 2024, 20:38:19 UTC
4 pl
4 pl
Tip revision: fff2918
integrator_leapfrog.c
/**
* @file integrator.c
* @brief Leap-frog integration scheme.
* @author Hanno Rein <hanno@hanno-rein.de>
* @details This file implements the leap-frog integration scheme.
* This scheme is second order accurate, symplectic and well suited for
* non-rotating coordinate systems. Note that the scheme is formally only
* first order accurate when velocity dependent forces are present.
*
* @section LICENSE
* Copyright (c) 2011 Hanno Rein, Shangfei Liu
*
* This file is part of rebound.
*
* rebound is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* rebound is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with rebound. If not, see <http://www.gnu.org/licenses/>.
*
*/
#include <stdio.h>
#include <stdlib.h>
#include "rebound.h"
// Leapfrog integrator (Drift-Kick-Drift)
// for non-rotating frame.
void reb_integrator_leapfrog_part1(struct reb_simulation* r){
r->gravity_ignore_terms = 0;
const unsigned int N = r->N;
struct reb_particle* restrict const particles = r->particles;
const double dt = r->dt;
#pragma omp parallel for schedule(guided)
for (unsigned int i=0;i<N;i++){
particles[i].x += 0.5* dt * particles[i].vx;
particles[i].y += 0.5* dt * particles[i].vy;
particles[i].z += 0.5* dt * particles[i].vz;
}
r->t+=dt/2.;
}
void reb_integrator_leapfrog_part2(struct reb_simulation* r){
const unsigned int N = r->N;
struct reb_particle* restrict const particles = r->particles;
const double dt = r->dt;
#pragma omp parallel for schedule(guided)
for (unsigned int i=0;i<N;i++){
particles[i].vx += dt * particles[i].ax;
particles[i].vy += dt * particles[i].ay;
particles[i].vz += dt * particles[i].az;
particles[i].x += 0.5* dt * particles[i].vx;
particles[i].y += 0.5* dt * particles[i].vy;
particles[i].z += 0.5* dt * particles[i].vz;
}
r->t+=dt/2.;
r->dt_last_done = r->dt;
}
void reb_integrator_leapfrog_synchronize(struct reb_simulation* r){
// Do nothing.
}
void reb_integrator_leapfrog_reset(struct reb_simulation* r){
// Do nothing.
}