DuQuad  v1.0
Quadratic Programming Optimizations
 All Data Structures Files Functions Variables Typedefs Macros
dfgm.c
Go to the documentation of this file.
1 /*
2  * DFGM.c
3  *
4  * Created on: Sep 19, 2014
5  * Author: sverre
6  */
7 
8 #include "dfgm.h"
9 
10 /* static functions declaration */
11 
12 static int32_t solve_DFGM();
13 static real_t dual_obj();
14 static void clean_up_inner_problem();
15 static void init_inner_problem();
16 
17 
18 /* public functions definition */
19 
20 int32_t DFGM(struct Struct_DFGM *s)
21 {
22  _SDEBUG("I'm in DFGM()\n");
23 
24  // Make z0 feasible
25  if(!s->info->ub_is_inf)
26  vector_min(s->prob->ub,s->prob->z0,s->prob->z0,N);
27  if(!s->info->lb_is_inf)
28  vector_max(s->prob->lb,s->prob->z0,s->prob->z0,N);
29 
30 
31  // Initialize the "inner" problem
32  struct Struct_FGM p_in;
33  init_inner_problem(s, &p_in);
34  if (solve_DFGM(s,&p_in) == -1){
35  ERROR("Error in solving solveDFGM()\n");
36  clean_up_inner_problem(&p_in);
37  return -1;
38  };
39  clean_up_inner_problem(&p_in);
40  _SDEBUG("End of DFGM()\n");
41  return 0;
42 }
43 
44 static int32_t solve_DFGM(struct Struct_DFGM *s, struct Struct_FGM *p_in)
45 {
46  _DEBUG2("problem_case: %d\n", s->info->problem_case);
47  // calculating alpha
48  real_t alpha = 1.0 / (2.0 * s->info->Ld);
49  _DEBUG("alpha %0.10f\n", alpha);
50 
51  // making b+lb_hat and b+ub_hat, because this is an constant operation
52  const uint32_t prob_case = s->info->problem_case;
53  switch (prob_case){
54  case 1:
55  vector_add(s->prob->b,s->prob->ub_hat,s->b_ub_hat,M);
56  vector_add(s->prob->b,s->prob->lb_hat,s->b_lb_hat,M);
57  break;
58  case 2:
59  vector_add(s->prob->b,s->prob->ub_hat,s->b_ub_hat,M);
60  break;
61  case 3:
62  vector_add(s->prob->b,s->prob->ub_hat,s->b_ub_hat,M);
63  break;
64  case 4:
65  vector_add(s->prob->b,s->prob->lb_hat,s->b_lb_hat,M);
66  break;
67  default:
68  ERROR("Error deciding the problem case\n");
69  return -1;
70  }
71 
72 
73  // Dynamc Arrays for ds and pf
74  struct Array ds_array, pf_array;
75  initArray(&ds_array, 100);
76  initArray(&pf_array, 100);
77 
78  // Declaration and initialization of other necessary variables
79  real_t pf = s->opt->eps_pf + 1;
80  uint32_t niter_feasible_ds = 0;
81  uint32_t niter_feasible_pf = 0;
82  uint32_t last_eps_ds = 0;
83  uint32_t niter = 1;
84  uint32_t niter_inner = 0;
85  real_t dual_value_new = 0.0;
86  real_t dual_value_diff = s->opt->eps_ds + 1;
87  s->res->exitflag = 1;
88  s->res->out->exitflag_inner = 1;
89  s->res->out->num_exceeded_max_niter_inner = 0;
90  // new:
91  real_t theta = 1.0;
92  real_t theta_new = 0.0;
93  real_t beta = 0.0;
94  real_t S = theta;
95 
96  // Clock
97  clock_t tic, toc, tic_in, toc_in;
98  tic = clock();
99 
100  // solve inner problem first time
101  // with: c_in = c + A'*y1 - A'*y2 = c + A'(y1-y2) = c + A'(lambda1-lambda2)
102  // NOTE: As long as lambda1 = lambda2 = y1 = y2 0, c_in = c
103  vector_copy(s->prob->c,p_in->c,N);
104  vector_copy(s->prob->z0,p_in->z0,N); // Initialize z0 ?????????????
105 
106  tic_in = clock();
107  niter_inner += FGM(p_in);
108  toc_in = clock();
109 
110  s->res->out->time_tot_inner = (real_t)(toc_in-tic_in);
111  s->time_inner_y = (real_t)(toc_in-tic_in);
112  if(p_in->exitflag == 2){
113  s->res->out->exitflag_inner = 2;
114  s->res->out->num_exceeded_max_niter_inner++;
115  }
116  vector_copy(p_in->zopt,s->z,N);
117  vector_copy(p_in->zopt,s->z_ds,N);
118  mtx_vec_mul(s->prob->A,s->z_ds,s->A_z_ds,M,N); // used in dual_obj
119  vector_copy(s->z,s->summ,N);
120  vector_scalar_mul(s->summ,(theta/S),s->summ,N);
121 
122  // find the value of the dual function (Lagrangian) first time
123  real_t dual_value = dual_obj(s,prob_case);
124 
125  /* ##################################################################################
126  START WHILE-LOOP
127  ################################################################################## */
128 
129  while (dual_value_diff > s->opt->eps_ds || pf > s->opt->eps_pf)
130  {
131 
132  if (niter > s->opt->maxiter_outer){
133  //printf("reached maximum number of iterations in DFGM\n");
134  niter_feasible_ds--;
135  niter_feasible_pf--;
136  s->res->exitflag = 2;
137  break;
138  }
139  if (dual_value_diff < s->opt->eps_ds){
140  niter_feasible_ds++;
141  last_eps_ds = 1;
142  }
143  else{
144  last_eps_ds = 0;
145  }
146  if (pf < s->opt->eps_pf)
147  niter_feasible_pf++;
148 
149  /* *********************************************************************
150  Finding the next lambda
151  ********************************************************************* */
152 
153  // lambda1 += alpha * (A*z - b_lb_hat)
154  // NOTE: re-use the 'lambda1' vector instead of creating 'lambda1_new'
155 
156  mtx_vec_mul(s->prob->A,s->z,s->A_z,M,N); // A*z
157 
158  switch (prob_case){
159  case 1:
160  // lambda1 = y1 + alpha*(A*z-b_ub_hat)
161  // NOTE: re-use the 'lambda1' vector instead of creating 'lambda1_new' etc
162  vector_sub(s->A_z,s->b_ub_hat,s->temp3_dim_M,M); // ANS - b_ub_hat
163  vector_scalar_mul(s->temp3_dim_M,alpha,s->temp3_dim_M,M); // ANS * alpha
164  vector_add(s->temp3_dim_M,s->y1,s->lambda1,M); // lambda1 = y1 + ANS
165  vector_max_with_zero(s->lambda1, M); // project lambda1 on R+
166 
167  // lambda2 = y2 + alpha*(-A*z+b_lb_hat)
168  vector_sub(s->b_lb_hat,s->A_z,s->temp3_dim_M,M); // -(A*z) + b_lb_hat
169  vector_scalar_mul(s->temp3_dim_M,alpha,s->temp3_dim_M,M); // ANS * alpha
170  vector_add(s->temp3_dim_M,s->y2,s->lambda2,M); // lambda2 = y2 + ANS
171  vector_max_with_zero(s->lambda2, M); // project lambda2 on R+
172  break;
173  case 2:
174  // lambda1 = y1 + alpha*(A*z-b_ub_hat)
175  vector_sub(s->A_z,s->b_ub_hat,s->temp3_dim_M,M); // ANS - b_ub_hat
176  vector_scalar_mul(s->temp3_dim_M,alpha,s->temp3_dim_M,M); // ANS * alpha
177  vector_add(s->temp3_dim_M,s->y1,s->lambda1,M); // lambda1 = y1 + ANS
178  // NOTE: No projection
179  break;
180  case 3:
181  // lambda1 = y1 + alpha*(A*z-b_ub_hat)
182  vector_sub(s->A_z,s->b_ub_hat,s->temp3_dim_M,M); // ANS - b_ub_hat
183  vector_scalar_mul(s->temp3_dim_M,alpha,s->temp3_dim_M,M); // ANS * alpha
184  vector_add(s->temp3_dim_M,s->y1,s->lambda1,M); // lambda1 = y1 + ANS
185  vector_max_with_zero(s->lambda1, M); // project lambda1 on R+
186  break;
187  case 4:
188  // lambda2 = y2 + alpha*(-A*z+b_lb_hat)
189  vector_sub(s->b_lb_hat,s->A_z,s->temp3_dim_M,M); // -(A*z) + b_lb_hat
190  vector_scalar_mul(s->temp3_dim_M,alpha,s->temp3_dim_M,M); // ANS * alpha
191  vector_add(s->temp3_dim_M,s->y2,s->lambda2,M); // lambda2 = y2 + ANS
192  vector_max_with_zero(s->lambda2, M); // project lambda2 on R+
193  break;
194  default:
195  ERROR("Error deciding the problem case\n");
196  return -1;
197  }
198 
199  /* *********************************************************************
200  Finding the next y
201  ********************************************************************* */
202 
203  // y = lambda + beta * (lambda-lambda_old)
204 
205  // calculating beta
206  theta_new = 0.5 * (1.0 + sqrt(1.0 + 4.0*(theta*theta)));
207  beta = (theta-1.0)/theta_new;
208 
209  switch (prob_case){
210  case 1:
211  // y1
212  vector_sub(s->lambda1,s->lambda1_old,s->temp2_dim_M,M); // lambda1 - lambda1_old
213  vector_scalar_mul(s->temp2_dim_M,beta,s->temp2_dim_M,M); // ANS * beta
214  vector_add(s->lambda1,s->temp2_dim_M,s->y1,M); // y1 = ANS + lambda1
215  //y2
216  vector_sub(s->lambda2,s->lambda2_old,s->temp2_dim_M,M); // lambda2 - lambda2_old
217  vector_scalar_mul(s->temp2_dim_M,beta,s->temp2_dim_M,M); // ANS * beta
218  vector_add(s->lambda2,s->temp2_dim_M,s->y2,M); // y1 = ANS + lambda2
219  break;
220  case 2:
221  // y1
222  vector_sub(s->lambda1,s->lambda1_old,s->temp2_dim_M,M); // lambda1 - lambda1_old
223  vector_scalar_mul(s->temp2_dim_M,beta,s->temp2_dim_M,M); // ANS * beta
224  vector_add(s->lambda1,s->temp2_dim_M,s->y1,M); // y1 = ANS + lambda1
225  break;
226  case 3:
227  // y1
228  vector_sub(s->lambda1,s->lambda1_old,s->temp2_dim_M,M); // lambda1 - lambda1_old
229  vector_scalar_mul(s->temp2_dim_M,beta,s->temp2_dim_M,M); // ANS * beta
230  vector_add(s->lambda1,s->temp2_dim_M,s->y1,M); // y1 = ANS + lambda1
231  break;
232  case 4:
233  //y2
234  vector_sub(s->lambda2,s->lambda2_old,s->temp2_dim_M,M); // lambda2 - lambda2_old
235  vector_scalar_mul(s->temp2_dim_M,beta,s->temp2_dim_M,M); // ANS * beta
236  vector_add(s->lambda2,s->temp2_dim_M,s->y2,M); // y1 = ANS + lambda2
237  break;
238  default:
239  ERROR("Error deciding the problem case\n");
240  return -1;
241  }
242 
243 
244  /* *********************************************************************
245  Solving the inner problem
246  ********************************************************************* */
247 
248  // First calculate the new c (c_hat) for the inner problem,
249  // c_hat = p_in->c = c + A' * (y1_new - y2_new)
250 
251  vector_copy(s->prob->c,p_in->c,N);
252  switch (prob_case){
253  case 1:
254  // p_in->c = c + A' * (y1_new - y2_new)
255  vector_sub(s->y1,s->y2,s->temp2_dim_M,M); // y1 - y2
256  mtx_vec_mul(s->prob->A_t,s->temp2_dim_M,s->temp1_dim_N,N,M); // ANS * A'
257  vector_add(p_in->c,s->temp1_dim_N,p_in->c,N); // p_in->c = ANS + c
258  break;
259  case 2:
260  // p_in->c = c + A' * y1
261  mtx_vec_mul(s->prob->A_t,s->y1,s->temp1_dim_N,N,M); // A' * y1
262  vector_add(p_in->c,s->temp1_dim_N,p_in->c,N); // p_in->c = ANS + c
263  break;
264  case 3:
265  // p_in->c = c + A' * y1
266  mtx_vec_mul(s->prob->A_t,s->y1,s->temp1_dim_N,N,M); // A' * y1
267  vector_add(p_in->c,s->temp1_dim_N,p_in->c,N); // p_in->c = ANS + c
268  break;
269  case 4:
270  // p_in->c = c - A' * y2
271  mtx_vec_mul(s->prob->A_t,s->y2,s->temp1_dim_N,N,M); // A' * y2
272  vector_sub(p_in->c,s->temp1_dim_N,p_in->c,N); // p_in->c = c - ANS
273  break;
274  default:
275  ERROR("Error deciding the problem case\n");
276  return -1;
277  }
278 
279  // Compute the new optimal z (s->z) from the inner problem
280  vector_copy(s->z, p_in->z0, N); // Warm start
281 
282  tic_in = clock();
283  niter_inner += FGM(p_in);
284  toc_in = clock();
285  s->res->out->time_tot_inner += (real_t)(toc_in-tic_in);
286  s->time_inner_y += (real_t)(toc_in-tic_in);
287  if(p_in->exitflag == 2){
288  s->res->out->exitflag_inner = 2;
289  s->res->out->num_exceeded_max_niter_inner++;
290  }
291  vector_copy(p_in->zopt,s->z,N);
292 
293  /* *********************************************************************
294  Finding dual_value_diff
295  ********************************************************************* */
296 
297  // calculate the difference between the new and previousdual function value (ds)
298  // NOTE: must calculate the 'inner problem' one more time, with lambda instead of y
299 
300  vector_copy(s->prob->c,p_in->c,N);
301 
302  switch (prob_case){
303  case 1:
304  // p_in->c = c + A' * (lambda1 - lambda2)
305  vector_sub(s->lambda1,s->lambda2,s->temp2_dim_M,M); // lambda1 - lambda2
306  mtx_vec_mul(s->prob->A_t,s->temp2_dim_M,s->temp1_dim_N,N,M); // ANS * A'
307  vector_add(p_in->c,s->temp1_dim_N,p_in->c,N); // p_in->c = ANS + c
308  break;
309  case 2:
310  // p_in->c = c + A' * lambda1
311  mtx_vec_mul(s->prob->A_t,s->lambda1,s->temp1_dim_N,N,M); // A' * lambda1
312  vector_add(p_in->c,s->temp1_dim_N,p_in->c,N); // p_in->c = ANS + c
313  break;
314  case 3:
315  // p_in->c = c + A' * lambda1
316  mtx_vec_mul(s->prob->A_t,s->lambda1,s->temp1_dim_N,N,M); // A' * lambda1
317  vector_add(p_in->c,s->temp1_dim_N,p_in->c,N); // p_in->c = ANS + c
318  break;
319  case 4:
320  // p_in->c = c - A' * lambda2
321  mtx_vec_mul(s->prob->A_t,s->lambda2,s->temp1_dim_N,N,M); // A' * lambda2
322  vector_sub(p_in->c,s->temp1_dim_N,p_in->c,N); // p_in->c = c - ANS
323  break;
324  default:
325  ERROR("Error deciding the problem case\n");
326  return -1;
327  }
328 
329  // Finding new z_ds
330  vector_copy(s->z_ds, p_in->z0, N); // Warm start
331  tic_in = clock();
332  niter_inner += FGM(p_in);
333  toc_in = clock();
334  s->res->out->time_tot_inner += (real_t)(toc_in-tic_in);
335 
336  if(p_in->exitflag == 2){
337  s->res->out->exitflag_inner = 2;
338  s->res->out->num_exceeded_max_niter_inner++;
339  }
340  vector_copy(p_in->zopt,s->z_ds,N);
341  mtx_vec_mul(s->prob->A,s->z_ds,s->A_z_ds,M,N); // used in dual_obj
342 
343  // Calculate dual_value_new
344  dual_value_new = dual_obj(s,prob_case); // z_ds is used in this function
345  dual_value_diff = abs_2(dual_value_new - dual_value);
346  dual_value = dual_value_new;
347 
348  /* *********************************************************************
349  Finding pf
350  ********************************************************************* */
351 
352  // pf = norm2( max([A*z_pf - b_ub_hat ; -A*z_pf + b_lb_hat]) , 0) )
353  // NOTE: Algorithm 3 (last) => we use z_pf = z_ds
354  // Algorithm 4 (average) => we use z_pf = z_avg = summ_k / (niter+1)
355 
356  if (s->opt->algorithm == 3) {
357  // *** LAST z in stopping criteria ***
358 
359  switch (prob_case){
360  case 1:
361  vector_sub(s->A_z_ds,s->b_ub_hat,s->pf_vec,M); // (A*z) - b_ub_hat
362  vector_sub(s->b_lb_hat,s->A_z_ds,&s->pf_vec[M],M); // b_ub_hat - (A*z)
363  vector_max_with_zero(s->pf_vec, s->info->pf_vec_length); // Projection: max(ANS,0)
364  break;
365  case 2:
366  vector_sub(s->A_z_ds,s->b_ub_hat,s->pf_vec,M); // (A*z) - b_ub_hat
367  // NOTE: No projection
368  break;
369  case 3:
370  vector_sub(s->A_z_ds,s->b_ub_hat,s->pf_vec,M); // (A*z) - b_ub_hat
371  vector_max_with_zero(s->pf_vec, s->info->pf_vec_length); // Projection: max(ANS,0)
372  break;
373  case 4:
374  vector_sub(s->b_lb_hat,s->A_z_ds,s->pf_vec,M); // (A*z) + b_ub_hat
375  vector_max_with_zero(s->pf_vec, s->info->pf_vec_length); // Projection: max(ANS,0)
376  break;
377  default:
378  ERROR("Error deciding the problem case\n");
379  return -1;
380  }
381  }
382  else if (s->opt->algorithm == 4) {
383  // *** AVERAGE z in pf stopping criteria ***
384 
385  // Calculating z_avg
386  S = S + theta_new;
387  vector_scalar_mul(s->z,theta_new,s->temp1_dim_N,N); // S * theta_new
388  vector_add(s->summ,s->temp1_dim_N,s->summ,N); // summ += ANS
389  vector_scalar_mul(s->summ,1.0/S,s->z_avg,N); //z_avg = summ / S
390  mtx_vec_mul(s->prob->A,s->z_avg,s->temp2_dim_M,M,N); // A * z_pf = A * z_avg
391  switch (prob_case){
392  case 1:
393  vector_sub(s->temp2_dim_M,s->b_ub_hat,s->pf_vec,M); // (A*z) - b_ub_hat
394  vector_sub(s->b_lb_hat,s->temp2_dim_M,&s->pf_vec[M],M); // b_ub_hat - (A*z)
395  vector_max_with_zero(s->pf_vec, s->info->pf_vec_length); // Projection: max(ANS,0)
396  break;
397  case 2:
398  vector_sub(s->temp2_dim_M,s->b_ub_hat,s->pf_vec,M); // (A*z) - b_ub_hat
399  // NOTE: No projection
400  break;
401  case 3:
402  vector_sub(s->temp2_dim_M,s->b_ub_hat,s->pf_vec,M); // (A*z) - b_ub_hat
403  vector_max_with_zero(s->pf_vec, s->info->pf_vec_length); // Projection: max(ANS,0)
404  break;
405  case 4:
406  vector_sub(s->b_lb_hat,s->temp2_dim_M,s->pf_vec,M); // (A*z) + b_ub_hat
407  vector_max_with_zero(s->pf_vec, s->info->pf_vec_length); // Projection: max(ANS,0)
408  break;
409  default:
410  ERROR("Error deciding the problem case\n");
411  return -1;
412  }
413  }
414  else {
415  ERROR("Choose algorithm 3 or 4 when running DFGM()\n");
416  return -1;
417  }
418 
419  // Take the norm
420  pf = vector_norm_2(s->pf_vec, s->info->pf_vec_length); // norm2
421 
422  // store the result of the stopping criteria
423  insertArray(&ds_array, dual_value_diff);
424  insertArray(&pf_array, pf);
425 
426  /* *********************************************************************
427  Update the variables
428  ********************************************************************* */
429 
430  vector_copy(s->lambda1,s->lambda1_old,M);
431  vector_copy(s->lambda2,s->lambda2_old,M);
432  theta = theta_new;
433  niter++;
434 
435  } /* #### END WHILE-LOOP #### */
436 
437  /* *********************************************************************
438  Setting the result
439  ********************************************************************* */
440 
441  // Clock
442  toc = clock();
443  s->res->out->time = (real_t)(toc - tic) / CLOCKS_PER_SEC;
444  s->res->out->time_tot_inner /= CLOCKS_PER_SEC;
445  s->time_inner_y /= CLOCKS_PER_SEC;
446  //printf("time inner y: %f\n", s->time_inner_y);
447 
448  if (s->opt->algorithm == 3) {
449  s->res->zopt = s->z_ds;
450  s->res->fopt = obj(s->z_ds, s->prob->H, s->prob->c, s->temp1_dim_N);
451  }
452  else if (s->opt->algorithm == 4) {
453  s->res->zopt = s->z_avg;
454  s->res->fopt = obj(s->z_avg, s->prob->H, s->prob->c, s->temp1_dim_N);
455  }
456  else {
457  ERROR("Choose algorithm 3 or 4 when running DFGM()\n");
458  return -1;
459  }
460 
461  s->res->lambda1 = s->lambda1;
462  s->res->lambda2 = s->lambda2;
463  s->res->out->iterations = --niter;
464  s->res->out->iterations_inner_tot = niter_inner;
465  s->res->out->niter_feasible_ds = ++niter_feasible_ds;
466  s->res->out->niter_feasible_pf = ++niter_feasible_pf;
467 
468  if (last_eps_ds == 1)
469  s->res->out->flag_last_satisfied = 2;
470  else
471  s->res->out->flag_last_satisfied = 1;
472 
473  // Copy ds_vector and pf_vector
474  s->res->out->ds_vector = (real_t *)realloc(s->res->out->ds_vector, (niter) * sizeof(real_t));
475  s->res->out->pf_vector = (real_t *)realloc(s->res->out->pf_vector, (niter) * sizeof(real_t));
476  vector_copy(ds_array.array,s->res->out->ds_vector,niter);
477  vector_copy(pf_array.array,s->res->out->pf_vector,niter);
478  freeArray(&ds_array);
479  freeArray(&pf_array);
480 
481  return 0;
482 }
483 
484 
485 /* Definition of static functions */
486 
487 static real_t dual_obj(struct Struct_DFGM *s, const uint32_t prob_case)
488 {
489  /*return: 0.5z'Hz + c'z + lambda1'*(A*z - b_ub_hat) + lambda2'*(-A*z + b_lb_hat) */
490 
491  // f_value = 0.5*z'Hz + c'z;
492  real_t f_value = obj(s->z_ds, s->prob->H, s->prob->c, s->temp1_dim_N);
493 
494  switch (prob_case){
495  case 1:
496  // f_value += lambda1'*(A*z-b-ub_hat)
497  vector_sub(s->A_z_ds,s->b_ub_hat,s->temp3_dim_M,M); // ANS - b_ub_hat
498  f_value += vector_mul(s->lambda1,s->temp3_dim_M,M);
499 
500  // f_value += lambda2'*(-A*z+b+lb_hat)
501  vector_sub(s->b_lb_hat,s->A_z_ds,s->temp3_dim_M,M); // b_lb_hat - (Az)
502  f_value += vector_mul(s->lambda2,s->temp3_dim_M,M);
503  break;
504  case 2:
505  // f_value += lambda1'*(A*z-b-ub_hat)
506  vector_sub(s->A_z_ds,s->b_ub_hat,s->temp3_dim_M,M); // ANS - b_ub_hat
507  f_value += vector_mul(s->lambda1,s->temp3_dim_M,M);
508  break;
509  case 3:
510  // f_value += lambda1'*(A*z-b-ub_hat)
511  vector_sub(s->A_z_ds,s->b_ub_hat,s->temp3_dim_M,M); // ANS - b_ub_hat
512  f_value += vector_mul(s->lambda1,s->temp3_dim_M,M);
513  break;
514  case 4:
515  // f_value += lambda2'*(-A*z+b+lb_hat)
516  vector_sub(s->b_lb_hat,s->A_z_ds,s->temp3_dim_M,M); // b_lb_hat - (Az)
517  f_value += vector_mul(s->lambda2,s->temp3_dim_M,M);
518  break;
519  default:
520  ERROR("Error deciding the problem case\n");
521  return -1;
522  }
523 
524  return f_value;
525 }
526 
527 static void init_inner_problem(const struct Struct_DFGM *s, struct Struct_FGM *p_in)
528 {
529  // Setting the inner problem to point to the outer problem except for c and z0
530  // Setting the inner problem to point to the same as the outer saves some memory
531  p_in->H = s->prob->H;
532  p_in->lb = s->prob->lb;
533  p_in->ub = s->prob->ub;
534  p_in->lb_is_inf = s->info->lb_is_inf;
535  p_in->ub_is_inf = s->info->ub_is_inf;
536  p_in->eigH_max = s->info->eigH_max;
537  p_in->eigH_min = s->info->eigH_min;
538 
539  // Make own allocation for c (This is changing for each outer iteration)
540  p_in->c = vector_alloc(N);
541  p_in->z0 = vector_alloc(N);
542 
543  // Allocate other necessary vectors
544  p_in->zopt = vector_alloc(N);
545  p_in->z = vector_alloc(N);
546  p_in->y = vector_alloc(N);
547  p_in->znew = vector_alloc(N);
548  p_in->ynew = vector_alloc(N);
549  p_in->temp1_dim_N = vector_alloc(N);
550 
551  // Initialize options
552  p_in->maxiter = s->opt->maxiter_inner;
553  p_in->eps = s->opt->eps_inner;
554 }
555 
556 static void clean_up_inner_problem(struct Struct_FGM *s)
557 {
558  free_pointer(s->c);
559  free_pointer(s->z0);
560  free_pointer(s->zopt);
561  free_pointer(s->z);
562  free_pointer(s->y);
563  free_pointer(s->znew);
564  free_pointer(s->ynew);
565  free_pointer(s->temp1_dim_N);
566 }
Struct_DFGM::b_ub_hat
real_t * b_ub_hat
Definition: dfgm.h:34
Struct_DFGM::temp3_dim_M
real_t * temp3_dim_M
Definition: dfgm.h:33
vector_add
void vector_add(const real_t *v1, const real_t *v2, real_t *res, const uint32_t length)
Definition: math_functions.c:97
Result::lambda2
real_t * lambda2
Definition: qp_structs.h:69
free_pointer
void free_pointer(real_t *pointer)
Definition: general_functions.c:16
Output::time_tot_inner
real_t time_tot_inner
Definition: qp_structs.h:54
insertArray
void insertArray(struct Array *a, real_t element)
Definition: general_functions.c:29
Struct_DFGM::time_inner_y
real_t time_inner_y
Definition: dfgm.h:49
abs_2
real_t abs_2(const real_t a)
Definition: math_functions.c:177
uint32_t
unsigned int uint32_t
Definition: typedefs.h:19
Output::time
real_t time
Definition: qp_structs.h:53
Problem::A_t
real_t * A_t
Definition: qp_structs.h:19
Struct_DFGM::A_z
real_t * A_z
Definition: dfgm.h:39
Options::algorithm
uint32_t algorithm
Definition: qp_structs.h:34
Struct_DFGM::z
real_t * z
Definition: dfgm.h:28
_DEBUG
#define _DEBUG(fmt, args...)
Definition: head.h:41
Struct_FGM::eigH_min
real_t eigH_min
Definition: fgm.h:37
obj
real_t obj(const real_t *z, const real_t *H, const real_t *c, real_t *temp)
Definition: math_functions.c:184
Struct_FGM::ub_is_inf
boolean ub_is_inf
Definition: fgm.h:35
vector_min
void vector_min(const real_t *v1, const real_t *v2, real_t *res, const uint32_t length)
Definition: math_functions.c:67
initArray
void initArray(struct Array *a, uint32_t initialSize)
Definition: general_functions.c:23
vector_alloc
real_t * vector_alloc(uint32_t size)
Definition: general_functions.c:10
Info::Ld
real_t Ld
Definition: qp_structs.h:45
Result::out
struct Output * out
Definition: qp_structs.h:70
Info::ub_is_inf
boolean ub_is_inf
Definition: qp_structs.h:40
Result::zopt
real_t * zopt
Definition: qp_structs.h:65
Output::iterations_inner_tot
uint32_t iterations_inner_tot
Definition: qp_structs.h:52
DFGM
int32_t DFGM(struct Struct_DFGM *s)
Definition: dfgm.c:20
Struct_DFGM::b_lb_hat
real_t * b_lb_hat
Definition: dfgm.h:35
Struct_FGM::ub
real_t * ub
Definition: fgm.h:25
Struct_DFGM::res
struct Result * res
Definition: dfgm.h:25
Problem::lb
real_t * lb
Definition: qp_structs.h:23
Struct_FGM::temp1_dim_N
real_t * temp1_dim_N
Definition: fgm.h:44
Struct_FGM::znew
real_t * znew
Definition: fgm.h:42
Output::ds_vector
real_t * ds_vector
Definition: qp_structs.h:60
Result::exitflag
uint32_t exitflag
Definition: qp_structs.h:67
Output::flag_last_satisfied
uint32_t flag_last_satisfied
Definition: qp_structs.h:55
Struct_FGM::eps
real_t eps
Definition: fgm.h:48
Struct_FGM::lb
real_t * lb
Definition: fgm.h:24
Struct_FGM::y
real_t * y
Definition: fgm.h:41
Struct_FGM::eigH_max
real_t eigH_max
Definition: fgm.h:36
Info::pf_vec_length
uint32_t pf_vec_length
Definition: qp_structs.h:47
Struct_DFGM::lambda2
real_t * lambda2
Definition: dfgm.h:30
Options::eps_ds
real_t eps_ds
Definition: qp_structs.h:31
Problem::H
real_t * H
Definition: qp_structs.h:16
Struct_DFGM::lambda1_old
real_t * lambda1_old
Definition: dfgm.h:42
Result::fopt
real_t fopt
Definition: qp_structs.h:66
Output::iterations
uint32_t iterations
Definition: qp_structs.h:51
Struct_FGM::H
real_t * H
Definition: fgm.h:22
Struct_FGM::lb_is_inf
boolean lb_is_inf
Definition: fgm.h:34
Struct_FGM::exitflag
uint32_t exitflag
Definition: fgm.h:31
Problem::A
real_t * A
Definition: qp_structs.h:18
Struct_DFGM::summ
real_t * summ
Definition: dfgm.h:37
Problem::ub
real_t * ub
Definition: qp_structs.h:24
Struct_DFGM::opt
struct Options * opt
Definition: dfgm.h:23
Options::maxiter_outer
uint32_t maxiter_outer
Definition: qp_structs.h:29
vector_sub
void vector_sub(const real_t *v1, const real_t *v2, real_t *res, const uint32_t length)
Definition: math_functions.c:89
Struct_DFGM
Definition: dfgm.h:20
Struct_DFGM::prob
struct Problem * prob
Definition: dfgm.h:22
Struct_DFGM::y2
real_t * y2
Definition: dfgm.h:46
FGM
uint32_t FGM(struct Struct_FGM *s)
Definition: fgm.c:12
int32_t
signed int int32_t
Definition: typedefs.h:15
vector_max_with_zero
void vector_max_with_zero(real_t *v, const uint32_t length)
Definition: math_functions.c:149
Struct_DFGM::lambda1
real_t * lambda1
Definition: dfgm.h:29
Struct_DFGM::temp1_dim_N
real_t * temp1_dim_N
Definition: dfgm.h:31
Output::niter_feasible_pf
uint32_t niter_feasible_pf
Definition: qp_structs.h:57
Result::lambda1
real_t * lambda1
Definition: qp_structs.h:68
vector_norm_2
real_t vector_norm_2(real_t *v, const uint32_t length)
Definition: math_functions.c:158
M
uint32_t M
Definition: head.h:34
Output::num_exceeded_max_niter_inner
uint32_t num_exceeded_max_niter_inner
Definition: qp_structs.h:59
Struct_DFGM::info
struct Info * info
Definition: dfgm.h:24
Struct_DFGM::temp2_dim_M
real_t * temp2_dim_M
Definition: dfgm.h:32
vector_copy
void vector_copy(const real_t *v1, real_t *v2, const uint32_t length)
Definition: math_functions.c:141
Options::eps_inner
real_t eps_inner
Definition: qp_structs.h:33
Struct_FGM::z
real_t * z
Definition: fgm.h:40
real_t
float64_t real_t
Definition: typedefs.h:25
Output::exitflag_inner
uint32_t exitflag_inner
Definition: qp_structs.h:58
Info::lb_is_inf
boolean lb_is_inf
Definition: qp_structs.h:39
Struct_DFGM::y1
real_t * y1
Definition: dfgm.h:44
vector_scalar_mul
void vector_scalar_mul(const real_t *v1, const real_t scalar, real_t *res, const uint32_t length)
Definition: math_functions.c:115
Struct_DFGM::z_avg
real_t * z_avg
Definition: dfgm.h:36
mtx_vec_mul
void mtx_vec_mul(const real_t *mtx, const real_t *v, real_t *res, const uint32_t rows, const uint32_t cols)
Definition: math_functions.c:23
Struct_FGM::ynew
real_t * ynew
Definition: fgm.h:43
_SDEBUG
#define _SDEBUG(fmt)
Definition: head.h:60
Struct_FGM::z0
real_t * z0
Definition: fgm.h:26
Options::maxiter_inner
uint32_t maxiter_inner
Definition: qp_structs.h:30
Struct_DFGM::A_z_ds
real_t * A_z_ds
Definition: dfgm.h:47
Struct_DFGM::z_ds
real_t * z_ds
Definition: dfgm.h:45
Problem::z0
real_t * z0
Definition: qp_structs.h:25
Struct_FGM::c
real_t * c
Definition: fgm.h:23
vector_max
void vector_max(const real_t *v1, const real_t *v2, real_t *res, const uint32_t length)
Definition: math_functions.c:78
Problem::ub_hat
real_t * ub_hat
Definition: qp_structs.h:22
Info::eigH_min
real_t eigH_min
Definition: qp_structs.h:44
Info::problem_case
uint32_t problem_case
Definition: qp_structs.h:46
ERROR
#define ERROR(fmt)
Definition: head.h:54
Struct_DFGM::pf_vec
real_t * pf_vec
Definition: dfgm.h:38
Output::pf_vector
real_t * pf_vector
Definition: qp_structs.h:61
_DEBUG2
#define _DEBUG2(fmt, args...)
Definition: head.h:46
Options::eps_pf
real_t eps_pf
Definition: qp_structs.h:32
Info::eigH_max
real_t eigH_max
Definition: qp_structs.h:43
Struct_FGM::maxiter
uint32_t maxiter
Definition: fgm.h:47
Output::niter_feasible_ds
uint32_t niter_feasible_ds
Definition: qp_structs.h:56
Struct_FGM
Definition: fgm.h:20
Problem::lb_hat
real_t * lb_hat
Definition: qp_structs.h:21
Struct_DFGM::lambda2_old
real_t * lambda2_old
Definition: dfgm.h:43
freeArray
void freeArray(struct Array *a)
Definition: general_functions.c:37
vector_mul
real_t vector_mul(const real_t *v1, const real_t *v2, const uint32_t length)
Definition: math_functions.c:105
Struct_FGM::zopt
real_t * zopt
Definition: fgm.h:29
Problem::b
real_t * b
Definition: qp_structs.h:20
N
uint32_t N
Definition: head.h:33
Problem::c
real_t * c
Definition: qp_structs.h:17