Mercurial > repos > public > sbplib_julia
comparison src/SbpOperators/quadrature/inverse_diagonal_quadrature.jl @ 504:21fba50cb5b0 feature/quadrature_as_outer_product
Use LazyOuterProduct to construct multi-dimensional quadratures. This change allwed to:
- Replace the types Quadrature and InverseQuadrature by functions returning outer products of the 1D operators.
- Avoid convoluted naming of types. DiagonalInnerProduct is now renamed to DiagonalQuadrature, similarly for InverseDiagonalInnerProduct.
author | Vidar Stiernström <vidar.stiernstrom@it.uu.se> |
---|---|
date | Sat, 07 Nov 2020 13:07:31 +0100 |
parents | src/SbpOperators/quadrature/inverse_diagonal_inner_product.jl@0d93d406c222 |
children | c2f991b819fc |
comparison
equal
deleted
inserted
replaced
503:fbbb3733650c | 504:21fba50cb5b0 |
---|---|
1 """ | |
2 inverse_diagonal_quadrature(g,quadrature_closure) | |
3 | |
4 Constructs the diagonal quadrature inverse operator `Qi` on a grid of `Dim` dimensions as | |
5 a `TensorMapping`. The one-dimensional operator is a InverseDiagonalQuadrature, while | |
6 the multi-dimensional operator is the outer-product of the one-dimensional operators | |
7 in each coordinate direction. | |
8 """ | |
9 function inverse_diagonal_quadrature(g::EquidistantGrid{Dim}, quadrature_closure) where Dim | |
10 Hi = InverseDiagonalQuadrature(restrict(g,1), quadrature_closure) | |
11 for i ∈ 2:Dim | |
12 Hi = Hi⊗InverseDiagonalQuadrature(restrict(g,i), quadrature_closure) | |
13 end | |
14 return Hi | |
15 end | |
16 export inverse_diagonal_quadrature | |
17 | |
18 | |
19 """ | |
20 InverseDiagonalQuadrature{Dim,T<:Real,M} <: TensorMapping{T,1,1} | |
21 | |
22 Implements the inverse diagonal inner product operator `Hi` of as a 1D TensorOperator | |
23 """ | |
24 struct InverseDiagonalQuadrature{T<:Real,M} <: TensorMapping{T,1,1} | |
25 h_inv::T | |
26 closure::NTuple{M,T} | |
27 size::Tuple{Int} | |
28 end | |
29 export InverseDiagonalQuadrature | |
30 | |
31 function InverseDiagonalQuadrature(g::EquidistantGrid{1}, quadrature_closure) | |
32 return InverseDiagonalQuadrature(inverse_spacing(g)[1], 1 ./ quadrature_closure, size(g)) | |
33 end | |
34 | |
35 | |
36 LazyTensors.range_size(Hi::InverseDiagonalQuadrature) = Hi.size | |
37 LazyTensors.domain_size(Hi::InverseDiagonalQuadrature) = Hi.size | |
38 | |
39 | |
40 function LazyTensors.apply(Hi::InverseDiagonalQuadrature{T}, v::AbstractVector{T}, I::Index{Lower}) where T | |
41 return @inbounds Hi.h_inv*Hi.closure[Int(I)]*v[Int(I)] | |
42 end | |
43 | |
44 function LazyTensors.apply(Hi::InverseDiagonalQuadrature{T}, v::AbstractVector{T}, I::Index{Upper}) where T | |
45 N = length(v); | |
46 return @inbounds Hi.h_inv*Hi.closure[N-Int(I)+1]*v[Int(I)] | |
47 end | |
48 | |
49 function LazyTensors.apply(Hi::InverseDiagonalQuadrature{T}, v::AbstractVector{T}, I::Index{Interior}) where T | |
50 return @inbounds Hi.h_inv*v[Int(I)] | |
51 end | |
52 | |
53 function LazyTensors.apply(Hi::InverseDiagonalQuadrature, v::AbstractVector{T}, index::Index{Unknown}) where T | |
54 N = length(v); | |
55 r = getregion(Int(index), closuresize(Hi), N) | |
56 i = Index(Int(index), r) | |
57 return LazyTensors.apply(Hi, v, i) | |
58 end | |
59 | |
60 LazyTensors.apply_transpose(Hi::InverseDiagonalQuadrature{T}, v::AbstractVector{T}, I::Index) where T = LazyTensors.apply(Hi,v,I) | |
61 | |
62 | |
63 closuresize(Hi::InverseDiagonalQuadrature{T,M}) where {T,M} = M | |
64 export closuresize |