# Maximize a six-dimensional function subject to joint positive-semidefiniteness constraints

I want to maximize

Abs[a1 b1] + Abs[a2 b2] + Abs[a3 b3]


subject to the joint constraints

9 (a1^2+a2^2)<=4&&18 (a1^2+a2^2)+9 (2+3 a2) a3^2<=8&&4 (b1^2+b2^2+b3^2)<=1


(The first constraint is circular and the last constraint, spherical in nature. The middle one is independent of the $$b$$'s.)

In lieu of an exact solution, a high-precision numerical one would be desired.

Conjecturally, the exact solution (to this quantum-information-related problem) has a denominator that is the product of powers of 2 and/or of 3.

To further expand, the constraints were obtained by requiring the joint positive-semidefiniteness of the $$3 \times 3$$ and $$4 \times 4$$ ("density") matrices

{{1/3 - a2/2, -((I a1)/2), (I a3)/2}, {(I a1)/2, 1/3 + a2/2, 0}, {-((I a3)/2), 0, 1/3}}


and

{{1/4, 0, b1/2, 0}, {0, 1/4, 1/2 (I b2 - b3), 0}, {b1/2, O1/2 (-I b2 - b3), 1/4, 0}, {0, 0, 0, 1/4}}

• Have you stated the problem correctly? If a2==-2/3, I think a3 is unbounded. NMaximize should handle problems of this sort. Commented Feb 1, 2020 at 13:33
• I did try NMaximize rather extensively--couldn't obtain any convincing convergence--despite various WorkingPrecision and MaxIterations choices. Not sure about your a2==-2/3 observation at this point. Oops--the first constraint is "circular" not "spherical". Will correct. Commented Feb 1, 2020 at 14:11
• As two others have observed, it is unbounded. Is one of the constraints slightly off? Commented Feb 1, 2020 at 16:35
• Do you have a typo in your matrices, specifically, the variable O1? Commented Feb 1, 2020 at 19:06
• I don't think your conditions ensure that your matrices are positive semi-definite. E.g. with b1==b2==b3==1/Sqrt[3] Commented Feb 1, 2020 at 19:17

\$Version

(* "12.0.0 for Mac OS X x86 (64-bit) (April 7, 2019)" *)

Clear["Global*"]


Assuming that all of the variables are real

sys1 = {Abs[a1 b1] + Abs[a2 b2] + Abs[a3 b3],
9 (a1^2 + a2^2) <= 4,
18 (a1^2 + a2^2) + 9 (2 + 3 a2) a3^2 <= 8,
4 (b1^2 + b2^2 + b3^2) <= 1} /.
Abs[z_] :> Sqrt[z^2];

var1 = Variables[Level[sys1, {-1}]]

(* {a1, a2, a3, b1, b2, b3} *)

NMaximize[sys1, var1, WorkingPrecision -> 20]

(* {0.33333333333333343737, {a1 -> -0.63260132112986447330,
a2 -> 0.21038063824695148785, a3 -> 4.7334038188327067403*10^-10,
b1 -> -0.47445099084739835497, b2 -> 0.15778547868521361589,
b3 -> -3.5500528641245300550*10^-10}} *)


Let a3 == 0 and b3 == 0

sol1 = {a3 -> 0, b3 -> 0};

sys2 = sys1 /. sol1

(* {Sqrt[a1^2 b1^2] + Sqrt[a2^2 b2^2], 9 (a1^2 + a2^2) <= 4,
18 (a1^2 + a2^2) <= 8, 4 (b1^2 + b2^2) <= 1} *)

var2 = Variables[Level[sys2, {-1}]]

(* {a1, a2, b1, b2} *)

sol2 = Maximize[sys2, var2]

(* {1/3, {a1 -> -(5/16), a2 -> -(Sqrt[799]/48), b1 -> -(15/64),
b2 -> -(Sqrt[799]/64)}} *)

sol = Join[sol1, sol2[[2]]] // Sort

(* {a1 -> -(5/16), a2 -> -(Sqrt[799]/48), a3 -> 0, b1 -> -(15/64),
b2 -> -(Sqrt[799]/64), b3 -> 0} *)

Abs[a1 b1] + Abs[a2 b2] + Abs[a3 b3] /. sol

(* 1/3 *)
`

EDIT: As pointed out in a comment by user64494, this is a local maximum.

• Wow--terrific looking! Seems like the denominator-power conjecture of mine is substantiated--that is, given the 1/3 answer. Commented Feb 1, 2020 at 15:07
• How about 9 (a1^2 + a2^2) <= 4 && 18 (a1^2 + a2^2) + 9 (2 + 3 a2) a3^2 <= 8 && 4 (b1^2 + b2^2 + b3^2) <= 1 /. {a1 -> 0, a2 -> -2/3, a3 -> 100, b1 -> 0, b2 -> 0, b3 -> 1/2} which produces True and Abs[a1 b1] + Abs[a2 b2] + Abs[a3 b3] /. {a1 -> 0, a2 -> -2/3, a3 -> 100, b1 -> 0, b2 -> 0, b3 -> 1/2} which produces 50? Commented Feb 1, 2020 at 15:18