Re: [petsc-users] Valgrind Errors with KSPSolve
Thanks Matt, I got it figured out. As with a lot of uninitialised value bugs it started a long way away (in a scalar value) and propagated through. What is slightly bizarre is that valgrind didn't pick it up where I define and first use it, even with --track-origins=yes. Thank you for your help. Phil On 26/10/2018 15:01, Matthew Knepley wrote: > On Fri, Oct 26, 2018 at 9:54 AM Phil Tooley > mailto:phil.too...@sheffield.ac.uk>> wrote: > > Hi All, > > Running valgrind over my code reveals a huge number (it stops > collecting > after 1e6) of uninitialised value errors occuring the KSPSolve > routine. > These don't occur with the ksp example scripts but I can't figure out > what I am doing wrong in my code. My code works as expected. The > same > errors occur with both version 3.9.4 and 3.10.0. > > I am developing in c++ and using smart pointers to allow > autodestruction > of Petsc objects, but I don't think this should have any side > effects. > > > From casual inspection, it sure seems like the matrix you are passing > KSP has unitialized > values. > > I used smart pointers for years, and work on > packages that had smart > pointers. In almost every case, we ended up ripping out the smart > pointers and going with > explicit destruction, the amount of work they save never outweighed > the problems with > detecting errors > > Thanks, > > Matt > > Thanks, > > Matt > > > I do the ksp setup as: > > KSP_unique m_ksp = create_unique_ksp(); > > perr = KSPCreate(m_comm, m_ksp.get());CHKERRABORT(m_comm, perr); > > perr = KSPSetOperators(*m_ksp, *normmat, *normmat);CHKERRABORT(m_comm, > perr); > perr = KSPSetUp(*m_ksp);CHKERRABORT(m_comm, perr); > > perr = KSPSetFromOptions(*m_ksp);CHKERRABORT(m_comm, perr); > perr = KSPSetUp(*m_ksp);CHKERRABORT(m_comm, perr); > > perr = KSPSolve(*m_ksp, *m_rhs, *m_delta);CHKERRABORT(m_comm, perr); > > Where normmat is a smart pointer to a previously created Mat > object and > m_rhs and m_delta are smart pointers to Vec objects which are > initialised using MatCreateVecs on normmat. > > I am using MPICH to provide a (hopefully) valgrind clean MPI > implementation. > > Any insights are appreciated. > > Many Thanks > > Phil Tooley > > Typical errors are as follows: > > ==7368== Conditional jump or move depends on uninitialised value(s) > ==7368== at 0x58894EE: PetscAbsScalar(double) (petscmath.h:389) > ==7368== by 0x5891702: MatLUFactorNumeric_SeqAIJ (aijfact.c:553) > ==7368== by 0x5406B83: MatLUFactorNumeric (matrix.c:3031) > ==7368== by 0x609F317: PCSetUp_ILU(_p_PC*) (ilu.c:176) > ==7368== by 0x6192D4F: PCSetUp (precon.c:923) > ==7368== by 0x62EAE56: KSPSetUp (itfunc.c:381) > ==7368== by 0x60190BB: PCSetUpOnBlocks_BJacobi_Singleblock(_p_PC*) > (bjacobi.c:621) > ==7368== by 0x619352C: PCSetUpOnBlocks (precon.c:954) > ==7368== by 0x62E8F0A: KSPSetUpOnBlocks (itfunc.c:213) > ==7368== by 0x62ECCB1: KSPSolve (itfunc.c:613) > ==7368== by 0x10FCDD: Elastic::innerstep(double) (in > /home/telemin/repos/pfire_petsc/bin/pfire) > ==7368== by 0x110600: Elastic::innerloop(int) (in > /home/telemin/repos/pfire_petsc/bin/pfire) > ==7368== > ==7368== Conditional jump or move depends on uninitialised value(s) > ==7368== at 0x58894EE: PetscAbsScalar(double) (petscmath.h:389) > ==7368== by 0x588A75F: MatPivotCheck_none(_p_Mat*, _p_Mat*, > MatFactorInfo const*, FactorShiftCtx*, int) (matimpl.h:704) > ==7368== by 0x588AD6E: MatPivotCheck(_p_Mat*, _p_Mat*, > MatFactorInfo > const*, FactorShiftCtx*, int) (matimpl.h:727) > ==7368== by 0x589177D: MatLUFactorNumeric_SeqAIJ (aijfact.c:558) > ==7368== by 0x5406B83: MatLUFactorNumeric (matrix.c:3031) > ==7368== by 0x609F317: PCSetUp_ILU(_p_PC*) (ilu.c:176) > ==7368== by 0x6192D4F: PCSetUp (precon.c:923) > ==7368== by 0x62EAE56: KSPSetUp (itfunc.c:381) > ==7368== by 0x60190BB: PCSetUpOnBlocks_BJacobi_Singleblock(_p_PC*) > (bjacobi.c:621) > ==7368== by 0x619352C: PCSetUpOnBlocks (precon.c:954) > ==7368== by 0x62E8F0A: KSPSetUpOnBlocks (itfunc.c:213) > ==7368== by 0x62ECCB1: KSPSolve (itfunc.c:613) > ==7368== > ==7368== Conditional jump or move depends on uninitialised value(s) > ==7368== at 0x588A76D: MatPivotCheck_none(_p_Mat*, _p_Mat*, > MatFactorInfo const*, FactorShiftCtx*, int) (matimpl.h:704) > ==7368== by 0x588AD6E: Mat
[petsc-users] Valgrind Errors with KSPSolve
Hi All, Running valgrind over my code reveals a huge number (it stops collecting after 1e6) of uninitialised value errors occuring the KSPSolve routine. These don't occur with the ksp example scripts but I can't figure out what I am doing wrong in my code. My code works as expected. The same errors occur with both version 3.9.4 and 3.10.0. I am developing in c++ and using smart pointers to allow autodestruction of Petsc objects, but I don't think this should have any side effects. I do the ksp setup as: KSP_unique m_ksp = create_unique_ksp(); perr = KSPCreate(m_comm, m_ksp.get());CHKERRABORT(m_comm, perr); perr = KSPSetOperators(*m_ksp, *normmat, *normmat);CHKERRABORT(m_comm, perr); perr = KSPSetUp(*m_ksp);CHKERRABORT(m_comm, perr); perr = KSPSetFromOptions(*m_ksp);CHKERRABORT(m_comm, perr); perr = KSPSetUp(*m_ksp);CHKERRABORT(m_comm, perr); perr = KSPSolve(*m_ksp, *m_rhs, *m_delta);CHKERRABORT(m_comm, perr); Where normmat is a smart pointer to a previously created Mat object and m_rhs and m_delta are smart pointers to Vec objects which are initialised using MatCreateVecs on normmat. I am using MPICH to provide a (hopefully) valgrind clean MPI implementation. Any insights are appreciated. Many Thanks Phil Tooley Typical errors are as follows: ==7368== Conditional jump or move depends on uninitialised value(s) ==7368== at 0x58894EE: PetscAbsScalar(double) (petscmath.h:389) ==7368== by 0x5891702: MatLUFactorNumeric_SeqAIJ (aijfact.c:553) ==7368== by 0x5406B83: MatLUFactorNumeric (matrix.c:3031) ==7368== by 0x609F317: PCSetUp_ILU(_p_PC*) (ilu.c:176) ==7368== by 0x6192D4F: PCSetUp (precon.c:923) ==7368== by 0x62EAE56: KSPSetUp (itfunc.c:381) ==7368== by 0x60190BB: PCSetUpOnBlocks_BJacobi_Singleblock(_p_PC*) (bjacobi.c:621) ==7368== by 0x619352C: PCSetUpOnBlocks (precon.c:954) ==7368== by 0x62E8F0A: KSPSetUpOnBlocks (itfunc.c:213) ==7368== by 0x62ECCB1: KSPSolve (itfunc.c:613) ==7368== by 0x10FCDD: Elastic::innerstep(double) (in /home/telemin/repos/pfire_petsc/bin/pfire) ==7368== by 0x110600: Elastic::innerloop(int) (in /home/telemin/repos/pfire_petsc/bin/pfire) ==7368== ==7368== Conditional jump or move depends on uninitialised value(s) ==7368== at 0x58894EE: PetscAbsScalar(double) (petscmath.h:389) ==7368== by 0x588A75F: MatPivotCheck_none(_p_Mat*, _p_Mat*, MatFactorInfo const*, FactorShiftCtx*, int) (matimpl.h:704) ==7368== by 0x588AD6E: MatPivotCheck(_p_Mat*, _p_Mat*, MatFactorInfo const*, FactorShiftCtx*, int) (matimpl.h:727) ==7368== by 0x589177D: MatLUFactorNumeric_SeqAIJ (aijfact.c:558) ==7368== by 0x5406B83: MatLUFactorNumeric (matrix.c:3031) ==7368== by 0x609F317: PCSetUp_ILU(_p_PC*) (ilu.c:176) ==7368== by 0x6192D4F: PCSetUp (precon.c:923) ==7368== by 0x62EAE56: KSPSetUp (itfunc.c:381) ==7368== by 0x60190BB: PCSetUpOnBlocks_BJacobi_Singleblock(_p_PC*) (bjacobi.c:621) ==7368== by 0x619352C: PCSetUpOnBlocks (precon.c:954) ==7368== by 0x62E8F0A: KSPSetUpOnBlocks (itfunc.c:213) ==7368== by 0x62ECCB1: KSPSolve (itfunc.c:613) ==7368== ==7368== Conditional jump or move depends on uninitialised value(s) ==7368== at 0x588A76D: MatPivotCheck_none(_p_Mat*, _p_Mat*, MatFactorInfo const*, FactorShiftCtx*, int) (matimpl.h:704) ==7368== by 0x588AD6E: MatPivotCheck(_p_Mat*, _p_Mat*, MatFactorInfo const*, FactorShiftCtx*, int) (matimpl.h:727) ==7368== by 0x589177D: MatLUFactorNumeric_SeqAIJ (aijfact.c:558) ==7368== by 0x5406B83: MatLUFactorNumeric (matrix.c:3031) ==7368== by 0x609F317: PCSetUp_ILU(_p_PC*) (ilu.c:176) ==7368== by 0x6192D4F: PCSetUp (precon.c:923) ==7368== by 0x62EAE56: KSPSetUp (itfunc.c:381) ==7368== by 0x60190BB: PCSetUpOnBlocks_BJacobi_Singleblock(_p_PC*) (bjacobi.c:621) ==7368== by 0x619352C: PCSetUpOnBlocks (precon.c:954) ==7368== by 0x62E8F0A: KSPSetUpOnBlocks (itfunc.c:213) ==7368== by 0x62ECCB1: KSPSolve (itfunc.c:613) ==7368== by 0x10FCDD: Elastic::innerstep(double) (in /home/telemin/repos/pfire_petsc/bin/pfire) ==7368== ==7368== Conditional jump or move depends on uninitialised value(s) ==7368== at 0x58894EE: PetscAbsScalar(double) (petscmath.h:389) ==7368== by 0x5889528: PetscIsNanScalar(double) (petscmath.h:674) ==7368== by 0x588A784: MatPivotCheck_none(_p_Mat*, _p_Mat*, MatFactorInfo const*, FactorShiftCtx*, int) (matimpl.h:704) ==7368== by 0x588AD6E: MatPivotCheck(_p_Mat*, _p_Mat*, MatFactorInfo const*, FactorShiftCtx*, int) (matimpl.h:727) ==7368== by 0x589177D: MatLUFactorNumeric_SeqAIJ (aijfact.c:558) ==7368== by 0x5406B83: MatLUFactorNumeric (matrix.c:3031) ==7368== by 0x609F317: PCSetUp_ILU(_p_PC*) (ilu.c:176) ==7368== by 0x6192D4F: PCSetUp (precon.c:923) ==7368== by 0x62EAE56: KSPSetUp (itfunc.c:381) ==7368== by 0x60190BB: PCSetUpOnBlocks_BJacobi_Singleblock(_p_PC*) (bjacobi.c:621) ==7368== by 0x619352C: PCSetUpOnBlocks (precon.c:954) ==7368== by 0x62E8F0A: KSPSetUpOnBlocks (itfunc.c:213) ==7368== ==736
[petsc-users] Better understanding VecScatter
Hi all, I am trying to understand how to create a custom scatter from petsc ordering in a local vector to natural ordering in a global vector. I have a 3D DMDA and local vector containing field data and am calculating the x y and z gradients into their own local vectors. I then need to scatter these gradients to different parts of a single vector so they are arranged [x_1 .. x_n y_1 .. y_n z_1 .. z_n] with natural ordering. (This is for use in a custom numerical scheme that I have been asked to parallelise) I can of course do this using multiple scatters but I need to perform this operation regularly and would like to express it as a single scatter operation if possible. I know I need to get the relevant AO from my DMDA and use it to construct appropriate ISs to build the scatter context but I am unsure how exactly to go about this. Can someone point me in the right direction please? Many Thanks Phil -- Phil Tooley Research Software Engineering University of Sheffield
Re: [petsc-users] DMDA with dimension of size 1
Hi Patrick, You are spot on, I was copying in a vector and not taking care to pad with 1s to length 3. Thanks Phil On 01/10/2018 15:57, Patrick Sanan wrote: > Whoops, sent that patch too fast (forgot to update SETERRQ3 to > SETERRQ4). Updated patch for maint attached. > > Am Mo., 1. Okt. 2018 um 16:55 Uhr schrieb Patrick Sanan > mailto:patrick.sa...@gmail.com>>: > > Meshes of size 1 should work. > > Looks like there is a bug in the code to produce this error > message (patch for maint attached). It's not outputting the "P" > (size in the 3rd dimension) component. > > This is just speculation without a full error message or code to > reproduce, but perhaps the size in the third dimension is an > uninitialized value which is triggering this warning with a 10 x > 10 x (huge garbage) x 1 (dof) mesh. > > Am Mo., 1. Okt. 2018 um 15:59 Uhr schrieb Phil Tooley > mailto:phil.too...@sheffield.ac.uk>>: > > Hi all, > > Is it valid to have a DMDA with one of the dimensions of size > 1. I was > hoping to avoid having to write explicit logic to handle the > case that I > am working on a 2D rather than a 3D image for my current > application. > Unfortunately when I try to construct such a DMDA I get an error: > > [0]PETSC ERROR: - Error Message > -- > [0]PETSC ERROR: Overflow in integer operation: > [0]PETSC ERROR: Mesh of 10 by 10 by 1 (dof) is too large for > 32 bit indices > > Is there a workaround other than "write everything twice"? > > Thanks > > Phil > > -- > Phil Tooley > Research Software Engineering > University of Sheffield > -- Phil Tooley Research Software Engineering University of Sheffield
[petsc-users] DMDA with dimension of size 1
Hi all, Is it valid to have a DMDA with one of the dimensions of size 1. I was hoping to avoid having to write explicit logic to handle the case that I am working on a 2D rather than a 3D image for my current application. Unfortunately when I try to construct such a DMDA I get an error: [0]PETSC ERROR: - Error Message -- [0]PETSC ERROR: Overflow in integer operation: [0]PETSC ERROR: Mesh of 10 by 10 by 1 (dof) is too large for 32 bit indices Is there a workaround other than "write everything twice"? Thanks Phil -- Phil Tooley Research Software Engineering University of Sheffield
Re: [petsc-users] Checking if a vector is a localvector of a given DMDA
Thanks both, I now have what I need. For now I am checking that the vector I am passed has the same local size, global size, and Comm as the vector provided by DMGetLocalVector, mostly because I already have a compatibility check function written. (I assume this requires a malloc and free behind the scenes) At some point I will likely change to explicitly checking for comm size of one and appropriate global and local sizes based on the DMDA properties instead, for now I want to get to an alpha version I can let people play with. Phil On 25/09/18 13:07, Dave May wrote: > > > On Tue, 25 Sep 2018 at 13:20, Matthew Knepley <mailto:knep...@gmail.com>> wrote: > > On Tue, Sep 25, 2018 at 7:03 AM Dave May <mailto:dave.mayhe...@gmail.com>> wrote: > > On Tue, 25 Sep 2018 at 11:49, Phil Tooley > <mailto:phil.too...@sheffield.ac.uk>> wrote: > > Hi all, > > Given a vector I know I can get an associated DM (if there > is one) by > calling VecGetDM, but I need to also be able to check that > > a) the vector is the localvector of that DM rather than > the global > > > Given the vector, you can check the communicator size via > PetscObjectGetComm() > > > > https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Sys/PetscObjectGetComm.html > and then MPI_Comm_size() > If the comm size 1, it is local vector. > > > In serial, both vectors have comm size 1. > > > Right - and the local and global sizes are the same. > > My point was to check the comm size first. If it's 1 then you have a > candidate for a local vector. Then you'd check the vec global size > matches the dmda local size. If the commsize is anything other than 1 > then it cannot be a local vector > > > Matt > > > You can check the size matches your local DMDA space by using > DMDAGetGhostCorners() > > > https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DMDA/DMDAGetGhostCorners.html > > and return the quantities m, n, and p. > > You also need to use DMDAGetInfo() > > > https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DMDA/DMDAGetInfo.html > > The important quantity you want returned is "dof" > > If m x n x p x dof matches the number returned by VecGetSize() > (assuming you know the vector is sequential) then you know the > local space will fit within your vector. > > > > > b) the DM is a DMDA rather than some other subclass > > > See Matt's answer > > > > I can't immediately see routines that do what I need, but > I am likely > missing something obvious. Is there a way to achieve the > above? > > Thanks > > Phil > > -- > Phil Tooley > Research Software Engineering > University of Sheffield > > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to > which their experiments lead. > -- Norbert Wiener > > https://www.cse.buffalo.edu/~knepley/ > <http://www.cse.buffalo.edu/%7Eknepley/> > -- Phil Tooley Research Software Engineering University of Sheffield
Re: [petsc-users] Checking if a vector is a localvector of a given DMDA
Hi Matt, Thanks for the rapid response. The reason I want to ensure I have a DMDA is that this is a custom gradient function, the results of which are meaningless for anything but a regular grid. I would rather raise an exception upfront with the user than allow the program to continue with garbage data. Phil On 25/09/18 11:55, Matthew Knepley wrote: > On Tue, Sep 25, 2018 at 6:49 AM Phil Tooley > mailto:phil.too...@sheffield.ac.uk>> wrote: > > Hi all, > > Given a vector I know I can get an associated DM (if there is one) by > calling VecGetDM, but I need to also be able to check that > > a) the vector is the localvector of that DM rather than the global > > > Get a local vector from the DM and check the sizes. > > > b) the DM is a DMDA rather than some other subclass > > > You can do this (PetscObjectTypeCompare), but it is not recommended. > Why would > you want to do that? > > Thanks, > > Matt > > > I can't immediately see routines that do what I need, but I am likely > missing something obvious. Is there a way to achieve the above? > > Thanks > > Phil > > -- > Phil Tooley > Research Software Engineering > University of Sheffield > > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which > their experiments lead. > -- Norbert Wiener > > https://www.cse.buffalo.edu/~knepley/ > <http://www.cse.buffalo.edu/%7Eknepley/> -- Phil Tooley Research Software Engineering University of Sheffield
[petsc-users] Checking if a vector is a localvector of a given DMDA
Hi all, Given a vector I know I can get an associated DM (if there is one) by calling VecGetDM, but I need to also be able to check that a) the vector is the localvector of that DM rather than the global b) the DM is a DMDA rather than some other subclass I can't immediately see routines that do what I need, but I am likely missing something obvious. Is there a way to achieve the above? Thanks Phil -- Phil Tooley Research Software Engineering University of Sheffield
[petsc-users] DMCreateXXXVector vs VecDuplicate
Hi all, Can I get come clarification regarding DMCreateXXXVector vs VecDuplicate. The manual says that I can obtain additional vectors with the latter, does this mean that I *must* only call DMCreateXXXVector once each for global and local? I ask because I want multiple instances of an image class all with the same distributed layout. It makes life much tidier if I can just pass the DMDA rather than also having to pass vectors to be copied. Cheers Phil -- Phil Tooley Research Software Engineering University of Sheffield
Re: [petsc-users] Indexing/using a 3D DMDA as a 1D vector
Ok, So I have figured this out a bit better now, facepalm moment as I finally realised that this is exactly what the DMDA exists to do... I need to create a vector using the DMCreate*Vector functions and can then insert values using the assembly functions. Then, if I understand correctly, in order to have the state I need to multiply with my operator I should make use of the DMDAGlobalToNaturalBegin() and End() functions to ensure the correct ordering. Do I then need to ensure I have reverted it to global ordering before I try to apply finite difference methods again? Thanks P On 12/09/18 14:27, Phil Tooley wrote: > I will preface this by saying I am new to PETSc and am still trying to > get my head around all of the layout mapping that is done. That means I > may well have fundamentally misunderstood something, but hopefully > someone will be able to to put me right. > > In my application I have some 3D pixel data which I want to manipulate > using finite difference methods and then transform by viewing as a 1-D > vector and multiplying by a large sparse matrix operator. > > I would assume that the correct way to do this is by creating a DMDA to > hold the image data and ghosting appropriately to apply my finite > difference operations. Then I had hoped that I could use some form of > application ordering to allow viewing the data as a vector that can be > multiplied with my operator matrix. This is where I have come unstuck, > I may just be missing something obivous but I can't figure out how to do > this. Can anyone point me in the correct direction please? > > Many Thanks > -- Phil Tooley Research Software Engineering University of Sheffield
[petsc-users] Indexing/using a 3D DMDA as a 1D vector
I will preface this by saying I am new to PETSc and am still trying to get my head around all of the layout mapping that is done. That means I may well have fundamentally misunderstood something, but hopefully someone will be able to to put me right. In my application I have some 3D pixel data which I want to manipulate using finite difference methods and then transform by viewing as a 1-D vector and multiplying by a large sparse matrix operator. I would assume that the correct way to do this is by creating a DMDA to hold the image data and ghosting appropriately to apply my finite difference operations. Then I had hoped that I could use some form of application ordering to allow viewing the data as a vector that can be multiplied with my operator matrix. This is where I have come unstuck, I may just be missing something obivous but I can't figure out how to do this. Can anyone point me in the correct direction please? Many Thanks -- Phil Tooley Research Software Engineering University of Sheffield