An Example

In this appendix we describe a simple 5-dimensional example. The integrands of interest are 5 expectations, $E\{x_i\}$, $i = 1, 2, \dots, 5$, with

$$E\{x_i\} = \frac{\int_0^1\int_0^{1-x_1}\cdots\int_0^{1-x_1-\... ...ots-x_4} \exp{(-\sum_{j=1}^5(jx_j)^2)}dx_5dx_4dx_3dx_2dx_1}.$$

Here is a Fortran program for this example.

MODULE INTEGRAND
  USE PRECISION_MODEL
  IMPLICIT NONE
  PRIVATE
  PUBLIC :: F
CONTAINS
  FUNCTION F( L, X ) RESULT(FUN)
    INTEGER,                       INTENT(IN) :: L
    REAL(kind=stnd), DIMENSION(:), INTENT(IN) :: X
    REAL(kind=stnd), DIMENSION(0:L-1)         :: FUN
    !
    REAL(kind=stnd)    :: S
    INTEGER            :: I
    INTEGER, PARAMETER :: N = 5
    S = 0
    DO I = 1, N
       S = S + ( I*X(I) )**2
    END DO
    FUN(0) = EXP(-S)
    FUN(1:N) = X(1:N)*FUN(0)
  END FUNCTION F
!
END MODULE INTEGRAND
!
PROGRAM Simplex_Example

USE Precision_Model
USE CUI                                 ! Cubpack User Interface
USE INTEGRAND

INTEGER, PARAMETER         :: n = 5,   & ! the dimension
                              nf = n + 1 ! number of integrand functions

INTEGER                             :: NEval, i, Inform
REAL(kind=stnd), DIMENSION(n,0:n,1) :: Simplex
REAL(kind=stnd), DIMENSION(n,0:n,2) :: TwoSimplices
REAL(kind=stnd), DIMENSION(0:n)     :: Value, AbsErr

Inform = -1                      ! Soft noisy errors            

Simplex = 0                      ! Build the unit simplex 
DO i = 1, n
   Simplex(i,i,1) = 1
END DO
TwoSimplices(:,:,1:1) = Simplex    ! Split it into 2 parts.
TwoSimplices(1,0,1) = 0.5_stnd     !  Therefore one has to change one
TwoSimplices(:,:,2:2) = Simplex    !    coordinate of the unit simplex.
TwoSimplices(1,1,2) = 0.5_stnd

!   The two subregions | are simplices    | |
!                      v                  v v
CALL CUBATR( n, nf, F, 2, TwoSimplices, (/1,1/), Value, AbsErr,    &
             Inform, NEval, EpsRel=5.0e-4_stnd, MaxPts=200000 )

Print "(""Expected values are ""    /5F12.8)",  Value(1:n)/Value(0) 
Print "(""with estimated errors < ""/5F12.8)", AbsErr(1:n)/Value(0) 
Print *,"The number of integrand evaluations used was ", NEval 

Inform = -1                      ! Soft noisy errors            
!    The single region | is a simplex    |
!                      v                 v
CALL CUBATR( n, nf, F, 1, Simplex,     (/1/), Value, AbsErr,      &
             Inform, NEval, EpsRel=5.0e-4_stnd, MaxPts=200000 )

Print "(""Expected values are ""    /5F12.8)",  Value(1:n)/Value(0) 
Print "(""with estimated errors < ""/5F12.8)", AbsErr(1:n)/Value(0) 
Print *,"The number of integrand evaluations used was ", NEval 

END PROGRAM Simplex_Example

In this test, CUBATR is called twice to illustrate its use with more than one subregion. On the first call, the original input simplex is split into two halves and the output gives the sum of results for both halves. For the second call, only the original simplex is used. Default values were used for all of the optional parameters except INFORM, EspRel and MaxPts. The following output results were produced.

Expected values are 
  0.22419120  0.17813751  0.14013325  0.11372875  0.09523524
with estimated errors < 
  0.00009069  0.00006152  0.00006977  0.00005160  0.00004529
 The number of integrand evaluations used was  178627
Expected values are 
  0.22419087  0.17813878  0.14013454  0.11372938  0.09523551
with estimated errors < 
  0.00006262  0.00006472  0.00006016  0.00005668  0.00004474
 The number of integrand evaluations used was  76735

Observe that the second call reached the requested accuracy much faster, which illustrates the effectiveness of the adaptive strategy used by CUBATR. The first call forces an initial subdivision which is not as efficient as the subdivision chosen by by CUBATR for the second call.