text: // Copyright (C) 2011-2012 by the Bem++ Authors
  //
  // Permission is hereby granted, free of charge, to any person obtaining a copy
  // of this software and associated documentation files (the "Software"), to deal
  // in the Software without restriction, including without limitation the rights
  // to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
  // copies of the Software, and to permit persons to whom the Software is
  // furnished to do so, subject to the following conditions:
  //
  // The above copyright notice and this permission notice shall be included in
  // all copies or substantial portions of the Software.
  //
  // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
  // IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
  // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
  // AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
  // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
  // OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
  // THE SOFTWARE.
  
  #ifndef fiber_element_pair_topology_hpp
  #define fiber_element_pair_topology_hpp
  
  #include "../common/common.hpp"
  
  #include "../common/armadillo_fwd.hpp"
  #include <cassert>
  #include <iostream>
  #include <boost/tuple/tuple_comparison.hpp>
  
  namespace Fiber
  {
  
  struct ElementPairTopology
  {
      ElementPairTopology() :
          type(Disjoint),
          testVertexCount(0), trialVertexCount(0),
          testSharedVertex0(-1), testSharedVertex1(-1),
          trialSharedVertex0(-1), trialSharedVertex1(-1)
      {}
  
      enum Type { Disjoint, SharedVertex, SharedEdge, Coincident };
      Type type;
      unsigned char testVertexCount;
      unsigned char trialVertexCount;
      signed char testSharedVertex0;
      signed char testSharedVertex1;
      signed char trialSharedVertex0;
      signed char trialSharedVertex1;
  
      bool operator<(const ElementPairTopology& other) const {
          using boost::tuples::make_tuple;
          return make_tuple(type, testVertexCount, trialVertexCount,
                            testSharedVertex0, testSharedVertex1,
                            trialSharedVertex0, trialSharedVertex1) <
                  make_tuple(other.type, other.testVertexCount, other.trialVertexCount,
                             other.testSharedVertex0, other.testSharedVertex1,
                             other.trialSharedVertex0, other.trialSharedVertex1);
      }
  
      bool operator==(const ElementPairTopology& other) const {
          return type == other.type &&
                  testVertexCount == other.testVertexCount &&
                  trialVertexCount == other.trialVertexCount &&
                  testSharedVertex0 == other.testSharedVertex0 &&
                  testSharedVertex1 == other.testSharedVertex1 &&
                  trialSharedVertex0 == other.trialSharedVertex0 &&
                  trialSharedVertex1 == other.trialSharedVertex1;
      }
  
      bool operator!=(const ElementPairTopology& other) const {
          return !operator==(other);
      }
  
      friend std::ostream&
      operator<< (std::ostream& dest, const ElementPairTopology& obj)
      {
          dest << obj.type << " "
               << (int)obj.testVertexCount << " "
               << (int)obj.trialVertexCount << " "
               << (int)obj.testSharedVertex0 << " "
               << (int)obj.testSharedVertex1 << " "
               << (int)obj.trialSharedVertex0 << " "
               << (int)obj.trialSharedVertex1;
          return dest;
      }
  };
  
  inline ElementPairTopology determineElementPairTopologyIn3D(
          const arma::Col<int>& testElementCornerIndices,
          const arma::Col<int>& trialElementCornerIndices)
  {
      ElementPairTopology topology;
  
      // Determine number of element corners
  #ifndef NDEBUG
      const int MIN_VERTEX_COUNT = 3;
  #endif
      const int MAX_VERTEX_COUNT = 4;
      topology.testVertexCount = testElementCornerIndices.n_rows;
      assert(MIN_VERTEX_COUNT <= topology.testVertexCount &&
             topology.testVertexCount <= MAX_VERTEX_COUNT);
      topology.trialVertexCount = trialElementCornerIndices.n_rows;
      assert(MIN_VERTEX_COUNT <= topology.trialVertexCount &&
             topology.trialVertexCount <= MAX_VERTEX_COUNT);
  
      // How many vertices coincide?
      int testSharedVertices[MAX_VERTEX_COUNT];
      int trialSharedVertices[MAX_VERTEX_COUNT];
      int hits = 0;
  
      for (int trialV = 0; trialV < topology.trialVertexCount; ++trialV)
          for (int testV = 0; testV < topology.testVertexCount; ++testV)
              if (testElementCornerIndices(testV) ==
                      trialElementCornerIndices(trialV))
              {
                  testSharedVertices[hits] = testV;
                  trialSharedVertices[hits] = trialV;
                  ++hits;
                  break;
              }
  
      if (hits == 0)
      {
          topology.type = ElementPairTopology::Disjoint;
      }
      else if (hits == 1)
      {
          topology.type = ElementPairTopology::SharedVertex;
          topology.testSharedVertex0 = testSharedVertices[0];
          topology.trialSharedVertex0 = trialSharedVertices[0];
      }
      else if (hits == 2) // && elementDim == 2)
      {
          topology.type = ElementPairTopology::SharedEdge;
          topology.testSharedVertex0 = testSharedVertices[0];
          topology.testSharedVertex1 = testSharedVertices[1];
          topology.trialSharedVertex0 = trialSharedVertices[0];
          topology.trialSharedVertex1 = trialSharedVertices[1];
      }
      // coincident
      else if (hits == topology.testVertexCount &&
               hits == topology.trialVertexCount)
      {
          // Note: we don't handle the case where elements are different, but
          // their vertices coincide
          for (int i = 0; i < hits; ++i)
              assert(testSharedVertices[i] == trialSharedVertices[i]);
          topology.type = ElementPairTopology::Coincident;
      }
      else
          throw std::runtime_error(
                  "Standard3DIntegrationManager::"
                  "selectTestKernelTrialQuadratureRules() :"
                  "Invalid element configuration");
  
      return topology;
  }
  
  } // namespace Fiber
  
  #endif
  