Mercurial > repos > public > sbplib_julia

comparison DiffOps/src/laplace.jl @ 246:3d83b4d78b55 boundary_conditions

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Add methods to Laplace type for creating boundary value and normal derivative operators
author Jonatan Werpers <jonatan@werpers.com>
date Thu, 27 Jun 2019 14:33:39 +0200
parents d9e262cb2e8d
children 89a101a63e7a
comparison
equal deleted inserted replaced
245:d9e262cb2e8d 246:3d83b4d78b55
66 66
67 struct Laplace{Dim,T<:Real,N,M,K} <: DiffOpCartesian{Dim} 67 struct Laplace{Dim,T<:Real,N,M,K} <: DiffOpCartesian{Dim}
68 grid::EquidistantGrid{Dim,T} 68 grid::EquidistantGrid{Dim,T}
69 a::T 69 a::T
70 op::D2{Float64,N,M,K} 70 op::D2{Float64,N,M,K}
71 # e::BoundaryValue
72 # d::NormalDerivative
73 end 71 end
72
73 boundary_value(L::Laplace, bId::CartesianBoundary) = BoundaryValue(L.op, L.grid, bId)
74 normal_derivative(L::Laplace, bId::CartesianBoundary) = NormalDerivative(L.op, L.grid, bId)
75 boundary_quadrature(L::Laplace, bId::CartesianBoundary) = throw(MethodError) # TODO: Implement this
74 76
75 function apply(L::Laplace{Dim}, v::AbstractArray{T,Dim} where T, I::CartesianIndex{Dim}) where Dim 77 function apply(L::Laplace{Dim}, v::AbstractArray{T,Dim} where T, I::CartesianIndex{Dim}) where Dim
76 error("not implemented") 78 error("not implemented")
77 end 79 end
78 80
101 103
102 104
103 struct Neumann{Bid<:BoundaryIdentifier} <: BoundaryCondition end 105 struct Neumann{Bid<:BoundaryIdentifier} <: BoundaryCondition end
104 106
105 function sat(L::Laplace{2,T}, bc::Neumann{Bid}, v::AbstractArray{T,2}, g::AbstractVector{T}, I::CartesianIndex{2}) where {T,Bid} 107 function sat(L::Laplace{2,T}, bc::Neumann{Bid}, v::AbstractArray{T,2}, g::AbstractVector{T}, I::CartesianIndex{2}) where {T,Bid}
106 e = BoundaryValue(L.op, L.grid, Bid()) 108 e = boundary_value(L.op, Bid())
107 d = NormalDerivative(L.op, L.grid, Bid()) 109 d = normal_derivative(L.op, Bid())
108 Hᵧ = BoundaryQuadrature(L.op, L.grid, Bid()) 110 Hᵧ = boundary_quadrature(L.op, Bid())
109 # TODO: Implement BoundaryQuadrature method
110 111
111 return -L.Hi*e*Hᵧ*(d'*v - g) 112 return -L.Hi*e*Hᵧ*(d'*v - g)
112 # Need to handle d'*v - g so that it is an AbstractArray that TensorMappings can act on 113 # Need to handle d'*v - g so that it is an AbstractArray that TensorMappings can act on
113 end 114 end
114 115
115 struct Dirichlet{Bid<:BoundaryIdentifier} <: BoundaryCondition 116 struct Dirichlet{Bid<:BoundaryIdentifier} <: BoundaryCondition
116 tau::Float64 117 tau::Float64
117 end 118 end
118 119
119 function sat(L::Laplace{2,T}, bc::Dirichlet{Bid}, v::AbstractArray{T,2}, g::AbstractVector{T}, i::CartesianIndex{2}) where {T,Bid} 120 function sat(L::Laplace{2,T}, bc::Dirichlet{Bid}, v::AbstractArray{T,2}, g::AbstractVector{T}, i::CartesianIndex{2}) where {T,Bid}
120 e = BoundaryValue(L.op, L.grid, Bid()) 121 e = boundary_value(L.op, Bid())
121 d = NormalDerivative(L.op, L.grid, Bid()) 122 d = normal_derivative(L.op, Bid())
122 Hᵧ = BoundaryQuadrature(L.op, L.grid, Bid()) 123 Hᵧ = boundary_quadrature(L.op, Bid())
123 # TODO: Implement BoundaryQuadrature method
124 124
125 return -L.Hi*(tau/h*e + d)*Hᵧ*(e'*v - g) 125 return -L.Hi*(tau/h*e + d)*Hᵧ*(e'*v - g)
126 # Need to handle scalar multiplication and addition of TensorMapping 126 # Need to handle scalar multiplication and addition of TensorMapping
127 end 127 end
128 128