[petsc-users] [EXTERNAL] Re: Example code of linear solve of a block matrix system in Fortran

[Satish Balay]

You can use https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatSetValuesBlocked.html

MatSetValues() will also work - but MatSetValuesBlocked() is more efficient wrt BAIJ

Satish

On Wed, 6 May 2020, Thomas S. Chyczewski wrote:

> Thanks Satish. I have seen that page and can create a block matrix. It's not clear to me how to fill it and use it in a Fortran code.
> 
> -----Original Message-----
> From: Satish Balay <balay at mcs.anl.gov> 
> Sent: Wednesday, May 6, 2020 3:43 PM
> To: Thomas S. Chyczewski <tsc109 at arl.psu.edu>
> Cc: petsc-users at mcs.anl.gov
> Subject: [EXTERNAL] Re: [petsc-users] Example code of linear solve of a block matrix system in Fortran
> 
> What you are looking for is:
> 
> https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateBAIJ.html
> 
> indoe is optimization for AIJ type - where you might not have a simple block size like 3 that you have.
> 
> BAIJ will perform better than AIJ/with-inode
> 
> To verify indoe usage - you can run the code with '-log_info' option - and do a 'grep inode'
> 
> Satish
> 
> On Wed, 6 May 2020, Thomas S. Chyczewski wrote:
> 
> > 
> > All,
> > 
> > I'm relatively new to PETSc and have relied pretty heavily on the example codes included in the distribution to figure out the finer points of using the PETSc library that I couldn't deduce from the manual. One thing I can't figure out is how to solve block matrix systems and I couldn't find an example in Fortran. I'm writing a 2D incompressible CFD solver so I have a 3x3 block Imax*Jmax system. The closest I've come to finding an example is ex19.c in the snes directory, but that is in c and for the nonlinear solver.
> > 
> > I have been able to run PETSc but unwrapping the block matrix into a monolithic system. But the manual says "Block matrices represent an important class of problems in numerical linear algebra and offer the possibility of far more efficient iterative solvers than just treating the entire matrix as black box." However, in the FAQs I saw a comment that PETSc scans the AIJ matrices for rows that have the same column layout and can deduce if it's a block system and use the more efficient solvers. I also saw in the archives for this email list a thread where it seems workaround for building fields in a Fortran code is discussed ("Back to struct in Fortran to represent field with dof > 1"), so I'm beginning to suspect building a block system in Fortran might not be straight forward.
> > 
> > All that being said, my questions:
> > 
> > Is there a significant advantage to building the block system as opposed to the analogous monolithic system if PETSc can figure out that it's a block system? Can you confirm that PETSc does figure this out?
> > If there is an advantage to loading the matrix as a block matrix, is there an example Fortran code that builds and solves a linear block system?
> > 
> > Thanks,
> > Tom C
> > 
> > 
>