From 2c186a49722bf462b9c8878f650281ca95e2fe88 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Jacob=20St=C3=B8ren?= <jacob.storen@ceetronSolutions.com>
Date: Thu, 11 Jun 2015 13:15:36 +0200
Subject: [PATCH] Added method to calculate characteristic element size in Fem
 parts

Needed for well path pipe radius #312
---
 .../GeoMech/GeoMechDataModel/RigFemPart.cpp   | 45 ++++++++++++++++++-
 .../GeoMech/GeoMechDataModel/RigFemPart.h     |  3 ++
 .../GeoMechDataModel/RigFemPartCollection.cpp | 15 +++++++
 .../GeoMechDataModel/RigFemPartCollection.h   |  1 +
 4 files changed, 63 insertions(+), 1 deletion(-)

diff --git a/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPart.cpp b/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPart.cpp
index 2995ed0146..699472a25e 100644
--- a/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPart.cpp
+++ b/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPart.cpp
@@ -25,7 +25,7 @@
 /// 
 //--------------------------------------------------------------------------------------------------
 RigFemPart::RigFemPart()
-    :m_elementPartId(-1)
+    :m_elementPartId(-1), m_characteristicElementSize(std::numeric_limits<float>::infinity())
 {
 
 }
@@ -278,3 +278,46 @@ cvf::Vec3f RigFemPart::faceNormal(int elmIdx, int faceIdx)
     return cvf::Vec3f::ZERO;
 }
 
+//--------------------------------------------------------------------------------------------------
+/// 
+//--------------------------------------------------------------------------------------------------
+float RigFemPart::characteristicElementSize()
+{
+    if (m_characteristicElementSize != std::numeric_limits<float>::infinity()) return m_characteristicElementSize;
+
+    // take 100 elements 
+    float elmIdxJump = elementCount() / 100.0f;
+    int elmIdxIncrement = elmIdxJump < 1 ? 1: (int) elmIdxJump;
+    int elmsToAverageCount = 0;
+    float sumMaxEdgeLength = 0;
+    for (int elmIdx = 0; elmIdx < elementCount(); elmIdx += elmIdxIncrement)
+    {
+        RigElementType eType = this->elementType(elmIdx);
+
+        if (eType == HEX8)
+        {
+            const int* elmentConn = this->connectivities(elmIdx);
+            cvf::Vec3f nodePos0 = this->nodes().coordinates[elmentConn[0]];
+            cvf::Vec3f nodePos1 = this->nodes().coordinates[elmentConn[1]];
+            cvf::Vec3f nodePos3 = this->nodes().coordinates[elmentConn[3]];
+            cvf::Vec3f nodePos4 = this->nodes().coordinates[elmentConn[4]];
+
+            float l1 = (nodePos1-nodePos0).length();
+            float l3 = (nodePos3-nodePos0).length();
+            float l4 = (nodePos4-nodePos0).length();
+
+            float maxLength = l1 > l3 ? l1: l3;
+            maxLength = maxLength > l4 ? maxLength: l4;
+
+            sumMaxEdgeLength += maxLength;
+            ++elmsToAverageCount;
+        }
+    }
+
+    CVF_ASSERT(elmsToAverageCount);
+
+    m_characteristicElementSize = sumMaxEdgeLength/elmsToAverageCount;
+
+    return m_characteristicElementSize;
+}
+
diff --git a/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPart.h b/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPart.h
index 1d9e7615a1..4b1355dd28 100644
--- a/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPart.h
+++ b/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPart.h
@@ -71,6 +71,7 @@ public:
     int                         neighborFace(int elementIndex, int faceIndex) const
                                 { return m_elmNeighbors[elementIndex].faceInNeighborElm[faceIndex]; }
 
+    float                       characteristicElementSize();
     const std::vector<int>&     possibleGridCornerElements() const { return m_possibleGridCornerElements; }
 
     cvf::Vec3f                  faceNormal(int elmentIndex, int faceIndex);
@@ -97,4 +98,6 @@ private:
     std::vector<  Neighbors > m_elmNeighbors;
     std::vector<int> m_possibleGridCornerElements;
 
+    float                       m_characteristicElementSize;
+
 };
diff --git a/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartCollection.cpp b/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartCollection.cpp
index b847cb4964..6e994d4e4d 100644
--- a/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartCollection.cpp
+++ b/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartCollection.cpp
@@ -83,3 +83,18 @@ size_t RigFemPartCollection::totalElementCount() const
 
     return elementCount;
 }
+
+//--------------------------------------------------------------------------------------------------
+/// 
+//--------------------------------------------------------------------------------------------------
+float RigFemPartCollection::characteristicElementSize()
+{
+    if (partCount())
+    {
+        return part(0)->characteristicElementSize();
+    }
+    else
+    {
+        return 0;
+    }
+}
diff --git a/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartCollection.h b/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartCollection.h
index 841b909bd9..0127143512 100644
--- a/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartCollection.h
+++ b/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartCollection.h
@@ -36,6 +36,7 @@ public:
     int partCount() const;
 
     size_t totalElementCount() const;
+    float characteristicElementSize();
 
 
 private: