Mercurial > repos > public > sbplib_julia
changeset 515:d55008f5e2f3 feature/boundary_ops
Fix the range of the BoundaryRestriction tensor mapping (the range is zero for the 1D operator). Add some documentation and todos.
author | Vidar Stiernström <vidar.stiernstrom@it.uu.se> |
---|---|
date | Tue, 24 Nov 2020 15:28:10 +0100 |
parents | 14e722e8607d |
children | f932af8c1e56 |
files | src/SbpOperators/boundaryops/boundary_restriction.jl test/testSbpOperators.jl |
diffstat | 2 files changed, 23 insertions(+), 15 deletions(-) [+] |
line wrap: on
line diff
--- a/src/SbpOperators/boundaryops/boundary_restriction.jl Mon Nov 23 21:19:08 2020 +0100 +++ b/src/SbpOperators/boundaryops/boundary_restriction.jl Tue Nov 24 15:28:10 2020 +0100 @@ -1,3 +1,8 @@ +""" + boundary_restriction(grid,closureStencil,boundary) + +Creates a BoundaryRestriction operator for the specified boundary +""" function boundary_restriction(grid::EquidistantGrid{1,T}, closureStencil::Stencil{T,M}, boundary::CartesianBoundary{1}) where {T,M} r = region(boundary) return e = BoundaryRestriction(grid, closureStencil, r()) @@ -16,15 +21,14 @@ I = IdentityMapping{T}(size(restrict(grid,1))) return I⊗e end - export boundary_restriction """ - BoundaryRestriction{T,N,R} <: TensorMapping{T,1,1} + BoundaryRestriction{T,N,R} <: TensorMapping{T,0,1} Implements the boundary operator `e` as a TensorMapping """ -struct BoundaryRestriction{T,M,R<:Region} <: TensorMapping{T,1,1} +struct BoundaryRestriction{T,M,R<:Region} <: TensorMapping{T,0,1} stencil::Stencil{T,M} size::NTuple{1,Int} end @@ -34,9 +38,13 @@ return BoundaryRestriction{T,M,typeof(region)}(closureStencil,size(grid)) end -LazyTensors.range_size(e::BoundaryRestriction) = (1,) +LazyTensors.range_size(e::BoundaryRestriction) = (0,) LazyTensors.domain_size(e::BoundaryRestriction) = e.size +# TODO: Currently not working. +# We need to handle getindex for LazyTensorMappingApplication such that we pass more #indices than the +# range size of the TensorMapping. Or we need to be able to handle the case where we dont pass any index, for +# 0-dimensional tensormappings. " Restricts a grid function v on a grid of size m to the scalar element v[1]" function LazyTensors.apply(e::BoundaryRestriction{T,M,Lower}, v::AbstractVector{T}, i::Index{Lower}) where {T,M} @boundscheck if Int(i)!=1 @@ -55,18 +63,18 @@ " Transpose of a restriction is an inflation or prolongation. Inflates the scalar (1-element) vector to a vector of size of the grid" -function LazyTensors.apply_transpose(e::BoundaryRestriction{T,M,Lower}, v::AbstractVector{T}, i) where {T,M} +function LazyTensors.apply_transpose(e::BoundaryRestriction{T,M,Lower}, v::AbstractArray{T,0}, i) where {T,M} @boundscheck if !(0 < Int(i) <= e.size[1]) throw(BoundsError()) end - return e.stencil[Int(i)-1]*v[1] + return e.stencil[Int(i)-1]*v[] end " Transpose of a restriction is an inflation or prolongation. Inflates the scalar (1-element) vector to a vector of size of the grid" -function LazyTensors.apply_transpose(e::BoundaryRestriction{T,M,Upper}, v::AbstractVector{T}, i) where {T,M} +function LazyTensors.apply_transpose(e::BoundaryRestriction{T,M,Upper}, v::AbstractArray{T,0}, i) where {T,M} @boundscheck if !(0 < Int(i) <= e.size[1]) throw(BoundsError()) end - return e.stencil[e.size[1] - Int(i)]*v[1] + return e.stencil[e.size[1] - Int(i)]*v[] end
--- a/test/testSbpOperators.jl Mon Nov 23 21:19:08 2020 +0100 +++ b/test/testSbpOperators.jl Tue Nov 24 15:28:10 2020 +0100 @@ -182,13 +182,13 @@ e_r = BoundaryRestriction(g,op.eClosure,Upper()) v = evalOn(g,x->1+x^2) - u = [3.124] #How to handle scalars having to be arrays? It's kind of ugly. + u = fill(3.124) e_l*v isa LazyTensorMappingApplication - @test (e_l*v)[Index{Lower}(1)] == v[1] - @test (e_r*v)[Index{Upper}(4)] == v[end] - @test e_l'*u == [u[1], 0, 0, 0] - @test e_r'*u == [0, 0, 0, u[1]] + @test_broken (e_l*v)[Index{Lower}(1)] == v[1] + @test_broken (e_r*v)[Index{Upper}(4)] == v[end] + @test e_l'*u == [u[], 0, 0, 0] + @test e_r'*u == [0, 0, 0, u[]] @test_throws BoundsError (e_l*v)[Index{Lower}(3)] @test_throws BoundsError (e_r*v)[Index{Upper}(3)] @@ -206,8 +206,8 @@ v[:,2] = [7, 8, 9, 10] v[:,1] = [10, 11, 12, 13] - @test_broken e_w isa TensorMapping{T,1,2} where T - @test_broken e_w' isa TensorMapping{T,2,1} where T + @test e_w isa TensorMapping{T,1,2} where T + @test e_w' isa TensorMapping{T,2,1} where T