How to remove empty row-column pairs from an asymmetric square matrix

Question

What is the best way to remove empty rows and columns of an asymmetric square matrix if and only if they are both empty for a given index. That is, it should remove row AND column i if and only if both row i AND column i are all zeros. The result will be a square matrix. This is useful for reducing the size of directed adjacency matrices for graphs or Markov models among other things.

This question is similar to this question, but requires a solution that works when we can't satisfy the condition: "when the i-th row is zeros, then i-th column is also zeros." I believe lacking this constraint makes the DeleteCases and Select-based methods unusable as they are offered. But there may be a way to adapt them to this more general case.

Currently my best method is the following:

First find the indices for which the respective matrix row AND column are zeros.

m = SomeSquareMatrix;
n = Length[m];
ids=Table[If[Total[m[[i]]]==0&&Total[m[[All,i]]]==0,{i},##&[]],{i,n}];

Note that Total works great for probabilities and edge weights that are always positive, but could be replaced by

ids=Table[If[m[[i]]==ConstantArray[0, n]&&Flatten[m[[All,i]]]==ConstantArray[0, n],{i},##&[]],{i,n}];

The ##&[] means that if the condition is false, return nothing. One you've identified the empty row-column combos, delete them:

m=Transpose[Delete[Transpose[Delete[m,ids]],ids]];

The last line is basically the same method used for this question and similar to the DeleteCases approach, so it might actually be the best way to go. However, I suspect there is a superior method that combines identifying and removing the empty pairs of rows and columns.

Answer 1

Generating a random, sparse square matrix.

m = n = 200000;
k = 1000000;
A = SparseArray[
   Rule[
    Transpose[{RandomInteger[{1, m}, {k}], 
      RandomInteger[{1, n}, {k}]}],
    RandomReal[{-1, 1}, {k}]
    ], {m, n}, 0.];

Detecting all rows and columns with nonzero entries and extracting the corresponding square matrix.

With[{B = SparseArray[A]},
   With[{ilist = Sort[DeleteDuplicates[Flatten[B["NonzeroPositions"]]]]},
    A[[ilist, ilist]]
    ]
   ]; // AbsoluteTiming

{0.076962, Null}