Python数据处理：应力应变曲线插值与分析

import pandas as pd
import numpy as np
from scipy.interpolate import interp1d

M0 = pd.read_csv('E:\shujushuchu\shuchu1.csv')  #导入四个参数
M1 = pd.read_csv('E:\shujushuchu\shuru1.csv')   #导入abaqus生成的应力应变曲线的点
B0=M0.values  
B1=M1.values  #将导入到参数转换为矩阵

a1 = int(np.size(B1, 0)/2)  #0表示行数，1表示列数，并赋值给a1
print(B1[1, 1])  #因为py语言里下标从0开始，检查数据的位置是否正确

B3 = np.zeros((a1, len(B1[0])))  #储存曲线点的矩阵
B4 = np.zeros((a1, len(B0[0])))  #储存预测参数的矩阵
                          #这两个矩阵根据前一步的计算结果修改这两个矩阵的大小，a1做为行元素的大小，而列元素大小要大于B1的shape
a2 = 0  #用于循环填充B3

for i in range(a1):  # a1值加1放在这儿，总循环，以行为单位
    B5 = B1[2*i:2*(i+1), :]
    B2 = np.zeros((2, len(B5[0])))  #创建一个空矩阵（该矩阵将用于储存应力应变曲线塑性阶段的点）
    C = B5[1, 2] / B5[0, 2]  #计算弹性模量
    a5 = np.size(B5, 1)  #计算B5矩阵有多少行
    a4 = 2   #用于之后循环跳出
    B7 = np.zeros((2, 3))   #创建空的B7矩阵便于储存条件屈服点（从弹性阶段进入塑性阶段）附近的三个点（用于寻找条件屈服点）
    for ii in range(2, a5):  

        if B5[1, ii - 1] <= C * (B5[0, ii - 1] - 0.002):  #寻找生成的数据中从弹性阶段进入塑性阶段瞬间的后一个点
            a4 = a4 + 1
            if a4 == 3:
                a6 = ii
                B7[0, 2] = B5[0, ii - 1]
                B7[1, 2] = B5[1, ii - 1]
            B2[0, a4 - 1] = B5[0, ii - 1]
            B2[1, a4 - 1] = B5[1, ii - 1]  #这两行循环储存条件屈服点后的实验数据于B2矩阵中


    B7[0, 0] = B5[0, a6 - 3]
    B7[0, 1] = B5[0, a6 - 2]
    B7[1, 0] = B5[1, a6 - 3]
    B7[1, 1] = B5[1, a6 - 2]    #该部分在B7中储存了条件屈服那一瞬间的后一个点和前两个点（用于后面的插值计算）
    B2[0, 1] = B5[0, a6 - 2]
    B2[1, 1] = B5[1, a6 - 2]    #该部分在B2中储存了条件屈服那一瞬间的后一个点
    a11 = (B7[0, 2] - B7[0, 1]) / 500   #将储存的三个点长度细分单位长度，用于插值计算获得中间更多的点
    f1 = interp1d(B7[0, :], B7[1, :], kind='linear')  #定义插值计算的数组、范围以及插值类型
    for i2 in np.arange(B7[0, 1], B7[0, 2], a11):  #利用for循环将数据细分并一一进行插值计算

        a12 = f1(i2)
        a13 = a12
        if a12 <= C * (i2 - 0.002):   #寻找条件屈服点的准确大概值
            B2[0, 1] = i2 - a11         #将条件屈服点应变储存进B2
            B2[1, 1] = a13 - C * a11      #将条件屈服点应力储存进B2
            break

    B6 = np.zeros((2, 10000))       #创建空白
    a7 = (np.max(B2[0, :]) - B2[0, 1]) / 500    #将塑性阶段细分单位长度用于插值计算
    a10 = 0    #用于循环递增储存插值后的细分的塑性阶段的点
    B2 = B2[:, np.nonzero(np.sum(B2, axis=0))[0]]    #去除矩阵末尾没储存数据的0值，保证插值计算的正常进行
    f2 = interp1d(B2[0, :], B2[1, :], kind='cubic')    #定义插值计算的数组、范围以及插值类型
    for i3 in np.arange(B2[0, 0], np.max(B2[0, :]) + a7, a7):    #利用for循环将数据细分并一一进行插值计算
        a20 = i3 + a7
        a22 = np.max(B2[0, :])
        if a20 > a22:
            a20 = a22
        a9 = f2(a20)
        a23 = i3
        if a23 > a22:
            a23 = a22
        a8 = f2(a23)       #为了寻找塑性失稳点，差值范围之前设置只到了原始数据的最大值处，插值时可能会超过这个范围而不能正常进行，所以将超过的部分等最大点处的数据
        derivative = (a9 - a8) / a7   #计算曲线的斜率
        if derivative - a9 >= 0:    #在斜率减缓到0之前，将所有的插值点都循环储存到B6矩阵里
            a10 += 1
            B6[0, a10 - 1] = i3
            B6[1, a10 - 1] = a8


    
    B6 = B6[:, np.nonzero(np.sum(B6, axis=0))[0]]     #去除矩阵末尾没储存数据的0值，保证插值计算的正常进行
    df = pd.DataFrame(B6)
    writer = pd.ExcelWriter('output3.xlsx', engine='openpyxl')
    df.to_excel(writer, sheet_name='Sheet1', index=False)
    writer.book.save('output3.xlsx')

    
    if np.max(B6[0, :]) > 0.47:      #将未能插值计算到应力极限点的数据去除掉（不能进行接下来的循环的话将不会储存进最终的矩阵里）
        continue
    else:
        a2 += 1
        for iiii in range(1, 5):
            B4[a2 - 1, iiii - 1] = B0[i - 1, iiii - 1]    #将参数储存进B4矩阵中

        B3[a2 - 1, 9] = np.min(B6[0, :])
        B3[a2 - 1, 10] = np.max(B6[0, :])    #将塑性阶段应变的最大值与最小值储存到B3矩阵中的第9列与第10列
        B6[0, :] = (B6[0, :] - B6[0, 0]) / (max(B6[0, :]) - B6[0, 0])   #将B6中的数据做归一化处理便于之后的插值计算（用于插值时取对应比例处的值）
        
                
        f = interp1d(B6[0, :], B6[1, :], kind='cubic')    #定义插值计算的数组、范围以及插值类型
        xx = np.array([0, 0.03, 0.09, 0.18, 0.3, 0.45, 0.63, 0.84, 1.0])    
        yy = f(xx)                                              #取得xx比例处的值

        for iii in range(1, 10):
            B3[a2 - 1, iii - 1] = yy[iii - 1]  #将按比例选取的9个值储存进B3矩阵当中
            
