<p>I wrote a code that follows the Barnes-Hut algorithm for gravitational dynamics. Everything looks pretty good, except that when I plot the total energy $$E_{tot}\equiv E_k+E_u=const$$ I get this result:</p>

<p><a href="https://i.sstatic.net/4C1rP.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/4C1rP.png" alt="enter image description here"></a></p>

<p>and you can see that although it seems constant for most of the time, it has some increasing transient in the beginning. I couldn't find a bug in the code that could be responsible for this, is there may be other explanation that can make sense of this?</p>

CC BY-SA 4.0
2018-08-15T13:38:10.073Z
Constants of motion not staying constant is a common problem with numerical integration codes. It needn't come from a very obvious bug, but merely from an inappropriate choice of algorithm or implementation.  I'd suggest you search the literature, as it's very unlikely you're the first to run into this.
Mark Olson
CC BY-SA 4.0
2018-08-15T21:15:01.197Z
The Barnes-Hut has major problems with discontinuity. If a star moves from an area into another, what will happen?
peterh
CC BY-SA 4.0
2018-08-16T04:07:17.600Z
This is a simulation question, not a programming question; don't migrate to Stack Overflow but consider [Computational Science SE](https://scicomp.stackexchange.com/). Anyway, what happens if you shorten the timestep?
Mike G
CC BY-SA 4.0
2018-08-16T06:18:04.853Z
What kind of integrator did you implement?
Mike G
CC BY-SA 4.0
2018-08-16T06:20:58.190Z
Runge-Kutte 4 order, shorten the time step didn't change anything.
אבנר יעקב
CC BY-SA 4.0
2018-08-16T13:50:44.007Z
A leapfrog or Verlet integrator may conserve energy better.
Mike G
<galaxy-cluster><n-body-simulations>
Energy conservation in Barnes-Hut algorithm
אבנר יעקב
# Energy conservation in Barnes-Hut algorithm

(Asked by: username_0 on 2018-08-15)

I wrote a code that follows the Barnes-Hut algorithm for gravitational dynamics. Everything looks pretty good, except that when I plot the total energy $$E\_{tot}\equiv E\_k+E\_u=const$$ I get this result:

[![enter image description here](https://i.sstatic.net/4C1rP.png)](https://i.sstatic.net/4C1rP.png)

and you can see that although it seems constant for most of the time, it has some increasing transient in the beginning. I couldn't find a bug in the code that could be responsible for this, is there may be other explanation that can make sense of this?

**Tags:** galaxy-cluster, n-body-simulations
> **Comments:**
> **username_1**: Constants of motion not staying constant is a common problem with numerical integration codes. It needn't come from a very obvious bug, but merely from an inappropriate choice of algorithm or implementation.  I'd suggest you search the literature, as it's very unlikely you're the first to run into this.  
> **username_2**: The Barnes-Hut has major problems with discontinuity. If a star moves from an area into another, what will happen?  
> **username_3**: This is a simulation question, not a programming question; don't migrate to Stack Overflow but consider [Computational Science SE](https://scicomp.stackexchange.com/). Anyway, what happens if you shorten the timestep?  
> **username_3**: What kind of integrator did you implement?  
> **username_0**: Runge-Kutte 4 order, shorten the time step didn't change anything.  
> **username_3**: A leapfrog or Verlet integrator may conserve energy better.  

---

## 0 answers

