Eric Chamberland <[email protected]> writes: > On 09/11/2014 12:52 PM, Matthew Knepley wrote: >> On Thu, Sep 11, 2014 at 11:47 AM, Eric Chamberland >> >> MatNest is absolutely the worst thing in the PETSc interface. You should > > Sounds so strange to me... :-)
MatNest is really supposed to be a backend optimization, not a public interface. Perhaps irresponsibly, we expose a few MatNest*() functions because those interfaces are convenient for a few legacy applications that made poor design/workflow decisions. There are a few other functions in PETSc with similar status, e.g., MatCreateMPIAIJWithSplitArrays(), but surely Matt would rate MatNest as the most popular and harmful such interface. > ... because we decided to use MatNest to create sub-matrices for linear > vs quadratic parts of a velocity field > (http://onlinelibrary.wiley.com/doi/10.1002/nla.757/abstract). We > choosed MatNest to minimize the memory used, and "overloaded" the > MatSetValues calls to be able to do the assembly of u-u blocks in all 4 > sub-matrices (doing global to local conversions with strides only). Better to work with local indices all along. The result is cleaner and faster. I understand that this option may not be available due to prior design choices. > This prevent us to do 4 different loops over all the elements (in other > words, it allows the assembly part of our code to not be aware that we > are doing the assembly into a matnest)! We just have to do the > numbering of the dofs correctly before everything. > >> never ever ever ever >> be calling MatNest directly. You should be assembling into one matrix >> from views. Then MatNest >> can be used for optimization in the background. > > By "views" you mean doing MatGetLocalSubMatrix and > MatRestoreLocalSubMatrix like in > http://www.mcs.anl.gov/petsc/petsc-current/src/snes/examples/tutorials/ex28.c.html > > ? > > but... If I do a velocity(u)-pressure(p) field problem, how do I fix a > block size of 3 for the u-u part and 1 for other parts without using a > MatNest (with sub-matrices with good options)? The isrow and iscol arguments should have a block size of 3. Then the Mat returned by MatGetLocalSubMatrix will inherit this block size such that MatSetValuesBlockedLocal behaves as desired.
pgpD5mtvZ_B2M.pgp
Description: PGP signature
