解讀python?cvxpy下SDP問題編程
python cvxpy下SDP問題編程
最近在做定位算法的復(fù)現(xiàn)問題,遇到了
Source Localization in Wireless Sensor Networks From Signal Time-of-Arrival Measurements
里面的一個(gè)半正定優(yōu)化算法,因此選用cvxpy庫實(shí)現(xiàn)。
官方文檔[cvxpy]的例程復(fù)現(xiàn)的算法都很簡單,因此對(duì)該問題的借鑒意義不大。
算法如下
對(duì)我而言,首先的難度就是拼接矩陣后的半正定約束條件,起初是另設(shè)立兩個(gè)矩陣變量,然后按部就班的增加限制條件。但最后求得的數(shù)據(jù)千奇百怪,與預(yù)測(cè)位置沒有任何關(guān)系。
后來不斷嘗試更改約束限制的表達(dá)形式,但均無效果。
后來輸出了每個(gè)變量的值查看,發(fā)現(xiàn)Q元素的物理意義為預(yù)測(cè)距離的平方,但是求出來的Q矩陣元素往往極大,因此擅自添加了一個(gè)約束條件,限制Q的最大元素在預(yù)測(cè)距離平方的量級(jí)上,完美解決問題。
附上代碼
class Program_t: def __init__(self,bt): self.BT = bt self.BT_x = [b[0] for b in self.BT] self.BT_y = [b[1] for b in self.BT] self.T=[b[2] for b in self.BT] self.number = len(bt) def LS_steps(self): num = len(self.T) up_control = 2*max(self.T)**2#限制最大元素量級(jí) Q = cp.Variable((num,num))#待求變量 tao = cp.Variable((num,1))#生成矩陣形式后面才可以拼接 y_ = cp.Variable((2,1)) y_s = cp.Variable((1,1))#矩陣形式用于拼接 yita = 0.000005*sum(self.T) / num#論文給出的參數(shù)選擇,可更改常數(shù) G = np.eye(num)-np.ones((num,num)) t = np.array([self.T]).T expr1 = cp.trace((cp.transpose(G)) @ G @ (Q- cp.multiply(2,t @ (cp.transpose(tao)))+t @ (cp.transpose(t)))) expr2 = yita*cp.sum(Q) expr = expr1+ expr2#目標(biāo)函數(shù) Q_ = cp.bmat([[Q,tao],[cp.transpose(tao),[[1]]]])#拼接矩陣 Y = cp.bmat([[np.eye(2),y_],[cp.transpose(y_),y_s]])#拼接矩陣 constraints = [Q_ >> 0, Y >> 0, cp.max(Q)<=up_control]#限制條件Q半正定,Y半正定,Q最大元素小于上限(這個(gè)約束非常重要,是我自己加上去的) for i in range(num): X = np.array([self.BT_x[i],self.BT_y[i],-1]).T constraints += [Q[i, i] == cp.transpose(X) @ Y @ X]#約束條件 for j in range(i+1,num): X_j = np.array([self.BT_x[j], self.BT_y[j],-1]).T constraints += [Q[i, j] >= cp.abs(cp.transpose(X) @ Y @ X_j)]#約束條件 obj = cp.Minimize(expr) prob = cp.Problem(obj, constraints) prob.solve() position = y_.value print(expr1.value)#輸出值 print(expr2.value) print(prob.value)#輸出值 print(prob.status)#輸出狀態(tài) print(position) return position
總結(jié)
1.理論算法與編程實(shí)現(xiàn)永遠(yuǎn)不等,不能輕易照搬,具體實(shí)現(xiàn)過程中要結(jié)合實(shí)際情況進(jìn)行考慮,當(dāng)求得的結(jié)果與預(yù)計(jì)相差很多時(shí),可以嘗試增加數(shù)值約束,因?yàn)橛?jì)算機(jī)仿真只是近似,不是理論上的完美條件。
2.編程實(shí)現(xiàn)調(diào)用庫時(shí),最好按照庫的標(biāo)準(zhǔn)寫,如本例中矩陣點(diǎn)乘可以用numpy 的dot或者cvxpy的@,以及轉(zhuǎn)置的.T和cp.transpose().但是dot有時(shí)會(huì)產(chǎn)生意想不到的情況,平白增加工作量。
3.復(fù)現(xiàn)算法時(shí)必須要對(duì)算法有深入理解,否則難以發(fā)現(xiàn)問題所在。
4.不要輕易懷疑工具包的問題,經(jīng)過大量使用的工具包一定比你的感覺可靠。
以上為個(gè)人經(jīng)驗(yàn),希望能給大家一個(gè)參考,也希望大家多多支持腳本之家。
相關(guān)文章
Python之time模塊的時(shí)間戳,時(shí)間字符串格式化與轉(zhuǎn)換方法(13位時(shí)間戳)
今天小編就為大家分享一篇Python之time模塊的時(shí)間戳,時(shí)間字符串格式化與轉(zhuǎn)換方法(13位時(shí)間戳),具有很好的參考價(jià)值,希望對(duì)大家有所幫助。一起跟隨小編過來看看吧2019-08-08PyMongo進(jìn)行MongoDB查詢和插入操作的高效使用示例
這篇文章主要為大家介紹了PyMongo進(jìn)行MongoDB查詢和插入操作的高效使用示例,有需要的朋友可以借鑒參考下,希望能夠有所幫助,祝大家多多進(jìn)步,早日升職加薪2023-11-11