mps.h

Go to the documentation of this file.

/*
 *
 *                 #####    #####   ######  ######  ###   ###
 *               ##   ##  ##   ##  ##      ##      ## ### ##
 *              ##   ##  ##   ##  ####    ####    ##  #  ##
 *             ##   ##  ##   ##  ##      ##      ##     ##
 *            ##   ##  ##   ##  ##      ##      ##     ##
 *            #####    #####   ##      ######  ##     ##
 *
 *
 *             OOFEM : Object Oriented Finite Element Code
 *
 *               Copyright (C) 1993 - 2025   Borek Patzak
 *
 *
 *
 *       Czech Technical University, Faculty of Civil Engineering,
 *   Department of Structural Mechanics, 166 29 Prague, Czech Republic
 *
 *  This library is free software; you can redistribute it and/or
 *  modify it under the terms of the GNU Lesser General Public
 *  License as published by the Free Software Foundation; either
 *  version 2.1 of the License, or (at your option) any later version.
 *
 *  This program 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
 *  Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public
 *  License along with this library; if not, write to the Free Software
 *  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
 */
 
//   ***********************************************************************************
//   *** CLASS Material Model based on Microprestress Solidification Theory & Status ***
//   ***********************************************************************************
 
#ifndef mps_h
#define mps_h
 
#include "kelvinChSolM.h"


#define _IFT_MPSMaterial_Name "mps"
#define _IFT_MPSMaterial_mode "mode"
#define _IFT_MPSMaterial_coupledanalysistype "coupledanalysistype"
#define _IFT_MPSMaterial_fc "fc"
#define _IFT_MPSMaterial_cc "cc"
#define _IFT_MPSMaterial_wc "w/c"
#define _IFT_MPSMaterial_wcr "wc" //for Python bindings
#define _IFT_MPSMaterial_ac "a/c"
#define _IFT_MPSMaterial_acr "ac" //for Python bindings
#define _IFT_MPSMaterial_q1 "q1"
#define _IFT_MPSMaterial_q2 "q2"
#define _IFT_MPSMaterial_q3 "q3"
#define _IFT_MPSMaterial_q4 "q4"
#define _IFT_MPSMaterial_lambda0 "lambda0"
// in input record t0 replaced by relMatAge
//#define _IFT_MPSMaterial_t0 "t0"
#define _IFT_MPSMaterial_ksh "ksh"
#define _IFT_MPSMaterial_qetor "qetor"
#define _IFT_MPSMaterial_qrtor "qrtor"
#define _IFT_MPSMaterial_qstor "qstor"
#define _IFT_MPSMaterial_alphae "alphae"
#define _IFT_MPSMaterial_alphar "alphar"
#define _IFT_MPSMaterial_alphas "alphas"
#define _IFT_MPSMaterial_mus "mus"
#define _IFT_MPSMaterial_ktm "ktm"
#define _IFT_MPSMaterial_ktc "ktc"
#define _IFT_MPSMaterial_khc "khc"
#define _IFT_MPSMaterial_stiffnessfactor "stiffnessfactor"
#define _IFT_MPSMaterial_p "p"
#define _IFT_MPSMaterial_p_tilde "p_tilde"
#define _IFT_MPSMaterial_k3 "k3"
#define _IFT_MPSMaterial_sh_a "sh_a"
#define _IFT_MPSMaterial_sh_hC "sh_hc"
#define _IFT_MPSMaterial_sh_n "sh_n"
#define _IFT_MPSMaterial_factor_ksh_h "factor_ksh_h"
#define _IFT_MPSMaterial_factor_ksh_fh "factor_ksh_f(h)"
#define _IFT_MPSMaterial_timedependent_ksh "timedependent_ksh"
#define _IFT_MPSMaterial_alpha_as "alpha_as"
#define _IFT_MPSMaterial_eps_cas0 "eps_cas0"
#define _IFT_MPSMaterial_B4_eps_au_infty "b4_eps_au_infty"
#define _IFT_MPSMaterial_B4_tau_au "b4_tau_au"
#define _IFT_MPSMaterial_B4_alpha "b4_alpha"
#define _IFT_MPSMaterial_B4_r_t "b4_r_t"
#define _IFT_MPSMaterial_B4_cem_type "b4_cem_type"
#define _IFT_MPSMaterial_temperInCelsius "temperincelsius"
#define _IFT_MPSMaterial_hydrationTimescaleTF "hydrationtimescaletf"
#define _IFT_MPSMaterial_autoShrinkageTF "autoshrinkagetf"
#define _IFT_MPSMaterial_tau_nano "tau_nano"
 
