Poster of Linux kernelThe best gift for a Linux geek
petsc_nonlinear_solver.C

petsc_nonlinear_solver.C

Section: C Library Functions (3) Updated: Thu Apr 7 2011
Local index Up
 

NAME

petsc_nonlinear_solver.C -  

SYNOPSIS


 

Typedefs


typedef int PetscErrorCode

typedef int PetscInt
 

Functions


PetscErrorCode __libmesh_petsc_snes_monitor (SNES, PetscInt its, PetscReal fnorm, void *)

PetscErrorCode __libmesh_petsc_snes_residual (SNES, Vec x, Vec r, void *ctx)

PetscErrorCode __libmesh_petsc_snes_jacobian (SNES, Vec x, Mat *jac, Mat *pc, MatStructure *msflag, void *ctx)
 

Typedef Documentation

 

typedef int PetscErrorCode

Definition at line 51 of file petsc_nonlinear_solver.C.  

typedef int PetscInt

Definition at line 52 of file petsc_nonlinear_solver.C.  

Function Documentation

 

PetscErrorCode __libmesh_petsc_snes_jacobian (SNES, Vecx, Mat *jac, Mat *pc, MatStructure *msflag, void *ctx)

Definition at line 124 of file petsc_nonlinear_solver.C.

References PetscMatrix< T >::close(), NonlinearSolver< T >::jacobian, NonlinearSolver< T >::matvec, PetscMatrix< T >::swap(), NonlinearSolver< T >::system(), and PetscMatrix< T >::zero().

Referenced by PetscNonlinearSolver< T >::solve().

  {
    int ierr=0;
    
    libmesh_assert (ctx != NULL);
    
    PetscNonlinearSolver<Number>* solver =
      static_cast<PetscNonlinearSolver<Number>*> (ctx);

    NonlinearImplicitSystem &sys = solver->system();
    
    PetscMatrix<Number> PC(*pc);
    PetscMatrix<Number> Jac(*jac);
    PetscVector<Number>& X_sys = *libmesh_cast_ptr<PetscVector<Number>*>(sys.solution.get());
    PetscMatrix<Number>& Jac_sys = *libmesh_cast_ptr<PetscMatrix<Number>*>(sys.matrix);
    PetscVector<Number> X_global(x);

    // Use the systems update() to get a good local version of the parallel solution
    X_global.swap(X_sys);
    Jac.swap(Jac_sys);

    sys.get_dof_map().enforce_constraints_exactly(sys);
    sys.update();

    X_global.swap(X_sys);
    Jac.swap(Jac_sys);

    PC.zero();

    if      (solver->jacobian != NULL) solver->jacobian (*sys.current_local_solution.get(), PC);
    else if (solver->matvec   != NULL) solver->matvec   (*sys.current_local_solution.get(), NULL, &PC);
    else libmesh_error();
    
    PC.close();
    Jac.close();
    X_global.close();
    
    *msflag = SAME_NONZERO_PATTERN;
    
    return ierr;
  }
 

PetscErrorCode __libmesh_petsc_snes_monitor (SNES, PetscIntits, PetscRealfnorm, void *)

Definition at line 58 of file petsc_nonlinear_solver.C.

Referenced by PetscNonlinearSolver< T >::init().

  {
    //int ierr=0;

    //if (its > 0)
      std::cout << '  NL step '
                << std::setw(2) << its
                << std::scientific
                << ', |residual|_2 = ' << fnorm
                << std::endl;

    //return ierr;
    return 0;
  }
 

PetscErrorCode __libmesh_petsc_snes_residual (SNES, Vecx, Vecr, void *ctx)

Definition at line 78 of file petsc_nonlinear_solver.C.

References NonlinearSolver< T >::matvec, NonlinearSolver< T >::residual, and NonlinearSolver< T >::system().

Referenced by PetscNonlinearSolver< T >::solve().

  {
    int ierr=0;

    libmesh_assert (x   != NULL);
    libmesh_assert (r   != NULL);
    libmesh_assert (ctx != NULL);
    
    PetscNonlinearSolver<Number>* solver =
      static_cast<PetscNonlinearSolver<Number>*> (ctx);
    
    NonlinearImplicitSystem &sys = solver->system();

    PetscVector<Number>& X_sys = *libmesh_cast_ptr<PetscVector<Number>*>(sys.solution.get());
    PetscVector<Number>& R_sys = *libmesh_cast_ptr<PetscVector<Number>*>(sys.rhs);
    PetscVector<Number> X_global(x), R(r);

    // Use the systems update() to get a good local version of the parallel solution
    X_global.swap(X_sys);
    R.swap(R_sys);

    sys.get_dof_map().enforce_constraints_exactly(sys);
    
    sys.update();

    //Swap back
    X_global.swap(X_sys);
    R.swap(R_sys);

    R.zero();

    if      (solver->residual != NULL) solver->residual (*sys.current_local_solution.get(), R);
    else if (solver->matvec   != NULL) solver->matvec   (*sys.current_local_solution.get(), &R, NULL);
    else libmesh_error();
    
    R.close();
    X_global.close();
    
    return ierr;
  }
 

Author

Generated automatically by Doxygen for libMesh from the source code.


 

Index

NAME
SYNOPSIS
Typedefs
Functions
Typedef Documentation
typedef int PetscErrorCode
typedef int PetscInt
Function Documentation
PetscErrorCode __libmesh_petsc_snes_jacobian (SNES, Vecx, Mat *jac, Mat *pc, MatStructure *msflag, void *ctx)
PetscErrorCode __libmesh_petsc_snes_monitor (SNES, PetscIntits, PetscRealfnorm, void *)
PetscErrorCode __libmesh_petsc_snes_residual (SNES, Vecx, Vecr, void *ctx)
Author

This document was created by man2html, using the manual pages.
Time: 21:52:24 GMT, April 16, 2011