## Description

This function evaluates the integral of the given function between the lower and upper limits using the weight and abscissa values specified in the rule data frame. The quadrature formula uses the weight function for generalized Hermite polynomials.

## Usage

 1 2 ghermite.h.quadrature(functn, rule, mu = 0, lower = -Inf, upper = Inf, weighted = TRUE, ...) 

## Arguments

 functn an R function which should take a numeric argument x and possibly some parameters. The function returns a numerical vector value for the given argument x. rule a data frame containing the order n generalized Hermite quadrature rule mu numeric value for the parameter for the generalized Hermite polynomials lower numeric value for the lower limit of the integral with a default value of - ∞ upper numeric value for the upper limit of the integral with a default value of + ∞ weighted a boolean value which if true causes the Hermite weight function to be included in the integrand ... other arguments passed to the give function

## Details

The rule argument corresponds to an order n generalized Hermite polynomial, weight function and interval ≤ft( { - ∞ ,∞ } \right) . The lower and upper limits of the integral must be infinite.

## Value

The value of definite integral evaluated using Gauss Hermite quadrature

## Author(s)

Frederick Novomestky [email protected]

## References

Abramowitz, M. and I. A. Stegun, 1968. Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables, Dover Publications, Inc., New York.

Press, W. H., S. A. Teukolsky, W. T. Vetterling, and B. P. Flannery, 1992. Numerical Recipes in C, Cambridge University Press, Cambridge, U.K.

Stroud, A. H., and D. Secrest, 1966. Gaussian Quadrature Formulas, Prentice-Hall, Englewood Cliffs, NJ.

hermite.h.quadrature.rules

## Examples

  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 ### ### this example evaluates the quadrature function for ### the generalized Hermite H polynomials. it computes the integral ### of the product for all pairs of orthogonal polynomials ### from order 0 to order 10. the results are compared to ### the diagonal matrix of the inner products for the ### polynomials. it also computes the integral of the product ### of all pairs of orthonormal polynomials from order 0 ### to order 10. the resultant matrix should be an identity matrix ### ### ### set the polynomial parameter ### mu <- 1 ### ### set the value for the maximum polynomial order ### n <- 10 ### ### maximum order plus 1 ### np1 <- n + 1 ### ### function to construct the polynomial products by column ### by.column.products <- function( c, p.list, p.p.list ) { ### ### function to construct the polynomial products by row ### by.row.products <- function( r, c, p.list ) { row.column.product <- p.list[[r]] * p.list[[c]] return (row.column.product ) } np1 <- length( p.list ) row.list <- lapply( 1:np1, by.row.products, c, p.list ) return( row.list ) } ### ### function construct the polynomial functions by column ### by.column.functions <- function( c, p.p.products ) { ### ### function to construct the polynomial functions by row ### by.row.functions <- function( r, c, p.p.products ) { row.column.function <- as.function( p.p.products[[r]][[c]] ) return( row.column.function ) } np1 <- length( p.p.products[[1]] ) row.list <- lapply( 1:np1, by.row.functions, c, p.p.products ) return( row.list ) } ### ### function to compute the integral of the polynomials by column ### by.column.integrals <- function( c, p.p.functions ) { ### ### function to compute the integral of the polynomials by row ### by.row.integrals <- function( r, c, p.p.functions ) { row.column.integral <- ghermite.h.quadrature( p.p.functions[[r]][[c]], order.np1.rule, mu ) return( row.column.integral ) } np1 <- length( p.p.functions[[1]] ) row.vector <- sapply( 1:np1, by.row.integrals, c, p.p.functions ) return( row.vector ) } ### ### construct a list of the generalized Hermite H orthogonal polynomials ### p.list <- ghermite.h.polynomials( n, mu ) ### ### construct the two dimensional list of pair products ### of polynomials ### p.p.products <- lapply( 1:np1, by.column.products, p.list ) ### ### compute the two dimensional list of functions ### corresponding to the polynomial products in ### the two dimensional list p.p.products ### p.p.functions <- lapply( 1:np1, by.column.functions, p.p.products ) ### ### get the rule table for the order np1 polynomial ### rules <- ghermite.h.quadrature.rules( np1, mu ) order.np1.rule <- rules[[np1]] ### ### construct the square symmetric matrix containing ### the definite integrals over the default limits ### corresponding to the two dimensional list of ### polynomial functions ### p.p.integrals <- sapply( 1:np1, by.column.integrals, p.p.functions ) ### ### construct the diagonal matrix with the inner products ### of the orthogonal polynomials on the diagonal ### p.p.inner.products <- diag( ghermite.h.inner.products( n, mu ) ) print( "Integral of cross products for the orthogonal polynomials " ) print( apply( p.p.integrals, 2, round, digits=5 ) ) print( apply( p.p.inner.products, 2, round, digits=5 ) ) ### ### construct a list of the Hermite H orthonormal polynomials ### p.list <- ghermite.h.polynomials( n, mu, TRUE ) ### ### construct the two dimensional list of pair products ### of polynomials ### p.p.products <- lapply( 1:np1, by.column.products, p.list ) ### ### compute the two dimensional list of functions ### corresponding to the polynomial products in ### the two dimensional list p.p.products ### p.p.functions <- lapply( 1:np1, by.column.functions, p.p.products ) ### ### get the rule table for the order np1 polynomial ### rules <- ghermite.h.quadrature.rules( np1, mu, TRUE ) order.np1.rule <- rules[[np1]] ### ### construct the square symmetric matrix containing ### the definite integrals over the default limits ### corresponding to the two dimensional list of ### polynomial functions ### p.p.integrals <- sapply( 1:np1, by.column.integrals, p.p.functions ) ### ### display the matrix of integrals ### print( "Integral of cross products for the orthonormal polynomials " ) print(apply( p.p.integrals, 2, round, digits=5 ) ) 

