[petsc-dev] Cleanup of dmhooks in snesdestroy (was: Re: Discussion about time-dependent optimization moved from PR)
Jed Brown
jed at jedbrown.org
Thu Oct 19 08:19:34 CDT 2017
Fixed here, along with some related issues.
https://bitbucket.org/petsc/petsc/pull-requests/780/jed-fix-dmcoarsenhookadd-identical/diff
The TS setting of hooks should probably be moved out of the individual
implementations so they can reliably be cleaned up if the user is
switching DMs around, but I don't have time to test a refactoring like
that and if I don't fix the pressing issue now, it'll become another
project that I start but don't finish...
Lawrence Mitchell <wencel at gmail.com> writes:
>> On 17 Oct 2017, at 17:40, Lawrence Mitchell <wencel at gmail.com> wrote:
>>
>>
>>> On 17 Oct 2017, at 17:35, Jed Brown <jed at jedbrown.org> wrote:
>>>
>>>>
>>>> So my initial concern was that if I do:
>>>>
>>>> SNESSetDM(snes, dm);
>>>> SNESSolve(snes);
>>>> SNESDestroy(&snes);
>>>>
>>>> Then dm still has stale pointers to snes (in the context of the
>>>> restrict hook).
>>>
>>> Thanks, this is fixable. Can someone make a test case for me?
>>
>> I'll try and cook one up tomorrow.
>
> Here we go.
>
> This just reuses a DA in a loop making two SNESes to solve the same problem.
>
> $ ./dmhookcleanup -da_refine 1 -pc_type mg -snes_fd_color -snes_monitor
> 0 SNES Function norm 4.375000000000e-01
> 1 SNES Function norm 1.086628858540e-05
> 2 SNES Function norm 7.356279464858e-12
> 0 SNES Function norm 4.375000000000e-01
> [0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
> [0]PETSC ERROR: Invalid argument
> [0]PETSC ERROR: Wrong type of object: Parameter # 2
> [0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.
> [0]PETSC ERROR: Petsc Development GIT revision: v3.7.6-3851-g17f6384 GIT Date: 2017-05-09 22:45:05 -0500
> [0]PETSC ERROR: ./dmhookcleanup on a arch-linux2-c-dbg named yam.doc.ic.ac.uk by lmitche1 Wed Oct 18 11:04:11 2017
> [0]PETSC ERROR: Configure options --download-chaco=1 --download-ctetgen=1 --download-exodusii=1 --download-hdf5=1 --download-hypre=1 --download-metis=1 --download-ml=1 --download-mumps=1 --download-netcdf=1 --download-parmetis=1 --download-ptscotch=1 --download-scalapack=1 --download-triangle=1 --download-eigen=1 --with-c2html=0 --with-debugging=1 --with-make-np=32 --with-shared-libraries=1 PETSC_ARCH=arch-linux2-c-dbg
> [0]PETSC ERROR: #1 MatRestrict() line 8044 in /data/lmitche1/src/deps/petsc/src/mat/interface/matrix.c
> [0]PETSC ERROR: #2 DMRestrictHook_SNESVecSol() line 528 in /data/lmitche1/src/deps/petsc/src/snes/interface/snes.c
> [0]PETSC ERROR: #3 DMRestrict() line 2421 in /data/lmitche1/src/deps/petsc/src/dm/interface/dm.c
> [0]PETSC ERROR: #4 PCSetUp_MG() line 729 in /data/lmitche1/src/deps/petsc/src/ksp/pc/impls/mg/mg.c
> [0]PETSC ERROR: #5 PCSetUp() line 924 in /data/lmitche1/src/deps/petsc/src/ksp/pc/interface/precon.c
> [0]PETSC ERROR: #6 KSPSetUp() line 378 in /data/lmitche1/src/deps/petsc/src/ksp/ksp/interface/itfunc.c
> [0]PETSC ERROR: #7 KSPSolve() line 609 in /data/lmitche1/src/deps/petsc/src/ksp/ksp/interface/itfunc.c
> [0]PETSC ERROR: #8 SNESSolve_NEWTONLS() line 224 in /data/lmitche1/src/deps/petsc/src/snes/impls/ls/ls.c
> [0]PETSC ERROR: #9 SNESSolve() line 4106 in /data/lmitche1/src/deps/petsc/src/snes/interface/snes.c
> [0]PETSC ERROR: #10 main() line 61 in /data/lmitche1/src/petsc-doodles/dmhookcleanup.c
> [0]PETSC ERROR: PETSc Option Table entries:
> [0]PETSC ERROR: -da_refine 1
> [0]PETSC ERROR: -pc_type mg
> [0]PETSC ERROR: -snes_fd_color
> [0]PETSC ERROR: -snes_monitor
> [0]PETSC ERROR: ----------------End of Error Message -------send entire error message to petsc-maint at mcs.anl.gov----------
> --------------------------------------------------------------------------
> MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD
> with errorcode 62.
>
> NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.
> You may or may not see output from other processes, depending on
> exactly when Open MPI kills them.
> --------------------------------------------------------------------------
>
> valgrind log:
>
> ==2215210== Memcheck, a memory error detector
> ==2215210== Copyright (C) 2002-2015, and GNU GPL'd, by Julian Seward et al.
> ==2215210== Using Valgrind-3.11.0 and LibVEX; rerun with -h for copyright info
> ==2215210== Command: ./dmhookcleanup -da_refine 1 -pc_type mg -snes_fd_color -snes_monitor
> ==2215210== Parent PID: 1523638
> ==2215210==
> ==2215210== Warning: noted but unhandled ioctl 0x30000001 with no size/direction hints.
> ==2215210== This could cause spurious value errors to appear.
> ==2215210== See README_MISSING_SYSCALL_OR_IOCTL for guidance on writing a proper wrapper.
> ==2215210== Warning: noted but unhandled ioctl 0x27 with no size/direction hints.
> ==2215210== This could cause spurious value errors to appear.
> ==2215210== See README_MISSING_SYSCALL_OR_IOCTL for guidance on writing a proper wrapper.
> ==2215210== Warning: noted but unhandled ioctl 0x7ff with no size/direction hints.
> ==2215210== This could cause spurious value errors to appear.
> ==2215210== See README_MISSING_SYSCALL_OR_IOCTL for guidance on writing a proper wrapper.
> ==2215210== Warning: noted but unhandled ioctl 0x25 with no size/direction hints.
> ==2215210== This could cause spurious value errors to appear.
> ==2215210== See README_MISSING_SYSCALL_OR_IOCTL for guidance on writing a proper wrapper.
> ==2215210== Warning: noted but unhandled ioctl 0x17 with no size/direction hints.
> ==2215210== This could cause spurious value errors to appear.
> ==2215210== See README_MISSING_SYSCALL_OR_IOCTL for guidance on writing a proper wrapper.
> ==2215210== Warning: set address range perms: large range [0x1000000000, 0x2100000000) (noaccess)
> ==2215210== Warning: set address range perms: large range [0x200000000, 0x400000000) (noaccess)
> ==2215210== Invalid read of size 8
> ==2215210== at 0x63E30DD: DMRestrictHook_SNESVecSol (snes.c:519)
> ==2215210== by 0x5A3C98B: DMRestrict (dm.c:2421)
> ==2215210== by 0x61A3C7A: PCSetUp_MG (mg.c:729)
> ==2215210== by 0x61E581C: PCSetUp (precon.c:924)
> ==2215210== by 0x6334EB8: KSPSetUp (itfunc.c:378)
> ==2215210== by 0x6336C76: KSPSolve (itfunc.c:609)
> ==2215210== by 0x6484F80: SNESSolve_NEWTONLS (ls.c:224)
> ==2215210== by 0x640A12D: SNESSolve (snes.c:4106)
> ==2215210== by 0x4014EA: main (dmhookcleanup.c:61)
> ==2215210== Address 0x150a5aa0 is 688 bytes inside a block of size 1,280 free'd
> ==2215210== at 0x4C2EDEB: free (in /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
> ==2215210== by 0x51565F6: PetscFreeAlign (mal.c:85)
> ==2215210== by 0x63FE427: SNESDestroy (snes.c:2896)
> ==2215210== by 0x401550: main (dmhookcleanup.c:62)
> ==2215210== Block was alloc'd at
> ==2215210== at 0x4C2FFC6: memalign (in /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
> ==2215210== by 0x515653C: PetscMallocAlign (mal.c:39)
> ==2215210== by 0x63EE313: SNESCreate (snes.c:1553)
> ==2215210== by 0x4013A7: main (dmhookcleanup.c:58)
> ==2215210==
> ==2215210==
> ==2215210== HEAP SUMMARY:
> ==2215210== in use at exit: 908,831 bytes in 9,597 blocks
> ==2215210== total heap usage: 52,827 allocs, 43,230 frees, 959,732,566 bytes allocated
> ==2215210==
> ==2215210== LEAK SUMMARY:
> ==2215210== definitely lost: 674 bytes in 9 blocks
> ==2215210== indirectly lost: 649,619 bytes in 9,354 blocks
> ==2215210== possibly lost: 2,800 bytes in 21 blocks
> ==2215210== still reachable: 255,738 bytes in 213 blocks
> ==2215210== suppressed: 0 bytes in 0 blocks
> ==2215210== Rerun with --leak-check=full to see details of leaked memory
> ==2215210==
> ==2215210== For counts of detected and suppressed errors, rerun with: -v
> ==2215210== ERROR SUMMARY: 2 errors from 1 contexts (suppressed: 0 from 0)
>
> Code:
>
> #include <petsc.h>
>
> static PetscErrorCode FormFunctionLocal(DMDALocalInfo *info, PetscScalar **xarr, PetscScalar **farr, void *ctx)
> {
> PetscInt i, j;
> PetscReal hx, hy, hxdhy, hydhx;
> PetscScalar u, ue, uw, un, us;
>
> PetscFunctionBeginUser;
>
> hx = 2.0/(PetscReal)(info->mx-1);
> hy = 2.0/(PetscReal)(info->my-1);
> hxdhy = hx/hy;
> hydhx = hy/hx;
>
> for (j=info->ys; j<info->ys+info->ym; j++) {
> for (i=info->xs; i<info->xs+info->xm; i++) {
> if (i == 0 || j == 0 || i == info->mx-1 || j == info->my-1) {
> farr[j][i] = 2.0*(hydhx + hxdhy)*xarr[j][i];
> } else {
> u = xarr[j][i];
> uw = xarr[j][i-1];
> ue = xarr[j][i+1];
> un = xarr[j-1][i];
> us = xarr[j+1][i];
>
> if (i-1 == 0) uw = 0.;
> if (i+1 == info->mx-1) ue = 0.;
> if (j-1 == 0) un = 0.;
> if (j+1 == info->my-1) us = 0.;
>
> farr[j][i] = (2.0*u - uw - ue)*hydhx + (2.0*u - un - us)*hxdhy - hx*hy;
> }
> }
> }
> PetscFunctionReturn(0);
> }
>
> int main(int argc,char **argv)
> {
> SNES snes;
> PetscErrorCode ierr;
> DM da;
> PetscInt i;
>
> ierr = PetscInitialize(&argc,&argv,NULL,NULL);
> if (ierr) return ierr;
>
> ierr = DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE,
> DMDA_STENCIL_STAR, 5, 5, PETSC_DECIDE, PETSC_DECIDE,
> 1, 1, 0, 0, &da);CHKERRQ(ierr);
> ierr = DMSetFromOptions(da);CHKERRQ(ierr);
> ierr = DMSetUp(da);CHKERRQ(ierr);
>
> ierr = DMDASNESSetFunctionLocal(da, INSERT_VALUES, (PetscErrorCode (*)(DMDALocalInfo*,void*,void*,void*))FormFunctionLocal,NULL);CHKERRQ(ierr);
>
> for (i = 0; i < 2; i++) {
> ierr = SNESCreate(PETSC_COMM_WORLD,&snes);CHKERRQ(ierr);
> ierr = SNESSetDM(snes,da);CHKERRQ(ierr);
> ierr = SNESSetFromOptions(snes);CHKERRQ(ierr);
> ierr = SNESSolve(snes,NULL,NULL);CHKERRQ(ierr);
> ierr = SNESDestroy(&snes);CHKERRQ(ierr);
> }
> ierr = DMDestroy(&da);CHKERRQ(ierr);
> ierr = PetscFinalize();
> return ierr;
> }
>
>
> Cheers,
>
> Lawrence
More information about the petsc-dev
mailing list