OilfieldToolsNet/ResInsight
4
0
Fork 0
mirror of https://github.com/OPM/ResInsight.git synced 2025-02-25 18:55:39 -06:00
Code Issues Packages Projects Releases Wiki Activity

#939 Added principal strains

Browse Source
This commit is contained in:
Jacob Støren
2016-10-31 10:24:17 +01:00
parent fc5ce8af1b
commit 91db241ee1
2 changed files with 215 additions and 115 deletions
Download Patch File Download Diff File Expand all files Collapse all files

411
ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartResultsCollection.cpp
View File

@@ -259,6 +259,9 @@ std::map<std::string, std::vector<std::string> > RigFemPartResultsCollection::sc
fieldCompNames["NE"].push_back("E12"); fieldCompNames["NE"].push_back("E12");
fieldCompNames["NE"].push_back("E13"); fieldCompNames["NE"].push_back("E13");
fieldCompNames["NE"].push_back("E23"); fieldCompNames["NE"].push_back("E23");
fieldCompNames["NE"].push_back("E1");
fieldCompNames["NE"].push_back("E2");
fieldCompNames["NE"].push_back("E3");
} }
else if (resPos == RIG_INTEGRATION_POINT) else if (resPos == RIG_INTEGRATION_POINT)
@@ -319,6 +322,9 @@ std::map<std::string, std::vector<std::string> > RigFemPartResultsCollection::sc
fieldCompNames["NE"].push_back("E12"); fieldCompNames["NE"].push_back("E12");
fieldCompNames["NE"].push_back("E13"); fieldCompNames["NE"].push_back("E13");
fieldCompNames["NE"].push_back("E23"); fieldCompNames["NE"].push_back("E23");
fieldCompNames["NE"].push_back("E1");
fieldCompNames["NE"].push_back("E2");
fieldCompNames["NE"].push_back("E3");
} }
else if(resPos == RIG_ELEMENT_NODAL_FACE) else if(resPos == RIG_ELEMENT_NODAL_FACE)
@@ -885,169 +891,16 @@ RigFemScalarResultFrames* RigFemPartResultsCollection::calculateSurfaceAlignedSt
//-------------------------------------------------------------------------------------------------- //--------------------------------------------------------------------------------------------------
/// ///
//-------------------------------------------------------------------------------------------------- //--------------------------------------------------------------------------------------------------
RigFemScalarResultFrames* RigFemPartResultsCollection::calculateDerivedResult(int partIndex, const RigFemResultAddress& resVarAddr) RigFemScalarResultFrames* RigFemPartResultsCollection::calculatePrincipalStressValues(int partIndex, const RigFemResultAddress &resVarAddr )
{ {
if (resVarAddr.isTimeLapse()) CVF_ASSERT(resVarAddr.componentName == "S1" || resVarAddr.componentName == "S2" || resVarAddr.componentName == "S3"
{
return calculateTimeLapseResult(partIndex, resVarAddr);
}
if(resVarAddr.resultPosType == RIG_ELEMENT_NODAL_FACE )
{
return calculateSurfaceAlignedStress(partIndex, resVarAddr);
}
if (resVarAddr.fieldName == "SE" && resVarAddr.componentName == "SFI")
{
return calculateSFI(partIndex, resVarAddr);
}
if(resVarAddr.fieldName == "NE" && resVarAddr.componentName == "EV")
{
return calculateVolumetricStrain(partIndex, resVarAddr);
}
if(resVarAddr.fieldName == "ST" && resVarAddr.componentName == "Q" )
{
return calculateDeviatoricStress(partIndex, resVarAddr);
}
if(resVarAddr.fieldName == "ST" && resVarAddr.componentName == "STM")
{
return calculateMeanStressSTM(partIndex, resVarAddr);
}
if(resVarAddr.fieldName == "SE" && resVarAddr.componentName == "SEM")
{
return calculateMeanStressSEM(partIndex, resVarAddr);
}
if (resVarAddr.fieldName == "S-Bar")
{
return calculateBarConvertedResult(partIndex, resVarAddr, "S");
}
if (resVarAddr.fieldName == "POR-Bar")
{
if (resVarAddr.resultPosType == RIG_NODAL)
return calculateBarConvertedResult(partIndex, resVarAddr, "POR");
else
return calculateEnIpPorBarResult(partIndex, resVarAddr);
}
if((resVarAddr.fieldName == "NE")
&& ( resVarAddr.componentName == "E11"
|| resVarAddr.componentName == "E22"
|| resVarAddr.componentName == "E33"
|| resVarAddr.componentName == "E12"
|| resVarAddr.componentName == "E13"
|| resVarAddr.componentName == "E23"))
{
caf::ProgressInfo frameCountProgress(this->frameCount() * 2, "");
frameCountProgress.setProgressDescription("Calculating " + QString::fromStdString(resVarAddr.fieldName + ": " + resVarAddr.componentName));
frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * srcDataFrames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(resVarAddr.resultPosType, "E", resVarAddr.componentName));
RigFemScalarResultFrames * dstDataFrames = m_femPartResults[partIndex]->createScalarResult(resVarAddr);
frameCountProgress.incrementProgress();
int frameCount = srcDataFrames->frameCount();
for (int fIdx = 0; fIdx < frameCount; ++fIdx)
{
const std::vector<float>& srcFrameData = srcDataFrames->frameData(fIdx);
std::vector<float>& dstFrameData = dstDataFrames->frameData(fIdx);
size_t valCount = srcFrameData.size();
dstFrameData.resize(valCount);
for (size_t vIdx = 0; vIdx < valCount; ++vIdx)
{
dstFrameData[vIdx] = -srcFrameData[vIdx];
}
frameCountProgress.incrementProgress();
}
return dstDataFrames;
}
if ((resVarAddr.fieldName == "SE")
&& ( resVarAddr.componentName == "S11"
|| resVarAddr.componentName == "S22"
|| resVarAddr.componentName == "S33"
|| resVarAddr.componentName == "S12"
|| resVarAddr.componentName == "S13"
|| resVarAddr.componentName == "S23" ))
{
caf::ProgressInfo frameCountProgress(this->frameCount() * 3, "");
frameCountProgress.setProgressDescription("Calculating " + QString::fromStdString(resVarAddr.fieldName + ": " + resVarAddr.componentName));
frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * srcDataFrames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(resVarAddr.resultPosType, "S-Bar", resVarAddr.componentName));
frameCountProgress.incrementProgress(); frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * srcPORDataFrames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(RIG_NODAL, "POR-Bar", ""));
RigFemScalarResultFrames * dstDataFrames = m_femPartResults[partIndex]->createScalarResult(resVarAddr);
frameCountProgress.incrementProgress();
const RigFemPart * femPart = m_femParts->part(partIndex);
float inf = std::numeric_limits<float>::infinity();
int frameCount = srcDataFrames->frameCount();
for (int fIdx = 0; fIdx < frameCount; ++fIdx)
{
const std::vector<float>& srcSFrameData = srcDataFrames->frameData(fIdx);
std::vector<float>& dstFrameData = dstDataFrames->frameData(fIdx);
size_t valCount = srcSFrameData.size();
dstFrameData.resize(valCount);
const std::vector<float>& srcPORFrameData = srcPORDataFrames->frameData(fIdx);
int elementCount = femPart->elementCount();
for (int elmIdx = 0; elmIdx < elementCount; ++elmIdx)
{
RigElementType elmType = femPart->elementType(elmIdx);
int elmNodeCount = RigFemTypes::elmentNodeCount(femPart->elementType(elmIdx));
if (elmType == HEX8P)
{
for (int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx)
{
size_t elmNodResIdx = femPart->elementNodeResultIdx(elmIdx, elmNodIdx);
dstFrameData[elmNodResIdx] = -srcSFrameData[elmNodResIdx];
}
}
else
{
for (int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx)
{
size_t elmNodResIdx = femPart->elementNodeResultIdx(elmIdx, elmNodIdx);
dstFrameData[elmNodResIdx] = inf;
}
}
}
frameCountProgress.incrementProgress();
}
return dstDataFrames;
}
if ( (resVarAddr.fieldName == "SE" || resVarAddr.fieldName == "ST" )
&& ( resVarAddr.componentName == "S1"
|| resVarAddr.componentName == "S2"
|| resVarAddr.componentName == "S3"
|| resVarAddr.componentName == "S1inc" || resVarAddr.componentName == "S1inc"
|| resVarAddr.componentName == "S1azi" || resVarAddr.componentName == "S1azi"
|| resVarAddr.componentName == "S2inc" || resVarAddr.componentName == "S2inc"
|| resVarAddr.componentName == "S2azi" || resVarAddr.componentName == "S2azi"
|| resVarAddr.componentName == "S3inc" || resVarAddr.componentName == "S3inc"
|| resVarAddr.componentName == "S3azi" ) || resVarAddr.componentName == "S3azi");
)
{
caf::ProgressInfo frameCountProgress(this->frameCount() * 7, ""); caf::ProgressInfo frameCountProgress(this->frameCount() * 7, "");
frameCountProgress.setProgressDescription("Calculating " + QString::fromStdString(resVarAddr.fieldName + ": " + resVarAddr.componentName)); frameCountProgress.setProgressDescription("Calculating " + QString::fromStdString(resVarAddr.fieldName + ": " + resVarAddr.componentName));
frameCountProgress.setNextProgressIncrement(this->frameCount()); frameCountProgress.setNextProgressIncrement(this->frameCount());
@@ -1165,6 +1018,250 @@ RigFemScalarResultFrames* RigFemPartResultsCollection::calculateDerivedResult(in
return requestedPrincipal; return requestedPrincipal;
} }
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigFemScalarResultFrames* RigFemPartResultsCollection::calculatePrincipalStrainValues(int partIndex, const RigFemResultAddress &resVarAddr)
{
CVF_ASSERT(resVarAddr.componentName == "E1" || resVarAddr.componentName == "E2" || resVarAddr.componentName == "E3");
caf::ProgressInfo frameCountProgress(this->frameCount() * 7, "");
frameCountProgress.setProgressDescription("Calculating " + QString::fromStdString(resVarAddr.fieldName + ": " + resVarAddr.componentName));
frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * s11Frames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(resVarAddr.resultPosType, resVarAddr.fieldName, "E11"));
frameCountProgress.incrementProgress(); frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * s22Frames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(resVarAddr.resultPosType, resVarAddr.fieldName, "E22"));
frameCountProgress.incrementProgress(); frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * s33Frames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(resVarAddr.resultPosType, resVarAddr.fieldName, "E33"));
frameCountProgress.incrementProgress(); frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * s12Frames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(resVarAddr.resultPosType, resVarAddr.fieldName, "E12"));
frameCountProgress.incrementProgress(); frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * s13Frames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(resVarAddr.resultPosType, resVarAddr.fieldName, "E13"));
frameCountProgress.incrementProgress(); frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * s23Frames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(resVarAddr.resultPosType, resVarAddr.fieldName, "E23"));
RigFemScalarResultFrames * s1Frames = m_femPartResults[partIndex]->createScalarResult(RigFemResultAddress(resVarAddr.resultPosType, resVarAddr.fieldName, "E1"));
RigFemScalarResultFrames * s2Frames = m_femPartResults[partIndex]->createScalarResult(RigFemResultAddress(resVarAddr.resultPosType, resVarAddr.fieldName, "E2"));
RigFemScalarResultFrames * s3Frames = m_femPartResults[partIndex]->createScalarResult(RigFemResultAddress(resVarAddr.resultPosType, resVarAddr.fieldName, "E3"));
frameCountProgress.incrementProgress();
int frameCount = s11Frames->frameCount();
for ( int fIdx = 0; fIdx < frameCount; ++fIdx )
{
const std::vector<float>& s11 = s11Frames->frameData(fIdx);
const std::vector<float>& s22 = s22Frames->frameData(fIdx);
const std::vector<float>& s33 = s33Frames->frameData(fIdx);
const std::vector<float>& s12 = s12Frames->frameData(fIdx);
const std::vector<float>& s13 = s13Frames->frameData(fIdx);
const std::vector<float>& s23 = s23Frames->frameData(fIdx);
std::vector<float>& s1 = s1Frames->frameData(fIdx);
std::vector<float>& s2 = s2Frames->frameData(fIdx);
std::vector<float>& s3 = s3Frames->frameData(fIdx);
size_t valCount = s11.size();
s1.resize(valCount);
s2.resize(valCount);
s3.resize(valCount);
for ( size_t vIdx = 0; vIdx < valCount; ++vIdx )
{
caf::Ten3f T(s11[vIdx], s22[vIdx], s33[vIdx], s12[vIdx], s23[vIdx], s13[vIdx]);
cvf::Vec3f principalDirs[3];
cvf::Vec3f principals = T.calculatePrincipals(principalDirs);
s1[vIdx] = principals[0];
s2[vIdx] = principals[1];
s3[vIdx] = principals[2];
}
frameCountProgress.incrementProgress();
}
RigFemScalarResultFrames* requestedPrincipal = this->findOrLoadScalarResult(partIndex, resVarAddr);
return requestedPrincipal;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigFemScalarResultFrames* RigFemPartResultsCollection::calculateDerivedResult(int partIndex, const RigFemResultAddress& resVarAddr)
{
if (resVarAddr.isTimeLapse())
{
return calculateTimeLapseResult(partIndex, resVarAddr);
}
if(resVarAddr.resultPosType == RIG_ELEMENT_NODAL_FACE )
{
return calculateSurfaceAlignedStress(partIndex, resVarAddr);
}
if (resVarAddr.fieldName == "SE" && resVarAddr.componentName == "SFI")
{
return calculateSFI(partIndex, resVarAddr);
}
if(resVarAddr.fieldName == "NE" && resVarAddr.componentName == "EV")
{
return calculateVolumetricStrain(partIndex, resVarAddr);
}
if(resVarAddr.fieldName == "ST" && resVarAddr.componentName == "Q" )
{
return calculateDeviatoricStress(partIndex, resVarAddr);
}
if(resVarAddr.fieldName == "ST" && resVarAddr.componentName == "STM")
{
return calculateMeanStressSTM(partIndex, resVarAddr);
}
if(resVarAddr.fieldName == "SE" && resVarAddr.componentName == "SEM")
{
return calculateMeanStressSEM(partIndex, resVarAddr);
}
if (resVarAddr.fieldName == "S-Bar")
{
return calculateBarConvertedResult(partIndex, resVarAddr, "S");
}
if (resVarAddr.fieldName == "POR-Bar")
{
if (resVarAddr.resultPosType == RIG_NODAL)
return calculateBarConvertedResult(partIndex, resVarAddr, "POR");
else
return calculateEnIpPorBarResult(partIndex, resVarAddr);
}
if((resVarAddr.fieldName == "NE")
&& ( resVarAddr.componentName == "E11"
|| resVarAddr.componentName == "E22"
|| resVarAddr.componentName == "E33"
|| resVarAddr.componentName == "E12"
|| resVarAddr.componentName == "E13"
|| resVarAddr.componentName == "E23"))
{
caf::ProgressInfo frameCountProgress(this->frameCount() * 2, "");
frameCountProgress.setProgressDescription("Calculating " + QString::fromStdString(resVarAddr.fieldName + ": " + resVarAddr.componentName));
frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * srcDataFrames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(resVarAddr.resultPosType, "E", resVarAddr.componentName));
RigFemScalarResultFrames * dstDataFrames = m_femPartResults[partIndex]->createScalarResult(resVarAddr);
frameCountProgress.incrementProgress();
int frameCount = srcDataFrames->frameCount();
for (int fIdx = 0; fIdx < frameCount; ++fIdx)
{
const std::vector<float>& srcFrameData = srcDataFrames->frameData(fIdx);
std::vector<float>& dstFrameData = dstDataFrames->frameData(fIdx);
size_t valCount = srcFrameData.size();
dstFrameData.resize(valCount);
for (size_t vIdx = 0; vIdx < valCount; ++vIdx)
{
dstFrameData[vIdx] = -srcFrameData[vIdx];
}
frameCountProgress.incrementProgress();
}
return dstDataFrames;
}
if ( (resVarAddr.fieldName == "NE")
&& ( resVarAddr.componentName == "E1"
|| resVarAddr.componentName == "E2"
|| resVarAddr.componentName == "E3") )
{
return calculatePrincipalStrainValues(partIndex, resVarAddr);
}
if ((resVarAddr.fieldName == "SE")
&& ( resVarAddr.componentName == "S11"
|| resVarAddr.componentName == "S22"
|| resVarAddr.componentName == "S33"
|| resVarAddr.componentName == "S12"
|| resVarAddr.componentName == "S13"
|| resVarAddr.componentName == "S23" ))
{
caf::ProgressInfo frameCountProgress(this->frameCount() * 3, "");
frameCountProgress.setProgressDescription("Calculating " + QString::fromStdString(resVarAddr.fieldName + ": " + resVarAddr.componentName));
frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * srcDataFrames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(resVarAddr.resultPosType, "S-Bar", resVarAddr.componentName));
frameCountProgress.incrementProgress(); frameCountProgress.setNextProgressIncrement(this->frameCount());
RigFemScalarResultFrames * srcPORDataFrames = this->findOrLoadScalarResult(partIndex, RigFemResultAddress(RIG_NODAL, "POR-Bar", ""));
RigFemScalarResultFrames * dstDataFrames = m_femPartResults[partIndex]->createScalarResult(resVarAddr);
frameCountProgress.incrementProgress();
const RigFemPart * femPart = m_femParts->part(partIndex);
float inf = std::numeric_limits<float>::infinity();
int frameCount = srcDataFrames->frameCount();
for (int fIdx = 0; fIdx < frameCount; ++fIdx)
{
const std::vector<float>& srcSFrameData = srcDataFrames->frameData(fIdx);
std::vector<float>& dstFrameData = dstDataFrames->frameData(fIdx);
size_t valCount = srcSFrameData.size();
dstFrameData.resize(valCount);
const std::vector<float>& srcPORFrameData = srcPORDataFrames->frameData(fIdx);
int elementCount = femPart->elementCount();
for (int elmIdx = 0; elmIdx < elementCount; ++elmIdx)
{
RigElementType elmType = femPart->elementType(elmIdx);
int elmNodeCount = RigFemTypes::elmentNodeCount(femPart->elementType(elmIdx));
if (elmType == HEX8P)
{
for (int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx)
{
size_t elmNodResIdx = femPart->elementNodeResultIdx(elmIdx, elmNodIdx);
dstFrameData[elmNodResIdx] = -srcSFrameData[elmNodResIdx];
}
}
else
{
for (int elmNodIdx = 0; elmNodIdx < elmNodeCount; ++elmNodIdx)
{
size_t elmNodResIdx = femPart->elementNodeResultIdx(elmIdx, elmNodIdx);
dstFrameData[elmNodResIdx] = inf;
}
}
}
frameCountProgress.incrementProgress();
}
return dstDataFrames;
}
if ( (resVarAddr.fieldName == "SE" || resVarAddr.fieldName == "ST" )
&& ( resVarAddr.componentName == "S1"
|| resVarAddr.componentName == "S2"
|| resVarAddr.componentName == "S3"
|| resVarAddr.componentName == "S1inc"
|| resVarAddr.componentName == "S1azi"
|| resVarAddr.componentName == "S2inc"
|| resVarAddr.componentName == "S2azi"
|| resVarAddr.componentName == "S3inc"
|| resVarAddr.componentName == "S3azi" )
)
{
return calculatePrincipalStressValues(partIndex, resVarAddr);
}
if ( resVarAddr.fieldName == "ST" if ( resVarAddr.fieldName == "ST"
&& ( resVarAddr.componentName == "S11" && ( resVarAddr.componentName == "S11"
|| resVarAddr.componentName == "S22" || resVarAddr.componentName == "S22"

3
ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartResultsCollection.h
View File

@@ -77,6 +77,7 @@ private:
RigFemScalarResultFrames* calculateDerivedResult(int partIndex, const RigFemResultAddress& resVarAddr); RigFemScalarResultFrames* calculateDerivedResult(int partIndex, const RigFemResultAddress& resVarAddr);
void calculateGammaFromFrames(int partIndex, void calculateGammaFromFrames(int partIndex,
const RigFemScalarResultFrames * totalStressComponentDataFrames, const RigFemScalarResultFrames * totalStressComponentDataFrames,
const RigFemScalarResultFrames * srcPORDataFrames, const RigFemScalarResultFrames * srcPORDataFrames,
@@ -92,6 +93,8 @@ private:
RigFemScalarResultFrames* calculateDeviatoricStress(int partIndex, const RigFemResultAddress& resVarAddr); RigFemScalarResultFrames* calculateDeviatoricStress(int partIndex, const RigFemResultAddress& resVarAddr);
RigFemScalarResultFrames* calculateVolumetricStrain(int partIndex, const RigFemResultAddress& resVarAddr); RigFemScalarResultFrames* calculateVolumetricStrain(int partIndex, const RigFemResultAddress& resVarAddr);
RigFemScalarResultFrames* calculateSurfaceAlignedStress(int partIndex, const RigFemResultAddress& resVarAddr); RigFemScalarResultFrames* calculateSurfaceAlignedStress(int partIndex, const RigFemResultAddress& resVarAddr);
RigFemScalarResultFrames* calculatePrincipalStressValues(int partIndex, const RigFemResultAddress &resVarAddr);
RigFemScalarResultFrames* calculatePrincipalStrainValues(int partIndex, const RigFemResultAddress &resVarAddr);
cvf::Collection<RigFemPartResults> m_femPartResults; cvf::Collection<RigFemPartResults> m_femPartResults;
cvf::ref<RifGeoMechReaderInterface> m_readerInterface; cvf::ref<RifGeoMechReaderInterface> m_readerInterface;