OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
rcm2.h
Go to the documentation of this file.
1 /*
2  *
3  * ##### ##### ###### ###### ### ###
4  * ## ## ## ## ## ## ## ### ##
5  * ## ## ## ## #### #### ## # ##
6  * ## ## ## ## ## ## ## ##
7  * ## ## ## ## ## ## ## ##
8  * ##### ##### ## ###### ## ##
9  *
10  *
11  * OOFEM : Object Oriented Finite Element Code
12  *
13  * Copyright (C) 1993 - 2013 Borek Patzak
14  *
15  *
16  *
17  * Czech Technical University, Faculty of Civil Engineering,
18  * Department of Structural Mechanics, 166 29 Prague, Czech Republic
19  *
20  * This library is free software; you can redistribute it and/or
21  * modify it under the terms of the GNU Lesser General Public
22  * License as published by the Free Software Foundation; either
23  * version 2.1 of the License, or (at your option) any later version.
24  *
25  * This program is distributed in the hope that it will be useful,
26  * but WITHOUT ANY WARRANTY; without even the implied warranty of
27  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
28  * Lesser General Public License for more details.
29  *
30  * You should have received a copy of the GNU Lesser General Public
31  * License along with this library; if not, write to the Free Software
32  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
33  */
34 
35 #ifndef rcm2_h
36 #define rcm2_h
37 
38 #include "material.h"
40 #include "../sm/Materials/structuralmaterial.h"
41 #include "../sm/Materials/structuralms.h"
42 #include "intarray.h"
43 
45 
46 #define _IFT_RCM2Material_gf "gf"
47 #define _IFT_RCM2Material_ft "ft"
48 
49 
50 namespace oofem {
51 // material contant's keys for give()
52 #define pscm_Ee 300
53 #define pscm_Et 301
54 #define pscm_Gf 302
55 #define pscm_Beta 303
56 #define pscm_G 304
57 #define pscm_Ft 305
58 
59 // crack statuses
60 #define pscm_NONE 0
61 #define pscm_OPEN 1
62 #define pscm_SOFTENING 2
63 #define pscm_RELOADING 3
64 #define pscm_UNLOADING 4
65 #define pscm_CLOSED 5
66 
67 #define pscm_NEW_CRACK 20
68 #define pscm_NEW_FULLY_OPEN_CRACK 21
69 #define pscm_REOPEN_CRACK 22
70 
71 #define rcm_SMALL_STRAIN 1.e-7
72 #define rcm2_BIGNUMBER 1.e8
73 
78 {
79 protected:
88 
90 
91  // FloatArray minEffStrainsForFullyOpenCrack;
92  // FloatArray *crackStrainVector, *crackStrainIncrementVector;
93  // charLengths and minEffStrainsForFullyOpenCrack are not temp,
94  // because they are set only if a new crack is created,
95  // if a new crack is created this is recognized based on value in
96  // tempCrackStatuses->at(i), so we need not create temp version
97  // of theese variables.
98  //
102 
103 public:
104  RCM2MaterialStatus(int n, Domain * d, GaussPoint * g);
105  virtual ~RCM2MaterialStatus();
106 
107  virtual void printOutputAt(FILE *file, TimeStep *tStep);
108 
113  void letPrincipalStrainVectorBe(FloatArray pv) { principalStrain = std :: move(pv); }
114  void letPrincipalStressVectorBe(FloatArray pv) { principalStress = std :: move(pv); }
115 
116  const IntArray & giveCrackMap() const { return crackMap; }
117  void letCrackMapBe(IntArray map) { crackMap = std :: move(map); }
118  virtual int isCrackActive(int i) const;
119  virtual int giveNumberOfActiveCracks() const;
120  virtual int giveNumberOfTempActiveCracks() const;
121  int giveTempAlreadyCrack() const { return this->giveNumberOfTempActiveCracks(); }
122 
123  //double giveMinCrackStrainsForFullyOpenCrack (int icrack) {return minEffStrainsForFullyOpenCrack.at(icrack);}
124  //void setMinCrackStrainsForFullyOpenCrack (int icrack, double val) {minEffStrainsForFullyOpenCrack.at(icrack) = val;}
125 
127  void letTempCrackDirsBe(FloatMatrix a) { tempCrackDirs = std :: move(a); }
128  //const FloatArray & giveTempMaxCrackStrain() { return tempMaxCrackStrains;}
129  double giveTempMaxCrackStrain(int icrack) { return tempMaxCrackStrains.at(icrack); }
130  void setTempMaxCrackStrain(int icrack, double val) { tempMaxCrackStrains.at(icrack) = val; }
132  int giveTempCrackStatus(int icrack) const { return tempCrackStatuses.at(icrack); }
133  void setTempCrackStatus(int icrack, int val) { tempCrackStatuses.at(icrack) = val; }
134 
136  double giveCrackStrain(int icrack) const { return crackStrainVector.at(icrack); }
138  void letCrackStrainVectorBe(FloatArray a) { crackStrainVector = std :: move(a); }
139  void letOldCrackStrainVectorBe(FloatArray a) { oldCrackStrainVector = std :: move(a); }
140 
141  double giveCharLength(int icrack) const {
142  if ( icrack ) {
143  return charLengths.at(icrack);
144  } else {
145  return 0.0;
146  }
147  }
148  void setCharLength(int icrack, double val) { charLengths.at(icrack) = val; }
149 
150  // query for non-tem variables (usefull for postprocessing)
151  const FloatMatrix &giveCrackDirs() { return crackDirs; }
153  int giveAlreadyCrack() const { return this->giveNumberOfActiveCracks(); }
154 
155  // definition
156  virtual const char *giveClassName() const { return "RCM2MaterialStatus"; }
157 
158  virtual void initTempStatus();
159  virtual void updateYourself(TimeStep *tStep);
160 
161  // saves current context(state) into stream
162  virtual contextIOResultType saveContext(DataStream &stream, ContextMode mode, void *obj = NULL);
163  virtual contextIOResultType restoreContext(DataStream &stream, ContextMode mode, void *obj = NULL);
164 };
165 
179 {
180 protected:
182  double Gf, Ft;
183  //double beta;
184 
185 public:
186  RCM2Material(int n, Domain * d);
187  virtual ~RCM2Material();
188 
189  // identification and auxiliary functions
190  virtual int hasNonLinearBehaviour() { return 1; }
191  virtual int hasMaterialModeCapability(MaterialMode mode);
192 
193  virtual const char *giveClassName() const { return "RCM2Material"; }
194 
196 
197  virtual double give(int aProperty, GaussPoint *gp);
198 
199  LinearElasticMaterial *giveLinearElasticMaterial() { return linearElasticMaterial; }
200 
201  virtual void give3dMaterialStiffnessMatrix(FloatMatrix &answer,
202  MatResponseMode mode,
203  GaussPoint *gp,
204  TimeStep *tStep);
205 
206  virtual void giveRealStressVector(FloatArray &answer, GaussPoint *gp,
207  const FloatArray &reducedStrain, TimeStep *tStep);
208 
209  virtual void giveRealStressVector_3d(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
210  { this->giveRealStressVector(answer, gp, reducedE, tStep); }
211  virtual void giveRealStressVector_PlaneStrain(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
212  { this->giveRealStressVector(answer, gp, reducedE, tStep); }
213  virtual void giveRealStressVector_PlaneStress(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
214  { this->giveRealStressVector(answer, gp, reducedE, tStep); }
215  virtual void giveRealStressVector_1d(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
216  { this->giveRealStressVector(answer, gp, reducedE, tStep); }
217  virtual void giveRealStressVector_2dBeamLayer(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
218  { this->giveRealStressVector(answer, gp, reducedE, tStep); }
219  virtual void giveRealStressVector_PlateLayer(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
220  { this->giveRealStressVector(answer, gp, reducedE, tStep); }
221 
222 
223  virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep);
224 
225  virtual MaterialStatus *CreateStatus(GaussPoint *gp) const { return new RCM2MaterialStatus(1, domain, gp); }
226 
227  virtual void giveThermalDilatationVector(FloatArray &answer, GaussPoint *gp, TimeStep *tStep) { linearElasticMaterial->giveThermalDilatationVector(answer, gp, tStep); }
228 
229 protected:
230 
231  virtual void checkForNewActiveCracks(IntArray &answer, GaussPoint *gp, const FloatArray &,
232  const FloatArray &, FloatArray &, const FloatArray &);
233  virtual void updateCrackStatus(GaussPoint *gp, const FloatArray &crackStrain);
234 
235  virtual void checkIfClosedCracks(GaussPoint *gp, FloatArray &crackStrainVector, IntArray &);
236  virtual int checkSizeLimit(GaussPoint *gp, double) { return 0; }
237  virtual double giveNormalCrackingStress(GaussPoint *gp, double eps_cr, int i) = 0;
238  virtual double giveMinCrackStrainsForFullyOpenCrack(GaussPoint *gp, int i) = 0;
239  virtual double computeStrength(GaussPoint *gp, double) = 0;
240  virtual void updateStatusForNewCrack(GaussPoint *, int, double);
241  virtual double giveCharacteristicElementLength(GaussPoint *gp, const FloatArray &crackPlaneNormal);
242  virtual double giveCrackingModulus(MatResponseMode rMode, GaussPoint *gp,
243  double effStrain, int i) { return 1.e20; }
244 
245  virtual void giveMaterialStiffnessMatrix(FloatMatrix &answer, MatResponseMode,
246  GaussPoint *gp,
247  TimeStep *tStep);
248 
249  void giveCrackedStiffnessMatrix(FloatMatrix &answer,
250  MatResponseMode rMode,
251  GaussPoint *gp,
252  TimeStep *tStep);
253  /*
254  * void computeTrialStressIncrement (FloatArray& answer, GaussPoint *gp,
255  * const FloatArray& strainIncrement, TimeStep* tStep);
256  */
257  FloatMatrix *GiveCrackTransformationMtrx(GaussPoint *gp, int i);
258  virtual void giveEffectiveMaterialStiffnessMatrix(FloatMatrix &answer,
259  MatResponseMode rMode,
260  GaussPoint *gp, TimeStep *tStep);
261 
262  void giveRealPrincipalStressVector3d(FloatArray &answer, GaussPoint *,
264  void giveNormalElasticStiffnessMatrix(FloatMatrix &answer,
265  bool reduce, MatResponseMode,
266  GaussPoint *, TimeStep *tStep,
267  const FloatMatrix &);
268  void updateActiveCrackMap(GaussPoint *gp, const IntArray *activatedCracks = NULL);
269  // Give3dMaterialStiffnessMatrix should return 3d material stiffness matrix
270  // taking into account possible failure or fracture of material
271  double giveResidualStrength() { return 0.01 * this->Ft; }
272 
273 
274  virtual void givePlaneStressStiffMtrx(FloatMatrix &answer, MatResponseMode mmode,
275  GaussPoint *gp,
276  TimeStep *tStep);
277  virtual void givePlaneStrainStiffMtrx(FloatMatrix &answer, MatResponseMode mmode,
278  GaussPoint *gp,
279  TimeStep *tStep);
280  virtual void give1dStressStiffMtrx(FloatMatrix &answer, MatResponseMode mmode,
281  GaussPoint *gp,
282  TimeStep *tStep);
283  virtual void give2dBeamLayerStiffMtrx(FloatMatrix &answer, MatResponseMode mmode,
284  GaussPoint *gp,
285  TimeStep *tStep);
286  virtual void givePlateLayerStiffMtrx(FloatMatrix &answer, MatResponseMode mmode,
287  GaussPoint *gp,
288  TimeStep *tStep);
289 };
290 } // end namespace oofem
291 #endif // rcm2_h
const FloatArray & getPrincipalStrainVector() const
Definition: rcm2.h:109
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
virtual void giveRealStressVector_3d(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
Default implementation relies on giveRealStressVector for second Piola-Kirchoff stress.
Definition: rcm2.h:209
virtual void giveRealStressVector_1d(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
Default implementation relies on giveRealStressVector_StressControl.
Definition: rcm2.h:215
virtual const char * giveClassName() const
Definition: rcm2.h:156
GaussPoint * gp
Associated integration point.
Class and object Domain.
Definition: domain.h:115
Domain * domain
Link to domain object, useful for communicating with other FEM components.
Definition: femcmpnn.h:82
void letTempCrackDirsBe(FloatMatrix a)
Definition: rcm2.h:127
virtual contextIOResultType restoreContext(DataStream &stream, ContextMode mode, void *obj=NULL)
Restores the receiver state previously written in stream.
Definition: rcm2.C:1195
virtual void giveRealStressVector_2dBeamLayer(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
Default implementation relies on giveRealStressVector_StressControl.
Definition: rcm2.h:217
virtual void updateYourself(TimeStep *tStep)
Update equilibrium history variables according to temp-variables.
Definition: rcm2.C:1110
The purpose of DataStream abstract class is to allow to store/restore context to different streams...
Definition: datastream.h:54
FloatMatrix tempCrackDirs
Definition: rcm2.h:87
const IntArray & giveCrackStatus()
Definition: rcm2.h:152
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Print receiver's output to given stream.
Definition: rcm2.C:1010
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
LinearElasticMaterial * linearElasticMaterial
Definition: rcm2.h:181
virtual void giveRealStressVector_PlaneStrain(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
Default implementation relies on giveRealStressVector_3d.
Definition: rcm2.h:211
LinearElasticMaterial * giveLinearElasticMaterial()
Definition: rcm2.h:199
This class implements a structural material status information.
Definition: structuralms.h:65
const IntArray & giveTempCrackStatus()
Definition: rcm2.h:131
double giveCrackStrain(int icrack) const
Definition: rcm2.h:136
FloatArray oldCrackStrainVector
Definition: rcm2.h:85
int giveTempCrackStatus(int icrack) const
Definition: rcm2.h:132
void setTempCrackStatus(int icrack, int val)
Definition: rcm2.h:133
void letCrackMapBe(IntArray map)
Definition: rcm2.h:117
virtual void giveRealStressVector_PlateLayer(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
Default implementation relies on giveRealStressVector_StressControl.
Definition: rcm2.h:219
MaterialMode
Type representing material mode of integration point.
Definition: materialmode.h:89
FloatArray charLengths
Definition: rcm2.h:89
const FloatArray & givePrevPrincStrainVector() const
Definition: rcm2.h:111
Class implementing an array of integers.
Definition: intarray.h:61
int & at(int i)
Coefficient access function.
Definition: intarray.h:103
MatResponseMode
Describes the character of characteristic material matrix.
virtual MaterialStatus * CreateStatus(GaussPoint *gp) const
Creates new copy of associated status and inserts it into given integration point.
Definition: rcm2.h:225
const FloatMatrix & giveCrackDirs()
Definition: rcm2.h:151
This class is a abstract base class for all linear elastic material models in a finite element proble...
int giveAlreadyCrack() const
Definition: rcm2.h:153
double giveResidualStrength()
Definition: rcm2.h:271
RCM2MaterialStatus(int n, Domain *d, GaussPoint *g)
Definition: rcm2.C:972
This class implements a Rotating Crack Model for fracture in smeared fashion ( only material stiffnes...
Definition: rcm2.h:178
const FloatArray & getPrincipalStressVector() const
Definition: rcm2.h:110
virtual int giveNumberOfActiveCracks() const
Definition: rcm2.C:1049
virtual void giveRealStressVector_PlaneStress(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
Default implementation relies on giveRealStressVector_StressControl.
Definition: rcm2.h:213
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Definition: matstatus.h:140
void setCharLength(int icrack, double val)
Definition: rcm2.h:148
FloatArray maxCrackStrains
Max crack strain reached.
Definition: rcm2.h:83
const IntArray & giveCrackMap() const
Definition: rcm2.h:116
FloatArray crackStrainVector
Components of crack strain vector.
Definition: rcm2.h:85
This class implements associated Material Status to SmearedCrackingMaterail.
Definition: rcm2.h:77
const FloatArray & givePrevPrincStressVector() const
Definition: rcm2.h:112
const FloatArray & giveCrackStrainVector() const
Definition: rcm2.h:135
virtual contextIOResultType saveContext(DataStream &stream, ContextMode mode, void *obj=NULL)
Stores receiver state to output stream.
Definition: rcm2.C:1131
int giveTempAlreadyCrack() const
Definition: rcm2.h:121
virtual ~RCM2MaterialStatus()
Definition: rcm2.C:986
void letCrackStrainVectorBe(FloatArray a)
Definition: rcm2.h:138
Abstract base class representing a material status information.
Definition: matstatus.h:84
Class representing vector of real numbers.
Definition: floatarray.h:82
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
IRResultType
Type defining the return values of InputRecord reading operations.
Definition: irresulttype.h:47
virtual void initTempStatus()
Initializes the temporary internal variables, describing the current state according to previously re...
Definition: rcm2.C:1086
const FloatArray & giveOldCrackStrainVector()
Definition: rcm2.h:137
virtual int giveNumberOfTempActiveCracks() const
Definition: rcm2.C:1067
virtual void giveThermalDilatationVector(FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
Returns a vector of coefficients of thermal dilatation in direction of each material principal (local...
Class representing the general Input Record.
Definition: inputrecord.h:101
void setTempMaxCrackStrain(int icrack, double val)
Definition: rcm2.h:130
IntArray tempCrackStatuses
Definition: rcm2.h:81
void letPrincipalStressVectorBe(FloatArray pv)
Definition: rcm2.h:114
virtual int checkSizeLimit(GaussPoint *gp, double)
Definition: rcm2.h:236
virtual void giveThermalDilatationVector(FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
Returns a vector of coefficients of thermal dilatation in direction of each material principal (local...
Definition: rcm2.h:227
FloatArray principalStrain
Definition: rcm2.h:99
FloatArray principalStress
Definition: rcm2.h:100
FloatArray tempMaxCrackStrains
Definition: rcm2.h:83
FloatArray oldPrincipalStrain
Definition: rcm2.h:99
FloatArray oldPrincipalStress
Definition: rcm2.h:100
long ContextMode
Context mode (mask), defining the type of information written/read to/from context.
Definition: contextmode.h:43
Abstract base class for all "structural" constitutive models.
FloatMatrix crackDirs
Storing direction of cracks in columwise format.
Definition: rcm2.h:87
virtual int hasNonLinearBehaviour()
Returns nonzero if receiver is non linear.
Definition: rcm2.h:190
const FloatMatrix & giveTempCrackDirs()
Definition: rcm2.h:126
the oofem namespace is to define a context or scope in which all oofem names are defined.
virtual int isCrackActive(int i) const
Definition: rcm2.C:992
double giveTempMaxCrackStrain(int icrack)
Definition: rcm2.h:129
IntArray crackStatuses
One value from (pscm_NONE, pscm_OPEN, pscm_SOFTENING, pscm_RELOADING, pscm_UNLOADING, pscm_CLOSED.
Definition: rcm2.h:81
double giveCharLength(int icrack) const
Definition: rcm2.h:141
Class representing integration point in finite element program.
Definition: gausspoint.h:93
void letOldCrackStrainVectorBe(FloatArray a)
Definition: rcm2.h:139
Class representing solution step.
Definition: timestep.h:80
virtual double giveCrackingModulus(MatResponseMode rMode, GaussPoint *gp, double effStrain, int i)
Definition: rcm2.h:242
void letPrincipalStrainVectorBe(FloatArray pv)
Definition: rcm2.h:113
virtual const char * giveClassName() const
Definition: rcm2.h:193

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:31 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011