在不知道函数的情况下计算限制在某个区域的曲线下的面积

问题描述 投票:0回答:1

我有图像中显示的蓝色和黑色图表的数据点(x 和 y 值)。我不知道蓝色和黑色图的功能。

enter image description here

如果我知道 x1 和 x2 坐标(如图所示),如何计算限制在 x1 和 x2 之间的 x 坐标的蓝线下方的面积。任何人都可以为我提供一个使用 numpy 和 scipy 库通过复合梯形 (numpy.trapz) 和辛普森 (scipy.integrate.simpson) 规则执行作业的示例。我是Python新手。我感谢您在这方面的帮助。

谢谢你。

python numpy scipy numerical-integration
1个回答
0
投票

首先我们加载您的文件:

import io
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import integrate

def load(file):
    with open(file, "r") as handler:
        raw = handler.readlines()[24:]
    data = pd.read_fwf(io.StringIO("".join(raw)), header=None, names=["x", "y"])
    return data

pol = load("pol_2000.xvg")
sol = load("sol_2000.xvg")

我们找到曲线交点(检查它们是否被同等采样后):

idx = np.argwhere(np.diff(np.sign(sol.y - pol.y))).flatten()
sel = slice(idx[0], idx[1] + 1)

交叉点发生在:

xl = sol.x[idx]
# 15      0.839916
# 239    12.978100

现在整合选择就足够了:

I = integrate.simpson(sol.y[sel], x=sol.x[sel])
# 3026.198371039349

以图形方式呈现如下:

fig, axe = plt.subplots()
axe.plot(sol.x, sol.y, label="Sol")
axe.plot(pol.x, pol.y, label="Pol")
axe.fill_between(sol.x[sel], sol.y[sel], alpha=0.4, label="A = %.2f" % I)
axe.legend()
axe.grid()

© www.soinside.com 2019 - 2024. All rights reserved.