Skip to content
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

MultiFlield Boundary Condition Problem #947

Closed
DaniStauber opened this issue Oct 17, 2023 · 3 comments
Closed

MultiFlield Boundary Condition Problem #947

DaniStauber opened this issue Oct 17, 2023 · 3 comments

Comments

@DaniStauber
Copy link

Hi there! Thanks for the amazing work that you been doing! Im using Gridap to try to solve the stationary Einstein equations in RG.
For now, i found this problem and i don't really understand it, and less how to solve it.
My toy problem is $$\Delta u= f_1$$ in a circle $\Omega$ with the boundary condition
$$\partial_{\theta} ^2 u = f_2$$ in $\partial \Omega$.
Using one exact solution $u= y^2 -x^2$ with $f_1 =0$ and $f_2 = 4(x^2 -y ^2)$.

So the first thing that i do was solve the boundary condition problem, in the Boundary FESpace.

Ω = Triangulation(model)
degree = 3
dΩ = Measure(Ω,degree)  
Γ= BoundaryTriangulation(model,tags="ext")#(model,tags="TotalBoundary")
Γ_model=Triangulation(Γ)
dΓ = Measure(Γ,degree) 
order=2
reffe_Γ = ReferenceFE(lagrangian,Float64,order)
V_Γ= TestFESpace(Γ_model,reffe_Γ;conformity=:H1,constraint=:zeromean)

and i got the solution with an absolute error $\sim e^{-12}$. The problem is when i tryied to mix the two problems. The weak formulation for the complete problem is $$\int_{\Omega}( \nabla v \cdot \nabla u) - \int_{\partial \Omega} ((\hat{r} \cdot \nabla u) v)+ \int_{\partial \Omega}(- (\hat{\theta} \cdot \nabla w )(\hat{\theta} \cdot \nabla t) - w r (\hat{r} \cdot \nabla t) ) + \int_{\partial \Omega} (v(u-t)) = \int_{\Omega} (-f_1 v) +\int_{\partial \Omega} (w f_2)$$
$t$ is the solution of the equation in the boundary and $w$ its test function, $u$ is the solution in the circle and $v$ its test function. The first two terms on the left is the weak formulation for the laplacian equation, the third terms is the weak formulation of my boundary condition, and the las term is the weak Dirichlet contion that means $u|_{\partial \Omega} = t$. On the rigth we have the sources.

On the code, i gave my FESpaces

V_Ω= TestFESpace(Ω,reffe;conformity=:H1)
U_Γ= TrialFESpace(V_Γ)
U_Ω= TrialFESpace(V_Ω)
V_T= MultiFieldFESpace([V_Ω, V_Γ]) #Test space
U_T = MultiFieldFESpace([U_Ω, U_Γ])

fuente(x)=0
partialT_θθ(x)=4. *(x[1]^2 -x[2]^2)/R 
T_rθ(x)=x[2]^2 -x[1]^2 
a_T((u,t),(v,w))= ∫(  ∇(u) ⋅ ∇(v) )*dΩ -∫( v * (n ⋅ ∇(u)) )*dΓ +∫(u*v)*dΓ - ∫( v * t )*dΓ+ ∫( -(θvec ⋅ ∇(t)) *(θvec ⋅ ∇(w)) -w * (rvec ⋅∇(t)) *R )*dΓ

b_T((v,w))=∫(- v * fuente )*dΩ +∫( R*(partialT_θθ*w) )*dΓ
op_T= AffineFEOperator(a_T,b_T,U_T,V_T)
trθ_T =solve(solver,op_T)

And i got this error.

AssertionError: A check failed

