<div dir="ltr">On Mon, May 20, 2013 at 11:13 AM, Jed Brown <span dir="ltr"><<a href="mailto:jedbrown@mcs.anl.gov" target="_blank">jedbrown@mcs.anl.gov</a>></span> wrote:<br><div class="gmail_extra"><div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="im">"Kaus, Boris" <<a href="mailto:kaus@uni-mainz.de">kaus@uni-mainz.de</a>> writes:<br>
<br>
> Hi,<br>
><br>
> Example 42 (src/ksp/ksp/examples/tutorials) with recursive fieldsplit used to work in petsc 3.2 and petsc 3.3p2 with the patch suggested below by Barry.<br>
<br>
</div>Damn,<br>
<br>
<a href="https://bitbucket.org/petsc/petsc/commits/ce780c64fa296d753b6ed81263eb2f3164d5f63f" target="_blank">https://bitbucket.org/petsc/petsc/commits/ce780c64fa296d753b6ed81263eb2f3164d5f63f</a><br>
<br>
The two reference commits are:<br>
<br>
<a href="https://bitbucket.org/petsc/petsc/commits/5847835f01a83a65f10e398c0972b3c7a1e1c5f4" target="_blank">https://bitbucket.org/petsc/petsc/commits/5847835f01a83a65f10e398c0972b3c7a1e1c5f4</a><br>
<a href="https://bitbucket.org/petsc/petsc/commits/4442dace165024532baeb2e4a78f8791bee44482" target="_blank">https://bitbucket.org/petsc/petsc/commits/4442dace165024532baeb2e4a78f8791bee44482</a><br>
<br>
We cannot just force the new block size into the input vector because<br>
that breaks other tests.  We can fix it with disgusting hacks:<br>
<br>
  x_bs = x->map->bs;<br>
<div class="im">  x->map->bs = jac->bs;<br>
</div>  VecStrideGatherAll(x,jac->x,INSERT_VALUES);<br>
  x->map->bs = x_bs;<br>
