67 real(kind=
rp),
public,
parameter ::
neko_eps = epsilon(1.0_rp)
73 real(kind=
rp),
public,
parameter ::
pi = 4._rp*atan(1._rp)
101 vcross,
vdot2,
vdot3,
vlsc3,
vlsc2,
add2,
add3,
add4,
sub2,
sub3, &
111 real(kind=
sp),
intent(in) :: x
112 real(kind=
sp),
intent(in) :: y
121 real(kind=
dp),
intent(in) :: x
122 real(kind=
dp),
intent(in) :: y
131 real(kind=
qp),
intent(in) :: x
132 real(kind=
qp),
intent(in) :: y
141 real(kind=
sp),
intent(in) :: x
142 real(kind=
sp),
intent(in) :: y
143 real(kind=
sp),
intent(in),
optional :: eps
145 if (
present(eps))
then
146 srelcmp = abs(x - y) .le. eps*abs(y)
155 real(kind=
dp),
intent(in) :: x
156 real(kind=
dp),
intent(in) :: y
157 real(kind=
dp),
intent(in),
optional :: eps
159 if (
present(eps))
then
160 drelcmp = abs(x - y) .le. eps*abs(y)
170 real(kind=
qp),
intent(in) :: x
171 real(kind=
qp),
intent(in) :: y
172 real(kind=
qp),
intent(in),
optional :: eps
174 if (
present(eps))
then
175 qrelcmp = abs(x - y)/abs(y) .lt. eps
184 integer,
intent(in) :: n
185 real(kind=
rp),
dimension(n),
intent(inout) :: a
195 integer,
intent(in) :: n
196 integer,
dimension(n),
intent(inout) :: a
206 integer,
intent(in) :: m, n, e
207 real(kind=
rp),
intent(inout) :: a(m,n)
217 integer,
intent(in) :: n
218 real(kind=
rp),
dimension(n),
intent(inout) :: a
228 integer,
intent(in) :: n
229 real(kind=
rp),
dimension(n),
intent(in) :: b
230 real(kind=
rp),
dimension(n),
intent(inout) :: a
247 integer,
intent(in) :: n, m
248 real(kind=
rp),
dimension(n),
intent(in) :: b
249 real(kind=
rp),
dimension(n),
intent(inout) :: a
250 integer,
dimension(0:m) :: mask
263 integer,
intent(in) :: size, mask_size
264 real(kind=
rp),
dimension(size),
intent(inout) :: a
265 real(kind=
rp),
intent(in) :: c
266 integer,
dimension(mask_size),
intent(in) :: mask
277 integer,
intent(in) :: n
278 real(kind=
rp),
dimension(n),
intent(inout) :: a
279 real(kind=
rp),
intent(in) :: c
289 integer,
intent(in) :: n
290 real(kind=
rp),
dimension(n),
intent(inout) :: a
291 real(kind=
rp),
intent(in) :: s
301 integer,
intent(in) :: n
302 real(kind=
rp),
dimension(n),
intent(inout) :: a
303 real(kind=
rp),
dimension(n),
intent(in) :: b
304 real(kind=
rp),
intent(in) :: s
314 integer,
intent(in) :: n
315 real(kind=
rp),
dimension(n),
intent(inout) :: a
316 real(kind=
rp),
intent(in) :: c
326 integer,
intent(in) :: n
327 real(kind=
rp),
dimension(n) :: a
334 call mpi_allreduce(tmp,
glsum, 1, &
341 integer,
intent(in) :: n
342 real(kind=
rp),
dimension(n) :: a
349 call mpi_allreduce(tmp,
glmax, 1, &
355 integer,
intent(in) :: n
356 integer,
dimension(n) :: a
363 call mpi_allreduce(tmp,
glimax, 1, &
369 integer,
intent(in) :: n
370 real(kind=
rp),
dimension(n) :: a
377 call mpi_allreduce(tmp,
glmin, 1, &
383 integer,
intent(in) :: n
384 integer,
dimension(n) :: a
391 call mpi_allreduce(tmp,
glimin, 1, &
400 integer,
intent(in) :: n
401 real(kind=
rp),
dimension(n),
intent(inout) :: a
413 real(kind=
rp),
intent(in) :: vec(n)
414 real(kind=
rp) :: tmax
417 tmax =
max(tmax, vec(i))
423 integer,
intent(in) :: n
424 real(kind=
rp),
intent(in) :: vec(n)
425 real(kind=
rp) :: tmin
429 tmin = min(tmin, vec(i))
435 integer,
intent(in) :: n
436 real(kind=
rp),
dimension(n),
intent(inout) :: a
447 integer,
intent(in) :: n
448 real(kind=
rp),
dimension(n),
intent(inout) :: a
449 real(kind=
rp),
dimension(n),
intent(in) :: b,c
460 integer,
intent(in) :: n
461 real(kind=
rp),
dimension(n),
intent(inout) :: a
462 real(kind=
rp),
dimension(n),
intent(in) :: b
473 subroutine vcross(u1, u2, u3, v1, v2, v3, w1, w2, w3, n)
474 integer,
intent(in) :: n
475 real(kind=
rp),
dimension(n),
intent(in) :: v1, v2, v3
476 real(kind=
rp),
dimension(n),
intent(in) :: w1, w2, w3
477 real(kind=
rp),
dimension(n),
intent(out) :: u1, u2, u3
481 u1(i) = v2(i)*w3(i) - v3(i)*w2(i)
482 u2(i) = v3(i)*w1(i) - v1(i)*w3(i)
483 u3(i) = v1(i)*w2(i) - v2(i)*w1(i)
490 subroutine vdot2(dot, u1, u2, v1, v2, n)
491 integer,
intent(in) :: n
492 real(kind=
rp),
dimension(n),
intent(in) :: u1, u2
493 real(kind=
rp),
dimension(n),
intent(in) :: v1, v2
494 real(kind=
rp),
dimension(n),
intent(out) :: dot
497 dot(i) = u1(i)*v1(i) + u2(i)*v2(i)
504 subroutine vdot3(dot, u1, u2, u3, v1, v2, v3, n)
505 integer,
intent(in) :: n
506 real(kind=
rp),
dimension(n),
intent(in) :: u1, u2, u3
507 real(kind=
rp),
dimension(n),
intent(in) :: v1, v2, v3
508 real(kind=
rp),
dimension(n),
intent(out) :: dot
512 dot(i) = u1(i)*v1(i) + u2(i)*v2(i) + u3(i)*v3(i)
518 function vlsc3(u, v, w, n)
result(s)
519 integer,
intent(in) :: n
520 real(kind=
rp),
dimension(n),
intent(in) :: u, v, w
526 s = s + u(i)*v(i)*w(i)
533 integer,
intent(in) :: n
534 real(kind=
rp),
dimension(n),
intent(in) :: u, v
547 integer,
intent(in) :: n
548 real(kind=
rp),
dimension(n),
intent(inout) :: a
549 real(kind=
rp),
dimension(n),
intent(in) :: b
560 integer,
intent(in) :: n
561 real(kind=
rp),
dimension(n),
intent(inout) :: a
562 real(kind=
rp),
dimension(n),
intent(in) :: b
563 real(kind=
rp),
dimension(n),
intent(in) :: c
574 integer,
intent(in) :: n
575 real(kind=
rp),
dimension(n),
intent(inout) :: d
576 real(kind=
rp),
dimension(n),
intent(inout) :: c
577 real(kind=
rp),
dimension(n),
intent(inout) :: b
578 real(kind=
rp),
dimension(n),
intent(out) :: a
582 a(i) = b(i) + c(i) + d(i)
589 integer,
intent(in) :: n
590 real(kind=
rp),
dimension(n),
intent(inout) :: a
591 real(kind=
rp),
dimension(n),
intent(inout) :: b
602 integer,
intent(in) :: n
603 real(kind=
rp),
dimension(n),
intent(inout) :: a
604 real(kind=
rp),
dimension(n),
intent(in) :: b
605 real(kind=
rp),
dimension(n),
intent(in) :: c
618 integer,
intent(in) :: n
619 real(kind=
rp),
dimension(n),
intent(inout) :: a
620 real(kind=
rp),
dimension(n),
intent(inout) :: b
621 real(kind=
rp),
intent(in) :: c1
625 a(i) = c1 * a(i) + b(i)
633 integer,
intent(in) :: n
634 real(kind=
rp),
dimension(n),
intent(inout) :: a
635 real(kind=
rp),
dimension(n),
intent(inout) :: b
636 real(kind=
rp),
intent(in) :: c1
640 a(i) = a(i) + c1 * b(i)
647 integer,
intent(in) :: n
648 real(kind=
rp),
dimension(n),
intent(inout) :: a
649 real(kind=
rp),
dimension(n),
intent(in) :: b
650 real(kind=
rp),
intent(in) :: c1
654 a(i) = a(i) + c1 * ( b(i) * b(i) )
661 integer,
intent(in) :: n
662 real(kind=
rp),
dimension(n),
intent(inout) :: a
663 real(kind=
rp),
dimension(n),
intent(in) :: b
664 real(kind=
rp),
intent(in) :: c
675 integer,
intent(in) :: n
676 real(kind=
rp),
dimension(n),
intent(inout) :: a
677 real(kind=
rp),
dimension(n),
intent(in) :: b
689 integer,
intent(in) :: n
690 real(kind=
rp),
dimension(n),
intent(inout) :: a
691 real(kind=
rp),
dimension(n),
intent(in) :: b
702 integer,
intent(in) :: n
703 real(kind=
rp),
dimension(n),
intent(inout) :: a
704 real(kind=
rp),
dimension(n),
intent(in) :: b
705 real(kind=
rp),
dimension(n),
intent(in) :: c
716 integer,
intent(in) :: n
717 real(kind=
rp),
dimension(n),
intent(inout) :: a
718 real(kind=
rp),
dimension(n),
intent(in) :: b
719 real(kind=
rp),
dimension(n),
intent(in) :: c
723 a(i) = a(i) - b(i) * c(i)
730 integer,
intent(in) :: n
731 real(kind=
rp),
dimension(n),
intent(inout) :: a
732 real(kind=
rp),
dimension(n),
intent(in) :: b
733 real(kind=
rp),
dimension(n),
intent(in) :: c
734 real(kind=
rp),
intent(in) :: c1, c2
738 a(i) = c1 * b(i) + c2 * c(i)
746 integer,
intent(in) :: n
747 real(kind=
rp),
dimension(n),
intent(inout) :: a
748 real(kind=
rp),
dimension(n),
intent(in) :: b
749 real(kind=
rp),
dimension(n),
intent(in) :: c
750 real(kind=
rp),
dimension(n),
intent(in) :: d
754 a(i) = a(i) - b(i) * c(i) * d(i)
761 integer,
intent(in) :: n
762 real(kind=
rp),
dimension(n),
intent(inout) :: a
763 real(kind=
rp),
dimension(n),
intent(in) :: b
764 real(kind=
rp),
dimension(n),
intent(in) :: c
768 a(i) = a(i) + b(i) * c(i)
775 integer,
intent(in) :: n
776 real(kind=
rp),
dimension(n),
intent(inout) :: a
777 real(kind=
rp),
dimension(n),
intent(in) :: b
778 real(kind=
rp),
dimension(n),
intent(in) :: c
779 real(kind=
rp),
dimension(n),
intent(in) :: d
783 a(i) = a(i) + b(i) * c(i) * d(i)
790 integer,
intent(in) :: n
791 real(kind=
rp),
dimension(n),
intent(inout) :: a
792 real(kind=
rp),
dimension(n),
intent(in) :: b
793 real(kind=
rp),
dimension(n),
intent(in) :: c
794 real(kind=
rp),
dimension(n),
intent(in) :: d
795 real(kind=
rp),
dimension(n),
intent(in) :: e
799 a(i) = b(i)*c(i)-d(i)*e(i)
806 integer,
intent(in) :: n
807 real(kind=
rp),
dimension(n),
intent(inout) :: a
808 real(kind=
rp),
dimension(n),
intent(in) :: b
809 real(kind=
rp),
dimension(n),
intent(in) :: c
810 real(kind=
rp),
intent(in) :: c1, c2
814 a(i) = b(i) + c1*(a(i)-c2*c(i))
821 integer,
intent(in) :: n
822 real(kind=
rp),
dimension(n),
intent(inout) :: a
823 real(kind=
rp),
dimension(n),
intent(in) :: b
824 real(kind=
rp),
dimension(n),
intent(in) :: c
825 real(kind=
rp),
intent(in) :: c1, c2
829 a(i) = a(i) + c1*b(i)+c2*c(i)
836 integer,
intent(in) :: n
837 real(kind=
rp),
dimension(n),
intent(in) :: a
838 real(kind=
rp),
dimension(n),
intent(in) :: b
844 tmp = tmp + a(i) * b(i)
847 call mpi_allreduce(tmp,
glsc2, 1, &
854 integer,
intent(in) :: n
855 real(kind=
rp),
dimension(n),
intent(in) :: a
856 real(kind=
rp),
dimension(n),
intent(in) :: b
857 real(kind=
rp),
dimension(n),
intent(in) :: c
863 tmp = tmp + a(i) * b(i) * c(i)
866 call mpi_allreduce(tmp,
glsc3, 1, &
871 integer,
intent(in) :: n
872 real(kind=
rp),
dimension(n),
intent(in) :: a
873 real(kind=
rp),
dimension(n),
intent(in) :: b
874 real(kind=
rp),
dimension(n),
intent(in) :: c
875 real(kind=
rp),
dimension(n),
intent(in) :: d
881 tmp = tmp + a(i) * b(i) * c(i) * d(i)
884 call mpi_allreduce(tmp,
glsc4, 1, &
895 integer,
intent(in) :: n
896 real(kind=
rp),
intent(inout) :: a(n)
897 integer,
intent(out) :: ind(n)
899 integer :: j, ir, i, ii, l
931 if ( a(j) .lt. a(j+1) ) j = j + 1
933 if (aa .lt. a(j))
then
953 integer,
intent(in) :: n
954 integer(i4),
intent(inout) :: a(n)
955 integer,
intent(out) :: ind(n)
957 integer :: j, ir, i, ii, l
988 if ( a(j) .lt. a(j + 1) ) j = j + 1
990 if (aa .lt. a(j))
then
1009 integer,
intent(in) :: n
1010 real(kind=
rp),
intent(inout) :: b(n)
1011 integer,
intent(in) :: ind(n)
1012 real(kind=
rp) :: temp(n)
1029 integer,
intent(in) :: n
1030 integer(i4),
intent(inout) :: b(n)
1031 integer,
intent(in) :: ind(n)
1032 integer(i4) :: temp(n)
1049 integer,
intent(in) :: n
1050 real(kind=
rp),
intent(inout) :: b(n)
1051 integer,
intent(in) :: ind(n)
1052 real(kind=
rp) :: temp(n)
1069 integer,
intent(in) :: n
1070 integer(i4),
intent(inout) :: b(n)
1071 integer,
intent(in) :: ind(n)
1072 integer(i4) :: temp(n)
1089 integer,
intent(in) :: n
1090 real(kind=
rp),
intent(inout) :: b(n)
1091 integer,
intent(inout) :: ind(n)
1092 real(kind=
rp) :: temp(n)
1093 integer :: tempind(n)
1099 tempind(jj) = ind(i)
1112 integer,
intent(in) :: n
1113 integer(i4),
intent(inout) :: b(n)
1114 integer,
intent(inout) :: ind(n)
1115 integer(i4) :: temp(n)
1116 integer :: tempind(n)
1122 tempind(jj) = ind(i)
type(mpi_comm) neko_comm
MPI communicator.
type(mpi_datatype) mpi_real_precision
MPI type for working precision of REAL types.
subroutine, public cmult(a, c, n)
Multiplication by constant c .
subroutine, public cmult2(a, b, c, n)
Multiplication by constant c .
subroutine, public row_zero(a, m, n, e)
Sets row e to 0 in matrix a.
subroutine, public invcol2(a, b, n)
Vector division .
real(kind=rp) function, public vlsc2(u, v, n)
Compute multiplication sum .
real(kind=rp), parameter, public pi
real(kind=rp) function, public glsc3(a, b, c, n)
Weighted inner product .
pure logical function qabscmp(x, y)
Return double precision absolute comparison .
subroutine, public ascol5(a, b, c, d, e, n)
Returns .
subroutine, public invers2(a, b, n)
Compute inverted vector .
subroutine, public cadd2(a, b, s, n)
Add a scalar to vector .
subroutine, public cadd(a, s, n)
Add a scalar to vector .
subroutine reorddp(b, ind, n)
reorder double precision array - inverse of swap
subroutine, public addsqr2s2(a, b, c1, n)
Returns .
real(kind=rp) function, public glsc4(a, b, c, d, n)
subroutine swapdp(b, ind, n)
sort double precision array acording to ind vector
subroutine flipvi4(b, ind, n)
Flip single integer vector b and ind.
subroutine, public add2s1(a, b, c1, n)
Vector addition with scalar multiplication (multiplication on first argument)
real(kind=rp) function, public glsc2(a, b, n)
Weighted inner product .
subroutine, public subcol3(a, b, c, n)
Returns .
subroutine, public rone(a, n)
Set all elements to one.
subroutine flipvdp(b, ind, n)
Flip double precision vector b and ind.
subroutine, public x_update(a, b, c, c1, c2, n)
Returns .
subroutine, public add3(a, b, c, n)
Vector addition .
subroutine swapi4(b, ind, n)
sort single integer array acording to ind vector
integer function, public glimin(a, n)
Min of an integer vector of length n.
real(kind=rp) function, public glsum(a, n)
Sum a vector of length n.
subroutine, public sub3(a, b, c, n)
Vector subtraction .
subroutine, public addcol4(a, b, c, d, n)
Returns .
subroutine, public add2(a, b, n)
Vector addition .
subroutine, public cfill(a, c, n)
Set all elements to a constant c .
subroutine, public invcol3(a, b, c, n)
Invert a vector .
pure logical function sabscmp(x, y)
Return single precision absolute comparison .
subroutine, public add3s2(a, b, c, c1, c2, n)
Returns .
subroutine, public subcol4(a, b, c, d, n)
Returns .
subroutine sorti4(a, ind, n)
Heap Sort for single integer arrays.
subroutine, public addcol3(a, b, c, n)
Returns .
pure logical function dabscmp(x, y)
Return double precision absolute comparison .
subroutine, public invcol1(a, n)
Invert a vector .
subroutine, public masked_copy(a, b, mask, n, m)
Copy a masked vector .
real(kind=rp), parameter, public neko_m_ln2
subroutine, public chsign(a, n)
Change sign of vector .
subroutine, public col2(a, b, n)
Vector multiplication .
subroutine, public izero(a, n)
Zero an integer vector.
real(kind=rp) function, public glmax(a, n)
Max of a vector of length n.
subroutine, public copy(a, b, n)
Copy a vector .
subroutine, public add4(a, b, c, d, n)
Vector addition .
subroutine, public col3(a, b, c, n)
Vector multiplication with 3 vectors .
pure logical function drelcmp(x, y, eps)
Return double precision relative comparison .
real(kind=rp), parameter, public neko_eps
Machine epsilon .
subroutine, public vdot3(dot, u1, u2, u3, v1, v2, v3, n)
Compute a dot product (3-d version) assuming vector components etc.
pure logical function qrelcmp(x, y, eps)
Return quad precision relative comparison .
subroutine, public rzero(a, n)
Zero a real vector.
subroutine, public vdot2(dot, u1, u2, v1, v2, n)
Compute a dot product (2-d version) assuming vector components etc.
real(kind=rp) function, public vlmin(vec, n)
minimun value of a vector of length n
real(kind=rp) function, public vlmax(vec, n)
maximum value of a vector of length n
integer function, public glimax(a, n)
Max of an integer vector of length n.
subroutine sortrp(a, ind, n)
Heap Sort for double precision arrays.
subroutine, public sub2(a, b, n)
Vector substraction .
subroutine, public cfill_mask(a, c, size, mask, mask_size)
Fill a constant to a masked vector. .
subroutine, public add2s2(a, b, c1, n)
Vector addition with scalar multiplication (multiplication on second argument)
real(kind=rp) function, public glmin(a, n)
Min of a vector of length n.
subroutine, public vcross(u1, u2, u3, v1, v2, v3, w1, w2, w3, n)
Compute a cross product assuming vector components etc.
pure logical function srelcmp(x, y, eps)
Return single precision relative comparison .
real(kind=rp) function, public vlsc3(u, v, w, n)
Compute multiplication sum .
subroutine reordi4(b, ind, n)
reorder single integer array - inverse of swap
subroutine, public p_update(a, b, c, c1, c2, n)
Returns .
integer, parameter, public qp
integer, parameter, public i4
integer, parameter, public dp
integer, parameter, public sp
integer, parameter, public rp
Global precision used in computations.