26 public :: filter_hyperdiff
28 interface filter_hyperdiff
29 module procedure filter_hyperdiff_2d
31 end interface filter_hyperdiff
49 subroutine filter_hyperdiff_2d( &
50 IA, IS, IE, JA, JS, JE, &
58 integer,
intent(in) :: IA, IS, IE
59 integer,
intent(in) :: JA, JS, JE
61 real(RP),
intent(inout) :: data(IA,JA)
63 integer,
intent(in) :: order
64 integer,
intent(in) :: nite
66 real(RP),
intent(in),
optional :: limiter_sign(IA,JA)
68 real(RP),
pointer :: p1(:,:)
69 real(RP),
pointer :: p2(:,:)
70 real(RP),
target :: work1(IA,JA)
71 real(RP),
target :: work2(IA,JA)
80 if ( mod(order,2) .ne. 0 )
then
81 log_error(
"FILTER_hyperdiff_2D", *)
"order must be even"
86 limiter =
present( limiter_sign )
91 call comm_vars8(
data(:,:), 1 )
92 call comm_wait (
data(:,:), 1, .false. )
97 work2(i,j) =
data(i,j)
105 do j = max(js,2), min(je,ja-1)
106 do i = max(is,2), min(ie,ia-1)
107 p2(i,j) = ( - p1(i+1,j) + p1(i,j)*2.0_rp - p1(i-1,j) &
108 - p1(i,j+1) + p1(i,j)*2.0_rp - p1(i,j-1) ) / 8.0_rp
112 do i = max(is,2), min(ie,ia-1)
113 p2(i,js) = ( - p1(i+1,js) + p1(i,js)*2.0_rp - p1(i-1,js) &
114 - p1(i,js+1) + p1(i,js) ) / 6.0_rp
117 do i = max(is,2), min(ie,ia-1)
118 p2(i,js-1) = p2(i,js)
122 do i = max(is,2), min(ie,ia-1)
123 p2(i,je) = ( - p1(i+1,je) + p1(i,je)*2.0_rp - p1(i-1,je) &
124 + p1(i,je) - p1(i,je-1) ) / 6.0_rp
127 do i = max(is,2), min(ie,ia-1)
128 p2(i,je+1) = p2(i,je)
132 do j = max(js,2), min(je,ja-1)
133 p2(is,j) = ( - p1(is+1,j) + p1(is,j) &
134 - p1(is,j+1) + p1(is,j)*2.0_rp - p1(is,j-1) ) / 6.0_rp
137 p2(is,js) = ( - p1(is+1,js) + p1(is,js) &
138 - p1(is,js+1) + p1(is,js) ) / 4.0_rp
141 p2(is,je) = ( - p1(is+1,je) + p1(is,je) &
142 + p1(is,je) - p1(is,je-1) ) / 4.0_rp
145 do j = max(js-1,1), min(je+1,ja)
146 p2(is-1,j) = p2(is,j)
150 do j = max(js,2), min(je,ja-1)
151 p2(ie,j) = ( + p1(ie,j) - p1(ie-1,j) &
152 - p1(ie,j+1) + p1(ie,j)*2.0_rp - p1(ie,j-1) ) / 6.0_rp
155 p2(ie,js) = ( + p1(ie,js) - p1(ie-1,js) &
156 - p1(ie,js+1) + p1(ie,js) ) / 4.0_rp
159 p2(ie,je) = ( + p1(ie,je) - p1(ie-1,je) &
160 + p1(ie,je) - p1(ie,je-1) ) / 4.0_rp
163 do j = max(js-1,1), min(je+1,ja)
164 p2(ie+1,j) = p2(ie,j)
168 call comm_vars8( p2(:,:), 1 )
169 call comm_wait ( p2(:,:), 1, .false. )
171 if ( mod(n,2) == 0 )
then
183 data(i,j) =
data(i,j) - p1(i,j)
191 data(i,j) = sign( max(
data(i,j) * limiter_sign(i,j), 0.0_rp ), limiter_sign(i,j) )
200 end subroutine filter_hyperdiff_2d
205 KA, KS, KE, IA, IS, IE, JA, JS, JE, &
213 integer,
intent(in) :: KA, KS, KE
214 integer,
intent(in) :: IA, IS, IE
215 integer,
intent(in) :: JA, JS, JE
217 real(RP),
intent(inout) :: data(KA,IA,JA)
219 integer,
intent(in) :: order
220 integer,
intent(in) :: nite
222 real(RP),
intent(in),
optional :: limiter_sign(KA,IA,JA)
224 real(RP) :: work_data(IA,JA)
225 real(RP) :: work_sign(IA,JA)
231 if (
present(limiter_sign) )
then
242 work_data(i,j) =
data(k,i,j)
249 work_sign(i,j) = limiter_sign(k,i,j)
253 call filter_hyperdiff_2d( ia, is, ie, ja, js, je, &
254 work_data(:,:), order, nite, &
255 limiter_sign = work_sign )
259 call filter_hyperdiff_2d( ia, is, ie, ja, js, je, &
260 work_data(:,:), order, nite )
267 data(k,i,j) = work_data(i,j)