#完成数据的处理

B3 = B3[:, np.nonzero(np.sum(B3, axis=0))[0]]
B4 = B4[:, np.nonzero(np.sum(B4, axis=0))[0]]  #去除矩阵中的多余0的值便于储存数据

df = pd.DataFrame(B3)
writer = pd.ExcelWriter('output1.xlsx', engine='openpyxl')
df.to_excel(writer, sheet_name='Sheet1', index=False)
writer.book.save('output1.xlsx')
df = pd.DataFrame(B4)
writer = pd.ExcelWriter('output2.xlsx', engine='openpyxl')
df.to_excel(writer, sheet_name='Sheet1', index=False)
writer.book.save('output2.xlsx')                          #将数据生成为表格

Here is an updated version of your code:

import pandas as pd
import numpy as np
from scipy.interpolate import interp1d

M0 = pd.read_csv('E:\shujushuchu\shuchu1.csv')  #导入四个参数
M1 = pd.read_csv('E:\shujushuchu\shuru1.csv')   #导入abaqus生成的应力应变曲线的点
B0=M0.values  
B1=M1.values  #将导入到参数转换为矩阵

a1 = int(np.size(B1, 0)/2)  #0表示行数，1表示列数，并赋值给a1
print(B1[1, 1])  #因为py语言里下标从0开始，检查数据的位置是否正确

B3 = np.zeros((a1, len(B1[0])))  #储存曲线点的矩阵
B4 = np.zeros((a1, len(B0[0])))  #储存预测参数的矩阵
                          #这两个矩阵根据前一步的计算结果修改这两个矩阵的大小，a1做为行元素的大小，而列元素大小要大于B1的shape
a2 = 0  #用于循环填充B3