Stacktrace:
  [1] macro expansion
    @ ~/.julia/packages/Gridap/MTWD6/src/Helpers/Macros.jl:60 [inlined]
  [2] get_triangulation(f::Gridap.MultiField.MultiFieldCellField{ReferenceDomain})
    @ Gridap.MultiField ~/.julia/packages/Gridap/MTWD6/src/MultiField/MultiFieldCellFields.jl:34
  [3] get_triangulation
    @ ~/.julia/packages/Gridap/MTWD6/src/MultiField/MultiFieldFEFunctions.jl:30 [inlined]
  [4] num_cells(a::Gridap.MultiField.MultiFieldFEFunction{Gridap.MultiField.MultiFieldCellField{ReferenceDomain}})
    @ Gridap.CellData ~/.julia/packages/Gridap/MTWD6/src/CellData/CellDataInterface.jl:59
  [5] show(io::IOContext{IOBuffer}, #unused#::MIME{Symbol("text/plain")}, f::Gridap.MultiField.MultiFieldFEFunction{Gridap.MultiField.MultiFieldCellField{ReferenceDomain}})
    @ Gridap.CellData ~/.julia/packages/Gridap/MTWD6/src/CellData/CellFields.jl:71
  [6] limitstringmime(mime::MIME{Symbol("text/plain")}, x::Gridap.MultiField.MultiFieldFEFunction{Gridap.MultiField.MultiFieldCellField{ReferenceDomain}}, forcetext::Bool)
    @ IJulia ~/.julia/packages/IJulia/Vo51o/src/inline.jl:43
  [7] limitstringmime
    @ ~/.julia/packages/IJulia/Vo51o/src/inline.jl:38 [inlined]
  [8] display_mimestring
    @ ~/.julia/packages/IJulia/Vo51o/src/display.jl:71 [inlined]
  [9] display_dict(x::Gridap.MultiField.MultiFieldFEFunction{Gridap.MultiField.MultiFieldCellField{ReferenceDomain}})
    @ IJulia ~/.julia/packages/IJulia/Vo51o/src/display.jl:102
 [10] #invokelatest#2
    @ ./essentials.jl:819 [inlined]
 [11] invokelatest
    @ ./essentials.jl:816 [inlined]
 [12] execute_request(socket::ZMQ.Socket, msg::IJulia.Msg)
    @ IJulia ~/.julia/packages/IJulia/Vo51o/src/execute_request.jl:112
 [13] #invokelatest#2
    @ ./essentials.jl:819 [inlined]
 [14] invokelatest
    @ ./essentials.jl:816 [inlined]
 [15] eventloop(socket::ZMQ.Socket)
    @ IJulia ~/.julia/packages/IJulia/Vo51o/src/eventloop.jl:8
 [16] (::IJulia.var"#15#18")()
    @ IJulia ./task.jl:514

I dont understand the error, and how to fix it. I tryied to solve the laplacian equation with weak boundari condition with the $ U_{\Omega}$ and $V_{\Omega}$ in this way, and it works with an absolute error of $e^{-10}$

a_Ts(u,v)= ∫(  ∇(u) ⋅ ∇(v) )*dΩ -∫( v * (n ⋅ ∇(u)) )*dΓ +∫(u*v)*dΓ 
b_Ts(v)=∫(- v * fuente )*dΩ +∫( v * exact )*dΓ

So the problems separately works fine. But when i tried to solve it together, that error comes out.
I provide here the toy code with the mesh in case that you want to reproduce the error
MWE_Issue.zip

I would appreciate that someone gives me an idea of what is going on, or what is the problem. Thanks!!!

@JordiManyer
Copy link
Member

JordiManyer commented Oct 17, 2023

@DaniStauber That seems a display error, related to the fact that the MultiFieldFESpace has two different triangulations. In that situation, get_triangulation is obviously not defined.
In any case, that is a display error. The variable should be just fine and everything is working. You can put a semicolon at the end of the statement if it bothers you.

@DaniStauber
Copy link
Author

Thank you for your answer!!. Yes, that thing works, thank you. But I really want to understand it, and knows if my problem is affected by the error

@JordiManyer
Copy link
Member

There is no error, just a display error.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants