NonlinearVariationalSolver does not pass linear solver parameters

Issue #10 resolved
Jan Blechta created an issue

NonlinearVariationalSolver has lu_solver and krylov_solver parameter sets in its parameters but it does not pass them to linear solver created by NewtonSolver.

Comments (29)

  1. Anders Logg (Chalmers)

    A fix has been implemented in logg/fix-issue-10 but testing is needed.

    This involves a change in the NewtonSolver interface whereby the variables linear_solver and preconditioner are changed into parameters. In addition, the lu_solver and krylov_solver parameters are moved from NonlinearVariationalProblem to NewtonSolver where they belong more naturally.

    The interface change might break some demos so I'll do some more testing before merging.

  2. Jan Blechta reporter

    @logg Take care not to break PETScSNESSolver. I think that @cmaurini or @pefarrell made a delivery of these parameters to PETScSNESSolver working so it would be silly to fix one and break another.

  3. Anders Logg (Chalmers)

    I see that now. Corrado and Patrick, could you move the linear solver parameters into PETScSNESSolver? That should make things simpler for you as well. You would not need to extract the parameters inside set_linear_solver_parameters and call that from NonlinearVariationalSolver.

    PS: Jan, how does one link to a user when posting on Bitbucket?

  4. Jan Blechta reporter

    Type at sign @ and then continue typing some civil name or user name and bitbucket will help you choose. Unfortunately it acts badly on usernames with dot, like Garth's one.

  5. Jan Blechta reporter

    @logg will be still those parameters accessible somehow when using high-level NonlinearVariationalSolver rather than NewtonSolver or PETScSNESSolver?

  6. Anders Logg (Chalmers)

    Yes @blechta, they will be accessible via the parameters of the Newton solver:

    nonlinear_variational_solver.parameters["newton_solver"]["linear_solver"]["..."]
    
  7. Jan Blechta reporter
    • changed status to open

    Maybe there is still something wrong. Check warnings below.

    >python demo_hyperelasticity.py 
    Solving nonlinear variational problem.
      *** Warning: Unable to set parameters for linear solver type "default".
      Newton iteration 0: r (abs) = 3.681e+00 (tol = 1.000e-10) r (rel) = 1.000e+00 (tol = 1.000e-09)
    ...
    
    >python demo_nonlinear-poisson.py 
    No Jacobian form specified for nonlinear variational problem.
    Differentiating residual form F to obtain Jacobian J = F'.
    Solving nonlinear variational problem.
      *** Warning: Unable to set parameters for linear solver type "default".
      Newton iteration 0: r (abs) = 5.745e+00 (tol = 1.000e-10) r (rel) = 1.000e+00 (tol = 1.000e-06)
      Solving linear system of size 1089 x 1089 using Epetra LU solver (Amesos_Klu).
    ...
    
    >python ConstrainedFunctionSpace.py 
    
    Testing constrained (periodic) FunctionSpace
    ------------------------------------------------
    ......Solving linear variational problem.
    No Jacobian form specified for nonlinear variational problem.
    Differentiating residual form F to obtain Jacobian J = F'.
    Solving nonlinear variational problem.
      *** Warning: Unable to set parameters for linear solver type "default".
      Newton iteration 0: r (abs) = 5.871e-02 (tol = 1.000e-10) r (rel) = 1.000e+00 (tol = 1.000e-09)
      Newton iteration 1: r (abs) = 1.547e-16 (tol = 1.000e-10) r (rel) = 2.635e-15 (tol = 1.000e-09)
      Newton solver finished in 1 iterations and 1 linear solver iterations.
    ..
    ----------------------------------------------------------------------
    Ran 8 tests in 0.348s
    
    OK
    

    And also in demo_auto-adaptive-navier-stokes.py and chapter_10.py test. Common denominator is solve(F == 0, ...).

  8. Anders Logg (Chalmers)

    Yes, this is a problem related to both NewtonSolver and LinearSolver needing to agree on which solver types are supported and what is meant by "default". I'll take a closer look.

  9. Anders Logg (Chalmers)

    I've fixed this now (again...). The naming of linear solvers (direct, iterative, krylov, lu, gmres etc) is somewhat of a jungle and the parameter sets (krylov_solver, lu_solver) must match and do not really correspond to the solvers. This should be discussed on the mailing list and cleaned up (after the release of 1.3).

  10. Jan Blechta reporter

    Yes, there are plenty of defficienscis, for example:

    When DOLFIN is not compiled with CHOLMOD, then LinearSolver('cholesky') redirects to CholmodCholeskySolver which implements a fallback to UmfpackLUSolver. It should be discussed if LinearSolver('cholesky') should call rather LUSolver with symmetric=True option.

    In PETScLUSolver::select_solver there is missing a logic for picking Cholesky-capable (as indicated by PETScLUSolver::_methods_cholesky) solver when symmetric=True.

  11. Chaffra Affouda

    Is there any parlicular reason _solver is only created when you call solve(...) on NewtonSolver() object and not in the constructor? I think this leads to seg faults.

    solver = NewtonSolver()
    solver.parameters['linear_solver'] = 'lu'
    ddm_lu = solver.linear_solver().parameters
    

    Any access to ddm_lu parameters will give seg fault.

  12. Jan Blechta reporter
    • changed status to open

    Probably, an error message should be raised in NewtonSolver::linear_solver() when no _solver has been created yet.

    Possibly the construction of the linear solver in NewtonSolver::solve could be moved to init method so that user can use approach linear solver before solve is called.

  13. Anders Logg (Chalmers)

    This seems to have already been fixed. Accessing the linear solver checks for the solver being instantiated:

    GenericLinearSolver& NewtonSolver::linear_solver() const
    {
        dolfin_assert(_solver);
        return *_solver;
    }
    

    Regardless of this, accessing the parameters of linear solver object won't work. The parameters must be set in the nested parameter set for the Newton solver. These are then passed to the linear solver when solve() is called. Changing the parameters for the linear solver object (if it had been instantiated) would have now effect as the parameters would be overwritten:

    _solver->update_parameters(parameters(_solver->parameter_type()));
    
  14. Jan Blechta reporter

    dolfin_assert is only consistency check in DOLFIN library which has no effect for production builds. Proper check when user might be doing wrong thing is if (!_solver) dolfin_error(...).

    Initialization on first use could be better.

  15. Anders Logg (Chalmers)

    Then if we want to be consistent, we would need to change a lot of dolfin_asserts throughout the code to if + sensible error messages.

    Do we want to go that route? Opinions?

    In this case, initialization on first use would trick the user into believing that the returned object can actually be used for setting for example parameters (which it can't even if it were initialized).

  16. Jan Blechta reporter

    To my understanding, asserts are for checking DOLFIN consistency, not user input.

    The proper way for tweaking linear solver should be NewtonSolver(solver, factory) ctor which is now broken - see #311. But it depends how #311 will be fixed - whether update_parameters will overwrite or not.

  17. Jan Blechta reporter

    Great!

    Does

    _solver->update_parameters(parameters(_solver->parameter_type()));
    

    rewrite parameters of solver passed to constructor of Newton solver?

  18. Log in to comment