Mercurial > repos > public > sbplib_julia

annotate src/SbpOperators/volumeops/derivatives/dissipation.jl @ 1096:9f0121e465a5 feature/dissipation_operators

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Docs for undivided_dissipation
author Jonatan Werpers <jonatan@werpers.com>
date Tue, 10 May 2022 21:10:31 +0200
parents 423a6442efc3
children 254934aac3f8
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
1096
9f0121e465a5 Docs for undivided_dissipation
Jonatan Werpers <jonatan@werpers.com>
parents: 1087
diff changeset
1 # REVIEW: Would it be more correct to
1087
423a6442efc3 REVIEW: Address documentation
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 1069
diff changeset
2 # call these undivided_differences instead of dissipation?
423a6442efc3 REVIEW: Address documentation
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 1069
diff changeset
3 # If I understand it correctly, this method simply provides
423a6442efc3 REVIEW: Address documentation
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 1069
diff changeset
4 # the operators required in order to compose a dissipation operator
423a6442efc3 REVIEW: Address documentation
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 1069
diff changeset
5 # and the dissipation operator are formed by a linear combination
423a6442efc3 REVIEW: Address documentation
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 1069
diff changeset
6 # of the products of Dᵀ and D for different orders.
1096
9f0121e465a5 Docs for undivided_dissipation
Jonatan Werpers <jonatan@werpers.com>
parents: 1087
diff changeset
7 """
9f0121e465a5 Docs for undivided_dissipation
Jonatan Werpers <jonatan@werpers.com>
parents: 1087
diff changeset
8 undivided_dissipation(g::EquidistantGrid, p, direction)
9f0121e465a5 Docs for undivided_dissipation
Jonatan Werpers <jonatan@werpers.com>
parents: 1087
diff changeset
9
9f0121e465a5 Docs for undivided_dissipation
Jonatan Werpers <jonatan@werpers.com>
parents: 1087
diff changeset
10 Create undivided difference operators approximating the `p`th derivative for
9f0121e465a5 Docs for undivided_dissipation
Jonatan Werpers <jonatan@werpers.com>
parents: 1087
diff changeset
11 building artificial dissipation. The operators and how they are used to create accurate artifical dissipation is described in
9f0121e465a5 Docs for undivided_dissipation
Jonatan Werpers <jonatan@werpers.com>
parents: 1087
diff changeset
12 "K. Mattsson, M. Svärd, and J. Nordström, “Stable and Accurate Artificial Dissipation,” Journal of Scientific Computing, vol. 21, no. 1, pp. 57–79, Aug. 2004"
9f0121e465a5 Docs for undivided_dissipation
Jonatan Werpers <jonatan@werpers.com>
parents: 1087
diff changeset
13 """
1069
c89c6b63c7f4 Fix typo
Jonatan Werpers <jonatan@werpers.com>
parents: 1035
diff changeset
14 function undivided_dissipation(g::EquidistantGrid, p, direction)
1034
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
15 T = eltype(g)
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
16 interior_weights = T.(dissipation_interior_weights(p))
1021
ee5a641a8277 Add some tests and start implementing dissipation()
Jonatan Werpers <jonatan@werpers.com>
parents: 1020
diff changeset
17
1034
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
18 D = stencil_operator_distinct_closures(
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
19 g,
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
20 dissipation_interior_stencil(interior_weights),
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
21 dissipation_lower_closure_stencils(interior_weights),
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
22 dissipation_upper_closure_stencils(interior_weights),
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
23 direction,
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
24 )
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
25 Dᵀ = stencil_operator_distinct_closures(
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
26 g,
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
27 dissipation_transpose_interior_stencil(interior_weights),
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
28 dissipation_transpose_lower_closure_stencils(interior_weights),
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
29 dissipation_transpose_upper_closure_stencils(interior_weights),
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
30 direction,
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
31 )
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
32
ed19c549c506 Implement undevided dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents: 1033
diff changeset
33 return D, Dᵀ
1021
ee5a641a8277 Add some tests and start implementing dissipation()
Jonatan Werpers <jonatan@werpers.com>
parents: 1020
diff changeset
34 end
ee5a641a8277 Add some tests and start implementing dissipation()
Jonatan Werpers <jonatan@werpers.com>
parents: 1020
diff changeset
35
1069
c89c6b63c7f4 Fix typo
Jonatan Werpers <jonatan@werpers.com>
parents: 1035
diff changeset
36 undivided_dissipation(g::EquidistantGrid{1}, p) = undivided_dissipation(g, p, 1)
1021
ee5a641a8277 Add some tests and start implementing dissipation()
Jonatan Werpers <jonatan@werpers.com>
parents: 1020
diff changeset
37
1020
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
38 function dissipation_interior_weights(p)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
39 if p == 0
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
40 return (1,)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
41 end
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
42
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
43 return (0, dissipation_interior_weights(p-1)...) .- (dissipation_interior_weights(p-1)..., 0)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
44 end
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
45
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
46 midpoint(weights) = length(weights)÷2 + 1
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
47 midpoint_transpose(weights) = length(weights)+1 - midpoint(weights)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
48
1035
ceda69b8f27a Add test for transpose equality and fix bugs found
Jonatan Werpers <jonatan@werpers.com>
parents: 1034
diff changeset
49 function dissipation_interior_stencil(weights)
ceda69b8f27a Add test for transpose equality and fix bugs found
Jonatan Werpers <jonatan@werpers.com>
parents: 1034
diff changeset
50 return Stencil(weights..., center=midpoint(weights))
ceda69b8f27a Add test for transpose equality and fix bugs found
Jonatan Werpers <jonatan@werpers.com>
parents: 1034
diff changeset
51 end
ceda69b8f27a Add test for transpose equality and fix bugs found
Jonatan Werpers <jonatan@werpers.com>
parents: 1034
diff changeset
52 function dissipation_transpose_interior_stencil(weights)
ceda69b8f27a Add test for transpose equality and fix bugs found
Jonatan Werpers <jonatan@werpers.com>
parents: 1034
diff changeset
53 if iseven(length(weights))
ceda69b8f27a Add test for transpose equality and fix bugs found
Jonatan Werpers <jonatan@werpers.com>
parents: 1034
diff changeset
54 weights = map(-, weights)
ceda69b8f27a Add test for transpose equality and fix bugs found
Jonatan Werpers <jonatan@werpers.com>
parents: 1034
diff changeset
55 end
ceda69b8f27a Add test for transpose equality and fix bugs found
Jonatan Werpers <jonatan@werpers.com>
parents: 1034
diff changeset
56
ceda69b8f27a Add test for transpose equality and fix bugs found
Jonatan Werpers <jonatan@werpers.com>
parents: 1034
diff changeset
57 return Stencil(weights..., center=midpoint_transpose(weights))
ceda69b8f27a Add test for transpose equality and fix bugs found
Jonatan Werpers <jonatan@werpers.com>
parents: 1034
diff changeset
58 end
1029
129262c8e897 Change signatures for interior stencil methods
Jonatan Werpers <jonatan@werpers.com>
parents: 1021
diff changeset
59
1020
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
60 dissipation_lower_closure_size(weights) = midpoint(weights) - 1
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
61 dissipation_upper_closure_size(weights) = length(weights) - midpoint(weights)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
62
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
63 dissipation_lower_closure_stencils(interior_weights) = ntuple(i->Stencil(interior_weights..., center=i ), dissipation_lower_closure_size(interior_weights))
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
64 dissipation_upper_closure_stencils(interior_weights) = ntuple(i->Stencil(interior_weights..., center=length(interior_weights)-dissipation_upper_closure_size(interior_weights)+i), dissipation_upper_closure_size(interior_weights))
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
65
1033
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
66 function dissipation_transpose_lower_closure_stencils(interior_weights)
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
67 closure = ntuple(i->dissipation_transpose_lower_closure_stencil(interior_weights, i), length(interior_weights))
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
68
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
69 N = maximum(s->length(s.weights), closure)
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
70 return right_pad.(closure, N)
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
71 end
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
72
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
73 function dissipation_transpose_upper_closure_stencils(interior_weights)
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
74 closure = reverse(ntuple(i->dissipation_transpose_upper_closure_stencil(interior_weights, i), length(interior_weights)))
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
75
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
76 N = maximum(s->length(s.weights), closure)
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
77 return left_pad.(closure, N)
0cb4c6b15d8e Make closures have the same number of weights in all stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1029
diff changeset
78 end
1020
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
79
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
80
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
81 function dissipation_transpose_lower_closure_stencil(interior_weights, i)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
82 w = ntuple(k->interior_weights[i], dissipation_lower_closure_size(interior_weights))
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
83
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
84 for k ∈ i:-1:1
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
85 w = (w..., interior_weights[k])
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
86 end
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
87
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
88 return Stencil(w..., center = i)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
89 end
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
90
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
91 function dissipation_transpose_upper_closure_stencil(interior_weights, i)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
92 j = length(interior_weights)+1-i
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
93 w = ntuple(k->interior_weights[j], dissipation_upper_closure_size(interior_weights))
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
94
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
95 for k ∈ j:1:length(interior_weights)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
96 w = (interior_weights[k], w...)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
97 end
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
98
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
99 return Stencil(w..., center = length(interior_weights)-midpoint(interior_weights)+1)
3f5137ce3aa1 Start adding helper functions for building dissipation operators
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
100 end