view SbpOperators/src/constantstenciloperator.jl @ 262:f1e90a92ad74 boundary_conditions

Add Quadrature and InverseQuadrature for Laplace as TensorMappings. Implement and test the 2D case. Fix implementation of apply_transpose for BoundaryQuadrature and add tests.
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Tue, 26 Nov 2019 08:28:26 -0800
parents 396eadb652bd
children ccef055233a2
line wrap: on
line source

export apply

abstract type ConstantStencilOperator end

# Apply for different regions Lower/Interior/Upper or Unknown region
@inline function apply(op::ConstantStencilOperator, h::Real, v::AbstractVector, i::Index{Lower})
    return @inbounds h*h*apply(op.closureStencils[Int(i)], v, Int(i))
end

@inline function apply(op::ConstantStencilOperator, h::Real, v::AbstractVector, i::Index{Interior})
    return @inbounds h*h*apply(op.innerStencil, v, Int(i))
end

@inline function apply(op::ConstantStencilOperator, h::Real, v::AbstractVector, i::Index{Upper})
    N = length(v)
    return @inbounds h*h*Int(op.parity)*apply_backwards(op.closureStencils[N-Int(i)+1], v, Int(i))
end

@inline function apply(op::ConstantStencilOperator, h::Real, v::AbstractVector, index::Index{Unknown})
    cSize = closuresize(op)
    N = length(v)

    i = Int(index)

    if 0 < i <= cSize
        return apply(op, h, v, Index{Lower}(i))
    elseif cSize < i <= N-cSize
        return apply(op, h, v, Index{Interior}(i))
    elseif N-cSize < i <= N
        return apply(op, h, v, Index{Upper}(i))
    else
        error("Bounds error") # TODO: Make this more standard
    end
end

# Wrapper functions for using regular indecies without specifying regions
@inline function apply(op::ConstantStencilOperator, h::Real, v::AbstractVector, i::Int)
    return apply(op, h, v, Index{Unknown}(i))
end