stdlib_intrinsics moduleThe stdlib_intrinsics module provides replacements for some of the well known intrinsic functions found in Fortran compilers for which either a faster and/or more accurate implementation is found which has also proven of interest to the Fortran community.
stdlib_sum functionThe stdlib_sum function can replace the intrinsic sum for real, complex or integer arrays. It follows a chunked implementation which maximizes vectorization potential as well as reducing the round-off error. This procedure is recommended when summing large (e..g, >2**10 elements) arrays, for repetitive summation of smaller arrays consider the classical sum.
res = stdlib_sum (x [,mask] )
res = stdlib_sum (x, dim [,mask] )
Experimental
Pure function.
x: N-D array of either real, complex or integer type. This argument is intent(in).
dim (optional): scalar of type integer with a value in the range from 1 to n, where n equals the rank of x.
mask (optional): N-D array of logical values, with the same shape as x. This argument is intent(in).
If dim is absent, the output is a scalar of the same type and kind as to that of x. Otherwise, an array of rank n-1, where n equals the rank of x, and a shape similar to that of x with dimension dim dropped is returned.
stdlib_sum_kahan functionThe stdlib_sum_kahan function can replace the intrinsic sum for real or complex arrays. It follows a chunked implementation which maximizes vectorization potential complemented by an elemental kernel based on the kahan summation strategy to reduce the round-off error:
elemental subroutine kahan_kernel_<kind>(a,s,c)
type(<kind>), intent(in) :: a
type(<kind>), intent(inout) :: s
type(<kind>), intent(inout) :: c
type(<kind>) :: t, y
y = a - c
t = s + y
c = (t - s) - y
s = t
end subroutine
res = stdlib_sum_kahan (x [,mask] )
res = stdlib_sum_kahan (x, dim [,mask] )
Experimental
Pure function.
x: 1D array of either real or complex type. This argument is intent(in).
dim (optional): scalar of type integer with a value in the range from 1 to n, where n equals the rank of x.
mask (optional): N-D array of logical values, with the same shape as x. This argument is intent(in).
If dim is absent, the output is a scalar of the same type and kind as to that of x. Otherwise, an array of rank n-1, where n equals the rank of x, and a shape similar to that of x with dimension dim dropped is returned.
program example_sum
use stdlib_kinds, only: sp
use stdlib_intrinsics, only: stdlib_sum, stdlib_sum_kahan
implicit none
real(sp), allocatable :: x(:)
real(sp) :: total_sum(3)
allocate( x(1000) )
call random_number(x)
total_sum(1) = sum(x) !> compiler intrinsic
total_sum(2) = stdlib_sum(x) !> chunked summation
total_sum(3) = stdlib_sum_kahan(x)!> chunked kahan summation
print *, total_sum(1:3)
end program example_sum
stdlib_dot_product functionThe stdlib_dot_product function can replace the intrinsic dot_product for 1D real, complex or integer arrays. It follows a chunked implementation which maximizes vectorization potential as well as reducing the round-off error. This procedure is recommended when crunching large arrays, for repetitive products of smaller arrays consider the classical dot_product.
res = stdlib_dot_product (x, y)
Experimental
Pure function.
x: 1D array of either real, complex or integer type. This argument is intent(in).
y: 1D array of the same type and kind as x. This argument is intent(in).
The output is a scalar of type and kind same as to that of x and y.
stdlib_dot_product_kahan functionThe stdlib_dot_product_kahan function can replace the intrinsic dot_product for 1D real or complex arrays. It follows a chunked implementation which maximizes vectorization potential, complemented by the same elemental kernel based on the kahan summation used for stdlib_sum to reduce the round-off error.
res = stdlib_dot_product_kahan (x, y)
Experimental
Pure function.
x: 1D array of either real or complex type. This argument is intent(in).
y: 1D array of the same type and kind as x. This argument is intent(in).
The output is a scalar of the same type and kind as to that of x and y.
program example_dot_product
use stdlib_kinds, only: sp
use stdlib_intrinsics, only: stdlib_dot_product, stdlib_dot_product_kahan
implicit none
real(sp), allocatable :: x(:), y(:)
real(sp) :: total_prod(3)
allocate( x(1000), y(1000) )
call random_number(x)
call random_number(y)
total_prod(1) = dot_product(x,y) !> compiler intrinsic
total_prod(2) = stdlib_dot_product(x,y) !> chunked summation over inner product
total_prod(3) = stdlib_dot_product_kahan(x,y) !> chunked kahan summation over inner product
print *, total_prod(1:3)
end program example_dot_product
stdlib_matmul functionThe extension of the intrinsic function matmul to handle more than 2 and less than or equal to 5 matrices, with error handling using linalg_state_type.
The optimal parenthesization to minimize the number of scalar multiplications is done using the Algorithm as outlined in Cormen, "Introduction to Algorithms", 4ed, ch-14, section-2.
The actual matrix multiplication is performed using the gemm interfaces.
It supports only real and complex matrices.
res = stdlib_matmul (m1, m2, m3, m4, m5, err)
Experimental
Function.
m1, m2: 2D arrays of the same kind and type. intent(in) arguments.
m3,m4,m5: 2D arrays of the same kind and type as the other matrices. intent(in), optional arguments.
err: type(linalg_state_type), intent(out), optional argument. Can be used for elegant error handling. It is assigned LINALG_VALUE_ERROR
in case the matrices are not of compatible sizes.
The output is a matrix of the appropriate size.
program example_matmul
use stdlib_intrinsics, only: stdlib_matmul
complex :: x(2, 2), y(2, 2), z(2, 2)
x = reshape([(0, 0), (1, 0), (1, 0), (0, 0)], [2, 2]) ! pauli x-matrix
y = reshape([(0, 0), (0, 1), (0, -1), (0, 0)], [2, 2]) ! pauli y-matrix
z = reshape([(1, 0), (0, 0), (0, 0), (-1, 0)], [2, 2]) ! pauli z-matrix
print *, stdlib_matmul(x, y) ! should be iota*z
print *, stdlib_matmul(y, z, x) ! should be iota*identity
print *, stdlib_matmul(x, x, z, y) ! should be -iota*x
print *, stdlib_matmul(x, x, z, y, y) ! should be z
end program example_matmul