gmpMaterialMDCreep.h

Go to the documentation of this file.

/************************************************************************
**
** Copyright (C) 2014 by Carlos Augusto Teixera Mendes
** All rights reserved.
**
** This file is part of the "GeMA" software. It's use should respect
** the terms in the license agreement that can be found together 
** with this source code.
** It is provided AS IS, with NO WARRANTY OF ANY KIND, 
** INCLUDING THE WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR 
** A PARTICULAR PURPOSE.
**
************************************************************************/

#ifndef _GEMA_PLUGIN_MECHANICALMATERIAL_MDCREEP_H_
#define _GEMA_PLUGIN_MECHANICALMATERIAL_MDCREEP_H_

#include "gmpMechanicalConfig.h"
#include <gmpFemPhysics.h>
#include <gmMathUtils.h>

#include <gmpMaterialDPCreep.h>
#include "gmpMechanicPoint.h"

class GMP_MECHANICAL_PHYSICS_API_EXPORT GmpMaterialMDCreep : public GmpMaterialDPCreep
{
protected:
  enum ElementPropertyIds
  {
    B1_ID = GmpMaterialDPCreep::ElementPropertyIds::NUM_PROPERTY_IDS, 
    B2_ID,                              
    Xq_ID,                              
    Xc_ID,                              
    Xm_ID,                              
    XKo_ID,                             
    alphaH_ID,          
    betaH_ID,                   
    deltaMin_ID,        

    NUM_PROPERTY_IDS,   
  };

  enum GaussAttributeIds
  {
    zeta_GA_ID = GmpMaterialDPCreep::GaussAttributeIds::NUM_GA_IDS,  
    zetaOld_GA_ID,    
    FT_GA_ID,         
    EcSS_GA_ID,       

    NUM_GA_IDS,  
  };
  
public:

  GmpMaterialMDCreep(int typeIndex, QString typeName, const GmLogCategory& logger)
    : GmpMaterialDPCreep(typeIndex, typeName, logger) {}

  virtual ~GmpMaterialMDCreep() {}

  static GmpFemPhysicsCommonMaterial* instance(GmSimulationData* simulation, int typeIndex, QString typeName, const GmLogCategory& logger)
  {
    Q_UNUSED(simulation);
    return new GmpMaterialMDCreep(typeIndex, typeName, logger);
  }
 
  // Returns a map with material associated properties. 
  virtual const QVariantMap* materialMetaDataMap()
  {
    S_TRACE();
    static QVariantMap m;
    if (m.isEmpty())
    {
      // read the material attribute map from elastic interface material
      m = *GmpMaterialDPCreep::materialMetaDataMap();
          
      // Adds material properties
      m["properties"] = m["properties"].value<GmpFemPhysicsCommon::ValueList>()                         
      <<GmpFemPhysicsCommon::ScalarValue(B1_ID, "B1", QObject::tr("Structure factor for DGL 1"), "1/s", true)                   
      <<GmpFemPhysicsCommon::ScalarValue(B2_ID, "B2", QObject::tr("Structure factor for DGL 2"), "1/s", true)                                   
      <<GmpFemPhysicsCommon::ScalarValue(Xq_ID, "Xq", QObject::tr("Stress constant"), "", true) 
      <<GmpFemPhysicsCommon::ScalarValue(Xc_ID, "Xc", QObject::tr("Theoretical constant"), "1/K", true)                 
      <<GmpFemPhysicsCommon::ScalarValue(Xm_ID, "Xm", QObject::tr("Theoretical power"), "", true)                       
      <<GmpFemPhysicsCommon::ScalarValue(XKo_ID, "XKo", QObject::tr("Transient parameter"), "", true)           
      <<GmpFemPhysicsCommon::ScalarValue(alphaH_ID, "alphaH", QObject::tr("Fitting parameter for hardening 1"), "", true)       
      <<GmpFemPhysicsCommon::ScalarValue(betaH_ID, "betaH", QObject::tr("Fitting parameter for hardening 2"), "", true)         
      <<GmpFemPhysicsCommon::ScalarValue(deltaMin_ID, "deltaMin", QObject::tr("Softening parameter"), "", true);        

      // Adds Gauss attributes
      m["gaussAttributes"] = m["gaussAttributes"].value<GmpFemPhysicsCommon::ValueList>()
        << GmpFemPhysicsCommon::ScalarValue(zeta_GA_ID, "zeta", QObject::tr("Isotropic hardening variable"), " ", true, false, 2, "12.4f")
        << GmpFemPhysicsCommon::HistoryValue(zetaOld_GA_ID, zeta_GA_ID, 1, true)
        << GmpFemPhysicsCommon::ScalarValue(FT_GA_ID, "Ft", QObject::tr("Hardening parameter"), " ", true, false, 1, "12.4f")
        << GmpFemPhysicsCommon::ScalarValue(EcSS_GA_ID, "Ecss", QObject::tr("Steady-state creep rate"), " ", true, false, 1, "12.4f");
    }
    return &m;
  }

