Python代码中速度数组未更新问题排查求助

1. Immediate Issue: Loop Never Executes #

The while loop condition s[-1] < St1Length is false from the start, so no iterations run. Here's why:

• St1Length = 400, ds=40, so St1Length//ds = 10.
• vel is initialized as np.zeros(10) (10 zeros).
• s = np.linspace(0, 400, len(vel)) creates an array where the last element is exactly 400, making s[-1] < 400 evaluate to False. No new velocity values are added because the loop never runs.

2. Secondary Issue: Incorrect Acceleration Calculation #

Even if the loop ran, the acceleration values (Ap and Aw) are normalized by gravity (g), which breaks the velocity update formula. The equation v² = u² + 2*a*ds requires a to be in m/s², not a dimensionless fraction of g.

Fixes: #

A. Correct Initialization

Start with initial position and velocity values, then build the arrays iteratively:

vel = np.array([0.0])  # Initial velocity (0 m/s)
s = np.array([0.0])    # Initial position (0 m)

B. Fix Acceleration Calculations

• For Ap: Remove division by g to get actual acceleration from power and drag:

  Ap = (PMax / ((m * vel[-1]) + 1e-9)) - (Fd / m)
  

• For Aw: Rewrite to calculate maximum non-wheelie acceleration in m/s² (derived from torque balance, no division by g):

  Aw = (b * g / h) - (Fd * ha) / (m * h)
  

Corrected Code #

import numpy as np
import matplotlib.pyplot as plt

g = 9.81  # Gravity
pa = 1.2  # Air density
m = 250  # mass of bike
h = 0.69  # centre of mass height
ha = 0.69  # height of centre pressure
w = 1.5  # wheelbase
b = 0.73  # longitudinal distance of centre of mass
CdA = 0.2  # drag coefficient
PMax = 180000  # max power
ux = 1.2  # Longitudinal friction coefficient
uy = 1.44  # Lateral friction coefficient

St1Length = 400
ds = 40

# Initialize with starting values instead of filling upfront
vel = np.array([0.0])
s = np.array([0.0])

while s[-1] < St1Length:
    Fd = 0.5 * pa * CdA * vel[-1]**2
    # Calculate actual acceleration from power (no division by g)
    Ap = (PMax / ((m * vel[-1]) + 1e-9)) - (Fd / m)
    yw = 0
    # Calculate maximum non-wheelie acceleration (no division by g)
    Aw = (b * g / h) - (Fd * ha) / (m * h)
    # Use minimum of the two accelerations (both in m/s²)
    Wheelie_Min = np.minimum(Ap, Aw)
    # Ensure acceleration isn't negative (avoid deceleration if drag exceeds power)
    Wheelie_Min = max(Wheelie_Min, 0)
    # Update velocity using correct acceleration units
    vel_new = np.sqrt(vel[-1]**2 + 2 * Wheelie_Min * ds)
    vel = np.append(vel, vel_new)
    
    s_new = s[-1] + ds
    s = np.append(s, s_new)

plt.plot(s, vel)
plt.xlabel('Distance (m)')
plt.ylabel('Velocity (m/s)')
plt.title('Bike Acceleration Over Distance')
plt.grid(True)
plt.show()

Explanation of Changes #

• Initialization: Starting with single-element arrays ensures the loop runs until s reaches 400m.
• Acceleration Units: Removing the division by g converts Ap and Aw to proper m/s² values, which aligns with the velocity update formula.
• Negative Acceleration Guard: Adding Wheelie_Min = max(Wheelie_Min, 0) prevents velocity from decreasing if drag ever exceeds available power.

When you run the corrected code, the velocity array will update with each iteration, and the chart will show a rising velocity curve limited by either power or wheelie constraints.

内容的提问来源于stack exchange，提问作者Karl Burnett