Mercurial > repos > public > sbplib_julia

annotate sbpD2.jl @ 85:8d505e9bc715 cell_based_test

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Merge with default
author Jonatan Werpers <jonatan@werpers.com>
date Fri, 25 Jan 2019 15:26:47 +0100
parents b795ec7f9ca0 48079bd39969
children c0f33eccd527 8324c82c2dfb
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
34
bb841977d198 Move stencil operator application to its own function
Jonatan Werpers <jonatan@werpers.com>
parents: 24
diff changeset
1 abstract type ConstantStencilOperator end
bb841977d198 Move stencil operator application to its own function
Jonatan Werpers <jonatan@werpers.com>
parents: 24
diff changeset
2
78
fbf7398f8154 Inline and inbounds everything
Ylva Rydin <ylva.rydin@telia.com>
parents: 68
diff changeset
3 @inline function apply(op::ConstantStencilOperator, h::Real, v::AbstractVector, i::Int)
34
bb841977d198 Move stencil operator application to its own function
Jonatan Werpers <jonatan@werpers.com>
parents: 24
diff changeset
4 cSize = closureSize(op)
56
27a8d3021a1c Convert apply functions to cell-based
Ylva Rydin <ylva.rydin@telia.com>
parents: 47
diff changeset
5 N = length(v)
34
bb841977d198 Move stencil operator application to its own function
Jonatan Werpers <jonatan@werpers.com>
parents: 24
diff changeset
6
56
27a8d3021a1c Convert apply functions to cell-based
Ylva Rydin <ylva.rydin@telia.com>
parents: 47
diff changeset
7 if i ∈ range(1; length=cSize)
72
4640839b1616 Make apply operator more efficient using @inbounds and @inline
Jonatan Werpers <jonatan@werpers.com>
parents: 68
diff changeset
8 @inbounds uᵢ = apply(op.closureStencils[i], v, i)/h^2
56
27a8d3021a1c Convert apply functions to cell-based
Ylva Rydin <ylva.rydin@telia.com>
parents: 47
diff changeset
9 elseif i ∈ range(N - cSize+1, length=cSize)
72
4640839b1616 Make apply operator more efficient using @inbounds and @inline
Jonatan Werpers <jonatan@werpers.com>
parents: 68
diff changeset
10 @inbounds uᵢ = Int(op.parity)*apply(flip(op.closureStencils[N-i+1]), v, i)/h^2
56
27a8d3021a1c Convert apply functions to cell-based
Ylva Rydin <ylva.rydin@telia.com>
parents: 47
diff changeset
11 else
72
4640839b1616 Make apply operator more efficient using @inbounds and @inline
Jonatan Werpers <jonatan@werpers.com>
parents: 68
diff changeset
12 @inbounds uᵢ = apply(op.innerStencil, v, i)/h^2
34
bb841977d198 Move stencil operator application to its own function
Jonatan Werpers <jonatan@werpers.com>
parents: 24
diff changeset
13 end
bb841977d198 Move stencil operator application to its own function
Jonatan Werpers <jonatan@werpers.com>
parents: 24
diff changeset
14
56
27a8d3021a1c Convert apply functions to cell-based
Ylva Rydin <ylva.rydin@telia.com>
parents: 47
diff changeset
15 return uᵢ
34
bb841977d198 Move stencil operator application to its own function
Jonatan Werpers <jonatan@werpers.com>
parents: 24
diff changeset
16 end
bb841977d198 Move stencil operator application to its own function
Jonatan Werpers <jonatan@werpers.com>
parents: 24
diff changeset
17
47
8c4cfa680e5a Use enum for parity
Jonatan Werpers <jonatan@werpers.com>
parents: 43
diff changeset
18 @enum Parity begin
8c4cfa680e5a Use enum for parity
Jonatan Werpers <jonatan@werpers.com>
parents: 43
diff changeset
19 odd = -1
8c4cfa680e5a Use enum for parity
Jonatan Werpers <jonatan@werpers.com>
parents: 43
diff changeset
20 even = 1
8c4cfa680e5a Use enum for parity
Jonatan Werpers <jonatan@werpers.com>
parents: 43
diff changeset
21 end
34
bb841977d198 Move stencil operator application to its own function
Jonatan Werpers <jonatan@werpers.com>
parents: 24
diff changeset
22
84
48079bd39969 Change to using tuples in stencils and ops
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
23 struct D2{T,N,M,K} <: ConstantStencilOperator
1
a286be5f8a32 Add type to sbpD2.jl
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
24 quadratureClosure::Vector{T}
84
48079bd39969 Change to using tuples in stencils and ops
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
25 innerStencil::Stencil{T,N}
48079bd39969 Change to using tuples in stencils and ops
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
26 closureStencils::NTuple{M, Stencil{T,K}}
1
a286be5f8a32 Add type to sbpD2.jl
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
27 eClosure::Vector{T}
a286be5f8a32 Add type to sbpD2.jl
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
28 dClosure::Vector{T}
47
8c4cfa680e5a Use enum for parity
Jonatan Werpers <jonatan@werpers.com>
parents: 43
diff changeset
29 parity::Parity
1
a286be5f8a32 Add type to sbpD2.jl
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
30 end
2
43be32298ae2 Add function to get closure size
Jonatan Werpers <jonatan@werpers.com>
parents: 1
diff changeset
31
43be32298ae2 Add function to get closure size
Jonatan Werpers <jonatan@werpers.com>
parents: 1
diff changeset
32 function closureSize(D::D2)::Int
43
ef060ab3b035 remove stride and remove some bugs
Ylva Rydin <ylva.rydin@telia.com>
parents: 40
diff changeset
33 return length(D.quadratureClosure)
24
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
34 end
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
35
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
36 function readOperator(D2fn, Hfn)
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
37 d = readSectionedFile(D2fn)
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
38 h = readSectionedFile(Hfn)
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
39
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
40 # Create inner stencil
35
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
41 innerStencilWeights = stringToVector(Float64, d["inner_stencil"][1])
24
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
42 width = length(innerStencilWeights)
35
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
43 r = (-div(width,2), div(width,2))
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
44
84
48079bd39969 Change to using tuples in stencils and ops
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
45 innerStencil = Stencil(r, Tuple(innerStencilWeights))
24
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
46
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
47 # Create boundary stencils
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
48 boundarySize = length(d["boundary_stencils"])
68
d485da6e3a77 Make D2 more type stable
Jonatan Werpers <jonatan@werpers.com>
parents: 56
diff changeset
49 closureStencils = Vector{typeof(innerStencil)}() # TBD: is the the right way to get the correct type?
35
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
50
24
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
51 for i ∈ 1:boundarySize
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
52 stencilWeights = stringToVector(Float64, d["boundary_stencils"][i])
35
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
53 width = length(stencilWeights)
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
54 r = (1-i,width-i)
84
48079bd39969 Change to using tuples in stencils and ops
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
55 closureStencils = (closureStencils..., Stencil(r, Tuple(stencilWeights)))
24
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
56 end
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
57
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
58 d2 = D2(
35
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
59 stringToVector(Float64, h["closure"][1]),
24
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
60 innerStencil,
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
61 closureStencils,
35
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
62 stringToVector(Float64, d["e"][1]),
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
63 stringToVector(Float64, d["d1"][1]),
34
bb841977d198 Move stencil operator application to its own function
Jonatan Werpers <jonatan@werpers.com>
parents: 24
diff changeset
64 even
24
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
65 )
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
66
35
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
67 return d2
24
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
68 end
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
69
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
70
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
71 function readSectionedFile(filename)::Dict{String, Vector{String}}
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
72 f = open(filename)
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
73 sections = Dict{String, Vector{String}}()
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
74 currentKey = ""
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
75
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
76 for ln ∈ eachline(f)
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
77 if ln == "" || ln[1] == '#' # Skip comments and empty lines
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
78 continue
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
79 end
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
80
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
81 if isletter(ln[1]) # Found start of new section
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
82 if ~haskey(sections, ln)
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
83 sections[ln] = Vector{String}()
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
84 end
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
85 currentKey = ln
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
86 continue
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
87 end
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
88
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
89 push!(sections[currentKey], ln)
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
90 end
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
91
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
92 return sections
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
93 end
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
94
35
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
95 function stringToVector(T::DataType, s::String)
91e662512e9a Fix bugs in sbpD2
Ylva Rydin <ylva.rydin@telia.com>
parents: 24
diff changeset
96 return T.(eval.(Meta.parse.(split(s))))
24
55fea1ceb6aa Start implementing reading 1D operator stencils from file into struct
Jonatan Werpers <jonatan@werpers.com>
parents: 8
diff changeset
97 end