diff options
| -rw-r--r-- | src/quadratures.f90 | 26 |
1 files changed, 19 insertions, 7 deletions
diff --git a/src/quadratures.f90 b/src/quadratures.f90 index 2eac963..43e6ac6 100644 --- a/src/quadratures.f90 +++ b/src/quadratures.f90 @@ -45,14 +45,10 @@ CONTAINS real(kind=8), intent(in) :: iend procedure(funint) :: fun integer(kind=4), intent(in) :: p - - real(kind=8) :: val, subinterval_width, qbeg, qend - real(kind=8), allocatable :: partval[:] + real(kind=8) :: val procedure(quadrature), pointer :: quad - integer(kind=8) :: min_i, max_i, i, subintervals_per_thread - integer(kind=4) :: im - + SELECT CASE (p) CASE (:-1) STOP "negative interpolationg polynomial order passed" @@ -66,6 +62,22 @@ CONTAINS STOP "Newton-Cotes quadratures only implemented for order < 3" END SELECT + val = integrate_generic(ibeg, iend, fun, quad) + END FUNCTION newton_cotes + + FUNCTION integrate_generic(ibeg, iend, fun, quad) result(val) + IMPLICIT none + real(kind=8), intent(in) :: ibeg + real(kind=8), intent(in) :: iend + procedure(funint) :: fun + procedure(quadrature) :: quad + + real(kind=8) :: val, subinterval_width, qbeg, qend + real(kind=8), allocatable :: partval[:] + + integer(kind=8) :: min_i, max_i, i, subintervals_per_thread + integer(kind=4) :: im + if (this_image() == 1) allocate(partval[*]) subintervals_per_thread = & @@ -98,7 +110,7 @@ CONTAINS val = partval[1] - END FUNCTION newton_cotes + END FUNCTION integrate_generic FUNCTION rectangle(qbeg, qend, fun) result(val) IMPLICIT none |