program main

!*****************************************************************************80
!
!! MAIN is the main program for QUAD.
!
!  Discussion:
!
!    QUAD is a sketch of a FORTRAN90 program to estimate the integral of
!    a function F(X) between A and B.
!
!    This program is only a sketch.  Your job is to finish it,
!    debug it, and get it to run.
!
!    Experiment to get a value of N that gives you a 'reasonable' accuracy.
!
!    NOW modify the program so that it will run under OpenMP.  Replace
!    CPU time by wall clock time.
!
!    Show that your program STILL gets the right answer.
!
!    Then determine the speedup when you run the OpenMP version on
!    1, 2, 4, 8 processors.
!
!
!  For an OpenMP version, you'll need 
!
!     use omp_lib

  integer n
  parameter ( n = 1000000 )

  double precision a
  double precision b
  double precision error
  double precision exact
  integer flops
  integer i
  double precision mflops
  double precision, parameter :: pi = 3.141592653589793D+00
  double precision total
  double precision wtime
  double precision wtime1
  double precision wtime2
  double precision x(n)

  a =  0.0
  b = 10.0
  exact = 0.49936338107645674464D+00

  write ( *, * ) ' '
  write ( *, * ) 'QUAD:'
  write ( *, * ) '  FORTRAN90 version'
  write ( *, * ) ' '
  write ( *, * ) '  Estimate the integral of f(x) from A to B.'
  write ( *, * ) '  f(x) = 50 / ( pi * ( 2500 * x * x + 1 ) ).'
  write ( *, * ) '  A = ', a
  write ( *, * ) '  B = ', b
  write ( *, * ) '  Exact integral is 0.49936338107645674464...'
!
!  Step 1:
!  Load the array X with evenly spaced values between A and B.
!
  do i = 1, n
    x(i) =  ( dble ( n - i     ) * a   &
            + dble (     i - 1 ) * b ) &
            / dble ( n     - 1 )
  end do
!
!  Step 2:
!
!  Use a DO loop to add to TOTAL the value of the function at each X.
!
!  Get your estimate by multiplying TOTAL by ( B - A ) and dividing by N.
!
!  Set FLOPS, the number of floating point operations, by multiplying N
!  by the number of floating point operations in the TOTAL update statement.
!
!  Add the OpenMP include file,
!
!  Insert an OpenMP directive before this loop of the form
!
!    !$omp parallel do private ( ... ) shared ( ... ) reduction ( + : ... )
!
!  Insert an OpenMP end directive of the form
!
!    !$omp end parallel do
!
  wtime1 = omp_get_wtime ( )

  total = 0.0D+00

  MISSING OPENMP DIRECTIVE
    MISSING DO LOOP BEGINNING
      MISSING UPDATE OF TOTAL
    MISSING END OF DO LOOP
  MISSING OPENMP END DIRECTIVE

  wtime2 = omp_get_wtime ( )

  total = ( b - a ) * total / dble ( n )

  flops = n * ?
!
!  Step 3:
!  Print quadrature estimate, Error, CPU time, MFLOPS rate
!
  error = abs ( total - exact )
  wtime = wtime2 - wtime1
  mflops = dble ( flops ) / 1000000.0 / wtime

  write ( *, * ) ' '
  write ( *, * ) '  Estimate = ', total
  write ( *, * ) '  Error    = ', error
  write ( *, * ) '  W time   = ', wtime
  write ( *, * ) '  FLOPS    = ', flops
  write ( *, * ) '  MFLOPS   = ', mflops

  stop
end