Neko  0.9.99
A portable framework for high-order spectral element flow simulations
field_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 !
60 module field_math
62  use num_types, only: rp
63  use field, only: field_t
64  use math, only: rzero, rone, copy, cmult, cadd, cfill, invcol1, vdot3, add2, &
73  implicit none
74  private
75 
83 
84 contains
85 
87  subroutine field_rzero(a, n)
88  integer, intent(in), optional :: n
89  type(field_t), intent(inout) :: a
90  integer :: size
91 
92  if (present(n)) then
93  size = n
94  else
95  size = a%size()
96  end if
97 
98  if (neko_bcknd_device .eq. 1) then
99  call device_rzero(a%x_d, size)
100  else
101  call rzero(a%x, size)
102  end if
103  end subroutine field_rzero
104 
106  subroutine field_rone(a, n)
107  integer, intent(in), optional :: n
108  type(field_t), intent(inout) :: a
109  integer :: size
110 
111  if (present(n)) then
112  size = n
113  else
114  size = a%size()
115  end if
116 
117  if (neko_bcknd_device .eq. 1) then
118  call device_rone(a%x_d, size)
119  else
120  call rone(a%x, size)
121  end if
122  end subroutine field_rone
123 
125  subroutine field_copy(a, b, n)
126  integer, intent(in), optional :: n
127  type(field_t), intent(in) :: b
128  type(field_t), intent(inout) :: a
129  integer :: size
130 
131  if (present(n)) then
132  size = n
133  else
134  size = a%size()
135  end if
136 
137  if (neko_bcknd_device .eq. 1) then
138  call device_copy(a%x_d, b%x_d, size)
139  else
140  call copy(a%x, b%x, size)
141  end if
142  end subroutine field_copy
143 
145  subroutine field_cmult(a, c, n)
146  integer, intent(in), optional :: n
147  type(field_t), intent(inout) :: a
148  real(kind=rp), intent(in) :: c
149  integer :: size
150 
151  if (present(n)) then
152  size = n
153  else
154  size = a%size()
155  end if
156 
157  if (neko_bcknd_device .eq. 1) then
158  call device_cmult(a%x_d, c, size)
159  else
160  call cmult(a%x, c, size)
161  end if
162  end subroutine field_cmult
163 
165  subroutine field_cadd(a, s, n)
166  integer, intent(in), optional :: n
167  type(field_t), intent(inout) :: a
168  real(kind=rp), intent(in) :: s
169  integer :: size
170 
171  if (present(n)) then
172  size = n
173  else
174  size = a%size()
175  end if
176 
177  if (neko_bcknd_device .eq. 1) then
178  call device_cadd(a%x_d, s, size)
179  else
180  call cadd(a%x, s, size)
181  end if
182  end subroutine field_cadd
183 
185  subroutine field_cfill(a, c, n)
186  integer, intent(in), optional :: n
187  type(field_t), intent(inout) :: a
188  real(kind=rp), intent(in) :: c
189  integer :: size
190 
191  if (present(n)) then
192  size = n
193  else
194  size = a%size()
195  end if
196 
197  if (neko_bcknd_device .eq. 1) then
198  call device_cfill(a%x_d, c, size)
199  else
200  call cfill(a%x, c, size)
201  end if
202  end subroutine field_cfill
203 
205  subroutine field_invcol1(a, n)
206  integer, intent(in), optional :: n
207  type(field_t), intent(inout) :: a
208  integer :: size
209 
210  if (present(n)) then
211  size = n
212  else
213  size = a%size()
214  end if
215 
216  if (neko_bcknd_device .eq. 1) then
217  call device_invcol1(a%x_d, size)
218  else
219  call invcol1(a%x, size)
220  end if
221 
222  end subroutine field_invcol1
223 
226  subroutine field_vdot3(dot, u1, u2, u3, v1, v2, v3, n)
227  integer, intent(in), optional :: n
228  type(field_t), intent(in) :: u1, u2, u3
229  type(field_t), intent(in) :: v1, v2, v3
230  type(field_t), intent(out) :: dot
231  integer :: size
232 
233  if (present(n)) then
234  size = n
235  else
236  size = dot%size()
237  end if
238 
239  if (neko_bcknd_device .eq. 1) then
240  call device_vdot3(dot%x_d, &
241  u1%x_d, u2%x_d, u3%x_d, &
242  v1%x_d, v2%x_d, v3%x_d, &
243  size)
244  else
245  call vdot3(dot%x, &
246  u1%x, u2%x, u3%x, &
247  v1%x, v2%x, v3%x, &
248  size)
249  end if
250 
251  end subroutine field_vdot3
252 
254  subroutine field_add2(a, b, n)
255  integer, intent(in), optional :: n
256  type(field_t), intent(inout) :: a
257  type(field_t), intent(in) :: b
258  integer :: size
259 
260  if (present(n)) then
261  size = n
262  else
263  size = a%size()
264  end if
265 
266  if (neko_bcknd_device .eq. 1) then
267  call device_add2(a%x_d, b%x_d, size)
268  else
269  call add2(a%x, b%x, size)
270  end if
271 
272  end subroutine field_add2
273 
275  subroutine field_add3(a, b, c, n)
276  integer, intent(in), optional :: n
277  type(field_t), intent(inout) :: a
278  type(field_t), intent(in) :: b, c
279  integer :: size
280 
281  if (present(n)) then
282  size = n
283  else
284  size = a%size()
285  end if
286 
287  if (neko_bcknd_device .eq. 1) then
288  call device_add3(a%x_d, b%x_d, c%x_d, size)
289  else
290  call add3(a%x, b%x, c%x, size)
291  end if
292 
293  end subroutine field_add3
294 
296  subroutine field_add4(a, b, c, d, n)
297  integer, intent(in), optional :: n
298  type(field_t), intent(inout) :: a
299  type(field_t), intent(in) :: b, c, d
300  integer :: size
301 
302  if (present(n)) then
303  size = n
304  else
305  size = a%size()
306  end if
307 
308  if (neko_bcknd_device .eq. 1) then
309  call device_add4(a%x_d, b%x_d, c%x_d, d%x_d, size)
310  else
311  call add4(a%x, b%x, c%x, d%x, size)
312  end if
313 
314  end subroutine field_add4
315 
317  subroutine field_sub2(a, b, n)
318  integer, intent(in), optional :: n
319  type(field_t), intent(inout) :: a
320  type(field_t), intent(inout) :: b
321  integer :: size
322 
323  if (present(n)) then
324  size = n
325  else
326  size = a%size()
327  end if
328 
329  if (neko_bcknd_device .eq. 1) then
330  call device_sub2(a%x_d, b%x_d, size)
331  else
332  call sub2(a%x, b%x, size)
333  end if
334 
335  end subroutine field_sub2
336 
338  subroutine field_sub3(a, b, c, n)
339  integer, intent(in), optional :: n
340  type(field_t), intent(inout) :: c
341  type(field_t), intent(inout) :: b
342  type(field_t), intent(out) :: a
343  integer :: size
344 
345  if (present(n)) then
346  size = n
347  else
348  size = a%size()
349  end if
350 
351  if (neko_bcknd_device .eq. 1) then
352  call device_sub3(a%x_d, b%x_d, c%x_d, size)
353  else
354  call sub3(a%x, b%x, c%x, size)
355  end if
356 
357  end subroutine field_sub3
358 
359 
362  subroutine field_add2s1(a, b, c1, n)
363  integer, intent(in), optional :: n
364  type(field_t), intent(inout) :: a
365  type(field_t), intent(inout) :: b
366  real(kind=rp), intent(in) :: c1
367  integer :: size
368 
369  if (present(n)) then
370  size = n
371  else
372  size = a%size()
373  end if
374 
375  if (neko_bcknd_device .eq. 1) then
376  call device_add2s1(a%x_d, b%x_d, c1, size)
377  else
378  call add2s1(a%x, b%x, c1, size)
379  end if
380 
381  end subroutine field_add2s1
382 
385  subroutine field_add2s2(a, b, c1, n)
386  integer, intent(in), optional :: n
387  type(field_t), intent(inout) :: a
388  type(field_t), intent(inout) :: b
389  real(kind=rp), intent(in) :: c1
390  integer :: size
391 
392  if (present(n)) then
393  size = n
394  else
395  size = a%size()
396  end if
397 
398  if (neko_bcknd_device .eq. 1) then
399  call device_add2s2(a%x_d, b%x_d, c1, size)
400  else
401  call add2s2(a%x, b%x, c1, size)
402  end if
403 
404  end subroutine field_add2s2
405 
407  subroutine field_addsqr2s2(a, b, c1, n)
408  integer, intent(in), optional :: n
409  type(field_t), intent(inout) :: a
410  type(field_t), intent(in) :: b
411  real(kind=rp), intent(in) :: c1
412  integer :: size
413 
414  if (present(n)) then
415  size = n
416  else
417  size = a%size()
418  end if
419 
420  if (neko_bcknd_device .eq. 1) then
421  call device_addsqr2s2(a%x_d, b%x_d, c1, size)
422  else
423  call addsqr2s2(a%x, b%x, c1, size)
424  end if
425 
426  end subroutine field_addsqr2s2
427 
429  subroutine field_cmult2(a, b, c, n)
430  integer, intent(in), optional :: n
431  type(field_t), intent(inout) :: a
432  type(field_t), intent(in) :: b
433  real(kind=rp), intent(in) :: c
434  integer :: size
435 
436  if (present(n)) then
437  size = n
438  else
439  size = a%size()
440  end if
441 
442  if (neko_bcknd_device .eq. 1) then
443  call device_cmult2(a%x_d, b%x_d, c, size)
444  else
445  call cmult2(a%x, b%x, c, size)
446  end if
447 
448  end subroutine field_cmult2
449 
451  subroutine field_invcol2(a, b, n)
452  integer, intent(in), optional :: n
453  type(field_t), intent(inout) :: a
454  type(field_t), intent(in) :: b
455  integer :: size
456 
457  if (present(n)) then
458  size = n
459  else
460  size = a%size()
461  end if
462 
463  if (neko_bcknd_device .eq. 1) then
464  call device_invcol2(a%x_d, b%x_d, size)
465  else
466  call invcol2(a%x, b%x, size)
467  end if
468 
469  end subroutine field_invcol2
470 
471 
473  subroutine field_col2(a, b, n)
474  integer, intent(in), optional :: n
475  type(field_t), intent(inout) :: a
476  type(field_t), intent(in) :: b
477  integer :: size
478 
479  if (present(n)) then
480  size = n
481  else
482  size = a%size()
483  end if
484 
485  if (neko_bcknd_device .eq. 1) then
486  call device_col2(a%x_d, b%x_d, size)
487  else
488  call col2(a%x, b%x, size)
489  end if
490 
491  end subroutine field_col2
492 
494  subroutine field_col3(a, b, c, n)
495  integer, intent(in), optional :: n
496  type(field_t), intent(inout) :: a
497  type(field_t), intent(in) :: b
498  type(field_t), intent(in) :: c
499  integer :: size
500 
501  if (present(n)) then
502  size = n
503  else
504  size = a%size()
505  end if
506 
507  if (neko_bcknd_device .eq. 1) then
508  call device_col3(a%x_d, b%x_d, c%x_d, size)
509  else
510  call col3(a%x, b%x, c%x, size)
511  end if
512 
513  end subroutine field_col3
514 
516  subroutine field_subcol3(a, b, c, n)
517  integer, intent(in), optional :: n
518  type(field_t), intent(inout) :: a
519  type(field_t), intent(in) :: b
520  type(field_t), intent(in) :: c
521  integer :: size
522 
523  if (present(n)) then
524  size = n
525  else
526  size = a%size()
527  end if
528 
529  if (neko_bcknd_device .eq. 1) then
530  call device_subcol3(a%x_d, b%x_d, c%x_d, size)
531  else
532  call subcol3(a%x, b%x, c%x, size)
533  end if
534 
535  end subroutine field_subcol3
536 
538  subroutine field_add3s2(a, b, c, c1, c2, n)
539  integer, intent(in), optional :: n
540  type(field_t), intent(inout) :: a
541  type(field_t), intent(in) :: b
542  type(field_t), intent(in) :: c
543  real(kind=rp), intent(in) :: c1, c2
544  integer :: size
545 
546  if (present(n)) then
547  size = n
548  else
549  size = a%size()
550  end if
551 
552  if (neko_bcknd_device .eq. 1) then
553  call device_add3s2(a%x_d, b%x_d, c%x_d, c1, c2, size)
554  else
555  call add3s2(a%x, b%x, c%x, c1, c2, size)
556  end if
557 
558  end subroutine field_add3s2
559 
561  subroutine field_addcol3(a, b, c, n)
562  integer, intent(in), optional :: n
563  type(field_t), intent(inout) :: a
564  type(field_t), intent(in) :: b
565  type(field_t), intent(in) :: c
566  integer :: size
567 
568  if (present(n)) then
569  size = n
570  else
571  size = a%size()
572  end if
573 
574  if (neko_bcknd_device .eq. 1) then
575  call device_addcol3(a%x_d, b%x_d, c%x_d, size)
576  else
577  call addcol3(a%x, b%x, c%x, size)
578  end if
579 
580  end subroutine field_addcol3
581 
583  subroutine field_addcol4(a, b, c, d, n)
584  integer, intent(in), optional :: n
585  type(field_t), intent(inout) :: a
586  type(field_t), intent(in) :: b
587  type(field_t), intent(in) :: c
588  type(field_t), intent(in) :: d
589  integer :: size
590 
591  if (present(n)) then
592  size = n
593  else
594  size = a%size()
595  end if
596 
597  if (neko_bcknd_device .eq. 1) then
598  call device_addcol4(a%x_d, b%x_d, c%x_d, d%x_d, size)
599  else
600  call addcol4(a%x, b%x, c%x, d%x, size)
601  end if
602 
603  end subroutine field_addcol4
604 
605  function field_glsum(a, n) result(sum)
606  integer, intent(in), optional :: n
607  type(field_t), intent(in) :: a
608  real(kind=rp) :: sum
609  integer :: size
610 
611  if (present(n)) then
612  size = n
613  else
614  size = a%size()
615  end if
616 
617  if (neko_bcknd_device .eq. 1) then
618  sum = device_glsum(a%x_d, size)
619  else
620  sum = glsum(a%x, size)
621  end if
622 
623  end function field_glsum
624 
625  function field_glsc2(a, b, n) result(norm)
626  integer, intent(in), optional :: n
627  type(field_t), intent(in) :: a, b
628  real(kind=rp) :: norm
629  integer :: size
630 
631  if (present(n)) then
632  size = n
633  else
634  size = a%size()
635  end if
636 
637  if (neko_bcknd_device .eq. 1) then
638  norm = device_glsc2(a%x_d, b%x_d, size)
639  else
640  norm = glsc2(a%x, b%x, size)
641  end if
642 
643  end function field_glsc2
644 
645  function field_glsc3(a, b, c, n) result(norm)
646  integer, intent(in), optional :: n
647  type(field_t), intent(in) :: a, b, c
648  real(kind=rp) :: norm
649  integer :: size
650 
651  if (present(n)) then
652  size = n
653  else
654  size = a%size()
655  end if
656 
657  if (neko_bcknd_device .eq. 1) then
658  norm = device_glsc3(a%x_d, b%x_d, c%x_d, size)
659  else
660  norm = glsc3(a%x, b%x, c%x, size)
661  end if
662 
663  end function field_glsc3
664 
665 end module field_math
subroutine, public device_add2(a_d, b_d, n)
Vector addition .
subroutine, public device_addcol3(a_d, b_d, c_d, n)
Returns .
subroutine, public device_col2(a_d, b_d, n)
Vector multiplication .
subroutine, public device_add2s1(a_d, b_d, c1, n)
subroutine, public device_rzero(a_d, n)
Zero a real vector.
subroutine, public device_rone(a_d, n)
Set all elements to one.
subroutine, public device_add2s2(a_d, b_d, c1, n)
Vector addition with scalar multiplication (multiplication on first argument)
subroutine, public device_invcol1(a_d, n)
Invert a vector .
subroutine, public device_col3(a_d, b_d, c_d, n)
Vector multiplication with 3 vectors .
subroutine, public device_add4(a_d, b_d, c_d, d_d, n)
subroutine, public device_cadd(a_d, c, n)
Add a scalar to vector .
subroutine, public device_vdot3(dot_d, u1_d, u2_d, u3_d, v1_d, v2_d, v3_d, n)
Compute a dot product (3-d version) assuming vector components etc.
subroutine, public device_cmult2(a_d, b_d, c, n)
Multiplication by constant c .
subroutine, public device_cmult(a_d, c, n)
Multiplication by constant c .
real(kind=rp) function, public device_glsc2(a_d, b_d, n)
Weighted inner product .
subroutine, public device_sub3(a_d, b_d, c_d, n)
Vector subtraction .
real(kind=rp) function, public device_glsc3(a_d, b_d, c_d, n)
Weighted inner product .
subroutine, public device_add3(a_d, b_d, c_d, n)
Vector addition .
real(kind=rp) function, public device_glsum(a_d, n)
Sum a vector of length n.
subroutine, public device_copy(a_d, b_d, n)
Copy a vector .
Definition: device_math.F90:76
subroutine, public device_add3s2(a_d, b_d, c_d, c1, c2, n)
Returns .
subroutine, public device_subcol3(a_d, b_d, c_d, n)
Returns .
subroutine, public device_sub2(a_d, b_d, n)
Vector substraction .
subroutine, public device_cfill(a_d, c, n)
Set all elements to a constant c .
subroutine, public device_addcol4(a_d, b_d, c_d, d_d, n)
Returns .
subroutine, public device_invcol2(a_d, b_d, n)
Vector division .
subroutine, public device_addsqr2s2(a_d, b_d, c1, n)
Returns .
subroutine, public field_add2s1(a, b, c1, n)
Vector addition with scalar multiplication (multiplication on first argument)
Definition: field_math.f90:363
subroutine, public field_invcol2(a, b, n)
Vector division .
Definition: field_math.f90:452
subroutine field_add4(a, b, c, d, n)
Vector addition .
Definition: field_math.f90:297
subroutine field_add3(a, b, c, n)
Vector addition .
Definition: field_math.f90:276
subroutine, public field_cadd(a, s, n)
Add a scalar to vector .
Definition: field_math.f90:166
subroutine, public field_invcol1(a, n)
Invert a vector .
Definition: field_math.f90:206
subroutine, public field_col2(a, b, n)
Vector multiplication .
Definition: field_math.f90:474
subroutine, public field_sub2(a, b, n)
Vector substraction .
Definition: field_math.f90:318
subroutine, public field_sub3(a, b, c, n)
Vector subtraction .
Definition: field_math.f90:339
subroutine, public field_cmult2(a, b, c, n)
Multiplication by constant c .
Definition: field_math.f90:430
real(kind=rp) function, public field_glsc2(a, b, n)
Definition: field_math.f90:626
subroutine, public field_cfill(a, c, n)
Set all elements to a constant c .
Definition: field_math.f90:186
subroutine, public field_add2s2(a, b, c1, n)
Vector addition with scalar multiplication (multiplication on second argument)
Definition: field_math.f90:386
real(kind=rp) function, public field_glsum(a, n)
Definition: field_math.f90:606
subroutine, public field_rzero(a, n)
Zero a real vector.
Definition: field_math.f90:88
subroutine, public field_addcol3(a, b, c, n)
Returns .
Definition: field_math.f90:562
subroutine, public field_rone(a, n)
Set all elements to one.
Definition: field_math.f90:107
subroutine, public field_add3s2(a, b, c, c1, c2, n)
Returns .
Definition: field_math.f90:539
subroutine, public field_add2(a, b, n)
Vector addition .
Definition: field_math.f90:255
subroutine, public field_col3(a, b, c, n)
Vector multiplication with 3 vectors .
Definition: field_math.f90:495
subroutine, public field_copy(a, b, n)
Copy a vector .
Definition: field_math.f90:126
subroutine, public field_subcol3(a, b, c, n)
Returns .
Definition: field_math.f90:517
subroutine, public field_cmult(a, c, n)
Multiplication by constant c .
Definition: field_math.f90:146
real(kind=rp) function, public field_glsc3(a, b, c, n)
Definition: field_math.f90:646
subroutine, public field_addcol4(a, b, c, d, n)
Returns .
Definition: field_math.f90:584
subroutine, public field_addsqr2s2(a, b, c1, n)
Returns .
Definition: field_math.f90:408
subroutine, public field_vdot3(dot, u1, u2, u3, v1, v2, v3, n)
Compute a dot product (3-d version) assuming vector components etc.
Definition: field_math.f90:227
Defines a field.
Definition: field.f90:34
Definition: math.f90:60
subroutine, public cmult(a, c, n)
Multiplication by constant c .
Definition: math.f90:311
subroutine, public cmult2(a, b, c, n)
Multiplication by constant c .
Definition: math.f90:701
subroutine, public invcol2(a, b, n)
Vector division .
Definition: math.f90:715
real(kind=rp) function, public glsc3(a, b, c, n)
Weighted inner product .
Definition: math.f90:895
subroutine, public cadd(a, s, n)
Add a scalar to vector .
Definition: math.f90:323
subroutine, public addsqr2s2(a, b, c1, n)
Returns .
Definition: math.f90:687
subroutine, public add2s1(a, b, c1, n)
Vector addition with scalar multiplication (multiplication on first argument)
Definition: math.f90:658
real(kind=rp) function, public glsc2(a, b, n)
Weighted inner product .
Definition: math.f90:876
subroutine, public subcol3(a, b, c, n)
Returns .
Definition: math.f90:756
subroutine, public rone(a, n)
Set all elements to one.
Definition: math.f90:228
subroutine, public add3(a, b, c, n)
Vector addition .
Definition: math.f90:600
real(kind=rp) function, public glsum(a, n)
Sum a vector of length n.
Definition: math.f90:360
subroutine, public sub3(a, b, c, n)
Vector subtraction .
Definition: math.f90:642
subroutine, public addcol4(a, b, c, d, n)
Returns .
Definition: math.f90:815
subroutine, public add2(a, b, n)
Vector addition .
Definition: math.f90:587
subroutine, public cfill(a, c, n)
Set all elements to a constant c .
Definition: math.f90:348
subroutine, public add3s2(a, b, c, c1, c2, n)
Returns .
Definition: math.f90:770
subroutine, public addcol3(a, b, c, n)
Returns .
Definition: math.f90:801
subroutine, public invcol1(a, n)
Invert a vector .
Definition: math.f90:475
subroutine, public col2(a, b, n)
Vector multiplication .
Definition: math.f90:729
subroutine, public copy(a, b, n)
Copy a vector .
Definition: math.f90:239
subroutine, public add4(a, b, c, d, n)
Vector addition .
Definition: math.f90:614
subroutine, public col3(a, b, c, n)
Vector multiplication with 3 vectors .
Definition: math.f90:742
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:545
subroutine, public rzero(a, n)
Zero a real vector.
Definition: math.f90:195
subroutine, public sub2(a, b, n)
Vector substraction .
Definition: math.f90:629
subroutine, public add2s2(a, b, c1, n)
Vector addition with scalar multiplication (multiplication on second argument)
Definition: math.f90:673
Build configurations.
Definition: neko_config.f90:34
integer, parameter neko_bcknd_device
Definition: neko_config.f90:44
integer, parameter, public rp
Global precision used in computations.
Definition: num_types.f90:12