1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch 2!! 3!! This program is free software; you can redistribute it and/or modify 4!! it under the terms of the GNU General Public License as published by 5!! the Free Software Foundation; either version 2, or (at your option) 6!! any later version. 7!! 8!! This program is distributed in the hope that it will be useful, 9!! but WITHOUT ANY WARRANTY; without even the implied warranty of 10!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 11!! GNU General Public License for more details. 12!! 13!! You should have received a copy of the GNU General Public License 14!! along with this program; if not, write to the Free Software 15!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 16!! 02110-1301, USA. 17!! 18 19#include "global.h" 20 21module stencil_cube_oct_m 22 use global_oct_m 23 use messages_oct_m 24 use stencil_oct_m 25 26 implicit none 27 28 private 29 public :: & 30 stencil_cube_size_lapl, & 31 stencil_cube_get_lapl, & 32 stencil_cube_polynomials_lapl, & 33 stencil_cube_size_grad, & 34 stencil_cube_get_grad, & 35 stencil_cube_polynomials_grad 36 37 38contains 39 40 ! --------------------------------------------------------- 41 integer function stencil_cube_size_lapl(dim, order) 42 integer, intent(in) :: dim 43 integer, intent(in) :: order 44 45 PUSH_SUB(stencil_cube_size_lapl) 46 47 stencil_cube_size_lapl = (2*order+1)**dim 48 49 POP_SUB(stencil_cube_size_lapl) 50 end function stencil_cube_size_lapl 51 52 53 ! --------------------------------------------------------- 54 subroutine stencil_cube_get_lapl(this, dim, order) 55 type(stencil_t), intent(out) :: this 56 integer, intent(in) :: dim 57 integer, intent(in) :: order 58 59 integer :: i, j, k, n 60 61 PUSH_SUB(stencil_cube_get_lapl) 62 63 call stencil_allocate(this, stencil_cube_size_lapl(dim, order)) 64 65 n = 1 66 select case(dim) 67 case(1) 68 do i = -order, order 69 this%points(1, n) = i 70 n = n + 1 71 end do 72 case(2) 73 do i = -order, order 74 do j = -order, order 75 this%points(1, n) = i 76 this%points(2, n) = j 77 n = n + 1 78 end do 79 end do 80 case(3) 81 do i = -order, order 82 do j = -order, order 83 do k = -order, order 84 this%points(1, n) = i 85 this%points(2, n) = j 86 this%points(3, n) = k 87 n = n + 1 88 end do 89 end do 90 end do 91 end select 92 93 call stencil_init_center(this) 94 95 POP_SUB(stencil_cube_get_lapl) 96 end subroutine stencil_cube_get_lapl 97 98 99 ! --------------------------------------------------------- 100 subroutine stencil_cube_polynomials_lapl(dim, order, pol) 101 integer, intent(in) :: dim 102 integer, intent(in) :: order 103 integer, intent(out) :: pol(:,:) !< pol(dim, order) 104 105 integer :: i, j, k, n 106 107 PUSH_SUB(stencil_cube_polynomials_lapl) 108 109 n = 1 110 select case(dim) 111 case(1) 112 do i = 0, 2*order 113 pol(1, n) = i 114 n = n + 1 115 end do 116 case(2) 117 do i = 0, 2*order 118 do j = 0, 2*order 119 pol(1, n) = i 120 pol(2, n) = j 121 n = n + 1 122 end do 123 end do 124 case(3) 125 do i = 0, 2*order 126 do j = 0, 2*order 127 do k = 0, 2*order 128 pol(1, n) = i 129 pol(2, n) = j 130 pol(3, n) = k 131 n = n + 1 132 end do 133 end do 134 end do 135 end select 136 137 POP_SUB(stencil_cube_polynomials_lapl) 138 end subroutine stencil_cube_polynomials_lapl 139 140 141 !> Now come the gradient routines. As this stencil is the same for 142 !! the laplacian and the gradient, these routines just call the 143 !! corresponding ones for the laplacian 144 145 ! --------------------------------------------------------- 146 integer function stencil_cube_size_grad(dim, order) 147 integer, intent(in) :: dim 148 integer, intent(in) :: order 149 150 PUSH_SUB(stencil_cube_size_grad) 151 152 stencil_cube_size_grad = stencil_cube_size_lapl(dim, order) 153 154 POP_SUB(stencil_cube_size_grad) 155 end function stencil_cube_size_grad 156 157 158 ! --------------------------------------------------------- 159 subroutine stencil_cube_get_grad(this, dim, order) 160 type(stencil_t), intent(out) :: this 161 integer, intent(in) :: dim 162 integer, intent(in) :: order 163 164 PUSH_SUB(stencil_cube_get_grad) 165 166 call stencil_cube_get_lapl(this, dim, order) 167 168 POP_SUB(stencil_cube_get_grad) 169 end subroutine stencil_cube_get_grad 170 171 172 ! --------------------------------------------------------- 173 subroutine stencil_cube_polynomials_grad(dim, order, pol) 174 integer, intent(in) :: dim 175 integer, intent(in) :: order 176 integer, intent(out) :: pol(:,:) !< pol(sb%dim, order) 177 178 PUSH_SUB(stencil_cube_polynomials_grad) 179 180 call stencil_cube_polynomials_lapl(dim, order, pol) 181 182 POP_SUB(stencil_cube_polynomials_grad) 183 end subroutine stencil_cube_polynomials_grad 184 185end module stencil_cube_oct_m 186 187!! Local Variables: 188!! mode: f90 189!! coding: utf-8 190!! End: 191