Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
time_interpolator.f90
Go to the documentation of this file.
1! Copyright (c) 2022-2025, 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!
35 use num_types, only : rp
36 use field, only : field_t
40 use math, only : add3s2, add4s3
41 use utils, only : neko_error
42 use, intrinsic :: iso_c_binding
43 use fast3d, only : fd_weights_full
44 implicit none
45 private
46
48 type, public :: time_interpolator_t
50 integer :: order
51 contains
53 procedure, pass(this) :: init => time_interpolator_init
55 procedure, pass(this) :: free => time_interpolator_free
57 procedure, pass(this) :: interpolate => time_interpolator_interpolate
59 procedure, pass(this) :: interpolate_scalar => time_interpolator_scalar
60
61 end type time_interpolator_t
62
63contains
64
67 subroutine time_interpolator_init(this, order)
68 class(time_interpolator_t), intent(inout) :: this
69 integer, intent(in), target :: order
70
72 call this%free()
73
75 this%order = order
76
77 end subroutine time_interpolator_init
78
79
81 subroutine time_interpolator_free(this)
82 class(time_interpolator_t), intent(inout) :: this
83
84 end subroutine time_interpolator_free
85
93 subroutine time_interpolator_interpolate(this, t, f, t_past, f_past, &
94 t_future, f_future)
95 class(time_interpolator_t), intent(inout) :: this
96 real(kind=rp), intent(inout) :: t, t_past, t_future
97 type(field_t), intent(inout) :: f, f_past, f_future
98 real(kind=rp) :: w_past, w_future !Weights for the interpolation
99 integer :: n
100
101 if (this%order .eq. 2) then
102
103 n = f%dof%size()
104 w_past = ( t_future - t ) / ( t_future - t_past )
105 w_future = ( t - t_past ) / ( t_future - t_past )
106
107 if (neko_bcknd_device .eq. 1) then
108 call device_add3s2(f%x_d, f_past%x_d, f_future%x_d, &
109 w_past, w_future, n)
110 else
111 call add3s2(f%x, f_past%x, f_future%x, w_past, w_future, n)
112 end if
113
114 else
115 call neko_error("Time interpolation of required order &
116 &is not implemented")
117 end if
118
119 end subroutine time_interpolator_interpolate
120
128 subroutine time_interpolator_scalar(this, t, f_interpolated, f_n, tlag, n)
129 class(time_interpolator_t), intent(inout) :: this
130 real(kind=rp), intent(in) :: t
131 integer, intent(in) :: n
132 type(field_series_t), intent(inout) :: f_n
133 type(field_t), intent(inout) :: f_interpolated
134 real(kind=rp), dimension(0:this%order), intent(in) :: tlag
135 type(c_ptr) :: f_interpolated_d
136 integer :: i, l
137
138 real(kind=rp), dimension(0:this%order - 1) :: wt
139 wt = 0
140
141 call fd_weights_full(t, tlag, this%order - 1, 0, wt) ! interpolation weights
142
143 if (neko_bcknd_device .eq. 1) then
144 call device_add4s3(f_interpolated%x_d, f_n%f%x_d, f_n%lf(1)%x_d, &
145 f_n%lf(2)%x_d, wt(0), wt(1), wt(2), n)
146 else
147 call add4s3(f_interpolated%x, f_n%f%x, f_n%lf(1)%x, f_n%lf(2)%x, &
148 wt(0), wt(1), wt(2), n)
149 end if
150
151 end subroutine time_interpolator_scalar
152
153end module time_interpolator
subroutine, public device_add4s3(a_d, b_d, c_d, d_d, c1, c2, c3, n, strm)
Returns .
subroutine, public device_add3s2(a_d, b_d, c_d, c1, c2, n, strm)
Returns .
Fast diagonalization methods from NEKTON.
Definition fast3d.f90:61
subroutine, public fd_weights_full(xi, x, n, m, c)
Compute finite-difference stencil weights for evaluating derivatives up to order at a point.
Definition fast3d.f90:105
Contains the field_serties_t type.
Defines a field.
Definition field.f90:34
Definition math.f90:60
subroutine, public add3s2(a, b, c, c1, c2, n)
Returns .
Definition math.f90:895
subroutine, public add4s3(a, b, c, d, c1, c2, c3, n)
Returns .
Definition math.f90:910
Build configurations.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Implements type time_interpolator_t.
subroutine time_interpolator_interpolate(this, t, f, t_past, f_past, t_future, f_future)
Interpolate a field at time t from fields at time t-dt and t+dt.
subroutine time_interpolator_free(this)
Destructor.
subroutine time_interpolator_init(this, order)
Constructor.
subroutine time_interpolator_scalar(this, t, f_interpolated, f_n, tlag, n)
Interpolate a scalar field at time t from known scalar fields at different time steps.
Utilities.
Definition utils.f90:35
Stores a series (sequence) of fields, logically connected to a base field, and arranged according to ...
Provides a tool to perform interpolation in time.