Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
math.f90
Go to the documentation of this file.
1! Copyright (c) 2008-2020, UCHICAGO ARGONNE, LLC.
2!
3! The UChicago Argonne, LLC as Operator of Argonne National
4! Laboratory holds copyright in the Software. The copyright holder
5! reserves all rights except those expressly granted to licensees,
6! and U.S. Government license rights.
7!
8! Redistribution and use in source and binary forms, with or without
9! modification, are permitted provided that the following conditions
10! are met:
11!
12! 1. Redistributions of source code must retain the above copyright
13! notice, this list of conditions and the disclaimer below.
14!
15! 2. Redistributions in binary form must reproduce the above copyright
16! notice, this list of conditions and the disclaimer (as noted below)
17! in the documentation and/or other materials provided with the
18! distribution.
19!
20! 3. Neither the name of ANL nor the names of its contributors
21! may be used to endorse or promote products derived from this software
22! without specific prior written permission.
23!
24! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
25! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
26! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
27! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL
28! UCHICAGO ARGONNE, LLC, THE U.S. DEPARTMENT OF
29! ENERGY OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
30! SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
31! TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
32! DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
33! THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
34! (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
35! OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36!
37! Additional BSD Notice
38! ---------------------
39! 1. This notice is required to be provided under our contract with
40! the U.S. Department of Energy (DOE). This work was produced at
41! Argonne National Laboratory under Contract
42! No. DE-AC02-06CH11357 with the DOE.
43!
44! 2. Neither the United States Government nor UCHICAGO ARGONNE,
45! LLC nor any of their employees, makes any warranty,
46! express or implied, or assumes any liability or responsibility for the
47! accuracy, completeness, or usefulness of any information, apparatus,
48! product, or process disclosed, or represents that its use would not
49! infringe privately-owned rights.
50!
51! 3. Also, reference herein to any specific commercial products, process,
52! or services by trade name, trademark, manufacturer or otherwise does
53! not necessarily constitute or imply its endorsement, recommendation,
54! or favoring by the United States Government or UCHICAGO ARGONNE LLC.
55! The views and opinions of authors expressed
56! herein do not necessarily state or reflect those of the United States
57! Government or UCHICAGO ARGONNE, LLC, and shall
58! not be used for advertising or product endorsement purposes.
59!
60module math
61 use num_types, only : rp, dp, sp, qp, i4, xp
63 use mpi_f08, only: mpi_min, mpi_max, mpi_sum, mpi_in_place, mpi_integer, &
64 mpi_allreduce
65 implicit none
66 private
67
69 real(kind=rp), public, parameter :: neko_eps = epsilon(1.0_rp)
70
72 real(kind=rp), public, parameter :: neko_m_ln2 = log(2.0_rp)
73
75 real(kind=rp), public, parameter :: pi = 4._rp*atan(1._rp)
76
77 interface abscmp
78 module procedure sabscmp, dabscmp, qabscmp
79 end interface abscmp
80
81 interface sort
82 module procedure sortrp, sorti4
83 end interface sort
84
85 interface swap
86 module procedure swapdp, swapi4
87 end interface swap
88
89 interface reord
90 module procedure reorddp, reordi4
91 end interface reord
92
93 interface flipv
94 module procedure flipvdp, flipvi4
95 end interface flipv
96
97 interface relcmp
98 module procedure srelcmp, drelcmp, qrelcmp
99 end interface relcmp
100
101 public :: abscmp, rzero, izero, row_zero, rone, copy, cmult, cadd, cfill, &
112
113contains
114
116 pure function sabscmp(x, y, tol)
117 real(kind=sp), intent(in) :: x
118 real(kind=sp), intent(in) :: y
119 real(kind=sp), intent(in), optional :: tol
120 logical :: sabscmp
121
122 if (present(tol)) then
123 sabscmp = abs(x - y) .lt. tol
124 else
125 sabscmp = abs(x - y) .lt. neko_eps
126 end if
127
128 end function sabscmp
129
131 pure function dabscmp(x, y, tol)
132 real(kind=dp), intent(in) :: x
133 real(kind=dp), intent(in) :: y
134 real(kind=dp), intent(in), optional :: tol
135 logical :: dabscmp
136
137 if (present(tol)) then
138 dabscmp = abs(x - y) .lt. tol
139 else
140 dabscmp = abs(x - y) .lt. neko_eps
141 end if
142
143 end function dabscmp
144
146 pure function qabscmp(x, y, tol)
147 real(kind=qp), intent(in) :: x
148 real(kind=qp), intent(in) :: y
149 real(kind=qp), intent(in), optional :: tol
150 logical :: qabscmp
151
152 if (present(tol)) then
153 qabscmp = abs(x - y) .lt. tol
154 else
155 qabscmp = abs(x - y) .lt. neko_eps
156 end if
157
158 end function qabscmp
159
161 pure function srelcmp(x, y, eps)
162 real(kind=sp), intent(in) :: x
163 real(kind=sp), intent(in) :: y
164 real(kind=sp), intent(in), optional :: eps
165 logical :: srelcmp
166 if (present(eps)) then
167 srelcmp = abs(x - y) .le. eps*abs(y)
168 else
169 srelcmp = abs(x - y) .le. neko_eps*abs(y)
170 end if
171
172 end function srelcmp
173
175 pure function drelcmp(x, y, eps)
176 real(kind=dp), intent(in) :: x
177 real(kind=dp), intent(in) :: y
178 real(kind=dp), intent(in), optional :: eps
179 logical :: drelcmp
180 if (present(eps)) then
181 drelcmp = abs(x - y) .le. eps*abs(y)
182 else
183 drelcmp = abs(x - y) .le. neko_eps*abs(y)
184 end if
185
186 end function drelcmp
187
188
190 pure function qrelcmp(x, y, eps)
191 real(kind=qp), intent(in) :: x
192 real(kind=qp), intent(in) :: y
193 real(kind=qp), intent(in), optional :: eps
194 logical :: qrelcmp
195 if (present(eps)) then
196 qrelcmp = abs(x - y)/abs(y) .lt. eps
197 else
198 qrelcmp = abs(x - y)/abs(y) .lt. neko_eps
199 end if
200
201 end function qrelcmp
202
204 subroutine rzero(a, n)
205 integer, intent(in) :: n
206 real(kind=rp), dimension(n), intent(inout) :: a
207 integer :: i
208
209 do i = 1, n
210 a(i) = 0.0_rp
211 end do
212 end subroutine rzero
213
215 subroutine izero(a, n)
216 integer, intent(in) :: n
217 integer, dimension(n), intent(inout) :: a
218 integer :: i
219
220 do i = 1, n
221 a(i) = 0
222 end do
223 end subroutine izero
224
226 subroutine row_zero(a, m, n, e)
227 integer, intent(in) :: m, n, e
228 real(kind=rp), intent(inout) :: a(m,n)
229 integer :: j
230
231 do j = 1,n
232 a(e,j) = 0.0_rp
233 end do
234 end subroutine row_zero
235
237 subroutine rone(a, n)
238 integer, intent(in) :: n
239 real(kind=rp), dimension(n), intent(inout) :: a
240 integer :: i
241
242 do i = 1, n
243 a(i) = 1.0_rp
244 end do
245 end subroutine rone
246
248 subroutine copy(a, b, n)
249 integer, intent(in) :: n
250 real(kind=rp), dimension(n), intent(in) :: b
251 real(kind=rp), dimension(n), intent(inout) :: a
252 integer :: i
253
254 do i = 1, n
255 a(i) = b(i)
256 end do
257
258 end subroutine copy
259
267 subroutine masked_copy_0(a, b, mask, n, n_mask)
268 integer, intent(in) :: n, n_mask
269 real(kind=rp), dimension(n), intent(in) :: b
270 real(kind=rp), dimension(n), intent(inout) :: a
271 integer, dimension(0:n_mask) :: mask
272 integer :: i, j
273
274 do i = 1, n_mask
275 j = mask(i)
276 a(j) = b(j)
277 end do
278
279 end subroutine masked_copy_0
280
288 subroutine masked_copy(a, b, mask, n, n_mask)
289 integer, intent(in) :: n, n_mask
290 real(kind=rp), dimension(n), intent(in) :: b
291 real(kind=rp), dimension(n), intent(inout) :: a
292 integer, dimension(n_mask) :: mask
293 integer :: i, j
294
295 do i = 1, n_mask
296 j = mask(i)
297 a(j) = b(j)
298 end do
299
300 end subroutine masked_copy
301
311 subroutine masked_gather_copy_0(a, b, mask, n, n_mask)
312 integer, intent(in) :: n, n_mask
313 real(kind=rp), dimension(n), intent(in) :: b
314 real(kind=rp), dimension(n_mask), intent(inout) :: a
315 integer, dimension(0:n_mask) :: mask
316 integer :: i, j
317
318 do i = 1, n_mask
319 j = mask(i)
320 a(i) = b(j)
321 end do
322
323 end subroutine masked_gather_copy_0
324
334 subroutine masked_gather_copy(a, b, mask, n, n_mask)
335 integer, intent(in) :: n, n_mask
336 real(kind=rp), dimension(n), intent(in) :: b
337 real(kind=rp), dimension(n_mask), intent(inout) :: a
338 integer, dimension(n_mask) :: mask
339 integer :: i, j
340
341 do i = 1, n_mask
342 j = mask(i)
343 a(i) = b(j)
344 end do
345
346 end subroutine masked_gather_copy
347
357 subroutine masked_scatter_copy_0(a, b, mask, n, n_mask)
358 integer, intent(in) :: n, n_mask
359 real(kind=rp), dimension(n_mask), intent(in) :: b
360 real(kind=rp), dimension(n), intent(inout) :: a
361 integer, dimension(0:n_mask) :: mask
362 integer :: i, j
363
364 do i = 1, n_mask
365 j = mask(i)
366 a(j) = b(i)
367 end do
368
369 end subroutine masked_scatter_copy_0
370
380 subroutine masked_scatter_copy(a, b, mask, n, n_mask)
381 integer, intent(in) :: n, n_mask
382 real(kind=rp), dimension(n_mask), intent(in) :: b
383 real(kind=rp), dimension(n), intent(inout) :: a
384 integer, dimension(n_mask) :: mask
385 integer :: i, j
386
387 do i = 1, n_mask
388 j = mask(i)
389 a(j) = b(i)
390 end do
391
392 end subroutine masked_scatter_copy
393
396 subroutine cfill_mask(a, c, n, mask, n_mask)
397 integer, intent(in) :: n, n_mask
398 real(kind=rp), dimension(n), intent(inout) :: a
399 real(kind=rp), intent(in) :: c
400 integer, dimension(n_mask), intent(in) :: mask
401 integer :: i
402
403 do i = 1, n_mask
404 a(mask(i)) = c
405 end do
406
407 end subroutine cfill_mask
408
410 subroutine cmult(a, c, n)
411 integer, intent(in) :: n
412 real(kind=rp), dimension(n), intent(inout) :: a
413 real(kind=rp), intent(in) :: c
414 integer :: i
415
416 do i = 1, n
417 a(i) = c * a(i)
418 end do
419 end subroutine cmult
420
422 subroutine cmult2(a, b, c, n)
423 integer, intent(in) :: n
424 real(kind=rp), dimension(n), intent(inout) :: a
425 real(kind=rp), dimension(n), intent(in) :: b
426 real(kind=rp), intent(in) :: c
427 integer :: i
428
429 do i = 1, n
430 a(i) = c * b(i)
431 end do
432
433 end subroutine cmult2
434
436 subroutine cdiv(a, c, n)
437 integer, intent(in) :: n
438 real(kind=rp), dimension(n), intent(inout) :: a
439 real(kind=rp), intent(in) :: c
440 integer :: i
441
442 do i = 1, n
443 a(i) = c / a(i)
444 end do
445 end subroutine cdiv
446
448 subroutine cdiv2(a, b, c, n)
449 integer, intent(in) :: n
450 real(kind=rp), dimension(n), intent(inout) :: a
451 real(kind=rp), dimension(n), intent(in) :: b
452 real(kind=rp), intent(in) :: c
453 integer :: i
454
455 do i = 1, n
456 a(i) = c / b(i)
457 end do
458 end subroutine cdiv2
459
461 subroutine cadd(a, s, n)
462 integer, intent(in) :: n
463 real(kind=rp), dimension(n), intent(inout) :: a
464 real(kind=rp), intent(in) :: s
465 integer :: i
466
467 do i = 1, n
468 a(i) = a(i) + s
469 end do
470 end subroutine cadd
471
473 subroutine cadd2(a, b, s, n)
474 integer, intent(in) :: n
475 real(kind=rp), dimension(n), intent(inout) :: a
476 real(kind=rp), dimension(n), intent(in) :: b
477 real(kind=rp), intent(in) :: s
478 integer :: i
479
480 do i = 1, n
481 a(i) = b(i) + s
482 end do
483 end subroutine cadd2
484
486 subroutine cfill(a, c, n)
487 integer, intent(in) :: n
488 real(kind=rp), dimension(n), intent(inout) :: a
489 real(kind=rp), intent(in) :: c
490 integer :: i
491
492 do i = 1, n
493 a(i) = c
494 end do
495 end subroutine cfill
496
498 function glsum(a, n)
499 integer, intent(in) :: n
500 real(kind=rp), dimension(n) :: a
501 real(kind=rp) :: glsum
502 real(kind=xp) :: tmp
503 integer :: i, ierr
504 tmp = 0.0_rp
505 do i = 1, n
506 tmp = tmp + a(i)
507 end do
508 call mpi_allreduce(mpi_in_place, tmp, 1, &
509 mpi_extra_precision, mpi_sum, neko_comm, ierr)
510 glsum = tmp
511
512 end function glsum
513
515 function glmax(a, n)
516 integer, intent(in) :: n
517 real(kind=rp), dimension(n) :: a
518 real(kind=rp) :: tmp, glmax
519 integer :: i, ierr
520
521 tmp = -huge(0.0_rp)
522 do i = 1, n
523 tmp = max(tmp,a(i))
524 end do
525 call mpi_allreduce(tmp, glmax, 1, &
526 mpi_real_precision, mpi_max, neko_comm, ierr)
527 end function glmax
528
530 function glimax(a, n)
531 integer, intent(in) :: n
532 integer, dimension(n) :: a
533 integer :: tmp, glimax
534 integer :: i, ierr
535
536 tmp = -huge(0)
537 do i = 1, n
538 tmp = max(tmp,a(i))
539 end do
540 call mpi_allreduce(tmp, glimax, 1, &
541 mpi_integer, mpi_max, neko_comm, ierr)
542 end function glimax
543
545 function glmin(a, n)
546 integer, intent(in) :: n
547 real(kind=rp), dimension(n) :: a
548 real(kind=rp) :: tmp, glmin
549 integer :: i, ierr
550
551 tmp = huge(0.0_rp)
552 do i = 1, n
553 tmp = min(tmp,a(i))
554 end do
555 call mpi_allreduce(tmp, glmin, 1, &
556 mpi_real_precision, mpi_min, neko_comm, ierr)
557 end function glmin
558
560 function glimin(a, n)
561 integer, intent(in) :: n
562 integer, dimension(n) :: a
563 integer :: tmp, glimin
564 integer :: i, ierr
565
566 tmp = huge(0)
567 do i = 1, n
568 tmp = min(tmp,a(i))
569 end do
570 call mpi_allreduce(tmp, glimin, 1, &
571 mpi_integer, mpi_min, neko_comm, ierr)
572 end function glimin
573
574
575
576
578 subroutine chsign(a, n)
579 integer, intent(in) :: n
580 real(kind=rp), dimension(n), intent(inout) :: a
581 integer :: i
582
583 do i = 1, n
584 a(i) = -a(i)
585 end do
586
587 end subroutine chsign
588
590 function vlmax(vec,n) result(tmax)
591 integer :: n, i
592 real(kind=rp), intent(in) :: vec(n)
593 real(kind=rp) :: tmax
594 tmax = real(-99d20, rp)
595 do i = 1, n
596 tmax = max(tmax, vec(i))
597 end do
598 end function vlmax
599
601 function vlmin(vec,n) result(tmin)
602 integer, intent(in) :: n
603 real(kind=rp), intent(in) :: vec(n)
604 real(kind=rp) :: tmin
605 integer :: i
606 tmin = real(99.0e20, rp)
607 do i = 1, n
608 tmin = min(tmin, vec(i))
609 end do
610 end function vlmin
611
613 subroutine invcol1(a, n)
614 integer, intent(in) :: n
615 real(kind=rp), dimension(n), intent(inout) :: a
616 integer :: i
617
618 do i = 1, n
619 a(i) = 1.0_xp / real(a(i),xp)
620 end do
621
622 end subroutine invcol1
623
625 subroutine invcol3(a, b, c, n)
626 integer, intent(in) :: n
627 real(kind=rp), dimension(n), intent(inout) :: a
628 real(kind=rp), dimension(n), intent(in) :: b,c
629 integer :: i
630
631 do i = 1, n
632 a(i) = real(b(i),xp) / c(i)
633 end do
634
635 end subroutine invcol3
636
638 subroutine invers2(a, b, n)
639 integer, intent(in) :: n
640 real(kind=rp), dimension(n), intent(inout) :: a
641 real(kind=rp), dimension(n), intent(in) :: b
642 integer :: i
643
644 do i = 1, n
645 a(i) = 1.0_xp / real(b(i),xp)
646 end do
647
648 end subroutine invers2
649
652 subroutine vcross(u1, u2, u3, v1, v2, v3, w1, w2, w3, n)
653 integer, intent(in) :: n
654 real(kind=rp), dimension(n), intent(in) :: v1, v2, v3
655 real(kind=rp), dimension(n), intent(in) :: w1, w2, w3
656 real(kind=rp), dimension(n), intent(out) :: u1, u2, u3
657 integer :: i
658
659 do i = 1, n
660 u1(i) = v2(i)*w3(i) - v3(i)*w2(i)
661 u2(i) = v3(i)*w1(i) - v1(i)*w3(i)
662 u3(i) = v1(i)*w2(i) - v2(i)*w1(i)
663 end do
664
665 end subroutine vcross
666
669 subroutine vdot2(dot, u1, u2, v1, v2, n)
670 integer, intent(in) :: n
671 real(kind=rp), dimension(n), intent(in) :: u1, u2
672 real(kind=rp), dimension(n), intent(in) :: v1, v2
673 real(kind=rp), dimension(n), intent(out) :: dot
674 integer :: i
675 do i = 1, n
676 dot(i) = u1(i)*v1(i) + u2(i)*v2(i)
677 end do
678
679 end subroutine vdot2
680
683 subroutine vdot3(dot, u1, u2, u3, v1, v2, v3, n)
684 integer, intent(in) :: n
685 real(kind=rp), dimension(n), intent(in) :: u1, u2, u3
686 real(kind=rp), dimension(n), intent(in) :: v1, v2, v3
687 real(kind=rp), dimension(n), intent(out) :: dot
688 integer :: i
689
690 do i = 1, n
691 dot(i) = u1(i)*v1(i) + u2(i)*v2(i) + u3(i)*v3(i)
692 end do
693
694 end subroutine vdot3
695
697 function vlsc3(u, v, w, n) result(s)
698 integer, intent(in) :: n
699 real(kind=rp), dimension(n), intent(in) :: u, v, w
700 real(kind=rp) :: s
701 integer :: i
702
703 s = 0.0_rp
704 do i = 1, n
705 s = s + u(i)*v(i)*w(i)
706 end do
707
708 end function vlsc3
709
711 function vlsc2(u, v, n) result(s)
712 integer, intent(in) :: n
713 real(kind=rp), dimension(n), intent(in) :: u, v
714 real(kind=rp) :: s
715 integer :: i
716
717 s = 0.0_rp
718 do i = 1, n
719 s = s + u(i)*v(i)
720 end do
721
722 end function vlsc2
723
725 subroutine add2(a, b, n)
726 integer, intent(in) :: n
727 real(kind=rp), dimension(n), intent(inout) :: a
728 real(kind=rp), dimension(n), intent(in) :: b
729 integer :: i
730
731 do i = 1, n
732 a(i) = a(i) + b(i)
733 end do
734
735 end subroutine add2
736
738 subroutine add3(a, b, c, n)
739 integer, intent(in) :: n
740 real(kind=rp), dimension(n), intent(inout) :: a
741 real(kind=rp), dimension(n), intent(in) :: b
742 real(kind=rp), dimension(n), intent(in) :: c
743 integer :: i
744
745 do i = 1, n
746 a(i) = b(i) + c(i)
747 end do
748
749 end subroutine add3
750
752 subroutine add4(a, b, c, d, n)
753 integer, intent(in) :: n
754 real(kind=rp), dimension(n), intent(out) :: a
755 real(kind=rp), dimension(n), intent(in) :: d
756 real(kind=rp), dimension(n), intent(in) :: c
757 real(kind=rp), dimension(n), intent(in) :: b
758 integer :: i
759
760 do i = 1, n
761 a(i) = b(i) + c(i) + d(i)
762 end do
763
764 end subroutine add4
765
767 subroutine sub2(a, b, n)
768 integer, intent(in) :: n
769 real(kind=rp), dimension(n), intent(inout) :: a
770 real(kind=rp), dimension(n), intent(in) :: b
771 integer :: i
772
773 do i = 1, n
774 a(i) = a(i) - b(i)
775 end do
776
777 end subroutine sub2
778
780 subroutine sub3(a, b, c, n)
781 integer, intent(in) :: n
782 real(kind=rp), dimension(n), intent(inout) :: a
783 real(kind=rp), dimension(n), intent(in) :: b
784 real(kind=rp), dimension(n), intent(in) :: c
785 integer :: i
786
787 do i = 1, n
788 a(i) = b(i) - c(i)
789 end do
790
791 end subroutine sub3
792
793
796 subroutine add2s1(a, b, c1, n)
797 integer, intent(in) :: n
798 real(kind=rp), dimension(n), intent(inout) :: a
799 real(kind=rp), dimension(n), intent(in) :: b
800 real(kind=rp), intent(in) :: c1
801 integer :: i
802
803 do i = 1, n
804 a(i) = c1 * a(i) + b(i)
805 end do
806
807 end subroutine add2s1
808
811 subroutine add2s2(a, b, c1, n)
812 integer, intent(in) :: n
813 real(kind=rp), dimension(n), intent(inout) :: a
814 real(kind=rp), dimension(n), intent(in) :: b
815 real(kind=rp), intent(in) :: c1
816 integer :: i
817
818 do i = 1, n
819 a(i) = a(i) + c1 * b(i)
820 end do
821
822 end subroutine add2s2
823
825 subroutine addsqr2s2(a, b, c1, n)
826 integer, intent(in) :: n
827 real(kind=rp), dimension(n), intent(inout) :: a
828 real(kind=rp), dimension(n), intent(in) :: b
829 real(kind=rp), intent(in) :: c1
830 integer :: i
831
832 do i = 1,n
833 a(i) = a(i) + c1 * ( b(i) * b(i) )
834 end do
835
836 end subroutine addsqr2s2
837
839 subroutine invcol2(a, b, n)
840 integer, intent(in) :: n
841 real(kind=rp), dimension(n), intent(inout) :: a
842 real(kind=rp), dimension(n), intent(in) :: b
843 integer :: i
844
845 do i = 1, n
846 a(i) = real(a(i),xp) /b(i)
847 end do
848
849 end subroutine invcol2
850
851
853 subroutine col2(a, b, n)
854 integer, intent(in) :: n
855 real(kind=rp), dimension(n), intent(inout) :: a
856 real(kind=rp), dimension(n), intent(in) :: b
857 integer :: i
858
859 do i = 1, n
860 a(i) = a(i) * b(i)
861 end do
862
863 end subroutine col2
864
866 subroutine col3(a, b, c, n)
867 integer, intent(in) :: n
868 real(kind=rp), dimension(n), intent(inout) :: a
869 real(kind=rp), dimension(n), intent(in) :: b
870 real(kind=rp), dimension(n), intent(in) :: c
871 integer :: i
872
873 do i = 1, n
874 a(i) = b(i) * c(i)
875 end do
876
877 end subroutine col3
878
880 subroutine subcol3(a, b, c, n)
881 integer, intent(in) :: n
882 real(kind=rp), dimension(n), intent(inout) :: a
883 real(kind=rp), dimension(n), intent(in) :: b
884 real(kind=rp), dimension(n), intent(in) :: c
885 integer :: i
886
887 do i = 1,n
888 a(i) = a(i) - b(i) * c(i)
889 end do
890
891 end subroutine subcol3
892
894 subroutine add3s2(a, b, c, c1, c2 ,n)
895 integer, intent(in) :: n
896 real(kind=rp), dimension(n), intent(inout) :: a
897 real(kind=rp), dimension(n), intent(in) :: b
898 real(kind=rp), dimension(n), intent(in) :: c
899 real(kind=rp), intent(in) :: c1, c2
900 integer :: i
901
902 do i = 1,n
903 a(i) = c1 * b(i) + c2 * c(i)
904 end do
905
906 end subroutine add3s2
907
909 subroutine add4s3(a, b, c, d, c1, c2, c3, n)
910 integer, intent(in) :: n
911 real(kind=rp), dimension(n), intent(inout) :: a
912 real(kind=rp), dimension(n), intent(in) :: b
913 real(kind=rp), dimension(n), intent(in) :: c
914 real(kind=rp), dimension(n), intent(in) :: d
915 real(kind=rp), intent(in) :: c1, c2, c3
916 integer :: i
917
918 do i = 1,n
919 a(i) = c1 * b(i) + c2 * c(i) + c3 * d(i)
920 end do
921
922 end subroutine add4s3
923
925 subroutine add5s4(a, b, c, d, e, c1, c2, c3, c4, n)
926 integer, intent(in) :: n
927 real(kind=rp), dimension(n), intent(inout) :: a
928 real(kind=rp), dimension(n), intent(in) :: b
929 real(kind=rp), dimension(n), intent(in) :: c
930 real(kind=rp), dimension(n), intent(in) :: d
931 real(kind=rp), dimension(n), intent(in) :: e
932 real(kind=rp), intent(in) :: c1, c2, c3, c4
933 integer :: i
934
935 do i = 1,n
936 a(i) = a(i) + c1 * b(i) + c2 * c(i) + c3 * d(i) + c4 * e(i)
937 end do
938
939 end subroutine add5s4
940
942 subroutine subcol4(a, b, c, d, n)
943 integer, intent(in) :: n
944 real(kind=rp), dimension(n), intent(inout) :: a
945 real(kind=rp), dimension(n), intent(in) :: b
946 real(kind=rp), dimension(n), intent(in) :: c
947 real(kind=rp), dimension(n), intent(in) :: d
948 integer :: i
949
950 do i = 1,n
951 a(i) = a(i) - b(i) * c(i) * d(i)
952 end do
953
954 end subroutine subcol4
955
957 subroutine addcol3(a, b, c, n)
958 integer, intent(in) :: n
959 real(kind=rp), dimension(n), intent(inout) :: a
960 real(kind=rp), dimension(n), intent(in) :: b
961 real(kind=rp), dimension(n), intent(in) :: c
962 integer :: i
963
964 do i = 1,n
965 a(i) = a(i) + b(i) * c(i)
966 end do
967
968 end subroutine addcol3
969
971 subroutine addcol4(a, b, c, d, n)
972 integer, intent(in) :: n
973 real(kind=rp), dimension(n), intent(inout) :: a
974 real(kind=rp), dimension(n), intent(in) :: b
975 real(kind=rp), dimension(n), intent(in) :: c
976 real(kind=rp), dimension(n), intent(in) :: d
977 integer :: i
978
979 do i = 1,n
980 a(i) = a(i) + b(i) * c(i) * d(i)
981 end do
982
983 end subroutine addcol4
984
986 subroutine addcol3s2(a, b, c, s, n)
987 integer, intent(in) :: n
988 real(kind=rp), dimension(n), intent(inout) :: a
989 real(kind=rp), dimension(n), intent(in) :: b
990 real(kind=rp), dimension(n), intent(in) :: c
991 real(kind=rp), intent(in) :: s
992 integer :: i
993
994 do i = 1,n
995 a(i) = a(i) + s * b(i) * c(i)
996 end do
997
998 end subroutine addcol3s2
999
1001 subroutine ascol5(a, b, c, d, e, n)
1002 integer, intent(in) :: n
1003 real(kind=rp), dimension(n), intent(inout) :: a
1004 real(kind=rp), dimension(n), intent(in) :: b
1005 real(kind=rp), dimension(n), intent(in) :: c
1006 real(kind=rp), dimension(n), intent(in) :: d
1007 real(kind=rp), dimension(n), intent(in) :: e
1008 integer :: i
1009
1010 do i = 1,n
1011 a(i) = b(i)*c(i)-d(i)*e(i)
1012 end do
1013
1014 end subroutine ascol5
1015
1017 subroutine p_update(a, b, c, c1, c2, n)
1018 integer, intent(in) :: n
1019 real(kind=rp), dimension(n), intent(inout) :: a
1020 real(kind=rp), dimension(n), intent(in) :: b
1021 real(kind=rp), dimension(n), intent(in) :: c
1022 real(kind=rp), intent(in) :: c1, c2
1023 integer :: i
1024
1025 do i = 1,n
1026 a(i) = b(i) + c1*(a(i)-c2*c(i))
1027 end do
1028
1029 end subroutine p_update
1030
1032 subroutine x_update(a, b, c, c1, c2, n)
1033 integer, intent(in) :: n
1034 real(kind=rp), dimension(n), intent(inout) :: a
1035 real(kind=rp), dimension(n), intent(in) :: b
1036 real(kind=rp), dimension(n), intent(in) :: c
1037 real(kind=rp), intent(in) :: c1, c2
1038 integer :: i
1039
1040 do i = 1,n
1041 a(i) = a(i) + c1*b(i)+c2*c(i)
1042 end do
1043
1044 end subroutine x_update
1045
1047 function glsc2(a, b, n)
1048 integer, intent(in) :: n
1049 real(kind=rp), dimension(n), intent(in) :: a
1050 real(kind=rp), dimension(n), intent(in) :: b
1051 real(kind=rp) :: glsc2
1052 real(kind=xp) :: tmp
1053 integer :: i, ierr
1054
1055 tmp = 0.0_xp
1056 do i = 1, n
1057 tmp = tmp + a(i) * b(i)
1058 end do
1059
1060 call mpi_allreduce(mpi_in_place, tmp, 1, &
1061 mpi_extra_precision, mpi_sum, neko_comm, ierr)
1062 glsc2 = tmp
1063 end function glsc2
1064
1066 function glsc3(a, b, c, n)
1067 integer, intent(in) :: n
1068 real(kind=rp), dimension(n), intent(in) :: a
1069 real(kind=rp), dimension(n), intent(in) :: b
1070 real(kind=rp), dimension(n), intent(in) :: c
1071 real(kind=rp) :: glsc3
1072 real(kind=xp) :: tmp
1073 integer :: i, ierr
1074
1075 tmp = 0.0_xp
1076 do i = 1, n
1077 tmp = tmp + a(i) * b(i) * c(i)
1078 end do
1079
1080 call mpi_allreduce(mpi_in_place, tmp, 1, &
1081 mpi_extra_precision, mpi_sum, neko_comm, ierr)
1082 glsc3 = tmp
1083
1084 end function glsc3
1085 function glsc4(a, b, c, d, n)
1086 integer, intent(in) :: n
1087 real(kind=rp), dimension(n), intent(in) :: a
1088 real(kind=rp), dimension(n), intent(in) :: b
1089 real(kind=rp), dimension(n), intent(in) :: c
1090 real(kind=rp), dimension(n), intent(in) :: d
1091 real(kind=rp) :: glsc4
1092 real(kind=xp) :: tmp
1093 integer :: i, ierr
1094
1095 tmp = 0.0_xp
1096 do i = 1, n
1097 tmp = tmp + a(i) * b(i) * c(i) * d(i)
1098 end do
1099
1100 call mpi_allreduce(mpi_in_place, tmp, 1, &
1101 mpi_extra_precision, mpi_sum, neko_comm, ierr)
1102 glsc4 = tmp
1103
1104 end function glsc4
1105
1108 function glsubnorm(a, b, n)
1109 integer, intent(in) :: n
1110 real(kind=rp), dimension(n), intent(in) :: a
1111 real(kind=rp), dimension(n), intent(in) :: b
1112 real(kind=rp) :: glsubnorm
1113 real(kind=xp) :: tmp
1114 integer :: i, ierr
1115
1116 tmp = 0.0_xp
1117 do i = 1, n
1118 tmp = tmp + (a(i) - b(i))**2
1119 end do
1120
1121 call mpi_allreduce(mpi_in_place, tmp, 1, &
1122 mpi_extra_precision, mpi_sum, neko_comm, ierr)
1123 glsubnorm = sqrt(tmp)
1124
1125 end function glsubnorm
1126
1132 subroutine sortrp(a, ind, n)
1133 integer, intent(in) :: n
1134 real(kind=rp), intent(inout) :: a(n)
1135 integer, intent(out) :: ind(n)
1136 real(kind=rp) :: aa
1137 integer :: j, ir, i, ii, l
1138
1139 do j = 1, n
1140 ind(j) = j
1141 end do
1142
1143 if (n .le. 1) return
1144
1145
1146 l = n/2+1
1147 ir = n
1148 do while (.true.)
1149 if (l .gt. 1) then
1150 l = l-1
1151 aa = a(l)
1152 ii = ind(l)
1153 else
1154 aa = a(ir)
1155 ii = ind(ir)
1156 a(ir) = a(1)
1157 ind(ir) = ind(1)
1158 ir = ir - 1
1159 if (ir .eq. 1) then
1160 a(1) = aa
1161 ind(1) = ii
1162 return
1163 end if
1164 end if
1165 i = l
1166 j = l+l
1167 do while (j .le. ir)
1168 if (j .lt. ir) then
1169 if ( a(j) .lt. a(j+1) ) j = j + 1
1170 end if
1171 if (aa .lt. a(j)) then
1172 a(i) = a(j)
1173 ind(i) = ind(j)
1174 i = j
1175 j = j+j
1176 else
1177 j = ir+1
1178 end if
1179 end do
1180 a(i) = aa
1181 ind(i) = ii
1182 end do
1183 end subroutine sortrp
1184
1190 subroutine sorti4(a, ind, n)
1191 integer, intent(in) :: n
1192 integer(i4), intent(inout) :: a(n)
1193 integer, intent(out) :: ind(n)
1194 integer(i4) :: aa
1195 integer :: j, ir, i, ii, l
1196
1197 do j = 1, n
1198 ind(j) = j
1199 end do
1200
1201 if (n .le. 1) return
1202
1203 l = n/2+1
1204 ir = n
1205 do while (.true.)
1206 if (l .gt. 1) then
1207 l = l - 1
1208 aa = a(l)
1209 ii = ind(l)
1210 else
1211 aa = a(ir)
1212 ii = ind(ir)
1213 a(ir) = a( 1)
1214 ind(ir) = ind( 1)
1215 ir = ir - 1
1216 if (ir .eq. 1) then
1217 a(1) = aa
1218 ind(1) = ii
1219 return
1220 end if
1221 end if
1222 i = l
1223 j = l + l
1224 do while (j .le. ir)
1225 if (j .lt. ir) then
1226 if ( a(j) .lt. a(j + 1) ) j = j + 1
1227 end if
1228 if (aa .lt. a(j)) then
1229 a(i) = a(j)
1230 ind(i) = ind(j)
1231 i = j
1232 j = j + j
1233 else
1234 j = ir + 1
1235 end if
1236 end do
1237 a(i) = aa
1238 ind(i) = ii
1239 end do
1240 end subroutine sorti4
1241
1246 subroutine swapdp(b, ind, n)
1247 integer, intent(in) :: n
1248 real(kind=rp), intent(inout) :: b(n)
1249 integer, intent(in) :: ind(n)
1250 real(kind=rp) :: temp(n)
1251 integer :: i, jj
1252
1253 do i = 1, n
1254 temp(i) = b(i)
1255 end do
1256 do i = 1, n
1257 jj = ind(i)
1258 b(i) = temp(jj)
1259 end do
1260 end subroutine swapdp
1261
1266 subroutine swapi4(b, ind, n)
1267 integer, intent(in) :: n
1268 integer(i4), intent(inout) :: b(n)
1269 integer, intent(in) :: ind(n)
1270 integer(i4) :: temp(n)
1271 integer :: i, jj
1272
1273 do i = 1, n
1274 temp(i) = b(i)
1275 end do
1276 do i = 1, n
1277 jj = ind(i)
1278 b(i) = temp(jj)
1279 end do
1280 end subroutine swapi4
1281
1286 subroutine reorddp(b, ind, n)
1287 integer, intent(in) :: n
1288 real(kind=rp), intent(inout) :: b(n)
1289 integer, intent(in) :: ind(n)
1290 real(kind=rp) :: temp(n)
1291 integer :: i, jj
1292
1293 do i = 1, n
1294 temp(i) = b(i)
1295 end do
1296 do i = 1, n
1297 jj = ind(i)
1298 b(jj) = temp(i)
1299 end do
1300 end subroutine reorddp
1301
1306 subroutine reordi4(b, ind, n)
1307 integer, intent(in) :: n
1308 integer(i4), intent(inout) :: b(n)
1309 integer, intent(in) :: ind(n)
1310 integer(i4) :: temp(n)
1311 integer :: i, jj
1312
1313 do i = 1, n
1314 temp(i) = b(i)
1315 end do
1316 do i = 1, n
1317 jj = ind(i)
1318 b(jj) = temp(i)
1319 end do
1320 end subroutine reordi4
1321
1326 subroutine flipvdp(b, ind, n)
1327 integer, intent(in) :: n
1328 real(kind=rp), intent(inout) :: b(n)
1329 integer, intent(inout) :: ind(n)
1330 real(kind=rp) :: temp(n)
1331 integer :: tempind(n)
1332 integer :: i, jj
1333
1334 do i = 1, n
1335 jj = n+1-i
1336 temp(jj) = b(i)
1337 tempind(jj) = ind(i)
1338 end do
1339 do i = 1,n
1340 b(i) = temp(i)
1341 ind(i) = tempind(i)
1342 end do
1343 end subroutine flipvdp
1344
1349 subroutine flipvi4(b, ind, n)
1350 integer, intent(in) :: n
1351 integer(i4), intent(inout) :: b(n)
1352 integer, intent(inout) :: ind(n)
1353 integer(i4) :: temp(n)
1354 integer :: tempind(n)
1355 integer :: i, jj
1356
1357 do i = 1, n
1358 jj = n+1-i
1359 temp(jj) = b(i)
1360 tempind(jj) = ind(i)
1361 end do
1362 do i = 1,n
1363 b(i) = temp(i)
1364 ind(i) = tempind(i)
1365 end do
1366 end subroutine flipvi4
1367
1371 subroutine absval(a, n)
1372 integer, intent(in) :: n
1373 real(kind=rp), dimension(n), intent(inout) :: a
1374 integer :: i
1375 do i = 1, n
1376 a(i) = abs(a(i))
1377 end do
1378 end subroutine absval
1379
1380 ! ========================================================================== !
1381 ! Point-wise operations
1382
1384 subroutine pwmax2(a, b, n)
1385 integer, intent(in) :: n
1386 real(kind=rp), dimension(n), intent(inout) :: a
1387 real(kind=rp), dimension(n), intent(in) :: b
1388 integer :: i
1389
1390 do i = 1, n
1391 a(i) = max(a(i), b(i))
1392 end do
1393 end subroutine pwmax2
1394
1396 subroutine pwmax3(a, b, c, n)
1397 integer, intent(in) :: n
1398 real(kind=rp), dimension(n), intent(inout) :: a
1399 real(kind=rp), dimension(n), intent(in) :: b, c
1400 integer :: i
1401
1402 do i = 1, n
1403 a(i) = max(b(i), c(i))
1404 end do
1405 end subroutine pwmax3
1406
1408 subroutine cpwmax2(a, b, n)
1409 integer, intent(in) :: n
1410 real(kind=rp), dimension(n), intent(inout) :: a
1411 real(kind=rp), intent(in) :: b
1412 integer :: i
1413
1414 do i = 1, n
1415 a(i) = max(a(i), b)
1416 end do
1417 end subroutine cpwmax2
1418
1420 subroutine cpwmax3(a, b, c, n)
1421 integer, intent(in) :: n
1422 real(kind=rp), dimension(n), intent(inout) :: a
1423 real(kind=rp), dimension(n), intent(in) :: b
1424 real(kind=rp), intent(in) :: c
1425 integer :: i
1426
1427 do i = 1, n
1428 a(i) = max(b(i), c)
1429 end do
1430 end subroutine cpwmax3
1431
1433 subroutine pwmin2(a, b, n)
1434 integer, intent(in) :: n
1435 real(kind=rp), dimension(n), intent(inout) :: a
1436 real(kind=rp), dimension(n), intent(in) :: b
1437 integer :: i
1438
1439 do i = 1, n
1440 a(i) = min(a(i), b(i))
1441 end do
1442 end subroutine pwmin2
1443
1445 subroutine pwmin3(a, b, c, n)
1446 integer, intent(in) :: n
1447 real(kind=rp), dimension(n), intent(inout) :: a
1448 real(kind=rp), dimension(n), intent(in) :: b, c
1449 integer :: i
1450
1451 do i = 1, n
1452 a(i) = min(b(i), c(i))
1453 end do
1454 end subroutine pwmin3
1455
1457 subroutine cpwmin2(a, b, n)
1458 integer, intent(in) :: n
1459 real(kind=rp), dimension(n), intent(inout) :: a
1460 real(kind=rp), intent(in) :: b
1461 integer :: i
1462
1463 do i = 1, n
1464 a(i) = min(a(i), b)
1465 end do
1466 end subroutine cpwmin2
1467
1469 subroutine cpwmin3(a, b, c, n)
1470 integer, intent(in) :: n
1471 real(kind=rp), dimension(n), intent(inout) :: a
1472 real(kind=rp), dimension(n), intent(in) :: b
1473 real(kind=rp), intent(in) :: c
1474 integer :: i
1475
1476 do i = 1, n
1477 a(i) = min(b(i), c)
1478 end do
1479 end subroutine cpwmin3
1480
1481 ! M33INV and M44INV by David G. Simpson pure function version from
1482 ! https://fortranwiki.org/fortran/show/Matrix+inversion
1483 ! Invert 3x3 matrix
1484 function matinv39(a11, a12, a13, a21, a22, a23, a31, a32, a33) &
1485 result(b)
1486 real(kind=rp), intent(in) :: a11, a12, a13, a21, a22, a23, a31, a32, a33
1487 real(xp) :: a(3,3) !! Matrix
1488 real(rp) :: b(3,3) !! Inverse matrix
1489 a(1,1) = a11
1490 a(1,2) = a12
1491 a(1,3) = a13
1492 a(2,1) = a21
1493 a(2,2) = a22
1494 a(2,3) = a23
1495 a(3,1) = a31
1496 a(3,2) = a32
1497 a(3,3) = a33
1498 b = matinv3(a)
1499 end function matinv39
1500
1505 function matinv3(A) result(B)
1506 !! Performs a direct calculation of the inverse of a 3×3 matrix.
1507 real(kind=xp), intent(in) :: a(3,3) !! Matrix
1508 real(kind=xp) :: b(3,3) !! Inverse matrix
1509 real(kind=xp) :: detinv
1510
1511 ! Calculate the inverse determinant of the matrix
1512 ! first index x,y,z, second r, s, t
1513 detinv = 1.0_xp/real(a(1,1)*a(2,2)*a(3,3) - a(1,1)*a(2,3)*a(3,2)&
1514 - a(1,2)*a(2,1)*a(3,3) + a(1,2)*a(2,3)*a(3,1)&
1515 + a(1,3)*a(2,1)*a(3,2) - a(1,3)*a(2,2)*a(3,1),xp)
1516 ! Calculate the inverse of the matrix
1517 ! first index r, s, t, second x, y, z
1518 b(1,1) = +detinv * (a(2,2)*a(3,3) - a(2,3)*a(3,2))
1519 b(2,1) = -detinv * (a(2,1)*a(3,3) - a(2,3)*a(3,1))
1520 b(3,1) = +detinv * (a(2,1)*a(3,2) - a(2,2)*a(3,1))
1521 b(1,2) = -detinv * (a(1,2)*a(3,3) - a(1,3)*a(3,2))
1522 b(2,2) = +detinv * (a(1,1)*a(3,3) - a(1,3)*a(3,1))
1523 b(3,2) = -detinv * (a(1,1)*a(3,2) - a(1,2)*a(3,1))
1524 b(1,3) = +detinv * (a(1,2)*a(2,3) - a(1,3)*a(2,2))
1525 b(2,3) = -detinv * (a(1,1)*a(2,3) - a(1,3)*a(2,1))
1526 b(3,3) = +detinv * (a(1,1)*a(2,2) - a(1,2)*a(2,1))
1527 end function matinv3
1528
1529end module math
double real
Definition comm.F90:1
type(mpi_datatype), public mpi_real_precision
MPI type for working precision of REAL types.
Definition comm.F90:50
type(mpi_comm), public neko_comm
MPI communicator.
Definition comm.F90:42
type(mpi_datatype), public mpi_extra_precision
Definition comm.F90:51
Object for handling masks in Neko.
Definition mask.f90:34
Definition math.f90:60
subroutine, public cmult(a, c, n)
Multiplication by constant c .
Definition math.f90:411
subroutine, public cmult2(a, b, c, n)
Multiplication by constant c .
Definition math.f90:423
subroutine, public row_zero(a, m, n, e)
Sets row e to 0 in matrix a.
Definition math.f90:227
subroutine, public invcol2(a, b, n)
Vector division .
Definition math.f90:840
real(kind=rp) function, public vlsc2(u, v, n)
Compute multiplication sum .
Definition math.f90:712
pure logical function sabscmp(x, y, tol)
Return single precision absolute comparison .
Definition math.f90:117
real(kind=rp), parameter, public pi
Definition math.f90:75
pure logical function qabscmp(x, y, tol)
Return double precision absolute comparison .
Definition math.f90:147
real(kind=rp) function, public glsc3(a, b, c, n)
Weighted inner product .
Definition math.f90:1067
subroutine, public ascol5(a, b, c, d, e, n)
Returns .
Definition math.f90:1002
subroutine, public addcol3s2(a, b, c, s, n)
Returns .
Definition math.f90:987
subroutine, public masked_scatter_copy(a, b, mask, n, n_mask)
Scatter a contigous vector to masked positions in a target array .
Definition math.f90:381
subroutine, public invers2(a, b, n)
Compute inverted vector .
Definition math.f90:639
subroutine, public cadd2(a, b, s, n)
Add a scalar to vector .
Definition math.f90:474
real(rp) function, dimension(3, 3), public matinv39(a11, a12, a13, a21, a22, a23, a31, a32, a33)
Definition math.f90:1486
subroutine, public cadd(a, s, n)
Add a scalar to vector .
Definition math.f90:462
subroutine, public masked_copy(a, b, mask, n, n_mask)
Copy a masked vector .
Definition math.f90:289
subroutine reorddp(b, ind, n)
reorder double precision array - inverse of swap
Definition math.f90:1287
subroutine, public addsqr2s2(a, b, c1, n)
Returns .
Definition math.f90:826
real(kind=rp) function, public glsc4(a, b, c, d, n)
Definition math.f90:1086
subroutine, public cdiv2(a, b, c, n)
Division of constant c by elements of a .
Definition math.f90:449
subroutine swapdp(b, ind, n)
sort double precision array acording to ind vector
Definition math.f90:1247
subroutine flipvi4(b, ind, n)
Flip single integer vector b and ind.
Definition math.f90:1350
subroutine, public add2s1(a, b, c1, n)
Vector addition with scalar multiplication (multiplication on first argument)
Definition math.f90:797
subroutine, public cpwmin2(a, b, n)
Point-wise minimum of scalar and vector .
Definition math.f90:1458
real(kind=rp) function, public glsc2(a, b, n)
Weighted inner product .
Definition math.f90:1048
subroutine, public masked_scatter_copy_0(a, b, mask, n, n_mask)
Scatter a contigous vector to masked positions in a target array .
Definition math.f90:358
subroutine, public subcol3(a, b, c, n)
Returns .
Definition math.f90:881
subroutine, public rone(a, n)
Set all elements to one.
Definition math.f90:238
subroutine flipvdp(b, ind, n)
Flip double precision vector b and ind.
Definition math.f90:1327
subroutine, public cpwmin3(a, b, c, n)
Point-wise minimum of scalar and vector .
Definition math.f90:1470
subroutine, public pwmax3(a, b, c, n)
Point-wise maximum of two vectors .
Definition math.f90:1397
subroutine, public masked_gather_copy(a, b, mask, n, n_mask)
Gather a masked vector to reduced contigous vector .
Definition math.f90:335
subroutine, public x_update(a, b, c, c1, c2, n)
Returns .
Definition math.f90:1033
subroutine, public add3(a, b, c, n)
Vector addition .
Definition math.f90:739
subroutine swapi4(b, ind, n)
sort single integer array acording to ind vector
Definition math.f90:1267
integer function, public glimin(a, n)
Min of an integer vector of length n.
Definition math.f90:561
real(kind=rp) function, public glsum(a, n)
Sum a vector of length n.
Definition math.f90:499
subroutine, public sub3(a, b, c, n)
Vector subtraction .
Definition math.f90:781
subroutine, public addcol4(a, b, c, d, n)
Returns .
Definition math.f90:972
subroutine, public add2(a, b, n)
Vector addition .
Definition math.f90:726
subroutine, public cfill(a, c, n)
Set all elements to a constant c .
Definition math.f90:487
subroutine, public absval(a, n)
Take the absolute value of an array.
Definition math.f90:1372
subroutine, public invcol3(a, b, c, n)
Invert a vector .
Definition math.f90:626
subroutine, public add3s2(a, b, c, c1, c2, n)
Returns .
Definition math.f90:895
real(kind=xp) function, dimension(3, 3), public matinv3(a)
Performs a direct calculation of the inverse of a 3×3 matrix. M33INV and M44INV by David G....
Definition math.f90:1506
subroutine, public pwmax2(a, b, n)
Point-wise maximum of two vectors .
Definition math.f90:1385
subroutine, public pwmin2(a, b, n)
Point-wise minimum of two vectors .
Definition math.f90:1434
subroutine, public masked_gather_copy_0(a, b, mask, n, n_mask)
Gather a masked vector to reduced contigous vector .
Definition math.f90:312
subroutine, public subcol4(a, b, c, d, n)
Returns .
Definition math.f90:943
subroutine sorti4(a, ind, n)
Heap Sort for single integer arrays.
Definition math.f90:1191
subroutine, public addcol3(a, b, c, n)
Returns .
Definition math.f90:958
subroutine, public invcol1(a, n)
Invert a vector .
Definition math.f90:614
subroutine, public cdiv(a, c, n)
Division of constant c by elements of a .
Definition math.f90:437
real(kind=rp), parameter, public neko_m_ln2
Definition math.f90:72
subroutine, public chsign(a, n)
Change sign of vector .
Definition math.f90:579
subroutine, public cpwmax3(a, b, c, n)
Point-wise maximum of scalar and vector .
Definition math.f90:1421
subroutine, public col2(a, b, n)
Vector multiplication .
Definition math.f90:854
subroutine, public izero(a, n)
Zero an integer vector.
Definition math.f90:216
real(kind=rp) function, public glmax(a, n)
Max of a vector of length n.
Definition math.f90:516
subroutine, public copy(a, b, n)
Copy a vector .
Definition math.f90:249
subroutine, public add4s3(a, b, c, d, c1, c2, c3, n)
Returns .
Definition math.f90:910
subroutine, public add4(a, b, c, d, n)
Vector addition .
Definition math.f90:753
pure logical function dabscmp(x, y, tol)
Return double precision absolute comparison .
Definition math.f90:132
subroutine, public col3(a, b, c, n)
Vector multiplication with 3 vectors .
Definition math.f90:867
subroutine, public add5s4(a, b, c, d, e, c1, c2, c3, c4, n)
Returns .
Definition math.f90:926
pure logical function drelcmp(x, y, eps)
Return double precision relative comparison .
Definition math.f90:176
real(kind=rp), parameter, public neko_eps
Machine epsilon .
Definition math.f90:69
subroutine, public vdot3(dot, u1, u2, u3, v1, v2, v3, n)
Compute a dot product (3-d version) assuming vector components etc.
Definition math.f90:684
pure logical function qrelcmp(x, y, eps)
Return quad precision relative comparison .
Definition math.f90:191
subroutine, public rzero(a, n)
Zero a real vector.
Definition math.f90:205
real(kind=rp) function, public glsubnorm(a, b, n)
Returns the norm of the difference of two vectors .
Definition math.f90:1109
subroutine, public vdot2(dot, u1, u2, v1, v2, n)
Compute a dot product (2-d version) assuming vector components etc.
Definition math.f90:670
subroutine, public cpwmax2(a, b, n)
Point-wise maximum of scalar and vector .
Definition math.f90:1409
real(kind=rp) function, public vlmin(vec, n)
minimun value of a vector of length n
Definition math.f90:602
subroutine, public cfill_mask(a, c, n, mask, n_mask)
Fill a constant to a masked vector. .
Definition math.f90:397
real(kind=rp) function, public vlmax(vec, n)
maximum value of a vector of length n
Definition math.f90:591
integer function, public glimax(a, n)
Max of an integer vector of length n.
Definition math.f90:531
subroutine sortrp(a, ind, n)
Heap Sort for double precision arrays.
Definition math.f90:1133
subroutine, public sub2(a, b, n)
Vector substraction .
Definition math.f90:768
subroutine, public pwmin3(a, b, c, n)
Point-wise minimum of two vectors .
Definition math.f90:1446
subroutine, public add2s2(a, b, c1, n)
Vector addition with scalar multiplication (multiplication on second argument)
Definition math.f90:812
real(kind=rp) function, public glmin(a, n)
Min of a vector of length n.
Definition math.f90:546
subroutine, public masked_copy_0(a, b, mask, n, n_mask)
Copy a masked vector .
Definition math.f90:268
subroutine, public vcross(u1, u2, u3, v1, v2, v3, w1, w2, w3, n)
Compute a cross product assuming vector components etc.
Definition math.f90:653
pure logical function srelcmp(x, y, eps)
Return single precision relative comparison .
Definition math.f90:162
real(kind=rp) function, public vlsc3(u, v, w, n)
Compute multiplication sum .
Definition math.f90:698
subroutine reordi4(b, ind, n)
reorder single integer array - inverse of swap
Definition math.f90:1307
subroutine, public p_update(a, b, c, c1, c2, n)
Returns .
Definition math.f90:1018
integer, parameter, public qp
Definition num_types.f90:10
integer, parameter, public i4
Definition num_types.f90:6
integer, parameter, public xp
Definition num_types.f90:14
integer, parameter, public dp
Definition num_types.f90:9
integer, parameter, public sp
Definition num_types.f90:8
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
#define max(a, b)
Definition tensor.cu:40