0
$\begingroup$

This is the mathematica code I used. My goal is trying to solve eq2==0 && eq3==0 simultaneously in the region "-200

The FindRoot doesn't gives me quick answers. So I tried to use Nminimize to do that. This is the code I use. All the input constants are real, but Nminimize gives me a complex value. The documentation also says that "By default, all variables are assumed to be real"

https://reference.wolfram.com/language/ref/NMinimize.html

Clear[u, H, A, J, c, beta, int1, int2, eq1, eq2, eq3, res, subdomain];

A = 1;

J = 1;
c = 1;

beta = 0.001;

int1 = Integrate[E^(c*beta*l*H*(3/2 x^2 - 1/2)), {x, 0, 1}];

int2 = Integrate[(3/2 x^2 - l/2)*E^(c*beta*l*H*(3/2 x^2 - 1/2)), {x, 
    0, 1}];

sum1 = Sum[l^(3/2)*E^(-u*l)*(int1), {l, 1, 400}];

sum2 = Sum[l^(5/2)*E^(-u*l)*(int2), {l, 1, 400}];

eq1 = A*beta^(-3/2)*E^(-J*beta - 1/2*c*beta*H^2);

eq2 = eq1*(sum1) - 1;

eq3 = eq1*(sum2) - H;

subdomain = -200 < H < 400 && 0 < u;

res1 = NMinimize[{Total[{eq2, eq3}^2], subdomain}, {u, H}]

res2 = Last[res1]; resval = First[res1];

{eq1, eq2, eq3, u, H, resval} /. res2

Internal`RealValuedNumericQ /@ {eq1, eq2, eq3, u, H, resval}
$\endgroup$
  • 1
    $\begingroup$ All the numbers I see have 0 imaginary parts. You are using an undocumented internal function here. Do you fully understand what it does, in particular that it returns false for 1. + 0. I as well as for any symbol? $\endgroup$ – Szabolcs Oct 11 '15 at 19:46
1
$\begingroup$

Use increased (arbitrary) precision to avoid precision issues

Clear[u, H, A, J, c, beta, int1, int2,
  eq1, eq2, eq3, res, subdomain];

A = 1; J = 1; c = 1;
beta = 0.001 // Rationalize;
int1 = Integrate[
   E^(c*beta*l*H*(3/2 x^2 - 1/2)),
   {x, 0, 1}];
int2 = Integrate[
   (3/2 x^2 - l/2)*
    E^(c*beta*l*H*(3/2 x^2 - 1/2)),
   {x, 0, 1}];
sum1 = Sum[
   l^(3/2)*E^(-u*l)*(int1),
   {l, 1, 400}];
sum2 = Sum[
   l^(5/2)*E^(-u*l)*(int2),
   {l, 1, 400}];
eq1 = A*beta^(-3/2)*E^(-J*beta - 1/2*c*beta*H^2);
eq2 = eq1*sum1 - 1;
eq3 = eq1*sum2 - H;
subdomain = -200 < H < 400 && 0 < u;

res1 = NMinimize[
   {Total[{eq2, eq3}^2], subdomain},
   {u, H},
   WorkingPrecision -> 30];
res2 = Last[res1]; resval = First[res1];
{eq1, eq2, eq3, u, H, resval} /. res2

(*  {0.010859851424065928700302970525, -1.7528412826622928835*10^-9, \
-7.6728992787359*10^-12, 0.220316517407452903602755920252, \
-172.530087602530220808736082091, 3.07251143558853376444449236742*10^-18}  *)
$\endgroup$
  • $\begingroup$ Thank you very much. This indeed gives the real value. $\endgroup$ – user56134 Oct 11 '15 at 22:25
  • $\begingroup$ When you said "Use increased (arbitrary) precision to avoid precision issues", do you refer to "Working Precision-> 30 " ? Could you elaborate more on that ? How could I eliminate complex value when I increase Working Precision ? $\endgroup$ – user56134 Oct 11 '15 at 22:31
  • $\begingroup$ That plus using Rationalize with definition of beta to convert it to an exact value (1/1000). $\endgroup$ – Bob Hanlon Oct 11 '15 at 22:35
  • 1
    $\begingroup$ The negligible imaginary artifacts arise using machine precision. With greater precision and associated precision tracking the algorithms eliminate the imaginary artifacts. You should also read the documentation for 'Chop`. $\endgroup$ – Bob Hanlon Oct 11 '15 at 23:25
  • $\begingroup$ Thank you so much for the help. I have struggled with this for a long time. Thanks again. $\endgroup$ – user56134 Oct 12 '15 at 0:09

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.