ResInsight/ApplicationCode/UserInterface/RiuFemTimeHistoryResultAccessor.cpp

/////////////////////////////////////////////////////////////////////////////////
//
//  Copyright (C) Statoil ASA
//  Copyright (C) Ceetron Solutions AS
//
//  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 <http://www.gnu.org/licenses/gpl.html>
//  for more details.
//
/////////////////////////////////////////////////////////////////////////////////

#include "RiuFemTimeHistoryResultAccessor.h"

#include "RigFemClosestResultIndexCalculator.h"
#include "RigFemPart.h"
#include "RigFemPartCollection.h"
#include "RigFemPartGrid.h"
#include "RigFemPartResultsCollection.h"
#include "RigFemTypes.h"
#include "RigGeoMechCaseData.h"
#include "RiuGeoMechXfTensorResultAccessor.h"

#include <cmath> // Needed for HUGE_VAL on Linux

//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RiuFemTimeHistoryResultAccessor::RiuFemTimeHistoryResultAccessor( RigGeoMechCaseData* geomData,
                                                                  RigFemResultAddress femResultAddress,
                                                                  size_t              gridIndex,
                                                                  int                 elementIndex,
                                                                  int                 face,
                                                                  const cvf::Vec3d&   intersectionPointInDomain )
    : m_geoMechCaseData( geomData )
    , m_femResultAddress( new RigFemResultAddress( femResultAddress ) )
    , m_gridIndex( gridIndex )
    , m_elementIndex( elementIndex )
    , m_face( face )
    , m_intersectionPointInDomain( intersectionPointInDomain )
    , m_hasIntersectionTriangle( false )
{
    computeTimeHistoryData();
}

//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RiuFemTimeHistoryResultAccessor::RiuFemTimeHistoryResultAccessor( RigGeoMechCaseData* geomData,
                                                                  RigFemResultAddress femResultAddress,
                                                                  size_t              gridIndex,
                                                                  int                 elementIndex,
                                                                  int                 face,
                                                                  const cvf::Vec3d&   intersectionPointInDomain,
                                                                  const std::array<cvf::Vec3f, 3>& intersectionTriangle )
    : m_geoMechCaseData( geomData )
    , m_femResultAddress( new RigFemResultAddress( femResultAddress ) )
    , m_gridIndex( gridIndex )
    , m_elementIndex( elementIndex )
    , m_face( face )
    , m_intersectionPointInDomain( intersectionPointInDomain )
    , m_hasIntersectionTriangle( true )
    , m_intersectionTriangle( intersectionTriangle )
{
    computeTimeHistoryData();
}

//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
QString RiuFemTimeHistoryResultAccessor::geometrySelectionText() const
{
    QString text;

    if ( m_geoMechCaseData )
    {
        RigFemPart* femPart   = m_geoMechCaseData->femParts()->part( m_gridIndex );
        int         elementId = femPart->elmId( m_elementIndex );
        text += QString( "Element : Id[%1]" ).arg( elementId );

        size_t i = 0;
        size_t j = 0;
        size_t k = 0;
        if ( m_geoMechCaseData->femParts()
                 ->part( m_gridIndex )
                 ->getOrCreateStructGrid()
                 ->ijkFromCellIndex( m_elementIndex, &i, &j, &k ) )
        {
            // Adjust to 1-based Eclipse indexing
            i++;
            j++;
            k++;

            cvf::Vec3d domainCoord = m_intersectionPointInDomain;
            text += QString( ", ijk[%1, %2, %3] " ).arg( i ).arg( j ).arg( k );

            QString formattedText;
            formattedText.sprintf( "Intersection point : [E: %.2f, N: %.2f, Depth: %.2f]",
                                   domainCoord.x(),
                                   domainCoord.y(),
                                   -domainCoord.z() );

            text += formattedText;
        }
    }

    return text;
}

//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::vector<double> RiuFemTimeHistoryResultAccessor::timeHistoryValues() const
{
    return m_timeHistoryValues;
}

//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
int RiuFemTimeHistoryResultAccessor::closestNodeId() const
{
    return m_closestNodeId;
}

//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RiuFemTimeHistoryResultAccessor::computeTimeHistoryData()
{
    m_timeHistoryValues.clear();

    RigFemClosestResultIndexCalculator closestCalc( m_geoMechCaseData->femParts()->part( m_gridIndex ),
                                                    m_femResultAddress->resultPosType,
                                                    m_elementIndex,
                                                    m_face,
                                                    m_intersectionPointInDomain );

    int scalarResultIndex = closestCalc.resultIndexToClosestResult();
    m_closestNodeId       = closestCalc.closestNodeId();

    RigFemPartResultsCollection* femPartResultsColl = m_geoMechCaseData->femPartResults();

    if ( m_femResultAddress->resultPosType == RIG_ELEMENT_NODAL_FACE && m_hasIntersectionTriangle )
    {
        int closestElmNodeResIndex = closestCalc.closestElementNodeResIdx();

        for ( int frameIdx = 0; frameIdx < femPartResultsColl->frameCount(); frameIdx++ )
        {
            RiuGeoMechXfTensorResultAccessor stressXfAccessor( femPartResultsColl, *m_femResultAddress, frameIdx );
            float scalarValue = stressXfAccessor.calculateElmNodeValue( m_intersectionTriangle, closestElmNodeResIndex );
            m_timeHistoryValues.push_back( scalarValue );
        }
    }
    else
    {
        if ( scalarResultIndex < 0 ) return;

        for ( int frameIdx = 0; frameIdx < femPartResultsColl->frameCount(); frameIdx++ )
        {
            const std::vector<float>& scalarResults =
                m_geoMechCaseData->femPartResults()->resultValues( *m_femResultAddress,
                                                                   static_cast<int>( m_gridIndex ),
                                                                   frameIdx );
            if ( scalarResults.size() )
            {
                float scalarValue = scalarResults[scalarResultIndex];

                m_timeHistoryValues.push_back( scalarValue );
            }
        }
    }
}