Skip to content

Commit

Permalink
Issue #724 while loop instead of for loop
Browse files Browse the repository at this point in the history
  • Loading branch information
jaeandersson committed Feb 20, 2018
1 parent 5e0b393 commit 8521cc3
Showing 1 changed file with 37 additions and 41 deletions.
78 changes: 37 additions & 41 deletions casadi/solvers/conic_as.cpp
Expand Up @@ -334,11 +334,8 @@ namespace casadi {
const casadi_int* a_row = A_.row();

// QP iterations

// QP iterations
casadi_int iter;
for (iter=0; iter<max_iter_; ++iter) {

casadi_int iter = 0;
while (true) {
// Debugging
if (verbose_) {
print("Current xk = \n");
Expand All @@ -351,39 +348,18 @@ namespace casadi {
print_vector(lam_ak, na_);
}

// Copy kkt to kktd
casadi_project(kkt, kkt_, kktd, kktd_, w);

// Loop over kktd entries (left two blocks of the transposed KKT)
for (casadi_int c=0; c<nx_; ++c) {
if (lam_xk[c]!=0) {
// Zero out column, set diagonal entry to 1
for (casadi_int k=kkt_colind[c]; k<kkt_colind[c+1]; ++k) {
kktd[k] = kkt_row[k]==c ? 1. : 0.;
}
}
}

// Loop over kktd entries (right two blocks of the transposed KKT)
for (casadi_int c=0; c<na_; ++c) {
if (lam_ak[c]==0) {
// Zero out column, set diagonal entry to -1
for (casadi_int k=kkt_colind[nx_+c]; k<kkt_colind[nx_+c+1]; ++k) {
kktd[k] = kkt_row[k]==nx_+c ? -1. : 0.;
}
}
}

// QR factorization
casadi_qr(kktd_, kktd, w, sp_v_, v, sp_r_, r, beta, get_ptr(prinv_), get_ptr(pc_));

// Evaluate gradient of the Lagrangian and constraint functions
casadi_copy(g, nx_, step);
casadi_mv(h, H_, xk, step, 0); // gradient of the objective
casadi_mv(a, A_, lam_ak, step, 1); // gradient of the Lagrangian
casadi_copy(gk, na_, step + nx_); // constraint evaluation

// Correct for active simple bounds
// Start new iteration
if (++iter==max_iter_) {
casadi_warning("Maximum number of iterations reached");
break;
}

// Calculate KKT residual: Correct for active simple bounds
for (i=0; i<nx_; ++i) {
if (lam_xk[i]!=0.) {
step[i] = xk[i];
Expand All @@ -392,7 +368,8 @@ namespace casadi {
}
}

// Correct for inactive constraints
// Calculate KKT residual: Correct for inactive constraints
casadi_copy(gk, na_, step + nx_); // constraint evaluation
for (i=0; i<na_; ++i) {
if (lam_ak[i]==0) {
step[nx_+i] = 0.; // -lam_ak[i]
Expand All @@ -408,10 +385,34 @@ namespace casadi {
print_vector(step, nx_ + na_);
}

// Negative residual
casadi_scal(nx_+na_, -1., step);
// Copy kkt to kktd
casadi_project(kkt, kkt_, kktd, kktd_, w);

// Loop over kktd entries (left two blocks of the transposed KKT)
for (casadi_int c=0; c<nx_; ++c) {
if (lam_xk[c]!=0) {
// Zero out column, set diagonal entry to 1
for (casadi_int k=kkt_colind[c]; k<kkt_colind[c+1]; ++k) {
kktd[k] = kkt_row[k]==c ? 1. : 0.;
}
}
}

// Loop over kktd entries (right two blocks of the transposed KKT)
for (casadi_int c=0; c<na_; ++c) {
if (lam_ak[c]==0) {
// Zero out column, set diagonal entry to -1
for (casadi_int k=kkt_colind[nx_+c]; k<kkt_colind[nx_+c+1]; ++k) {
kktd[k] = kkt_row[k]==nx_+c ? -1. : 0.;
}
}
}

// QR factorization
casadi_qr(kktd_, kktd, w, sp_v_, v, sp_r_, r, beta, get_ptr(prinv_), get_ptr(pc_));

// Solve to get primal-dual step
casadi_scal(nx_+na_, -1., step);
casadi_qr_solve(step, 1, 1, sp_v_, v, sp_r_, r, beta,
get_ptr(prinv_), get_ptr(pc_), w);

Expand Down Expand Up @@ -693,11 +694,6 @@ namespace casadi {
break;
}

// Maximum number of iterations reached
if (iter==max_iter_) {
casadi_warning("Maximum number of iterations reached");
}

// Calculate optimal cost
if (f) *f = fk;

Expand Down

0 comments on commit 8521cc3

Please sign in to comment.