[petsc-dev] Valgrind defect: memory leak with PetscCommDuplicate?
Brendan Kochunas
bkochuna at umich.edu
Fri Feb 6 13:17:07 CST 2015
Hi Barry/Boyce,
I did a bit more digging and was able to reproduce this defect using a
slightly modified version of ex2f.F located at:
http://www.mcs.anl.gov/petsc/petsc-3.3/src/ksp/ksp/examples/tutorials/ex2f.F
I have pasted this modified version below. Essentially the modification
replaces PETSC_COMM_WORLD with a duplicate comm of MPI_COMM_SELF in each
of the Mat/Vec/KSPCreate calls. (tested with 1 mpi proc)
I'll note that the defect was NOT reported if a duplicate of
PETSC_COMM_WORLD was used or if MPI_COMM_WORLD was used or a duplicate
of it was used.
It should be worthwhile to test this modified example against mpich to
confirm that the OpenMPI is the culprit in this case.
I have had similar experiences between MPICH and OpenMPI, so I would not
be surprised in the least if OpenMPI is the source of the leak.
Thanks,
-Brendan
!
! Description: Solves a linear system in parallel with KSP (Fortran code).
! Also shows how to set a user-defined monitoring routine.
!
! Program usage: mpiexec -n <procs> ex2f [-help] [all PETSc options]
!
!/*T
! Concepts: KSP^basic parallel example
! Concepts: KSP^setting a user-defined monitoring routine
! Processors: n
!T*/
!
! -----------------------------------------------------------------------
program main
implicit none
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Include files
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!
! This program uses CPP for preprocessing, as indicated by the use of
! PETSc include files in the directory petsc/include/finclude. This
! convention enables use of the CPP preprocessor, which allows the use
! of the #include statements that define PETSc objects and variables.
!
! Use of the conventional Fortran include statements is also supported
! In this case, the PETsc include files are located in the directory
! petsc/include/foldinclude.
!
! Since one must be very careful to include each file no more than once
! in a Fortran routine, application programmers must exlicitly list
! each file needed for the various PETSc components within their
! program (unlike the C/C++ interface).
!
! See the Fortran section of the PETSc users manual for details.
!
! The following include statements are required for KSP Fortran programs:
! petscsys.h - base PETSc routines
! petscvec.h - vectors
! petscmat.h - matrices
! petscpc.h - preconditioners
! petscksp.h - Krylov subspace methods
! Additional include statements may be needed if using additional
! PETSc routines in a Fortran program, e.g.,
! petscviewer.h - viewers
! petscis.h - index sets
!
#include <finclude/petscsys.h>
#include <finclude/petscvec.h>
#include <finclude/petscmat.h>
#include <finclude/petscpc.h>
#include <finclude/petscksp.h>
!
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Variable declarations
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!
! Variables:
! ksp - linear solver context
! ksp - Krylov subspace method context
! pc - preconditioner context
! x, b, u - approx solution, right-hand-side, exact solution vectors
! A - matrix that defines linear system
! its - iterations for convergence
! norm - norm of error in solution
! rctx - random number generator context
!
! Note that vectors are declared as PETSc "Vec" objects. These vectors
! are mathematical objects that contain more than just an array of
! double precision numbers. I.e., vectors in PETSc are not just
! double precision x(*).
! However, local vector data can be easily accessed via VecGetArray().
! See the Fortran section of the PETSc users manual for details.
!
integer testComm
double precision norm
PetscInt i,j,II,JJ,m,n,its
PetscInt Istart,Iend,ione
PetscErrorCode ierr
PetscMPIInt rank,size
PetscBool flg
PetscScalar v,one,neg_one
Vec x,b,u
Mat A
KSP ksp
PetscRandom rctx
! These variables are not currently used.
! PC pc
! PCType ptype
! double precision tol
! Note: Any user-defined Fortran routines (such as MyKSPMonitor)
! MUST be declared as external.
external MyKSPMonitor,MyKSPConverged
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Beginning of program
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
m = 3
n = 3
one = 1.0
neg_one = -1.0
ione = 1
call PetscOptionsGetInt(PETSC_NULL_CHARACTER,'-m',m,flg,ierr)
call PetscOptionsGetInt(PETSC_NULL_CHARACTER,'-n',n,flg,ierr)
call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
call MPI_Comm_size(PETSC_COMM_WORLD,size,ierr)
! call MPI_Comm_dup(MPI_COMM_SELF,testComm,ierr)
testComm=MPI_COMM_SELF
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Compute the matrix and right-hand-side vector that define
! the linear system, Ax = b.
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Create parallel matrix, specifying only its global dimensions.
! When using MatCreate(), the matrix format can be specified at
! runtime. Also, the parallel partitioning of the matrix is
! determined by PETSc at runtime.
call MatCreate(testComm,A,ierr)
call MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m*n,m*n,ierr)
call MatSetFromOptions(A,ierr)
call MatSetUp(A,ierr)
! Currently, all PETSc parallel matrix formats are partitioned by
! contiguous chunks of rows across the processors. Determine which
! rows of the matrix are locally owned.
call MatGetOwnershipRange(A,Istart,Iend,ierr)
! Set matrix elements for the 2-D, five-point stencil in parallel.
! - Each processor needs to insert only elements that it owns
! locally (but any non-local elements will be sent to the
! appropriate processor during matrix assembly).
! - Always specify global row and columns of matrix entries.
! - Note that MatSetValues() uses 0-based row and column numbers
! in Fortran as well as in C.
! Note: this uses the less common natural ordering that orders first
! all the unknowns for x = h then for x = 2h etc; Hence you see JH =
II +- n
! instead of JJ = II +- m as you might expect. The more standard
ordering
! would first do all variables for y = h, then y = 2h etc.
do 10, II=Istart,Iend-1
v = -1.0
i = II/n
j = II - i*n
if (i.gt.0) then
JJ = II - n
call MatSetValues(A,ione,II,ione,JJ,v,INSERT_VALUES,ierr)
endif
if (i.lt.m-1) then
JJ = II + n
call MatSetValues(A,ione,II,ione,JJ,v,INSERT_VALUES,ierr)
endif
if (j.gt.0) then
JJ = II - 1
call MatSetValues(A,ione,II,ione,JJ,v,INSERT_VALUES,ierr)
endif
if (j.lt.n-1) then
JJ = II + 1
call MatSetValues(A,ione,II,ione,JJ,v,INSERT_VALUES,ierr)
endif
v = 4.0
call MatSetValues(A,ione,II,ione,II,v,INSERT_VALUES,ierr)
10 continue
! Assemble matrix, using the 2-step process:
! MatAssemblyBegin(), MatAssemblyEnd()
! Computations can be done while messages are in transition,
! by placing code between these two statements.
call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
! Create parallel vectors.
! - Here, the parallel partitioning of the vector is determined by
! PETSc at runtime. We could also specify the local dimensions
! if desired -- or use the more general routine VecCreate().
! - When solving a linear system, the vectors and matrices MUST
! be partitioned accordingly. PETSc automatically generates
! appropriately partitioned matrices and vectors when MatCreate()
! and VecCreate() are used with the same communicator.
! - Note: We form 1 vector from scratch and then duplicate as needed.
call VecCreateMPI(testComm,PETSC_DECIDE,m*n,u,ierr)
call VecSetFromOptions(u,ierr)
call VecDuplicate(u,b,ierr)
call VecDuplicate(b,x,ierr)
! Set exact solution; then compute right-hand-side vector.
! By default we use an exact solution of a vector with all
! elements of 1.0; Alternatively, using the runtime option
! -random_sol forms a solution vector with random components.
call PetscOptionsHasName(PETSC_NULL_CHARACTER, &
& "-random_exact_sol",flg,ierr)
if (flg) then
call PetscRandomCreate(testComm,rctx,ierr)
call PetscRandomSetFromOptions(rctx,ierr)
call VecSetRandom(u,rctx,ierr)
call PetscRandomDestroy(rctx,ierr)
else
call VecSet(u,one,ierr)
endif
call MatMult(A,u,b,ierr)
! View the exact solution vector if desired
call PetscOptionsHasName(PETSC_NULL_CHARACTER, &
& "-view_exact_sol",flg,ierr)
if (flg) then
call VecView(u,PETSC_VIEWER_STDOUT_WORLD,ierr)
endif
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Create the linear solver and set various options
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Create linear solver context
call KSPCreate(testComm,ksp,ierr)
! Set operators. Here the matrix that defines the linear system
! also serves as the preconditioning matrix.
call KSPSetOperators(ksp,A,A,DIFFERENT_NONZERO_PATTERN,ierr)
! Set linear solver defaults for this problem (optional).
! - By extracting the KSP and PC contexts from the KSP context,
! we can then directly directly call any KSP and PC routines
! to set various options.
! - The following four statements are optional; all of these
! parameters could alternatively be specified at runtime via
! KSPSetFromOptions(). All of these defaults can be
! overridden at runtime, as indicated below.
! We comment out this section of code since the Jacobi
! preconditioner is not a good general default.
! call KSPGetPC(ksp,pc,ierr)
! ptype = PCJACOBI
! call PCSetType(pc,ptype,ierr)
! tol = 1.e-7
! call KSPSetTolerances(ksp,tol,PETSC_DEFAULT_DOUBLE_PRECISION,
! & PETSC_DEFAULT_DOUBLE_PRECISION,PETSC_DEFAULT_INTEGER,ierr)
! Set user-defined monitoring routine if desired
call PetscOptionsHasName(PETSC_NULL_CHARACTER,'-my_ksp_monitor', &
& flg,ierr)
if (flg) then
call KSPMonitorSet(ksp,MyKSPMonitor,PETSC_NULL_OBJECT, &
& PETSC_NULL_FUNCTION,ierr)
endif
! Set runtime options, e.g.,
! -ksp_type <type> -pc_type <type> -ksp_monitor -ksp_rtol <rtol>
! These options will override those specified above as long as
! KSPSetFromOptions() is called _after_ any other customization
! routines.
call KSPSetFromOptions(ksp,ierr)
! Set convergence test routine if desired
call PetscOptionsHasName(PETSC_NULL_CHARACTER, &
& '-my_ksp_convergence',flg,ierr)
if (flg) then
call KSPSetConvergenceTest(ksp,MyKSPConverged, &
& PETSC_NULL_OBJECT,PETSC_NULL_FUNCTION,ierr)
endif
!
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Solve the linear system
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
call KSPSolve(ksp,b,x,ierr)
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Check solution and clean up
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Check the error
call VecAXPY(x,neg_one,u,ierr)
call VecNorm(x,NORM_2,norm,ierr)
call KSPGetIterationNumber(ksp,its,ierr)
if (rank .eq. 0) then
if (norm .gt. 1.e-12) then
write(6,100) norm,its
else
write(6,110) its
endif
endif
100 format('Norm of error ',e11.4,' iterations ',i5)
110 format('Norm of error < 1.e-12,iterations ',i5)
! Free work space. All PETSc objects should be destroyed when they
! are no longer needed.
call KSPDestroy(ksp,ierr)
call VecDestroy(u,ierr)
call VecDestroy(x,ierr)
call VecDestroy(b,ierr)
call MatDestroy(A,ierr)
! Always call PetscFinalize() before exiting a program. This routine
! - finalizes the PETSc libraries as well as MPI
! - provides summary and diagnostic information if certain runtime
! options are chosen (e.g., -log_summary). See PetscFinalize()
! manpage for more information.
call PetscFinalize(ierr)
end
! --------------------------------------------------------------
!
! MyKSPMonitor - This is a user-defined routine for monitoring
! the KSP iterative solvers.
!
! Input Parameters:
! ksp - iterative context
! n - iteration number
! rnorm - 2-norm (preconditioned) residual value (may be estimated)
! dummy - optional user-defined monitor context (unused here)
!
subroutine MyKSPMonitor(ksp,n,rnorm,dummy,ierr)
implicit none
#include <finclude/petscsys.h>
#include <finclude/petscvec.h>
#include <finclude/petscksp.h>
KSP ksp
Vec x
PetscErrorCode ierr
PetscInt n,dummy
PetscMPIInt rank
double precision rnorm
! Build the solution vector
call KSPBuildSolution(ksp,PETSC_NULL_OBJECT,x,ierr)
! Write the solution vector and residual norm to stdout
! - Note that the parallel viewer PETSC_VIEWER_STDOUT_WORLD
! handles data from multiple processors so that the
! output is not jumbled.
call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
if (rank .eq. 0) write(6,100) n
call VecView(x,PETSC_VIEWER_STDOUT_WORLD,ierr)
if (rank .eq. 0) write(6,200) n,rnorm
100 format('iteration ',i5,' solution vector:')
200 format('iteration ',i5,' residual norm ',e11.4)
ierr = 0
end
! --------------------------------------------------------------
!
! MyKSPConverged - This is a user-defined routine for testing
! convergence of the KSP iterative solvers.
!
! Input Parameters:
! ksp - iterative context
! n - iteration number
! rnorm - 2-norm (preconditioned) residual value (may be estimated)
! dummy - optional user-defined monitor context (unused here)
!
subroutine MyKSPConverged(ksp,n,rnorm,flag,dummy,ierr)
implicit none
#include <finclude/petscsys.h>
#include <finclude/petscvec.h>
#include <finclude/petscksp.h>
KSP ksp
PetscErrorCode ierr
PetscInt n,dummy
KSPConvergedReason flag
double precision rnorm
if (rnorm .le. .05) then
flag = 1
else
flag = 0
endif
ierr = 0
end
On 2/6/2015 2:13 PM, Boyce Griffith wrote:
>> On Feb 6, 2015, at 2:04 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>>
>>
>> If you don't see any memory leaks directly from PETSc routines then it is a problem in OpenMPI.
> OpenMPI intentionally does stuff that is not Valgrind clean:
>
> http://www.open-mpi.org/faq/?category=debugging#valgrind_clean
>
> I've personally had mixed success with the OpenMPI-provided suppression file in the past, but I haven't tried it in a while.
>
> -- Boyce
>
>> We've found that MPICH seems to have less memory leaks than OpenMPI.
>>
>> Barry
>>
>>> On Feb 6, 2015, at 8:52 AM, Brendan Kochunas <bkochuna at umich.edu> wrote:
>>>
>>> Hi we are trying to clear valgrind defects from our code and presently we are valgrind is reporting memory leaks like the following:
>>>
>>> ==2884== at 0x4A07EB7: malloc (vg_replace_malloc.c:296)
>>> ==2884== by 0x8519874: set_value.isra.0.part.1 (in /gcc-4.6.1/toolset/openmpi-1.4.3/lib/libmpi.so.0.0.2)
>>> ==2884== by 0x8547E4D: PMPI_Attr_put (in /gcc-4.6.1/toolset/openmpi-1.4.3/lib/libmpi.so.0.0.2)
>>> ==2884== by 0x113A153: PetscCommDuplicate
>>> ==2884== by 0x113BFA3: PetscHeaderCreate_Private
>>> ==2884== by 0x129ADC6: MatCreate
>>> ==2884== by 0x123C355: MatMPIAIJSetPreallocation_MPIAIJ
>>> ==2884== by 0x124DAB4: MatMPIAIJSetPreallocation
>>> ==2884== by 0x12549B0: MatSetUp_MPIAIJ
>>> ==2884== by 0x1191186: MatSetUp
>>> ==2884== by 0x10D9333: matsetup_
>>>
>>> and...
>>>
>>> ==2884== at 0x4A07EB7: malloc (vg_replace_malloc.c:296)
>>> ==2884== by 0x8519874: set_value.isra.0.part.1 (in /gcc-4.6.1/toolset/openmpi-1.4.3/lib/libmpi.so.0.0.2)
>>> ==2884== by 0x8547E4D: PMPI_Attr_put (in /gcc-4.6.1/toolset/openmpi-1.4.3/lib/libmpi.so.0.0.2)
>>> ==2884== by 0x113A22D: PetscCommDuplicate
>>> ==2884== by 0x113BFA3: PetscHeaderCreate_Private
>>> ==2884== by 0x10ED42B: KSPCreate
>>> ==2884== by 0x10DA55C: kspcreate_
>>>
>>> Is the development team aware of any memory leaks that may be coming from PetscCommDuplicate as it may be used in the call stack shown above?
>>>
>>> The version of PETSc we are linking with is 3.3 patch 4. And this is built on OpenMPI v. 1.4.3.
>>>
>>> We are trying to determine if the leak is due to
>>> 1. Our codes usage of PETSc
>>> 2. The actual PETSc library
>>> 3. PETSc's usage of MPI
>>> 3. The OpenMPI library that PETSc was built against
>>>
>>> Any assistance in being able to point to the culprit or suggestions for particular tests (e.g. a PETSc example) worth trying to identify the root issue would be appreciated.
>>>
>>> Thanks in advance!
>>> -Brendan
>>>
More information about the petsc-dev
mailing list