|
OOFEM
2.1
|
00001 /* 00002 * 00003 * ##### ##### ###### ###### ### ### 00004 * ## ## ## ## ## ## ## ### ## 00005 * ## ## ## ## #### #### ## # ## 00006 * ## ## ## ## ## ## ## ## 00007 * ## ## ## ## ## ## ## ## 00008 * ##### ##### ## ###### ## ## 00009 * 00010 * 00011 * OOFEM : Object Oriented Finite Element Code 00012 * 00013 * Copyright (C) 1993 - 2013 Borek Patzak 00014 * 00015 * 00016 * 00017 * Czech Technical University, Faculty of Civil Engineering, 00018 * Department of Structural Mechanics, 166 29 Prague, Czech Republic 00019 * 00020 * This program is free software; you can redistribute it and/or modify 00021 * it under the terms of the GNU General Public License as published by 00022 * the Free Software Foundation; either version 2 of the License, or 00023 * (at your option) any later version. 00024 * 00025 * This program is distributed in the hope that it will be useful, 00026 * but WITHOUT ANY WARRANTY; without even the implied warranty of 00027 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00028 * GNU General Public License for more details. 00029 * 00030 * You should have received a copy of the GNU General Public License 00031 * along with this program; if not, write to the Free Software 00032 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. 00033 */ 00034 00035 #ifndef rcm2_h 00036 #define rcm2_h 00037 00038 #include "material.h" 00039 #include "linearelasticmaterial.h" 00040 #include "structuralmaterial.h" 00041 #include "structuralms.h" 00042 00043 namespace oofem { 00044 // material contant's keys for give() 00045 #define pscm_Ee 300 00046 #define pscm_Et 301 00047 #define pscm_Gf 302 00048 #define pscm_Beta 303 00049 #define pscm_G 304 00050 #define pscm_Ft 305 00051 00052 // crack statuses 00053 #define pscm_NONE 0 00054 #define pscm_OPEN 1 00055 #define pscm_SOFTENING 2 00056 #define pscm_RELOADING 3 00057 #define pscm_UNLOADING 4 00058 #define pscm_CLOSED 5 00059 00060 #define pscm_NEW_CRACK 20 00061 #define pscm_NEW_FULLY_OPEN_CRACK 21 00062 #define pscm_REOPEN_CRACK 22 00063 00064 #define rcm_SMALL_STRAIN 1.e-7 00065 #define rcm2_BIGNUMBER 1.e8 00066 00070 class RCM2MaterialStatus : public StructuralMaterialStatus 00071 { 00072 protected: 00074 IntArray crackStatuses, tempCrackStatuses; 00076 FloatArray maxCrackStrains, tempMaxCrackStrains; 00078 FloatArray crackStrainVector, oldCrackStrainVector; 00080 FloatMatrix crackDirs, tempCrackDirs; 00081 00082 FloatArray charLengths; 00083 00084 // FloatArray minEffStrainsForFullyOpenCrack; 00085 // FloatArray *crackStrainVector, *crackStrainIncrementVector; 00086 // charLengths and minEffStrainsForFullyOpenCrack are not temp, 00087 // because they are set only if a new crack is created, 00088 // if a new crack is created this is recognized based on value in 00089 // tempCrackStatuses->at(i), so we need not create temp version 00090 // of theese variables. 00091 // 00092 FloatArray principalStrain, oldPrincipalStrain; 00093 FloatArray principalStress, oldPrincipalStress; 00094 IntArray crackMap; 00095 00096 public: 00097 RCM2MaterialStatus(int n, Domain *d, GaussPoint *g); 00098 virtual ~RCM2MaterialStatus(); 00099 00100 virtual void printOutputAt(FILE *file, TimeStep *tStep); 00101 00102 void getPrincipalStrainVector(FloatArray &answer) const { answer = principalStrain; } 00103 void getPrincipalStressVector(FloatArray &answer) const { answer = principalStress; } 00104 void givePrevPrincStrainVector(FloatArray &answer) const { answer = oldPrincipalStrain; } 00105 void givePrevPrincStressVector(FloatArray &answer) const { answer = oldPrincipalStress; } 00106 void letPrincipalStrainVectorBe(const FloatArray &pv) { principalStrain = pv; } 00107 void letPrincipalStressVectorBe(const FloatArray &pv) { principalStress = pv; } 00108 00109 void giveCrackMap(IntArray &answer) const { answer = crackMap; } 00110 void letCrackMapBe(IntArray &map) { crackMap = map; } 00111 virtual int isCrackActive(int i) const; 00112 virtual int giveNumberOfActiveCracks() const; 00113 virtual int giveNumberOfTempActiveCracks() const; 00114 int giveTempAlreadyCrack() const { return this->giveNumberOfTempActiveCracks(); } 00115 00116 //double giveMinCrackStrainsForFullyOpenCrack (int icrack) {return minEffStrainsForFullyOpenCrack.at(icrack);} 00117 //void setMinCrackStrainsForFullyOpenCrack (int icrack, double val) {minEffStrainsForFullyOpenCrack.at(icrack) = val;} 00118 00119 void giveTempCrackDirs(FloatMatrix &answer) { answer = tempCrackDirs; } 00120 void letTempCrackDirsBe(const FloatMatrix &a) { tempCrackDirs = a; } 00121 //void giveTempMaxCrackStrain(FloatArray& answer) {answer = tempMaxCrackStrains;} 00122 double giveTempMaxCrackStrain(int icrack) { return tempMaxCrackStrains.at(icrack); } 00123 void setTempMaxCrackStrain(int icrack, double val) { tempMaxCrackStrains.at(icrack) = val; } 00124 void giveTempCrackStatus(IntArray &answer) { answer = tempCrackStatuses; } 00125 int giveTempCrackStatus(int icrack) const { return tempCrackStatuses.at(icrack); } 00126 void setTempCrackStatus(int icrack, int val) { tempCrackStatuses.at(icrack) = val; } 00127 00128 void giveCrackStrainVector(FloatArray &answer) { answer = crackStrainVector; } 00129 double giveCrackStrain(int icrack) const { return crackStrainVector.at(icrack); } 00130 void giveOldCrackStrainVector(FloatArray &answer) { answer = oldCrackStrainVector; } 00131 void letCrackStrainVectorBe(const FloatArray &a) { crackStrainVector = a; } 00132 void letOldCrackStrainVectorBe(const FloatArray &a) { oldCrackStrainVector = a; } 00133 00134 double giveCharLength(int icrack) const { if ( icrack ) { return charLengths.at(icrack); } else { return 0.0; } } 00135 void setCharLength(int icrack, double val) { charLengths.at(icrack) = val; } 00136 00137 // query for non-tem variables (usefull for postprocessing) 00138 void giveCrackDirs(FloatMatrix &answer) { answer = crackDirs; } 00139 void giveCrackStatus(IntArray &answer) { answer = crackStatuses; } 00140 int giveAlreadyCrack() const { return this->giveNumberOfActiveCracks(); } 00141 00142 // definition 00143 virtual const char *giveClassName() const { return "RCM2MaterialStatus"; } 00144 virtual classType giveClassID() const { return RCMMaterialStatusClass; } 00145 00146 virtual void initTempStatus(); 00147 virtual void updateYourself(TimeStep *tStep); 00148 00149 // saves current context(state) into stream 00150 virtual contextIOResultType saveContext(DataStream *stream, ContextMode mode, void *obj = NULL); 00151 virtual contextIOResultType restoreContext(DataStream *stream, ContextMode mode, void *obj = NULL); 00152 }; 00153 00166 class RCM2Material : public StructuralMaterial 00167 { 00168 protected: 00169 LinearElasticMaterial *linearElasticMaterial; 00170 double Gf, Ft; 00171 //double beta; 00172 00173 public: 00174 RCM2Material(int n, Domain *d); 00175 virtual ~RCM2Material(); 00176 00177 // identification and auxiliary functions 00178 virtual int hasNonLinearBehaviour() { return 1; } 00179 virtual int hasMaterialModeCapability(MaterialMode mode); 00180 00181 virtual const char *giveClassName() const { return "RCM2Material"; } 00182 virtual classType giveClassID() const { return RCMMaterialClass; } 00183 00184 virtual IRResultType initializeFrom(InputRecord *ir); 00185 00186 virtual double give(int aProperty, GaussPoint *gp); 00187 00188 LinearElasticMaterial *giveLinearElasticMaterial() { return linearElasticMaterial; } 00189 00190 virtual void give3dMaterialStiffnessMatrix(FloatMatrix &answer, 00191 MatResponseForm form, MatResponseMode mode, 00192 GaussPoint *gp, 00193 TimeStep *tStep); 00194 00195 virtual void giveRealStressVector(FloatArray &answer, MatResponseForm form, GaussPoint *gp, 00196 const FloatArray &reducedStrain, TimeStep *tStep); 00197 00198 virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep); 00199 virtual int giveIntVarCompFullIndx(IntArray &answer, InternalStateType type, MaterialMode mmode); 00200 virtual InternalStateValueType giveIPValueType(InternalStateType type); 00201 virtual int giveIPValueSize(InternalStateType type, GaussPoint *gp); 00202 00203 virtual MaterialStatus *CreateStatus(GaussPoint *gp) const { return new RCM2MaterialStatus(1, domain, gp); } 00204 00205 protected: 00206 00207 virtual void initTempStatus(GaussPoint *gp); 00208 virtual void checkForNewActiveCracks(IntArray &answer, GaussPoint *gp, const FloatArray &, 00209 const FloatArray &, FloatArray &, const FloatArray &); 00210 virtual void updateCrackStatus(GaussPoint *gp, const FloatArray &crackStrain); 00211 00212 virtual void checkIfClosedCracks(GaussPoint *gp, FloatArray &crackStrainVector, IntArray &); 00213 virtual int checkSizeLimit(GaussPoint *gp, double) { return 0; } 00214 virtual double giveNormalCrackingStress(GaussPoint *gp, double eps_cr, int i) = 0; 00215 virtual double giveMinCrackStrainsForFullyOpenCrack(GaussPoint *gp, int i) = 0; 00216 virtual double computeStrength(GaussPoint *gp, double) = 0; 00217 virtual void updateStatusForNewCrack(GaussPoint *, int, double); 00218 virtual double giveCharacteristicElementLenght(GaussPoint *gp, const FloatArray &); 00219 virtual double giveCrackingModulus(MatResponseMode rMode, GaussPoint *gp, 00220 double effStrain, int i) { return 1.e20; } 00221 00222 virtual void giveMaterialStiffnessMatrix(FloatMatrix & answer, MatResponseForm, MatResponseMode, 00223 GaussPoint * gp, 00224 TimeStep * tStep); 00225 00226 void giveCrackedStiffnessMatrix(FloatMatrix &answer, 00227 MatResponseMode rMode, 00228 GaussPoint *gp, 00229 TimeStep *tStep); 00230 /* 00231 * void computeTrialStressIncrement (FloatArray& answer, GaussPoint *gp, 00232 * const FloatArray& strainIncrement, TimeStep* tStep); 00233 */ 00234 FloatMatrix *GiveCrackTransformationMtrx(GaussPoint *gp, int i); 00235 virtual void giveEffectiveMaterialStiffnessMatrix(FloatMatrix &answer, MatResponseForm form, 00236 MatResponseMode rMode, 00237 GaussPoint *gp, TimeStep *tStep); 00238 00239 void giveRealPrincipalStressVector3d(FloatArray &answer, GaussPoint *, 00240 FloatArray &, FloatMatrix &, TimeStep *); 00241 void giveNormalElasticStiffnessMatrix(FloatMatrix & answer, 00242 MatResponseForm, MatResponseMode, 00243 GaussPoint *, TimeStep * tStep, 00244 const FloatMatrix &); 00245 void updateActiveCrackMap(GaussPoint *gp, const IntArray *activatedCracks = NULL); 00246 // Give3dMaterialStiffnessMatrix should return 3d material stiffness matrix 00247 // taking into account possible failure or fracture of material 00248 double giveResidualStrength() { return 0.01 * this->Ft; } 00249 00250 00251 virtual void givePlaneStressStiffMtrx(FloatMatrix &answer, MatResponseForm form, MatResponseMode mmode, 00252 GaussPoint *gp, 00253 TimeStep *tStep); 00254 virtual void givePlaneStrainStiffMtrx(FloatMatrix &answer, MatResponseForm form, MatResponseMode mmode, 00255 GaussPoint *gp, 00256 TimeStep *tStep); 00257 virtual void give1dStressStiffMtrx(FloatMatrix &answer, MatResponseForm form, MatResponseMode mmode, 00258 GaussPoint *gp, 00259 TimeStep *tStep); 00260 virtual void give2dBeamLayerStiffMtrx(FloatMatrix &answer, MatResponseForm form, MatResponseMode mmode, 00261 GaussPoint *gp, 00262 TimeStep *tStep); 00263 virtual void give2dPlateLayerStiffMtrx(FloatMatrix &answer, MatResponseForm form, MatResponseMode mmode, 00264 GaussPoint *gp, 00265 TimeStep *tStep); 00266 virtual void give3dShellLayerStiffMtrx(FloatMatrix &answer, MatResponseForm form, MatResponseMode mmode, 00267 GaussPoint *gp, 00268 TimeStep *tStep); 00269 }; 00270 } // end namespace oofem 00271 #endif // rcm2_h