Neko  0.9.99
A portable framework for high-order spectral element flow simulations
wall_model_bc.f90
Go to the documentation of this file.
1 ! Copyright (c) 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 num_types, only : rp
37  use bc, only : bc_t
38  use, intrinsic :: iso_c_binding, only : c_ptr
39  use utils, only : neko_error, nonlinear_index
40  use coefs, only : coef_t
41  use wall_model, only : wall_model_t, wall_model_factory
42  use rough_log_law, only : rough_log_law_t
43  use spalding, only : spalding_t
44  use shear_stress, only : shear_stress_t
45  use json_module, only : json_file
46  implicit none
47  private
48 
52  type, public, extends(shear_stress_t) :: wall_model_bc_t
54  class(wall_model_t), allocatable :: wall_model
55  contains
56  procedure, pass(this) :: apply_scalar => wall_model_bc_apply_scalar
57  procedure, pass(this) :: apply_vector => wall_model_bc_apply_vector
58  procedure, pass(this) :: apply_scalar_dev => &
60  procedure, pass(this) :: apply_vector_dev => &
62  procedure, pass(this) :: init_wall_model_bc => &
64  end type wall_model_bc_t
65 
66  contains
67 
69  subroutine wall_model_bc_apply_scalar(this, x, n, t, tstep)
70  class(wall_model_bc_t), intent(inout) :: this
71  integer, intent(in) :: n
72  real(kind=rp), intent(inout), dimension(n) :: x
73  real(kind=rp), intent(in), optional :: t
74  integer, intent(in), optional :: tstep
75 
76  call neko_error("The wall model bc is not applicable to scalar fields.")
77 
78  end subroutine wall_model_bc_apply_scalar
79 
87  subroutine wall_model_bc_apply_vector(this, x, y, z, n, t, tstep)
88  class(wall_model_bc_t), intent(inout) :: this
89  integer, intent(in) :: n
90  real(kind=rp), intent(inout), dimension(n) :: x
91  real(kind=rp), intent(inout), dimension(n) :: y
92  real(kind=rp), intent(inout), dimension(n) :: z
93  real(kind=rp), intent(in), optional :: t
94  integer, intent(in), optional :: tstep
95  integer :: i, m, k, fid
96  real(kind=rp) :: magtau
97 
98  ! Compute the wall stress using the wall model.
99  call this%wall_model%compute(t, tstep)
100 
101  ! Populate the 3D wall stress field for post-processing.
102  do i = 1, this%msk(0)
103  magtau = sqrt(this%wall_model%tau_x(i)**2 + this%wall_model%tau_y(i)**2&
104  + this%wall_model%tau_z(i)**2)
105 
106  ! Mark sampling nodes with a -1 for debugging
107  this%wall_model%tau_field%x(this%wall_model%ind_r(i), &
108  this%wall_model%ind_s(i), &
109  this%wall_model%ind_t(i), &
110  this%wall_model%ind_e(i)) = -1.0_rp
111  this%wall_model%tau_field%x(this%msk(i),1,1,1) = magtau
112  end do
113 
114  ! Set the computed stress for application by the underlying Neumann
115  ! boundary conditions.
116  call this%set_stress(this%wall_model%tau_x, this%wall_model%tau_y, &
117  this%wall_model%tau_z)
118 
119  ! Add the stress as a forcing to the right hand side arrays
120  call this%shear_stress_t%apply_vector(x, y, z, n, t, tstep)
121 
122  end subroutine wall_model_bc_apply_vector
123 
126  subroutine wall_model_bc_apply_scalar_dev(this, x_d, t, tstep)
127  class(wall_model_bc_t), intent(inout), target :: this
128  type(c_ptr) :: x_d
129  real(kind=rp), intent(in), optional :: t
130  integer, intent(in), optional :: tstep
131 
132  call neko_error("wall_model_bc bc not implemented on the device")
133 
134  end subroutine wall_model_bc_apply_scalar_dev
135 
138  subroutine wall_model_bc_apply_vector_dev(this, x_d, y_d, z_d, t, tstep)
139  class(wall_model_bc_t), intent(inout), target :: this
140  type(c_ptr) :: x_d
141  type(c_ptr) :: y_d
142  type(c_ptr) :: z_d
143  real(kind=rp), intent(in), optional :: t
144  integer, intent(in), optional :: tstep
145 
146  call neko_error("wall_model_bc bc not implemented on the device")
147 
148  end subroutine wall_model_bc_apply_vector_dev
149 
152  subroutine wall_model_bc_init_wall_model_bc(this, json, nu)
153  class(wall_model_bc_t), intent(inout) :: this
154  type(json_file), intent(inout) :: json
155  real(kind=rp), intent(in) :: nu
156 
157  call this%shear_stress_t%init_shear_stress(this%coef)
158 
159  call wall_model_factory(this%wall_model, this%coef, this%msk, &
160  this%facet, nu, json)
161 
162  end subroutine wall_model_bc_init_wall_model_bc
163 
164  end module wall_model_bc
__device__ void nonlinear_index(const int idx, const int lx, int *index)
Defines a boundary condition.
Definition: bc.f90:34
Coefficients.
Definition: coef.f90:34
integer, parameter, public rp
Global precision used in computations.
Definition: num_types.f90:12
Implements rough_log_law_t.
Defines a shear stress boundary condition for a vector field. Maintainer: Timofey Mukha.
Implements spalding_t.
Definition: spalding.f90:35
Utilities.
Definition: utils.f90:35
Defines the wall_model_bc_t type. Maintainer: Timofey Mukha.
subroutine wall_model_bc_apply_scalar(this, x, n, t, tstep)
Apply shear stress for a scalar field x.
subroutine wall_model_bc_init_wall_model_bc(this, json, nu)
Constructor.
subroutine wall_model_bc_apply_scalar_dev(this, x_d, t, tstep)
Boundary condition apply for a generic wall_model_bc condition to a vector x (device version)
subroutine wall_model_bc_apply_vector(this, x, y, z, n, t, tstep)
Apply the boundary condition to the right-hand side.
subroutine wall_model_bc_apply_vector_dev(this, x_d, y_d, z_d, t, tstep)
Boundary condition apply for a generic wall_model_bc condition to vectors x, y and z (device version)
Implements wall_model_t.
Definition: wall_model.f90:35
Base type for a boundary condition.
Definition: bc.f90:51
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition: coef.f90:55
Wall model based on the log-law for a rough wall. The formula defining the law is ....
A shear stress boundary condition.
Wall model based on Spalding's law of the wall. Reference: http://dx.doi.org/10.1115/1....
Definition: spalding.f90:53
Base abstract type for wall-stress models for wall-modelled LES.
Definition: wall_model.f90:54
A shear stress boundary condition, computing the stress values using a wall model.