<html><head><meta http-equiv="content-type" content="text/html; charset=us-ascii"></head><body style="overflow-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;"><div><br></div>It is in the loop below that it errors due to an unassembled vector. T<br><div><br></div><div><div>static PetscErrorCode VecCopy_Nest(Vec x, Vec y)</div><div>{</div><div> Vec_Nest *bx = (Vec_Nest *)x->data;</div><div> Vec_Nest *by = (Vec_Nest *)y->data;</div><div> PetscInt i;</div><div><br></div><div> PetscFunctionBegin;</div><div> VecNestCheckCompatible2(x, 1, y, 2);</div><div> for (i = 0; i < bx->nb; i++) PetscCall(VecCopy(bx->v[i], by->v[i]));</div><div> PetscFunctionReturn(PETSC_SUCCESS);</div><div>}</div><div><br></div><div>thus it is one of the individual vectors that make up the full vectors that is marked as not assembled</div><div><br></div><div>In ex70.c I see </div><div><br></div><div><div>PetscErrorCode StokesExactSolution(Stokes *s)</div><div>{</div><div> PetscInt row, start, end, i, j;</div><div> PetscScalar val;</div><div> Vec y0, y1;</div><div><br></div><div> PetscFunctionBeginUser;</div><div> /* velocity part */</div><div> PetscCall(VecGetSubVector(s->y, s->isg[0], &y0));</div><div> PetscCall(VecGetOwnershipRange(y0, &start, &end));</div><div> for (row = start; row < end; row++) {</div><div> PetscCall(StokesGetPosition(s, row, &i, &j));</div><div> if (row < s->nx * s->ny) {</div><div> val = StokesExactVelocityX(j * s->hy + s->hy / 2);</div><div> } else {</div><div> val = 0;</div><div> }</div><div> PetscCall(VecSetValue(y0, row, val, INSERT_VALUES));</div><div> }</div><div> PetscCall(VecRestoreSubVector(s->y, s->isg[0], &y0));</div><div><br></div><div>VecSetValues() is called on y0 but VecAssemblyBegin/End is never called hence that vector is not assembled and will result in the error reported. </div><div>All the vectors in this example on which VecSetValue[s] is called need to have VecAssemblyBegin/End. </div><div><br></div><div> Barry</div><div><br></div><div>For sequential vectors where entries set with VecSetValues() are inserted directly in the array it is possible to "get away" with not calling VecAssemblyBegin/End but</div><div>on should never do that so even in the sequential case this example is technically buggy</div><div><br></div><div><br></div></div><blockquote type="cite"><div>On Jul 15, 2023, at 7:08 PM, Rochan Upadhyay <rochan.upadhyay@gmail.com> wrote:</div><br class="Apple-interchange-newline"><div><div dir="ltr"><div>Dear Petsc developers,</div><div>
<div>Is the use of VecNest format supported for field-split type preconditioner or for use with MatNest formats ?</div><div>A quick change to example snes/tutorials/ex70.c seems to suggest that it is not possible.</div><div><br></div><div>For example the changes :</div><div>Vec subx[2], subb[2]; (declared in Stokes struct)<br>In routine StokesSetupVectors(Stokes *s) :<br>MatCreateVecs(s->subA[0], &s->subx[0], &s->subb[0]);<br>MatCreateVecs(s->subA[3], &s->subx[1], &s->subb[1]);<br>//<br>IS is_rows[2];<br>MatNestGetISs(s->A, is_rows, NULL);<br>//<br>VecCreateNest(PETSC_COMM_WORLD, 2, is_rows, s->subx, &s->x);<br>VecCreateNest(PETSC_COMM_WORLD, 2, is_rows, s->subb, &s->b);<br>VecAssemblyBegin(s->x);<br>VecAssemblyEnd(s->x);<br>VecAssemblyBegin(s->b);<br>VecAssemblyEnd(s->b);<br>....<br>The serial case runs fine but the parallel with np = 2 fails with the error message copied below. Could you find out what I am doing wrong ? In th error message, which vector(s) needs assembling ?<br>I would normally not use VecNest but many finite element packages that allow hooking up with petsc seem to have their block vectors of Petsc VecType VecNest</div><div>Any help/insights would be appreciated.</div><div>Regards,</div><div>Rochan</div><div><br></div><div>[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------<br>[0]PETSC ERROR: Object is in wrong state<br>[0]PETSC ERROR: Not for unassembled vector, did you call VecAssemblyBegin()/VecAssemblyEnd()?<br>[0]PETSC ERROR: See <a href="https://petsc.org/release/faq/">https://petsc.org/release/faq/</a> for trouble shooting.<br>[0]PETSC ERROR: Petsc Release Version 3.19.3, Jun 29, 2023<br>[0]PETSC ERROR: ./ex70 on a named localhost.localdomain by user Sat Jul 15 17:51:01 2023<br>[0]PETSC ERROR: Configure options --with-clean=yes --with-debugging=no --with-mpi=1 --with-mpi-dir=/home/user/software_installs/mpich-4.0.2/ --prefix=/home/user/software_installs/petsc-3.19.3/ --PETSC_ARCH=linux-gnu-opt --PETSC_DIR=/home/user/software_downloads/petsc-3.19.3/ --with-cxx-dialect=C++11 --with-hypre=1 --download-hypre=yes --with-suitesparse=1 --download-suitesparse=yes --with-metis=1 --download-metis=yes --with-shared-libraries=1 --with-zlib=1 --download-zlib=no --with-superlu_dist=1 --download-superlu_dist=yes --with-mumps=1 --download-mumps=yes --with-parmetis=1 --download-parmetis=yes --with-hdf5=yes --download-hdf5=1 --with-chaco=1 --download-chaco=yes --with-scalapack=1 --download-scalapack=yes --with-p4est=1 --download-p4est=yes --download-fblaslapack=1 --with-strumpack=1 --download-strumpack=yes<br>[0]PETSC ERROR: #1 VecCopy() at /home/user/software_downloads/petsc-3.19.3/src/vec/vec/interface/vector.c:1684<br>[0]PETSC ERROR: #2 VecCopy_Nest() at /home/user/software_downloads/petsc-3.19.3/src/vec/vec/impls/nest/vecnest.c:58<br>[0]PETSC ERROR: #3 VecCopy() at /home/user/software_downloads/petsc-3.19.3/src/vec/vec/interface/vector.c:1723<br>[0]PETSC ERROR: #4 KSPInitialResidual() at /home/user/software_downloads/petsc-3.19.3/src/ksp/ksp/interface/itres.c:60<br>[0]PETSC ERROR: #5 KSPSolve_GMRES() at /home/user/software_downloads/petsc-3.19.3/src/ksp/ksp/impls/gmres/gmres.c:226<br>[0]PETSC ERROR: #6 KSPSolve_Private() at /home/user/software_downloads/petsc-3.19.3/src/ksp/ksp/interface/itfunc.c:898<br>[0]PETSC ERROR: #7 KSPSolve() at /home/user/software_downloads/petsc-3.19.3/src/ksp/ksp/interface/itfunc.c:1070<br>[0]PETSC ERROR: #8 main() at /home/user/Projects/Misc/Petsc/SNES/ex70.c:684<br>[0]PETSC ERROR: PETSc Option Table entries:<br>[0]PETSC ERROR: -ksp_monitor (source: command line)<br>[0]PETSC ERROR: -nx 4 (source: command line)<br>[0]PETSC ERROR: -ny 8 (source: command line)<br>[0]PETSC ERROR: ----------------End of Error Message -------send entire error message to petsc-maint@mcs.anl.gov----------<br>Abort(73) on node 0 (rank 0 in comm 16): application called MPI_Abort(MPI_COMM_SELF, 73) - process 0<br></div>
</div></div>
</div></blockquote></div><br></body></html>