Mercurial > repos > public > sbplib_julia

comparison src/SbpOperators/volumeops/laplace/laplace.jl @ 1607:7216448d0c5a feature/boundary_conditions

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
REVIEW: Suggest deduplication of positivity decompostion code
author Jonatan Werpers <jonatan@werpers.com>
date Sun, 09 Jun 2024 00:02:40 +0200
parents 93b86625fcfd
children 8315c456e3b4
comparison
equal deleted inserted replaced
1606:93b86625fcfd 1607:7216448d0c5a
92 92
93 penalty_tensor = -H⁻¹∘e'∘Hᵧ 93 penalty_tensor = -H⁻¹∘e'∘Hᵧ
94 return penalty_tensor, d 94 return penalty_tensor, d
95 end 95 end
96 96
97
98 function positivity_decomposition(Δ::Laplace, g::Grid, bc::DirichletCondition; H_tuning, R_tuning)
99 Nτ_H, τ_R = positivity_limits(Δ,g,bc)
100 return H_tuning*Nτ_H + R_tuning*τ_R
101 end
102
103
97 # TODO: We should consider implementing a proper BoundaryIdentifier for EquidistantGrid and then 104 # TODO: We should consider implementing a proper BoundaryIdentifier for EquidistantGrid and then
98 # change bc::BoundaryCondition to id::BoundaryIdentifier 105 # change bc::BoundaryCondition to id::BoundaryIdentifier
99 106 function positivity_limits(Δ::Laplace, g::EquidistantGrid, bc::DirichletCondition)
100 function positivity_decomposition(Δ::Laplace, g::EquidistantGrid, bc::BoundaryCondition; H_tuning, R_tuning)
101 pos_prop = positivity_properties(Δ) 107 pos_prop = positivity_properties(Δ)
102 h = spacing(g) 108 h = spacing(g)
103 θ_H = pos_prop.theta_H 109 θ_H = pos_prop.theta_H
104 τ_H = H_tuning*ndims(g)/(h*θ_H) 110 τ_H = 1/(h*θ_H)
105 θ_R = pos_prop.theta_R 111 θ_R = pos_prop.theta_R
106 τ_R = R_tuning/(h*θ_R) 112 τ_R = 1/(h*θ_R)
107 B = τ_H + τ_R 113 return τ_H, τ_R
108 return B
109 end 114 end
110 115
111 function positivity_decomposition(Δ::Laplace, g::TensorGrid, bc::BoundaryCondition; H_tuning, R_tuning) 116 function positivity_limits(Δ::Laplace, g::TensorGrid, bc::DirichletCondition)
112 pos_prop = positivity_properties(Δ) 117 τ_H, τ_R = positivity_limits(Δ, g.grids[grid_id(boundary(bc))], bc)
113 h = spacing(g.grids[grid_id(boundary(bc))]) # grid spacing of the 1D grid normal to the boundary 118 return τ_H*ndims(g), τ_R
114 θ_H = pos_prop.theta_H
115 τ_H = H_tuning*ndims(g)/(h*θ_H)
116 θ_R = pos_prop.theta_R
117 τ_R = R_tuning/(h*θ_R)
118 B = τ_H + τ_R
119 return B
120 end 119 end
120
121 121
122 function positivity_properties(Δ::Laplace) 122 function positivity_properties(Δ::Laplace)
123 D2_pos_prop = parse_named_tuple(Δ.stencil_set["D2"]["positivity"]) 123 D2_pos_prop = parse_named_tuple(Δ.stencil_set["D2"]["positivity"])
124 H_closure = parse_tuple(Δ.stencil_set["H"]["closure"]) 124 H_closure = parse_tuple(Δ.stencil_set["H"]["closure"])
125 return merge(D2_pos_prop, (theta_H = H_closure[1],)) 125 return merge(D2_pos_prop, (theta_H = H_closure[1],))