[petsc-users] MATBDIAG in > 3.1

[Parker, Andrew (UK Filton)]

I have tried to use the following, this was before I had both you're
emails (Barry, Matthew) ;-)

1) First I set this up

      PetscInt size = 7;
      MatCreate(MPI_COMM_SELF,&LHS);
      MatSetSizes(LHS,size,size,size,size);
      MatSetType(LHS,MATSEQDENSE);  //**************NOTE THIS LINE
      MatAssemblyBegin(LHS,MAT_FINAL_ASSEMBLY);
      MatAssemblyEnd(LHS,MAT_FINAL_ASSEMBLY);

      VecCreate(PETSC_COMM_SELF,&RHS);
      VecSetSizes(RHS,PETSC_DECIDE, size);
      VecSetFromOptions(RHS);
      VecAssemblyBegin(RHS);
      VecAssemblyEnd(RHS);

      VecCreate(PETSC_COMM_SELF,&Updates);
      VecSetSizes(RHS,PETSC_DECIDE, size);
      VecSetFromOptions(Updates);
      VecAssemblyBegin(Updates);
      VecAssemblyEnd(Updates);
      
      KSPCreate(PETSC_COMM_SELF,&Solver);
      KSPSetOperators(Solver, LHS, LHS, SAME_NONZERO_PATTERN);
      KSPSetType(Solver, KSPPREONLY);  //**************NOTE THIS LINE
      KSPGetPC(Solver, &prec);
      PCSetType(prec, PCLU);  //**************NOTE THIS LINE
      KSPSetFromOptions(Solver);
      KSPSetUp(Solver);

For each cell I then do:

2) I copy out the appropriate data from global structures into LHS and
RHS 

3) KSPSolve(Solver, RHS, Updates);

4) Stick the local Updates back into my global data structure

Could I have your thoughts on whether the above is ok?  In terms of
efficiency and the correct thing to do?  Would I be correct in assuming
that Barry's method and mine are identical with the exception that I've
got some copies going on to set up these small 7x7 and 1x7 matrices and
vectors, but Barry's method does this inline by doing a preonly, lu
solve DIRECT solve on each 7x7??

Many thanks,
Andy

-----Original Message-----
From: petsc-users-bounces at mcs.anl.gov
[mailto:petsc-users-bounces at mcs.anl.gov] On Behalf Of Barry Smith
Sent: 01 December 2011 14:09
To: PETSc users list
Subject: Re: [petsc-users] MATBDIAG in > 3.1


                    *** WARNING ***

  This message has originated outside your organisation,
  either from an external partner or the Global Internet. 
      Keep this in mind if you answer this message.
 

On Dec 1, 2011, at 7:05 AM, Matthew Knepley wrote:

> You can get this effect using PCBJACOBI, setting the number of local
blocks, and using -sub_pc_type lu -sub_ksp_type preonly.

   Do not do this. It will give bad performance since it explicitly
pulls out each block with MatGetSubMatrices() and solves each small
matrix as a separate PETSc matrix.

    Instead use the BAIJ matrix format with a block size of 7 and then
use a PCPBJACOBI (note the P in the name for point block, each point
block is that little 7 by 7 block that it solves efficiently using a
tuned kernel) with -ksp_type preonly.  This should be very fast and does
not do any any neccessary scatters or communication.

    Barry

> 
>    Matt
> 
> On Thu, Dec 1, 2011 at 4:28 AM, Parker, Andrew (UK Filton)
<Andrew.Parker2 at baesystems.com> wrote:
> Hi,
> 
>  
> 
> I've been using PETSc for a conventional full parallel implicit
solver, currently in MPIBAIJ format.  However, I'd like to understand
the performance verse memory and speed trade-offs of just having a block
diagonal system.  In this case I do not scatter to the off diagonal
blocks.  Each block is 7x7 in size and fully dense.  In my view this
should be a point implicit solution of the system.  I am looking for the
best way to solve this in PETSc, clearly because of the lack of coupling
between neighbouring cells by only having a block diagonal it's actually
a 7x7 dense linear solve, but to test my theory at present I'd prefer to
just not scatter to off diagonals, and use:
> 
>  
> 
> MATBDIAG with my current ksp solver and take it from there.  BDIAG
seems to have been removed as off 3.0:
> 
>  
> 
> http://www.mcs.anl.gov/petsc/documentation/changes/300.html
> 
>  
> 
> Can anybody help me with this, your thoughts on the above and what I
should use to solve this reduced system would be appreciated??
> 
> 
> Cheers,
> 
> Andy
> 
> 
> ********************************************************************
> This email and any attachments are confidential to the intended 
> recipient and may also be privileged. If you are not the intended 
> recipient please delete it from your system and notify the sender.
> You should not copy it or use it for any purpose nor disclose or 
> distribute its contents to any other person.
> ********************************************************************
> 
> 
> 
> 
> --
> What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which
their experiments lead.
> -- Norbert Wiener