stencil_star.F90

Go to the documentation of this file.

!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
!!
!! This program is free software; you can redistribute it and/or modify
!! it under the terms of the GNU General Public License as published by
!! the Free Software Foundation; either version 2, or (at your option)
!! any later version.
!!
!! This program is distributed in the hope that it will be useful,
!! but WITHOUT ANY WARRANTY; without even the implied warranty of
!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!! GNU General Public License for more details.
!!
!! You should have received a copy of the GNU General Public License
!! along with this program; if not, write to the Free Software
!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
!! 02110-1301, USA.
!!
 
#include "global.h"
 
module stencil_star_oct_m
  use debug_oct_m
  use global_oct_m
  use math_oct_m
  use messages_oct_m
  use nl_operator_oct_m
  use profiling_oct_m
  use stencil_oct_m
 
  implicit none
 
  private
  public ::                        &
    stencil_star_size_lapl,        &
    stencil_star_get_lapl,         &
    stencil_star_polynomials_lapl, &
    stencil_star_coeff_lapl,       &
    stencil_star_size_grad,        &
    stencil_star_get_grad,         &
    stencil_star_polynomials_grad, &
    stencil_star_coeff_grad
 
contains
 
  ! ---------------------------------------------------------
  integer function stencil_star_size_lapl(dim, order)
    integer, intent(in) :: dim
    integer, intent(in) :: order
 
    push_sub(stencil_star_size_lapl)
 
    stencil_star_size_lapl = 2*dim*order + 1
 
    pop_sub(stencil_star_size_lapl)
  end function stencil_star_size_lapl
 
 
  ! ---------------------------------------------------------
  subroutine stencil_star_get_lapl(this, dim, order)
    type(stencil_t), intent(inout) :: this
    integer,         intent(in)    :: dim
    integer,         intent(in)    :: order
 
    integer :: ii, jj, nn
    logical :: got_center
 
    push_sub(stencil_star_get_lapl)
 
    call stencil_allocate(this, dim, stencil_star_size_lapl(dim, order))
 
    got_center = .false.
 
    nn = 0
    do ii = 1, dim
      do jj = -order, order
 
        if (jj == 0) then
          if (got_center) then
            cycle
          else
            got_center = .true.
          end if
        end if
 
        nn = nn + 1
        this%points(ii, nn) = jj
      end do
    end do
 
    call stencil_init_center(this)
 
    pop_sub(stencil_star_get_lapl)
  end subroutine stencil_star_get_lapl
 
 
  ! ---------------------------------------------------------
  subroutine stencil_star_polynomials_lapl(dim, order, pol)
    integer, intent(in)  :: dim
    integer, intent(in)  :: order
    integer, intent(out) :: pol(:,:)
 
    integer :: i, j, n
 
    push_sub(stencil_star_polynomials_lapl)
 
    n = 1
    pol(:,:) = 0
    do i = 1, dim
      do j = 1, 2*order
        n = n + 1
        pol(i, n) = j
      end do
    end do
 
    pop_sub(stencil_star_polynomials_lapl)
  end subroutine stencil_star_polynomials_lapl
 
 
  ! ---------------------------------------------------------
  subroutine stencil_star_coeff_lapl(dim, order, h, lapl)
    integer,                intent(in)    :: dim
    integer,                intent(in)    :: order
    real(real64),           intent(in)    :: h(:)
    type(nl_operator_t),    intent(inout) :: lapl
 
    integer :: k, i, j, morder
    real(real64), allocatable :: cc(:,:,:)
 
    push_sub(stencil_star_coeff_lapl)
 
    assert(order >= 1)
 
    morder = 2*order
    safe_allocate(cc(0:morder, 0:morder, 0:2))
    call weights(2, morder, cc)
    lapl%w(1,:) = cc(0, morder, 2)*sum(1/h(1:dim)**2)
 
    k = 1
    do i = 1, dim
      do j = -order, -1
        k = k + 1
        lapl%w(k,:) = cc(-2*j-1, morder, 2) / h(i)**2
      end do
 
      do j = 1, order
        k = k + 1
        lapl%w(k,:) = cc( 2*j,   morder, 2) / h(i)**2
      end do
    end do
 
    safe_deallocate_a(cc)
 
    pop_sub(stencil_star_coeff_lapl)
  end subroutine stencil_star_coeff_lapl
 
 
 
  ! ---------------------------------------------------------
  integer function stencil_star_size_grad(order)
    integer, intent(in) :: order
 
    push_sub(stencil_star_size_grad)
 
    stencil_star_size_grad = 2*order + 1
 
    pop_sub(stencil_star_size_grad)
  end function stencil_star_size_grad
 
 
  ! ---------------------------------------------------------
  subroutine stencil_star_get_grad(this, dim, dir, order)
    type(stencil_t), intent(inout) :: this
    integer,         intent(in)    :: dim
    integer,         intent(in)    :: dir
    integer,         intent(in)    :: order
 
    integer :: i, n
 
    push_sub(stencil_star_get_grad)
 
    call stencil_allocate(this, dim, stencil_star_size_grad(order))
 
    n = 1
    do i = -order, order
      this%points(dir, n) = i
      n = n + 1
    end do
 
    call stencil_init_center(this)
 
    pop_sub(stencil_star_get_grad)
  end subroutine stencil_star_get_grad
 
 
  ! ---------------------------------------------------------
  subroutine stencil_star_polynomials_grad(dir, order, pol)
    integer, intent(in)  :: dir
    integer, intent(in)  :: order
    integer, intent(out) :: pol(:,:)
 
    integer :: j
 
    push_sub(stencil_star_polynomials_grad)
 
    pol(:,:) = 0
    do j = 0, 2*order
      pol(dir, j+1) = j
    end do
 
    pop_sub(stencil_star_polynomials_grad)
  end subroutine stencil_star_polynomials_grad
 
 
  ! ---------------------------------------------------------
  subroutine stencil_star_coeff_grad(order, h, grad)
    integer,             intent(in)    :: order
    real(real64),        intent(in)    :: h
    type(nl_operator_t), intent(inout) :: grad
 
    integer :: j, k, morder
    real(real64), allocatable :: cc(:,:,:)
 
    push_sub(stencil_star_coeff_grad)
 
    assert(order >= 1)
 
    morder = 2*order
    safe_allocate(cc(0:morder, 0:morder, 0:1))
    call weights(1, morder, cc)
 
    k = 1
    do j = -order, -1
      grad%w(k,:) = cc(-2*j-1, morder, 1) / h
      k = k + 1
    end do
 
    grad%w(k,:) = cc(0, morder, 1) / h
 
    do j = 1, order
      k = k + 1
      grad%w(k,:) = cc(2*j, morder, 1) / h
    end do
 
    safe_deallocate_a(cc)
 
    pop_sub(stencil_star_coeff_grad)
  end subroutine stencil_star_coeff_grad
 
end module stencil_star_oct_m
 
!! Local Variables:
!! mode: f90
!! coding: utf-8
!! End: