Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
tree_amg_utils.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!
35 use num_types, only : rp
36 use tree_amg, only : tamg_hierarchy_t
37 use gather_scatter, only : gs_op_add
38 implicit none
39 private
40
42
43contains
44
51 subroutine tamg_sample_matrix_val(val, amg, lvl, i, j)
52 real(kind=rp), intent(inout) :: val
53 type(tamg_hierarchy_t), intent(inout) :: amg
54 integer, intent(in) :: lvl, i, j
55 real(kind=rp), allocatable :: u(:), v(:)
56 integer :: n
57 n = amg%lvl(lvl+1)%fine_lvl_dofs
58 allocate( u(n) )
59 allocate( v(n) )
60 u = 0.0_rp
61 v = 0.0_rp
62
63 u(i) = 1.0_rp
64
66 if (lvl.eq.0) call amg%gs_h%op(u, n, gs_op_add)
67
68 call amg%matvec(v, u, lvl)
69
70 val = v(j)
71 deallocate( u )
72 deallocate( v )
73 end subroutine tamg_sample_matrix_val
74
77 subroutine tamg_print_matrix(amg, lvl)
78 type(tamg_hierarchy_t), intent(inout) :: amg
79 integer, intent(in) :: lvl
80 integer :: n, i, j
81 real(kind=rp) :: val
82 real(kind=rp), allocatable :: u(:), v(:)
83
84 n = amg%lvl(lvl+1)%fine_lvl_dofs
85 allocate( u(n) )
86 allocate( v(n) )
87
88 do i = 1, n
89 do j = 1, n
90 u = 0.0_rp
91 v = 0.0_rp
92
93 u(i) = 1.0_rp
95 if (lvl.eq.0) call amg%gs_h%op(u, n, gs_op_add)
96
97 call amg%matvec(v, u, lvl)
98
99 val = v(j)
100
101 print *, j, i, val
102 end do
103 end do
104
105 deallocate( u )
106 deallocate( v )
107 end subroutine tamg_print_matrix
108
111 subroutine tamg_print_restriction_matrix(amg, lvl)
112 type(tamg_hierarchy_t), intent(inout) :: amg
113 integer, intent(in) :: lvl
114 integer :: n, m, i, j
115 real(kind=rp) :: val
116 real(kind=rp), allocatable :: u(:), v(:)
117
118 n = amg%lvl(lvl+1)%fine_lvl_dofs
119 m = amg%lvl(lvl+2)%fine_lvl_dofs
120 allocate( u(n) )
121 allocate( v(m) )
122
123 do i = 1, n
124 do j = 1, m
125 u = 0.0_rp
126 v = 0.0_rp
127
128 u(i) = 1.0_rp
130 if (lvl.eq.0) call amg%gs_h%op(u, n, gs_op_add)
131
132 call amg%interp_f2c(v, u, lvl+1)
133
134 val = v(j)
135
136 print *, j, i, val
137 end do
138 end do
139
140 deallocate( u )
141 deallocate( v )
142 end subroutine tamg_print_restriction_matrix
143
147 type(tamg_hierarchy_t), intent(inout) :: amg
148 integer, intent(in) :: lvl
149 integer :: n, m, i, j
150 real(kind=rp) :: val
151 real(kind=rp), allocatable :: u(:), v(:)
152
153 n = amg%lvl(lvl+1)%fine_lvl_dofs
154 m = amg%lvl(lvl+2)%fine_lvl_dofs
155 allocate( u(m) )
156 allocate( v(n) )
157
158 do i = 1, m
159 do j = 1, n
160 u = 0.0_rp
161 v = 0.0_rp
162
163 u(i) = 1.0_rp
164 !---!> Get the ghosts to agree on the value of 1.0
165 !---if (lvl.eq.0) call amg%gs_h%op(u, n, GS_OP_ADD)
166
167 call amg%interp_c2f(v, u, lvl+1)
168
169 val = v(j)
170
171 print *, j, i, val
172 end do
173 end do
174
175 deallocate( u )
176 deallocate( v )
177 end subroutine tamg_print_prolongation_matrix
178
179
180end module tree_amg_utils
Gather-scatter.
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Implements utilities for the TreeAMG hierarchy structure.
subroutine tamg_print_restriction_matrix(amg, lvl)
never call this function also this does not account for duplicates/ghosts in dofs
subroutine tamg_print_matrix(amg, lvl)
never call this function also this does not account for duplicates/ghosts in dofs
subroutine tamg_print_prolongation_matrix(amg, lvl)
never call this function also this does not account for duplicates/ghosts in dofs
subroutine, public tamg_sample_matrix_val(val, amg, lvl, i, j)
Sample the values in a matix (expensive, use with caution)
Implements the base type for TreeAMG hierarchy structure.
Definition tree_amg.f90:34
Type for a TreeAMG hierarchy.
Definition tree_amg.f90:94