这不应使用迭代拟合;该问题是精确确定的，并且应该通过分析计算：

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from scipy.optimize import curve_fit

def quadratic(x: np.ndarray, a: float, b: float, c: float) -> np.ndarray:
    return a*x**2 + b*x + c

def estimate(x: np.ndarray, y: np.ndarray, p: float, q: float) -> np.ndarray:
    """Calculate parabolic parameters exactly through both roots and a point close to the vertex"""
    iv = y.argmin()  # index of point closest to vertex
    vx = x[iv]       # x-value near to vertex
    vy = y[iv]       # y-value near to vertex

    return np.linalg.solve(
        a=(
            (  p*p,  p, 1),  # p root
            (  q*q,  q, 1),  # q root
            (vx*vx, vx, 1),  # near the vertex
        ),
        b=(0, 0, vy),
    )

def analytic_integral(a: float, b: float, c: float, p: float, q: float) -> float:
    """Definite integral of quadratic(a, b, c) from p to q"""
    return (
        a*(q**3 - p**3)/3
        + 0.5*b*(q**2 - p**2)
        + c*(q - p)
    )

def solve(p: float, q: float, integ: float) -> tuple[float, float, float]:
    # to fix the integral
    b = integ/(
        0.5*(q**2 - p**2)
        - (
            q*p*(q - p)
            + (q**3 - p**3)/3
        )/(q + p)
    )
    # to fix the first and second roots
    a = -b/(q + p)
    c = a*p*q
    return a, b, c

def main() -> None:
    x = np.array((258.6669469, 258.831, 259.129, 259.428, 259.726,
                  260.025, 260.324, 260.622, 261.15238824))
    y = np.array((0, -0.062, -0.139, -0.181, -0.197, -0.193, -0.17, -0.129, 0))

    # empirical integral, to be preserved
    integ = np.trapz(y, x)

    # assume that the first and last data points are the roots
    p = x[0]
    q = x[-1]

    # estimate through three points, disregarding integral
    popt_est = estimate(x, y, p, q)
    # generic unconstrained fit
    popt_approx, _ = curve_fit(f=quadratic, xdata=x, ydata=y, p0=popt_est)

    df = pd.DataFrame(
        {
            'estimated': popt_est,
            'approx': popt_approx,
            'exact': solve(p, q, integ),  # through both roots, preserving integral
        },
        index=pd.Index(name='parameter', data=('a', 'b', 'c')),
    ).T
    params = df.values.T

    dev = quadratic(x[:, np.newaxis], *params) - y[:, np.newaxis]
    df['integral'] = analytic_integral(*params, p, q)
    df['integral_err'] = df.integral - integ
    df['root1_err'] = quadratic(p, *params)
    df['root2_err'] = quadratic(q, *params)
    df['sqerr'] = (dev*dev).sum(axis=0)
    pd.set_option('display.float_format', '{:.6e}'.format)
    print(df.T)

    xlong = np.linspace(p, q, 1000)
    fig, ax = plt.subplots()
    ax.scatter(x, y, label='Empirical data')

    for title, param in df[['a', 'b', 'c']].iterrows():
        ax.plot(xlong, quadratic(xlong, *param), label=title)
    ax.legend()
    plt.show()

if __name__ == '__main__':
    main()

个字符
x1c 0d1x的数据