# Solution of 20-variable Constrained NMaximize with Method->"DifferentialEvolution" does not fully satisfy the constraints

I want to find the maximum of

d2=    (Sqrt[Abs[a^2 b^2]] + Sqrt[Abs[a^2 b^2]] + Sqrt[Abs[a^2 b^2]] + Sqrt[Abs[a^2 b^2]] + Sqrt[Abs[a^2 b^2]] + Sqrt[Abs[a^2 b^2]] + (1/Sqrt)(\[Sqrt]Abs[a^2 b^2 + a^2 b^2 + a^2 b^2 + a^2 b^2 - \[Sqrt](((a b + a b)^2 + (a b - a b)^2) ((a b - a b)^2 + (a b +a b)^2))]) + (1/Sqrt)(\[Sqrt]Abs[a^2 b^2 + a^2 b^2 + a^2 b^2 + a^2 b^2 + \[Sqrt](((a b + a b)^2 + (a b - a b)^2) ((a b - a b)^2 + (a b +a b)^2))]))^2


subject to a very large set (LeafCount of 27947) of constraints designated "vcf" in the notebook

https://www.wolframcloud.com/obj/slater/Published/MaximizeNorms

The command

G1 = NMaximize[{d2, vcf}, Join[Array[a, 10], Array[b, 10]]]


yields the result

{0.1565, {a -> -0.217171, a -> 0.228215, a -> 0.242746, a -> 0.185082, a -> -0.194495, a -> -0.239223, a -> -0.25139, a -> -0.236464, a -> 0.157249,  a -> 0.152502, b -> -0.187483, b -> 0.197017, b -> 0.217313, b -> 0.201103, b -> 0.19137,b -> 0.269263, b -> -0.256231, b -> 0.197663, b -> 0.0729522, b -> 0.17673}}


which does in fact satisfy the constraint vcf.

However, the command

G2 = NMaximize[{d2, vcf}, Join[Array[a, 10], Array[b, 10]], Method -> "DifferentialEvolution"]


yields the larger (and more plausible/attractive) result

{0.186769, {a -> -0.268552, a -> -0.255555, a -> 0.227544,a -> -0.18653, a -> 0.196016, a -> -0.203221, a -> 0.193385, a -> -0.316546, a -> -0.162042, a -> -0.153076, b -> -0.237133, b -> -0.249411,  b -> 0.215763, b -> 0.20219, b -> -0.192236, b -> -0.193805, b -> 0.203839, b -> -0.31299, b -> -0.026066, b -> -0.201624}}


but does not fully satisfy the constraint vcf.

