We are having some difficulties solving two coupled multi-variate non-linear PDEs using Gridap. The various tutorials haven’t quite filled the gaps in our knowledge.
We have two PDEs that we need to solve together. Their residuals are:
res_PF(s,ϕ) = ∫( Gc(tags)*ls*∇(ϕ)⋅ ∇(s) + 2*ψPlusPrev_in*s*ϕ + (Gc(tags)/ls)*s*ϕ )*dΩ - ∫( (Gc(tags)/ls)*ϕ )*dΩ
res_Disp(u,v) = ∫( (ε(v) ⊙ (σ_mod∘(ε(u),ε(uh_in),sh_in)) ) )*dΩ
Note that uh_in and sh_in are the fields at the end of the last increment.
(We have Jacobians too, but with the automatic differentation, let’s not complicate matters yet!)
As per the tutorials, to solve one non-linear PDE, we can use
op_Disp = FEOperator(res_Disp, U_Disp, V0_Disp)
nls_Disp = NLSolver(show_trace=true, method=:newton, linesearch=BackTracking(), ftol = ftol, iterations = 4)
solver_Disp = FESolver(nls_Disp)
uh_out, = solve!(uh_in, solver_Disp, op_Disp)
where U_Disp and V0_Disp are the trial and test spaces, respectively.
Referring to Tutorial 11 on FSI, it suggests that we need to add our residuals, and solve like this:
residual((u,s),(v,ϕ)) = ∫( (ε(v) ⊙ (σ_mod∘(ε(u),ε(uh_in),sh_in)) ) )*dΩ + ∫( Gc(tags)*ls*∇(ϕ)⋅ ∇(s) + 2*ψPlusPrev_in*s*ϕ + (Gc(tags)/ls)*s*ϕ )*dΩ - ∫( (Gc(tags)/ls)*ϕ )*dΩ
op = FEOperator(residual, MultiFieldFESpace([U_Disp; U_PF]), MultiFieldFESpace([V0_Disp; V0_PF]))
nls = NLSolver(show_trace=true, method=:newton, linesearch=BackTracking(), ftol = ftol, iterations = 5)
solver = FESolver(nls)
out = solve(solver, op)
This at least runs; however, it prompts two questions:
-
Is it correct to add residuals in this instance? I was expecting to solve the two residual functions in a 2D Newton’s method. If not, what is the correct way, please?
-
If it is the correct way, how can we specify the two initial fields (i.e. from the end of the last increment)? We can’t get any variation of the following to work:
out = solve([uh_in, sh_in], solver, op)
Many thanks!