 penm.lbx = [-Inf, -Inf, -1, 0, -3, -6, -77]; 
%   penm.ubx = [10, Inf, Inf, 4, 66, 0, Inf];
%
% * [optional] *box constraints treatment*, by default a                
% penalty/barrier function is used to include box constraints into      
% the Augmented Lagrangian which might cause that the variables         
% are feasible only up to the given tolerance (see option settings).         
% Alternatively, it is possible to apply barrier function to all or     
% some of them instead and achieve stri