为什么Python（NumPy）比MATLAB慢100倍以上呢？

Question

几年来，我每天都在使用MATLAB，但最近我决定改用Python，因为它有额外的功能（例如机器学习库等）。你知道吗

我首先使用NumPy将一个简单的MATLAB代码翻译成Python，然后在Spyder环境中运行它。结果是惊人的： Python慢了100多！下面我附上了两个代码。对于给定的CSV文件（gmr），Python在25秒内执行了代码，而MATLAB在几毫秒内执行了代码。更重要的是，我想对大量CSV文件（例如200）重复相同的计算。令人吃惊的是，MATLAB需要大约50秒来处理200个CSV文件，而Python需要超过1小时。你知道吗

你能建议一些修改吗？我做错什么了？你知道吗

我尝试了几种不同的方法来制定代码，但不幸的是，它们都不起作用。你知道吗

MATLAB：

clear variables; close all ; clc ;

tic;
damp=0.05;
T=[0.01:0.01:3]';
Spectra = cell(num_records,2);

name='rec_1_out_scaled.csv';
gmr=load(name);
dt = gmr(3,1) - gmr(2,1);
cd (main_Dir)
Sa=zeros(length(T),1);
% Newmark's Direct Integration Method
for j=1:length(T)
    Sa(j)=newmark_dim(gmr(:,2),T(j),damp,dt);
end
Spectra(i,:) = {name,Sa};

toc;


function Sa = newmark_dim(gacc,T,zeta,dt)

% newmark coefficients
beta = 1/4;
gamma = 1/2;

% natural circular frequency of SDOF system
wn = 2*pi/T;

% initialization
npun = length(gacc);
y = zeros(npun,1);
yp = zeros(npun,1);
ypp = zeros(npun,1);
ypp(1) = -gacc(1)-2*wn*zeta*yp0-wn^2*y0;

% Integration coefficients
keff = wn^2 + 1/(beta*dt^2) + gamma*2*wn*zeta/(beta*dt);
a1 = 1/(beta*dt^2)+gamma*2*wn*zeta/(beta*dt);
a2 = 1/(beta*dt)+2*wn*zeta*(gamma/beta-1);
a3 = (1/(2*beta)-1)+2*wn*zeta*dt*(gamma/(2*beta)-1);

for i=1:npun-1
    y(i+1)   = (-gacc(i+1)+a1*y(i)+a2*yp(i)+a3*ypp(i))/keff;
    ypp(i+1) = (y(i+1)-y(i)-dt*yp(i)-dt^2*ypp(i)/2)/(beta*dt^2) + ypp(i);
    yp(i+1)  = yp(i)+dt*ypp(i)+dt*gamma*(ypp(i+1)-ypp(i));
end

Sd = max(abs(y));
Sa = Sd*wn^2;

return

PYTHON：

import time
import numpy as np
import math
import pandas as pd

start = time.time()

# Initialize
name = "rec_1_out_scaled.csv"
T = np.arange(0.01,3.01,0.01)
zeta = 0.05 # 5% damping
extension = 'csv'
num_records = 200
Spectra = np.zeros((num_records,T.size,2)) # Multidimentional array

"""
Newmark's Direct Integration Method
"""
def Newmark_DIM(gmr,T,zeta,dt):
    # newmark coefficients for constant acceleration
    beta = 1/4
    gamma = 1/2

    # natural circular frequency of the SDOF system
    wn = 2*math.pi/T

    # initialization
    npun =  len(gmr)
    y = np.zeros(npun)
    yp = np.zeros(npun)
    ypp = np.zeros(npun)
    ypp[0] = - gmr[0] - 2 * wn * zeta* yp[0] - wn**2 * y[0]

    # Integration coefficients
    keff = wn**2 + 1/(beta * dt**2) + gamma * 2 * wn * zeta /(beta * dt) 
    a1 = 1/(beta * dt**2) + gamma * 2 * wn * zeta/(beta * dt)
    a2 = 1/(beta * dt) + 2 * wn * zeta * (gamma/beta - 1)
    a3 = (1/(2 * beta) - 1) + 2 * wn * zeta * dt * (gamma/(2 * beta) - 1)    

    # Solves ODE
    for i in range(npun-1):
        y[i+1]   =  (-gmr[i+1] + a1 * y[i] + a2 * yp[i] + a3 * ypp[i])/keff
        ypp[i+1] = (y[i+1] -y[i] - dt * yp[i] - dt**2 * ypp[i]/2) / (beta * dt**2) + ypp[i]
        yp[i+1]  = yp[i] + dt * ypp[i] + dt * gamma * (ypp[i+1] - ypp[i])

    Sd = np.amax(np.abs(y))
    Sa = Sd * wn**2

    return (Sa)

# Calculates Spectrum
gmr = pd.read_csv(name, header=None)
gmr = gmr.values # makes it numpy array
dt = gmr[2,0] - gmr[1,0]
Sa = np.zeros(len(T))
for j in range(len(T)):
   Sa[j] = Newmark_DIM(gmr[:,1],T[j],zeta,dt)
Spectra[0][:,0] = T
Spectra[0][:,1] = Sa

print("Finished")
end = time.time()
print(end - start)