代码如下，该怎么读取优化出来的个变量的结果
import pandas as pd
import math
import matplotlib.pyplot as plt
import numpy as np
import scipy.optimize as optimize
import pickle
from joblib import Parallel, delayed
from sko.PSO import PSO
class crm_model():
def __init__(self,inj,pro,pre,time,dt):
self.inj = inj
self.pro = pro
self.pre = pre
self.time = time
self.dt = dt
def init_para(self,n1,n2):
tao = np.random.rand(n2)
gain = np.random.rand(n1,n2)
p_j = np.random.rand(n2)
We = np.random.rand(n2)*10000
dt=1
x0 = np.hstack((tao,gain.reshape(1,-1)[0],p_j,We))
return x0
def boundary(self,n1,n2,num):
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,0),)*n2
if num == 0:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,0),)*n2
elif num == 1:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,np.inf),)*n2
elif num == 2:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,np.inf),)*n2
elif num == 3:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,0),)*n2
elif num == 4:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,0),)*n2
elif num == 5:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,0),)*n2
elif num == 6:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,np.inf),)*n2
elif num == 7:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,0),)*n2
elif num == 8:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,np.inf),)*n2
elif num == 9:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,np.inf),)*n2
elif num == 10:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,np.inf),)*n2
elif num == 11:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,np.inf),)*n2
elif num == 12:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,0),)*n2
elif num == 13:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,np.inf),)*n2
elif num == 14:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,np.inf),)*n2
else:
raise('the num is wrong')
bound1 = list(bound1)
return bound1
def rever_para(self,x0,n1,n2):
tao = x0[:n2]
gain = x0[n2:n2+n2*n1]
gain = gain.reshape([n1,n2])
p_j = x0[n2+n2*n1:n2+n2*n1+n2]
We = x0[n2+n2*n1+n2:n2+n2*n1+n2+n2]
return tao,gain,p_j,We
def cal_q_hat(self,x0):
inj = self.inj
pro = self.pro
pre = self.pre
time = self.time
dt = self.dt
q_hat = np.zeros([1,len(pro)])
q_hat2 = 0
tao,gain,p_j,We = self.rever_para(x0,1,1)
for j in range(1): #生产井数
for tn in range(len(pro )): #时间天数
q_hat1 = pro[j][0]*math.exp(-(time[tn]-time[0])/tao[j])
for k in range(1,tn+1):
summ = 0
#summ = summ + gain1[j]*inj[1][k]+gain2[j]*inj[2][k]+gain3[j]*inj[3][k]
for i in range(len(inj)):
summ = summ + gain[i][j]*inj[i][k]
S1 = summ - p_j[j]*tao[j]*(pre[j][k]-pre[j][k-1])/dt + We[j]
q_hat2 += math.exp(-(time[tn]-time[k])/tao[j]) * (1-math.exp(dt)/tao[j])*S1
q_hat[j][tn] = q_hat1 + q_hat2
return q_hat
def mse(self,x0):
inj = self.inj
pro = self.pro
pre = self.pre
time = self.time
q_hat = self.cal_q_hat(x0)
error = sum(sum((q_hat-pro)**2))
return error
def fit(self):
boundry
lb = []
ub = []
for i,j in boundry:
lb.append(i)
ub.append(j)
result = optimize.minimize(self.mse,x0,method='SLSQP',bounds=boundry,options={'maxiter':5,'disp':True})
# pso = PSO(func=self.mse,n_dim=4,lb = lb,ub = ub)
a = 1
return result
opts_perwell = [opts(r['x']) for r in result]
## read data
inj = pd.read_csv('./data/injection.csv').values
pro = pd.read_csv('./data/production.csv').values
pre = pd.read_csv('./data/well_pressure.csv').values
time = pd.read_csv('./data/time.csv').values
n=len(time)
#dp=np.zeros((proj.shape[1],n))
inj=inj.transpose()
pro=pro.transpose()
pre = pre.transpose()
dt = 1
res = []
crm = crm_model(inj, pro, pre, time, dt)
for inj_i in inj:
ls = []
for num,(pro_j,pre_j) in enumerate(zip(pro,pre)):
inj_i = np.expand_dims(inj_i,0)
pro_j = np.expand_dims(pro_j,0)
pre_j = np.expand_dims(pre_j,0)
crm = crm_model(inj_i, pro_j, pre_j, time, dt)
x0 = crm.init_para(len(inj_i),len(pro_j))
boundry = crm.boundary(len(inj_i),len(pro_j),num)
crm.cal_q_hat(x0)
result= crm.fit()
ls.append(result)
a = 1
res.append(ls)
print(res)
看不懂，但是你boundary这个函数写的。。。，一共16种情况写16个if...elif...else，大量重复。。。，建议你改一改，然后加一些备注，哪些是干什么的，要做什么，因为是真的看不懂
if num in [0,3,4,5,7,12]:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,0),)*n2
elif num in [1,2,6,8,9,10,11,13,14]:
bound1 = ((0,np.inf),)*n2 + ((0,1),)*n2+ ((0,np.inf),)*n2+ ((0,np.inf),)*n2
else:
raise('the num is wrong')
fit函数return result之后的操作是无效的，把 opts_perwell = [opts(r['x']) for r in result] 放到return之前
您好，我是有问必答小助手，您的问题已经有小伙伴解答了，您看下是否解决，可以追评进行沟通哦~
如果有您比较满意的答案 / 帮您提供解决思路的答案，可以点击【采纳】按钮，给回答的小伙伴一些鼓励哦～～
ps:问答VIP仅需29元，即可享受5次/月 有问必答服务，了解详情>>>https://vip.csdn.net/askvip?utm_source=1146287632