3
私は粒子のシステムを持っています。彼らは加速、速度と位置を持っています。この粒子系はどのようにしてエネルギーが増加し続けていますか？
粒子が壁に当たったとき、その速度が反転します。 2個の粒子が互いに近づくと、それらはこの力によって互いに反発：
F=1/r^2
又は
F_x=delta(x)/r^3
F_y=delta(y)/r^3
システムが実行中であった場合、私は全ての粒子の合計速度が増加している感じていました。どちらが奇妙ですか？粒子はそのエネルギーを に与えるべきです。したがって、システムの総エネルギーは一定に保たれなければなりません。
システムの運動エネルギーは、私は、システム全体の総エネルギーを監視維持しcout介して印刷
E_k=Sigma v^2
に等しく、それは増加し続けることを確認します。それはエネルギーの保守性に反する。コード内で私はどこで間違いを犯していますか？
#include <vector>
#include <cstdlib>
#include <cmath>
#include <iostream>
constexpr size_t N=1000;
struct Point
{
double x, y;
double v_x, v_y;
double a_x, a_y;
};
Point points[N];
void next_frame()
{
double energy=0.0;
// calculate forces
for(size_t i = 0; i < N; ++i)
{
double fx=0.0,fy=0.0;
for(size_t j = 0; j < N; ++j)
{
if(i!=j)
{
double dx=points[i].x-points[j].x;
double dy=points[i].y-points[j].y;
double r2=dx*dx+dy*dy;
if(r2>0.01 && r2<100.0) // avoid nan and also unnecessary computation
{
// F=1/r^2
double r=sqrt(r2);
fx+=dx/(r*r*r);
fy+=dy/(r*r*r);
}
}
}
points[i].a_x=0.01*fx;
points[i].a_y=0.01*fy;
energy+=points[i].v_x*points[i].v_x+points[i].v_y*points[i].v_y;
}
std::cout<<energy<<std::endl;
for(size_t i = 0; i < N; ++i)
{
// integrations
points[i].v_x += points[i].a_x;
points[i].v_y += points[i].a_y;
points[i].x += points[i].v_x;
points[i].y += points[i].v_y;
// wall
if(points[i].x < -50.0)
points[i].v_x = +std::abs(points[i].v_x);
else if(points[i].x > +50.0)
points[i].v_x = -std::abs(points[i].v_x);
if(points[i].y < -50.0)
points[i].v_y = +std::abs(points[i].v_y);
else if(points[i].y > +50.0)
points[i].v_y = -std::abs(points[i].v_y);
}
}
int main(int argc, char **argv)
{
// initialize particles
for(size_t i = 0; i < N; ++i)
{
Point p;
p.x = -50 + ((rand() % 1000)/1000.0)*100.0;
p.y = -50 + ((rand() % 1000)/1000.0)*100.0;
p.a_x=0.0;
p.a_y=0.0;
p.v_x=0.001*((rand() % 1000)/1000.0-0.5);
p.v_y=0.001*((rand() % 1000)/1000.0-0.5);
points[i]=p;
}
while(1)
{
next_frame();
}
return 0;
}
これは、反復回数を超えるエネルギープロファイルです：
この質問のタグを変更することは避けてください。
私が物理的なフォーラムでこの質問をすると、彼らは物理的ではなくプログラミングの問題であると私に教えてくれるでしょう。
あなたは神粒子を紛失しています：P –
あなたの計算が正しいと仮定すると、数値エラーが蓄積されている可能性があります。 – CodeMonkey
また、離散化手法に問題がある可能性があります。解像度を上げたり下げたりして、解像度が結果を変えるかどうかを調べることができます。 – CodeMonkey