#! /usr/bin/env python
#
def r8_dawson ( x ):

#*****************************************************************************80
#
## R8_DAWSON evaluates Dawson's integral of an R8 argument.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    25 April 2016
#
#  Author:
#
#    Original FORTRAN77 version by Wayne Fullerton.
#    Python 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 value of Dawson's integral at X.
#
  import numpy as np
  from r8_csevl import r8_csevl
  from r8_inits import r8_inits
  from machine import r8_mach

  daw2cs = np.array ( [ \
      -0.56886544105215527114160533733674E-01, \
      -0.31811346996168131279322878048822, \
      +0.20873845413642236789741580198858, \
      -0.12475409913779131214073498314784, \
      +0.67869305186676777092847516423676E-01, \
      -0.33659144895270939503068230966587E-01, \
      +0.15260781271987971743682460381640E-01, \
      -0.63483709625962148230586094788535E-02, \
      +0.24326740920748520596865966109343E-02, \
      -0.86219541491065032038526983549637E-03, \
      +0.28376573336321625302857636538295E-03, \
      -0.87057549874170423699396581464335E-04, \
      +0.24986849985481658331800044137276E-04, \
      -0.67319286764160294344603050339520E-05, \
      +0.17078578785573543710504524047844E-05, \
      -0.40917551226475381271896592490038E-06, \
      +0.92828292216755773260751785312273E-07, \
      -0.19991403610147617829845096332198E-07, \
      +0.40963490644082195241210487868917E-08, \
      -0.80032409540993168075706781753561E-09, \
      +0.14938503128761465059143225550110E-09, \
      -0.26687999885622329284924651063339E-10, \
      +0.45712216985159458151405617724103E-11, \
      -0.75187305222043565872243727326771E-12, \
      +0.11893100052629681879029828987302E-12, \
      -0.18116907933852346973490318263084E-13, \
      +0.26611733684358969193001612199626E-14, \
      -0.37738863052129419795444109905930E-15, \
      +0.51727953789087172679680082229329E-16, \
      -0.68603684084077500979419564670102E-17, \
      +0.88123751354161071806469337321745E-18, \
      -0.10974248249996606292106299624652E-18, \
      +0.13261199326367178513595545891635E-19, \
      -0.15562732768137380785488776571562E-20, \
      +0.17751425583655720607833415570773E-21, \
      -0.19695006967006578384953608765439E-22, \
      +0.21270074896998699661924010120533E-23, \
      -0.22375398124627973794182113962666E-24, \
      +0.22942768578582348946971383125333E-25, \
      -0.22943788846552928693329592319999E-26, \
      +0.22391702100592453618342297600000E-27, \
      -0.21338230616608897703678225066666E-28, \
      +0.19866196585123531518028458666666E-29, \
      -0.18079295866694391771955199999999E-30, \
      +0.16090686015283030305450666666666E-31 ] )

  dawacs = np.array ( [ \
      +0.1690485637765703755422637438849E-01, \
      +0.8683252278406957990536107850768E-02, \
      +0.2424864042417715453277703459889E-03, \
      +0.1261182399572690001651949240377E-04, \
      +0.1066453314636176955705691125906E-05, \
      +0.1358159794790727611348424505728E-06, \
      +0.2171042356577298398904312744743E-07, \
      +0.2867010501805295270343676804813E-08, \
      -0.1901336393035820112282492378024E-09, \
      -0.3097780484395201125532065774268E-09, \
      -0.1029414876057509247398132286413E-09, \
      -0.6260356459459576150417587283121E-11, \
      +0.8563132497446451216262303166276E-11, \
      +0.3033045148075659292976266276257E-11, \
      -0.2523618306809291372630886938826E-12, \
      -0.4210604795440664513175461934510E-12, \
      -0.4431140826646238312143429452036E-13, \
      +0.4911210272841205205940037065117E-13, \
      +0.1235856242283903407076477954739E-13, \
      -0.5788733199016569246955765071069E-14, \
      -0.2282723294807358620978183957030E-14, \
      +0.7637149411014126476312362917590E-15, \
      +0.3851546883566811728777594002095E-15, \
      -0.1199932056928290592803237283045E-15, \
      -0.6313439150094572347334270285250E-16, \
      +0.2239559965972975375254912790237E-16, \
      +0.9987925830076495995132891200749E-17, \
      -0.4681068274322495334536246507252E-17, \
      -0.1436303644349721337241628751534E-17, \
      +0.1020822731410541112977908032130E-17, \
      +0.1538908873136092072837389822372E-18, \
      -0.2189157877645793888894790926056E-18, \
      +0.2156879197938651750392359152517E-20, \
      +0.4370219827442449851134792557395E-19, \
      -0.8234581460977207241098927905177E-20, \
      -0.7498648721256466222903202835420E-20, \
      +0.3282536720735671610957612930039E-20, \
      +0.8858064309503921116076561515151E-21, \
      -0.9185087111727002988094460531485E-21, \
      +0.2978962223788748988314166045791E-22, \
      +0.1972132136618471883159505468041E-21, \
      -0.5974775596362906638089584995117E-22, \
      -0.2834410031503850965443825182441E-22, \
      +0.2209560791131554514777150489012E-22, \
      -0.5439955741897144300079480307711E-25, \
      -0.5213549243294848668017136696470E-23, \
      +0.1702350556813114199065671499076E-23, \
      +0.6917400860836148343022185660197E-24, \
      -0.6540941793002752512239445125802E-24, \
      +0.6093576580439328960371824654636E-25, \
      +0.1408070432905187461501945080272E-24, \
      -0.6785886121054846331167674943755E-25, \
      -0.9799732036214295711741583102225E-26, \
      +0.2121244903099041332598960939160E-25, \
      -0.5954455022548790938238802154487E-26, \
      -0.3093088861875470177838847232049E-26, \
      +0.2854389216344524682400691986104E-26, \
      -0.3951289447379305566023477271811E-27, \
      -0.5906000648607628478116840894453E-27, \
      +0.3670236964668687003647889980609E-27, \
      -0.4839958238042276256598303038941E-29, \
      -0.9799265984210443869597404017022E-28, \
      +0.4684773732612130606158908804300E-28, \
      +0.5030877696993461051647667603155E-29, \
      -0.1547395051706028239247552068295E-28, \
      +0.6112180185086419243976005662714E-29, \
      +0.1357913399124811650343602736158E-29, \
      -0.2417687752768673088385304299044E-29, \
      +0.8369074582074298945292887587291E-30, \
      +0.2665413042788979165838319401566E-30, \
      -0.3811653692354890336935691003712E-30, \
      +0.1230054721884951464371706872585E-30, \
      +0.4622506399041493508805536929983E-31, \
      -0.6120087296881677722911435593001E-31, \
      +0.1966024640193164686956230217896E-31 ] )

  dawcs = np.array ( [ \
      -0.6351734375145949201065127736293E-02, \
      -0.2294071479677386939899824125866, \
      +0.2213050093908476441683979161786E-01, \
      -0.1549265453892985046743057753375E-02, \
      +0.8497327715684917456777542948066E-04, \
      -0.3828266270972014924994099521309E-05, \
      +0.1462854806250163197757148949539E-06, \
      -0.4851982381825991798846715425114E-08, \
      +0.1421463577759139790347568183304E-09, \
      -0.3728836087920596525335493054088E-11, \
      +0.8854942961778203370194565231369E-13, \
      -0.1920757131350206355421648417493E-14, \
      +0.3834325867246327588241074439253E-16, \
      -0.7089154168175881633584099327999E-18, \
      +0.1220552135889457674416901120000E-19, \
      -0.1966204826605348760299451733333E-21, \
      +0.2975845541376597189113173333333E-23, \
      -0.4247069514800596951039999999999E-25, \
      +0.5734270767391742798506666666666E-27, \
      -0.7345836823178450261333333333333E-29, \
      +0.8951937667516552533333333333333E-31 ] )

  eps = r8_mach ( 3 )
  ntdaw  = r8_inits ( dawcs,  21, 0.1 * eps )
  ntdaw2 = r8_inits ( daw2cs, 45, 0.1 * eps )
  ntdawa = r8_inits ( dawacs, 75, 0.1 * eps )
  xsml = np.sqrt ( 1.5 * eps )
  xbig = np.sqrt ( 0.5 / eps )
  xmax = np.exp ( min ( - np.log ( 2.0 * r8_mach ( 1 ) ), \
      np.log ( r8_mach ( 2 ) ) ) - 0.01 )

  y = abs ( x )

  if ( y <= xsml ):
    value = x
  elif ( y <= 1.0 ):
    value = x * ( 0.75 + r8_csevl ( 2.0 * y * y - 1.0, dawcs, ntdaw ) )
  elif ( y <= 4.0 ):
    value = x * ( 0.25 + r8_csevl ( 0.125 * y * y - 1.0, daw2cs, ntdaw2 ) )
  elif ( y < xbig ):
    value = ( 0.5 + r8_csevl ( 32.0 / y / y - 1.0, dawacs, ntdawa ) ) / x
  elif ( y <= xmax ):
    value = 0.5 / x
  else:
    value = 0.0

  return value

def r8_dawson_test ( ):

#*****************************************************************************80
#
## R8_DAWSON_TEST tests R8_DAWSON.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    25 April 2016
#
#  Author:
#
#    John Burkardt
#
  import platform
  from dawson_values import dawson_values

  print ( '' )
  print ( 'R8_DAWSON_TEST:' )
  print ( '  Python version: %s' % ( platform.python_version ( ) ) )
  print ( '  R8_DAWSON evaluates Dawson\'s integral.' )
  print ( '' )
  print ( '             X      DAWSON(X)  R8_DAWSON(X)         Diff' )
  print ( '' )

  n_data = 0

  while ( True ):

    n_data, x, fx1 = dawson_values ( n_data )

    if ( n_data == 0 ):
      break

    fx2 = r8_dawson ( x )

    print ( '  %14.4f  %14.6g  %14.6g  %14.6g' % ( x, fx1, fx2, abs ( fx1 - fx2 ) ) )
#
#  Terminate.
#
  print ( '' )
  print ( 'R8_DAWSON_TEST:' )
  print ( '  Normal end of execution.' )
  return

if ( __name__ == '__main__' ):
  from timestamp import timestamp
  timestamp ( )
  r8_dawson_test ( )
  timestamp ( )