Computation of hydraulic balance-ODE fails with non-linear solver

You sometimes need to manually guard against nans in the variational form. Using a PETSc SNES solver with backtracking ("bt") line search can also improve robustness. For example, the following modifications run without error:

#f1 = Constant(0.5)
f1 = Constant(0.01)

#F = v*u.dx(0)*dx+v*sqrt(u)*dx-v*f1*dx
def safeSqrt(x):
    return conditional(gt(x,DOLFIN_EPS),sqrt(x),Constant(0.0))
F = v*u.dx(0)*dx+v*safeSqrt(u)*dx-v*f1*dx

#solve(F == 0, u, bc,solver_parameters={"newton_solver": {
#        "relative_tolerance": 1e-7, "maximum_iterations": 50,
#        "convergence_criterion": "residual",
#        "linear_solver":"gmres", "preconditioner": "jacobi",
#        "krylov_solver":{"nonzero_initial_guess":True}}})
solve(F==0,u,bc,
      solver_parameters
      ={"nonlinear_solver":"snes",
        "snes_solver":{"line_search":"bt"}})
1 Like