2

周期的境界条件または反射境界条件のいずれかを使用して、レナード・ジョーンズ流体の分子動力学シミュレーションを 2d ボックスに記述しようとしています。シミュレーションは反射境界でうまく動作しているように見えますが、何らかの理由で周期ボックス内の粒子は、数百回の積分ステップの後、余分な速度を得始め、最終的には粒子がボックス内に残りません。

時間ステップの設定が小さすぎると暴走速度が発生する可能性があることを理解しています。ただし、同じ速度-ベルレット アルゴリズムを使用してシミュレーションを実行すると、特に反射境界が非常にうまく機能するため、これが暴走速度の理由になる可能性があることがわかりません。この 2 つのケースでは、粒子にかかる力または粒子の分離を評価するために使用される方法に誤りがある可能性があると思いますが、自分で見つけることはできません。

力を評価するために使用しているコードは次のようになります。

public static Vector3d LJForce(Particle3d a, Particle3d b,
        Vector3d particleSep) {

    //Define R (for simplicity)
    Vector3d R = particleSep; //(b-a)

    return new Vector3d(Vector3d.multVector3d(
            -24*EPSILON*(2*power12(SIGMA/R.getMag())/R.getMag()
                    -power6(SIGMA/R.getMag())/R.getMag()), R));
}

これは、各粒子に作用するレナード・ジョーンズ力を与えると考えられています。粒子の分離は、境界条件によって決まります。反射の場合は次のとおりです。

@Override
public Vector3d getParticleSeparation(Particle3d a, Particle3d b) {
    return Particle3d.sepParticle3d(a, b);
}

と組み合わせて使用​​ されます

@Override
public void checkBoundaries(Particle3d a)   {
    Vector3d position = a.getPosition();
    Vector3d velocity = a.getVelocity();

    /*
     * We want to check the x,y and z coordinates. If
     * they break the boundaries we need to reflect 
     * the position in that coordinate and invert the
     * velocity in that axis.
     */
    //X coordinate
    if(position.getX()<0.0) {
        position.setX(Math.abs(position.getX()));
        velocity.setX(-velocity.getX());    
    }
    else if(position.getX()>getWidth()) {
        double howFar = position.getX()-getWidth();
        position.setX(getWidth()-howFar);
        velocity.setX(-velocity.getX());
    }   
    else    {}
    //Y coordinate
    if(position.getY()<0.0) {
        position.setY(Math.abs(position.getY()));
        velocity.setY(-velocity.getY());    
    }
    else if(position.getY()>getWidth()) {
        double howFar = position.getY()-getWidth();
        position.setY(getWidth()-howFar);
        velocity.setY(-velocity.getY());
    }   
    else    {}
    //Z coordinate
    if(position.getZ()<0.0) {
        position.setZ(Math.abs(position.getZ()));
        velocity.setZ(-velocity.getZ());    
    }
    else if(position.getZ()>getWidth()) {
        double howFar = position.getZ()-getWidth();
        position.setZ(getWidth()-howFar);
        velocity.setZ(-velocity.getZ());
    }   
    else    {}

    a.setPosition(position);
    a.setVelocity(velocity);
}//End checkBoundaries(i)

ボックス内に粒子を保持します。一方、定期的に...

@Override
public Vector3d getParticleSeparation(Particle3d a,Particle3d b) {
    Vector3d sep = Particle3d.sepParticle3d(a, b);
    if(sep.getX()>=(double)getWidth()/2)    {
        sep.setX(sep.getX()-getWidth());
    }   else{}
    if(sep.getY()>=(double)getWidth()/2)    {
        sep.setY(sep.getY()-getWidth());
    }   else{}
    if(sep.getZ()>=(double)getWidth()/2)    {
        sep.setZ(sep.getZ()-getWidth());
    }   else{}

    return sep;
}

粒子分離を与え、

@Override
public void checkBoundaries(Particle3d a) {
    Vector3d position = new Vector3d();

    position.setX((getWidth()+a.getPosition().getX())%getWidth());
    position.setY((getWidth()+a.getPosition().getY())%getWidth());
    position.setZ((getWidth()+a.getPosition().getZ())%getWidth());

    a.setPosition(position);
}

境界をチェックします。どちらの分離方法も呼び出します

//Relative seperation
public static Vector3d sepParticle3d( Particle3d a, Particle3d b) {
    return new Vector3d(Vector3d.subVector3d(b.getPosition(),a.getPosition()));
}

これは、デカルト空間における 2 つの粒子のベクトル分離です。周期的な境界条件でシミュレーションが失敗するのは、コーディングに誤りがあるか、またはその他の理由がありますか?

ありがとう

4

1 に答える 1

2

もうお気づきかもしれませんが、周期境界の checkBoundaries コードに誤りがあると思います。私は間違っているかもしれませんが、float/double では % 演算子を使用できないと思います。コードはそうあるべきだと思います（https://en.wikipedia.org/wiki/Periodic_boundary_conditionsも参照）

Vector3d op = a.getPosition(); //old position
double w = getWidth();
position.setX( op.getX() - w*((int)(2*op.getX()-1)) );
position.setX( op.getY() - w*((int)(2*op.getY()-1)) );
position.setX( op.getZ() - w*((int)(2*op.getZ()-1)) );
于 2012-06-06T13:47:09.857 に答える