主页 / 问题 / 1287310
Accepted
Fire13nyu
Asked:2022-05-29 02:19:20 +0800 CST2022-05-29 02:19:20 +0800 CST

Werlet集成方法的实现

  • 772

再会！我想在不使用速度的情况下实现 Verlet 集成方法（https://en.wikipedia.org/wiki/Verlet_integration），但我偶然发现了拥有当前坐标和以前坐标的问题。我正在介绍一个具有三个实例（一维、二维和多维案例）的坐标类：

class Coord coord where  
    infixl 6 +++
    (+++) :: coord -> coord -> coord
    infixl 7 ***
    (***) :: Double -> coord -> coord
    infixl 6 -.-
    (-.-) :: coord -> coord -> coord
    (-.-) coord1 coord2 = coord1 +++ (-1.0)***coord2

instance Coord Double where
    coord1 +++ coord2 = coord1 + coord2
    coord1 *** coord2 = coord1 * coord2

instance Coord (Double, Double) where
    (coord11,coord12) +++ (coord21,coord22) = (coord11+coord21,coord12+coord22)
    num *** (coord11,coord12) = (num*coord11,num*coord12)

instance Coord [Double] where 
    list1 +++ list2 = zipWith (+) list1 list2
    num *** list = map (* num) list

然后我用输入参数做积分函数本身：加速度（函数），初始坐标，初始速度和时间步长，并在每一步之后返回一个无限的点位置列表

verlet :: Coord coord => (coord -> coord) -> coord -> coord -> Double -> [coord]
verlet a x0 v0 dt = 
    let x1 = x0 +++ (dt***v0) +++ (0.5*(dt**2))***(a x0)
        verlet' xPrev xNow = let xNext = 2***xNow -.- xPrev +++ (dt**2))***(a xNow)
                            in (xPrev,xNow) : verlet' xNow xNext
    --- in (,) : verlet'  - вот тут проблема

如果有的话，我是一个真正的初学者，所以，很可能，这很容易解决，但它对我不起作用，所以我寻求帮助）

haskell

2 个回答

  1. Best Answer

    您已经拥有x0，这是给您的，并且x1您根据 Wikipedia 的说明计算得出。然后，按照相同的说明，只需将它们用作 的前两个参数verlet'：

    verlet :: Coord coord => (coord -> coord) -> coord -> coord -> Double -> [coord]
verlet a x0 v0 dt = 
    let x1 = x0 +++ (dt***v0) +++ (0.5*(dt**2))***(a x0)
        verlet' xPrev xNow = let xNext = 2***xNow -.- xPrev +++ (dt**2))***(a xNow)
                             in (xPrev,xNow) : verlet' xNow xNext
    in verlet' x0 x1

    添加：

    从您的评论来看，您需要一个单独的坐标列表作为输出，而不是它们的对。如果是这样，那么只需让函数verlet'完全返回此列表：

    ...
in xPrev : verlet' xNow xNext
...
    • 1

  2. verlet :: Coord coord => (coord -> coord) -> coord -> coord -> Double -> [coord] 奇怪的系统...从好的方面来说，它是这样实现的：

    struct Body
{
    Vec3d pos { 0.0, 0.0, 0.0 };
    Vec3d vel { 2.0, 0.0, 0.0 }; // 2m/s along x-axis
    Vec3d acc { 0.0, 0.0, 0.0 }; // no acceleration at first
    double mass = 1.0; // 1kg
    double drag = 0.1; // rho*C*Area - simplified drag for this example

    /**
     * Update pos and vel using "Velocity Verlet" integration
     * @param dt DeltaTime / time step [eg: 0.01]
     */
    void update(double dt)
    {
        Vec3d new_pos = pos + vel*dt + acc*(dt*dt*0.5);
        Vec3d new_acc = apply_forces(); // only needed if acceleration is not constant
        Vec3d new_vel = vel + (acc+new_acc)*(dt*0.5);
        pos = new_pos;
        vel = new_vel;
        acc = new_acc;
    }

    Vec3d apply_forces() const
    {
        Vec3d grav_acc = Vec3d{0.0, 0.0, -9.81 }; // 9.81m/s^2 down in the Z-axis
        Vec3d drag_force = 0.5 * drag * (vel * abs(vel)); // D = 0.5 * (rho * C * Area * vel^2)
        Vec3d drag_acc = drag_force / mass; // a = F/m
        return grav_acc - drag_acc;
    }
};

    也许会有所帮助）

    • -1

相关问题