fix for variable h

zjwegert · Dec 9, 2024 · e810988 · e810988
1 parent bd474f8
commit e810988
Show file tree

Hide file tree

Showing 3 changed files with 6 additions and 5 deletions.
diff --git a/scripts/Embedded/Examples/fsi/gmsh/GMSH_TO-5-Brinkmann_stokes_P1-P1_Ersatz_elast_fsi.jl b/scripts/Embedded/Examples/fsi/gmsh/GMSH_TO-5-Brinkmann_stokes_P1-P1_Ersatz_elast_fsi.jl
@@ -29,6 +29,7 @@ writevtk(model,path*"model")
 
 Ω_act = Triangulation(model)
 hₕ = CellField(get_element_sizes(model),Ω_act)
+hmin = GridapTopOpt._get_minimum_element_diameter(hₕ)
 
 # Cut the background model
 reffe_scalar = ReferenceFE(lagrangian,Float64,1)
@@ -107,7 +108,7 @@ cl = H # Characteristic length
 u0_max = maximum(abs,get_dirichlet_dof_values(U))
 μ = ρ*cl*u0_max/Re # Viscosity
 # Stabilization parameters
-γ(h) = 100.0*μ/h
+γ(h) = 1e5*μ/h
 β1(h) = 1/3*(h^2/4μ)
 
 # Terms
@@ -118,7 +119,7 @@ c_Ω(p,q) = (β1 ∘ hₕ)*∇(p)⋅∇(q)
 
 a_fluid((u,p),(v,q)) =
   ∫( a_Ω(u,v)+b_Ω(u,q)+b_Ω(v,p)-c_Ω(p,q)) * Ω.dΩf +
-  ∫( a_Ω(u,v)+b_Ω(u,q)+b_Ω(v,p)-c_Ω(p,q)+γ*u⋅v ) * Ω.dΩs
+  ∫( a_Ω(u,v)+b_Ω(u,q)+b_Ω(v,p)-c_Ω(p,q)+(γ ∘ hₕ)*(u⋅v)) * Ω.dΩs
 
 ## Structure
 # Stabilization and material parameters
@@ -166,7 +167,7 @@ vel_ext = VelocityExtension(a_hilb,U_reg,V_reg)
 ## Optimiser
 converged(m) = GridapTopOpt.default_al_converged(
   m;
-  L_tol = 0.5h,
+  L_tol = 0.5*hmin,
   C_tol = 0.01vf
 )
 function has_oscillations(m,os_it)

diff --git a/src/LevelSetEvolution/UnfittedEvolution/CutFEMEvolve.jl b/src/LevelSetEvolution/UnfittedEvolution/CutFEMEvolve.jl
@@ -58,7 +58,7 @@ function get_transient_operator(φh,velh,s::CutFEMEvolve)
   v_norm = maximum(abs,get_free_dof_values(velh))
   β(vh,∇φ) = vh/(ϵ + v_norm) * ∇φ/(ϵ + norm(∇φ))
 
-  stiffness(t,u,v) = ∫(((β ∘ (velh,∇(φh))) ⋅ ∇(u)) * v)dΩ_bg + ∫(γg*h^2*jump(∇(u) ⋅ n_Γg)*jump(∇(v) ⋅ n_Γg))dΓg
+  stiffness(t,u,v) = ∫(((β ∘ (velh,∇(φh))) ⋅ ∇(u)) * v)dΩ_bg + ∫((γg*h*h)*jump(∇(u) ⋅ n_Γg)*jump(∇(v) ⋅ n_Γg))dΓg
   mass(t, ∂ₜu, v) = ∫(∂ₜu * v)dΩ_bg
   forcing(t,v) = ∫(0v)dΩ_bg + ∫(0*jump(∇(v) ⋅ n_Γg))dΓg
   # Second term is added to address the following issue:

diff --git a/src/LevelSetEvolution/UnfittedEvolution/StabilisedReinit.jl b/src/LevelSetEvolution/UnfittedEvolution/StabilisedReinit.jl
@@ -135,7 +135,7 @@ function get_residual_and_jacobian(s::StabilisedReinit{InteriorPenalty},φh)
   dΓ = Ωs.dΓ
 
   W(u,∇u) = sign(u) * ∇u / (ϵ + norm(∇u))
-  a(w,u,v) = ∫(v*(W ∘ (w,∇(w))) ⋅ ∇(u))dΩ_bg + ∫(h^2*jump(∇(u)) ⋅ jump(∇(v)))dΛ + ∫((γd/h)*v*u)dΓ
+  a(w,u,v) = ∫(v*(W ∘ (w,∇(w))) ⋅ ∇(u))dΩ_bg + ∫(h*h*jump(∇(u)) ⋅ jump(∇(v)))dΛ + ∫((γd/h)*v*u)dΓ
   b(w,v) = ∫((sign ∘ w)*v)dΩ_bg
   res(u,v) = a(u,u,v) - b(u,v)
   jac(u,du,v) = a(u,du,v)