司守奎《数学建模算法与应用》第二版第五章插值与拟合答案

最新推荐文章于 2024-07-23 14:36:35 发布

我玩游戏像芥雏子

最新推荐文章于 2024-07-23 14:36:35 发布

阅读量344

点赞数 3

文章标签： python

本文链接：https://blog.csdn.net/qq_55069195/article/details/138258374

版权

司守奎《数学建模算法与应用》第二版第五章插值与拟合答案

文章目录

- 5.1
- 5.2
- 5.3
- 5.4

5.1

在这里插入图片描述

import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d

# 解决中文显示问题
plt.rcParams["font.sans-serif"] = ["KaiTi"]  # 指定默认字体
plt.rcParams["axes.unicode_minus"] = False  # 解决保存图像是负号'-'显示为方块的问题

if __name__=="__main__":
    # 定义多项式的系数
    coefficients = [1,-6,5,-3]

    # 定义x范围
    x_min = -5
    x_max = 5
    num_points = 20

    # 在定义的范围内生成x值
    x_values = np.linspace(x_min, x_max, num_points)

    # 计算多项式的值
    y_values_true = np.polyval(coefficients, x_values)

    # 添加(0,1)之间的随机噪声模拟数据
    noise = np.random.rand(num_points)
    y_values_noisy = y_values_true + noise

    # 创建样条插值函数
    f_3 = interp1d(x_values, y_values_noisy, kind='cubic') # 三次样条插值
    f_2 = interp1d(x_values, y_values_noisy, kind='quadratic') # 二次样条插值

    # 估算新的数据点
    y_values_3 = f_3(x_values)  # 通过三次样条插值函数估算新数据点的值
    y_values_2 = f_2(x_values) # 通过二次样条插值函数估算新数据点的值

    # 绘制图形
    plt.plot(x_values,y_values_true,marker='o',markersize=4,linestyle='--',linewidth=3,label='真实数据',color='blue')
    plt.plot(x_values,y_values_noisy, marker='^',markersize=4,linestyle='--',linewidth=1,label='添加噪声后', color='red')
    plt.plot(x_values,y_values_3, marker='+',markersize=4,linestyle=':',linewidth=1,label='三次插值', color='black')
    plt.plot(x_values,y_values_2, marker='s',markersize=4,linestyle='-.',linewidth=1,label='二次插值', color='green')

    # 添加图例
    plt.legend()

    # 显示图形
    plt.grid(True)
    plt.show()

绘制的图形如下

在这里插入图片描述

5.2

在这里插入图片描述

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from scipy.interpolate import interp2d

if __name__=="__main__":
    # 原始数据点
    y = np.array([0,400,800,1200,1600,2000,2400,2800,3200,3600,4000,4400,4800])
    x=np.array([0,400,800,1200,1600,2000,2400,2800,3200,3600,4000,4400,4800,5200,5600])
    z=np.array([[370,470,550,600,670,690,670,620,580,450,400,300,100,150,250],
                [510,620,730,800,850,870,850,780,720,650,500,200,300,350,320],
                [650,760,880,970,1020,1050,1020,830,800,700,300,500,550,480,350],
                [740,880,1080,1130,1250,1280,1230,1040,900,500,700,780,750,650,550],
                [830,980,1180,1320,1450,1420,400,1300,700,900,850,810,380,780,750],
                [880,1060,1230,1390,1500,1500,1400,900,1100,1060,950,870,900,936,950],
                [910,1090,1270,1500,1200,1100,1350,1450,1200,1150,1010,880,1000,1050,1100],
                [950,1190,1370,1500,1200,1100,1550,1600,1550,1380,1070,900,1050,1150,1200],
                [1430,1450,1460,1500,1550,1600,1550,1600,1600,1600,1550,1500,1500,1550,1550],
                [1420,1430,1450,1480,1500,1550,1510,1430,1300,1200,980,850,750,550,500],
                [1380,1410,1430,1450,1470,1320,1280,1200,1080,940,780,620,460,370,350],
                [1370,1390,1410,1430,1440,1140,1110,1050,950,820,690,540,380,300,210],
                [1350,1370,1390,1400,1410,960,940,880,800,690,570,430,290,210,150]])

    # 创建二维插值函数
    interp_func = interp2d(x, y, z, kind='cubic')

    # 生成间隔为50的新的坐标点
    new_x = np.linspace(0, 5600, 113)
    new_y = np.linspace(0, 4800, 97)

    # 生成网格坐标矩阵
    X, Y = np.meshgrid(new_x, new_y)

    # 计算新坐标点对应的高程值
    new_z = interp_func(new_x, new_y)

    # 绘制二维平面的等高线
    plt.contour(X,Y,new_z)
    plt.xlabel('X')
    plt.ylabel('Y')
    plt.title('Contour Plot')
    plt.colorbar(label='Z') # 为z轴添加颜色条

    # 绘制三位平面的等高线
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')# 创建3D图形对象
    # surface = ax.plot_surface(X, Y, new_z, cmap='viridis')# 绘制3D曲面
    contour = ax.contour3D(X, Y, new_z, 50,cmap='viridis')# 绘制3D等高线图
    ax.set_xlabel('X')
    ax.set_ylabel('Y')
    ax.set_zlabel('Z')
    ax.set_title('3D Surface and Contour Plot')
    # fig.colorbar(surface, ax=ax, shrink=0.5, label='Z')# 添加颜色条

    # 显示图像
    plt.show()

