OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
oofem::MPSDamMaterial Class Reference

This class extends the material model based on MPS theory (microprestress-solidification) for concrete creep and shrinkage by a simple isotropic damage model to take into account cracking in tension. More...

#include <mpsdammat.h>

+ Inheritance diagram for oofem::MPSDamMaterial:
+ Collaboration diagram for oofem::MPSDamMaterial:

Public Member Functions

 MPSDamMaterial (int n, Domain *d)
 
virtual ~MPSDamMaterial ()
 
virtual int hasNonLinearBehaviour ()
 Returns nonzero if receiver is non linear. More...
 
virtual int hasMaterialModeCapability (MaterialMode mode)
 Tests if material supports material mode. More...
 
virtual const char * giveInputRecordName () const
 
virtual const char * giveClassName () const
 
virtual IRResultType initializeFrom (InputRecord *ir)
 Initializes receiver according to object description stored in input record. More...
 
virtual int giveIPValue (FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
 Returns the integration point corresponding value in Reduced form. More...
 
virtual void giveRealStressVector (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)
 Computes the real stress vector for given total strain and integration point. More...
 
double givee0 (GaussPoint *gp)
 
double givegf (GaussPoint *gp)
 
void initDamaged (double kappa, FloatArray &totalStrainVector, GaussPoint *gp, TimeStep *tStep)
 Abstract service allowing to perform some initialization, when damage first appear. More...
 
void initDamagedFib (GaussPoint *gp, TimeStep *tStep)
 
void computeDamage (double &omega, double kappa, GaussPoint *gp)
 Computes the value of damage omega, based on given value of equivalent strain. More...
 
void computeDamageForCohesiveCrack (double &omega, double kappa, GaussPoint *gp)
 computes the value of damage parameter omega, based on a given value of equivalent strain, using iterations to achieve objectivity, based on the crack band concept (effective element size used) More...
 
virtual MaterialStatusCreateStatus (GaussPoint *gp) const
 Creates new copy of associated status and inserts it into given integration point. More...
 
virtual double computeTensileStrength (double equivalentTime)
 
virtual double computeFractureEnergy (double equivalentTime)
 
virtual void give3dMaterialStiffnessMatrix (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Computes full 3d material stiffness matrix at given integration point, time, respecting load history in integration point. More...
 
virtual void givePlaneStressStiffMtrx (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Method for computing plane stress stiffness matrix of receiver. More...
 
virtual void givePlaneStrainStiffMtrx (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Method for computing plane strain stiffness matrix of receiver. More...
 
virtual void give1dStressStiffMtrx (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Method for computing 1d stiffness matrix of receiver. More...
 
- Public Member Functions inherited from oofem::MPSMaterial
 MPSMaterial (int n, Domain *d)
 
virtual ~MPSMaterial ()
 
virtual void giveShrinkageStrainVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep, ValueModeType mode)
 Computes, for the given integration point, the strain vector induced by stress-independent shrinkage. More...
 
virtual double computeCreepFunction (double t, double t_prime, GaussPoint *gp, TimeStep *tStep)
 Evaluation of the basic creep compliance function - can be used to compute elastic modulus in derived damage material. More...
 
- Public Member Functions inherited from oofem::KelvinChainSolidMaterial
 KelvinChainSolidMaterial (int n, Domain *d)
 
virtual ~KelvinChainSolidMaterial ()
 
void computeHiddenVars (GaussPoint *gp, TimeStep *tStep)
 
- Public Member Functions inherited from oofem::RheoChainMaterial
 RheoChainMaterial (int n, Domain *d)
 
virtual ~RheoChainMaterial ()
 
virtual void giveRealStressVector_3d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 Default implementation relies on giveRealStressVector for second Piola-Kirchoff stress. More...
 
virtual void giveRealStressVector_PlaneStrain (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 Default implementation relies on giveRealStressVector_3d. More...
 
virtual void giveRealStressVector_PlaneStress (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 Default implementation relies on giveRealStressVector_StressControl. More...
 
virtual void giveRealStressVector_1d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 Default implementation relies on giveRealStressVector_StressControl. More...
 
virtual void giveRealStressVector_2dBeamLayer (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 Default implementation relies on giveRealStressVector_StressControl. More...
 
virtual void giveRealStressVector_PlateLayer (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 Default implementation relies on giveRealStressVector_StressControl. More...
 
virtual void giveThermalDilatationVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
 Returns a vector of coefficients of thermal dilatation in direction of each material principal (local) axis. More...
 
virtual int hasCastingTimeSupport ()
 Tests if material supports casting time. More...
 
virtual contextIOResultType saveIPContext (DataStream &stream, ContextMode mode, GaussPoint *gp)
 Stores integration point state to output stream. More...
 
virtual contextIOResultType restoreIPContext (DataStream &stream, ContextMode mode, GaussPoint *gp)
 Reads integration point state to output stream. More...
 
virtual void give2dLatticeStiffMtrx (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 Method for computing 2d lattice stiffness matrix of receiver. More...
 
virtual void give3dLatticeStiffMtrx (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 Method for computing 3d lattice stiffness matrix of receiver. More...
 
virtual void computeStressIndependentStrainVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep, ValueModeType mode)
 Computes reduced strain vector in given integration point, generated by internal processes in material, which are independent on loading in particular integration point. More...
 
double giveAlphaOne () const
 
double giveAlphaTwo () const
 
virtual bool isActivated (TimeStep *tStep)
 
- Public Member Functions inherited from oofem::StructuralMaterial
 StructuralMaterial (int n, Domain *d)
 Constructor. More...
 
virtual ~StructuralMaterial ()
 Destructor. More...
 
virtual void giveInputRecord (DynamicInputRecord &input)
 Setups the input record string of receiver. More...
 
virtual void giveStiffnessMatrix (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Computes the stiffness matrix for giveRealStressVector of receiver in given integration point, respecting its history. More...
 
virtual void giveRealStressVector_StressControl (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, const IntArray &strainControl, TimeStep *tStep)
 Iteratively calls giveRealStressVector_3d to find the stress controlled equal to zero· More...
 
virtual void giveRealStressVector_ShellStressControl (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, const IntArray &strainControl, TimeStep *tStep)
 
virtual void giveRealStressVector_Warping (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 Default implementation relies on giveRealStressVector_StressControl. More...
 
virtual void giveRealStressVector_Fiber (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 Default implementation relies on giveRealStressVector_StressControl. More...
 
virtual void giveRealStressVector_Lattice2d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 
virtual void giveRealStressVector_Lattice3d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 
virtual void giveRealStressVector_2dPlateSubSoil (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 Default implementation is not provided. More...
 
virtual void giveRealStressVector_3dBeamSubSoil (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
 
virtual void giveEshelbyStressVector_PlaneStrain (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedF, TimeStep *tStep)
 Prototype for computation of Eshelby stress. More...
 
void give_dPdF_from (const FloatMatrix &dSdE, FloatMatrix &answer, GaussPoint *gp)
 
void convert_dSdE_2_dPdF (FloatMatrix &answer, const FloatMatrix &dSdE, const FloatArray &S, const FloatArray &F, MaterialMode matMode)
 
double giveReferenceTemperature ()
 Returns the reference temperature of receiver. More...
 
virtual void computeStressIndependentStrainVector_3d (FloatArray &answer, GaussPoint *gp, TimeStep *tStep, ValueModeType mode)
 
virtual void give3dMaterialStiffnessMatrix_dPdF (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 
virtual void give3dMaterialStiffnessMatrix_dCde (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 
void giveStressDependentPartOfStrainVector (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrainVector, TimeStep *tStep, ValueModeType mode)
 Method for subtracting from reduced space strain vector its stress-independent parts (caused by temperature, shrinkage, creep and possibly by other phenomena). More...
 
void giveStressDependentPartOfStrainVector_3d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrainVector, TimeStep *tStep, ValueModeType mode)
 
virtual int setIPValue (const FloatArray &value, GaussPoint *gp, InternalStateType type)
 Sets the value of a certain variable at a given integration point to the given value. More...
 
virtual void give2dBeamLayerStiffMtrx (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 Method for computing 2d beam layer stiffness matrix of receiver. More...
 
virtual void givePlateLayerStiffMtrx (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 Method for computing 2d plate layer stiffness matrix of receiver. More...
 
virtual void giveFiberStiffMtrx (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 Method for computing 1d fiber stiffness matrix of receiver. More...
 
virtual void give2dPlateSubSoilStiffMtrx (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 Method for computing stiffness matrix of plate subsoil model. More...
 
virtual void give3dBeamSubSoilStiffMtrx (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 Method for computing stiffness matrix of beam3d subsoil model. More...
 
virtual void giveFirstPKStressVector_3d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedF, TimeStep *tStep)
 Default implementation relies on giveRealStressVector for second Piola-Kirchoff stress. More...
 
virtual void giveFirstPKStressVector_PlaneStrain (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedF, TimeStep *tStep)
 Default implementation relies on giveFirstPKStressVector_3d. More...
 
virtual void giveFirstPKStressVector_PlaneStress (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedF, TimeStep *tStep)
 Default implementation relies on giveFirstPKStressVector_3d. More...
 
virtual void giveFirstPKStressVector_1d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedF, TimeStep *tStep)
 Default implementation relies on giveFirstPKStressVector_3d. More...
 
virtual void giveCauchyStressVector_3d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedF, TimeStep *tStep)
 
virtual void giveCauchyStressVector_PlaneStrain (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedF, TimeStep *tStep)
 
virtual void giveCauchyStressVector_PlaneStress (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedF, TimeStep *tStep)
 
virtual void giveCauchyStressVector_1d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedF, TimeStep *tStep)
 
virtual void givePlaneStressStiffMtrx_dPdF (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 
virtual void givePlaneStressStiffMtrx_dCde (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 
virtual void givePlaneStrainStiffMtrx_dPdF (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 
virtual void givePlaneStrainStiffMtrx_dCde (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 
virtual void give1dStressStiffMtrx_dPdF (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 
virtual void give1dStressStiffMtrx_dCde (FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
 
- Public Member Functions inherited from oofem::Material
 Material (int n, Domain *d)
 Constructor. More...
 
virtual ~Material ()
 Destructor. More...
 
virtual bool isCharacteristicMtrxSymmetric (MatResponseMode rMode)
 Returns true if stiffness matrix of receiver is symmetric Default implementation returns true. More...
 
virtual double give (int aProperty, GaussPoint *gp)
 Returns the value of material property 'aProperty'. More...
 
virtual bool hasProperty (int aProperty, GaussPoint *gp)
 Returns true if 'aProperty' exists on material. More...
 
virtual void modifyProperty (int aProperty, double value, GaussPoint *gp)
 Modify 'aProperty', which already exists on material. More...
 
double giveCastingTime ()
 
virtual void printYourself ()
 Prints receiver state on stdout. Useful for debugging. More...
 
virtual int checkConsistency ()
 Allows programmer to test some internal data, before computation begins. More...
 
virtual int initMaterial (Element *element)
 Optional function to call specific procedures when initializing a material. More...
 
virtual MaterialStatusgiveStatus (GaussPoint *gp) const
 Returns material status of receiver in given integration point. More...
 
virtual int packUnknowns (DataStream &buff, TimeStep *tStep, GaussPoint *ip)
 Pack all necessary data of integration point (according to element parallel_mode) into given communication buffer. More...
 
virtual int unpackAndUpdateUnknowns (DataStream &buff, TimeStep *tStep, GaussPoint *ip)
 Unpack and updates all necessary data of given integration point (according to element parallel_mode) into given communication buffer. More...
 
virtual int estimatePackSize (DataStream &buff, GaussPoint *ip)
 Estimates the necessary pack size to hold all packed data of receiver. More...
 
virtual double predictRelativeComputationalCost (GaussPoint *gp)
 Returns the weight representing relative computational cost of receiver The reference material model is linear isotropic material - its weight is set to 1.0 The other material models should compare to this reference model. More...
 
virtual double predictRelativeRedistributionCost (GaussPoint *gp)
 Returns the relative redistribution cost of the receiver. More...
 
virtual void initTempStatus (GaussPoint *gp)
 Initializes temporary variables stored in integration point status at the beginning of new time step. More...
 
- Public Member Functions inherited from oofem::FEMComponent
 FEMComponent (int n, Domain *d)
 Regular constructor, creates component with given number and belonging to given domain. More...
 
virtual ~FEMComponent ()
 Virtual destructor. More...
 
DomaingiveDomain () const
 
virtual void setDomain (Domain *d)
 Sets associated Domain. More...
 
int giveNumber () const
 
void setNumber (int num)
 Sets number of receiver. More...
 
virtual void updateLocalNumbering (EntityRenumberingFunctor &f)
 Local renumbering support. More...
 
virtual contextIOResultType saveContext (DataStream &stream, ContextMode mode, void *obj=NULL)
 Stores receiver state to output stream. More...
 
virtual contextIOResultType restoreContext (DataStream &stream, ContextMode mode, void *obj=NULL)
 Restores the receiver state previously written in stream. More...
 
virtual void printOutputAt (FILE *file, TimeStep *tStep)
 Prints output of receiver to stream, for given time step. More...
 
virtual InterfacegiveInterface (InterfaceType t)
 Interface requesting service. More...
 
std::string errorInfo (const char *func) const
 Returns string for prepending output (used by error reporting macros). More...
 

Protected Types

enum  SofteningType { ST_Exponential_Cohesive_Crack, ST_Linear_Cohesive_Crack, ST_Disable_Damage }
 Type characterizing the formula for the damage law. More...
 
- Protected Types inherited from oofem::MPSMaterial
enum  coupledAnalysisType { Basic, MPS_full, MPS_humidity, MPS_temperature }
 

Protected Attributes

bool timeDepFracturing
 
double fib_s
 
double fib_fcm28
 
bool isotropic
 
double E
 dummy Young's modulus More...
 
double maxOmega
 Maximum limit on omega. The purpose is elimination of a too compliant material which may cause convergence problems. Set to something like 0.99 if needed. More...
 
double ft
 Equivalent strain at stress peak (or a similar parameter). More...
 
double const_gf
 Determines the softening -> corresponds to the initial fracture energy. More...
 
int checkSnapBack
 Check possible snap back flag. More...
 
SofteningType softType
 Parameter specifying the type of softening (damage law). More...
 
ElementCharSizeMethod ecsMethod
 Method used for evaluation of characteristic element size. More...
 
double ft28
 28-day value of tensile strength. Used only with "timedepfracturing" More...
 
double gf28
 28-day value of fracture energy. Used only with "timedepfracturing" More...
 
- Protected Attributes inherited from oofem::MPSMaterial
double t0
 age when temperature or humidity starts to change More...
 
double q1
 compliances of the B3 model More...
 
double q2
 
double q3
 
double q4
 
double lambda0
 constant equal to one day in time units of analysis (eg. 86400 if the analysis runs in seconds) More...
 
enum oofem::MPSMaterial::coupledAnalysisType CoupledAnalysis
 
double EspringVal
 
double kSh
 additional parameters for sorption isotherm (used to compute relative humidity from water content) More...
 
double muS
 fluidity parameter used in viscosity evolution equation More...
 
double k3
 
double kTm
 kTm replaces ln(h) on RHS of the differential equation describing evolution of MPS More...
 
double kTc
 parameter reducing creep effects of thermal cycling (replaces kTm in such case) More...
 
double ct
 parameter reducing creep effects of thermal cycling More...
 
double roomTemperature
 reference room temperature for MPS algorithm [K] More...
 
double QEtoR
 activation energies More...
 
double QRtoR
 
double QStoR
 
double alphaE
 parameters that control the effect of humidity on rates of hydration, creep and microprestress relaxation More...
 
double alphaR
 
double alphaS
 
double p
 exponent in the microprestress/viscosity governing equation More...
 
double sh_a
 parameters for nonlinear shrinkage function More...
 
double sh_hC
 
double sh_n
 
double eps_cas0
 parameter for autogenous shrinkage according to fib MC 2010 More...
 
double b4_eps_au_infty
 parameters for autogenous shrinkage according to B4 model More...
 
double b4_tau_au
 
double b4_alpha
 
double b4_r_t
 
double stiffnessFactor
 scaling factor 1. for Pa, 1.e6 for MPa - only for empirical formulas - q1-q4 and ft and gf More...
 
double temperScaleDifference
 0 for Kelvin, 273.15 for Celsius More...
 
- Protected Attributes inherited from oofem::RheoChainMaterial
double talpha
 thermal dilatation coeff. More...
 
int nUnits
 Number of (Maxwell or Kelvin) units in the rheologic chain. More...
 
double relMatAge
 Physical age of the material at simulation time = 0. More...
 
bool lattice
 
double nu
 Poisson's ratio (assumed to be constant, unaffected by creep). More...
 
double alphaOne
 Parameters for the lattice model. More...
 
double alphaTwo
 
double EparValTime
 Time for which the partial moduli of individual units have been evaluated. More...
 
double begOfTimeOfInterest
 Time from which the model should give a good approximation. Optional field. Default value is 0.1 [day]. More...
 
double endOfTimeOfInterest
 Time (age???) up to which the model should give a good approximation. More...
 
LinearElasticMateriallinearElasticMaterial
 Associated linearElasticMaterial, with E = 1. More...
 
FloatArray EparVal
 Partial moduli of individual units. More...
 
FloatArray charTimes
 Characteristic times of individual units (relaxation or retardation times). More...
 
FloatArray discreteTimeScale
 Times at which the errors are evaluated if the least-square method is used. More...
 
double timeFactor
 Scaling factor transforming the simulation time units into days (gives the number of simulation time units in one day, e.g. More...
 
int preCastingTimeMat
 Stiffness at time less than casting time - optional parameter, negative by default. More...
 
- Protected Attributes inherited from oofem::StructuralMaterial
double referenceTemperature
 Reference temperature (temperature, when material has been built into structure). More...
 
- Protected Attributes inherited from oofem::Material
Dictionary propertyDictionary
 Property dictionary. More...
 
double castingTime
 Casting time. More...
 
- Protected Attributes inherited from oofem::FEMComponent
int number
 Component number. More...
 
Domaindomain
 Link to domain object, useful for communicating with other FEM components. More...
 

Additional Inherited Members

- Static Public Member Functions inherited from oofem::StructuralMaterial
static int giveSymVI (int ind1, int ind2)
 
static int giveVI (int ind1, int ind2)
 
static int giveVoigtVectorMask (IntArray &answer, MaterialMode mmode)
 Returns a mask of the vector indicies corresponding to components in a general (non-symmetric) second order tensor of some stress/strain/deformation measure that performes work. More...
 
static int giveVoigtSymVectorMask (IntArray &answer, MaterialMode mmode)
 The same as giveVoigtVectorMask but returns a mask corresponding to a symmetric second order tensor. More...
 
static void giveInvertedVoigtVectorMask (IntArray &answer, MaterialMode mmode)
 Gives the inverted version of giveVoigtVectorMask. More...
 
static int giveSizeOfVoigtVector (MaterialMode mmode)
 Returns the size of reduced stress/strain vector according to given mode. More...
 
static int giveSizeOfVoigtSymVector (MaterialMode mmode)
 Returns the size of symmetric part of a reduced stress/strain vector according to given mode. More...
 
static void giveFullVectorForm (FloatArray &answer, const FloatArray &strainVector, MaterialMode matMode)
 Converts the reduced symmetric Voigt vector (2nd order tensor) to full form. More...
 
static void giveFullVectorFormF (FloatArray &answer, const FloatArray &strainVector, MaterialMode matMode)
 Converts the reduced deformation gradient Voigt vector (2nd order tensor). More...
 
static void giveFullSymVectorForm (FloatArray &answer, const FloatArray &vec, MaterialMode matMode)
 Converts the reduced unsymmetric Voigt vector (2nd order tensor) to full form. More...
 
static void giveReducedVectorForm (FloatArray &answer, const FloatArray &vec, MaterialMode matMode)
 Converts the full symmetric Voigt vector (2nd order tensor) to reduced form. More...
 
static void giveReducedSymVectorForm (FloatArray &answer, const FloatArray &vec, MaterialMode matMode)
 Converts the full unsymmetric Voigt vector (2nd order tensor) to reduced form. More...
 
static void giveFullSymMatrixForm (FloatMatrix &answer, const FloatMatrix &red, MaterialMode matMode)
 Converts the full unsymmetric Voigt matrix (4th order tensor) to reduced form. More...
 
static void giveReducedMatrixForm (FloatMatrix &answer, const FloatMatrix &full, MaterialMode matMode)
 Converts the full symmetric Voigt matrix (4th order tensor) to reduced form. More...
 
static void giveReducedSymMatrixForm (FloatMatrix &answer, const FloatMatrix &full, MaterialMode matMode)
 Converts the full unsymmetric Voigt matrix (4th order tensor) to reduced form. More...
 
static void transformStrainVectorTo (FloatArray &answer, const FloatMatrix &base, const FloatArray &strainVector, bool transpose=false)
 Transforms 3d strain vector into another coordinate system. More...
 
static void transformStressVectorTo (FloatArray &answer, const FloatMatrix &base, const FloatArray &stressVector, bool transpose=false)
 Transforms 3d stress vector into another coordinate system. More...
 
static double computeVonMisesStress (const FloatArray *currentStress)
 Computes equivalent of von Mises stress. More...
 
static void giveStrainVectorTranformationMtrx (FloatMatrix &answer, const FloatMatrix &base, bool transpose=false)
 Computes 3d strain vector transformation matrix from standard vector transformation matrix. More...
 
static void give2DStrainVectorTranformationMtrx (FloatMatrix &answer, const FloatMatrix &base, bool transpose=false)
 Computes 2d strain vector transformation matrix from standard vector transformation matrix. More...
 
static void giveStressVectorTranformationMtrx (FloatMatrix &answer, const FloatMatrix &base, bool transpose=false)
 Computes 3d stress vector transformation matrix from standard vector transformation matrix. More...
 
static void givePlaneStressVectorTranformationMtrx (FloatMatrix &answer, const FloatMatrix &base, bool transpose=false)
 Computes 2d stress vector transformation matrix from standard vector transformation matrix. More...
 
static void sortPrincDirAndValCloseTo (FloatArray *pVal, FloatMatrix *pDir, FloatMatrix *toPDir)
 Method for sorting newly computed principal values (pVal) and corresponding principal directions (pDir) to be closed to some (often previous) principal directions (toPDir). More...
 
static void computePrincipalValues (FloatArray &answer, const FloatArray &s, stressStrainPrincMode mode)
 Common functions for convenience. More...
 
static void computePrincipalValDir (FloatArray &answer, FloatMatrix &dir, const FloatArray &s, stressStrainPrincMode mode)
 Computes principal values and directions of stress or strain vector. More...
 
static double computeDeviatoricVolumetricSplit (FloatArray &dev, const FloatArray &s)
 Computes split of receiver into deviatoric and volumetric part. More...
 
static void computeDeviatoricVolumetricSum (FloatArray &s, const FloatArray &dev, double mean)
 
static void applyDeviatoricElasticCompliance (FloatArray &strain, const FloatArray &stress, double EModulus, double nu)
 
static void applyDeviatoricElasticCompliance (FloatArray &strain, const FloatArray &stress, double GModulus)
 
static void applyDeviatoricElasticStiffness (FloatArray &stress, const FloatArray &strain, double EModulus, double nu)
 
static void applyDeviatoricElasticStiffness (FloatArray &stress, const FloatArray &strain, double GModulus)
 
static void applyElasticStiffness (FloatArray &stress, const FloatArray &strain, double EModulus, double nu)
 
static void applyElasticCompliance (FloatArray &strain, const FloatArray &stress, double EModulus, double nu)
 
static double computeStressNorm (const FloatArray &stress)
 
static double computeFirstInvariant (const FloatArray &s)
 
static double computeSecondStressInvariant (const FloatArray &s)
 
static double computeThirdStressInvariant (const FloatArray &s)
 
static double computeFirstCoordinate (const FloatArray &s)
 
static double computeSecondCoordinate (const FloatArray &s)
 
static double computeThirdCoordinate (const FloatArray &s)
 
- Static Public Attributes inherited from oofem::StructuralMaterial
static std::vector< std::vector< int > > vIindex
 Voigt index map. More...
 
static std::vector< std::vector< int > > svIndex
 Symmetric Voigt index map. More...
 
- Protected Member Functions inherited from oofem::MPSMaterial
void predictParametersFrom (double, double, double, double)
 
virtual void computeCharTimes ()
 Evaluation of characteristic times. More...
 
virtual void computeCharCoefficients (FloatArray &answer, double, GaussPoint *gp, TimeStep *tStep)
 Evaluation of characteristic moduli of the non-aging Kelvin chain. More...
 
virtual double giveEModulus (GaussPoint *gp, TimeStep *tStep)
 Evaluation of the incremental modulus. More...
 
virtual double computeSolidifiedVolume (GaussPoint *gp, TimeStep *tStep)
 Evaluation of the relative volume of the solidified material. More...
 
virtual double computeBetaMu (GaussPoint *gp, TimeStep *tStep, int Mu)
 factors for exponential algorithm More...
 
virtual double computeLambdaMu (GaussPoint *gp, TimeStep *tStep, int Mu)
 
double computeFlowTermViscosity (GaussPoint *gp, TimeStep *tStep)
 Evaluation of the flow term viscosity. More...
 
double giveInitViscosity (TimeStep *tStep)
 Returns initial value of the flow term viscosity. More...
 
virtual void giveEigenStrainVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep, ValueModeType mode)
 Computes, for the given integration point, the strain vector induced by the stress history (typically creep strain). More...
 
virtual int hasIncrementalShrinkageFormulation ()
 If only incremental shrinkage strain formulation is provided, then total shrinkage strain must be tracked in status in order to be able to compute total value. More...
 
void computePointShrinkageStrainVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
 Evaluation of the shrinkageStrainVector - shrinkage is fully dependent on humidity rate in given GP. More...
 
void computeFibAutogenousShrinkageStrainVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
 Evaluation of the autogenousShrinkageStrainVector according to fib MC 2010 - autogenous shrinkage is fully dependent on the equivalent age at given GP. More...
 
void computeB4AutogenousShrinkageStrainVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
 Evaluation of the autogenousShrinkageStrainVector according to Bazant's B4 model. In the model the evolution depends on temperature adjusted age, here on equivalent age (additional humidity influence) More...
 
double giveHumidity (GaussPoint *gp, TimeStep *tStep, int option)
 Gives value of humidity at given GP and timestep option = 0 ... More...
 
double giveTemperature (GaussPoint *gp, TimeStep *tStep, int option)
 Gives value of temperature at given GP and timestep option = 0 ... More...
 
double computePsiR (GaussPoint *gp, TimeStep *tStep, int option)
 Evaluation of the factor transforming real time to reduced time (effect on the flow term) option = 0 ... More...
 
double computePsiS (GaussPoint *gp, TimeStep *tStep)
 Evaluation of the factor transforming real time to reduced time (effect on the evolution of microprestress) More...
 
double computePsiE (GaussPoint *gp, TimeStep *tStep)
 Evaluation of the factor transforming real time to equivalent time (effect on the solidified volume) More...
 
double computeEquivalentTime (GaussPoint *gp, TimeStep *tStep, int option)
 Computes equivalent time at given time step and GP. More...
 
- Protected Member Functions inherited from oofem::RheoChainMaterial
const FloatArraygiveDiscreteTimes ()
 
void computeDiscreteRelaxationFunction (FloatArray &answer, const FloatArray &tSteps, double t0, double tr, GaussPoint *gp, TimeStep *tSte)
 Evaluation of the relaxation function at given times. More...
 
void giveUnitComplianceMatrix (FloatMatrix &answer, GaussPoint *gp, TimeStep *tStep)
 Evaluation of elastic compliance matrix for unit Young's modulus. More...
 
void giveUnitStiffnessMatrix (FloatMatrix &answer, GaussPoint *gp, TimeStep *tStep)
 Evaluation of elastic stiffness matrix for unit Young's modulus. More...
 
virtual void updateEparModuli (double tPrime, GaussPoint *gp, TimeStep *tStep)
 Update of partial moduli of individual chain units. More...
 
double giveEparModulus (int iChain)
 Access to partial modulus of a given unit. More...
 
double giveCharTime (int) const
 Access to the characteristic time of a given unit. More...
 
virtual double giveCharTimeExponent (int i)
 Exponent to be used with the char time of a given unit, usually = 1.0. More...
 
LinearElasticMaterialgiveLinearElasticMaterial ()
 Access to the underlying linear elastic material with unit Young's modulus. More...
 
double giveEndOfTimeOfInterest ()
 Access to the time up to which the response should be accurate. More...
 
void computeTrueStressIndependentStrainVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep, ValueModeType mode)
 Computes, for the given integration point, the strain vector induced by stress-independent internal processes in the material. More...
 
- Static Protected Member Functions inherited from oofem::RheoChainMaterial
static void generateLogTimeScale (FloatArray &answer, double from, double to, int nsteps)
 Generates discrete times starting from time "from" to time "to" uniformly distributed in log time scale. More...
 

Detailed Description

This class extends the material model based on MPS theory (microprestress-solidification) for concrete creep and shrinkage by a simple isotropic damage model to take into account cracking in tension.

Most of the methods are simply copied from "idm1" and "isodamagemodel" files. Since the damage is not the main feature of this model the methods were simplified to make the code brief. The following simplifications are assumed:

  • only two softening laws are implemented: linear and exponential
  • the equivalent strain measure uses Rankine (standard) criterion (max. principal strain)
  • the damage parameters (tensile strength, fracture energy) can be either specified or they can be computed according to recommendations from fib Model Code 2010
  • on contrary to idm1, the crack vector specifies the direction of crack opening and not the crack orientation
  • snapback is checked by default
  • element length is computed using ecsMethod = ECSM_Projection

Definition at line 182 of file mpsdammat.h.

Member Enumeration Documentation

Type characterizing the formula for the damage law.

For example, linear softening can be specified with fracturing strain or crack opening.

Enumerator
ST_Exponential_Cohesive_Crack 
ST_Linear_Cohesive_Crack 
ST_Disable_Damage 

Definition at line 214 of file mpsdammat.h.

Constructor & Destructor Documentation

oofem::MPSDamMaterial::MPSDamMaterial ( int  n,
Domain d 
)
virtual oofem::MPSDamMaterial::~MPSDamMaterial ( )
inlinevirtual

Definition at line 232 of file mpsdammat.h.

Member Function Documentation

void oofem::MPSDamMaterial::computeDamage ( double &  omega,
double  kappa,
GaussPoint gp 
)

Computes the value of damage omega, based on given value of equivalent strain.

Parameters
[out]omegaContains result.
kappaEquivalent strain measure.
strainTotal strain in full form.
gpIntegration point.

Definition at line 644 of file mpsdammat.C.

References computeDamageForCohesiveCrack(), softType, and ST_Disable_Damage.

Referenced by giveRealStressVector().

void oofem::MPSDamMaterial::computeDamageForCohesiveCrack ( double &  omega,
double  kappa,
GaussPoint gp 
)

computes the value of damage parameter omega, based on a given value of equivalent strain, using iterations to achieve objectivity, based on the crack band concept (effective element size used)

Parameters
[out]omegaContains the resulting damage.
kappaEquivalent strain measure.
gpIntegration point.

Definition at line 654 of file mpsdammat.C.

References checkSnapBack, E, oofem::MPSDamMaterialStatus::giveCharLength(), givee0(), oofem::GaussPoint::giveElement(), givegf(), oofem::FEMComponent::giveNumber(), oofem::Material::giveStatus(), MPSDAMMAT_ITERATION_LIMIT, OOFEM_ERROR, OOFEM_WARNING, oofem::MPSDamMaterialStatus::setResidualTensileStrength(), softType, ST_Exponential_Cohesive_Crack, and ST_Linear_Cohesive_Crack.

Referenced by computeDamage().

double oofem::MPSDamMaterial::computeFractureEnergy ( double  equivalentTime)
virtual
double oofem::MPSDamMaterial::computeTensileStrength ( double  equivalentTime)
virtual
MaterialStatus * oofem::MPSDamMaterial::CreateStatus ( GaussPoint gp) const
virtual

Creates new copy of associated status and inserts it into given integration point.

Parameters
gpIntegration point where newly created status will be stored.
Returns
Reference to new status.

Reimplemented from oofem::MPSMaterial.

Definition at line 807 of file mpsdammat.C.

References oofem::FEMComponent::giveDomain(), and oofem::RheoChainMaterial::nUnits.

void oofem::MPSDamMaterial::give1dStressStiffMtrx ( FloatMatrix answer,
MatResponseMode  mmode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Method for computing 1d stiffness matrix of receiver.

Default implementation computes 3d stiffness matrix using give3dMaterialStiffnessMatrix and reduces it to 1d stiffness using reduce method described above. However, this reduction is quite time consuming and if it is possible, it is recommended to overload this method and provide direct method for computing particular stiffness matrix.

Parameters
answerStiffness matrix.
mmodeMaterial response mode.
gpIntegration point, which load history is used.
tStepTime step (most models are able to respond only when tStep is current time step).

Reimplemented from oofem::RheoChainMaterial.

Definition at line 873 of file mpsdammat.C.

References oofem::RheoChainMaterial::give1dStressStiffMtrx(), oofem::Material::giveStatus(), oofem::MPSDamMaterialStatus::giveTempDamage(), isotropic, maxOmega, oofem::min(), and oofem::FloatMatrix::times().

void oofem::MPSDamMaterial::give3dMaterialStiffnessMatrix ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Computes full 3d material stiffness matrix at given integration point, time, respecting load history in integration point.

Parameters
answerComputed results.
modeMaterial response mode.
gpIntegration point.
tStepTime step (most models are able to respond only when tStep is current time step).

Reimplemented from oofem::RheoChainMaterial.

Definition at line 817 of file mpsdammat.C.

References oofem::RheoChainMaterial::give3dMaterialStiffnessMatrix(), oofem::Material::giveStatus(), oofem::MPSDamMaterialStatus::giveTempDamage(), isotropic, maxOmega, oofem::min(), and oofem::FloatMatrix::times().

virtual const char* oofem::MPSDamMaterial::giveClassName ( ) const
inlinevirtual
double oofem::MPSDamMaterial::givee0 ( GaussPoint gp)
double oofem::MPSDamMaterial::givegf ( GaussPoint gp)
virtual const char* oofem::MPSDamMaterial::giveInputRecordName ( ) const
inlinevirtual
Returns
Input record name of the receiver.

Reimplemented from oofem::MPSMaterial.

Definition at line 238 of file mpsdammat.h.

References _IFT_MPSDamMaterial_Name.

int oofem::MPSDamMaterial::giveIPValue ( FloatArray answer,
GaussPoint gp,
InternalStateType  type,
TimeStep tStep 
)
virtual

Returns the integration point corresponding value in Reduced form.

Parameters
answerContain corresponding ip value, zero sized if not available.
gpIntegration point to which the value refers.
typeDetermines the type of internal variable.
tStepDetermines the time step.
Returns
Nonzero if the assignment can be done, zero if this type of variable is not supported.

Reimplemented from oofem::MPSMaterial.

Definition at line 892 of file mpsdammat.C.

References oofem::FloatArray::at(), oofem::Material::castingTime, computeTensileStrength(), ft, oofem::MPSDamMaterialStatus::giveCharLength(), oofem::MPSDamMaterialStatus::giveCrackVector(), oofem::MPSDamMaterialStatus::giveCrackWidth(), oofem::MPSDamMaterialStatus::giveDamage(), oofem::MPSMaterialStatus::giveEquivalentTime(), oofem::MPSMaterial::giveIPValue(), oofem::StructuralMaterial::giveIPValue(), oofem::MPSDamMaterialStatus::giveResidualTensileStrength(), oofem::Material::giveStatus(), oofem::MPSDamMaterialStatus::giveTempDamage(), oofem::FloatArray::resize(), and oofem::FloatArray::zero().

void oofem::MPSDamMaterial::givePlaneStrainStiffMtrx ( FloatMatrix answer,
MatResponseMode  mmode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Method for computing plane strain stiffness matrix of receiver.

Default implementation computes 3d stiffness matrix using give3dMaterialStiffnessMatrix and reduces it to plane strain stiffness using reduce method described above. However, this reduction is quite time consuming and if it is possible, it is recommended to overload this method and provide direct method for computing particular stiffness matrix. Note: as already described, if zero strain component is imposed (Plane strain, ..) this condition must be taken into account in geometrical relations, and corresponding component has to be included in reduced vector. (So plane strain conditions are $ \epsilon_z = \gamma_{xz} = \gamma_{yz} = 0 $, but relations for $ \epsilon_z$ and $\sigma_z$ are included).

Parameters
answerStiffness matrix.
mmodeMaterial response mode.
gpIntegration point, which load history is used.
tStepTime step (most models are able to respond only when tStep is current time step).

Reimplemented from oofem::RheoChainMaterial.

Definition at line 854 of file mpsdammat.C.

References oofem::RheoChainMaterial::givePlaneStrainStiffMtrx(), oofem::Material::giveStatus(), oofem::MPSDamMaterialStatus::giveTempDamage(), isotropic, maxOmega, oofem::min(), and oofem::FloatMatrix::times().

void oofem::MPSDamMaterial::givePlaneStressStiffMtrx ( FloatMatrix answer,
MatResponseMode  mmode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Method for computing plane stress stiffness matrix of receiver.

Default implementation computes 3d stiffness matrix using give3dMaterialStiffnessMatrix and reduces it to plane stress stiffness using reduce method described above. However, this reduction is quite time consuming and if it is possible, it is recommended to overload this method and provide direct method for computing particular stiffness matrix.

Parameters
answerStiffness matrix.
mmodeMaterial response mode.
gpIntegration point, which load history is used.
tStepTime step (most models are able to respond only when tStep is current time step).

Reimplemented from oofem::RheoChainMaterial.

Definition at line 836 of file mpsdammat.C.

References oofem::RheoChainMaterial::givePlaneStressStiffMtrx(), oofem::Material::giveStatus(), oofem::MPSDamMaterialStatus::giveTempDamage(), isotropic, maxOmega, oofem::min(), and oofem::FloatMatrix::times().

void oofem::MPSDamMaterial::giveRealStressVector ( FloatArray answer,
GaussPoint gp,
const FloatArray reducedStrain,
TimeStep tStep 
)
virtual

Computes the real stress vector for given total strain and integration point.

The total strain is defined as strain computed directly from displacement field at given time. The stress independent parts (temperature, eigenstrains) are subtracted in constitutive driver. The service should use previously reached equilibrium history variables. Also it should update temporary history variables in status according to newly reached state. The temporary history variables are moved into equilibrium ones after global structure equilibrium has been reached by iteration process.

Parameters
answerStress vector in reduced form. For large deformations it is treated as the second Piola-Kirchoff stress.
gpIntegration point.
reducedStrainStrain vector in reduced form. For large deformations it is treated as the Green-Lagrange strain.
tStepCurrent time step (most models are able to respond only when tStep is current time step).
Todo:
Move this to StructuralCrossSection ?

Reimplemented from oofem::MPSMaterial.

Definition at line 326 of file mpsdammat.C.

References oofem::FloatArray::add(), oofem::FloatArray::at(), oofem::FloatMatrix::at(), oofem::MPSMaterial::computeCreepFunction(), computeDamage(), oofem::StressVector::computePrincipalValDir(), oofem::StructuralMaterial::computePrincipalValues(), oofem::StressStrainBaseVector::convertFromFullForm(), E, oofem::MPSDamMaterialStatus::giveCharLength(), oofem::MPSDamMaterialStatus::giveDamage(), givee0(), givegf(), oofem::MPSDamMaterialStatus::giveKappa(), oofem::GaussPoint::giveMaterialMode(), oofem::StructuralMaterial::givePlaneStressVectorTranformationMtrx(), oofem::MPSMaterial::giveRealStressVector(), oofem::FloatArray::giveSize(), oofem::StructuralMaterial::giveSizeOfVoigtSymVector(), oofem::Material::giveStatus(), oofem::StructuralMaterial::giveStressVectorTranformationMtrx(), oofem::MPSMaterialStatus::giveTempAutogenousShrinkageStrain(), oofem::MPSMaterialStatus::giveTempDryingShrinkageStrain(), oofem::RheoChainMaterialStatus::giveTempThermalStrain(), initDamaged(), initDamagedFib(), oofem::RheoChainMaterial::isActivated(), isotropic, oofem::StructuralMaterialStatus::letTempStrainVectorBe(), oofem::StructuralMaterialStatus::letTempStressVectorBe(), oofem::MPSDamMaterialStatus::letTempViscoelasticStressVectorBe(), OOFEM_ERROR, oofem::principal_strain, oofem::FloatArray::resize(), oofem::FloatArray::resizeWithValues(), oofem::FloatArray::rotatedWith(), oofem::MPSDamMaterialStatus::setCrackWidth(), oofem::MPSDamMaterialStatus::setResidualTensileStrength(), oofem::MPSDamMaterialStatus::setTempDamage(), oofem::MPSDamMaterialStatus::setTempKappa(), softType, ST_Exponential_Cohesive_Crack, ST_Linear_Cohesive_Crack, timeDepFracturing, oofem::FloatArray::times(), and oofem::FloatArray::zero().

int oofem::MPSDamMaterial::hasMaterialModeCapability ( MaterialMode  mode)
virtual

Tests if material supports material mode.

Parameters
modeRequired material mode.
Returns
Nonzero if supported, zero otherwise.

Reimplemented from oofem::RheoChainMaterial.

Definition at line 239 of file mpsdammat.C.

virtual int oofem::MPSDamMaterial::hasNonLinearBehaviour ( )
inlinevirtual

Returns nonzero if receiver is non linear.

Reimplemented from oofem::KelvinChainSolidMaterial.

Definition at line 234 of file mpsdammat.h.

void oofem::MPSDamMaterial::initDamaged ( double  kappa,
FloatArray totalStrainVector,
GaussPoint gp,
TimeStep tStep 
)
IRResultType oofem::MPSDamMaterial::initializeFrom ( InputRecord ir)
virtual

Initializes receiver according to object description stored in input record.

This function is called immediately after creating object using constructor. Input record can be imagined as data record in component database belonging to receiver. Receiver may use value-name extracting functions to extract particular field from record.

See also
IR_GIVE_FIELD
IR_GIVE_OPTIONAL_FIELD
Parameters
irInput record to initialize from.
Returns
IRResultType

if sortion parameters are provided then it is assumed that the transport analysis uses moisture content/ratio; if not then it exports relative humidity

flag - if true, external fields and reference temperature are in Celsius

auxiliary parameters for autogenous shrinkage according to B4 model

Reimplemented from oofem::MPSMaterial.

Definition at line 250 of file mpsdammat.C.

References _IFT_MPSDamMaterial_checkSnapBack, _IFT_MPSDamMaterial_damageLaw, _IFT_MPSDamMaterial_fib_s, _IFT_MPSDamMaterial_ft, _IFT_MPSDamMaterial_ft28, _IFT_MPSDamMaterial_gf, _IFT_MPSDamMaterial_gf28, _IFT_MPSDamMaterial_isotropic, _IFT_MPSDamMaterial_timedepfracturing, _IFT_MPSMaterial_fc, checkSnapBack, const_gf, fib_fcm28, fib_s, ft, ft28, gf28, oofem::InputRecord::hasField(), oofem::MPSMaterial::initializeFrom(), IR_GIVE_FIELD, IR_GIVE_OPTIONAL_FIELD, oofem::IRRT_OK, isotropic, OOFEM_ERROR, softType, ST_Disable_Damage, ST_Exponential_Cohesive_Crack, ST_Linear_Cohesive_Crack, and timeDepFracturing.

Member Data Documentation

int oofem::MPSDamMaterial::checkSnapBack
protected

Check possible snap back flag.

Definition at line 217 of file mpsdammat.h.

Referenced by computeDamageForCohesiveCrack(), initDamaged(), initializeFrom(), and MPSDamMaterial().

double oofem::MPSDamMaterial::const_gf
protected

Determines the softening -> corresponds to the initial fracture energy.

For a linear law, it is the area under the stress/strain curve. For an exponential law, it is the area bounded by the elastic range and a tangent to the softening part of the curve at the peak stress. For a bilinear law, gf corresponds to area bounded by elasticity and the first linear softening line projected to zero stress.

Definition at line 209 of file mpsdammat.h.

Referenced by givegf(), initializeFrom(), and MPSDamMaterial().

double oofem::MPSDamMaterial::E
protected

dummy Young's modulus

Definition at line 192 of file mpsdammat.h.

Referenced by computeDamageForCohesiveCrack(), givee0(), giveRealStressVector(), initDamaged(), initDamagedFib(), and MPSDamMaterial().

ElementCharSizeMethod oofem::MPSDamMaterial::ecsMethod
protected

Method used for evaluation of characteristic element size.

Definition at line 223 of file mpsdammat.h.

Referenced by initDamaged(), and MPSDamMaterial().

double oofem::MPSDamMaterial::fib_fcm28
protected

Definition at line 188 of file mpsdammat.h.

Referenced by computeFractureEnergy(), computeTensileStrength(), and initializeFrom().

double oofem::MPSDamMaterial::fib_s
protected

Definition at line 187 of file mpsdammat.h.

Referenced by computeTensileStrength(), and initializeFrom().

double oofem::MPSDamMaterial::ft
protected

Equivalent strain at stress peak (or a similar parameter).

constant tensile strength

Definition at line 201 of file mpsdammat.h.

Referenced by givee0(), giveIPValue(), and initializeFrom().

double oofem::MPSDamMaterial::ft28
protected

28-day value of tensile strength. Used only with "timedepfracturing"

Definition at line 226 of file mpsdammat.h.

Referenced by computeTensileStrength(), and initializeFrom().

double oofem::MPSDamMaterial::gf28
protected

28-day value of fracture energy. Used only with "timedepfracturing"

Definition at line 228 of file mpsdammat.h.

Referenced by computeFractureEnergy(), and initializeFrom().

bool oofem::MPSDamMaterial::isotropic
protected
double oofem::MPSDamMaterial::maxOmega
protected

Maximum limit on omega. The purpose is elimination of a too compliant material which may cause convergence problems. Set to something like 0.99 if needed.

Definition at line 195 of file mpsdammat.h.

Referenced by give1dStressStiffMtrx(), give3dMaterialStiffnessMatrix(), givePlaneStrainStiffMtrx(), givePlaneStressStiffMtrx(), and MPSDamMaterial().

SofteningType oofem::MPSDamMaterial::softType
protected

Parameter specifying the type of softening (damage law).

Definition at line 220 of file mpsdammat.h.

Referenced by computeDamage(), computeDamageForCohesiveCrack(), giveRealStressVector(), initDamaged(), initializeFrom(), and MPSDamMaterial().

bool oofem::MPSDamMaterial::timeDepFracturing
protected

Definition at line 186 of file mpsdammat.h.

Referenced by givee0(), givegf(), giveRealStressVector(), initDamaged(), and initializeFrom().


The documentation for this class was generated from the following files:

This page is part of the OOFEM documentation. Copyright (c) 2011 Borek Patzak
Project e-mail: info@oofem.org
Generated at Tue Jan 2 2018 20:07:38 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011