Laplace's method - higher order terms

\(TraditionalForm
    \`We  consider  here  the  limit  x → 
      ∞  of  integrals 
        \(∫\_a\%b\( f(t)\) \(e\^\(x \( ϕ(t)\)\)\) 
            dt\), where\ we\ have\ assumed\)\(TraditionalForm\`	
    i.\ \ \ \(ϕ(x)\) takes  its  maximum  at  x = 
      c  ϵ[a, b] \((at  which \(\( ϕ'\)'\)[x] < 0)\), and\)\(TraditionalForm\`	ii.\ \ \(f(c)\) ≠ 0.\)

Test example

\(As\ a\ test\ case\ 
              \((arising\ from\ the\ gamma - function; \ \ 
                Γ(x)\  = 
                  x\^x\ \(∫\_0\%∞ 1\/t\ 
                        e\^\(\(-\ x\)\ \((t - ln\ t)\)\)\ 
                        dt\)\ )\)\ \ \ we\ choose

            f(t) = \(1\/t\ \ and\ \ \(ϕ(t)\) = \(-t\) + ln\ \(t.\)\)\) The maximun of ϕ occurs for  c =1. We also specify how many terms we want in the generated expansion (here nterm  =  
                4. The computer cost  goes up rapidly with the value
                   of nterm,  
              and values higher than 5 or 6 may take an excessive 
                time):

f[t_] := 1/t ;  ϕ[t_] := -t + Log[t] ;  c  =  1 ;   
    nterm  =  4;

Algorithm

\(TraditionalForm
    \`\(The\ algorithm\ generalizes\ to\ arbitrary\ order\ the\ description\ 
          in\ Bender\ &\)\ Orszag, \ pp\ 272 - 273.\)\(TraditionalForm
    \`We\ start\ by\ Taylor\ expanding\ \(f(t)\).\ In\ this\ and\ next\ 
      statement, \ 
    note\ that\ for\ each\ extra\ term\ in\ the\ final\ result, \ 
    we\ need\ two\ more\ terms\ in\ the\ expansion\ for\ \(f(t)\)\ and\ six\ 
      more\ in\ the\ expansion\ of\ \(e\^\(ϕ(t)\) : \)\)

\(\(ft = 	Series[f[t], {t, c, 2*nterm - 2}]; \)\)

Next,  
    we express  (e^(x ϕ(t))) as  
      (e^(x (ϕ[c] + (t - c)^2 ϕ''[c]/2))) * 
      (a  Taylor expansion in  t) :

ex[t_] := x (ϕ[c] + (t - c)^2 ϕ''[c]/2); 
    ϕs  =  
      Series[Exp[x ϕ[t] - ex[t]], {t, c, 6*(nterm - 1)}];

We multiply together the two expansions,  
    and substitute  s  =  (Sqrt[x]) (t - c) :

Normal[ft]*Normal[ϕs]  /.  t -> c + s*x^(-1/2);

\(TraditionalForm
    \`\(The\ integrand\ is\ now\ in\ the\ form\ of\ a\ generalization\ of\ B&
          \) O\ equation\ \((6.4  .34)\).\ The\ next\ step\ is\ \ very\ \ 
        computing - intensive.\ To\ save\ some\ time, \ 
    we\ \ exploit\ that\ all\ the\ integrals\ involved\ have\ a\ relatively\ 
        simple\ closed - form, \ namely\ \ for\ \ i\ \ even, \ 
    ∫\_\(-∞\)\%∞\( s\^i\) 
          e\^\(\(s\^2\) \(ϕ'\)' \((c)\)/2\)\ ds\  = \ 
      \(\(\(i!\) \@\(2  π\)\)\/\(\(2\^\(i/2\)\) \(\((i\/2)\)!\) 
                \((\(-\(ϕ'\)'\) \((c)\))\)\^\(\((i + 1)\)/2\)\)\).\ \ 
          Using\ this\ result\ 
        \((rather\ than\ having\ Mathematica\ spend\ time\ to\ re - 
            evaluate\ these\ integrals\ repeatedly)\), \ we\ get\)

Sum[
      Coefficient[%, s, i]*Sqrt[2*Pi]*
        i!/(
            2^(i/2)*(i/2)!*
              ((-ϕ''[c])^((i + 1)/2))), {
        i, 0, 8*(nterm - 1), 2}];

\(TraditionalForm
    \`It\ now\ only\ remains\ to\ get\ rid\ of\ high\ terms\ that\ are\ not\ 
        exact.\ The\ O[x]\ and\ Normal[\ ]\  - \ 
      commands\ work\ only\ for\ positive\ powers\ of\ x, \ 
    so\ we\ temporarily\ introduce\ \ ix\  = \ 1/x, \ truncate, \ 
    and\ then\ return\ to\ x\ \(again : \)\)

Normal[Expand[%  /. x -> 1/ix] + O[ix]^nterm];

%  /.  ix -> 1/x;

\(TraditionalForm
    \`The\ final\ step\ before\ presenting\ the\ asymptotic\ expansion\ is\ to
      \ multiply\ the\ result\ of\ the\ integration\ \ by\ the\ leading\ 
      factor\ in\ \((6.4  .34)\)\ 
      \(\((which\ we\ have\ omitted\ up\ to\ this\ point)\) : \)\)

Exp[x*ϕ[c]]/Sqrt[x]*%