LTSpice压敏电阻参数拟合问题：minimize算法无法输出正确结果

我编写了脚本，以9个随机参数为初始值，填入LTSpice的.asc配置文件进行压敏电阻（Varistor）仿真，将仿真结果传入scipy minimize函数迭代优化，最终将优化结果与实测数据对比绘图。脚本运行正常，但minimize始终无法输出正确结果，仅随初始参数的调整接近或偏离目标值，且无法计算出准确的初始参数。目前已尝试更换优化方法、调整初始参数，但问题仍未解决，相关实现代码如下：

import PyLTSpice as lt
from PyLTSpice.LTSpice_RawRead import RawRead
import matplotlib.pyplot as plt
import openpyxl
import numpy as np
import lmfit as lm
from sklearn.metrics import mean_squared_error as mse
from scipy.stats import linregress

sheet1_row = []
sheet2_row = []

def getInput(string, default):
    #val = input(string)
    return default
    if val != "":
        return val
    else:
        return default

def createWorkbook(path):
    return 

def getData(worksheet, productcode, productmode):
    data = []
    rows = worksheet.max_row
    cols = worksheet.max_column
    for j in range(1, rows):
        if worksheet.cell(row = j, column = 3).value == productcode and worksheet.cell(row = j, column = 4).value == productmode:
            for i in range(5, cols):
                if worksheet.cell(row = j, column = i).value is not None:
                    data.append(float(worksheet.cell(row = j, column = i).value))
            break

    if data == []:
        print("Product not found")
    return data

global_voltage = []

def array_to_string():
    string = ""

    for nbr in global_voltage:
        string = string + str(nbr) + " "
    return string

def inputRonRoff(output_filename, a1, a2, a3, a4, a5, a6, a7, a8,a9):
    file = open(r"./res/D_ideal_3.asc", 'r').read()
    new_file = (file.replace("Ron=a1", f"Ron={a1}")
    .replace("Roff=a2", f"Roff={a2}")
    .replace("Vfwd=a3", f"Vfwd={a3}")
    .replace("Vrev=a4", f"Vrev={a4}")
    .replace("Ilimit=a5", f"Ilimit={a5}")
    .replace("RevIlimit=a6", f"RevIlimit={a6}")
    .replace("Epsilon=a7", f"Epsilon={a7}")
    .replace("RevEpsilon=a8", f"RevEpsilon={a8}")
    .replace("Rrev=a9", f"Rrev={a9}")
    .replace("1 2 3 4 5 6 7 8 9", array_to_string()))
    print("String is " + array_to_string()) 
    output_filepath = "./" + output_filename + ".asc"
    open(output_filepath, 'w').write(new_file)

def Simulation(filename):
    filepath = "./" + filename + ".asc"
    sim = lt.SimCommander(filepath)
    sim.run()
    sim.wait_completion()
    
def Get_Results_From_Raw(output_filename):
    raw_filepath = "./" + output_filename + "_1.raw" # or + 
    print(raw_filepath)
    raw = RawRead(raw_filepath)
    voltage = raw.get_trace("V(n001)").data #v1
    current = raw.get_trace("I(D1)").data
    return voltage, current

def run_simulation(r1, r2, r3, r4, r5, r6, r7, r8,r9):
    output_filename = "test"
    inputRonRoff(output_filename, r1, r2, r3, r4, r5, r6, r7, r8,r9)
    Simulation(output_filename)
    return 
    
def calculate_loss(v_params, meas_current):
    r1 = v_params["r1"].value
    r2 = v_params["r2"].value
    r3 = v_params["r3"].value
    r4 = v_params["r4"].value
    r5 = v_params["r5"].value
    r6 = v_params["r6"].value
    r7 = v_params["r7"].value
    r8 = v_params["r8"].value
    r9 = v_params["r9"].value
    run_simulation(np.log(r1), np.log(r2), np.log(r3), np.log(r4), np.log(r5), np.log(r6), np.log(r7), np.log(r8), np.log(r9))
    calc_voltage, calc_current = Get_Results_From_Raw("test")
    loss = mse(meas_current, calc_current)
    return loss

#sperdoret
def init_parameters():
    v_params_init = [1,1,1,1,1,1,1,1,1]
    r1i,r2i,r3i,r4i,r5i,r6i,r7i,r8i,r9i = v_params_init
    return [r1i,r2i,r3i,r4i,r5i,r6i,r7i,r8i,r9i]

