1!! Copyright (C) 2008 X. Andrade 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_oct_m 22 use global_oct_m 23 use messages_oct_m 24 use profiling_oct_m 25 26 implicit none 27 28 private 29 public :: & 30 stencil_t, & 31 stencil_allocate, & 32 stencil_copy, & 33 stencil_end, & 34 stencil_init_center, & 35 stencil_union 36 37 type stargeneral_arms_t 38 ! Components are public by default 39 integer :: arms(1:3,1:3) 40 integer :: narms 41 end type stargeneral_arms_t 42 43 44 type stencil_t 45 ! Components are public by default 46 integer :: center 47 integer :: size 48 integer :: npoly 49 integer, pointer :: points(:, :) 50 51 ! The stargeneral arms 52 type(stargeneral_arms_t) :: stargeneral 53 end type stencil_t 54 55contains 56 57 !------------------------------------------------------- 58 subroutine stencil_allocate(this, size) 59 type(stencil_t), intent(out) :: this 60 integer, intent(in) :: size 61 62 PUSH_SUB(stencil_allocate) 63 64 this%size = size 65 this%npoly = size 66 67 SAFE_ALLOCATE(this%points(1:MAX_DIM, 1:size)) 68 69 this%points = 0 70 71 POP_SUB(stencil_allocate) 72 end subroutine stencil_allocate 73 74 !------------------------------------------------------- 75 subroutine stencil_copy(input, output) 76 type(stencil_t), intent(in) :: input 77 type(stencil_t), intent(out) :: output 78 79 PUSH_SUB(stencil_copy) 80 81 call stencil_allocate(output, input%size) 82 output%points(1:MAX_DIM, 1:output%size) = input%points(1:MAX_DIM, 1:output%size) 83 output%center = input%center 84 output%npoly = input%npoly 85 86 POP_SUB(stencil_copy) 87 end subroutine stencil_copy 88 89 90 !------------------------------------------------------- 91 subroutine stencil_end(this) 92 type(stencil_t), intent(inout) :: this 93 94 PUSH_SUB(stencil_end) 95 96 SAFE_DEALLOCATE_P(this%points) 97 98 POP_SUB(stencil_end) 99 end subroutine stencil_end 100 101 102 !------------------------------------------------------- 103 subroutine stencil_union(dim, st1, st2, stu) 104 integer, intent(in) :: dim 105 type(stencil_t), intent(inout) :: st1 106 type(stencil_t), intent(inout) :: st2 107 type(stencil_t), intent(inout) :: stu 108 109 integer :: idir, ii, jj, nstu 110 logical :: not_in_st1 111 112 PUSH_SUB(stencil_union) 113 114 call stencil_allocate(stu, st1%size + st2%size) 115 116 ! copy the first stencil 117 do ii = 1, st1%size 118 do idir = 1, dim 119 stu%points(idir, ii) = st1%points(idir, ii) 120 end do 121 end do 122 123 nstu = st1%size 124 125 do ii = 1, st2%size 126 127 not_in_st1 = .true. 128 129 ! check whether that point was already in the stencil 130 do jj = 1, st1%size 131 if(all(st1%points(1:dim, jj) == st2%points(1:dim, ii))) then 132 not_in_st1 = .false. 133 exit 134 end if 135 end do 136 137 if(not_in_st1) then !add it 138 nstu = nstu + 1 139 stu%points(1:dim, nstu) = st2%points(1:dim, ii) 140 end if 141 142 end do 143 144 stu%points(dim + 1:MAX_DIM, 1:nstu) = 0 145 146 stu%size = nstu 147 148 ! this is not defined for a union, which could be any combination. The 149 ! weights should already be known here 150 stu%npoly = -1 151 152 call stencil_init_center(stu) 153 154 POP_SUB(stencil_union) 155 end subroutine stencil_union 156 157 158 !------------------------------------------------------- 159 subroutine stencil_init_center(this) 160 type(stencil_t), intent(inout) :: this 161 162 integer :: ii 163 164 PUSH_SUB(stencil_init_center) 165 166 this%center = -1 167 168 do ii = 1, this%size 169 if(all(this%points(1:MAX_DIM, ii) == 0)) this%center = ii 170 end do 171 172 POP_SUB(stencil_init_center) 173 end subroutine stencil_init_center 174 175end module stencil_oct_m 176 177!! Local Variables: 178!! mode: f90 179!! coding: utf-8 180!! End: 181