DISPATCH
time_slices_mod.f90
1 !===============================================================================
2 !> Compute time derivative from a sequence of time slices, using Lagrange
3 !> interpolation
4 !===============================================================================
6  USE io_mod
7  USE io_unit_mod
8  USE trace_mod
9  USE patch_mod
10  USE kinds_mod
11  USE lagrange_mod
12  implicit none
13  private
14  type, public:: time_slices_t
15  integer:: order=3
16  contains
17  procedure:: init
18  procedure, nopass:: interpolate
19  procedure, nopass:: derivative1
20  procedure, nopass:: derivative2
21  end type
22  integer, save:: verbose=0
23  integer, save:: order=3
24  integer, save:: id_debug=0
25  type(time_slices_t), public:: time_slices
26 CONTAINS
27 
28 !===============================================================================
29 !===============================================================================
30 SUBROUTINE init (self, nt)
31  class(time_slices_t):: self
32  integer:: nt
33  !.............................................................................
34  integer:: iostat
35  namelist /time_slices_params/ verbose, order, id_debug
36  logical, save:: first_time=.true.
37  !-----------------------------------------------------------------------------
38  call trace%begin ('time_slices_t%init')
39  !$omp critical (input_cr)
40  if (first_time) then
41  first_time = .false.
42  rewind(io_unit%input)
43  read (io_unit%input, time_slices_params, iostat=iostat)
44  order = min(order,nt-2)
45  write (io%output, time_slices_params)
46  write (io_unit%nml, time_slices_params)
47  flush (io_unit%nml)
48  call lagrange%test
49  end if
50  !$omp end critical (input_cr)
51  self%order = order
52  call trace%end()
53 END SUBROUTINE init
54 
55 !===============================================================================
56 !> Time interpolation. mem() is assumed to contain nt time slices of a 3-D
57 !> variable, with times consistent with the times produced by patch%timeslots()
58 !===============================================================================
59 SUBROUTINE interpolate (patch, mem, buffer)
60  class(patch_t):: patch
61  integer:: iv
62  real(kind=KindScalarVar), dimension(:,:,:,:):: mem
63  real(kind=4), dimension(:,:,:):: buffer
64  !.............................................................................
65  real(8):: times(patch%nt-1)
66  integer:: iit(patch%nt-1)
67  real(8), dimension(order+1):: w
68  integer:: i1, i2
69  !-----------------------------------------------------------------------------
70  if (verbose > 0) write(io_unit%log,*) patch%id, 'time_slices_t%interpolate'
71  call patch%timeslots (iit, times)
72  call lagrange%sequence_weights (patch%out_next, times, i1, i2, w, order)
73  call weighted_sum (patch, times, w, i1, i2, mem, buffer)
74 END SUBROUTINE interpolate
75 
76 !===============================================================================
77 !> First time derivate
78 !===============================================================================
79 SUBROUTINE derivative1 (patch, mem, buffer)
80  class(patch_t):: patch
81  integer:: iv
82  real(kind=KindScalarVar), dimension(:,:,:,:):: mem
83  real(kind=4), dimension(:,:,:):: buffer
84  !.............................................................................
85  real(8):: times(patch%nt-1)
86  real(8), dimension(order+1):: w
87  integer:: i1, i2
88  !-----------------------------------------------------------------------------
89  if (verbose > 0) write(io_unit%log,*) patch%id, 'time_slices_t%derivative1'
90  times = patch%t(patch%iit(1:patch%nt-1))
91  call lagrange%deriv_sequence_weights (patch%out_next, times, i1, i2, w, order)
92  call weighted_sum (patch, times, w, i1, i2, mem, buffer)
93 END SUBROUTINE derivative1
94 
95 !===============================================================================
96 !> Second time derivate
97 !===============================================================================
98 SUBROUTINE derivative2 (patch, mem, buffer)
99  class(patch_t):: patch
100  integer:: iv
101  real(kind=KindScalarVar), dimension(:,:,:,:):: mem
102  real(kind=4), dimension(:,:,:):: buffer
103  !.............................................................................
104  real(8):: times(patch%nt-1)
105  real(8), dimension(order+1):: w
106  integer:: i1, i2
107  !-----------------------------------------------------------------------------
108  if (verbose > 0) write(io_unit%log,*) patch%id, 'time_slices_t%derivative2'
109  times = patch%t(patch%iit(1:patch%nt-1))
110  call lagrange%deriv2_sequence_weights (patch%out_next, times, i1, i2, w, order)
111  call weighted_sum (patch, times, w, i1, i2, mem, buffer)
112 END SUBROUTINE derivative2
113 
114 !===============================================================================
115 !> Weighted sum of time slice values
116 !===============================================================================
117 SUBROUTINE weighted_sum (patch, times, w, i1, i2, mem, buffer)
118  class(patch_t):: patch
119  real(8), dimension(:):: times, w
120  integer:: i1, i2
121  real(kind=KindScalarVar), dimension(:,:,:,:):: mem
122  real(kind=4), dimension(:,:,:):: buffer
123  !.............................................................................
124  integer:: i
125  !-----------------------------------------------------------------------------
126  if (verbose>0 .or. patch%id==id_debug) then
127  write (io_unit%output,'(a,2i4,2x,10i4)') 'indices:', i1, i2, patch%iit(i1:i2)
128  write (io_unit%output,'(a,f12.6,2x,10f14.6)') ' times:', patch%out_next, times(i1:i2)
129  write (io_unit%output,'(a,f12.6,2x,10f14.3)') 'weights:', sum(w), w
130  end if
131  buffer = 0.0
132  do i=i1,i2
133  buffer = buffer + w(1+i-i1)*mem(:,:,:,i)
134  end do
135 END SUBROUTINE weighted_sum
136 
137 END MODULE time_slices_mod
Compute time derivative from a sequence of time slices, using Lagrange interpolation.
Template module for patches, which adds pointers to memory and mesh, and number of dimensions and var...
Definition: patch_mod.f90:6
Module for Lagrange interpolation.
Definition: lagrange_mod.f90:4
Definition: io_mod.f90:4