Home > source > @penlab > lmltm_update.m

# lmltm_update

## PURPOSE

Lagrangian Multipliers (for Matrix constraints) update

## SYNOPSIS

function [] = lmltm_update(obj)

## DESCRIPTION

     Lagrangian Multipliers (for Matrix constraints) update

## CROSS-REFERENCE INFORMATION

This function calls:

This function is called by:

## SUBFUNCTIONS

- function [M] = m_stabilize(M)

## SOURCE CODE

    0001 % Lagrangian Multipliers (for Matrix constraints) update
    0002 function [] = lmltm_update(obj)
    0003
    0004   % obviously only pen/bar
    0005   % I have in the object: UYbox, UYboxnew, UA, UAnew - cell array with indices
    0006   % matching Yboxindphi & Aindphi
    0007
    0008   mu2 = obj.allopts.mmlt_update;
    0009
    0010   % compute new candidates in *new; this will be done in eval*() as a by-product
    0011
    0012   % matrix variable - pen/bar
    0013   for k=obj.Yboxindphi
    0014     pkx=obj.PYbox(k);            %2*p(sdpdata.Ng+k);
    0015     %Akx=evalax(x,k,sdpdata);
    0016     % mapping xall first !!! TODO take obj.Y{...}
    0017     Ykx = obj.Y{obj.Yboxmap(k)};
    0018     Akx=obj.Yboxshift(k)*speye(size(Ykx)) + obj.Yboxmlt(k)*Ykx;
    0019     umatk=obj.UYbox{k};
    0020
    0021     %Z=(pkx*speye(size(Akx))-Akx);
    0022     Z=(pkx*speye(size(Akx))-Akx);
    0023     invZ=inv(Z);
    0024     pZUZ=pkx^2*invZ*umatk*invZ;
    0025     obj.UYboxnew{k}=pZUZ;
    0026   end
    0027
    0028   % matrix constraints - pen/bar
    0029   for k=obj.Aindphi
    0030     pkx=obj.PA(k);  % I used to use 2*         !!!!!!!!
    0031     % TODO need to map the matrix first! - is it correct???
    0032     kuser=obj.Amap(k);
    0033     [Akuserx, obj.userdata] = obj.mconfun(obj.x, obj.Y, kuser, obj.userdata);
    0034     Akx = obj.Ashift(k)*speye(size(Akuserx)) + obj.Amlt(k) .* Akuserx;
    0035     umatk=obj.UA{k};
    0036
    0037     %Z=(pkx*speye(size(Akx))-Akx);
    0038     Z=(pkx*speye(size(Akx))-Akx);
    0039     invZ=inv(Z);
    0040     pZUZ=pkx^2*invZ*umatk*invZ;
    0041     obj.UAnew{k}=pZUZ;
    0042   end
    0043
    0044   % update as in update_multipliers() @ pbmfnc/Pennon1.0
    0045   % [in SDP/BMI there is still active one more option]
    0046   %  rStabilize = 1.0E-16*sqrt((double)(Umat->mat[i]->m)*diagnorm(Umat->mat[i], nSparse));
    0047   %  m_stabilize(Umat->mat[i],  rStabilize, nSparse);
    0048   %    rLambda = 1-mu2;
    0049   %    m_ConvComb(rLambda, Umattmp->mat[i], Umat->mat[i], Umat->mat[i], nSparse);
    0050   %  rStabilize = 1.0E-16*sqrt((double)(Umat->mat[i]->m)*diagnorm(Umat->mat[i], nSparse));
    0051   %  m_stabilize(Umat->mat[i],  rStabilize, nSparse);
    0052   % Note that in the code Umat & Umattmp are swapped (Umat is the new one in
    0053   % this part)
    0054
    0055   % matrix variable - pen/bar
    0056   for k=obj.Yboxindphi
    0057     umatk = obj.UYboxnew{k};
    0058     umatk = m_stabilize(umatk);
    0059     %lambda=1-mu2;
    0060     %umatk = (1-lambda)*old + lambda*new;
    0061     umatk = mu2*obj.UYbox{k} + (1-mu2)*umatk;
    0062     obj.UYbox{k} = m_stabilize(umatk);
    0063     % & push it back but keep the 'old ones'??
    0064     % are they used behind update_multipliers??
    0065   end
    0066
    0067   % matrix constraints - pen/bar
    0068   for k=obj.Aindphi
    0069     umatk = obj.UAnew{k};
    0070     umatk = m_stabilize(umatk);
    0071     %lambda=1-mu2;
    0072     %umatk = (1-lambda)*old + lambda*new;
    0073     umatk = mu2*obj.UA{k} + (1-mu2)*umatk;
    0074     obj.UA{k} = m_stabilize(umatk);
    0075     % & push it back but keep the 'old ones'??
    0076     % are they used behind update_multipliers??
    0077   end
    0078
    0079 end
    0080
    0081 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    0082 % stabilize M
    0083 % only on diag elements: m_ii = |m_ii| + r
    0084 % where r is given by the norm of the diag
    0085 % how to access the diagonal of a dense matrix?
    0086 %   M(1:n+1:end) = newdiag
    0087 %   iii=1:n+1:n^2; M(iii)=...   [slightly slower]
    0088 %   or normal for loop
    0089 function [M] = m_stabilize(M)
    0090   [n m] = size(M);
    0091   Mdiag = diag(M);
    0092   rStabilize = 1.0e-16*realsqrt(n*norm(Mdiag,2));
    0093   M(1:n+1:end) = abs(Mdiag) + rStabilize;
    0094 end
    0095
    0096