[petsc-dev] [Fortran] MatLoad() error: MatLoad is not supported for type: mpiaij!

[Jed Brown]

On Fri, 3 Sep 2010 11:57:59 -0500, Barry Smith <bsmith at mcs.anl.gov> wrote:
> 
>   Jed,
> 
>     Any sorting algorithm that gives good performance is what we need. The particular algorithm class doesn't matter.

The issue is the possibility/probability of hitting a bad-case for
quicksort.  Robust quicksorts usually do something different if they
detect that pivot choice is going poorly (like do a lot of work to
choose a good pivot or switch to heapsort).  Heapsort usually costs
slightly more than quicksort, but I think most users would be willing to
pay that in exchange for guaranteed O(n log n) performance (as far as I
know, sorting is never an especially performance-sensitive kernel in
PETSc, but it's a big deal if 1 second turns into hours occasionally
(e.g. when the mesh is changed)).

>     Thanks for looking at this and improving the results dramatically. BTW: if maintaining the PETScStack and that stuff is slowing down the debug version a lot because of the recursion you could take the PetscFunctionBegin()/Return stuff from the private function that recursives.

I did, but for a different reason, see the third bullet:

changeset:   16901:fb7c48e8b485
user:        Jed Brown <jed at 59A2.org>
date:        Fri Sep 03 18:44:43 2010 +0200
files:       src/sys/utils/sorti.c
description:
Speed up PetscSortInt

* Use median of first,middle,last for quicksort pivot selection.

* Use less intrusive partitioning algorithm (many fewer swaps for
  almost-sorted array).

* Don't use PetscFunctionBegin/Return for the private recursive
  function.  It costs some, but more importantly, the depth of recursion
  is data-dependent, so it can blow PETSc's stack, thus making
  PETSc-provided stack traces useless.  The function could never fail
  anyway, so this does not lose any error-handling capability.