def init_minimizer(r1i,r2i,r3i,r4i,r5i,r6i,r7i,r8i,r9i):
    v_params = lm.Parameters()
    v_params["r1"] = lm.Parameter(name="r1", value=r1i, min=1e-5, max=1, vary=True)
    v_params["r2"] = lm.Parameter(name="r2", value=r2i, min=3, max=3e5, vary=False)
    v_params["r3"] = lm.Parameter(name="r3", value=r3i, min=36, max=100, vary=False)
    v_params["r4"] = lm.Parameter(name="r4", value=r4i, min=-1, max=-36, vary=False)
    v_params["r5"] = lm.Parameter(name="r5", value=r5i, min=-1e24, max=1e24, vary=True)
    v_params["r6"] = lm.Parameter(name="r6", value=r6i, min=-1e24, max=1e24, vary=True)
    v_params["r7"] = lm.Parameter(name="r7", value=r7i, min=-1e24, max=1e24, vary=False)
    v_params["r8"] = lm.Parameter(name="r8", value=r8i, min=-1e24, max=1e24, vary=False)
    v_params["r9"] = lm.Parameter(name="r9", value=r9i, min=1e1, max=1e5, vary=True)
    return v_params
    
def solver(v_params, meas_current, meas_voltage):    
    result = lm.Minimizer(calculate_loss, v_params, fcn_args=([meas_current])).minimize(method='least_squares')
    return result

def solver_result(result):
    r1 = float('{:.3e}'.format(((result.params["r1"].value))))
    r2 = float('{:.3e}'.format(((result.params["r2"].value))))
    r3 = float('{:.3e}'.format(((result.params["r3"].value))))
    r4 = float('{:.3e}'.format(((result.params["r4"].value))))
    r5 = float('{:.3e}'.format(((result.params["r5"].value))))
    r6 = float('{:.3e}'.format(((result.params["r6"].value))))
    r7 = float('{:.3e}'.format(((result.params["r7"].value))))
    r8 = float('{:.3e}'.format(((result.params["r8"].value))))
    r9 = float('{:.3e}'.format(((result.params["r9"].value))))

    print("r1 is " + str(r1))
    print("r2 is " + str(r2))
    print("r3 is " + str(r3))
    print("r4 is " + str(r4))
    print("r5 is " + str(r5))
    print("r6 is " + str(r6))
    print("r7 is " + str(r7))
    print("r8 is " + str(r8))
    print("r9 is " + str(r9))

    run_simulation(r1, r2, r3, r4, r5, r6, r7, r8,r9)
    calc_voltage, calc_current = Get_Results_From_Raw("test")
    
    return calc_voltage, calc_current


def plot_result(c, v, current, voltage):
    plt.loglog()
    plt.plot(c, v, color="red", label="Final Result") # output is red
    plt.plot(current, voltage, color="blue", label="Sample")
    plt.xlabel("Voltage")
    plt.ylabel("Current")
    plt.legend(loc="upper left")
    plt.show()

def calc_Ron(current, voltage):
    delta_current = current[1] - current[0]
    delta_voltage = voltage[1] - voltage[0]
    Ron = delta_voltage / delta_current if delta_current != 0 else 1
    return Ron

def calc_Roff(current, voltage):
    delta_current = current[1] - current[0]
    delta_voltage = voltage[1] - voltage[0]
    return (delta_voltage / delta_current)

def calc_Vfwd(current, voltage):
    return voltage[5]

#MAIN FUNCTION
def main():
    #productcode = getInput("Enter Product Code: ", "820412001")820413001
    productcode = getInput("Enter Product Code: ", "820412001")
    productmode = getInput("Enter Product Mode: ", "combined")
    workbook = openpyxl.load_workbook("./Varistors.xlsx", data_only=True)
    current = getData(workbook.worksheets[1], productcode, productmode)
    voltage = getData(workbook.worksheets[2], productcode, productmode)
    global global_voltage
    global_voltage = voltage
    result = solver( init_minimizer(calc_Ron(current, voltage),
                                    calc_Roff(current, voltage),
                                    36.3, 
                                    -29.7, 0, 0, 0.01, 0.01, 2.5e4), current, voltage)
    v, c = solver_result(result)
    plot_result(c, v, current, voltage)

if __name__ == "__main__":
    main()

附上对比曲线图：

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