Using a combination of Jacobi and KSP scaling, you almost never need a diagonal matrix in solves. However, if you really do, I would just use MatDiagonScale():
http://www-unix.mcs.anl.gov/petsc/petsc-as/snapshots/petsc-current/docs/manualpages/Mat/MatDiagonalScale.html#MatDiagonalScale and if you really need to, you can wrap that operation in a MatShell. Matt On Mon, Mar 31, 2008 at 4:11 PM, Shi Jin <jinzishuai at yahoo.com> wrote: > > Hi, > > I have a related question to the one above: what is the best way to create a > purely diagonal matrix? > I am using > ierr = MatCreateMPIAIJ(comm,localM,localN,PETSC_DETERMINE,PETSC_DETERMINE,1, > PETSC_NULL, 0, PETSC_NULL, &lumpedMP);CHKERRQ(ierr); > But I would think there might be a even better (faster/less memory) way? > > I thought about the BDiag Matrix. I tried > ierr=MatCreateMPIBDiag(comm,localM,globalSize,globalSize,PETSC_NULL,1, > PETSC_NULL,PETSC_NULL,&lumpedM);CHKERRQ(ierr); > but not sure if it is done correctly. I am a little suspicious about those > pointers I set as NULL. > > If the above two are both correct, is there a better one? > > Thank you very much. > -- > Shi Jin, PhD > > > ________________________________ > Special deal for Yahoo! users & friends - No Cost. Get a month of > Blockbuster Total Access now -- 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
