Quadratic Objective with two Quadratic Constraints

Question

Krishnendu K 2023-11-8

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2044402-quadratic-objective-with-two-quadratic-constraints

评论： Torsten 2023-11-9

Hii,

My programming knowledge is limited. I would like to solve a quadratic objective function with two quadratic constraints. I read and understood the concepts and examples given in https://in.mathworks.com/help/optim/ug/linear-or-quadratic-problem-with-quadratic-constraints.html?s_tid=mwa_osa_a, but I am having difficulty in implementing this for my problem.

Matrices A, B and C be

symmetric positive semi-definite matrices. How to find the

subject to the constraints

and

? The documentation deals with a single quadratic constraint only. But my problem has two quadratic constraints. How can I code the second constraint as per the documentation in the link above?

Thank You..

3 个评论
显示 1更早的评论隐藏 1更早的评论

Krishnendu K 2023-11-9

Is there any condition can I check whether the constrained region have an intersection?

I have generated the matrices A, B and C and attached the .mat files. Is there any possible way of optimisation?

Matt J 2023-11-9

编辑：Matt J 2023-11-9

Even if they have an intersection, I don't see how that solves the problem. Do you have an initial guess close enough to the global solution that disconnected solution regions will be avoided?

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Torsten 2023-11-8

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2044402-quadratic-objective-with-two-quadratic-constraints#answer_1348787

编辑：Torsten 2023-11-8

在 MATLAB Online 中打开

What about

function [y,yeq,grady,gradyeq] = quadconstr(x,B,C)
  y = [];
  yeq(1) = x.'*B*x - 1;
  yeq(2) = x.'*C*x - 1;
  if nargout > 2
    grady = [];
    gradyeq(:,1) = 2*B*x;    % Assumes B is symmetric, otherwise (B+B.')*x
    gradyeq(:,2) = 2*C*x;    % Assumes C is symmetric, otherwise (C+C.')*x
  end
end

10 个评论
显示 8更早的评论隐藏 8更早的评论

Krishnendu K 2023-11-9

编辑：Torsten 2023-11-9

在 MATLAB Online 中打开

Still there is problem with the initial values of x. Can you pls check whether the code works fine for the attached A, B and C matrices of n=6 size?

Amat = load("Matrix_A.mat");
A = Amat.P1
A = 
   0.1667 + 0.0000i   0.0625 + 0.0361i  -0.0208 - 0.0361i  -0.0000 - 0.0000i  -0.0208 + 0.0361i   0.0625 - 0.0361i
   0.0625 - 0.0361i   0.9167 + 0.0000i   0.2500 + 0.1443i  -0.2708 - 0.0361i   0.0000 - 0.0000i  -0.0833 + 0.1443i
  -0.0208 + 0.0361i   0.2500 - 0.1443i   0.9167 + 0.0000i   0.0625 + 0.0361i  -0.3333 - 0.1443i   0.0000 - 0.0000i
  -0.0000 + 0.0000i  -0.2708 + 0.0361i   0.0625 - 0.0361i   0.6667 + 0.0000i   0.0625 + 0.0361i  -0.0208 - 0.0361i
  -0.0208 - 0.0361i   0.0000 + 0.0000i  -0.3333 + 0.1443i   0.0625 - 0.0361i   0.9167 + 0.0000i   0.2500 + 0.1443i
   0.0625 + 0.0361i  -0.0833 - 0.1443i   0.0000 + 0.0000i  -0.0208 + 0.0361i   0.2500 - 0.1443i   0.4167 + 0.0000i
Bmat = load("Matrix_B.mat");
B = Bmat.Lambda_n
B = 6×6
    1.0000         0   -0.5000         0         0         0
         0    1.0000         0   -0.5000         0         0
   -0.5000         0    1.0000         0   -0.5000         0
         0   -0.5000         0    1.0000         0   -0.5000
         0         0   -0.5000         0    1.0000         0
         0         0         0   -0.5000         0    1.0000
Cmat = load("Matrix_C.mat");
C = Cmat.Gamma_n
C = 
   0.3333 + 0.0000i  -0.1250 - 0.0722i  -0.0417 - 0.0722i   0.0000 + 0.0000i  -0.0417 + 0.0722i  -0.1250 + 0.0722i
  -0.1250 + 0.0722i   1.3333 + 0.0000i   0.2500 + 0.1443i  -0.0417 - 0.0722i   0.0000 + 0.0000i   0.0833 - 0.1443i
  -0.0417 + 0.0722i   0.2500 - 0.1443i   1.3333 + 0.0000i  -0.1250 - 0.0722i   0.0833 + 0.1443i   0.0000 + 0.0000i
   0.0000 - 0.0000i  -0.0417 + 0.0722i  -0.1250 + 0.0722i   0.3333 + 0.0000i  -0.1250 - 0.0722i  -0.0417 - 0.0722i
  -0.0417 - 0.0722i   0.0000 - 0.0000i   0.0833 - 0.1443i  -0.1250 + 0.0722i   1.3333 + 0.0000i   0.2500 + 0.1443i
  -0.1250 - 0.0722i   0.0833 + 0.1443i   0.0000 - 0.0000i  -0.0417 + 0.0722i   0.2500 - 0.1443i   1.3333 + 0.0000i
N =size(A,1);
fun = @(x)quadobj(x,A);
nonlconstr = @(x)quadconstr(x,B,C);
x0 = zeros(N,1); % Column vector
options = optimoptions(@fmincon,'Algorithm','interior-point',...
    'SpecifyObjectiveGradient',true,'SpecifyConstraintGradient',true);
[x,fval,eflag,output,lambda] = fmincon(fun,x0,...
    [],[],[],[],[],[],nonlconstr,options);
Converged to an infeasible point.

fmincon stopped because the size of the current step is less than
the value of the step size tolerance but constraints are not
satisfied to within the value of the constraint tolerance.


Consider enabling the interior point method feasibility mode.
function [y,grady] = quadobj(x,A)
  y = -x.'*A*x ;
  if nargout > 1
    grady = -2*A*x ;
  end
end
function [y,yeq,grady,gradyeq] = quadconstr(x,B,C)
  y = [];
  yeq(1) = x.'*B*x - 1;
  yeq(2) = x.'*C*x - 1;
  if nargout > 2
    grady = [];
    gradyeq(:,1) = 2*B*x;    % Assumes B is symmetric, otherwise (B+B.')*x
    gradyeq(:,2) = 2*C*x;    % Assumes C is symmetric, otherwise (C+C.')*x
  end
end

Krishnendu K 2023-11-9

What I meant is I tried to solve it by using eigen value decomposition. But the area of optimisation is entirely new to me. When I apply these constraints the problem cannot be solved by eigen value decomposition method and I have to rely entirely on convex optimisation which is new to me..

Torsten 2023-11-9

What I meant is I tried to solve it by using eigen value decomposition.

You formulated an optimization problem. Obviously, the problem is not well-posed because the objective function is not real-valued.

What is the underlying problem that you tried to solved via this optimization formulation ? (You again talk of "I tried to solve it", but you don't explain what "it" is).

请先登录，再进行评论。