#! /usr/bin/env python
#
def hermite_probabilist_set ( n ):

#*****************************************************************************80
#
## HERMITE_PROBABILIST_SET: probabilist Hermite quadrature.
#
#  Discussion:
#
#    The integral:
#
#      integral ( -oo < x < +oo ) f(x) * rho(x) dx
#
#    The weight:
#
#      rho(x) = exp ( - x * x / 2 ) / sqrt ( 2 * pi ) dx
#
#    The quadrature rule:
#
#      sum ( 1 <= i <= n ) w(i) * f ( x(i) )
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    18 June 2015
#
#  Author:
#
#    John Burkardt
#
#  Reference:
#
#    Stephen Wolfram,
#    The Mathematica Book,
#    Fourth Edition,
#    Cambridge University Press, 1999,
#    ISBN: 0-521-64314-7,
#    LC: QA76.95.W65.
#
#  Parameters:
#
#    Input, integer N, the order.
#
#    Output, real X(N), the abscissas.
#
#    Output, real W(N), the weights.
#
  import numpy as np
  from sys import exit

  if ( n == 1 ):

    x = np.array ( [ \
           0.0 ] )

    w = np.array ( [ \
           1.00000000000000000000000000000 ] )

  elif ( n == 2 ):

    x = np.array ( [ \
           -1.00000000000000000000000000000, \
           +1.00000000000000000000000000000 ] )

    w = np.array ( [ \
           0.50000000000000000000000000000, \
           0.50000000000000000000000000000 ] )

  elif ( n == 3 ):

    x = np.array ( [ \
           -1.73205080756887729352744634151, \
            0.0, \
           +1.73205080756887729352744634151 ] )

    w = np.array ( [ \
           0.166666666666666666666666666667, \
           0.666666666666666666666666666667, \
           0.166666666666666666666666666667 ] )

  elif ( n == 4 ):

    x = np.array ( [ \
           -2.33441421833897723931751226721, \
           -0.741963784302725857648513596726, \
           +0.741963784302725857648513596726, \
           +2.33441421833897723931751226721 ] )

    w = np.array ( [ \
           0.0458758547680684918168929937745, \
           0.454124145231931508183107006225, \
           0.454124145231931508183107006225, \
           0.0458758547680684918168929937745 ] )

  elif ( n == 5 ):

    x = np.array ( [ \
           -2.85697001387280565416230426401, \
           -1.35562617997426586583052129087, \
            0.0, \
           +1.35562617997426586583052129087, \
           +2.85697001387280565416230426401 ] )

    w = np.array ( [ \
           0.0112574113277206889333702151856, \
           0.222075922005612644399963118148, \
           0.533333333333333333333333333333, \
           0.222075922005612644399963118148, \
           0.0112574113277206889333702151856 ] )

  elif ( n == 6 ):

    x = np.array ( [ \
           -3.32425743355211895236183546247, \
           -1.88917587775371067550566789858, \
           -0.616706590192594152193686099399, \
           +0.616706590192594152193686099399, \
           +1.88917587775371067550566789858, \
           +3.32425743355211895236183546247 ] )

    w = np.array ( [ \
           0.00255578440205624643060629074383, \
           0.0886157460419145274808558830057, \
           0.40882846955602922608853782625, \
           0.40882846955602922608853782625, \
           0.0886157460419145274808558830057, \
           0.00255578440205624643060629074383 ] )

  elif ( n == 7 ):

    x = np.array ( [ \
           -3.75043971772574225630392202571, \
           -2.36675941073454128861885646856, \
           -1.15440539473996812723959775884, \
            0.0, \
           +1.15440539473996812723959775884, \
           +2.36675941073454128861885646856, \
           +3.75043971772574225630392202571 ] )

    w = np.array ( [ \
           0.000548268855972217791621570532802, \
           0.0307571239675864970396450057164, \
           0.240123178605012713740161995179, \
           0.457142857142857142857142857143, \
           0.240123178605012713740161995179, \
           0.0307571239675864970396450057164, \
           0.000548268855972217791621570532802 ] )

  elif ( n == 8 ):

    x = np.array ( [ \
           -4.14454718612589433206019783917, \
           -2.80248586128754169911301080618, \
           -1.63651904243510799922544657297, \
           -0.539079811351375108072461918694, \
           +0.539079811351375108072461918694, \
           +1.63651904243510799922544657297, \
           +2.80248586128754169911301080618, \
           +4.14454718612589433206019783917 ] )

    w = np.array ( [ \
           0.00011261453837536777039380201687, \
           0.00963522012078826718691913771988, \
           0.117239907661759015117137525962, \
           0.373012257679077349925549534301, \
           0.373012257679077349925549534301, \
           0.117239907661759015117137525962, \
           0.00963522012078826718691913771988, \
           0.00011261453837536777039380201687 ] )

  elif ( n == 9 ):

    x = np.array ( [ \
           -4.51274586339978266756667884317, \
           -3.20542900285646994336567590292, \
           -2.07684797867783010652215614374, \
           -1.02325566378913252482814822581, \
            0.0, \
           +1.02325566378913252482814822581, \
           +2.07684797867783010652215614374, \
           +3.20542900285646994336567590292, \
           +4.51274586339978266756667884317 ] )

    w = np.array ( [ \
           0.0000223458440077465836484639907118, \
           0.00278914132123176862881344575164, \
           0.0499164067652178740433414693826, \
           0.2440975028949394361410220177, \
           0.406349206349206349206349206349, \
           0.2440975028949394361410220177, \
           0.0499164067652178740433414693826, \
           0.00278914132123176862881344575164, \
           0.0000223458440077465836484639907118 ] )

  elif ( n == 10 ):

    x = np.array ( [ \
            -4.85946282833231215015516494660, \
            -3.58182348355192692277623675546, \
            -2.48432584163895458087625118368, \
            -1.46598909439115818325066466416, \
            -0.484935707515497653046233483105, \
            +0.484935707515497653046233483105, \
            +1.46598909439115818325066466416, \
            +2.48432584163895458087625118368, \
            +3.58182348355192692277623675546, \
            +4.85946282833231215015516494660 ] )

    w = np.array ( [ \
            0.0000043106526307182867322209547262, \
            0.000758070934312217670069636036508, \
            0.0191115805007702856047383687629, \
            0.135483702980267735563431657727, \
            0.344642334932019042875028116518, \
            0.344642334932019042875028116518, \
            0.135483702980267735563431657727, \
            0.0191115805007702856047383687629, \
            0.000758070934312217670069636036508, \
            0.0000043106526307182867322209547262 ] )

  else:

    print ( '' )
    print ( 'HERMITE_PROBABILIST_SET - Fatal error!' )
    print ( '  Illegal value of N = %d' % ( n ) )
    print ( '  Legal values are 1:10.' )
    exit ( 'HERMITE_PROBABILIST_SET - Fatal error!' )

  return x, w

def hermite_probabilist_set_test ( ):

#*****************************************************************************80
#
## HERMITE_PROBABILIST_SET_TEST tests HERMITE_PROBABILIST_SET.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license. 
#
#  Modified:
#
#    18 June 2015
#
#  Author:
#
#    John Burkardt
#
  import platform

  print ( '' )
  print ( 'HERMITE_PROBABILIST_SET_TEST' )
  print ( '  Python version: %s' % ( platform.python_version ( ) ) )
  print ( '  HERMITE_PROBABILIST_SET sets a Hermite quadrature rule.' )
  print ( '  The integration interval is ( -oo, +oo ).' )
  print ( '  The weight is exp ( - x * x / 2 ) / sqrt ( 2 * pi ).' )

  print ( '' )
  print ( '  Index       X             W' )
  print ( '' )

  for n in range ( 1, 11 ):

    x, w = hermite_probabilist_set ( n )

    print ( '' )

    for i in range ( 0, n ):
      print ( '  %2d  %24.16g  %24.16g' % ( i, x[i], w[i] ) )
#
#  Terminate.
#
  print ( '' )
  print ( 'HERMITE_PROBABILIST_SET_TEST:' )
  print ( '  Normal end of execution.' )
  return

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