0
$\begingroup$

I have the following problem:

I have the following function of 3 parameters A, B, mu: $$ g(A,B,\mu) = \int_{BZ} dk_1 dk_2\sqrt{(\mu+ F (\cos k_1 + \cos k_2 + \cos (k_1+k_2)))^2 - A^2 (\sin k_1 + \sin k_2 - \sin (k_1+k_2))^2} - 1.28 \mu + 3A^2 - 3 F^2 $$ where $BZ$ denotes the regular hexagon

H =  Polygon[ 4 \[Pi]/ 
 3 {{1, 0}, {1/2, Sqrt[3]/2}, {-1/2, Sqrt[3]/2}, {-1, 
  0}, {-1/2, -Sqrt[3]/2}, {1/2, -Sqrt[3]/2}}];

The expression under integral becomes purely imaginary for some choices of $A,F,\mu$. My aim is to find the local maximum of this function (where it has real values). I know, that it should be at the point $(A,F,\mu) = (0.181,-0.0528,0.4025)$. But using

 FindMaximum 
gives an error (becouse the procedure gets gomplex values of the function). Is it possible to tackle this problem?

My code is the following:

Ns = 100;  J = 1 ; kappa = 0.28;
H =  Polygon[ 
   4 \[Pi]/ 
     3 {{1, 0}, {1/2, Sqrt[3]/2}, {-1/2, Sqrt[3]/2}, {-1, 
      0}, {-1/2, -Sqrt[3]/2}, {1/2, -Sqrt[3]/2}}];