  // Returns the stresses according to the material behavior adopted
  virtual double fillCreepStrainRate(const GmElement* e, const GmpMechanicPoint* mp, const GmVector* coord, const GmVector& Time, double MISES) const
  {
    S_TRACE();
    assert(e); 
    // Variable definition
    double eE, eNU, eG, TEMP, Rg, A1, n1, Q1, A2, n2, Q2, B1, B2, So, Xq;
    double Xc, Xm, XKo, alphaH, betaH, deltaMin, zeta, dzeta, creepTrLim, deltaMai;
    double dTime, creepRate, creepSSRate, creepSSRateDCL, creepSSRateUMC, creepSSRateDGL, creepTrFunction;
    
    // Reads elastic parameters 
    eE = elasticModulus(e, coord, mp->_index);
    eNU = poissonRatio(e, coord, mp->_index);
    eG = GmpMechanicUtils::shearModulus(eE, eNU);
    
    // Reads time parameters 
    dTime = Time(0)-Time(1);
    
    // Reads steady-state creep parameters
    Rg = universalGasConstant(e, coord, mp->_index);  
    A1 = structureFactorDCL(e, coord, mp->_index);
    n1 = stressPowerDCL(e, coord, mp->_index);
    Q1 = thermalActivEnergyDCL(e, coord, mp->_index);
    A2 = structureFactorUMC(e, coord, mp->_index);
    n2 = stressPowerUMC(e, coord, mp->_index);
    Q2 = thermalActivEnergyUMC(e, coord, mp->_index);
    B1 = structureFactorDGL1(e, coord, mp->_index);
    B2 = structureFactorDGL2(e, coord, mp->_index);
    So = thresholdDevStress(e, coord, mp->_index);
    Xq = stressConstant(e, coord, mp->_index);

    // Reads the transient creep parameters
    Xc = theoreticalConstant(e, coord, mp->_index);
    Xm = theoreticalPower(e, coord, mp->_index);
    XKo = transientParameter(e, coord, mp->_index);
    alphaH = fittingParamHard1(e, coord, mp->_index);
    betaH = fittingParamHard2(e, coord, mp->_index);
    deltaMin = softParameter(e, coord, mp->_index);

    // Reads the current value of the internal isotropic hardening variable
    GmGaussAccessor* zetaAcc = gaussAttrAc(zeta_GA_ID);
    GmGaussAccessor* zetaOldAcc = gaussAttrAc(zetaOld_GA_ID);
    zeta = zetaOldAcc->scalarValueAt(e, mp->_index, coord);

    // No creep for Misses equal to zero
    if (MISES == 0.0)
    {
      return 0.0;
    }

    // Reads the temperature at integration point
    if (nodeAttrAc(T_NA_ID)!= NULL)
    {
      TEMP = fillTemperatureFromNodalAttr(e, coord);
    }
    else
    {
      TEMP = temperatureProperty(e, coord, mp->_index);
    }    
    
    // Calculations
    creepSSRateDCL = A1 * exp(-1.0 * Q1/(Rg*TEMP)) * pow(MISES/eG , n1);
    creepSSRateUMC = A2 * exp(-1.0 * Q2/(Rg*TEMP)) * pow(MISES/eG , n2);
    creepSSRateDGL = 0.0;
    
    if (MISES > So)
    {
      creepSSRateDGL = (B1 * exp(-1.0 * Q1/(Rg * TEMP)) + B2 * exp(-1.0 * Q2/(Rg * TEMP))) * sinh(Xq * (MISES-So)/eG);
    }
    
    // Stady-state creep rate
    creepSSRate = creepSSRateDCL + creepSSRateUMC + creepSSRateDGL;     

    // Transient creep limit
    creepTrLim = XKo * pow((MISES / eG), Xm) * exp(Xc * TEMP);

    // Hardening parameter
    deltaMai = alphaH + betaH * log10(MISES / eG);

    // Evaluates the transient function Ft
    if (zeta <= creepTrLim)
    {
      creepTrFunction = exp(deltaMai*pow((1.0 - zeta / creepTrLim), 2.0));
    }
    else // (zeta > creepTrLim)
    {
      creepTrFunction = exp(-1.0*deltaMin*pow((1.0 - zeta / creepTrLim), 2.0));
    }
 
    // Evaluates and saves the internal isotropic hardening variable for the next increment
    dzeta = (creepTrFunction - 1.0)*(creepSSRate)*dTime;
    zeta = zeta + dzeta;            
    zetaAcc->setScalarValue(e, mp->_index, zeta);
    // Saves the transient function Ft 
    gaussAttrAc(FT_GA_ID)->setScalarValue(e, mp->_index, creepTrFunction);
    // Saves the steady-state creep rate
    gaussAttrAc(EcSS_GA_ID)->setScalarValue(e, mp->_index, creepSSRate);

    // Calculates the creep strain rate
    creepRate = creepTrFunction * creepSSRate;
        
    return creepRate;
  }

