How to use the casadi.MX.zeros function in casadi

To help you get started, we’ve selected a few casadi examples, based on popular ways it is used in public projects.

Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.

github helgeanl / GP-MPC / gp_mpc / gp_functions.py View on Github external
X: Training data matrix with inputs of size NxNx, with Nx number of
            inputs to the GP.
        Y: Training data matrix with outpyts of size (N x Ny).
        hyper: Array with hyperparameters [ell_1 .. ell_Nx sf sn].
        inputmean: Mean from the last GP iteration of size (1 x Nx)
        inputcov: Covariance matrix from the last GP iteration of size (Nx x Nx)

    # Returns
        mean: Array of the output mean of size (Ny x 1).
        covariance: Covariance matrix of size (Ny x Ny).
    """

    hyper = ca.log(hyper)
    Ny     = len(invK)
    N, Nx     = ca.MX.size(X)
    mean  = ca.MX.zeros(Ny, 1)
    beta  = ca.MX.zeros(N, Ny)
    log_k = ca.MX.zeros(N, Ny)
    v     = X - ca.repmat(inputmean, N, 1)

    covariance = ca.MX.zeros(Ny, Ny)

    #TODO: Fix that LinsolQr don't work with the extended graph?
    A = ca.SX.sym('A', inputcov.shape)
    [Q, R2] = ca.qr(A)
    determinant = ca.Function('determinant', [A], [ca.exp(ca.trace(ca.log(R2)))])

    for a in range(Ny):
        beta[:, a] = ca.mtimes(invK[a], Y[:, a])
        iLambda   = ca.diag(ca.exp(-2 * hyper[a, :Nx]))
        R  = inputcov + ca.diag(ca.exp(2 * hyper[a, :Nx]))
        iR = ca.mtimes(iLambda, (ca.MX.eye(Nx) - ca.solve((ca.MX.eye(Nx)
github meco-group / omg-tools / omgtools / problems / dualdecomposition.py View on Github external
if sym.name() in child.symbol_dict:
                        name = child.symbol_dict[sym.name()][1]
                        v = x_i[label][name]
                        ind = self.q_i[child][name]
                        sym2 = MX.zeros(sym.size())
                        sym2[ind] = v
                        sym2 = reshape(sym2, sym.shape)
                        c = substitute(c, sym, sym2)
                        break
                for nghb in self.q_ij.keys():
                    for label, child in nghb.group.items():
                        if sym.name() in child.symbol_dict:
                            name = child.symbol_dict[sym.name()][1]
                            v = z_ij[nghb.label, label, name]
                            ind = self.q_ij[nghb][child][name]
                            sym2 = MX.zeros(sym.size())
                            sym2[ind] = v
                            sym2 = reshape(sym2, sym.shape)
                            c = substitute(c, sym, sym2)
                            break
                for name, s in self.par_global.items():
                    if s.name() == sym.name():
                        c = substitute(c, sym, par[name])
            lb, ub = con[1], con[2]
            self.define_constraint(c, lb, ub)
        # construct problem
        prob, buildtime = self.father_updx.construct_problem(
            self.options, str(self._index), problem)
        self.problem_upd_xz = prob
        self.father_updx.init_transformations(self.problem.init_primal_transform,
                                         self.problem.init_dual_transform)
        self.init_var_dd()
github helgeanl / GP-MPC / GP / prediction.py View on Github external
def gp_casadi(invK, hyp, X, Y, z):
    E = len(invK)
    n = ca.MX.size(X[:, 1])[0]
    D = ca.MX.size(X[1, :])[1]

    mean  = ca.MX.zeros(E, 1)
    var  = ca.MX.zeros(E, 1)
    for a in range(E):
        ell = ca.MX(hyp[a, 0:D])
        sf2 = ca.MX(hyp[a, D]**2)
        m = 0 #hyp[a, D + 2]
        kss = covSEard2(z, z, ell, sf2)
        ks = ca.MX.zeros(n, 1)

        for i in range(n):
            ks[i] = covSEard2(X[i, :] - m, z - m, ell, sf2)
        #ks = repmat()
        ksK = ca.mtimes(ks.T, invK[a])

        mean[a] = ca.mtimes(ksK, Y[:, a] - m) + m
        var[a] = kss - ca.mtimes(ksK, ks)

    return mean, var
github helgeanl / GP-MPC / GP / prediction.py View on Github external
def calc_NLL_casadi(hyper, X, Y):
    """ Objective function """
    # Calculate NLL
    n, D = ca.MX.size(X)

    ell = hyper[:D]
    sf2 = hyper[D]**2
    lik = hyper[D + 1]**2
    m   = hyper[D + 2]
    K   = ca.MX.zeros(n, n)

    for i in range(n):
        for j in range(n):
            K[i, j] = covSEard2(X[i, :] - m, X[j, :] - m, ell, sf2)

    K = K + lik * ca.MX.eye(n)
    K = (K + K.T) * 0.5   # Make sure matrix is symmentric

    A = ca.SX.sym('A', ca.MX.size(K))
    #L = ca.chol(A)      # Should check for PD !!!
    cholesky = ca.Function('Cholesky', [A], [ca.chol(A)])
    L = cholesky(K).T

    logK = 2 * ca.sum1(ca.MX.log(ca.fabs(ca.diag(L))))

    invLy = ca.solve(L, Y - m)
github helgeanl / GP-MPC / GP / prediction.py View on Github external
def gp_taylor_approx(invK, X, F, hyper, D, inputmean, inputcov):
    #hyper = ca.MX.log(hyper)
    E     = len(invK)
    n     = ca.MX.size(F[:, 1])[0]
    mean  = ca.MX.zeros(E, 1)
    var = ca.MX.zeros(D, 1)
    v     = X - ca.repmat(inputmean, n, 1)
    covar = ca.MX.zeros(E, E)
    covariance = ca.MX.zeros(E, E)
    d_mean = ca.MX.zeros(E, 1)
    dd_var = ca.MX.zeros(E, E)

    for a in range(E):
        ell = hyper[a, :D]
        w = 1 / ell**2
        sf2 = ca.MX(hyper[a, D]**2)
        iK = ca.MX(invK[a])
        alpha = ca.mtimes(iK, F[:, a])
        kss = sf2
        ks = ca.MX.zeros(n, 1)
github helgeanl / GP-MPC / GP / noisyGP.py View on Github external
def GP_noisy_input(invK, X, F, hyper, D, inputmean, inputcov):
    hyper = ca.MX.log(hyper)
    E     = len(invK)
    n     = ca.MX.size(F[:, 1])[0]
    mean  = ca.MX.zeros(E, 1)
    beta  = ca.MX.zeros(n, E)
    log_k = ca.MX.zeros(n, E)
    v     = X - ca.repmat(inputmean, n, 1)

    #invK = MX(invK)
    covariance = ca.MX.zeros(E, E)

    A = ca.SX.sym('A', inputcov.shape)
    [Q, R2] = ca.qr(A)
    determinant = ca.Function('determinant', [A], [ca.exp(ca.trace(ca.log(R2)))])

    for a in range(E):
        beta[:, a] = ca.mtimes(invK[a], F[:, a])
        iLambda   = ca.diag(ca.exp(-2 * hyper[a, :D]))
        R  = inputcov + ca.diag(ca.exp(2 * hyper[a, :D]))
        iR = ca.mtimes(iLambda, (ca.MX.eye(D) - ca.solve((ca.MX.eye(D) + ca.mtimes(inputcov, iLambda)), (ca.mtimes(inputcov, iLambda)))))
github helgeanl / GP-MPC / GP / prediction.py View on Github external
def gp_taylor_approx(invK, X, F, hyper, D, inputmean, inputcov):
    #hyper = ca.MX.log(hyper)
    E     = len(invK)
    n     = ca.MX.size(F[:, 1])[0]
    mean  = ca.MX.zeros(E, 1)
    var = ca.MX.zeros(D, 1)
    v     = X - ca.repmat(inputmean, n, 1)
    covar = ca.MX.zeros(E, E)
    covariance = ca.MX.zeros(E, E)
    d_mean = ca.MX.zeros(E, 1)
    dd_var = ca.MX.zeros(E, E)

    for a in range(E):
        ell = hyper[a, :D]
        w = 1 / ell**2
        sf2 = ca.MX(hyper[a, D]**2)
        iK = ca.MX(invK[a])
        alpha = ca.mtimes(iK, F[:, a])
        kss = sf2
        ks = ca.MX.zeros(n, 1)

        for i in range(n):
            ks[i] = covSEard2(X[i, :], inputmean, ell, sf2)
github meco-group / omg-tools / omgtools / problems / dualdecomposition.py View on Github external
for child in q_j.keys():
                    for name in q_j[child].keys():
                        z = z_ij[str(nghb), child.label, name]
                        l = l_ij[str(nghb), child.label, name]
                        obj -= mtimes(l.T, z)
            self.define_objective(obj)
        # construct constraints
        for con in self.constraints:
            c = con[0]
            for sym in symvar(c):
                for label, child in self.group.items():
                    if sym.name() in child.symbol_dict:
                        name = child.symbol_dict[sym.name()][1]
                        v = x_i[label][name]
                        ind = self.q_i[child][name]
                        sym2 = MX.zeros(sym.size())
                        sym2[ind] = v
                        sym2 = reshape(sym2, sym.shape)
                        c = substitute(c, sym, sym2)
                        break
                for nghb in self.q_ij.keys():
                    for label, child in nghb.group.items():
                        if sym.name() in child.symbol_dict:
                            name = child.symbol_dict[sym.name()][1]
                            v = z_ij[nghb.label, label, name]
                            ind = self.q_ij[nghb][child][name]
                            sym2 = MX.zeros(sym.size())
                            sym2[ind] = v
                            sym2 = reshape(sym2, sym.shape)
                            c = substitute(c, sym, sym2)
                            break
                for name, s in self.par_i.items():
github helgeanl / GP-MPC / GP / noisyGP.py View on Github external
def GP_noisy_input(invK, X, F, hyper, D, inputmean, inputcov):
    hyper = ca.MX.log(hyper)
    E     = len(invK)
    n     = ca.MX.size(F[:, 1])[0]
    mean  = ca.MX.zeros(E, 1)
    beta  = ca.MX.zeros(n, E)
    log_k = ca.MX.zeros(n, E)
    v     = X - ca.repmat(inputmean, n, 1)

    #invK = MX(invK)
    covariance = ca.MX.zeros(E, E)

    A = ca.SX.sym('A', inputcov.shape)
    [Q, R2] = ca.qr(A)
    determinant = ca.Function('determinant', [A], [ca.exp(ca.trace(ca.log(R2)))])

    for a in range(E):
        beta[:, a] = ca.mtimes(invK[a], F[:, a])
        iLambda   = ca.diag(ca.exp(-2 * hyper[a, :D]))
        R  = inputcov + ca.diag(ca.exp(2 * hyper[a, :D]))
        iR = ca.mtimes(iLambda, (ca.MX.eye(D) - ca.solve((ca.MX.eye(D) + ca.mtimes(inputcov, iLambda)), (ca.mtimes(inputcov, iLambda)))))
        T  = ca.mtimes(v, iR)
        c  = ca.exp(2 * hyper[a, D]) / ca.sqrt(determinant(R)) * ca.exp(ca.sum2(hyper[a, :D]))