I have no objections to a single mat type and would be happy to refactor it that way. But before I do that, I think some discussion on the subject is useful.
AIJ matrices in PETSc are structured similarly to the LMVM matrices. They have a base type, MATAIJ. They have major subtypes like MATSEQAIJ and MATMPIAIJ. And then they have a third layer, minor subtypes, like MATSEQAIJMKL. Furthermore, the minor subtypes also have their own special constructors like MatCreateSeqAIJMKL. Most importantly though, there is no such thing as a MatAIJType or MatSeqAIJType. All of the subtypes, on both second and third layers, are directly defined as MatType. LMVM is doing the same thing. In fact it’s explicitly modeled after the AIJ code structure. It is indeed missing the *additional* ability to get and set subtypes, but that’s not a replacement to the existing code. It’s an additional functionality. In that sense, creating a MatLMVMType and preventing subtypes from being created directly in the top level Mat interface would go against what has been previously done for other Mat types in PETSc. Jed is right, by the way, about packing optimization specific things here. However, the modularity benefits are significant. Furthermore, good and bad Broyden methods are actually used almost exclusively for general nonlinear system solutions instead of optimization, because they cannot produce symmetric approximations necessary for Hessians. Implementing these things as Mat objects enables both SNES and TAO to make use of them (although currently only TAO uses them). They are not true matrices. They encode nonlinear algorithm information in matrix-like formats. That’s unusual, I agree, but hardcoding them to the algorithms themselves results in massive amounts of duplicate code whenever anyone wants to use these types of things in other places besides the specific algorithm they’ve been written into. The two basic Broyden methods can indeed bet combined into one object easily, but that’s not as true for other subtypes. Mathematical formulations differ significantly between BFGS, DFP, symmetric Broyden and SR1 methods. They can be combined on paper, because they’re all symmetric Broyden-class of updates. However, doing so causes BFGS, DFP and SR1 to have inflated memory footprints and additional algebra operations that they don’t actually need. Eliminating those in a single object requires either a lot of ugly conditionals/switches, or playing with function pointers. Doing the latter basically gets us 99% of the way to separating them into different objects though, which is what the current implementation does. —— Alp Dener Argonne National Laboratory http://www.mcs.anl.gov/person/alp-dener On Sep 11, 2018, at 10:33 AM, Munson, Todd <[email protected]<mailto:[email protected]>> wrote: I would be happy with one matrix type and removing the convenience functions from the code. I am surprised that the dark lord did not mention it when he approved the pull request... Matt: the difference in the subtypes is that they are using different formulas that result in different approximations. BFGS is for symmetric problems and is a rank-2 update that maintains positive definiteness, SR1 is for symmetric problems and is a rank-1 update that can be indefinite. The nonsymmetric equivalents (BadBroyden and GoodBroyden) are also different and bad Broyden is typically better than good Broyden. Performance on a particular problem varies greatly depending on both the specific problem and the approximation type selected. Recall to get performance, you can either make the computational kernels run faster, change the algorithm, or both. Changing the matrix approximations changes the algorithm. Todd. On Sep 11, 2018, at 10:12 AM, Matthew Knepley <[email protected]<mailto:[email protected]>> wrote: On Tue, Sep 11, 2018 at 9:47 AM Dener, Alp <[email protected]<mailto:[email protected]>> wrote: 1 base type and 8 subtypes. If there’s a better convention/structure to do this in PETSc, I’d be happy to get a refactoring done ASAP this week so that it’s cleaner in the release. I guess what Lisandro is asking is, what backend differences are there to make so many type necessary. Looking at the few I recognize, it seems like they encode the sum of rank 1 matrices (BFGS, SR1). Why would we have separate types for this two things? Thanks, Matt — Alp Dener Argonne National Laboratory https://mcs.anl.gov/person/alp-Dener On Sep 11, 2018, at 5:10 AM, Lisandro Dalcin <[email protected]> wrote: So now we have 9 new, top-level, public matrix types for LMVM... Really? #define MATLMVM "lmvm" #define MATLMVMDFP "lmvmdfp" #define MATLMVMBFGS "lmvmbfgs" #define MATLMVMSR1 "lmvmsr1" #define MATLMVMBRDN "lmvmbrdn" #define MATLMVMBADBRDN "lmvmbadbrdn" #define MATLMVMSYMBRDN "lmvmsymbrdn" #define MATLMVMSYMBADBRDN "lmvmsymbadbrdn" #define MATLMVMDIAGBRDN "lmvmdiagbrdn" -- Lisandro Dalcin ============ Research Scientist Computer, Electrical and Mathematical Sciences & Engineering (CEMSE) Extreme Computing Research Center (ECRC) King Abdullah University of Science and Technology (KAUST) http://ecrc.kaust.edu.sa/ 4700 King Abdullah University of Science and Technology al-Khawarizmi Bldg (Bldg 1), Office # 0109 Thuwal 23955-6900, Kingdom of Saudi Arabia http://www.kaust.edu.sa Office Phone: +966 12 808-0459 -- 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 https://www.cse.buffalo.edu/~knepley/