namespace oofem {
class MPSMaterialStatus : public KelvinChainSolidMaterialStatus
{
protected:
    double hum = -1.;
    double hum_increment = -1.;
    double h_min = 1.;
    double T = -1.;
    double T_increment = -1.;
    double T_max = 0.;
    double equivalentTime = 0.;
    double equivalentTimeTemp = 0.;
    double flowTermViscosity = 0.;
    double flowTermViscosityTemp = -1.;
    bool storedEmodulusFlag = false;
    double storedEmodulus = -1.;

    double humNano = -1.;
    double humNanoTemp = -1.;
 
#ifdef keep_track_of_strains
    double dryingShrinkageStrain = 0.;
    double tempDryingShrinkageStrain = 0.;
    double autogenousShrinkageStrain = 0.;
    double tempAutogenousShrinkageStrain = 0.;
    FloatArray creepStrain;
    FloatArray creepStrainIncrement;
#endif
 
 
public:
    MPSMaterialStatus(GaussPoint *g, int nunits);
 
    void initTempStatus() override;
    void updateYourself(TimeStep *tStep) override;
 
    void saveContext(DataStream &stream, ContextMode mode) override;
    void restoreContext(DataStream &stream, ContextMode mode) override;

    double giveHum() { return hum; }
    void setHum(double src) { hum = src; }
    double giveHumIncrement() { return hum_increment; }
    void setHumIncrement(double src) { hum_increment = src; }
 

    double giveHumNano() { return humNano; }
    void setHumNanoTemp(double src) { humNanoTemp = src; }

    double giveHumMin() { return h_min; }
    void setHumMin(double src) { h_min = src; }
 

    double giveT() { return T; }
    void setT(double src) { T = src; }

    double giveTIncrement() { return T_increment; }
    void setTIncrement(double src) { T_increment = src; }

    double giveTmax() { return T_max; }
    void setTmax(double src) { T_max = src; }

    double giveEquivalentTime() { return equivalentTime; }
    void setEquivalentTimeTemp(double src) { equivalentTimeTemp = src; }

    double giveFlowTermViscosity() { return flowTermViscosity; }
    double giveFlowTermViscosityTemp() { return flowTermViscosityTemp; }
    void setFlowTermViscosityTemp(double src) { flowTermViscosityTemp = src; }

    void storeEmodulus(double src) { storedEmodulus = src; }
    void setEmodulusFlag(bool src) { storedEmodulusFlag = src; }
    double giveStoredEmodulus(void) { return storedEmodulus; }
    bool giveStoredEmodulusFlag(void) { return storedEmodulusFlag; }
 
#ifdef keep_track_of_strains
    void setTempDryingShrinkageStrain(double src) { tempDryingShrinkageStrain = src; }
    double giveTempDryingShrinkageStrain(void) { return tempDryingShrinkageStrain; }
    double giveDryingShrinkageStrain(void) { return dryingShrinkageStrain; }
 
    void setTempAutogenousShrinkageStrain(double src) { tempAutogenousShrinkageStrain = src; }
    double giveTempAutogenousShrinkageStrain(void) { return tempAutogenousShrinkageStrain; }
    double giveAutogenousShrinkageStrain(void) { return autogenousShrinkageStrain; }
 
    void setCreepStrainIncrement(FloatArray src) { creepStrainIncrement  = std::move(src); }
    const FloatArray &giveCreepStrain() const { return creepStrain; }
#endif
 
    // definition
    const char *giveClassName() const override { return "MPSMaterialStatus"; }
};
 

class MPSMaterial : public KelvinChainSolidMaterial
{
protected:
    double t0 = 0.;
    double q1 = 0., q2 = 0., q3 = 0., q4 = 0.;
    double lambda0 = 0.;
 
 
    enum coupledAnalysisType { Basic, MPS_full, MPS_humidity, MPS_temperature } CoupledAnalysis;
 
    mutable double EspringVal = 0.; // elastic modulus of the aging spring (first member of Kelvin chain if retardation spectrum is used)

