matlab把每次循环的结果保存idx,来自Matlab的循环Python
我開始用Python編寫代碼,我來自Matlab背景。我嘗試執(zhí)行的for循環(huán)出現(xiàn)問題。在
這是我在Matlab中的for循環(huán)ix = indoor(1);
idx = indoor(2)-indoor(1);
%Initialize X apply I.C
X = [ix;idx];
for k=(1:1:287)
X(:,k+1) = Abest*X(:,k) + Bbest*outdoor(k+1) + B1best* (cbest4/cbest1);
end
在這個代碼中,Abest是2x2矩陣,Bbest是2x1矩陣,outdoor是288x1向量,B1best是2x1矩陣。矩陣是通過使用matrix expidential命令從函數(shù)中找到的。c4和c1是在常量之前定義的術(shù)語。在
在Python中,我已經(jīng)能夠讓matrix exponential命令在我的函數(shù)中工作,但是我不能讓循環(huán)工作。在
^{pr2}$
同樣,Ae,Be,Be1都是和Matlab一樣大小的矩陣。室外向量也是一樣。在
我想盡一切辦法讓它成功。。。唯一對我有用的是Xo = np.zeros(())
#Initial COnditions
ix = np.array(indoor[0])
idx = np.array(indoor[1]-indoor[0])
Xo = np.array([[ix],[idx]])
#Range for the for loop
num1 = range(0,1)
for k in num1:
Xo = Ae*Xo[k] + Be*outdoor[k+1] + Be1*(c4/c1)
但現(xiàn)在只給我兩點。如果我改變范圍,我會得到一個錯誤。我假設(shè)這段代碼可以工作,因為我最初的Xo只有兩個狀態(tài),所以k會經(jīng)歷這兩個狀態(tài),但這不是我想要的。在
如果有人能幫我,那將是非常有幫助的!如果我犯了一些代碼錯誤,那是因為我沒有很好地理解python中的“For循環(huán)”,當(dāng)涉及到數(shù)據(jù)分析時,讓它遍歷行并增加列。謝謝你抽出時間。在
如有要求,這里是我的完整代碼:import scipy.io as sc
import math as m
import numpy as np
import matplotlib.pyplot as plt
import sys
from scipy.linalg import expm, sinm, cosm
import pandas as pd
df = pd.read_excel('datatemp.xlsx')
outdoor = np.array(df[['Outdoor']])
indoor = np.array(df[['Indoor']])
###########################. FUNCTION DEFINE. #################################################
#Progress bar
def progress(count, total, status=''):
percents = round(100.0 * count / float(total), 1)
sys.stdout.write(' %s%s ...%s\r' % ( percents, '%', status))
sys.stdout.flush()
#Define Matrix for Model
def Matrixbuild(c1,c2,c3):
A = np.array([[0,1],[-c3/c1,-c2/c1]])
B = np.array([[0],[1/c1]])
B1 = np.array([[1],[0]])
C = np.zeros((2,2))
D = np.zeros((2,2))
F = np.array([[0,1,0,1],[-c3/c1,-c2/c1,1/c1,0],[0,0,0,0],[0,0,0,0]])
R = np.array(expm(F))
Ae = np.array([[R.item(0),R.item(1)],[R.item(4),R.item(5)]])
Be = np.array([[R.item(2)],[R.item(6)]])
Be1 = np.array([[R.item(3)],[R.item(7)]])
return Ae,Be,Be1;
###########################. Data. #################################################
#USED FOR JUST TRYING WITHOUT ACTUAL DATA
# outdoor = np.array([5.8115,4.394,5.094,5.1123,5.1224])
# indoor = np.array([15.595,15.2429,15.0867,14.9982,14.8993])
###########################. Model Define. #################################################
Xo = np.zeros((2,288))
ix = np.array(indoor[0])
idx = np.array(indoor[1])
err_min = m.inf
c1spam = np.linspace(0.05,0.001,30)
c2spam = np.linspace(6.2,6.5,30)
c3spam = np.linspace(7.1,7.45,30)
totalspam = len(c1spam)*len(c2spam)*len(c3spam)
ind = 0
for c1 in c1spam:
for c2 in c2spam:
for c3 in c3spam:
c4 = 1.1
#MatrixBuild Function
result = Matrixbuild(c1,c2,c3)
Ae,Be,Be1 = result
Xo = np.array([ix,idx])
Datarange = range(0,len(outdoor)-1,1)
for k in Datarange:
Xo[:,k+1] = np.matmul(Ae,Xo[:,k]) + np.matmul(Be,outdoor[k+1]) + Be1*(c4/c1)
ind = ind + 1
print(Xo)
err = np.linalg.norm(Xo[0,range(0,287)]-indoor.T)
if err
err_min = err
cbest = np.array([[c1],[c2],[c3],[c4]])
progress(ind,totalspam,status='Done')
# print(X)
# print(err)
# print(cbest)
###########################. Model with Cbest Values. #################################################
c1 = cbest[0]
c2 = cbest[1]
c3 = cbest[2]
result2 = Matrixbuild(c1,c2,c3)
AeBest,BeBest,Be1Best = result2
Xo = np.array([ix,idx])
Datarange = np.arange(0,len(outdoor)-1)
for k in Datarange:
Xo[:,k+1] = np.matmul(AeBestb,Xo[:,k]) + np.matmul(BeBest,outdoor[k+1]) + Be1Best*(c4/c1)
err = np.linalg.norm(Xo[0,range(0,287)]-indoor.T)
print(cbest)
print(err)
###########################. Plots. #################################################
plt.figure(0)
time = np.linspace(1,2,2)
plt.scatter(time,X[0],s=15,c="blue")
plt.scatter(time,indoor[0:2],s=15,c="red")
plt.show()
我的錯誤再次出現(xiàn)在for循環(huán)的行中
對于數(shù)據(jù)范圍中的k:Xo[:,k+1] = np.matmul(Ae,Xo[k]) + np.matmul(Be,outdoor[k+1]) + Be1*(c4/c1)
我想用np.材料對于矩陣乘法,但即使沒有它,它也不起作用。在
如果對我的代碼還有其他問題,請詢問。從本質(zhì)上講,我試圖通過一個基本的二階常系數(shù)模型,找到最適合我的數(shù)據(jù)的c1,c2,c3系數(shù),這就是室內(nèi)溫度。在
總結(jié)
以上是生活随笔為你收集整理的matlab把每次循环的结果保存idx,来自Matlab的循环Python的全部內(nèi)容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: php事件和行为,Yii框架组件和事件行
- 下一篇: 举一个通货膨胀的例子