Gridap Solvers

[JCHostos]

Hi everyone,

I am currently studying the following example to use GridapSolvers.jl in the solution of Navier-Stokes equations.

https://gridap.github.io/GridapSolvers.jl/stable/Examples/NavierStokes/

I was wondering if there is an appropriate iterative solver for the velocity block (solver_u), instead of using the direct solver LUSolver, in this simple sequential example.

.......Π_Qh = LocalProjectionMap(divergence,Q,qdegree)
graddiv(u,v,dΩ) = ∫(α*(∇⋅v)⋅Π_Qh(u))dΩ

conv(u,∇u) = (∇u')⋅u
dconv(du,∇du,u,∇u) = conv(u,∇du)+conv(du,∇u)
c(u,v,dΩ) = ∫(v⊙(conv∘(u,∇(u))))dΩ
dc(u,du,dv,dΩ) = ∫(dv⊙(dconv∘(du,∇(du),u,∇(u))))dΩ

lap(u,v,dΩ) = ∫(ν*∇(v)⊙∇(u))dΩ
rhs(v,dΩ) = ∫(v⋅f)dΩ

jac_u(u,du,dv,dΩ) = lap(du,dv,dΩ) + dc(u,du,dv,dΩ) + graddiv(du,dv,dΩ)
jac((u,p),(du,dp),(dv,dq),dΩ) = jac_u(u,du,dv,dΩ) - ∫(divergence(dv)*dp)dΩ - ∫(divergence(du)*dq)dΩ

res_u(u,v,dΩ) = lap(u,v,dΩ) + c(u,v,dΩ) + graddiv(u,v,dΩ) - rhs(v,dΩ)
res((u,p),(v,q),dΩ) = res_u(u,v,dΩ) - ∫(divergence(v)*p)dΩ - ∫(divergence(u)*q)dΩ

Ω = Triangulation(model)
dΩ = Measure(Ω,qdegree)
jac_h(x,dx,dy) = jac(x,dx,dy,dΩ)
res_h(x,dy) = res(x,dy,dΩ)
op = FEOperator(res_h,jac_h,X,Y)

solver_u = LUSolver()
solver_p = CGSolver(JacobiLinearSolver();maxiter=20,atol=1e-14,rtol=1.e-6,verbose=i_am_main(parts))
solver_p.log.depth = 4

bblocks = [NonlinearSystemBlock() LinearSystemBlock();
LinearSystemBlock() BiformBlock((p,q) -> ∫(-(1.0/α)pq)dΩ,Q,Q)]
coeffs = [1.0 1.0;
0.0 1.0]
P = BlockTriangularSolver(bblocks,[solver_u,solver_p],coeffs,:upper)
solver = FGMRESSolver(20,P;atol=1e-11,rtol=1.e-8,verbose=i_am_main(parts))
solver.log.depth = 2

nlsolver = NewtonSolver(solver;maxiter=20,atol=1e-10,rtol=1.e-12,verbose=i_am_main(parts))
xh = solve(nlsolver,op); ........

[amartinhuertas]

have you taken a look at https://gridap.github.io/GridapSolvers.jl/stable/Examples/NavierStokesGMG/?

[JCHostos]

Yes, I actually did, but I realized that It uses mh = CartesianModelHierarchy(parts,np_per_level,domain,nc;add_labels! = add_labels!), which I understand is for models made with the Gridap Discretization tool, and not for an external mesh generated with Gmsh. am I wrong?.

Another thing is that it seems that the Navier-Stokes GMG version seems to be intended for parallel distributed computations, and I want to star with a more simple serial version (since I am still a begginer using Gridap, although I have already understood really well and solved my own problems starting from the Gridap tutorials).

Thank you very much

[JordiManyer]

You can refine models coming from gmsh and create the model hierarchy by hand, or rely on https://github.com/gridap/GridapP4est.jl (i.e like this) .

[JordiManyer]

Just to clarify:

• If you do not use the Augmented Lagrangian term, you can use any solver that works for a Laplacian matrix. That gives you many options in PETSc, for instance.
• If you use the AL term, however, you need a solver that is robust wrt the kernel that the divergence operator introduces. GMG is the best candiate, in my opinion, but there are others (that we do not have implemented).

[JCHostos]

Thank you Jordi,

It has been really useful.

Regards