-
Notifications
You must be signed in to change notification settings - Fork 2
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Starting to draft primal HHO #20
base: main
Are you sure you want to change the base?
Conversation
Codecov ReportAttention: Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #20 +/- ##
==========================================
+ Coverage 81.80% 83.93% +2.13%
==========================================
Files 14 18 +4
Lines 2022 2595 +573
==========================================
+ Hits 1654 2178 +524
- Misses 368 417 +49 ☔ View full report in Codecov by Sentry. |
draft_hho.jl
Outdated
a(u,v)=∫(∇(v)⋅∇(u))dΩ | ||
l( (uK,u∂K), v)=∫((∇(v)⋅nK)*u∂K)d∂K+∫(-Δ(v)*uK)dΩ | ||
|
||
op=LocalAffineFEOperator((a,l),UKR,VKR) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
# Define HHO local reconstruction operator
m(ur,v) = ∫∇(ur)*∇(v)dΩ
l( (uK,u∂K), v)=∫((∇(v)⋅nK)*u∂K)d∂K+∫(-Δ(v)*uK)dΩ
P=LocalProjector((m,l),UKR,VKR)
# Global formulation
UK∂K = MultiFieldFESpace([UK, U∂K])
a(u,v)=∫(∇(P(v))⋅∇(P(u)))dΩ
function f(v)
vK, v∂K = UK∂K
f(v) = ∫(f*vK)dΩ
end
# as HDG
draft_hho.jl
Outdated
|
||
a(u,v)=∫(∇(v)⋅∇(u))dΩ | ||
l( (uK,u∂K), v)=∫((∇(v)⋅nK)*u∂K)d∂K+∫(-Δ(v)*uK)dΩ | ||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
# Stabilisation
# Method 1: Orthogonal basis using change of basis with M^{-1}
refferecᵤ = ReferenceFE(orthogonal,Float64,order+1;space=:P)
\mu_T = h^{-2} # k^4
m(u,v)=∫(∇(v)⋅∇(u))dΩ + \mu_T ∫(v⋅u)dΩ
l( (uK,u∂K), v)=∫((∇(v)⋅nK)*u∂K)d∂K+∫(-Δ(v)*uK)dΩ + \mu_T ∫(v⋅u∂K)dΩ
# Optimisation
# \pi(k,u) that makes 0 dofs higher than k
# Method 2: Define subspaces of orthogonal spaces
# orthogonal basis from 1 to order+1 (excluding 0)
refferecᵤ = ReferenceFE(orthogonal,Float64,1:order+1;space=:P)
refferecᵤ_0 = ReferenceFE(orthogonal,Float64,0;space=:P)
VKR_1 = ...
VKR_0 = ...
VKR = MultiFieldFESpace([UKR_0,UKR_1])
# Define HHO local reconstruction operator
function m(u,v )
v0, v1 = v
m(u,v) = ∫∇(u)*∇(v1)dΩ + ∫(u*v0) dΩ
end
function l(u,v )
v0, v1 = v
l( (uK,u∂K), v)=∫((∇(v1)⋅nK)*u∂K)d∂K+∫(-Δ(v1)*uK)dΩ + ∫(uK*v0) dΩ
end
P=LocalProjector((m,l),UKR,VKR)
m(ur,v) = ∫u_mean*vdΩ
l(uK, v)=∫v*uKdΩ
\pi_0=LocalProjector((m,l),UKR,VKR)
# Global formulation
UK∂K = MultiFieldFESpace([UK, U∂K])
a(u,v)=∫(∇(P(v))⋅∇(P(u)))dΩ
function f(v)
vK, v∂K = UK∂K
f(v) = ∫(f*vK)dΩ
end
#
m(u,v)=∫(∇(v)⋅∇(u))dΩ
l( (uK,u∂K), v)=∫((∇(v)⋅nK)*u∂K)d∂K+∫(-Δ(v)*uK)dΩ
∫(v⋅u∂K)dΩ
draft_hho.jl
Outdated
l( (uK,u∂K), v)=∫((∇(v)⋅nK)*u∂K)d∂K+∫(-Δ(v)*uK)dΩ | ||
|
||
op=LocalAffineFEOperator((a,l),UKR,VKR) | ||
UKR_basis=op(get_trial_fe_basis.((UK,U∂K))...) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
# Stab in Badia et al 2.11
# Define spaces for the projectors
# Cell term
m(ur,v) = ∫u*vdΩ
l( (uK,u∂K), v)=∫(P((uK,u∂K)) - uK)*vdΩ
\PiK =LocalProjector((m,l),UKR,VKR)
# Face term
m(ur,v) = ∫u*vd∂K
l( (uK,u∂K), v)=∫(P((uK,u∂K)) - u∂K)*vd∂K
\Pi∂K =LocalProjector((m,l),UKR,VKR)
# stab term in the global operator
a(u,v)=∫(∇(P(v))⋅∇(P(u)))dΩ + ∫\Pi∂K(u)*\Pi∂K(v)d∂K + ∫\PiK(u)*\PiK(v)dΩ
(i.e., cell-wise + cell-boundary-wise L2 projection of bulk FE functions)
I took a k+1-order polynomial from a H1-conforming space, reduced to the k-order hybrid space, and reconstructed it again to the k+1-th order non-conforming space. The first and the last match.
TO-DEBUG: Currently it does not solve a problem with manufactured solution in the solution FE space.
* Added stabilization
It does not work for order=2 nor order=4, though, :-(((. No idea why not at this point. Also, the code is currently quite slow, I do not know why yet, but clearly needs to be investigated.
ProjectionFEOperator
form associated to ProjectionFEOperator. As of now, the implementation seems to be working for any polynomial order with analytical solution living in the FE space * Developed functions required to evaluate the residual of the HH0 method.
* Adding tests with manufactured solutions in FE space up to order 3.
WIP ...