<div dir="ltr"><br><div class="gmail_extra"><br><div class="gmail_quote">On 4 August 2016 at 10:10, Patrick Sanan <span dir="ltr"><<a href="mailto:patrick.sanan@gmail.com" target="_blank">patrick.sanan@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">I have a patch that I got from Dave that he got from Jed which seems<br>
to be related to this. I'll make a PR.<br></blockquote><div><br></div><div>Jed wrote this variant of the VTK viewer so please mark him as a reviewer for my bug fix.<br></div><br> <blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div class="HOEnZb"><div class="h5"><br>
On Wed, Aug 3, 2016 at 8:18 PM, Mohammad Mirzadeh <<a href="mailto:mirzadeh@gmail.com">mirzadeh@gmail.com</a>> wrote:<br>
> OK so I just ran the example under valgrind, and if I use two VecViews, it<br>
> complains about following leak:<br>
><br>
> ==66838== 24,802 (544 direct, 24,258 indirect) bytes in 1 blocks are<br>
> definitely lost in loss record 924 of 926<br>
> ==66838==    at 0x100009EBB: malloc (in<br>
> /usr/local/Cellar/valgrind/3.11.0/lib/valgrind/vgpreload_memcheck-amd64-darwin.so)<br>
> ==66838==    by 0x10005E638: PetscMallocAlign (in<br>
> /usr/local/Cellar/petsc/3.7.2/real/lib/libpetsc.3.7.2.dylib)<br>
> ==66838==    by 0x100405F00: DMCreate_DA (in<br>
> /usr/local/Cellar/petsc/3.7.2/real/lib/libpetsc.3.7.2.dylib)<br>
> ==66838==    by 0x1003CFFA4: DMSetType (in<br>
> /usr/local/Cellar/petsc/3.7.2/real/lib/libpetsc.3.7.2.dylib)<br>
> ==66838==    by 0x100405B7F: DMDACreate (in<br>
> /usr/local/Cellar/petsc/3.7.2/real/lib/libpetsc.3.7.2.dylib)<br>
> ==66838==    by 0x1003F825F: DMDACreate2d (in<br>
> /usr/local/Cellar/petsc/3.7.2/real/lib/libpetsc.3.7.2.dylib)<br>
> ==66838==    by 0x100001D89: main (main_test.cpp:7)<br>
><br>
> By I am destroying the dm ... also I dont get this when using a single<br>
> VecView. As a bonus info, PETSC_VIEWER_STDOUT_WORLD is just fine, so this<br>
> looks like it is definitely vtk related.<br></div></div></blockquote><div><br></div><div>Mohammad,<br></div><div><br>I can confirm that this VTK functionality bleeds memory if you write more than 1 vector to disk.<br></div><div><br></div><div>Cheers<br></div><div>  Dave<br></div><div><br></div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="HOEnZb"><div class="h5">
><br>
> On Wed, Aug 3, 2016 at 2:05 PM, Mohammad Mirzadeh <<a href="mailto:mirzadeh@gmail.com">mirzadeh@gmail.com</a>><br>
> wrote:<br>
>><br>
>> On Wed, Aug 3, 2016 at 10:59 AM, Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>><br>
>> wrote:<br>
>>><br>
>>> On Tue, Aug 2, 2016 at 12:40 PM, Mohammad Mirzadeh <<a href="mailto:mirzadeh@gmail.com">mirzadeh@gmail.com</a>><br>
>>> wrote:<br>
>>>><br>
>>>> I often use the memory usage information in log_view as a way to check<br>
>>>> on memory leaks and so far it has worked perfect. However, I had long<br>
>>>> noticed a false-positive report in memory leak for Viewers, i.e. destruction<br>
>>>> count is always one less than creation.<br>
>>><br>
>>><br>
>>> Yes, I believe that is the Viewer being used to print this information.<br>
>><br>
>><br>
>> That makes sense.<br>
>>><br>
>>><br>
>>>><br>
>>>> Today, I noticed what seems to be a second one. If you use VecView to<br>
>>>> write the same DA to vtk, i.e. call VecView(A, vtk); twice, it also report a<br>
>>>> memory leak for vectors, vecscatters, dm, etc. I am calling this a<br>
>>>> false-positive since the code is valgrind-clean.<br>
>>>><br>
>>>> Is this known/expected?<br>
>>><br>
>>><br>
>>> The VTK viewers have to hold everything they output until they are<br>
>>> destroyed since the format does not allow immediate writing.<br>
>>> I think the VTK viewer is not destroyed at the time of this output. Can<br>
>>> you make a small example that does this?<br>
>><br>
>><br>
>> Here's a small example that illustrates the issues<br>
>><br>
>> #include <petsc.h><br>
>><br>
>><br>
>> int main(int argc, char *argv[]) {<br>
>><br>
>>   PetscInitialize(&argc, &argv, NULL, NULL);<br>
>><br>
>><br>
>>   DM dm;<br>
>><br>
>>   DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE,<br>
>> DMDA_STENCIL_BOX,<br>
>><br>
>>                -10, -10, PETSC_DECIDE, PETSC_DECIDE, 1, 1, NULL, NULL,<br>
>><br>
>>                &dm);<br>
>><br>
>> //  DMDASetUniformCoordinates(dm, -1, 1, -1, 1, 0, 0);<br>
>><br>
>><br>
>>   Vec sol;<br>
>><br>
>>   DMCreateGlobalVector(dm, &sol);<br>
>><br>
>>   VecSet(sol, 0);<br>
>><br>
>><br>
>>   PetscViewer vtk;<br>
>><br>
>>   PetscViewerVTKOpen(PETSC_COMM_WORLD, "test.vts", FILE_MODE_WRITE, &vtk);<br>
>><br>
>>   VecView(sol, vtk);<br>
>><br>
>> //  VecView(sol, vtk);<br>
>><br>
>>   PetscViewerDestroy(&vtk);<br>
>><br>
>><br>
>>   DMDestroy(&dm);<br>
>><br>
>>   VecDestroy(&sol);<br>
>><br>
>><br>
>>   PetscFinalize();<br>
>><br>
>>   return 0;<br>
>><br>
>> }<br>
>><br>
>><br>
>> If you uncomment the second VecView you get reports for leaks in<br>
>> VecScatter and dm. If you also uncomment the DMDASetUniformCoordinates, and<br>
>> use both VecViews, you also get a leak report for Vecs ... its quite bizarre<br>
>> ...<br>
>><br>
>>><br>
>>> I have switched to HDF5 and XDMF due to the limitations of VTK format.<br>
>>><br>
>><br>
>> I had used XDMF + raw binary in the past and was satisfied with the<br>
>> result. Do you write a single XDMF as a "post-processing" step when the<br>
>> simulation is finished? If I remember correctly preview could not open xmf<br>
>> files as time-series.<br>
>>><br>
>>>   Thanks,<br>
>>><br>
>>>      Matt<br>
>>><br>
>>>><br>
>>>> Here's the relevant bit from log_view:<br>
>>>><br>
>>>> --- Event Stage 0: Main Stage<br>
>>>><br>
>>>>               Vector     8              7       250992     0.<br>
>>>>       Vector Scatter     2              0            0     0.<br>
>>>>     Distributed Mesh     2              0            0     0.<br>
>>>> Star Forest Bipartite Graph     4              0            0     0.<br>
>>>>      Discrete System     2              0            0     0.<br>
>>>>            Index Set     4              4        83136     0.<br>
>>>>    IS L to G Mapping     2              0            0     0.<br>
>>>>               Viewer     2              1          784     0.<br>
>>>><br>
>>>> ========================================================================================================================<br>
>>>><br>
>>><br>
>>><br>
>>><br>
>>> --<br>
>>> What most experimenters take for granted before they begin their<br>
>>> experiments is infinitely more interesting than any results to which their<br>
>>> experiments lead.<br>
>>> -- Norbert Wiener<br>
>><br>
>><br>
><br>
</div></div></blockquote></div><br></div></div>