Optimization for large XML read

Grid/serialisation/XmlIO.h

@@ -184,11 +184,29 @@ namespace Grid
       std::cout << GridLogWarning << "XML: cannot open node '" << s << "'";
       std::cout << std::endl;
     } else {
-      for(unsigned int i = 0; node_.child("elem"); )
+      // Consumed children are *removed* rather than renamed to "elem-done".
+      //
+      // pugixml's node_.child("elem") scans the child list from the start.
+      // Leaving consumed nodes in place under a new name meant iteration i
+      // had to skip i dead siblings, and the lookup happens three times per
+      // element (loop condition, read, rename) -- so reading an N-element
+      // vector cost O(N^2) node visits.  Measured on this reader: 11 s for
+      // 32k elements, 278 s for 131k, i.e. hours for the ~500k elements a
+      // 512x512 complex matrix produces.  Removing the node instead keeps
+      // the next "elem" as the first child, making each lookup O(1).
+      //
+      // Counting first also lets the vector be sized once instead of grown
+      // one element at a time.
+      size_t n = 0;
+      for (pugi::xml_node e = node_.child("elem"); e; e = e.next_sibling("elem"))
       {
-        output.resize(i + 1);
-        read("elem", output[i++]);
-        node_.child("elem").set_name("elem-done");
+        n++;
+      }
+      output.resize(n);
+      for (size_t i = 0; i < n; ++i)
+      {
+        read("elem", output[i]);
+        node_.remove_child(node_.child("elem"));
       }
       pop();
     }
@@ -225,7 +243,11 @@ namespace Grid
       {
         NumElements--;
         read("elem", x);
-        node_.child("elem").set_name("elem-done");
+        // Remove rather than rename, for the same O(N^2) reason as the
+        // std::vector reader above.  Only "elem" nodes are removed; the
+        // tensor-level renaming below is left alone because push(sName)
+        // relies on it to find the next unvisited tensor block.
+        node_.remove_child(node_.child("elem"));
         int i;
         for( i = Rank - 1 ; i != -1 && ++MyIndex[i] == dim[i] ; i-- )
           MyIndex[i] = 0;