(I believe that there is an exact value--possibly rational in nature [not too much larger than the 0.186769]--for this [quantum-information-theoretic] problem, but it certainly seems to be too demanding a task to obtain it. So, I'd like to "believe" the 0.186769 result, but the evidence for it is clearly less than 100% convincing.)

There is a second expression

d1=a^2 a^2 a^2 a^2 a^2 a^2 b^2 b^2 b^2 b^2 b^2 b^2 (a a b b - a a b b)^2

which it is also of interest to maximize with respect to the same constraints.

These constraints are so complicated that it is very easy to get caught in local maxima. Maybe RandomSearch can find a way out.

I got a much higher maximum than 0.186, namely 3.7354. (First I did some simplifications). Verification yields "True".

d2 = d2 //
Simplify[#, Thread[Join[Array[a, 10], Array[b, 10]] \[Element] Reals]] &

vcfl = vcf[] // Rationalize //
Simplify[#, Thread[Join[Array[a, 10], Array[b, 10]] \[Element] Reals]] &;

(H1 = NMaximize[{d2, vcfl}, Join[Array[a, 10], Array[b, 10]],
Method -> "RandomSearch"]) // Timing

(*   {209.609, {3.7354, {a[1.] -> -0.240179, a[2.] -> -0.228551,
a[3.] -> 0.196117, a[4.] -> 0.222086, a[5.] -> 0.211336,
a[6.] -> 0.211402, a[7.] -> 0.222153, a[8.] -> 0.145827,
a[9.] -> -0.195623, a[10.] -> -0.281385, b[1.] -> -0.836558,
b[2.] -> -1.1664, b[3.] -> 0.971168, b[4.] -> -1.37882,
b[5.] -> 0.96628, b[6.] -> -1.12641, b[7.] -> -1.2148,
b[8.] -> 0.811886, b[9.] -> -0.481526, b[10.] -> 0.812585}}}   *)

vcfl /. H1[]

(*   True   *)


May be you can play with "MaxIterations" to get an even higher maximum.

• Due to the awkwardness of using the Wolfram cloud to post the notebook, I broke giving the constraint "vcf" into three steps. First, I gave an expression "vcf", then "vcf2", then I did the command vcf =Join[vcf,vcf2]. I note that the norm of the solution 20-vector is 3.25719, and some of the components are greater than 1 in absolute value. While, I wait to fully confirm this answer, I'll post one of my own (of much smaller magnitude)--that I arrived at several hours ago--before seeing this one. My (primitive) laptop seems to be taking a long time for the RandomSearch step. Mar 3 '20 at 10:13
• Important! Important! My computation gives me that the constraint "vcf" is NOT satisfied by the answer given--I get False, not True. Near the end of the posted Mathematica notebook, the answer 0.1565 is correctly given, so the constraint seems to be working properly there. I "apologize" for the hugeness of the constraint and possibly ensuing difficulties in using it. Mar 3 '20 at 10:38
• Also, in the quantum-information-theoretic ("entanglement constraint") context arxiv.org/abs/1708.05336 of the problem, it seems to be a question of whether an answer greater than 1 (like 3.7354) is in fact possible. Mar 3 '20 at 10:57
• Well, I downloaded a copy of the notebook from the wolframcloud, and executed the commands one-by-one and got the indicated answers--so it does seem kosher. The LeafCount for the constraint set in this notebook is 27948, while in my other (not wolframcloud) executions it comes out as 27947--but it doesn't seem to influence the computations. Mar 3 '20 at 11:14
• Well, I now fully understand why this 3.7354 does NOT properly maximize the expression d2. This is because by using vcf[] and not vcf, the LeafCount was reduced from 27948 to 13890, and all the constraints pertaining to the ten b[i] variables were omitted. Thus the analysis above is UNCONSTRAINED with respect to the ten b[i] variables. When I created the wolframcloud object, due to the size of vcf, it seemed I needed to copy-and-paste in two steps, which I then joined. The first insertion pertained to the ten a[i] variables, and the second insertion to the ten b[i] variables. Mar 3 '20 at 14:44

Well, here's something of a (better/greater) "answer" that I arrived at a number of hours before the quite surprising one of Akku14. If his/her answer does stand up, it makes this one look rather "irrelevant". (But a series of comments of mine to his/her answer now show that by employing vcf[], and not vcf, the ten b[i] variables wound up being unconstrained.)

I tried reducing the problem from a 20-variable one to a 10-variable one, by the command Do[b[i] = a[i], {i, 1, 10}]. (It also works "symmetrically" if one uses Do[a[i] = b[i], {i, 1, 10}].) Then, the use of "DifferentialEvolution" gave me a solution (NOW satisfying the constraints) of

{0.184342, {a -> -0.251841, a -> 0.239653, a -> -0.241901, a -> -0.22769, a -> 0.21667, a -> -0.185068, a -> 0.194481, a -> -0.281275, a -> -0.0274505, a -> -0.0274505}},

but still inferior to 0.186769.

However, the command NMaximize[{d2, vcf}, Array[a, 10]] gave

{0.194815, {a -> -0.255072, a -> -0.268045, a -> 0.228979, a -> 0.186145, a -> 0.195611, a -> -0.194521, a -> 0.204414, a -> -0.315437, a -> -0.149244, a -> -0.149244}}

with 0.194815 > 0.186769, also satisfying the constraints. (If I explicitly use Method->"NelderMead", I get the same result, so it seems that this was the default Method employed.)

This line of reasoning does not explain why the original use of "DifferentialEvolution" appeared flawed (in not fully satisfying the constraints), but in giving a greater-valued answer, seemed to make it a less pressing question.

Also, note that in moving from the 20-variable problem to these 10-variable problems, the [now somewhat notorious] set ("vcf") of constraints can be preliminary simplified.

So, can a maximum greater than 0.194815 be found?

Also, the maximization of the indicated

d1=a^2 a^2 a^2 a^2 a^2 a^2 b^2 b^2 b^2 b^2 b^2 b^2 (a a b b - a a b b)^2


is of interest. The maximum appears (from other analyses) to be greater than $$2 \cdot 10^{-21}$$.