On Mon, 21 Apr 2008, tribur at vision.ee.ethz.ch wrote: > Dear all, > > Sorry for switching from Schur to Hypre and back, but I'm trying two > approaches at the same time to find the optimal solution for our > convection-diffusion/Stokes problems: a) solving the global stiffness matrix > directly and in parallel using Petsc and a suitable preconditioner (???) and > b) applying first non-overlapping domain decomposition and than solving the > Schur complement system. > > Being concerned with b in the moment, I managed to set up and solve the global > Schur system using MATDENSE. The solving works well with, e.g., gmres+jacobi, > but the assembling of the global Schur matrix takes too long.
Hmm - with dense - if you have some other efficient way of assembling the matrix - you can specify this directly to MatCreateMPIDense() - [or use MatGetArray() - and set the values directly into this array] > Therefore, I'm > trying to use the matrix in unassembled form using MatShell. Not very > successfully, however: > > 1) When I use KSPGMRES, I got the error > [1]PETSC ERROR: MatMult() line 1632 in src/mat/interface/matrix.c > [1]PETSC ERROR: PCApplyBAorAB() line 584 in src/ksp/pc/interface/precon.c > [1]PETSC ERROR: GMREScycle() line 159 in src/ksp/ksp/impls/gmres/gmres.c > [1]PETSC ERROR: KSPSolve_GMRES() line 241 in src/ksp/ksp/impls/gmres/gmres.c > [1]PETSC ERROR: KSPSolve() line 379 in src/ksp/ksp/interface/itfunc.c Which version of PETSc is this? I can't place the line numbers correctly with latest petsc-2.3.3. [Can you send the complete error trace?] > 2) Using KSPBICG, it iterates without error message, but the result is wrong > (norm of residual 1.42768 instead of something like 1.0e-10), although my > Mat-functions PETSC_SchurMatMult and PETSC_SchurMatMultTranspose seem to be > correct. I tested the latter comparing the vectors y1 and y2 computed by, > e.g., MatMult(S,x,y1) and PETS_SchurMatMult(S,x,y2). Norm(y1-y2) was < e-15 > for both functions. Not sure what the problem could be. Can you confirm that the code is valgrind clean? It could explain the issue 1 aswell. With mpich2 you can do the following on linux: mpiexec -np 2 valgrind --tool=memcheck ./executable Satish > > > Could you please have a look at my code snippet below? > > Thank you very much! > Kathrin > > > > PS: My Code: > > Vec gtot, x; > ... > Mat Stot; IS is; > ISCreateGeneral(PETSC_COMM_SELF, NPb, &uBId_global[0], &is); > localData ctx; > ctx.NPb = NPb; //size of local Schur system S > ctx.Sloc = &S[0]; > ctx.is = is; > MatCreateShell(PETSC_COMM_WORLD,m,n,NPb_tot,NPb_tot,&ctx,&Stot); > MatShellSetOperation(Stot,MATOP_MULT,(void(*)(void)) PETSC_SchurMatMult); > MatShellSetOperation(Stot,MATOP_MULT_TRANSPOSE,(void(*)(void))PETSC_SchurMatMultTranspose); > KSP ksp; > KSPCreate(PETSC_COMM_WORLD,&ksp); > PC prec; > KSPSetOperators(ksp,Stot,Stot,DIFFERENT_NONZERO_PATTERN); > KSPGetPC(ksp,&prec); > PCSetType(prec, PCNONE); > KSPSetType(ksp, KSPBICG); > KSPSetTolerances(ksp, 1.e-10, 1.e-50,PETSC_DEFAULT,PETSC_DEFAULT); > KSPSolve(ksp,gtot,x); > ... > > >
