Mercurial > repos > public > sbplib_julia
diff SbpOperators/src/laplace/laplace.jl @ 300:b00eea62c78e
Create 1D tensor mapping for diagonal norm quadratures, and make the multi-dimensional quadrature use those. Move Qudrature from laplace.jl into Quadrature.jl
author | Vidar Stiernström <vidar.stiernstrom@it.uu.se> |
---|---|
date | Tue, 23 Jun 2020 17:32:54 +0200 |
parents | f63232aeb1c6 |
children | 6fa2ba769ae3 |
line wrap: on
line diff
--- a/SbpOperators/src/laplace/laplace.jl Tue Jun 23 17:31:50 2020 +0200 +++ b/SbpOperators/src/laplace/laplace.jl Tue Jun 23 17:32:54 2020 +0200 @@ -2,8 +2,8 @@ Laplace{Dim,T<:Real,N,M,K} <: TensorOperator{T,Dim} Implements the Laplace operator `L` in Dim dimensions as a tensor operator -The multi-dimensional tensor operator simply consists of a tuple of the 1D -Laplace tensor operator as defined by ConstantLaplaceOperator. +The multi-dimensional tensor operator consists of a tuple of 1D SecondDerivative +tensor operators. """ struct Laplace{Dim,T<:Real,N,M,K} <: TensorOperator{T,Dim} D2::NTuple(Dim,SecondDerivative{T,N,M,K}) @@ -17,20 +17,23 @@ error("not implemented") end +function LazyTensors.apply_transpose(L::Laplace{Dim,T}, v::AbstractArray{T,Dim}, I::NTuple{Dim,Index}) where {T,Dim} = LazyTensors.apply(L, v, I) + # u = L*v function LazyTensors.apply(L::Laplace{1,T}, v::AbstractVector{T}, I::NTuple{1,Index}) where T - return apply(L.D2[1],v,I) + @inbounds u = apply(L.D2[1],v,I) + return u end @inline function LazyTensors.apply(L::Laplace{2,T}, v::AbstractArray{T,2}, I::NTuple{2,Index}) where T # 2nd x-derivative @inbounds vx = view(v, :, Int(I[2])) - @inbounds uᵢ = apply(L.D2[1], vx , (I[1],)) #Tuple conversion here is ugly. How to do it? Should we use indexing here? + @inbounds uᵢ = apply(L.D2[1], vx , I[1]) # 2nd y-derivative @inbounds vy = view(v, Int(I[1]), :) - @inbounds uᵢ += apply(L.D2[2], vy , (I[2],)) #Tuple conversion here is ugly. How to do it? + @inbounds uᵢ += apply(L.D2[2], vy , I[2]) return uᵢ end @@ -42,31 +45,6 @@ boundary_quadrature(L::Laplace, bId::CartesianBoundary) = BoundaryQuadrature(L.op, L.grid, bId) export quadrature -# At the moment the grid property is used all over. It could possibly be removed if we implement all the 1D operators as TensorMappings -""" - Quadrature{Dim,T<:Real,N,M,K} <: TensorMapping{T,Dim,Dim} - -Implements the quadrature operator `H` of Dim dimension as a TensorMapping -""" -struct Quadrature{Dim,T<:Real,N,M,K} <: TensorOperator{T,Dim} - op::D2{T,N,M,K} - grid::EquidistantGrid{Dim,T} -end -export Quadrature - -LazyTensors.domain_size(H::Quadrature{Dim}, range_size::NTuple{Dim,Integer}) where Dim = range_size - -@inline function LazyTensors.apply(H::Quadrature{2,T}, v::AbstractArray{T,2}, I::NTuple{2,Index}) where T - N = size(H.grid) - # Quadrature in x direction - @inbounds q = apply_quadrature(H.op, spacing(H.grid)[1], v[I] , I[1], N[1]) - # Quadrature in y-direction - @inbounds q = apply_quadrature(H.op, spacing(H.grid)[2], q, I[2], N[2]) - return q -end - -LazyTensors.apply_transpose(H::Quadrature{2,T}, v::AbstractArray{T,2}, I::NTuple{2,Index}) where T = LazyTensors.apply(H,v,I) - """ InverseQuadrature{Dim,T<:Real,N,M,K} <: TensorMapping{T,Dim,Dim}