1
$\begingroup$

I was trying to minimize a quartic polyonomial equation with a constraint using Mathematica but I got some conflicting results and I wanted to ask to you. Am I making a conceptual mistake or a coding mistake?

First I wanted to minimize using the Minimize function:

Minimize[pol1[1.9, 0.1], 
  c0r^2 + c0i^2 + c1r^2 + c1i^2 == 1, {c0r, c0i, c1r, c1i}];

As you can see the function called pol1 and the constraint is a unit sphere. Function is a quartic polynomial with funny coefficients:

 Chop[N[pol1[1.9, 0.1]]]
2.23232 c0i^2 - 0.0995481 c0i^4 + 2.23232 c0r^2 - 
 0.199096 c0i^2 c0r^2 - 0.0995481 c0r^4 - 0.433493 c0i c1i + 
 0.126205 c0i^3 c1i + 2.07557 c0r c1i - 1.69557 c0i^2 c0r c1i + 
 0.126205 c0i c0r^2 c1i - 1.69557 c0r^3 c1i + 7.66189 c1i^2 - 
 0.763494 c0i^2 c1i^2 + 1.0748 c0i c0r c1i^2 - 7.94349 c0r^2 c1i^2 + 
 0.458615 c0i c1i^3 - 6.16151 c0r c1i^3 - 1.31455 c1i^4 - 
 2.07557 c0i c1r + 1.69557 c0i^3 c1r - 0.433493 c0r c1r + 
 0.126205 c0i^2 c0r c1r + 1.69557 c0i c0r^2 c1r + 
 0.126205 c0r^3 c1r - 1.0748 c0i^2 c1i c1r + 14.36 c0i c0r c1i c1r + 
 1.0748 c0r^2 c1i c1r + 6.16151 c0i c1i^2 c1r + 
 0.458615 c0r c1i^2 c1r + 7.66189 c1r^2 - 7.94349 c0i^2 c1r^2 - 
 1.0748 c0i c0r c1r^2 - 0.763494 c0r^2 c1r^2 + 
 0.458615 c0i c1i c1r^2 - 6.16151 c0r c1i c1r^2 - 
 2.6291 c1i^2 c1r^2 + 6.16151 c0i c1r^3 + 0.458615 c0r c1r^3 - 
 1.31455 c1r^4

It gave me a real solution.

But then I thought maybe I should do it in the proper way: use Lagrange multiplier method, take derivatives with respect to each parameters, solve the equations and then find the solution that has the global minimum:

    eqs = pol1[1.9, 0.1] + \[Mu]*(c0r^2 + c0i^2 + c1r^2 + c1i^2 - 1)
eq1 = N[D[eqs, c0r]];
eq2 = N[D[eqs, c0i]];
eq3 = N[D[eqs, c1r]];
eq4 = N[D[eqs, c1i]];
eq5 = N[D[eqs, \[Mu]]];

NSolve[eq1 == 0 && eq2 == 0 && eq3 == 0 && eq4 == 0 && eq5 == 0, {c0r,
   c0i, c1r, c1i, \[Mu]}]

Now this second method doesn't even return a real solution set. Where is the mistake?

I will be grateful if you help.

Mathematica version 11.01.0

