Neko  0.9.0
A portable framework for high-order spectral element flow simulations
scalar_residual_device.F90
Go to the documentation of this file.
1 ! Copyright (c) 2022-2024, The Neko Authors
2 ! All rights reserved.
3 !
4 ! Redistribution and use in source and binary forms, with or without
5 ! modification, are permitted provided that the following conditions
6 ! are met:
7 !
8 ! * Redistributions of source code must retain the above copyright
9 ! notice, this list of conditions and the following disclaimer.
10 !
11 ! * Redistributions in binary form must reproduce the above
12 ! copyright notice, this list of conditions and the following
13 ! disclaimer in the documentation and/or other materials provided
14 ! with the distribution.
15 !
16 ! * Neither the name of the authors nor the names of its
17 ! contributors may be used to endorse or promote products derived
18 ! from this software without specific prior written permission.
19 !
20 ! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21 ! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22 ! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
23 ! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
24 ! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
25 ! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
26 ! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
27 ! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28 ! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29 ! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30 ! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31 ! POSSIBILITY OF SUCH DAMAGE.
32 !
36  use ax_product, only : ax_t
37  use field, only : field_t
38  use coefs, only : coef_t
39  use space, only : space_t
40  use mesh, only : mesh_t
41  use num_types, only : rp
42  use, intrinsic :: iso_c_binding
43  implicit none
44  private
45 
46  type, public, extends(scalar_residual_t) :: scalar_residual_device_t
47  contains
48  procedure, nopass :: compute => scalar_residual_device_compute
50 
51 #ifdef HAVE_HIP
52  interface
53  subroutine scalar_residual_update_hip(s_res_d, f_s_d, n) &
54  bind(c, name = 'scalar_residual_update_hip')
55  use, intrinsic :: iso_c_binding
56  implicit none
57  type(c_ptr), value :: s_res_d
58  type(c_ptr), value :: f_s_d
59  integer(c_int) :: n
60  end subroutine scalar_residual_update_hip
61  end interface
62 #elif HAVE_CUDA
63 
64  interface
65  subroutine scalar_residual_update_cuda(s_res_d, f_s_d, n) &
66  bind(c, name = 'scalar_residual_update_cuda')
67  use, intrinsic :: iso_c_binding
68  implicit none
69  type(c_ptr), value :: s_res_d
70  type(c_ptr), value :: f_s_d
71  integer(c_int) :: n
72  end subroutine scalar_residual_update_cuda
73  end interface
74 #elif HAVE_OPENCL
75 
76  interface
77  subroutine scalar_residual_update_opencl(s_res_d, f_s_d, n) &
78  bind(c, name = 'scalar_residual_update_opencl')
79  use, intrinsic :: iso_c_binding
80  implicit none
81  type(c_ptr), value :: s_res_d
82  type(c_ptr), value :: f_s_d
83  integer(c_int) :: n
84  end subroutine scalar_residual_update_opencl
85  end interface
86 #endif
87 
88 
89 contains
90 
91 
92  subroutine scalar_residual_device_compute(Ax, s, s_res, f_Xh, c_Xh, msh, Xh, &
93  lambda, rhocp, bd, dt, n)
94  class(ax_t), intent(in) :: Ax
95  type(mesh_t), intent(inout) :: msh
96  type(space_t), intent(inout) :: Xh
97  type(field_t), intent(inout) :: s
98  type(field_t), intent(inout) :: s_res
99  type(field_t), intent(inout) :: f_Xh
100  type(coef_t), intent(inout) :: c_Xh
101  type(field_t), intent(in) :: lambda
102  real(kind=rp), intent(in) :: rhocp
103  real(kind=rp), intent(in) :: bd
104  real(kind=rp), intent(in) :: dt
105  integer, intent(in) :: n
106 
107  call device_copy(c_xh%h1_d, lambda%x_d, n)
108  call device_cfill(c_xh%h2_d, rhocp * (bd / dt), n)
109  c_xh%ifh2 = .true.
110 
111  call ax%compute(s_res%x, s%x, c_xh, msh, xh)
112 
113 #ifdef HAVE_HIP
114  call scalar_residual_update_hip(s_res%x_d, f_xh%x_d, n)
115 #elif HAVE_CUDA
116  call scalar_residual_update_cuda(s_res%x_d, f_xh%x_d, n)
117 #elif HAVE_OPENCL
118  call scalar_residual_update_opencl(s_res%x_d, f_xh%x_d, n)
119 #endif
120 
121  end subroutine scalar_residual_device_compute
122 
123 end module scalar_residual_device
Defines a Matrix-vector product.
Definition: ax.f90:34
Coefficients.
Definition: coef.f90:34
subroutine, public device_copy(a_d, b_d, n)
Copy a vector .
Definition: device_math.F90:76
subroutine, public device_cfill(a_d, c, n)
Set all elements to a constant c .
Defines a field.
Definition: field.f90:34
Defines a mesh.
Definition: mesh.f90:34
integer, parameter, public rp
Global precision used in computations.
Definition: num_types.f90:12
subroutine scalar_residual_device_compute(Ax, s, s_res, f_Xh, c_Xh, msh, Xh, lambda, rhocp, bd, dt, n)
Defines the residual for the scalar transport equation.
Defines a function space.
Definition: space.f90:34
void scalar_residual_update_opencl(void *s_res, void *f_s, int *n)
void scalar_residual_update_cuda(void *s_res, void *f_s, int *n)
Base type for a matrix-vector product providing .
Definition: ax.f90:43
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition: coef.f90:55
Abstract type to compute scalar residual.
The function space for the SEM solution fields.
Definition: space.f90:62