for i in range(a1):  # a1值加1放在这儿，总循环，以行为单位
    B5 = B1[2*i:2*(i+1), :]
    B2 = np.zeros((2, len(B5[0])))  #创建一个空矩阵（该矩阵将用于储存应力应变曲线塑性阶段的点）
    C = B5[1, 2] / B5[0, 2]  #计算弹性模量
    a5 = np.size(B5, 1)  #计算B5矩阵有多少行
    a4 = 2   #用于之后循环跳出
    B7 = np.zeros((2, 3))   #创建空的B7矩阵便于储存条件屈服点（从弹性阶段进入塑性阶段）附近的三个点（用于寻找条件屈服点）
    for ii in range(2, a5):  

        if B5[1, ii - 1] <= C * (B5[0, ii - 1] - 0.002):  #寻找生成的数据中从弹性阶段进入塑性阶段瞬间的后一个点
            a4 = a4 + 1
            if a4 == 3:
                a6 = ii
                B7[0, 2] = B5[0, ii - 1]
                B7[1, 2] = B5[1, ii - 1]
            B2[0, a4 - 1] = B5[0, ii - 1]
            B2[1, a4 - 1] = B5[1, ii - 1]  #这两行循环储存条件屈服点后的实验数据于B2矩阵中


    B7[0, 0] = B5[0, a6 - 3]
    B7[0, 1] = B5[0, a6 - 2]
    B7[1, 0] = B5[1, a6 - 3]
    B7[1, 1] = B5[1, a6 - 2]    #该部分在B7中储存了条件屈服那一瞬间的后一个点和前两个点（用于后面的插值计算）
    B2[0, 1] = B5[0, a6 - 2]
    B2[1, 1] = B5[1, a6 - 2]    #该部分在B2中储存了条件屈服那一瞬间的后一个点
    a11 = (B7[0, 2] - B7[0, 1]) / 500   #将储存的三个点长度细分单位长度，用于插值计算获得中间更多的点
    f1 = interp1d(B7[0, :], B7[1, :], kind='linear')  #定义插值计算的数组、范围以及插值类型
    for i2 in np.arange(B7[0, 1], B7[0, 2], a11):  #利用for循环将数据细分并一一进行插值计算

        a12 = f1(i2)
        a13 = a12
        if a12 <= C * (i2 - 0.002):   #寻找条件屈服点的准确大概值
            B2[0, 1] = i2 - a11         #将条件屈服点应变储存进B2
            B2[1, 1] = a13 - C * a11      #将条件屈服点应力储存进B2
            break

    B6 = np.zeros((2, 10000))       #创建空白矩阵，用于插值计算条件屈服点后的点
    a7 = (np.max(B2[0, :]) - B2[0, 1]) / 500    #将塑性阶段细分单位长度用于插值计算
    a10 = 0    #用于循环递增储存插值后的细分的塑性阶段的点
    B2 = B2[:, np.nonzero(np.sum(B2, axis=0))[0]]    #去除矩阵末尾没储存数据的0值，保证插值计算的正常进行
    f2 = interp1d(B2[0, :], B2[1, :], kind='cubic')    #定义插值计算的数组、范围以及插值类型
    for i3 in np.arange(B2[0, 0], np.max(B2[0, :]) + a7, a7):    #利用for循环将数据细分并一一进行插值计算
        a20 = i3 + a7
        a22 = np.max(B2[0, :])
        if a20 > a22:
            a20 = a22
        a9 = f2(a20)
        a23 = i3
        if a23 > a22:
            a23 = a22
        a8 = f2(a23)       #为了寻找塑性失稳点，差值范围之前设置只到了原始数据的最大值处，插值时可能会超过这个范围而不能正常进行，所以将超过的部分等最大点处的数据
        derivative = (a9 - a8) / a7   #计算曲线的斜率
        if derivative - a9 >= 0:    #在斜率减缓到0之前，将所有的插值点都循环储存到B6矩阵里
            a10 += 1
            B6[0, a10 - 1] = i3
            B6[1, a10 - 1] = a8


    
    B6 = B6[:, np.nonzero(np.sum(B6, axis=0))[0]]     #去除矩阵末尾没储存数据的0值，保证插值计算的正常进行
    df = pd.DataFrame(B6)
    writer = pd.ExcelWriter('output3.xlsx', engine='openpyxl')
    df.to_excel(writer, sheet_name='Sheet1', index=False)
    writer.book.save('output3.xlsx')

    
    if np.max(B6[0, :]) > 0.47:      #将未能插值计算到应力极限点的数据去除掉（不能进行接下来的循环的话将不会储存进最终的矩阵里）
        continue
    else:
        a2 += 1
        for iiii in range(1, 5):
            B4[a2 - 1, iiii - 1] = B0[i - 1, iiii - 1]    #将参数储存进B4矩阵中

        B3[a2 - 1, 9] = np.min(B6[0, :])
        B3[a2 - 1, 10] = np.max(B6[0, :])    #将塑性阶段应变的最大值与最小值储存到B3矩阵中的第9列与第10列
        B6[0, :] = (B6[0, :] - B6[0, 0]) / (max(B6[0, :]) - B6[0, 0])   #将B6中的数据做归一化处理便于之后的插值计算（用于插值时取对应比例处的值）
        
                
        f = interp1d(B6[0, :], B6[1, :], kind='cubic')    #定义插值计算的数组、范围以及插值类型
        xx = np.array([0, 0.03, 0.09, 0.18, 0.3, 0.45, 0.63, 0.84, 1.0])    
        yy = f(xx)                                              #取得xx比例处的值

        for iii in range(1, 10):
            B3[a2 - 1, iii - 1] = yy[iii - 1]  #将按比例选取的9个值储存进B3矩阵当中
            
#完成数据的处理

B3 = B3[:, np.nonzero(np.sum(B3, axis=0))[0]]
B4 = B4[:, np.nonzero(np.sum(B4, axis=0))[0]]  #去除矩阵中的多余0的值便于储存数据

df = pd.DataFrame(B3)
writer = pd.ExcelWriter('output1.xlsx', engine='openpyxl')
df.to_excel(writer, sheet_name='Sheet1', index=False)
writer.book.save('output1.xlsx')
df = pd.DataFrame(B4)
writer = pd.ExcelWriter('output2.xlsx', engine='openpyxl')
df.to_excel(writer, sheet_name='Sheet1', index=False)
writer.book.save('output2.xlsx')                          #将数据生成为表格