$\endgroup$
  • 3
    $\begingroup$ Please provide the definition of pol1[...] $\endgroup$ – Ulrich Neumann Aug 7 '18 at 12:48
  • 1
    $\begingroup$ I get real solutions. $\endgroup$ – Daniel Lichtblau Aug 7 '18 at 13:18
  • $\begingroup$ So you got real solutions for the NSolve case right? When you plug them into the original equation do you get 1.97153? Thanks for the comment. $\endgroup$ – Mad_Maximus Aug 7 '18 at 14:12
  • $\begingroup$ (1) Use @ sign in front of name to notify when you comment. (2) Yes I think I got something like that for one solution. Another gave a smaller value, which I could replicate with NMinimize and `Method->"DifferentialEvolution". $\endgroup$ – Daniel Lichtblau Aug 7 '18 at 22:58
  • $\begingroup$ Re 1.97, sometimes exact solvers, such as Minimize, have trouble with approximate numbers (floats) -- and you've only given 6 of the 16 digits of each coefficient, although those rounding errors do not seem numerically significant here. You can improve the result converting them to arbitrary precision numbers: Minimize[SetPrecision[poly, 16], c0r^2 + c0i^2 + c1r^2 + c1i^2 == 1, {c0r, c0i, c1r, c1i}] yields 1.55 (in V11.3). [Aside from giving the exact polynomial you used, it would be considerate to others to indicate in the question why you ask about 1.97153.] $\endgroup$ – Michael E2 Aug 7 '18 at 23:35
1
$\begingroup$

You only need NMinimize to solve your problem!

J=2.23232 c0i^2 - 0.0995481 c0i^4 + 2.23232 c0r^2 - 
0.199096 c0i^2 c0r^2 - 0.0995481 c0r^4 - 0.433493 c0i c1i + 
0.126205 c0i^3 c1i + 2.07557 c0r c1i - 1.69557 c0i^2 c0r c1i + 
0.126205 c0i c0r^2 c1i - 1.69557 c0r^3 c1i + 7.66189 c1i^2 - 
0.763494 c0i^2 c1i^2 + 1.0748 c0i c0r c1i^2 - 7.94349 c0r^2 c1i^2 + 
0.458615 c0i c1i^3 - 6.16151 c0r c1i^3 - 1.31455 c1i^4 - 
2.07557 c0i c1r + 1.69557 c0i^3 c1r - 0.433493 c0r c1r + 
0.126205 c0i^2 c0r c1r + 1.69557 c0i c0r^2 c1r + 
0.126205 c0r^3 c1r - 1.0748 c0i^2 c1i c1r + 14.36 c0i c0r c1i c1r + 
1.0748 c0r^2 c1i c1r + 6.16151 c0i c1i^2 c1r + 
0.458615 c0r c1i^2 c1r + 7.66189 c1r^2 - 7.94349 c0i^2 c1r^2 - 
1.0748 c0i c0r c1r^2 - 0.763494 c0r^2 c1r^2 + 
0.458615 c0i c1i c1r^2 - 6.16151 c0r c1i c1r^2 - 
2.6291 c1i^2 c1r^2 + 6.16151 c0i c1r^3 + 0.458615 c0r c1r^3 -1.31455 c1r^4

sol=NMinimize[J, c0r^2 + c0i^2 + c1r^2 + c1i^2 == 1, {c0r, c0i, c1r, c1i}]
(*{1.97153, {c0r -> 0.948685, c0i -> -0.118812, c1r -> 0.00523234,c1i -> -0.293008}}*)

NSolve doesn't find a solution!

eqn = Map[# == 0 &, D[J, {{c0i, c0r, c1i, c1r}}]];
NSolve[ Join[eqn , {c0r^2 + c0i^2 + c1r^2 + c1i^2 == 1}], {c0i, c0r,c1i, c1r}, Reals]
(*{}*)
$\endgroup$
  • $\begingroup$ Yes I did but why the other method is not working? That is confusing me. It is the proper mathematical method. $\endgroup$ – Mad_Maximus Aug 7 '18 at 13:09
  • $\begingroup$ Thank you very much with the effort but there wasn't supposed to be a difference, if there is then there is a problem with Mathematica. $\endgroup$ – Mad_Maximus Aug 7 '18 at 13:34
  • $\begingroup$ Sorry I didn't see Lagrangeparameter... $\endgroup$ – Ulrich Neumann Aug 7 '18 at 13:38
1
$\begingroup$

Not really an answer, just showing we can solve the system in question.

poly = 2.23232 c0i^2 - 0.0995481 c0i^4 + 2.23232 c0r^2 - 
   0.199096 c0i^2 c0r^2 - 0.0995481 c0r^4 - 0.433493 c0i c1i + 
   0.126205 c0i^3 c1i + 2.07557 c0r c1i - 1.69557 c0i^2 c0r c1i + 
   0.126205 c0i c0r^2 c1i - 1.69557 c0r^3 c1i + 7.66189 c1i^2 - 
   0.763494 c0i^2 c1i^2 + 1.0748 c0i c0r c1i^2 - 
   7.94349 c0r^2 c1i^2 + 0.458615 c0i c1i^3 - 6.16151 c0r c1i^3 - 
   1.31455 c1i^4 - 2.07557 c0i c1r + 1.69557 c0i^3 c1r - 
   0.433493 c0r c1r + 0.126205 c0i^2 c0r c1r + 
   1.69557 c0i c0r^2 c1r + 0.126205 c0r^3 c1r - 
   1.0748 c0i^2 c1i c1r + 14.36 c0i c0r c1i c1r + 
   1.0748 c0r^2 c1i c1r + 6.16151 c0i c1i^2 c1r + 
   0.458615 c0r c1i^2 c1r + 7.66189 c1r^2 - 7.94349 c0i^2 c1r^2 - 
   1.0748 c0i c0r c1r^2 - 0.763494 c0r^2 c1r^2 + 
   0.458615 c0i c1i c1r^2 - 6.16151 c0r c1i c1r^2 - 
   2.6291 c1i^2 c1r^2 + 6.16151 c0i c1r^3 + 0.458615 c0r c1r^3 - 
   1.31455 c1r^4;

vars = Variables[poly];
lagrangian = poly + mu*(vars.vars - 1);
derivs = Grad[lagrangian, Variables[lagrangian]];

solns = NSolve[derivs];
realsolns = Select[solns, FreeQ[#, Complex] &];

Check we made the sum of squares vanish.

In[303]:= (vars.vars - 1) /. realsolns

(* Out[303]= {-1.11022302463*10^-16, 7.77156117238*10^-16, 
 6.66133814775*10^-16, -1.11022302463*10^-16, 0., \
-6.10622663544*10^-16, -2.38524477947*10^-16, -1.82145964978*10^-16, \
-6.39245600897*10^-16, 
 2.08166817117*10^-15, -1.2490009027*10^-16, -2.00811589579*10^-14, \
-1.89015469942*10^-14, -3.88578058619*10^-16, 0., 1.38777878078*10^-16} *)

Check resulting polynomial values.

In[304]:= poly /. realsolns

(* Out[304]= {6.70664718833, 6.70664718833, 6.70664718833, \
6.70664718833, 1.55255414283, 1.55255414283, 2.14374447723, \
2.14374447723, 2.14374443673, 1.9715274388, 1.9715274388, \
1.9715274388, 1.9715274388, 1.97152732129, 1.97152732129, \
1.97152732129} *)

We can get that smallest value using optimization directly.

{min, vals} = 
 Minimize[pol1, 
  c0r^2 + c0i^2 + c1r^2 + c1i^2 == 1, {c0r, c0i, c1r, c1i}, 
  Method -> "DifferentialEvolution"]

(* Out[305]= {1.55255311155, {c0r -> 0.13967123044, 
  c0i -> 0.781692557246, c1r -> -0.602194200876, 
  c1i -> 0.0825305977305}} *)

So what is different? Just the version. For NSolve some of the version 10 releases had trouble getting explicitly real values. Use of Chop at some modest level might suffice to recover reasonable results.

$\endgroup$
  • $\begingroup$ By the way I wrote it in another comment but when I solve it I get 1.97... but then I used the the method Differential evolution and I get your result. $\endgroup$ – Mad_Maximus Aug 8 '18 at 6:44

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.