diff --git a/ApplicationLibCode/Commands/WellPathCommands/RicCreateWellTargetsPickEventHandler.cpp b/ApplicationLibCode/Commands/WellPathCommands/RicCreateWellTargetsPickEventHandler.cpp index 1c60cf8a69..6925f3e03e 100644 --- a/ApplicationLibCode/Commands/WellPathCommands/RicCreateWellTargetsPickEventHandler.cpp +++ b/ApplicationLibCode/Commands/WellPathCommands/RicCreateWellTargetsPickEventHandler.cpp @@ -241,7 +241,7 @@ cvf::Vec3d RicCreateWellTargetsPickEventHandler::findHexElementIntersection( gsl RigFemPart* femPart = geoMechView->femParts()->part( femPartIndex ); RigElementType elType = femPart->elementType( elementIndex ); - if ( elType == HEX8 || elType == HEX8P ) + if ( RigFemTypes::is8NodeElement( elType ) ) { cellIndex = elementIndex; const RigFemPartGrid* femGrid = femPart->getOrCreateStructGrid(); diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.cpp index 0e44fd634d..29258b93ba 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.cpp @@ -92,6 +92,14 @@ int RigFemPart::elementCount() const return static_cast( m_elementId.size() ); } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +int RigFemPart::nodeCount() const +{ + return static_cast( m_nodes.nodeIds.size() ); +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -442,7 +450,7 @@ float RigFemPart::characteristicElementSize() const { if ( m_characteristicElementSize != std::numeric_limits::infinity() ) return m_characteristicElementSize; - std::vector elementPriority = { HEX8P, HEX8 }; + std::vector elementPriority = { RigElementType::HEX8P, RigElementType::HEX8 }; for ( auto elmType : elementPriority ) { @@ -604,7 +612,7 @@ size_t RigFemPart::resultValueIdxFromResultPosType( RigFemResultPosEnum resultPo bool RigFemPart::isHexahedron( size_t elementIdx ) const { RigElementType elType = elementType( elementIdx ); - return elType == HEX8 || elType == HEX8P; + return RigFemTypes::is8NodeElement( elType ); } //-------------------------------------------------------------------------------------------------- diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.h b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.h index 06aa53b17c..11c83aa732 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.h +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.h @@ -58,6 +58,7 @@ public: void appendElement( RigElementType elmType, int elementId, const int* connectivities ); int elementCount() const; + int nodeCount() const; int allConnectivitiesCount() const; int elmId( size_t elementIdx ) const; diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartGrid.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartGrid.cpp index 52cdcef5ff..5e03d21da6 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartGrid.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartGrid.cpp @@ -200,7 +200,7 @@ void RigFemPartGrid::generateStructGridData() RigElementType elementType = m_femPart->elementType( elmIdx ); size_t i, j, k; bool validIndex = ijkFromCellIndex( elmIdx, &i, &j, &k ); - if ( elementType == HEX8P && validIndex ) + if ( ( elementType == RigElementType::HEX8P ) && validIndex ) { if ( i < min.x() ) min.x() = i; if ( j < min.y() ) min.y() = j; @@ -254,7 +254,7 @@ cvf::Vec3i RigFemPartGrid::findMainIJKFaces( int elementIndex ) const // Record three independent main direction vectors for the element, and what face they are created from cvf::Vec3f mainElmDirections[3]; int mainElmDirOriginFaces[3]; - if ( eType == HEX8 || eType == HEX8P ) + if ( RigFemTypes::is8NodeElement( eType ) ) { mainElmDirections[0] = normals[0] - normals[1]; // To get a better "average" direction vector mainElmDirections[1] = normals[2] - normals[3]; diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorEnIpPorBar.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorEnIpPorBar.cpp index c2a643ba4b..7d421f643a 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorEnIpPorBar.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorEnIpPorBar.cpp @@ -92,7 +92,7 @@ RigFemScalarResultFrames* RigFemPartResultCalculatorEnIpPorBar::calculate( int p { RigElementType elmType = femPart->elementType( elmIdx ); - if ( elmType == HEX8P ) + if ( elmType == RigElementType::HEX8P ) { int elmNodeCount = RigFemTypes::elementNodeCount( elmType ); for ( int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx ) diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorGamma.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorGamma.cpp index eb17f6201d..e8f706e52e 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorGamma.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorGamma.cpp @@ -135,7 +135,7 @@ void RigFemPartResultCalculatorGamma::calculateGammaFromFrames( int int elmNodeCount = RigFemTypes::elementNodeCount( femPart->elementType( elmIdx ) ); - if ( elmType == HEX8P ) + if ( elmType == RigElementType::HEX8P ) { for ( int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx ) { diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorInitialPorosity.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorInitialPorosity.cpp index e7feee3f55..bd64955ccc 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorInitialPorosity.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorInitialPorosity.cpp @@ -94,7 +94,7 @@ RigFemScalarResultFrames* RigFemPartResultCalculatorInitialPorosity::calculate( int elmNodeCount = RigFemTypes::elementNodeCount( femPart->elementType( elmIdx ) ); - if ( elmType == HEX8P ) + if ( elmType == RigElementType::HEX8P ) { for ( int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx ) { diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNodalGradients.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNodalGradients.cpp index acce939115..7bf0874e6f 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNodalGradients.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNodalGradients.cpp @@ -124,7 +124,7 @@ RigFemScalarResultFrames* RigFemPartResultCalculatorNodalGradients::calculate( i for ( int elmIdx : elements ) { RigElementType elmType = femPart->elementType( elmIdx ); - if ( elmType == HEX8P ) + if ( elmType == RigElementType::HEX8P ) { // Find the corner coordinates and values for the node std::array hexCorners; diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalSE.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalSE.cpp index de914efa9f..6964643ad7 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalSE.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalSE.cpp @@ -94,7 +94,7 @@ RigFemScalarResultFrames* RigFemPartResultCalculatorNormalSE::calculate( int par int elmNodeCount = RigFemTypes::elementNodeCount( femPart->elementType( elmIdx ) ); - if ( elmType == HEX8P ) + if ( elmType == RigElementType::HEX8P ) { for ( int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx ) { diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalST.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalST.cpp index 20ce822237..010ae24675 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalST.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalST.cpp @@ -102,7 +102,7 @@ RigFemScalarResultFrames* RigFemPartResultCalculatorNormalST::calculate( int par int elmNodeCount = RigFemTypes::elementNodeCount( femPart->elementType( elmIdx ) ); - if ( elmType == HEX8P ) + if ( elmType == RigElementType::HEX8P ) { for ( int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx ) { diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalized.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalized.cpp index 40b830484e..89ad4320b5 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalized.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorNormalized.cpp @@ -123,7 +123,7 @@ RigFemScalarResultFrames* RigFemPartResultCalculatorNormalized::calculate( int p for ( int elmIdx = 0; elmIdx < femPart->elementCount(); ++elmIdx ) { RigElementType elmType = femPart->elementType( elmIdx ); - if ( elmType != HEX8 && elmType != HEX8P ) continue; + if ( !RigFemTypes::is8NodeElement( elmType ) ) continue; bool porRegion = false; for ( int elmLocalNodeIdx = 0; elmLocalNodeIdx < 8; ++elmLocalNodeIdx ) diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorPoreCompressibility.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorPoreCompressibility.cpp index 09a941f5f3..e4cca5e513 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorPoreCompressibility.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorPoreCompressibility.cpp @@ -170,7 +170,7 @@ RigFemScalarResultFrames* RigFemPartResultCalculatorPoreCompressibility::calcula int elmNodeCount = RigFemTypes::elementNodeCount( femPart->elementType( elmIdx ) ); - if ( elmType == HEX8P ) + if ( elmType == RigElementType::HEX8P ) { for ( int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx ) { diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorPorosityPermeability.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorPorosityPermeability.cpp index e97eb56752..56afae59e1 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorPorosityPermeability.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorPorosityPermeability.cpp @@ -165,7 +165,7 @@ RigFemScalarResultFrames* RigFemPartResultCalculatorPorosityPermeability::calcul int elmNodeCount = RigFemTypes::elementNodeCount( femPart->elementType( elmIdx ) ); - if ( elmType == HEX8P ) + if ( elmType == RigElementType::HEX8P ) { for ( int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx ) { diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorShearSE.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorShearSE.cpp index 4a74500d4c..d3d75b9c55 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorShearSE.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorShearSE.cpp @@ -96,7 +96,7 @@ RigFemScalarResultFrames* RigFemPartResultCalculatorShearSE::calculate( int part int elmNodeCount = RigFemTypes::elementNodeCount( femPart->elementType( elmIdx ) ); - if ( elmType == HEX8P ) + if ( elmType == RigElementType::HEX8P ) { for ( int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx ) { diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorStressGradients.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorStressGradients.cpp index 2bfdf03486..6d2f611982 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorStressGradients.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultCalculatorStressGradients.cpp @@ -133,7 +133,7 @@ RigFemScalarResultFrames* RigFemPartResultCalculatorStressGradients::calculate( const int* cornerIndices = femPart->connectivities( elmIdx ); RigElementType elmType = femPart->elementType( elmIdx ); - if ( elmType != HEX8P && elmType != HEX8 ) continue; + if ( !RigFemTypes::is8NodeElement( elmType ) ) continue; // Find the corner coordinates for element std::array hexCorners; diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultsCollection.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultsCollection.cpp index 278078cc6c..f378652ca6 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultsCollection.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPartResultsCollection.cpp @@ -443,6 +443,7 @@ RigFemScalarResultFrames* RigFemPartResultsCollection::findOrLoadScalarResult( i frames = calculateDerivedResult( partIndex, resVarAddr ); if ( frames ) return frames; + // is it available as an imported property if ( resVarAddr.resultPosType == RIG_ELEMENT ) { std::map> elementProperties = @@ -463,11 +464,6 @@ RigFemScalarResultFrames* RigFemPartResultsCollection::findOrLoadScalarResult( i { return frames; } - else - { - // Create a dummy empty result - return m_femPartResults[partIndex]->createScalarResult( resVarAddr ); - } } // We need to read the data as bulk fields, and populate the correct scalar caches @@ -504,6 +500,9 @@ RigFemScalarResultFrames* RigFemPartResultsCollection::findOrLoadScalarResult( i case RIG_NODAL: m_readerInterface->readNodeField( resVarAddr.fieldName, partIndex, stepIndex, frameIndex, &componentDataVectors ); break; + case RIG_ELEMENT: + m_readerInterface->readElementField( resVarAddr.fieldName, partIndex, stepIndex, frameIndex, &componentDataVectors ); + break; case RIG_ELEMENT_NODAL: m_readerInterface->readElementNodeField( resVarAddr.fieldName, partIndex, stepIndex, frameIndex, &componentDataVectors ); break; @@ -573,14 +572,20 @@ std::map> RigFemPartResultsCollection::sca if ( resPos == RIG_NODAL ) { fieldCompNames = m_readerInterface->scalarNodeFieldAndComponentNames(); - if ( fieldCompNames.contains( "U" ) ) fieldCompNames["U"].push_back( "U_LENGTH" ); fieldCompNames[RigFemAddressDefines::porBar()]; - fieldCompNames[FIELD_NAME_COMPACTION]; + + if ( m_readerInterface->populateDerivedResultNames() ) + { + if ( fieldCompNames.contains( "U" ) ) fieldCompNames["U"].push_back( "U_LENGTH" ); + fieldCompNames[FIELD_NAME_COMPACTION]; + } } else if ( resPos == RIG_ELEMENT_NODAL ) { fieldCompNames = m_readerInterface->scalarElementNodeFieldAndComponentNames(); + if ( !m_readerInterface->populateDerivedResultNames() ) return fieldCompNames; + fieldCompNames["SE"].push_back( "SM" ); fieldCompNames["SE"].push_back( "SFI" ); fieldCompNames["SE"].push_back( "DSM" ); @@ -677,6 +682,8 @@ std::map> RigFemPartResultsCollection::sca { fieldCompNames = m_readerInterface->scalarIntegrationPointFieldAndComponentNames(); + if ( !m_readerInterface->populateDerivedResultNames() ) return fieldCompNames; + fieldCompNames["SE"].push_back( "SM" ); fieldCompNames["SE"].push_back( "SFI" ); fieldCompNames["SE"].push_back( "DSM" ); @@ -779,6 +786,8 @@ std::map> RigFemPartResultsCollection::sca } else if ( resPos == RIG_ELEMENT_NODAL_FACE ) { + if ( !m_readerInterface->populateDerivedResultNames() ) return fieldCompNames; + fieldCompNames["Plane"].push_back( "Pinc" ); fieldCompNames["Plane"].push_back( "Pazi" ); @@ -799,9 +808,11 @@ std::map> RigFemPartResultsCollection::sca } else if ( resPos == RIG_ELEMENT ) { + fieldCompNames = m_readerInterface->scalarElementFieldAndComponentNames(); + for ( const std::string& field : m_elementPropertyReader->scalarElementFields() ) { - fieldCompNames[field]; + fieldCompNames[field] = {}; } } else if ( resPos == RIG_WELLPATH_DERIVED ) @@ -883,6 +894,9 @@ std::vector RigFemPartResultsCollection::getResAddrToCompon case RIG_NODAL: fieldAndComponentNames = m_readerInterface->scalarNodeFieldAndComponentNames(); break; + case RIG_ELEMENT: + fieldAndComponentNames = m_readerInterface->scalarElementFieldAndComponentNames(); + break; case RIG_ELEMENT_NODAL: fieldAndComponentNames = m_readerInterface->scalarElementNodeFieldAndComponentNames(); break; diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemTypes.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemTypes.cpp index 5cb8003311..fb5520bff5 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemTypes.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemTypes.cpp @@ -30,7 +30,7 @@ int RigFemTypes::elementNodeCount( RigElementType elmType ) { static int elementTypeCounts[3] = { 8, 8, 4 }; - return elementTypeCounts[elmType]; + return elementTypeCounts[(unsigned int)elmType]; } //-------------------------------------------------------------------------------------------------- @@ -40,7 +40,7 @@ int RigFemTypes::elementFaceCount( RigElementType elmType ) { const static int elementFaceCounts[3] = { 6, 6, 1 }; - return elementFaceCounts[elmType]; + return elementFaceCounts[(unsigned int)elmType]; } //-------------------------------------------------------------------------------------------------- @@ -63,12 +63,12 @@ const int* RigFemTypes::localElmNodeIndicesForFace( RigElementType elmType, int switch ( elmType ) { - case HEX8: - case HEX8P: + case RigElementType::HEX8: + case RigElementType::HEX8P: ( *faceNodeCount ) = 4; return HEX8_Faces[faceIdx]; break; - case CAX4: + case RigElementType::CAX4: ( *faceNodeCount ) = 4; return CAX4_Faces; break; @@ -86,11 +86,11 @@ int RigFemTypes::oppositeFace( RigElementType elmType, int faceIdx ) switch ( elmType ) { - case HEX8: - case HEX8P: + case RigElementType::HEX8: + case RigElementType::HEX8P: return HEX8_OppositeFaces[faceIdx]; break; - case CAX4: + case RigElementType::CAX4: return faceIdx; break; default: @@ -111,11 +111,11 @@ const int* RigFemTypes::localElmNodeToIntegrationPointMapping( RigElementType el switch ( elmType ) { - case HEX8: - case HEX8P: + case RigElementType::HEX8: + case RigElementType::HEX8P: return HEX8_Mapping; break; - case CAX4: + case RigElementType::CAX4: return HEX8_Mapping; // First four is identical to HEX8 break; default: @@ -129,22 +129,22 @@ const int* RigFemTypes::localElmNodeToIntegrationPointMapping( RigElementType el //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -QString RigFemTypes::elementTypeText( RigElementType elemType ) +std::string RigFemTypes::elementTypeText( RigElementType elemType ) { - QString txt = "UNKNOWN_ELM_TYPE"; + std::string txt = "UNKNOWN_ELM_TYPE"; switch ( elemType ) { - case HEX8: + case RigElementType::HEX8: txt = "HEX8"; break; - case HEX8P: + case RigElementType::HEX8P: txt = "HEX8P"; break; - case CAX4: + case RigElementType::CAX4: txt = "CAX4"; break; - case UNKNOWN_ELM_TYPE: + case RigElementType::UNKNOWN_ELM_TYPE: break; default: break; @@ -152,3 +152,45 @@ QString RigFemTypes::elementTypeText( RigElementType elemType ) return txt; } + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::map RigFemTypes::femTypeMap() +{ + std::map typeMap; + typeMap["C3D8R"] = RigElementType::HEX8; + typeMap["C3D8"] = RigElementType::HEX8; + typeMap["C3D8P"] = RigElementType::HEX8P; + typeMap["CAX4"] = RigElementType::CAX4; + typeMap["C3D20RT"] = RigElementType::HEX8; + typeMap["C3D8RT"] = RigElementType::HEX8; + typeMap["C3D8T"] = RigElementType::HEX8; + + return typeMap; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RigElementType RigFemTypes::toRigElementType( const std::string odbTypeName ) +{ + static std::map odbElmTypeToRigElmTypeMap = femTypeMap(); + + std::map::iterator it = odbElmTypeToRigElmTypeMap.find( odbTypeName ); + + if ( it == odbElmTypeToRigElmTypeMap.end() ) + { + return RigElementType::UNKNOWN_ELM_TYPE; + } + + return it->second; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RigFemTypes::is8NodeElement( RigElementType elmType ) +{ + return ( elmType == RigElementType::HEX8 ) || ( elmType == RigElementType::HEX8P ); +} diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemTypes.h b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemTypes.h index 2a053923a0..f7ffde404c 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemTypes.h +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemTypes.h @@ -19,9 +19,10 @@ #pragma once -class QString; +#include +#include -enum RigElementType +enum class RigElementType { HEX8, HEX8P, @@ -32,10 +33,15 @@ enum RigElementType class RigFemTypes { public: - static int elementNodeCount( RigElementType elmType ); - static int elementFaceCount( RigElementType elmType ); - static const int* localElmNodeIndicesForFace( RigElementType elmType, int faceIdx, int* faceNodeCount ); - static int oppositeFace( RigElementType elmType, int faceIdx ); - static const int* localElmNodeToIntegrationPointMapping( RigElementType elmType ); - static QString elementTypeText( RigElementType elemType ); + static int elementNodeCount( RigElementType elmType ); + static int elementFaceCount( RigElementType elmType ); + static const int* localElmNodeIndicesForFace( RigElementType elmType, int faceIdx, int* faceNodeCount ); + static int oppositeFace( RigElementType elmType, int faceIdx ); + static const int* localElmNodeToIntegrationPointMapping( RigElementType elmType ); + static std::string elementTypeText( RigElementType elemType ); + static RigElementType toRigElementType( const std::string odbTypeName ); + static bool is8NodeElement( RigElementType elmType ); + +private: + static std::map femTypeMap(); }; diff --git a/ApplicationLibCode/GeoMech/OdbReader/CMakeLists.txt b/ApplicationLibCode/GeoMech/OdbReader/CMakeLists.txt index f3165f1c96..8ba372c02b 100644 --- a/ApplicationLibCode/GeoMech/OdbReader/CMakeLists.txt +++ b/ApplicationLibCode/GeoMech/OdbReader/CMakeLists.txt @@ -46,6 +46,8 @@ add_library( RifInpReader.cpp RifGeoMechReaderInterface.h RifGeoMechReaderInterface.cpp + RifInpIncludeReader.h + RifInpIncludeReader.cpp OdbSetup.cmake ) diff --git a/ApplicationLibCode/GeoMech/OdbReader/RifGeoMechReaderInterface.h b/ApplicationLibCode/GeoMech/OdbReader/RifGeoMechReaderInterface.h index ed983ecd1e..aa814273da 100644 --- a/ApplicationLibCode/GeoMech/OdbReader/RifGeoMechReaderInterface.h +++ b/ApplicationLibCode/GeoMech/OdbReader/RifGeoMechReaderInterface.h @@ -51,6 +51,7 @@ public: virtual std::vector elementSet( int partIndex, std::string partName, int setIndex ) = 0; virtual std::map> scalarNodeFieldAndComponentNames() = 0; + virtual std::map> scalarElementFieldAndComponentNames() = 0; virtual std::map> scalarElementNodeFieldAndComponentNames() = 0; virtual std::map> scalarIntegrationPointFieldAndComponentNames() = 0; @@ -61,6 +62,11 @@ public: int stepIndex, int frameIndex, std::vector*>* resultValues ) = 0; + virtual void readElementField( const std::string& fieldName, + int partIndex, + int stepIndex, + int frameIndex, + std::vector*>* resultValues ) = 0; virtual void readElementNodeField( const std::string& fieldName, int partIndex, int stepIndex, @@ -72,6 +78,8 @@ public: int frameIndex, std::vector*>* resultValues ) = 0; + virtual bool populateDerivedResultNames() const = 0; + void setTimeStepFilter( const std::vector& fileTimeStepIndices, bool readOnlyLastFrame ); bool isTimeStepIncludedByFilter( int timeStepIndex ) const; int timeStepIndexOnFile( int timeStepIndex ) const; diff --git a/ApplicationLibCode/GeoMech/OdbReader/RifInpIncludeReader.cpp b/ApplicationLibCode/GeoMech/OdbReader/RifInpIncludeReader.cpp new file mode 100644 index 0000000000..cea3a88fa9 --- /dev/null +++ b/ApplicationLibCode/GeoMech/OdbReader/RifInpIncludeReader.cpp @@ -0,0 +1,117 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2023- Equinor ASA +// +// ResInsight is free software: you can redistribute it and/or modify +// it under the terms of the GNU General Public License as published by +// the Free Software Foundation, either version 3 of the License, or +// (at your option) any later version. +// +// ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY +// WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. +// +// See the GNU General Public License at +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// + +#pragma once + +#include "RifInpIncludeReader.h" + +#include "RiaStdStringTools.h" + +#include +#include + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RifInpIncludeReader::RifInpIncludeReader() +{ +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RifInpIncludeReader::~RifInpIncludeReader() +{ + if ( isOpen() ) close(); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RifInpIncludeReader::openFile( const std::string& fileName ) +{ + m_stream.open( fileName ); + return m_stream.good(); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RifInpIncludeReader::isOpen() const +{ + return m_stream.is_open(); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RifInpIncludeReader::close() +{ + m_stream.close(); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RifInpIncludeReader::readData( int columnIndex, const std::map& parts, std::map>& data ) +{ + std::map nameToPartId; + for ( auto p : parts ) + { + nameToPartId[p.second] = p.first; + } + + std::string line; + + while ( true ) + { + std::getline( m_stream, line ); + + if ( m_stream ) + { + if ( line.starts_with( '*' ) ) continue; + if ( line.starts_with( ',' ) ) continue; + + // is the requested column present? + auto columns = RiaStdStringTools::splitString( line, ',' ); + if ( columnIndex >= columns.size() ) continue; + + // split part/set/node/element in first column + auto partNode = RiaStdStringTools::splitString( columns[0], '.' ); + if ( partNode.size() != 3 ) continue; + auto& partName = partNode[0]; + int nodeOrElmIndex = RiaStdStringTools::toInt( partNode[2] ) - 1; + if ( nodeOrElmIndex < 0 ) continue; + + // is it a valid part name? + if ( nameToPartId.count( partName ) == 0 ) continue; + int partId = nameToPartId[partName]; + + // is the index as expected? + if ( nodeOrElmIndex >= (int)data[partId].size() ) continue; + + // is it a valid value? + double value = std::numeric_limits::infinity(); + RiaStdStringTools::toDouble( RiaStdStringTools::trimString( columns[columnIndex] ), value ); + + data[partId][nodeOrElmIndex] = value; + } + + if ( m_stream.eof() ) break; + } +} diff --git a/ApplicationLibCode/GeoMech/OdbReader/RifInpIncludeReader.h b/ApplicationLibCode/GeoMech/OdbReader/RifInpIncludeReader.h new file mode 100644 index 0000000000..a66e762b3a --- /dev/null +++ b/ApplicationLibCode/GeoMech/OdbReader/RifInpIncludeReader.h @@ -0,0 +1,45 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2023- Equinor ASA +// +// ResInsight is free software: you can redistribute it and/or modify +// it under the terms of the GNU General Public License as published by +// the Free Software Foundation, either version 3 of the License, or +// (at your option) any later version. +// +// ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY +// WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. +// +// See the GNU General Public License at +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// + +#pragma once + +#include +#include +#include +#include + +//================================================================================================== +// +//================================================================================================== +class RifInpIncludeReader +{ +public: + RifInpIncludeReader(); + ~RifInpIncludeReader(); + + bool openFile( const std::string& fileName ); + bool isOpen() const; + + void readData( int column, const std::map& parts, std::map>& data ); + +private: + void close(); + +private: + std::ifstream m_stream; +}; diff --git a/ApplicationLibCode/GeoMech/OdbReader/RifInpReader.cpp b/ApplicationLibCode/GeoMech/OdbReader/RifInpReader.cpp index 4908f9e95b..5f02b782c8 100644 --- a/ApplicationLibCode/GeoMech/OdbReader/RifInpReader.cpp +++ b/ApplicationLibCode/GeoMech/OdbReader/RifInpReader.cpp @@ -18,6 +18,8 @@ #include "RifInpReader.h" +#include "RifInpIncludeReader.h" + #include "RigFemPart.h" #include "RigFemPartCollection.h" #include "RigFemTypes.h" @@ -38,6 +40,7 @@ /// //-------------------------------------------------------------------------------------------------- RifInpReader::RifInpReader() + : m_enableIncludes( true ) { } @@ -48,6 +51,14 @@ RifInpReader::~RifInpReader() { } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RifInpReader::enableIncludes( bool enable ) +{ + m_enableIncludes = enable; +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -56,13 +67,28 @@ void RifInpReader::close() m_stream.close(); } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RifInpReader::populateDerivedResultNames() const +{ + return false; +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- bool RifInpReader::openFile( const std::string& fileName, std::string* errorMessage ) { m_stream.open( fileName ); - return m_stream.good(); + bool bOK = m_stream.good(); + + if ( bOK ) + { + m_inputPath = std::filesystem::path( fileName ).parent_path(); + } + + return bOK; } //-------------------------------------------------------------------------------------------------- @@ -86,9 +112,23 @@ bool RifInpReader::readFemParts( RigFemPartCollection* femParts ) std::map>>> elements; std::map>>> elementSets; - read( m_stream, parts, nodes, elements, elementSets ); + auto elementType = read( m_stream, parts, nodes, elements, elementSets, m_stepNames, m_enableIncludes, m_includeEntries ); - RiaLogging::debug( QString( "Read FEM parts: %1 nodes: %2 elements: %3" ).arg( parts.size() ).arg( nodes.size() ).arg( elements.size() ) ); + for ( int i = 0; i < m_includeEntries.size(); i++ ) + { + m_includeEntries[i].fileName = ( m_inputPath / m_includeEntries[i].fileName ).string(); + } + + RiaLogging::debug( QString( "Read FEM parts: %1, steps: %2, element type: %3" ) + .arg( parts.size() ) + .arg( m_stepNames.size() ) + .arg( QString::fromStdString( RigFemTypes::elementTypeText( elementType ) ) ) ); + + if ( !RigFemTypes::is8NodeElement( elementType ) ) + { + RiaLogging::error( QString( "Unsupported element type." ) ); + return false; + } caf::ProgressInfo modelProgress( parts.size() * 2, "Reading Inp Parts" ); @@ -135,9 +175,7 @@ bool RifInpReader::readFemParts( RigFemPartCollection* femParts ) elementIdToIdxMap[elmId] = elmIdx; - // TODO: handle more types? - RigElementType elmType = RigElementType::HEX8; - int nodeCount = RigFemTypes::elementNodeCount( elmType ); + int nodeCount = RigFemTypes::elementNodeCount( elementType ); indexBasedConnectivities.resize( nodeCount ); for ( int lnIdx = 0; lnIdx < nodeCount; ++lnIdx ) @@ -145,12 +183,12 @@ bool RifInpReader::readFemParts( RigFemPartCollection* femParts ) indexBasedConnectivities[lnIdx] = nodeIdToIdxMap[nodesInElement[lnIdx]]; } - femPart->appendElement( elmType, elmId, indexBasedConnectivities.data() ); + femPart->appendElement( elementType, elmId, indexBasedConnectivities.data() ); } // read element sets - auto elementSetsForPart = elementSets[partId]; - for ( auto [setName, elementSet] : elementSetsForPart ) + auto& elementSetsForPart = elementSets[partId]; + for ( auto& [setName, elementSet] : elementSetsForPart ) { femPart->addElementSet( setName, elementSet ); } @@ -161,22 +199,33 @@ bool RifInpReader::readFemParts( RigFemPartCollection* femParts ) modelProgress.incrementProgress(); } + readScalarData( femParts, parts, m_includeEntries ); + return true; } //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -void RifInpReader::read( std::istream& stream, - std::map& parts, - std::map>>& nodes, - std::map>>>& elements, - std::map>>>& elementSets ) +RigElementType RifInpReader::read( std::istream& stream, + std::map& parts, + std::map>>& nodes, + std::map>>>& elements, + std::map>>>& elementSets, + std::vector& stepNames, + bool enableIncludes, + std::vector& includeEntries ) { std::string line; + std::string prevline; + + RigElementType elementType = RigElementType::UNKNOWN_ELM_TYPE; int partId = 0; std::string partName; + int stepId = -1; + std::string stepName; + int timeSteps = 0; while ( true ) { @@ -208,6 +257,8 @@ void RifInpReader::read( std::istream& // "Element" section. else if ( uppercasedLine.starts_with( "*ELEMENT," ) ) { + auto nodeType = parseLabel( line, "type" ); + elementType = RigFemTypes::toRigElementType( nodeType ); skipComments( stream ); elements[partId] = readElements( stream ); } @@ -219,12 +270,70 @@ void RifInpReader::read( std::istream& auto elementSet = isGenerateSet ? readElementSetGenerate( stream ) : readElementSet( stream ); elementSets[partId].push_back( { setName, elementSet } ); } + else if ( uppercasedLine.starts_with( "*STEP" ) ) + { + stepName = parseLabel( line, "name" ); + stepNames.push_back( stepName ); + stepId = stepNames.size() - 1; + } + else if ( uppercasedLine.starts_with( "*END STEP" ) ) + { + stepId = -1; + stepName = ""; + } + else if ( enableIncludes && uppercasedLine.starts_with( "*INCLUDE" ) ) + { + auto filename = parseLabel( line, "input" ); + RigFemResultPosEnum resultType = RigFemResultPosEnum::RIG_ELEMENT; + std::string propertyName( "" ); + int columnIndex = 1; + if ( prevline.starts_with( "*BOUNDARY" ) ) + { + propertyName = "POR"; + resultType = RigFemResultPosEnum::RIG_NODAL; + columnIndex = 3; + } + else if ( prevline.starts_with( "*TEMPERATURE" ) ) + { + propertyName = "TEMP"; + resultType = RigFemResultPosEnum::RIG_NODAL; + } + else if ( prevline.starts_with( "*INITIAL" ) ) + { + auto label = parseLabel( prevline, "type" ); + if ( label == "RATIO" ) propertyName = "RATIO"; + resultType = RigFemResultPosEnum::RIG_NODAL; + } + if ( propertyName.empty() ) + { + std::string uppercasedFilename = RiaStdStringTools::toUpper( filename ); + + if ( uppercasedFilename.find( "DENSITY" ) != std::string::npos ) + { + propertyName = "DENSITY"; + } + else if ( uppercasedFilename.find( "ELASTICS" ) != std::string::npos ) + { + includeEntries.push_back( RifInpIncludeEntry( "MODULUS", RigFemResultPosEnum::RIG_ELEMENT, stepId, filename, 1 ) ); + + propertyName = "RATIO"; + columnIndex = 2; + } + } + if ( !propertyName.empty() ) + { + includeEntries.push_back( RifInpIncludeEntry( propertyName, resultType, stepId, filename, columnIndex ) ); + } + } + prevline = uppercasedLine; continue; } if ( stream.eof() ) break; } + + return elementType; } //-------------------------------------------------------------------------------------------------- @@ -431,11 +540,7 @@ void RifInpReader::skipComments( std::istream& stream ) //-------------------------------------------------------------------------------------------------- std::vector RifInpReader::allStepNames() const { - // TODO: read step names from file. - std::vector stepNames; - stepNames.push_back( "step 1" ); - - return stepNames; + return m_stepNames; } //-------------------------------------------------------------------------------------------------- @@ -443,7 +548,7 @@ std::vector RifInpReader::allStepNames() const //-------------------------------------------------------------------------------------------------- std::vector RifInpReader::filteredStepNames() const { - // TODO: add filtering. + // no filter supported return RifInpReader::allStepNames(); } @@ -452,9 +557,8 @@ std::vector RifInpReader::filteredStepNames() const //-------------------------------------------------------------------------------------------------- std::vector RifInpReader::frameTimes( int stepIndex ) const { - // TODO: get from file? - std::vector frameValues; - frameValues.push_back( 1.0 ); + // only one frame from INP file + std::vector frameValues( { 1.0 } ); return frameValues; } @@ -463,8 +567,6 @@ std::vector RifInpReader::frameTimes( int stepIndex ) const //-------------------------------------------------------------------------------------------------- int RifInpReader::frameCount( int stepIndex ) const { - if ( shouldReadOnlyLastFrame() ) return 1; - return frameTimes( stepIndex ).size(); } @@ -473,9 +575,7 @@ int RifInpReader::frameCount( int stepIndex ) const //-------------------------------------------------------------------------------------------------- std::vector RifInpReader::elementSetNames( int partIndex, std::string partName ) { - // TODO: not implemented yet if ( partIndex >= m_partElementSetNames.size() ) return {}; - return m_partElementSetNames.at( partIndex ); } @@ -489,12 +589,34 @@ std::vector RifInpReader::elementSet( int partIndex, std::string partNam return elementIndexes; } -// //-------------------------------------------------------------------------------------------------- -// /// -// //-------------------------------------------------------------------------------------------------- +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- std::map> RifInpReader::scalarNodeFieldAndComponentNames() { - return {}; + std::map> retVal; + + for ( auto& entry : m_propertyPartDataNodes ) + { + retVal[entry.first] = {}; + } + + return retVal; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::map> RifInpReader::scalarElementFieldAndComponentNames() +{ + std::map> retVal; + + for ( auto& entry : m_propertyPartDataElements ) + { + retVal[entry.first] = {}; + } + + return retVal; } //-------------------------------------------------------------------------------------------------- @@ -521,6 +643,41 @@ void RifInpReader::readDisplacements( int partIndex, int stepIndex, int frameInd CVF_ASSERT( displacements ); } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RifInpReader::readField( RigFemResultPosEnum resultType, + const std::string& fieldName, + int partIndex, + int stepIndex, + std::vector*>* resultValues ) +{ + CVF_ASSERT( resultValues ); + + auto dataMap = propertyDataMap( resultType ); + if ( dataMap == nullptr ) return; + + if ( dataMap->count( fieldName ) == 0 ) return; + + // is there only a static result? Use it for all steps. + if ( ( *dataMap )[fieldName].count( stepIndex ) == 0 ) + { + if ( ( *dataMap )[fieldName].count( -1 ) == 0 ) return; + stepIndex = -1; + } + if ( ( *dataMap )[fieldName][stepIndex].count( partIndex ) == 0 ) return; + + auto dataSize = ( *dataMap )[fieldName][stepIndex][partIndex].size(); + + ( *resultValues )[0]->resize( dataSize ); + + std::vector* singleComponentValues = ( *resultValues )[0]; + for ( size_t i = 0; i < dataSize; i++ ) + { + ( *singleComponentValues )[i] = (float)( *dataMap )[fieldName][stepIndex][partIndex][i]; + } +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -530,7 +687,19 @@ void RifInpReader::readNodeField( const std::string& fieldName, int frameIndex, std::vector*>* resultValues ) { - CVF_ASSERT( resultValues ); + readField( RigFemResultPosEnum::RIG_NODAL, fieldName, partIndex, stepIndex, resultValues ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RifInpReader::readElementField( const std::string& fieldName, + int partIndex, + int stepIndex, + int frameIndex, + std::vector*>* resultValues ) +{ + readField( RigFemResultPosEnum::RIG_ELEMENT, fieldName, partIndex, stepIndex, resultValues ); } //-------------------------------------------------------------------------------------------------- @@ -542,7 +711,7 @@ void RifInpReader::readElementNodeField( const std::string& field int frameIndex, std::vector*>* resultValues ) { - CVF_ASSERT( resultValues ); + readField( RigFemResultPosEnum::RIG_ELEMENT_NODAL, fieldName, partIndex, stepIndex, resultValues ); } //-------------------------------------------------------------------------------------------------- @@ -554,5 +723,61 @@ void RifInpReader::readIntegrationPointField( const std::string& int frameIndex, std::vector*>* resultValues ) { - CVF_ASSERT( resultValues ); + readField( RigFemResultPosEnum::RIG_INTEGRATION_POINT, fieldName, partIndex, stepIndex, resultValues ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RifInpReader::readScalarData( RigFemPartCollection* femParts, + std::map& parts, + std::vector& includeEntries ) +{ + for ( auto& entry : includeEntries ) + { + auto map = propertyDataMap( entry.resultType ); + if ( map == nullptr ) continue; + if ( map->count( entry.propertyName ) == 0 ) + { + ( *map )[entry.propertyName] = {}; + } + + int stepId = entry.stepId; + + if ( ( *map )[entry.propertyName].count( stepId ) == 0 ) + { + ( *map )[entry.propertyName][stepId] = {}; + } + + for ( int partId = 0; partId < femParts->partCount(); partId++ ) + { + ( *map )[entry.propertyName][stepId][partId] = {}; + size_t dataSize = 0; + if ( entry.resultType == RigFemResultPosEnum::RIG_NODAL ) + { + dataSize = femParts->part( partId )->nodeCount(); + } + if ( entry.resultType == RigFemResultPosEnum::RIG_ELEMENT ) + { + dataSize = femParts->part( partId )->elementCount(); + } + ( *map )[entry.propertyName][stepId][partId].resize( dataSize, 0.0 ); + } + RifInpIncludeReader reader; + if ( !reader.openFile( entry.fileName ) ) continue; + reader.readData( entry.columnIndex, parts, ( *map )[entry.propertyName][stepId] ); + } + + return; +} + +//-------------------------------------------------------------------------------------------------- +/// Map keys: result name / time step / part id +//-------------------------------------------------------------------------------------------------- +std::map>>>* RifInpReader::propertyDataMap( RigFemResultPosEnum resultType ) +{ + if ( resultType == RigFemResultPosEnum::RIG_ELEMENT ) return &m_propertyPartDataElements; + if ( resultType == RigFemResultPosEnum::RIG_NODAL ) return &m_propertyPartDataNodes; + + return nullptr; } diff --git a/ApplicationLibCode/GeoMech/OdbReader/RifInpReader.h b/ApplicationLibCode/GeoMech/OdbReader/RifInpReader.h index b28c2fcf9e..d1d94b0ed1 100644 --- a/ApplicationLibCode/GeoMech/OdbReader/RifInpReader.h +++ b/ApplicationLibCode/GeoMech/OdbReader/RifInpReader.h @@ -20,16 +20,39 @@ #include "RifGeoMechReaderInterface.h" +#include "RigFemResultPosEnum.h" +#include "RigFemTypes.h" + +#include #include #include #include +#include class RigFemPartCollection; +struct RifInpIncludeEntry +{ +public: + RifInpIncludeEntry( std::string propertyName, RigFemResultPosEnum resultType, int stepId, std::string fileName, int columnIndex ) + { + this->propertyName = propertyName; + this->stepId = stepId; + this->fileName = fileName; + this->resultType = resultType; + this->columnIndex = columnIndex; + } + +public: + std::string propertyName; + int stepId; + std::string fileName; + RigFemResultPosEnum resultType; + int columnIndex; +}; + //================================================================================================== // -// Data interface base class -// //================================================================================================== class RifInpReader : public RifGeoMechReaderInterface { @@ -37,6 +60,8 @@ public: RifInpReader(); ~RifInpReader() override; + void enableIncludes( bool enable ); + bool openFile( const std::string& fileName, std::string* errorMessage ) override; bool isOpen() const override; bool readFemParts( RigFemPartCollection* geoMechCase ) override; @@ -49,6 +74,7 @@ public: std::vector elementSet( int partIndex, std::string partName, int setIndex ) override; std::map> scalarNodeFieldAndComponentNames() override; + std::map> scalarElementFieldAndComponentNames() override; std::map> scalarElementNodeFieldAndComponentNames() override; std::map> scalarIntegrationPointFieldAndComponentNames() override; @@ -59,6 +85,11 @@ public: int stepIndex, int frameIndex, std::vector*>* resultValues ) override; + void readElementField( const std::string& fieldName, + int partIndex, + int stepIndex, + int frameIndex, + std::vector*>* resultValues ) override; void readElementNodeField( const std::string& fieldName, int partIndex, int stepIndex, @@ -70,9 +101,21 @@ public: int frameIndex, std::vector*>* resultValues ) override; + bool populateDerivedResultNames() const override; + private: void close(); + void readField( RigFemResultPosEnum resultType, + const std::string& fieldName, + int partIndex, + int stepIndex, + std::vector*>* resultValues ); + + void readScalarData( RigFemPartCollection* femParts, std::map& parts, std::vector& includeEntries ); + + std::map>>>* propertyDataMap( RigFemResultPosEnum resultType ); + static void skipComments( std::istream& stream ); static std::string parseLabel( const std::string& line, const std::string& labelName ); static std::vector> readNodes( std::istream& stream ); @@ -80,14 +123,22 @@ private: static std::vector readElementSet( std::istream& stream ); static std::vector readElementSetGenerate( std::istream& stream ); - static void read( std::istream& stream, - std::map& parts, - std::map>>& nodes, - std::map>>>& elements, - std::map>>>& elementSets ); + static RigElementType read( std::istream& stream, + std::map& parts, + std::map>>& nodes, + std::map>>>& elements, + std::map>>>& elementSets, + std::vector& stepNames, + bool enableIncludes, + std::vector& includeEntries ); private: - std::map> m_partElementSetNames; - - std::ifstream m_stream; + bool m_enableIncludes; + std::map> m_partElementSetNames; + std::vector m_stepNames; + std::vector m_includeEntries; + std::ifstream m_stream; + std::filesystem::path m_inputPath; + std::map>>> m_propertyPartDataNodes; + std::map>>> m_propertyPartDataElements; }; diff --git a/ApplicationLibCode/GeoMech/OdbReader/RifOdbReader.cpp b/ApplicationLibCode/GeoMech/OdbReader/RifOdbReader.cpp index 02a5a853b1..cf5bd13f84 100644 --- a/ApplicationLibCode/GeoMech/OdbReader/RifOdbReader.cpp +++ b/ApplicationLibCode/GeoMech/OdbReader/RifOdbReader.cpp @@ -101,43 +101,6 @@ private: size_t RifOdbReader::sm_instanceCount = 0; -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::map initFemTypeMap() -{ - std::map typeMap; - typeMap["C3D8R"] = HEX8; - typeMap["C3D8"] = HEX8; - typeMap["C3D8P"] = HEX8P; - typeMap["CAX4"] = CAX4; - typeMap["C3D20RT"] = HEX8; - typeMap["C3D8RT"] = HEX8; - typeMap["C3D8T"] = HEX8; - - return typeMap; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -RigElementType toRigElementType( const odb_String& odbTypeName ) -{ - static std::map odbElmTypeToRigElmTypeMap = initFemTypeMap(); - - std::map::iterator it = odbElmTypeToRigElmTypeMap.find( odbTypeName.cStr() ); - - if ( it == odbElmTypeToRigElmTypeMap.end() ) - { -#if 0 - std::cout << "Unsupported element type :" << odbElm.type().cStr() << std::endl; -#endif - return UNKNOWN_ELM_TYPE; - } - - return it->second; -} - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -176,6 +139,14 @@ void RifOdbReader::close() } } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RifOdbReader::populateDerivedResultNames() const +{ + return true; +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -380,8 +351,8 @@ bool RifOdbReader::readFemParts( RigFemPartCollection* femParts ) elementIdToIdxMap[odbElm.label()] = elmIdx; - RigElementType elmType = toRigElementType( odbElm.type() ); - if ( elmType == UNKNOWN_ELM_TYPE ) continue; + RigElementType elmType = RigFemTypes::toRigElementType( std::string( odbElm.type().cStr() ) ); + if ( elmType == RigElementType::UNKNOWN_ELM_TYPE ) continue; int nodeCount = 0; const int* idBasedConnectivities = odbElm.connectivity( nodeCount ); @@ -594,6 +565,15 @@ std::map> RifOdbReader::scalarNodeFieldAnd return fieldAndComponentNames( NODAL ); } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::map> RifOdbReader::scalarElementFieldAndComponentNames() +{ + // not supported, yet + return {}; +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -864,6 +844,20 @@ void RifOdbReader::readNodeField( const std::string& fieldName, } } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RifOdbReader::readElementField( const std::string& fieldName, + int partIndex, + int stepIndex, + int frameIndex, + std::vector*>* resultValues ) +{ + CVF_ASSERT( resultValues ); + + // Not supported, yet +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -986,7 +980,7 @@ void RifOdbReader::readIntegrationPointField( const std::string& int* elementLabels = bulkData.elementLabels(); float* data = bulkDataGetter.data(); - RigElementType eType = toRigElementType( bulkData.baseElementType() ); + RigElementType eType = RigFemTypes::toRigElementType( bulkData.baseElementType().cStr() ); const int* elmNodeToIpResultMapping = RigFemTypes::localElmNodeToIntegrationPointMapping( eType ); for ( int elem = 0; elem < elemCount; elem++ ) diff --git a/ApplicationLibCode/GeoMech/OdbReader/RifOdbReader.h b/ApplicationLibCode/GeoMech/OdbReader/RifOdbReader.h index 02193b43cf..1ca07b766c 100644 --- a/ApplicationLibCode/GeoMech/OdbReader/RifOdbReader.h +++ b/ApplicationLibCode/GeoMech/OdbReader/RifOdbReader.h @@ -53,6 +53,7 @@ public: std::vector elementSet( int partIndex, std::string partName, int setIndex ) override; std::map> scalarNodeFieldAndComponentNames() override; + std::map> scalarElementFieldAndComponentNames() override; std::map> scalarElementNodeFieldAndComponentNames() override; std::map> scalarIntegrationPointFieldAndComponentNames() override; @@ -63,6 +64,11 @@ public: int stepIndex, int frameIndex, std::vector*>* resultValues ) override; + void readElementField( const std::string& fieldName, + int partIndex, + int stepIndex, + int frameIndex, + std::vector*>* resultValues ) override; void readElementNodeField( const std::string& fieldName, int partIndex, int stepIndex, @@ -74,6 +80,8 @@ public: int frameIndex, std::vector*>* resultValues ) override; + bool populateDerivedResultNames() const override; + private: enum ResultPosition { diff --git a/ApplicationLibCode/ModelVisualization/Intersections/RivFemIntersectionGrid.cpp b/ApplicationLibCode/ModelVisualization/Intersections/RivFemIntersectionGrid.cpp index b9219fa01d..854fb7d564 100644 --- a/ApplicationLibCode/ModelVisualization/Intersections/RivFemIntersectionGrid.cpp +++ b/ApplicationLibCode/ModelVisualization/Intersections/RivFemIntersectionGrid.cpp @@ -110,7 +110,7 @@ void RivFemIntersectionGrid::cellCornerIndices( size_t globalCellIndex, size_t c auto [part, elementIdx] = m_femParts->partAndElementIndex( globalCellIndex ); RigElementType elmType = part->elementType( elementIdx ); - if ( elmType != HEX8 && elmType != HEX8P ) return; + if ( !RigFemTypes::is8NodeElement( elmType ) ) return; int elmIdx = static_cast( elementIdx ); const int partId = part->elementPartId(); diff --git a/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechContourMapProjection.cpp b/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechContourMapProjection.cpp index 3fad4a4409..1575ea7963 100644 --- a/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechContourMapProjection.cpp +++ b/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechContourMapProjection.cpp @@ -542,7 +542,7 @@ std::vector RimGeoMechContourMapProjection::gridCellValues( RigFemResult for ( size_t globalCellIdx = 0; globalCellIdx < static_cast( m_femPart->elementCount() ); ++globalCellIdx ) { RigElementType elmType = m_femPart->elementType( globalCellIdx ); - if ( elmType != HEX8 && elmType != HEX8P ) continue; + if ( !RigFemTypes::is8NodeElement( elmType ) ) continue; if ( resAddr.resultPosType == RIG_ELEMENT ) { diff --git a/ApplicationLibCode/ReservoirDataModel/RigCaseToCaseCellMapperTools.cpp b/ApplicationLibCode/ReservoirDataModel/RigCaseToCaseCellMapperTools.cpp index 910c0c2caf..7c334efb48 100644 --- a/ApplicationLibCode/ReservoirDataModel/RigCaseToCaseCellMapperTools.cpp +++ b/ApplicationLibCode/ReservoirDataModel/RigCaseToCaseCellMapperTools.cpp @@ -277,7 +277,7 @@ int RigCaseToCaseCellMapperTools::quadVxClosestToXYOfPoint( const cvf::Vec3d poi bool RigCaseToCaseCellMapperTools::elementCorners( const RigFemPart* femPart, int elmIdx, cvf::Vec3d elmCorners[8] ) { RigElementType elmType = femPart->elementType( elmIdx ); - if ( elmType != HEX8 && elmType != HEX8P ) return false; + if ( !RigFemTypes::is8NodeElement( elmType ) ) return false; const std::vector& nodeCoords = femPart->nodes().coordinates; const int* cornerIndices = femPart->connectivities( elmIdx ); @@ -300,7 +300,8 @@ bool RigCaseToCaseCellMapperTools::elementCorners( const RigFemPart* femPart, in int RigCaseToCaseCellMapperTools::findMatchingPOSKFaceIdx( const cvf::Vec3d baseCell[8], bool isBaseCellNormalsOutwards, const cvf::Vec3d c2[8] ) { int faceNodeCount; - const int* posKFace = RigFemTypes::localElmNodeIndicesForFace( HEX8, (int)( cvf::StructGridInterface::POS_K ), &faceNodeCount ); + const int* posKFace = + RigFemTypes::localElmNodeIndicesForFace( RigElementType::HEX8, (int)( cvf::StructGridInterface::POS_K ), &faceNodeCount ); double sign = isBaseCellNormalsOutwards ? 1.0 : -1.0; @@ -311,7 +312,7 @@ int RigCaseToCaseCellMapperTools::findMatchingPOSKFaceIdx( const cvf::Vec3d base int bestFace = -1; for ( int faceIdx = 5; faceIdx >= 0; --faceIdx ) // Backwards. might hit earlier more often { - const int* face = RigFemTypes::localElmNodeIndicesForFace( HEX8, faceIdx, &faceNodeCount ); + const int* face = RigFemTypes::localElmNodeIndicesForFace( RigElementType::HEX8, faceIdx, &faceNodeCount ); cvf::Vec3d normal = ( c2[face[2]] - c2[face[0]] ) ^ ( c2[face[3]] - c2[face[1]] ); normal.normalize(); double sqDiff = ( posKnormal - normal ).lengthSquared(); @@ -367,11 +368,13 @@ void RigCaseToCaseCellMapperTools::rotateCellTopologicallyToMatchBaseCell( const tmpFemCorners[6] = cell[6]; tmpFemCorners[7] = cell[7]; - int femShallowZFaceIdx = RigFemTypes::oppositeFace( HEX8, femDeepZFaceIdx ); + int femShallowZFaceIdx = RigFemTypes::oppositeFace( RigElementType::HEX8, femDeepZFaceIdx ); int faceNodeCount; - const int* localElmNodeIndicesForPOSKFace = RigFemTypes::localElmNodeIndicesForFace( HEX8, femDeepZFaceIdx, &faceNodeCount ); - const int* localElmNodeIndicesForNEGKFace = RigFemTypes::localElmNodeIndicesForFace( HEX8, femShallowZFaceIdx, &faceNodeCount ); + const int* localElmNodeIndicesForPOSKFace = + RigFemTypes::localElmNodeIndicesForFace( RigElementType::HEX8, femDeepZFaceIdx, &faceNodeCount ); + const int* localElmNodeIndicesForNEGKFace = + RigFemTypes::localElmNodeIndicesForFace( RigElementType::HEX8, femShallowZFaceIdx, &faceNodeCount ); cell[0] = tmpFemCorners[localElmNodeIndicesForNEGKFace[0]]; cell[1] = tmpFemCorners[localElmNodeIndicesForNEGKFace[1]]; diff --git a/ApplicationLibCode/ReservoirDataModel/RigGeoMechWellLogExtractor.cpp b/ApplicationLibCode/ReservoirDataModel/RigGeoMechWellLogExtractor.cpp index 96dd3a4a09..9c2fd324db 100644 --- a/ApplicationLibCode/ReservoirDataModel/RigGeoMechWellLogExtractor.cpp +++ b/ApplicationLibCode/ReservoirDataModel/RigGeoMechWellLogExtractor.cpp @@ -874,7 +874,7 @@ T RigGeoMechWellLogExtractor::interpolateGridResultValue( RigFemResultPosEnum size_t elmIdx = intersectedCellsGlobIdx()[intersectionIdx]; RigElementType elmType = femPart->elementType( elmIdx ); - if ( elmType != HEX8 && elmType != HEX8P ) return T(); + if ( !RigFemTypes::is8NodeElement( elmType ) ) return T(); if ( resultPosType == RIG_FORMATION_NAMES ) { @@ -986,7 +986,7 @@ void RigGeoMechWellLogExtractor::calculateIntersection() for ( size_t ccIdx = 0; ccIdx < closeCells.size(); ++ccIdx ) { RigElementType elmType = femPart->elementType( closeCells[ccIdx] ); - if ( elmType != HEX8 && elmType != HEX8P ) continue; + if ( elmType != RigElementType::HEX8 && elmType != RigElementType::HEX8P ) continue; const int* cornerIndices = femPart->connectivities( closeCells[ccIdx] ); @@ -1249,7 +1249,7 @@ std::vector RigGeoMechWellLogExtractor::interpolateInterfaceValues( RigFemRes { size_t elmIdx = intersectedCellsGlobIdx()[intersectionIdx]; RigElementType elmType = femPart->elementType( elmIdx ); - if ( elmType != HEX8 && elmType != HEX8P ) continue; + if ( !RigFemTypes::is8NodeElement( elmType ) ) continue; interpolatedInterfaceValues[intersectionIdx] = interpolateGridResultValue( nativeAddr.resultPosType, unscaledResultValues, intersectionIdx ); diff --git a/ApplicationLibCode/UserInterface/RiuFemResultTextBuilder.cpp b/ApplicationLibCode/UserInterface/RiuFemResultTextBuilder.cpp index 8ef22bbfbe..0f1f0b8281 100644 --- a/ApplicationLibCode/UserInterface/RiuFemResultTextBuilder.cpp +++ b/ApplicationLibCode/UserInterface/RiuFemResultTextBuilder.cpp @@ -135,7 +135,8 @@ QString RiuFemResultTextBuilder::geometrySelectionText( QString itemSeparator ) int elementId = femPart->elmId( m_cellIndex ); auto elementType = femPart->elementType( m_cellIndex ); - text += QString( "Element : Id[%1], Type[%2]" ).arg( elementId ).arg( RigFemTypes::elementTypeText( elementType ) ); + text += + QString( "Element : Id[%1], Type[%2]" ).arg( elementId ).arg( QString::fromStdString( RigFemTypes::elementTypeText( elementType ) ) ); size_t i = 0; size_t j = 0; diff --git a/ApplicationLibCode/UserInterface/RiuMohrsCirclePlot.cpp b/ApplicationLibCode/UserInterface/RiuMohrsCirclePlot.cpp index a4daa17ed7..bf7fbbb581 100644 --- a/ApplicationLibCode/UserInterface/RiuMohrsCirclePlot.cpp +++ b/ApplicationLibCode/UserInterface/RiuMohrsCirclePlot.cpp @@ -371,7 +371,7 @@ bool RiuMohrsCirclePlot::addOrUpdateCurves( const RimGeoMechResultDefinition* ge { RigFemPart* femPart = geomResDef->ownerCaseData()->femParts()->part( gridIndex ); - if ( femPart->elementType( elmIndex ) != HEX8P ) return false; + if ( femPart->elementType( elmIndex ) != RigElementType::HEX8P ) return false; RigFemPartResultsCollection* resultCollection = geomResDef->geoMechCase()->geoMechData()->femPartResults();