0
$\begingroup$

I would like to find the solution of an equation containing a slow function

FindRoot[func[x, onsideS, hoppingS, 0.01] - (1/E) == 0, {x, 1,10}]

for different values of onsideS and hoppingS. func is given by

func[tau_?NumericQ, onsideS_?NumericQ, hoppingS_?NumericQ, eta_?NumericQ] := 
    func[tau, onsideS, hoppingS, eta] =
    Min[Exp[eta*tau]*Abs[NIntegrate[Exp[-I*x*tau]*(1/(2 Pi))*((1/(x - onsideS + I*eta)) 
    + ((1/(x - onsideS + I*eta))^2)*(hoppingS^2)*localGF[x,eta]*
    (1/(1 - (hoppingS^2) * (1/(x - onsideS + I*eta))*localGF[x, eta]))), {x, -Infinity, Infinity},  
    Method -> {"LevinRule", "Kernel" -> Exp[-I*x*tau]}], 1]

with

localGF[energy_?NumericQ, eta_?NumericQ] := localGF[energy, eta] =  
    NIntegrate[ldos[x]/(energy - x + I*eta), {x, -hbw, hbw}, 
    Method -> {"GaussKronrodRule", "Points" -> 25,"SymbolicProcessing" -> 0}]
ldos[energy_?NumericQ] := (1/(Pi))*Sqrt[1 - (energy/2)^2]
hbw = 2;

However Mathematica always returns a bunch of error messages containing of

NIntegrate: Integrand is not a Levin function

FindRoot: Failed to converge to the requested accuracy or precision within 100 iterations

and FindRoot never gives a valid result, even though calling func manually returns something useful. The same happens when trying to Plot this function. I tried using NSolve aswell

NSolve[func[x, onsideS, hoppingS, 0.01] - 1/E == 0, x, Method -> {Automatic, "SymbolicProcessing" -> 0}]

and it works atleast for certain values (i.e. onsideS = 0 and hoppingS = 1) quite fast, but evaluation time increases very quickly when going to other values (i.e. onsideS = 0 and hoppingS = 0.05).

Why does this happen and how could I fix this?

$\endgroup$
  • $\begingroup$ what value have hbw? $\endgroup$ – Mariusz Iwaniuk Jul 7 '18 at 17:52
  • $\begingroup$ @MariuszIwaniuk Sorry, I forgot adding it into the question. hbw = 2; $\endgroup$ – StirriX Jul 7 '18 at 18:17
  • $\begingroup$ The error message you get suggests that localGF has not been passed into the FindRoot properly, I can't see why immediately though. $\endgroup$ – KraZug Jul 8 '18 at 6:03
2
$\begingroup$

This is not a complete answer, but some hints.

First, do as much integrations, you can, analytically.

hbw = 2;

ldos[energy_] = (1/(Pi))*Sqrt[1 - (energy/2)^2] // Simplify;

localGF[energy_, eta_] = 
   Integrate[ldos[x]/(energy - x + I*eta), {x, -hbw, hbw}, 
         Assumptions -> eta > 0 && energy \[Element] Reals]

(*   -(1/2) (-1 + Sqrt[1 - 4/(energy + I eta)^2]) (energy + I eta)   *)

Then regard the integrand for NIntegrate at different parameters.

integrand1[tau_, onsideS_, hoppingS_, eta_, x_] = 
    Exp[-I*x*tau]*(1/(2 Pi))*((1/(x - onsideS + 
     I*eta)) + ((1/(x - onsideS + I*eta))^2)*(hoppingS^2)*
    localGF[x, 
     eta]*(1/(1 - (hoppingS^2)*(1/(x - onsideS + I*eta))*
       localGF[x, eta]))) // Simplify;

Manipulate[
  int = integrand1[tau, onsideS, hoppingS, 1/100, x]; {Plot[
   Re@int, {x, -10, 10}], Plot[Im@int, {x, -10, 10}]}, {tau, 1, 
   10}, {onsideS, 0, 4}, {hoppingS, 0, 2}]

enter image description here

For some parameter combinations, NIntegrate works quite good

NIntegrate[integrand1[1, 0, 1, 1/100, x], {x, -Infinity, Infinity}]

(*   2.77556*10^-17 - 0.57145 I   *)

NIntegrate[integrand1[2, 0, 1, 1/100, x], {x, -Infinity, Infinity}]

(*   0. + 0.0321347 I   *)

For others, like shown in picture, there is a non integrable singularity at zero or other x-values, you have to use Method -> "PrincipalValue" an explicitly give the x-value of it in integration range.

NIntegrate[integrand1[1, 0, 0, 1/100, x], {x, -Infinity, 0, Infinity},
   Method -> "PrincipalValue"]

(*   0. - 0.99005 I   *)

For other parameter values, integral failes to converge.

NIntegrate[
  integrand1[1, 0, 2/10, 1/100, x], {x, -Infinity, 0, Infinity}, 
  Method -> "PrincipalValue"]

It's a tricky job to further examine integrand behavior and positions of singularities.

$\endgroup$

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.