15 #include "../base/base_uses.f90"
37 INTEGER,
INTENT(in) :: itype, nd
38 INTEGER,
INTENT(out) :: nrange
39 REAL(kind=
dp),
DIMENSION(0:nd),
INTENT(out) :: a, x
41 INTEGER :: i, i_all, m, ni, nt
42 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: y
43 REAL(kind=
dp),
DIMENSION(:),
POINTER :: cg, cgt, ch, cht
54 ALLOCATE (y(0:nd), stat=i_all)
56 WRITE (*, *)
' scaling_function: problem of memory allocation'
68 CALL back_trans(nd, nt, x, y, m, ch, cg)
69 CALL dcopy(nt, y, 1, x, 1)
77 a(i) = 1._dp*i*ni/nd - (.5_dp*ni - 1._dp)
81 DEALLOCATE (ch, cg, cgt, cht)
92 SUBROUTINE wavelet_function(itype, nd, a, x)
95 INTEGER,
INTENT(in) :: itype, nd
96 REAL(kind=
dp),
DIMENSION(0:nd),
INTENT(out) :: a, x
98 INTEGER :: i, i_all, m, ni, nt
99 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: y
100 REAL(kind=
dp),
DIMENSION(:),
POINTER :: cg, cgt, ch, cht
109 ALLOCATE (y(0:nd), stat=i_all)
111 WRITE (*, *)
' wavelet_function: problem of memory allocation'
119 x(nt + nt/2 - 1) = 1._dp
123 CALL back_trans(nd, nt, x, y, m, ch, cg)
124 CALL dcopy(nd, y, 1, x, 1)
132 a(i) = 1._dp*i*ni/nd - (.5_dp*ni - .5_dp)
136 DEALLOCATE (ch, cg, cgt, cht)
139 END SUBROUTINE wavelet_function
151 INTEGER,
INTENT(in) :: itype, n_iter, n_range
152 REAL(kind=
dp),
INTENT(inout) :: kernel_scf(-n_range:n_range)
153 REAL(kind=
dp),
INTENT(out) :: kern_1_scf(-n_range:n_range)
156 REAL(kind=
dp),
DIMENSION(:),
POINTER :: cg, cgt, ch, cht
161 CALL scf_recurs(n_iter, n_range, kernel_scf, kern_1_scf, m, ch)
162 DEALLOCATE (ch, cg, cgt, cht)
171 SUBROUTINE zero(n, x)
172 INTEGER,
INTENT(in) :: n
173 REAL(kind=
dp),
INTENT(out) :: x(n)
196 SUBROUTINE for_trans(nd, nt, x, y, m, cgt, cht)
197 INTEGER,
INTENT(in) :: nd, nt
198 REAL(kind=
dp),
INTENT(in) :: x(0:nd - 1)
199 REAL(kind=
dp),
INTENT(out) :: y(0:nd - 1)
201 REAL(kind=
dp),
DIMENSION(:),
POINTER :: cgt, cht
219 IF (ind .GE. nt)
THEN
226 y(i) = y(i) + cht(j)*x(ind)
227 y(nt/2 + i) = y(nt/2 + i) + cgt(j)*x(ind)
232 END SUBROUTINE for_trans
244 SUBROUTINE back_trans(nd, nt, x, y, m, ch, cg)
250 INTEGER,
INTENT(in) :: nd, nt
251 REAL(kind=
dp),
INTENT(in) :: x(0:nd - 1)
252 REAL(kind=
dp),
INTENT(out) :: y(0:nd - 1)
254 REAL(kind=
dp),
DIMENSION(:),
POINTER :: ch, cg
273 IF (ind .GE. nt/2)
THEN
280 y(2*i + 0) = y(2*i + 0) + ch(2*j - 0)*x(ind) + cg(2*j - 0)*x(ind + nt/2)
281 y(2*i + 1) = y(2*i + 1) + ch(2*j + 1)*x(ind) + cg(2*j + 1)*x(ind + nt/2)
286 END SUBROUTINE back_trans
296 SUBROUTINE ftest(m, ch, cg, cgt, cht)
298 REAL(kind=
dp),
DIMENSION(:),
POINTER :: ch, cg, cgt, cht
300 CHARACTER(len=*),
PARAMETER :: fmt22 =
"(a,i3,i4,4(e17.10))"
303 REAL(kind=
dp) :: eps, t1, t2, t3, t4
316 IF (l - 2*i .GE. -m .AND. l - 2*i .LE. m .AND. &
317 l - 2*j .GE. -m .AND. l - 2*j .LE. m)
THEN
318 t1 = t1 + ch(l - 2*i)*cht(l - 2*j)
319 t2 = t2 + cg(l - 2*i)*cgt(l - 2*j)
320 t3 = t3 + ch(l - 2*i)*cgt(l - 2*j)
321 t4 = t4 + cht(l - 2*i)*cg(l - 2*j)
326 IF (abs(t1 - 1._dp) .GT. eps .OR. abs(t2 - 1._dp) .GT. eps .OR. &
327 abs(t3) .GT. eps .OR. abs(t4) .GT. eps)
THEN
328 WRITE (*, fmt22)
'Orthogonality ERROR', i, j, t1, t2, t3, t4
331 IF (abs(t1) .GT. eps .OR. abs(t2) .GT. eps .OR. &
332 abs(t3) .GT. eps .OR. abs(t4) .GT. eps)
THEN
333 WRITE (*, fmt22)
'Orthogonality ERROR', i, j, t1, t2, t3, t4
339 WRITE (*, *)
'FILTER TEST PASSED'
353 SUBROUTINE scf_recurs(n_iter, n_range, kernel_scf, kern_1_scf, m, ch)
354 INTEGER,
INTENT(in) :: n_iter, n_range
355 REAL(kind=
dp),
INTENT(inout) :: kernel_scf(-n_range:n_range)
356 REAL(kind=
dp),
INTENT(out) :: kern_1_scf(-n_range:n_range)
358 REAL(kind=
dp),
DIMENSION(:),
POINTER :: ch
360 INTEGER :: i, i_iter, ind, j
361 REAL(kind=
dp) :: kern, kern_tot
365 loop_iter_scf:
DO i_iter = 1, n_iter
366 kern_1_scf(:) = kernel_scf(:)
367 kernel_scf(:) = 0._dp
368 loop_iter_i:
DO i = 0, n_range
372 IF (abs(ind) > n_range)
THEN
375 kern = kern_1_scf(ind)
377 kern_tot = kern_tot + ch(j)*kern
379 IF (kern_tot == 0._dp)
THEN
383 kernel_scf(i) = 0.5_dp*kern_tot
384 kernel_scf(-i) = kernel_scf(i)
388 END SUBROUTINE scf_recurs
Defines the basic variable types.
integer, parameter, public dp
Filters for interpolating scaling functions .
subroutine, public lazy_arrays(itype, m, ch, cg, cgt, cht)
...
Creates the wavelet kernel for the wavelet based poisson solver.
subroutine, public scf_recursion(itype, n_iter, n_range, kernel_scf, kern_1_scf)
Do iterations to go from p0gauss to pgauss order interpolating scaling function.
subroutine, public scaling_function(itype, nd, nrange, a, x)
Calculate the values of a scaling function in real uniform grid.