绘制出的图形如图所示

在这里插入图片描述

5.3

在这里插入图片描述

import numpy as np
from scipy.optimize import curve_fit

if __name__ == "__main__":
    # 准备经验数据
    x = np.arange(1,9)
    y = np.array([15.3,20.5,27.4,36.6,49.1,65.6,87.87,117.6])

    # 定义拟合模型
    def empirical_model(x, a, b):
        return a * np.exp(b*x)

    # 使用最小二乘法拟合数据
    params, covariance = curve_fit(empirical_model, x, y)

    # 提取拟合参数
    a, b = params

    # 输出拟合结果
    print("拟合参数 a:", round(a,1))
    print("拟合参数 b:", round(b,1))

    # 输出经验公式
    print("经验公式：y ={}*exp({}*x)".format(round(a,1),round(b,1)))
# output:拟合参数 a: 11.4
# 拟合参数 b: 0.3
# 经验公式：y =11.4*exp(0.3*x)

5.4

在这里插入图片描述

将数据事先录入 $E x ce l$ 表格中， 5-4.xlsx 。

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d

# 解决中文显示问题
plt.rcParams["font.sans-serif"] = ["KaiTi"]  # 指定默认字体
plt.rcParams["axes.unicode_minus"] = False  # 解决保存图像是负号'-'显示为方块的问题

if __name__=="__main__":
    E = 30.24
    D = 57*E
    S = D**2*np.pi/4
    L1 = 514*1000
    L2 = 677.6*1000

    df = pd.DataFrame(pd.read_excel('5-4.xlsx'))
    time = np.round(df['时间'].values/3600, decimals=2)
    volume = df['水位'].values*0.01*E*S

    # 先求导求出水流量后进行插值
    tmp = np.logical_not(pd.isna(df['水位']).values) # 泵水时刻的索引
    volume = volume[tmp] # 剔除掉泵水的水位数据
    time = time[tmp] # 剔除掉泵水的时间数据
    velcity = -np.gradient(volume, time) # 计算剔除掉泵水后的水流量
    f = interp1d(time, velcity, kind='cubic') # 计算时间-水流量的三次样条插值函数
    x_values = np.linspace(0,24,100) # 在定义的范围内生成x值
    y_values = f(x_values) # 对生成的x值进行插值
    plt.subplot(2,2,1)
    plt.plot(time,velcity,marker='o',markersize=4,linestyle='--',linewidth=3,label='时间(h)-水流量(g/h)',color='blue')
    plt.title('时间(h)-水流量(g/h)')
    plt.subplot(2,2,2)
    plt.plot(x_values,y_values,marker='o',markersize=4,linestyle='--',linewidth=3,label='时间(h)-水流量(g/h)',color='blue')
    plt.title('时间(h)-水流量(g/h)')

    # 先进行插值后进行求导求出水流量
    f1 = interp1d(time, volume, kind='cubic') # 计算时间-水量的三次样条插值函数
    x_values = np.linspace(0,24,100) # 在定义的范围内生成x值
    y_values_1 = (f1(x_values)) # 对生成的x值进行插值
    velcity1 = -(np.gradient(y_values_1,x_values)) # 进行插值后进行求导求出水流量
    plt.subplot(2,2,3)
    plt.plot(time,volume,marker='o',markersize=4,linewidth=3,label='时间(h)-水量(g)',color='red')
    plt.title('时间(h)-水量(g)')
    plt.subplot(2,2,4)
    plt.plot(x_values,velcity1,marker='o',markersize=4,linewidth=3,label='时间(h)-水量(g)',color='red')
    plt.title('时间(h)-水流量(g/h)')

    plt.tight_layout()  # 调整子图布局，使其紧凑显示
    plt.show()