<div dir="ltr"><div dir="ltr">On Tue, Jan 14, 2020 at 6:58 AM Y. Shidi <<a href="mailto:ys453@cam.ac.uk">ys453@cam.ac.uk</a>> wrote:<br></div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">Thanks Matt,<br>
<br>
Instead of constructing the nested vector by VecCreateNest(), I use<br>
VecRestoreSubVector() to solve this issue.<br>
<br>
But I got problems for field splitting method,<br>
I use the following options:<br>
#PETSc Option Table entries:<br>
-fieldsplit_0_ksp_converged_reason true<br>
-fieldsplit_0_ksp_rtol 1e-12<br>
-fieldsplit_0_ksp_type cg<br>
-fieldsplit_0_pc_type ilu<br>
-fieldsplit_1_ksp_converged_reason true<br>
-fieldsplit_1_ksp_rtol 1e-12<br>
-fieldsplit_1_ksp_type cg<br>
-fieldsplit_1_pc_type hypre<br>
-ksp_converged_reason<br>
-ksp_rtol 1e-12<br>
-ksp_type fgmres<br>
-ksp_view<br>
-pc_fieldsplit_schur_fact_type upper<br>
-pc_fieldsplit_schur_precondition seflp<br>
-pc_type fieldsplit<br>
#End of PETSc Option Table entries<br>
<br>
There are some options that were not used.<br>
<br>
WARNING! There are options you set that were not used!<br>
WARNING! could be spelling mistake, etc!<br>
Option left: name:-fieldsplit_0_ksp_converged_reason value: true<br>
Option left: name:-fieldsplit_1_ksp_converged_reason value: true<br>
Option left: name:-pc_fieldsplit_schur_fact_type value: upper<br>
Option left: name:-pc_fieldsplit_schur_precondition value: seflp<br></blockquote><div><br></div><div>PCFIELDSPLIT has to know how to divide your vector into fields. VecNest provides that information</div><div>since it is divided into the nested pieces. Other ways you can provide this:</div><div><br></div><div>  1) Use a DM to describe the data layout. This is the best in my opinion.</div><div><br></div><div>  2) Fix the VecNest. If an overwrite is happening, its best to track it down.</div><div><br></div><div>  3) Provide an IS to PCFIELDSPLIT using PCFieldSplitSetIS().</div><div><br></div><div>  4) Split by block component</div><div><br></div><div>  Thanks,</div><div><br></div><div>    Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
Kind Regards,<br>
Shidi<br>
<br>
On 2020-01-14 11:53, Matthew Knepley wrote:<br>
> It says that memory is being overwritten somewhere. You can track this<br>
> down using valgrind,<br>
> as it suggests in the error message.<br>
> <br>
>   Thanks,<br>
> <br>
>      Matt<br>
> <br>
> On Tue, Jan 14, 2020 at 5:44 AM Y. Shidi <<a href="mailto:ys453@cam.ac.uk" target="_blank">ys453@cam.ac.uk</a>> wrote:<br>
> <br>
>> Dear developers,<br>
>> <br>
>> I have a 2x2 nested matrix and the corresponding nested vector.<br>
>> When I running the code with field splitting, it gets the following<br>
>> errors:<br>
>> <br>
>> [0]PETSC ERROR: PetscTrFreeDefault() called from<br>
>> VecRestoreArray_Nest()<br>
>> line 678 in<br>
>> /home/ys453/Sources/petsc/src/vec/vec/impls/nest/vecnest.c<br>
>> [0]PETSC ERROR: Block at address 0x3f95f60 is corrupted; cannot<br>
>> free;<br>
>> may be block not allocated with PetscMalloc()<br>
>> [0]PETSC ERROR: --------------------- Error Message<br>
>> --------------------------------------------------------------<br>
>> [0]PETSC ERROR: Memory corruption:<br>
>> <br>
> <a href="http://www.mcs.anl.gov/petsc/documentation/installation.html#valgrind" rel="noreferrer" target="_blank">http://www.mcs.anl.gov/petsc/documentation/installation.html#valgrind</a><br>
>> [0]PETSC ERROR: Bad location or corrupted memory<br>
>> [0]PETSC ERROR: See<br>
>> <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" rel="noreferrer" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a><br>
>> for trouble shooting.<br>
>> [0]PETSC ERROR: Petsc Release Version 3.9.3, unknown<br>
>> [0]PETSC ERROR: 2DPetscSpuriousTest on a arch-linux2-c-debug named<br>
>> merlin by ys453 Tue Jan 14 10:36:53 2020<br>
>> [0]PETSC ERROR: Configure options --download-scalapack<br>
>> --download-mumps<br>
>> --download-parmetis --download-metis --download-ptscotch<br>
>> --download-superlu_dist --download-hypre<br>
>> [0]PETSC ERROR: #1 PetscTrFreeDefault() line 269 in<br>
>> /home/ys453/Sources/petsc/src/sys/memory/mtr.c<br>
>> [0]PETSC ERROR: #2 VecRestoreArray_Nest() line 678 in<br>
>> /home/ys453/Sources/petsc/src/vec/vec/impls/nest/vecnest.c<br>
>> [0]PETSC ERROR: #3 VecRestoreArrayRead() line 1835 in<br>
>> /home/ys453/Sources/petsc/src/vec/vec/interface/rvector.c<br>
>> [0]PETSC ERROR: #4 VecRestoreArrayPair() line 511 in<br>
>> /home/ys453/Sources/petsc/include/petscvec.h<br>
>> [0]PETSC ERROR: #5 VecScatterBegin_SSToSS() line 671 in<br>
>> /home/ys453/Sources/petsc/src/vec/vscat/impls/vscat.c<br>
>> [0]PETSC ERROR: #6 VecScatterBegin() line 1779 in<br>
>> /home/ys453/Sources/petsc/src/vec/vscat/impls/vscat.c<br>
>> [0]PETSC ERROR: #7 PCApply_FieldSplit() line 1010 in<br>
>> /home/ys453/Sources/petsc/src/ksp/pc/impls/fieldsplit/fieldsplit.c<br>
>> [0]PETSC ERROR: #8 PCApply() line 457 in<br>
>> /home/ys453/Sources/petsc/src/ksp/pc/interface/precon.c<br>
>> [0]PETSC ERROR: #9 KSP_PCApply() line 276 in<br>
>> /home/ys453/Sources/petsc/include/petsc/private/kspimpl.h<br>
>> [0]PETSC ERROR: #10 KSPFGMRESCycle() line 166 in<br>
>> /home/ys453/Sources/petsc/src/ksp/ksp/impls/gmres/fgmres/fgmres.c<br>
>> [0]PETSC ERROR: #11 KSPSolve_FGMRES() line 291 in<br>
>> /home/ys453/Sources/petsc/src/ksp/ksp/impls/gmres/fgmres/fgmres.c<br>
>> [0]PETSC ERROR: #12 KSPSolve() line 669 in<br>
>> /home/ys453/Sources/petsc/src/ksp/ksp/interface/itfunc.c<br>
>> <br>
>> I am not sure why it happens.<br>
>> <br>
>> Thank you for your time.<br>
>> <br>
>> Kind Regards,<br>
>> Shidi<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<br>
> their experiments lead.<br>
> -- Norbert Wiener<br>
> <br>
> <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a> [1]<br>
> <br>
> <br>
> Links:<br>
> ------<br>
> [1] <a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a><br>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>