### Example output

Loading required package: polynom
[1] "Integral of cross products for the orthogonal polynomials "
[,1]    [,2]     [,3]     [,4]     [,5]     [,6]     [,7]
[1,] 0.88623 0.00000  0.00000   0.0000    0.000     0.00      0.0
[2,] 0.00000 5.31736  0.00000   0.0000    0.000     0.00      0.0
[3,] 0.00000 0.00000 21.26945   0.0000    0.000     0.00      0.0
[4,] 0.00000 0.00000  0.00000 212.6945    0.000     0.00      0.0
[5,] 0.00000 0.00000  0.00000   0.0000 1701.556     0.00      0.0
[6,] 0.00000 0.00000  0.00000   0.0000    0.000 23821.78      0.0
[7,] 0.00000 0.00000  0.00000   0.0000    0.000     0.00 285861.4
[8,] 0.00000 0.00000  0.00000   0.0000    0.000     0.00      0.0
[9,] 0.00000 0.00000  0.00000   0.0000    0.000     0.00      0.0
[10,] 0.00000 0.00000  0.00000   0.0000    0.000     0.00      0.0
[11,] 0.00000 0.00000  0.00000   0.0000    0.000     0.00      0.0
[,8]          [,9]        [,10]         [,11]
[1,]       0.00000  0.000000e+00 0.000000e+00  0.000000e+00
[2,]       0.00000  0.000000e+00 0.000000e+00  0.000000e+00
[3,]       0.00000  0.000000e+00 0.000000e+00  0.000000e+00
[4,]       0.00000  0.000000e+00 0.000000e+00  0.000000e+00
[5,]       0.00000  0.000000e+00 0.000000e+00  0.000000e+00
[6,]       0.00000  0.000000e+00 0.000000e+00  0.000000e+00
[7,]       0.00000  0.000000e+00 0.000000e+00  0.000000e+00
[8,] 5145504.42733  0.000000e+00 0.000000e+00 -2.000000e-05
[9,]       0.00000  8.232807e+07 0.000000e+00 -8.000000e-05
[10,]       0.00000  0.000000e+00 1.811218e+09  1.010000e-03
[11,]      -0.00002 -8.000000e-05 1.010000e-03  3.622435e+10
[,1]    [,2]     [,3]     [,4]     [,5]     [,6]     [,7]    [,8]
[1,] 0.88623 0.00000  0.00000   0.0000    0.000     0.00      0.0       0
[2,] 0.00000 5.31736  0.00000   0.0000    0.000     0.00      0.0       0
[3,] 0.00000 0.00000 21.26945   0.0000    0.000     0.00      0.0       0
[4,] 0.00000 0.00000  0.00000 212.6945    0.000     0.00      0.0       0
[5,] 0.00000 0.00000  0.00000   0.0000 1701.556     0.00      0.0       0
[6,] 0.00000 0.00000  0.00000   0.0000    0.000 23821.78      0.0       0
[7,] 0.00000 0.00000  0.00000   0.0000    0.000     0.00 285861.4       0
[8,] 0.00000 0.00000  0.00000   0.0000    0.000     0.00      0.0 5145504
[9,] 0.00000 0.00000  0.00000   0.0000    0.000     0.00      0.0       0
[10,] 0.00000 0.00000  0.00000   0.0000    0.000     0.00      0.0       0
[11,] 0.00000 0.00000  0.00000   0.0000    0.000     0.00      0.0       0
[,9]      [,10]       [,11]
[1,]        0          0           0
[2,]        0          0           0
[3,]        0          0           0
[4,]        0          0           0
[5,]        0          0           0
[6,]        0          0           0
[7,]        0          0           0
[8,]        0          0           0
[9,] 82328071          0           0
[10,]        0 1811217558           0
[11,]        0          0 36224351168
[1] "Integral of cross products for the orthonormal polynomials "
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
[1,]    1    0    0    0    0    0    0    0    0     0     0
[2,]    0    1    0    0    0    0    0    0    0     0     0
[3,]    0    0    1    0    0    0    0    0    0     0     0
[4,]    0    0    0    1    0    0    0    0    0     0     0
[5,]    0    0    0    0    1    0    0    0    0     0     0
[6,]    0    0    0    0    0    1    0    0    0     0     0
[7,]    0    0    0    0    0    0    1    0    0     0     0
[8,]    0    0    0    0    0    0    0    1    0     0     0
[9,]    0    0    0    0    0    0    0    0    1     0     0
[10,]    0    0    0    0    0    0    0    0    0     1     0
[11,]    0    0    0    0    0    0    0    0    0     0     1


gaussquad documentation built on May 30, 2017, 8:06 a.m.