pts = RandomPoint[H, Ns^2];
ptsKspace = {#[[1]], -1/2 #[[1]] + (Sqrt[3]/2) #[[2]]} & /@ pts;

f[A_, F_, mu_, k1_, k2_] = 
  Sqrt[(mu + F (Cos[k1] + Cos[k2] + Cos[k1 + k2]))^2 - 
    A^2 (Sin[k1] + Sin[k2] - Sin[k1 + k2])^2];

En3 [A_, F_, mu_] := Module[{re, im},
   re = Sqrt[3]/(8 \[Pi]^2) Integrate[
       Re[f[A, F, mu, k1, -1/2 k1 + Sqrt[3]/2 k2] ], {k1, 
         k2} \[Element] H] - mu (1 + kappa) + 3  Abs[A]^2  - 
     3  Abs[F]^2;
   im = Sqrt[3]/(8 \[Pi]^2) Integrate[
      Im[f[A, F, mu, k1, -1/2 k1 + Sqrt[3]/2 k2] ], {k1, 
        k2} \[Element] H] ;
   Return[re + I im]; 
];

FindMaximum[En3[x, y, z], {x, 0.181}, {y, -0.05}, {z, 0.45}]

And I ge an error. The firs lines with random points from hexagon are for possible Monte carlo calculation of integral.

$\endgroup$
7
  • 3
    $\begingroup$ Please show how exactly you used FindMaximum, i.e. post the code by editing the question. $\endgroup$
    – corey979
    Mar 14, 2017 at 19:21
  • 1
    $\begingroup$ Also include the code for the definition of the function. $\endgroup$
    – Bob Hanlon
    Mar 14, 2017 at 19:41
  • $\begingroup$ Recode the objective so that it replaces complex values by a suitably large real value. $\endgroup$ Mar 14, 2017 at 22:56
  • $\begingroup$ @Daniel Lichtblau Thank you for the answer - I tried to set the value -100000 for non - zero imaginary part, but then I got an error : function value -Null is not a real number at {x,y,z} = \ {0.181,-0.05,0.45} $\endgroup$
    – HBoson
    Mar 15, 2017 at 10:19
  • $\begingroup$ Instead of Return[re + I im]; I will recommend e + I im, that is, no Return and no semicolon. See if that gets rid of the Null problem. $\endgroup$ Mar 15, 2017 at 14:46

1 Answer 1

2
$\begingroup$

I redefined so that the numeric functions only evaluate for explicitly numeric input, and the objective returns a negative value well below the max for cases where there is a nontrivial imaginary part.

ns = 100; J = 1; kappa = 0.28;
hh = Polygon[
   4 \[Pi]/3 {{1, 0}, {1/2, Sqrt[3]/2}, {-1/2, Sqrt[3]/2}, {-1, 
      0}, {-1/2, -Sqrt[3]/2}, {1/2, -Sqrt[3]/2}}];
pts = RandomPoint[hh, ns^2];
ptsKspace = {#[[1]], -1/2 #[[1]] + (Sqrt[3]/2) #[[2]]} & /@ pts;

f[a_?NumberQ, f_?NumberQ, mu_?NumberQ, k1_?NumberQ, 
  k2_?NumberQ] := 
 Sqrt[(mu + f (Cos[k1] + Cos[k2] + Cos[k1 + k2]))^2 - 
   a^2 (Sin[k1] + Sin[k2] - Sin[k1 + k2])^2]

en3[aA_?NumberQ, fF_?NumberQ, mu_?NumberQ] := 
 Module[{re, im}, 
  re = Sqrt[3]/(8 \[Pi]^2) NIntegrate[
      Re[f[aA, fF, mu, k1, -1/2 k1 + Sqrt[3]/2 k2]], {k1, 
        k2} \[Element] hh] - mu (1 + kappa) + 3 Abs[aA]^2 - 
    3 Abs[fF]^2;
  im = NIntegrate[
    Im[f[aA, fF, mu, k1, -1/2 k1 + Sqrt[3]/2 k2]], {k1, 
      k2} \[Element] hh];
  If[Abs[im] > .0001, -100, re]
  ]

This is slow and gives warnings but it also delivers a result after several minutes of deep cogitation.

In[21]:= FindMaximum[en3[x, y, z], {x, 0.181}, {y, -0.05}, {z, 0.45}]

(* During evaluation of In[21]:= NIntegrate::izero: Integral and error estimates are 0 on all integration subregions. Try increasing the value of the MinRecursion option. If value of integral may be 0, specify a finite value for the AccuracyGoal option.

During evaluation of In[21]:= NIntegrate::izero: Integral and error estimates are 0 on all integration subregions. Try increasing the value of the MinRecursion option. If value of integral may be 0, specify a finite value for the AccuracyGoal option.

During evaluation of In[21]:= NIntegrate::izero: Integral and error estimates are 0 on all integration subregions. Try increasing the value of the MinRecursion option. If value of integral may be 0, specify a finite value for the AccuracyGoal option.

During evaluation of In[21]:= General::stop: Further output of NIntegrate::izero will be suppressed during this calculation.

During evaluation of In[21]:= NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small.

During evaluation of In[21]:= NIntegrate::eincr: The global error of the strategy GlobalAdaptive has increased more than 2000 times. The global error is expected to decrease monotonically after a number of integrand evaluations. Suspect one of the following: the working precision is insufficient for the specified precision goal; the integrand is highly oscillatory or it is not a (piecewise) smooth function; or the true value of the integral is 0. Increasing the value of the GlobalAdaptive option MaxErrorIncreases might lead to a convergent numerical integration. NIntegrate obtained 10.077065638905038` and 0.00025052670223458276` for the integral and error estimates.

During evaluation of In[21]:= NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small.

During evaluation of In[21]:= NIntegrate::eincr: The global error of the strategy GlobalAdaptive has increased more than 2000 times. The global error is expected to decrease monotonically after a number of integrand evaluations. Suspect one of the following: the working precision is insufficient for the specified precision goal; the integrand is highly oscillatory or it is not a (piecewise) smooth function; or the true value of the integral is 0. Increasing the value of the GlobalAdaptive option MaxErrorIncreases might lead to a convergent numerical integration. NIntegrate obtained 7.2069039646135185` and 0.0002487749558483572` for the integral and error estimates.

During evaluation of In[21]:= NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small.

During evaluation of In[21]:= General::stop: Further output of NIntegrate::slwcon will be suppressed during this calculation.

During evaluation of In[21]:= NIntegrate::eincr: The global error of the strategy GlobalAdaptive has increased more than 2000 times. The global error is expected to decrease monotonically after a number of integrand evaluations. Suspect one of the following: the working precision is insufficient for the specified precision goal; the integrand is highly oscillatory or it is not a (piecewise) smooth function; or the true value of the integral is 0. Increasing the value of the GlobalAdaptive option MaxErrorIncreases might lead to a convergent numerical integration. NIntegrate obtained 10.077065638905038` and 0.00025052670223458276` for the integral and error estimates.

During evaluation of In[21]:= General::stop: Further output of NIntegrate::eincr will be suppressed during this calculation.

During evaluation of In[21]:= FindMaximum::lstol: The line search decreased the step size to within the tolerance specified by AccuracyGoal and PrecisionGoal but was unable to find a sufficient increase in the function. You may need more than MachinePrecision digits of working precision to meet these tolerances.

Out[21]= {-0.0860541069032, {x -> 0.247531713662, 
  y -> -0.103777453975, z -> 0.487104944188}} *)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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