  double structureFactorDGL1(const GmElement* e, const GmVector* coord, int ip) const
  {
    return propertyAc(B1_ID)->scalarValueAt(e, coord, ip);
  }
  
  double structureFactorDGL2(const GmElement* e, const GmVector* coord, int ip) const
  {
    return propertyAc(B2_ID)->scalarValueAt(e, coord, ip);
  }

  double stressConstant(const GmElement* e, const GmVector* coord, int ip) const
  {
    return propertyAc(Xq_ID)->scalarValueAt(e, coord, ip);
  }

  virtual double theoreticalConstant(const GmElement* e, const GmVector* coord, int ip) const
  {
    return propertyAc(Xc_ID)->scalarValueAt(e, coord, ip);
  }
  
  double theoreticalPower(const GmElement* e, const GmVector* coord, int ip) const
  {
    return propertyAc(Xm_ID)->scalarValueAt(e, coord, ip);
  }
  
  double transientParameter(const GmElement* e, const GmVector* coord, int ip) const
  {
    return propertyAc(XKo_ID)->scalarValueAt(e, coord, ip);
  }
  
  double fittingParamHard1(const GmElement* e, const GmVector* coord, int ip) const
  {
    return propertyAc(alphaH_ID)->scalarValueAt(e, coord, ip);
  }
  
  double fittingParamHard2(const GmElement* e, const GmVector* coord, int ip) const
  {
    return propertyAc(betaH_ID)->scalarValueAt(e, coord, ip);
  }
  
  double softParameter(const GmElement* e, const GmVector* coord, int ip) const
  {
    return propertyAc(deltaMin_ID)->scalarValueAt(e, coord, ip);
  }
    
};

#endif