如何求数值时间序列的多重极大值

Question

我正在用欧拉方法数值求解一阶微分方程。我从时间t=0开始建立解，直到某个任意时间以0.05的时间步长进行。结果解决方案的一个示例如下图所示

我想找出我们在图中看到的最大值（以及它们出现的次数），并将它们存储在字典中。如果在整个时间范围内只有一个最大值，我可以使用这个代码

y=[xinitial,viinitial,ainitial]            
t=0
maximum=-20000
maxai={}
h=0.05
ai=-2.1
for i in range(0,3701):
   dydt=computederivs(y)
   y = euler(y,dydt,h)
   t+=h
   if y[0]>maximum:
      maximum = y[0]
 maxai[ai]=maximum

由于我有多个局部极大值，我必须在时间t中移动时检测它们，通过某种方式检查函数在爬升几个时间步后何时下降。我还需要将最大值存储在一个列表中，该列表是字典键的值。 我在想象这是一个足够普遍的任务，必须有众所周知的方法来做它或多或少简单？你知道吗

Answer 1

您应该寻找的事件是当一阶导数v改变其符号时。对于最大值，则测试二阶导数a是否为负。你知道吗

 yp, y = y, euler(y,dydt,h)
 if py[1]*y[1] < 0 and y[2] < 0:
      # do what you need to do to list the maximum

可以使用线性插值来获得更精确的最大位置坐标。你知道吗

---

如果你想得到严肃的结果，可以采用更高阶的ODE积分方法。你知道吗

Answer 2

我建议使用来自scipy.signal模块的find_peaks。该函数采用一维数组，通过对相邻值的简单比较找到所有的局部极大值。可选地，可以通过指定峰值特性的条件来选择这些峰值的子集。

下面是一段代码片段，让您开始：

import matplotlib.pyplot as plt
import numpy as np
from scipy.signal import find_peaks

Fs = 8000
f = 5
sample = 8000
x = np.arange(sample)
y = np.sin(2 * np.pi * f * x / Fs)
peaks = find_peaks(y)

plt.scatter(peaks[0], np.ones(f), c='red')
plt.plot(x, y)
plt.xlabel('sample(n)')
plt.ylabel('voltage(V)')
plt.show()

你的最高标准：

print(peaks[0])
[ 400 2000 3600 5200 6800]