UXO_TEM_Widget.py 文件源码

def updatePolarizationsQP(self,r0,UB):

        # Set operator and solution array
        Hp = self.computeHp(r0=r0)
        Brx = self.computeBrx(r0=r0)
        P = self.computeP(Hp,Brx)
        dunc = self.dunc
        dobs = self.dobs

        K = np.shape(dobs)[1]
        q = np.zeros((6,K))

        # Bounds
        ub = UB*np.ones(6)
        e = matrix(np.r_[np.zeros(9),ub])

        # Constraints
        C1 = np.r_[np.c_[-1,0,0,0,0,0],np.c_[0,0,0,-1,0,0],np.c_[0,0,0,0,0,-1]]
        C2 = np.r_[np.c_[-0.5,1,0,-0.5,0,0],np.c_[-0.5,0,1,0,0,-0.5],np.c_[0,0,0,-0.5,1,-0.5]]
        C3 = np.r_[np.c_[-0.5,-1,0,-0.5,0,0],np.c_[-0.5,0,-1,0,0,-0.5],np.c_[0,0,0,-0.5,-1,-0.5]]

        # G = sp.kron(sp.diags([-np.ones(K),np.ones(K)],[0,1],shape=(K-1,K)),sp.diags(np.ones(6)))
        I = np.diag(np.ones(6))
        G = np.r_[C1,C2,C3,I]
        G = matrix(G)

        solvers.options['show_progress'] = False

        for kk in range(0,K):

            LHS = P/np.c_[dunc[:,kk],dunc[:,kk],dunc[:,kk],dunc[:,kk],dunc[:,kk],dunc[:,kk]]
            RHS = dobs[:,kk]/dunc[:,kk]

            A = matrix(np.dot(LHS.T,LHS))
            b = matrix(np.dot(LHS.T,-RHS))

            sol = solvers.qp(A, b, G=G, h=e)

            q[:,kk] = np.reshape(sol['x'],(6))

        self.q = q