<br>
or we can replace our use of VecStrideGatherAll and VecStrideScatterAll<br>
with something else (or add an argument to the VecStride functions so<br>
that we can pass the block size---traditional users would pass<br>
PETSC_DECIDE).<br>
<br>
I don't consider making PetscLayout intentionally mutable is acceptable<br>
because it is shared many different ways.</blockquote><div><br></div><div style>I think we have to allow an argument to the strided operations, but the arg</div><div style>should be a factor of the bs of the Vec I think.</div>
<div style><br></div><div style>   Matt</div><div style> <br></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">
> mpiexec -np 4 ./ex42 \<br>
> -stokes_ksp_type gcr \<br>
> -stokes_ksp_rtol 1.0e-6 \<br>
> -stokes_pc_type fieldsplit \<br>
> -stokes_pc_fieldsplit_type SCHUR \<br>
> -stokes_pc_fieldsplit_schur_factorization_type UPPER \<br>
> -stokes_fieldsplit_u_ksp_type fgmres \<br>
> -stokes_fieldsplit_u_ksp_rtol 1e-3 \<br>
> -stokes_fieldsplit_u_pc_type fieldsplit \<br>
> -stokes_fieldsplit_u_fieldsplit_ksp_type preonly \<br>
> -stokes_fieldsplit_u_fieldsplit_pc_type ml \<br>
> -stokes_fieldsplit_u_pc_fieldsplit_block_size 3 \<br>
> -stokes_fieldsplit_u_pc_fieldsplit_type ADDITIVE \<br>
> -stokes_fieldsplit_p_ksp_type preonly \<br>
> -stokes_fieldsplit_p_pc_type jacobi \<br>
> -stokes_ksp_monitor_blocks \<br>
> -mx 16 \<br>
> -model 3<br>
><br>
> It no longer works in  petsc 3.4.0. Is this something that can be fixed and potentially added to 3.4.1?<br>
> Our production code relies on similar functionality.<br>
><br>
> thanks a lot!<br>
><br>
> Boris<br>
><br>
><br>
><br>
> On Jun 28, 2012, at 8:55 PM, Barry Smith <<a href="mailto:bsmith@mcs.anl.gov">bsmith@mcs.anl.gov</a>> wrote:<br>
><br>
>>  Anton,<br>
>><br>
>>   This came about because we are now being much more pedantic about the blocksizes of PETSc objects and not allowing them to be causally changed when they shouldn't be.<br>
>><br>
>>   You can resolve this problem by editing the file src/ksp/pc/impls/fieldsplit/fieldsplit.c locate the function<br>
>><br>
>> #undef __FUNCT__<br>
>> #define __FUNCT__ "PCApply_FieldSplit"<br>
>> static PetscErrorCode PCApply_FieldSplit(PC pc,Vec x,Vec y)<br>
>> {<br>
>>  PC_FieldSplit     *jac = (PC_FieldSplit*)pc->data;<br>
>>  PetscErrorCode    ierr;<br>
>>  PC_FieldSplitLink ilink = jac->head;<br>
>>  PetscInt          cnt,bs;<br>
>><br>
>>  PetscFunctionBegin;<br>
>><br>
>> and add the two lines right here<br>
>><br>
>>  x->map->bs = jac->bs;<br>
>>  y->map->bs = jac->bs;<br>
>><br>
>><br>
>> then run make cmake in that directory.<br>
>><br>
>> To resolve this permanently we will need to figure out how to insure those inner vectors are created with the correct block size. Are you willing to share your code with <a href="mailto:petsc-maint@mcs.anl.gov">petsc-maint@mcs.anl.gov</a> so that we can reproduce the problem and fix it properly for the long run? (The problem is in PETSc not in your code).<br>

>><br>
>>   Barry<br>
>><br>
>><br>
>><br>
>> On Jun 28, 2012, at 10:44 AM, Anton Popov wrote:<br>
>><br>
>>> Dear petsc team,<br>
>>><br>
>>> I'm trying to use fieldsplit preconditioner for the velocity block in the Stokes system which is also preconditioned by<br>
>>> fieldsplit (kind of recursive).<br>
>>><br>
>>> Running example 42 from src/ksp/ksp/examples/tutorials with petsc-3.2, as follows:<br>
>>><br>
>>> mpiexec -np 4 ./ex42 \<br>
>>> -stokes_ksp_type gcr \<br>
>>> -stokes_ksp_rtol 1.0e-6 \<br>
>>> -stokes_pc_type fieldsplit \<br>
>>> -stokes_pc_fieldsplit_type SCHUR \<br>
>>> -stokes_pc_fieldsplit_schur_factorization_type UPPER \<br>
>>> -stokes_fieldsplit_u_ksp_type fgmres \<br>
>>> -stokes_fieldsplit_u_ksp_rtol 1e-3 \<br>
>>> -stokes_fieldsplit_u_pc_type fieldsplit \<br>
>>> -stokes_fieldsplit_u_fieldsplit_ksp_type preonly \<br>
>>> -stokes_fieldsplit_u_fieldsplit_pc_type ml \<br>
>>> -stokes_fieldsplit_u_pc_fieldsplit_block_size 3 \<br>
>>> -stokes_fieldsplit_u_pc_fieldsplit_type ADDITIVE \<br>
>>> -stokes_fieldsplit_p_ksp_type preonly \<br>
>>> -stokes_fieldsplit_p_pc_type jacobi \<br>
>>> -stokes_ksp_monitor_blocks \<br>
>>> -mx 16 \<br>
>>> -model 3<br>
>>><br>
>>> gives nicely looking output.<br>
>>><br>
>>> But! Repeating the same exercise with petsc-3.3, like this (actually, there is only one difference: factorization -> fact):<br>
>>><br>
>>> mpiexec -np 4 ./ex42 \<br>
>>> -stokes_ksp_type gcr \<br>
>>> -stokes_ksp_rtol 1.0e-6 \<br>
>>> -stokes_pc_type fieldsplit \<br>
>>> -stokes_pc_fieldsplit_type SCHUR \<br>
>>> -stokes-pc_fieldsplit_schur_fact_type UPPER \<br>
>>> -stokes_fieldsplit_u_ksp_type fgmres \<br>
>>> -stokes_fieldsplit_u_ksp_rtol 1e-3 \<br>
>>> -stokes_fieldsplit_u_pc_type fieldsplit \<br>
>>> -stokes_fieldsplit_u_fieldsplit_ksp_type preonly \<br>
>>> -stokes_fieldsplit_u_fieldsplit_pc_type ml \<br>
>>> -stokes_fieldsplit_u_pc_fieldsplit_block_size 3 \<br>
>>> -stokes_fieldsplit_u_pc_fieldsplit_type ADDITIVE \<br>
>>> -stokes_fieldsplit_p_ksp_type preonly \<br>
>>> -stokes_fieldsplit_p_pc_type jacobi \<br>
>>> -stokes_ksp_monitor_blocks \<br>
>>> -mx 16 \<br>
>>> -model 3<br>
>>><br>
>>> curses me hardly by claiming:<br>
>>><br>
>>> [0]PETSC ERROR: --------------------- Error Message ------------------------------------<br>
>>> [0]PETSC ERROR: Object is in wrong state!<br>
>>> [0]PETSC ERROR: Blocksize of x vector 1 does not match fieldsplit blocksize 3!<br>
>>> [0]PETSC ERROR: ------------------------------------------------------------------------<br>
>>> [0]PETSC ERROR: Petsc Release Version 3.3.0, Patch 0, Tue Jun  5 14:20:42 CDT 2012<br>
>>> [0]PETSC ERROR: See docs/changes/index.html for recent updates.<br>
>>> [0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.<br>
>>> [0]PETSC ERROR: See docs/index.html for manual pages.<br>
>>> [0]PETSC ERROR: ------------------------------------------------------------------------<br>
>>> [0]PETSC ERROR: ./ex42 on a int32-deb named <a href="http://mac11-005.geo.uni-mainz.de" target="_blank">mac11-005.geo.uni-mainz.de</a> by anton Thu Jun 28 17:06:53 2012<br>
>>> [0]PETSC ERROR: Libraries linked from /Users/anton/LIB/petsc-3.3-p0/int32-debug/lib<br>
>>> [0]PETSC ERROR: Configure run at Tue Jun 12 15:32:21 2012<br>
>>> [0]PETSC ERROR: Configure options PETSC_DIR=/Users/anton/LIB/petsc-3.3-p0 PETSC_ARCH=int32-debug --download-f-blas-lapack=1 --with-debugging=1 --COPTFLAGS="-g -O0" --FOPTFLAGS="-g -O0" --CXXOPTFLAGS="-g -O0" --with-c++-support=1 --with-fortran=1 --with-fortran-kernels=1 --with-large-file-io=1 --with-mpi-compilers=1 --with-cc=mpicc --with-cxx=mpicxx --with-fc=mpif90 --download-ml=1 --download-hypre=1 --download-blacs=1 --download-scalapack=1 --download-metis=1 --download-parmetis=1 --download-mumps=1 --download-superlu_dist=1<br>

>>> [0]PETSC ERROR: ------------------------------------------------------------------------<br>
>>> [0]PETSC ERROR: PCApply_FieldSplit() line 726 in /Users/anton/LIB/petsc-3.3-p0/src/ksp/pc/impls/fieldsplit/fieldsplit.c<br>
>>> [0]PETSC ERROR: PCApply() line 384 in /Users/anton/LIB/petsc-3.3-p0/src/ksp/pc/interface/precon.c<br>
>>> [0]PETSC ERROR: KSPFGMRESCycle() line 169 in /Users/anton/LIB/petsc-3.3-p0/src/ksp/ksp/impls/gmres/fgmres/fgmres.c<br>
>>> [0]PETSC ERROR: KSPSolve_FGMRES() line 294 in /Users/anton/LIB/petsc-3.3-p0/src/ksp/ksp/impls/gmres/fgmres/fgmres.c<br>
>>> [0]PETSC ERROR: KSPSolve() line 446 in /Users/anton/LIB/petsc-3.3-p0/src/ksp/ksp/interface/itfunc.c<br>
>>> [0]PETSC ERROR: PCApply_FieldSplit_Schur() line 693 in /Users/anton/LIB/petsc-3.3-p0/src/ksp/pc/impls/fieldsplit/fieldsplit.c<br>
>>> [0]PETSC ERROR: PCApply() line 384 in /Users/anton/LIB/petsc-3.3-p0/src/ksp/pc/interface/precon.c<br>
>>> [0]PETSC ERROR: KSPSolve_GCR_cycle() line 47 in /Users/anton/LIB/petsc-3.3-p0/src/ksp/ksp/impls/gcr/gcr.c<br>
>>> [0]PETSC ERROR: KSPSolve_GCR() line 117 in /Users/anton/LIB/petsc-3.3-p0/src/ksp/ksp/impls/gcr/gcr.c<br>
>>> [0]PETSC ERROR: KSPSolve() line 446 in /Users/anton/LIB/petsc-3.3-p0/src/ksp/ksp/interface/itfunc.c<br>
>>> [0]PETSC ERROR: solve_stokes_3d_coupled() line 2045 in src/ksp/ksp/examples/tutorials/ex42.c<br>
>>> [0]PETSC ERROR: main() line 2106 in src/ksp/ksp/examples/tutorials/ex42.c<br>
>>> application called MPI_Abort(MPI_COMM_WORLD, 1) - process 0<br>
>>><br>
>>> Similar error appeared in our code after upgrading to petsc-3.3, and we're using similar functionality and options as I posted above.<br>
>>><br>
>>> Please explain this issue. An advice how to get rid of the error is also appreciated.<br>
>>><br>
>>> Thanks a lot<br>
>>><br>
>>> Anton<br>
>><br>
><br>
><br>
><br>
> -----------------------------------------------------------------------------<br>
> Boris J.P. Kaus<br>
><br>
> Institute of Geosciences,<br>
> Geocycles Research Center &<br>
> Center for Computational Sciences.<br>
> University of Mainz, Mainz, Germany<br>
> Office:       00-285<br>
> Tel:          +49.6131.392.4527<br>
><br>
> <a href="http://www.geophysik.uni-mainz.de" target="_blank">http://www.geophysik.uni-mainz.de</a><br>
> -----------------------------------------------------------------------------<br>
</div></div></blockquote></div><br><br clear="all"><div><br></div>-- <br>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>