    //double w_h, n, a; //constant (obtained from experiments) A [Pedersen, 1990]
 
    // MPS theory parameters
    double kSh = 0.;
    double muS = 0., k3 = 0.;
    double kTm = 0.;
    double kTc = 0.;
    double khc = 0.;
    // commented - unused parameter? PH
    //double ct = 0.;
    double roomTemperature = 0.;
    double QEtoR = 0., QRtoR = 0., QStoR = 0.; //[K]
    double alphaE = 0., alphaR = 0., alphaS = 0.; //[-]
    double p = 0.;
    double sh_a = 0., sh_hC = 0., sh_n = 0.;
    bool timeDependent_ksh = false;
    FloatArray ksh_h, ksh_fh;
    double eps_cas0 = 0.;
    double b4_eps_au_infty = 0., b4_tau_au = 0., b4_alpha = 0., b4_r_t = 0.;

    double stiffnessFactor = 0.;

    double temperScaleDifference = 0.;

    int hydrationTimescaleTF;

    int autoShrinkageTF;

    double tau_nano = 0.;
 
 
 
public:
    MPSMaterial(int n, Domain *d) : KelvinChainSolidMaterial(n, d) { }
 
    const char *giveInputRecordName() const override { return _IFT_MPSMaterial_Name; }
    const char *giveClassName() const override { return "MPSMaterial"; }
 
    void initializeFrom(InputRecord &ir) override;
 
    int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep) override;
 
    void giveRealStressVector(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep) const override;
    //void updateYourself(GaussPoint *gp, TimeStep *tStep) override;
 
    void giveShrinkageStrainVector(FloatArray &answer, GaussPoint *gp, TimeStep *tStep, ValueModeType mode) const override;

    double computeCreepFunction(double t, double t_prime, GaussPoint *gp, TimeStep *tStep) const override;
 
    double giveEquivalentTime(GaussPoint *gp, TimeStep *tStep) const override;
 
    std::unique_ptr<MaterialStatus> CreateStatus(GaussPoint *gp) const override;
 
protected:
    void predictParametersFrom(double, double, double, double);
 
    void computeCharTimes() override;

    FloatArray computeCharCoefficients(double tPrime, GaussPoint *gp, TimeStep *tStep) const override;
 
    double giveEModulus(GaussPoint *gp, TimeStep *tStep) const override;
 
    double computeSolidifiedVolume(GaussPoint *gp, TimeStep *tStep) const override;
 
    double computeBetaMu(GaussPoint *gp, TimeStep *tStep, int Mu) const override;
    double computeLambdaMu(GaussPoint *gp, TimeStep *tStep, int Mu) const override;

    double computeFlowTermViscosity(GaussPoint *gp, TimeStep *tStep) const;

    double giveInitViscosity(TimeStep *tStep) const;
 
    void giveEigenStrainVector(FloatArray &answer, GaussPoint *gp, TimeStep *tStep, ValueModeType mode) const override;
 
    bool hasIncrementalShrinkageFormulation() const override { return true; }

    void computePointShrinkageStrainVector(FloatArray &answer, GaussPoint *gp, TimeStep *tStep) const;

    void computeFibAutogenousShrinkageStrainVector(FloatArray &answer, GaussPoint *gp, TimeStep *tStep) const;

    void computeB4AutogenousShrinkageStrainVector(FloatArray &answer, GaussPoint *gp, TimeStep *tStep) const;

    void computeAutogenousShrinkageDefinedByTF(FloatArray &answer, GaussPoint *gp, TimeStep *tStep) const;
 
    //double inverse_sorption_isotherm(double w) const;

    double giveHumidity(GaussPoint *gp, TimeStep *tStep, int option) const;

    double giveHumidityNano(GaussPoint *gp, TimeStep *tStep, int option) const;

    double giveTemperature(GaussPoint *gp, TimeStep *tStep, int option) const;

    double computePsiR(GaussPoint *gp, TimeStep *tStep, int option) const;

    double computePsiS(GaussPoint *gp, TimeStep *tStep) const;

    double computePsiE(GaussPoint *gp, TimeStep *tStep) const;

    double computeEquivalentTime(GaussPoint *gp, TimeStep *tStep, int option) const;
 
    friend class RankineMPSmat;
};
} // end namespace oofem
#endif // mps_h