<div>Sorry about the last email...I hit sent it too early. </div><div><br></div><div>The full calling sequence is:</div><div><br></div><div><div><div><br class="Apple-interchange-newline"> call KSPCreate(SUMB_COMM_WORLD, ksp, PETScIerr)</div>
<div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> useAD = .False.</div><div> useTranspose = .True.</div><div> usePC = .True.</div><div> call setupStateResidualMatrix(drdwpret,useAD,usePC,useTranspose)</div>
<div><br></div><div> ! dRdwT has already been created</div><div> call KSPSetOperators(ksp, dRdWT, dRdwpreT, &</div><div> DIFFERENT_NONZERO_PATTERN, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div>
<div><br></div><div> call KSPSetFromOptions(ksp, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call KSPSetType(ksp, KSPFGMRES, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div>
<div><br></div><div> call KSPGMRESSetRestart(ksp, adjRestart, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call KSPGMRESSetCGSRefinementType(ksp, &</div><div> KSP_GMRES_CGS_REFINE_IFNEEDED, PETScIerr)</div>
<div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call kspgetpc(ksp, pc, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call pcsettype(pc, PCMG, PETScIerr)</div>
<div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call PetscOptionsSetValue('-pc_mg_levels','2',PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div></div><div><div><br>
</div><div> ! Set command line options </div>
<div> call PCSetFromOptions(pc, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div class="im"><div><br></div><div> call PCMGSetGalerkin(pc, PETSC_TRUE, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div>
<div><br></div></div><div> ! Create the restriction operator between finest and one below: </div>
<div> call createRestrictOperator(RL1, 1)</div><div><br></div><div class="im"> call PCMGSetRestriction(pc, 1, RL1, PETScierr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call PCSetup(pc, PETScierr)</div>
<div> call EChk(PETScIerr,__FILE__,__LINE__)</div></div><div class="HOEnZb"><div class="h5"></div></div></div><div><br></div><div> The error occurs in PCSetup() which calls PCSetup_MG() and the issue is in KSPGetVecs() on line 644 in mg.c. Just to clarify, this isn't a memory problem/segfault. PETSc returns control cleanly and I my code exits when the error is detected. Therefore a backtrace or valgrind doesn't yield any additional information. ( I Ran both). The terminal output is:</div>
<div><br></div><div><div>[0]PETSC ERROR: ------------------------------------------------------------------------</div><div>[0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, probably memory access out of range</div>
<div>[0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger</div><div>[0]PETSC ERROR: or see <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[0]PETSC">http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[0]PETSC</a> ERROR: or try <a href="http://valgrind.org">http://valgrind.org</a> on GNU/linux and Apple Mac OS X to find memory corruption errors</div>
<div>[0]PETSC ERROR: likely location of problem given in stack below</div><div>[0]PETSC ERROR: --------------------- Stack Frames ------------------------------------</div><div> =================================================================</div>
<div>PETSc or MPI Error. Error Code 73. Detected on Proc 0</div><div>Error at line: 130 in file: setupPETScKsp.F90</div><div> =================================================================</div><div>[0]PETSC ERROR: Note: The EXACT line numbers in the stack are not available,</div>
<div>[0]PETSC ERROR: INSTEAD the line number of the start of the function</div><div>[0]PETSC ERROR: is given.</div><div>[0]PETSC ERROR: [0] MatGetVecs line 8142 /nfs/mica/home/kenway/Downloads/petsc-3.3/src/mat/interface/matrix.c</div>
<div>[0]PETSC ERROR: [0] KSPGetVecs line 774 /nfs/mica/home/kenway/Downloads/petsc-3.3/src/ksp/ksp/interface/iterativ.c</div><div>[0]PETSC ERROR: [0] PCSetUp_MG line 508 /nfs/mica/home/kenway/Downloads/petsc-3.3/src/ksp/pc/impls/mg/mg.c</div>
<div>[0]PETSC ERROR: [0] PCSetUp line 810 /nfs/mica/home/kenway/Downloads/petsc-3.3/src/ksp/pc/interface/precon.c</div></div><div><br></div>I had a look at ex42.c but it uses a DM which I am not doing. From the example there is the following:<div>
<br></div><div><div>ierr = PCMGSetLevels(pc,nlevels,PETSC_NULL);CHKERRQ(ierr);</div><div> ierr = PCMGSetType(pc,PC_MG_MULTIPLICATIVE);CHKERRQ(ierr);</div><div> ierr = PCMGSetGalerkin(pc,PETSC_TRUE);CHKERRQ(ierr);</div><div>
<br></div><div> for (k=1; k<nlevels; k++){</div><div> ierr = DMCreateInterpolation(da_list[k-1],da_list[k],&R,PETSC_NULL);CHKERRQ(ierr);</div><div> ierr = PCMGSetInterpolation(pc,k,R);CHKERRQ(ierr);</div><div>
ierr = MatDestroy(&R);CHKERRQ(ierr);</div><div> }</div><div><br></div><div>Is there anything special about "R" in this case? As I mentioned before, I know what the issue is, the coarse grid operators are not set at all. From what I can tell, they should be set in the routine PCSetup_MG in mg.c, but the code is never run unless pc->dm is True. I guess I'll just create the operators myself. Even if I modify the code to make it enter the if statement at line 577</div>
<div><br></div><div><div> if (mg->galerkin == 1) {</div></div><div><div> Mat B;</div><div> /* currently only handle case where mat and pmat are the same on coarser levels */</div><div> ierr = KSPGetOperators(mglevels[n-1]->smoothd,&dA,&dB,&uflag);CHKERRQ(ierr);</div>
<div> if (!pc->setupcalled) {</div><div> for (i=n-2; i>-1; i--) {</div><div> ierr = MatPtAP(dB,mglevels[i+1]->interpolate,MAT_INITIAL_MATRIX,1.0,&B);CHKERRQ(ierr);</div><div> ierr = KSPSetOperators(mglevels[i]->smoothd,B,B,uflag);CHKERRQ(ierr);</div>
<div> if (i != n-2) {ierr = PetscObjectDereference((PetscObject)dB);CHKERRQ(ierr);}</div><div> dB = B;</div><div> }</div><div> if (n > 1) {ierr = PetscObjectDereference((PetscObject)dB);CHKERRQ(ierr);}</div>
<div> } else {</div><div> for (i=n-2; i>-1; i--) {</div><div> ierr = KSPGetOperators(mglevels[i]->smoothd,PETSC_NULL,&B,PETSC_NULL);CHKERRQ(ierr);</div><div> ierr = MatPtAP(dB,mglevels[i+1]->interpolate,MAT_REUSE_MATRIX,1.0,&B);CHKERRQ(ierr);</div>
<div> ierr = KSPSetOperators(mglevels[i]->smoothd,B,B,uflag);CHKERRQ(ierr);</div><div> dB = B;</div><div> }</div><div> }</div></div><div><br></div><div>which seems like where it should be done....the MatPTAP will not run. I'm using a mpibaij matrix for 'dRdwPreT' and my restriction matrix 'RL1' is a mpiaij matrix. When I try to run MatRART myself with the code:</div>
<div><br></div><div><div> call MatRART(dRdwT, RL1, MAT_INITIAL_MATRIX, 1.0, newMat, PETScIERR)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div>I receive error code 56: /* no support for requested operation */ (this is a clean exit from PETSc)</div>
</div><div><br></div><div>In conclusion, I have two issues: 1) The coarse grid operators are not being set automatically in mg.c. 2) Using MatRART() with a mpibaij matrix and a mpiaij matrix is unsupported. </div><div><br>
</div><div>Is there a PCMG example that doesn't use a DM to generate the restriction/prolongation/matrix information I can have a look at?</div><div><br></div><div>Thanks,</div><div>Gaetan</div><div><br></div><div><br>
</div><div><br></div><div><br></div><br><div class="gmail_quote">On Sat, Sep 1, 2012 at 4:54 PM, Gaetan Kenway <span dir="ltr"><<a href="mailto:kenway@utias.utoronto.ca" target="_blank">kenway@utias.utoronto.ca</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">The full calling sequence is:<div><div><br></div><div><br></div><div> call KSPCreate(SUMB_COMM_WORLD, ksp, PETScIerr)</div>
<div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> useAD = .False.</div><div>
useTranspose = .True.</div><div> usePC = .True.</div><div> call setupStateResidualMatrix(drdwpret,useAD,usePC,useTranspose)</div><div><br></div><div> call KSPSetOperators(ksp, dRdWT, dRdwpreT, &</div>
<div> DIFFERENT_NONZERO_PATTERN, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call KSPSetFromOptions(ksp, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div>
<div><br></div><div> call KSPSetType(ksp, KSPFGMRES, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call KSPGMRESSetRestart(ksp, adjRestart, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div>
<div><br></div><div> call KSPGMRESSetCGSRefinementType(ksp, &</div><div> KSP_GMRES_CGS_REFINE_IFNEEDED, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call kspgetpc(ksp, pc, PETScIerr)</div>
<div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call pcsettype(pc, PCMG, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call PetscOptionsSetValue('-pc_mg_levels','2',PETScIerr)</div>
<div> call EChk(PETScIerr,__FILE__,__LINE__)</div></div><div><div><br></div><div> ! Set command line options </div>
<div> call PCSetFromOptions(pc, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div class="im"><div><br></div><div> call PCMGSetGalerkin(pc, PETSC_TRUE, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div>
<div><br></div></div><div> ! Create the restriction operator between finest and one below: </div>
<div> call createRestrictOperator(RL1, 1)</div><div><br></div><div class="im"><div> call PCMGSetRestriction(pc, 1, RL1, PETScierr)</div></div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div> call PCSetup(pc, PETScierr)</div>
<div> call EChk(PETScIerr,__FILE__,__LINE__)</div></div><div class="HOEnZb"><div class="h5"><div><br></div><div><br><br><div class="gmail_quote">On Sat, Sep 1, 2012 at 3:27 PM, Jed Brown <span dir="ltr"><<a href="mailto:five9a2@gmail.com" target="_blank">five9a2@gmail.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><p>Can you send the whole calling sequence and the full error message? It's not clear from your email how to reproduce. You might also look at src/ksp/ksp/examples/tutorials/ex42.c which uses PCMGSetGalerkin. Did you set the fine grid operators before calling PCSetUp()? The coarse grid operators should have been constructed earlier in PCSetUp_MG().</p>
<div><div>
<div class="gmail_quote">On Sep 1, 2012 10:05 AM, "Gaetan Kenway" <<a href="mailto:kenway@utias.utoronto.ca" target="_blank">kenway@utias.utoronto.ca</a>> wrote:<br type="attribution"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
I believe I partially tracked down my problem. The issue arises at line 644 in mg.c with a call to KSPGetVecs(). This is used to get the coarse grid vectors for the RHS. The issue is that up until this point, the actual operator for mglevels[i]->smoothd have not actually been set and therefore it cannot get the vectors. <div>
<br></div><div>I am using </div><div> call PCMGSetGalerkin(pc, PETSC_TRUE, PETScIerr)</div><div> call EChk(PETScIerr,__FILE__,__LINE__)</div><div><br></div><div>with the intention of having the coarse matrices approximately automatically. However, from looking over the code, it appears when it is used in this manner, the coarse grid operators never get set. From what I can tell, it is only possible to create these matrices automatically if (pc->dm) is True which I think means that the preconditioner is derived from a distributed matrix. </div>
<div><br></div><div>I've set my own restriction operator, but mg->galerkin has a value of -1 in mg.c and none of the KSPSetOperators() are executed. </div><div><br></div><div>I've also tried to generate the matrices myself using MatRARt() or MatPtAP, but I receive error code 56 "<strong><font color="#228B22"> </font><font color="#B22222">/* no support for requested operation */" . </font></strong>Is this op only implemented from matrices derived from from DM objects?</div>
<div><br></div><div>Any other suggestions for creating the coarse grid operators manually?</div><div><br></div><div>Thank you,</div><div><br></div><div>Gaetan</div><div><br><br><div class="gmail_quote">On Sat, Sep 1, 2012 at 12:29 AM, Jed Brown <span dir="ltr"><<a href="mailto:jedbrown@mcs.anl.gov" target="_blank">jedbrown@mcs.anl.gov</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="gmail_quote"><div><div>On Fri, Aug 31, 2012 at 11:22 PM, Gaetan Kenway <span dir="ltr"><<a href="mailto:kenway@utias.utoronto.ca" target="_blank">kenway@utias.utoronto.ca</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Hi Again<div><br></div><div>I also tried petsc-3.2 version and I still get the same backtrace. </div><div><div><br></div><div>If its not possible to figure out where the communicator segfault is coming from its not a huge deal...I've just set the option using PetscOptionsSetValue() and then use PCSetFromOptions() to pull it back out. That seems to work fine. </div>
<div><br></div><div>Even avoiding the above problem, with the PetscOptionsSetValue I'm still receiving an error code 73 when I run the multigrid solver. I've included the backtrace output below but its not a lot of help since the code exited cleaning using my error checking procedure</div>
<div><br></div><div><div> =================================================================</div><div>PETSc or MPI Error. Error Code 73. Detected on Proc 0</div><div>Error at line: 122 in file: solveADjointTransposePETSc.F90</div>
<div> =================================================================</div><div><br></div><div>Program received signal SIGSEGV, Segmentation fault.</div><div>0xb68c98f0 in PMPI_Abort () from /usr/local/lib/libmpi.so.0</div>
<div>(gdb) bt</div><div>#0 0xb68c98f0 in PMPI_Abort () from /usr/local/lib/libmpi.so.0</div><div>#1 0xb519f190 in pmpi_abort__ () from /usr/local/lib/libmpi_f77.so.0</div><div>#2 0xb44c9e4c in echk (ierr=@0x49, file=..., line=@0x7a, .tmp.FILE.len_V$eb=30) at terminate.f90:154</div>
<div>#3 0xb44bd68f in solveadjointtransposepetsc () at solveADjointTransposePETSc.F90:122</div><div>#4 0xb44138a9 in f2py_rout_sumb_solveadjointtransposepetsc () from /tmp/tmpKYF_DT/sumb.so</div><div>#5 0xb440fd35 in fortran_call () from /tmp/tmpKYF_DT/sumb.so</div>
<div>#6 0x0805fd6a in PyObject_Call ()</div><div>#7 0x080dd5b0 in PyEval_EvalFrameEx ()</div><div>#8 0x080dfbb2 in PyEval_EvalCodeEx ()</div><div>#9 0x08168f1f in ?? ()</div><div>#10 0x0805fd6a in PyObject_Call ()</div>
<div>#11 0x080dcbeb in PyEval_EvalFrameEx ()</div><div>#12 0x080dfbb2 in PyEval_EvalCodeEx ()</div><div>#13 0x080de145 in PyEval_EvalFrameEx ()</div><div>#14 0x080dfbb2 in PyEval_EvalCodeEx ()</div><div>#15 0x080dfca7 in PyEval_EvalCode ()</div>
<div>#16 0x080fd956 in PyRun_FileExFlags ()</div><div>#17 0x080fdbb2 in PyRun_SimpleFileExFlags ()</div><div>#18 0x0805b6d3 in Py_Main ()</div><div>#19 0x0805a8ab in main ()</div></div></div></blockquote><div><br></div></div>
</div><div>
This stack doesn't involve PETSc at all.</div><div><div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div><div><br></div><div>Valgrid was clean right up until the end where I get the normal error message:</div>
<div><br></div><div><div><div>[0]PETSC ERROR: ------------------------------------------------------------------------</div><div>[0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, probably memory access out of range</div>
<div>[0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger</div><div>[0]PETSC ERROR: or see <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[0]PETSC" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[0]PETSC</a> ERROR: or try <a href="http://valgrind.org" target="_blank">http://valgrind.org</a> on GNU/linux and Apple Mac OS X to find memory corruption errors</div>
<div>[0]PETSC ERROR: likely location of problem given in stack below</div><div>[0]PETSC ERROR: --------------------- Stack Frames ------------------------------------</div><div>[0]PETSC ERROR: Note: The EXACT line numbers in the stack are not available,</div>
<div>[0]PETSC ERROR: INSTEAD the line number of the start of the function</div><div>[0]PETSC ERROR: is given.</div></div><div><div>[0]PETSC ERROR: [0] MatGetVecs line 8142 /nfs/mica/home/kenway/Downloads/petsc-3.3/src/mat/interface/matrix.c</div>
<div>[0]PETSC ERROR: [0] KSPGetVecs line 774 /nfs/mica/home/kenway/Downloads/petsc-3.3/src/ksp/ksp/interface/iterativ.c</div><div>[0]PETSC ERROR: [0] PCSetUp_MG line 508 /nfs/mica/home/kenway/Downloads/petsc-3.3/src/ksp/pc/impls/mg/mg.c</div>
<div>[0]PETSC ERROR: [0] PCSetUp line 810 /nfs/mica/home/kenway/Downloads/petsc-3.3/src/ksp/pc/interface/precon.c</div><div>[0]PETSC ERROR: [0] KSPSetUp line 182 /nfs/mica/home/kenway/Downloads/petsc-3.3/src/ksp/ksp/interface/itfunc.c</div>
<div>[0]PETSC ERROR: [0] KSPSolve line 351 /nfs/mica/home/kenway/Downloads/petsc-3.3/src/ksp/ksp/interface/itfunc.c</div></div></div><div><br></div><div>I will try it on a different system tomorrow to see if I have any more luck.</div>
<div><br></div><div>Thanks,</div><div><br></div><div>Gaetan</div><div><div>
<div><br></div><div><br></div><br><div class="gmail_quote">On Fri, Aug 31, 2012 at 11:08 PM, Jed Brown <span dir="ltr"><<a href="mailto:jedbrown@mcs.anl.gov" target="_blank">jedbrown@mcs.anl.gov</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div><div class="gmail_quote">On Fri, Aug 31, 2012 at 10:06 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div>True, but the backtrace also shows that comm = 0x0 on the call to KSPCreate(), which<div>leads me to believe that your petsc4py has not initialized PETSc, and therefor not</div><div>initialized PETSC_COMM_WORLD.</div>
</div></blockquote></div><br></div><div>Could catch, maybe PetscFunctionBegin() should check that PETSc has been initialized.</div>
</blockquote></div><br></div></div></div>
</blockquote></div></div></div><br>
</blockquote></div><br></div>
</blockquote></div>
</div></div></blockquote></div><br></div>
</div></div></blockquote></div><br></div>