function value = r8_aie ( x )

%*****************************************************************************80
%
%% R8_AIE evaluates the exponentially scaled Airy function Ai of an R8 argument.
%
%  Discussion:
%
%    if X <= 0,
%      R8_AIE ( X ) = R8_AI ( X )
%    else
%      R8_AIE ( X ) = R8_AI ( X ) * exp ( 2/3 * X^(3/2) )
%
%    Thanks to Aleksandra Piper for pointing out a correction involving a
%    missing assignment to SQRTX, 27 January 2012.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    27 January 2012
%
%  Author:
%
%    Original FORTRAN77 version by Wayne Fullerton.
%    MATLAB version by John Burkardt.
%
%  Reference:
%
%    Wayne Fullerton,
%    Portable Special Function Routines,
%    in Portability of Numerical Software,
%    edited by Wayne Cowell,
%    Lecture Notes in Computer Science, Volume 57,
%    Springer 1977,
%    ISBN: 978-3-540-08446-4,
%    LC: QA297.W65.
%
%  Parameters:
%
%    Input, real X, the argument.
%
%    Output, real VALUE, the exponentially scaled Airy
%    function Ai of X.
%
  persistent aifcs
  persistent aigcs
  persistent aip1cs
  persistent aip2cs
  persistent naif
  persistent naig
  persistent naip1
  persistent naip2
  persistent x32sml
  persistent x3sml
  persistent xbig

  if ( isempty ( naif ) )

    aifcs = [ ...
      -0.37971358496669997496197089469414E-01, ...
      +0.59191888537263638574319728013777E-01, ...
      +0.98629280577279975365603891044060E-03, ...
      +0.68488438190765667554854830182412E-05, ...
      +0.25942025962194713019489279081403E-07, ...
      +0.61766127740813750329445749697236E-10, ...
      +0.10092454172466117901429556224601E-12, ...
      +0.12014792511179938141288033225333E-15, ...
      +0.10882945588716991878525295466666E-18, ...
      +0.77513772196684887039238400000000E-22, ...
      +0.44548112037175638391466666666666E-25, ...
      +0.21092845231692343466666666666666E-28, ...
      +0.83701735910741333333333333333333E-32 ]';
    aigcs = [ ...
      +0.18152365581161273011556209957864E-01, ...
      +0.21572563166010755534030638819968E-01, ...
      +0.25678356987483249659052428090133E-03, ...
      +0.14265214119792403898829496921721E-05, ...
      +0.45721149200180426070434097558191E-08, ...
      +0.95251708435647098607392278840592E-11, ...
      +0.13925634605771399051150420686190E-13, ...
      +0.15070999142762379592306991138666E-16, ...
      +0.12559148312567778822703205333333E-19, ...
      +0.83063073770821340343829333333333E-23, ...
      +0.44657538493718567445333333333333E-26, ...
      +0.19900855034518869333333333333333E-29, ...
      +0.74702885256533333333333333333333E-33 ]';
    aip1cs = [ ...
      -0.2146951858910538455460863467778E-01, ...
      -0.7535382535043301166219720865565E-02, ...
      +0.5971527949026380852035388881994E-03, ...
      -0.7283251254207610648502368291548E-04, ...
      +0.1110297130739299666517381821140E-04, ...
      -0.1950386152284405710346930314033E-05, ...
      +0.3786973885159515193885319670057E-06, ...
      -0.7929675297350978279039072879154E-07, ...
      +0.1762247638674256075568420122202E-07, ...
      -0.4110767539667195045029896593893E-08, ...
      +0.9984770057857892247183414107544E-09, ...
      -0.2510093251387122211349867730034E-09, ...
      +0.6500501929860695409272038601725E-10, ...
      -0.1727818405393616515478877107366E-10, ...
      +0.4699378842824512578362292872307E-11, ...
      -0.1304675656297743914491241246272E-11, ...
      +0.3689698478462678810473948382282E-12, ...
      -0.1061087206646806173650359679035E-12, ...
      +0.3098414384878187438660210070110E-13, ...
      -0.9174908079824139307833423547851E-14, ...
      +0.2752049140347210895693579062271E-14, ...
      -0.8353750115922046558091393301880E-15, ...
      +0.2563931129357934947568636168612E-15, ...
      -0.7950633762598854983273747289822E-16, ...
      +0.2489283634603069977437281175644E-16, ...
      -0.7864326933928735569664626221296E-17, ...
      +0.2505687311439975672324470645019E-17, ...
      -0.8047420364163909524537958682241E-18, ...
      +0.2604097118952053964443401104392E-18, ...
      -0.8486954164056412259482488834184E-19, ...
      +0.2784706882142337843359429186027E-19, ...
      -0.9195858953498612913687224151354E-20, ...
      +0.3055304318374238742247668225583E-20, ...
      -0.1021035455479477875902177048439E-20, ...
      +0.3431118190743757844000555680836E-21, ...
      -0.1159129341797749513376922463109E-21, ...
      +0.3935772844200255610836268229154E-22, ...
      -0.1342880980296717611956718989038E-22, ...
      +0.4603287883520002741659190305314E-23, ...
      -0.1585043927004064227810772499387E-23, ...
      +0.5481275667729675908925523755008E-24, ...
      -0.1903349371855047259064017948945E-24, ...
      +0.6635682302374008716777612115968E-25, ...
      -0.2322311650026314307975200986453E-25, ...
      +0.8157640113429179313142743695359E-26, ...
      -0.2875824240632900490057489929557E-26, ...
      +0.1017329450942901435079714319018E-26, ...
      -0.3610879108742216446575703490559E-27, ...
      +0.1285788540363993421256640342698E-27, ...
      -0.4592901037378547425160693022719E-28, ...
      +0.1645597033820713725812102485333E-28, ...
      -0.5913421299843501842087920271360E-29, ...
      +0.2131057006604993303479369509546E-29, ...
      -0.7701158157787598216982761745066E-30, ...
      +0.2790533307968930417581783777280E-30, ...
      -0.1013807715111284006452241367039E-30, ...
      +0.3692580158719624093658286216533E-31 ]';
    aip2cs = [ ...
      -0.174314496929375513390355844011E-02, ...
      -0.167893854325541671632190613480E-02, ...
      +0.359653403352166035885983858114E-04, ...
      -0.138081860273922835457399383100E-05, ...
      +0.741122807731505298848699095233E-07, ...
      -0.500238203900133013130422866325E-08, ...
      +0.400693917417184240675446866355E-09, ...
      -0.367331242795905044199318496207E-10, ...
      +0.376034439592373852439592002918E-11, ...
      -0.422321332718747538026564938968E-12, ...
      +0.513509454033657070919618754120E-13, ...
      -0.669095850390477595651681356676E-14, ...
      +0.926667545641290648239550724382E-15, ...
      -0.135514382416070576333397356591E-15, ...
      +0.208115496312830995299006549335E-16, ...
      -0.334116499159176856871277570256E-17, ...
      +0.558578584585924316868032946585E-18, ...
      -0.969219040152365247518658209109E-19, ...
      +0.174045700128893206465696557738E-19, ...
      -0.322640979731130400247846333098E-20, ...
      +0.616074471106625258533259618986E-21, ...
      -0.120936347982490059076420676266E-21, ...
      +0.243632763310138108261570095786E-22, ...
      -0.502914221497457468943403144533E-23, ...
      +0.106224175543635689495470626133E-23, ...
      -0.229284284895989241509856324266E-24, ...
      +0.505181733929503744986884778666E-25, ...
      -0.113498123714412404979793920000E-25, ...
      +0.259765565985606980698374144000E-26, ...
      -0.605124621542939506172231679999E-27, ...
      +0.143359777966772800720295253333E-27, ...
      -0.345147757060899986280721066666E-28, ...
      +0.843875190213646740427025066666E-29, ...
      -0.209396142298188169434453333333E-29, ...
      +0.527008873478945503182848000000E-30, ...
      -0.134457433014553385789030399999E-30, ...
      +0.347570964526601147340117333333E-31 ]';

    eta = 0.1 * r8_mach ( 3 );
    naif = r8_inits ( aifcs, 13, eta );
    naig = r8_inits ( aigcs, 13, eta );
    naip1 = r8_inits ( aip1cs, 57, eta );
    naip2 = r8_inits ( aip2cs, 37, eta );
    x3sml = eta^0.3333;
    x32sml = 1.3104 * x3sml * x3sml;
    xbig = r8_mach ( 2 )^0.6666;

  end

  if ( x < - 1.0 )
    [ xm, theta ] = r8_aimp ( x );
    value = xm * cos ( theta );
  elseif ( 0.0 <= x && x <= x32sml )
    z = 0.0;
    value = 0.3750 + ( r8_csevl ( z, aifcs, naif ) ...
      - x * ( 0.25 + r8_csevl ( z, aigcs, naig ) ) );
  elseif ( abs ( x ) <= x3sml )
    z = 0.0;
    value = 0.3750 + ( r8_csevl ( z, aifcs, naif ) ...
      - x * ( 0.25 + r8_csevl ( z, aigcs, naig ) ) );
    value = value * exp ( 2.0 * x * sqrt ( x ) / 3.0 );
  elseif ( x <= 1.0 )
    z = x * x * x;
    value = 0.3750 + ( r8_csevl ( z, aifcs, naif ) ...
      - x * ( 0.25 + r8_csevl ( z, aigcs, naig ) ) );
    value = value * exp ( 2.0 * x * sqrt ( x ) / 3.0 );
  elseif ( x <= 4.0 )
    sqrtx = sqrt ( x );
    z = ( 16.0 / ( x * sqrtx ) - 9.0 ) / 7.0;
    value = ( 0.28125 + r8_csevl ( z, aip1cs, naip1 ) ) ...
      / sqrt ( sqrtx );
  elseif ( x < xbig )
    sqrtx = sqrt ( x );
    z = 16.0 / ( x * sqrtx ) - 1.0;
    value = ( 0.28125 + r8_csevl ( z, aip2cs, naip2 ) ) ...
      / sqrt ( sqrtx );
  else
    sqrtx = sqrt ( x );
    z = - 1.0;
    value = ( 0.28125 + r8_csevl ( z, aip2cs, naip2 ) ) ...
      / sqrt ( sqrtx );
  end

  return
end