When I use TFOCS to solve a constrained optimization problem, sometimes the solution does not satisfy the constraints. Is this a bug with TFOCS?

This is not a bug in TFOCS, but rather it is a consequence of the fact that TFOCS can dualize constraints. By putting constraints into the dual function, the problem because tractable but the consequence is that the primal variable does not need to satisfy the constraint exactly. All we know is that if you run the program forever and get the exact dual variable, then the corresponding exact primal variable will satisfy the constraint.

In practice, if the primal constraint is not satisfied, it can usually be fixed by simply asking TFOCS to take more iterations. Below is an example that shows this behavior and illustrate some tools you can use.

Consider the following basis pursuit problem. If youâre not familiar with this problem, the only important thing to know is that we are asking for a variable *x* that satisfies \|Ax-b\| \le \epsilon.

```
reset( RandStream.getDefaultStream );
A = randn(20,40);
b = randn(20,1);
epsilon = 1e-2;
mu = 1e-3;
x0 = [];
z0 = [];
opts = struct( 'tol', 1e-5 );
x = solver_sBPDN( A, b, epsilon, mu, x0, z0, opts);
```

and then we ask if the constraints are satisfied:

`norm( A*x - b )`

. The answer is `0.3519`

, which is clearly not less than epsilon (which was .01).

The problem is that TFOCS did not run for enough iterations. Letâs try adding the following code:

```
opts = struct('printStopCrit',true,'tol',1e-7);
opts.errFcn = @(f,p,d) norm( A*d-b );
x = solver_sBPDN( A, b, epsilon, mu, x0, z0,opts);
```

The most obvious change is that we are solving to a higher precision (1e-7 instead of 1e-5). We have also requested `printStopCrit`

. TFOCS can use different stopping criteria, and each criteria is a test of the form `currentValue < tol`

. Asking for `printStopCrit`

just tells TFOCS to print `currentValue`

to the output so the user gets an idea of whether they want to increase or decrease the tolerance.

The final difference is that we are supplying an error function. This takes three arguments: the objective function value at the iterate, the main variable at the iterate, and the dual variable at the iterate. The only tricky thing to remember is that the âmainâ variable refers to the dual variable (since in the âSCDâ mode, we solve the dual problem), so the dual of the dual variable is the primal. In the above code, we label these variables by `f,p,d`

respectively. If we want to check the value of \|Ax-b\| at every iteration, then we just compute `norm( A*d - b )`

. The output of this will be displayed in the âerrorâ column.

Letâs see the output of this new code.

```
Auslender & Teboulle's single-projection method
Iter Objective |dx|/|x| step errors stopping criteria
----+----------------------------------+----------+-------------------
100 | +9.12944e-01 2.40e-03 2.01e-05 | 2.26e+00 | 2.40e-03
200 | +1.52453e+00 4.73e-03 6.10e-05 | 2.57e+00 | 4.73e-03
300 | +1.93661e+00 2.23e-03 1.69e-05 | 9.37e+00 | 2.23e-03
400 | +2.04548e+00 1.31e-03 1.87e-05 | 1.76e+00 | 1.31e-03
500 | +2.05480e+00 3.19e-04 2.48e-05 | 4.65e-01 | 3.19e-04
600 | +2.06953e+00 4.09e-04 1.78e-05 | 4.23e-01 | 4.09e-04
700 | +2.09176e+00 5.56e-04 1.52e-05 | 3.72e-01 | 5.56e-04
800 | +2.12275e+00 1.08e-03 3.33e-05 | 3.71e-01 | 1.08e-03
900 | +2.12776e+00 3.39e-04 4.72e-05 | 4.21e-01 | 3.39e-04
1000| +2.14052e+00 4.23e-04 2.35e-05 | 6.14e-01 | 4.23e-04
...
```

We see that the error function is decreasing, but still not below epsilon. Given the value of the stopping criteria at iteration 1000, we can see that a tolerance of just 1e-4 or 1e-5 is probably too loose.

The end of the output is:

```
...
5300| +2.36216e+00 6.21e-05 2.81e-05 | 1.55e-01 | 6.21e-05
5400| +2.36402e+00 1.34e-04 3.47e-05 | 1.55e-01 | 1.34e-04
5500| +2.36707e+00 1.65e-04 2.44e-05 | 1.55e-01 | 1.65e-04
5600| +2.36821e+00 6.81e-05 2.84e-05 | 2.23e-01 | 6.81e-05
5700| +2.36857e+00 4.16e-05 2.92e-05 | 9.92e-02 | 4.16e-05
5800| +2.36895e+00 3.07e-05 2.38e-05 | 1.60e-01 | 3.07e-05
5900| +2.36896e+00 3.99e-06 1.43e-05 | 6.62e-02 | 3.99e-06
6000| +2.36896e+00 2.50e-06 2.40e-05 | 1.23e-02 | 2.50e-06
6037| +2.36896e+00 5.56e-08 5.24e-06*| 1.03e-02 | 5.56e-08
Finished: Step size tolerance reached
```

After 5000 iterations, the constraint is still not satisfied. Around 5900 iterations, the variable enters a different type of region, and suddenly the convergence is very rapid. By 6037 iterations, the inequality constraint is basically satisfied (the norm of `A*x-b`

is 1.03e-2, and the desired norm was 1.0e-2).

Of course, another method to ensure the constraints are satisfied is to increase convergence speed by playing around with the smoothing parameter `mu`

, the initial values of `x0`

, and continuation settings.