OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
beam2d.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 beam2d_h
36 #define beam2d_h
37 
38 #include "../sm/Elements/Beams/beambaseelement.h"
39 #include "../sm/CrossSections/layeredcrosssection.h"
40 #include "dofmanager.h"
41 
43 
44 #define _IFT_Beam2d_Name "beam2d"
45 #define _IFT_Beam2d_dofstocondense "dofstocondense"
46 
47 
48 namespace oofem {
49 class FEI2dLineLin;
50 class FEI2dLineHermite;
51 
63 {
64 protected:
65  double kappa, pitch, length;
75 
78 
79 public:
80  Beam2d(int n, Domain *aDomain);
81  virtual ~Beam2d();
82 
83  virtual void computeConsistentMassMatrix(FloatMatrix &answer, TimeStep *tStep, double &mass, const double *ipDensity = NULL);
84  virtual void computeInitialStressMatrix(FloatMatrix &answer, TimeStep *tStep);
85  virtual void computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep);
86  virtual int giveLocalCoordinateSystem(FloatMatrix &answer);
87  virtual void giveInternalForcesVector(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord = 0);
88  virtual void giveEndForcesVector(FloatArray &answer, TimeStep *tStep);
89 
90  virtual int testElementExtension(ElementExtension ext) { return ( ext == Element_EdgeLoadSupport ); }
91 
93 
94  virtual FEInterpolation *giveInterpolation() const;
95  virtual FEInterpolation *giveInterpolation(DofIDItem id) const { return NULL; }
96 
97  virtual int computeNumberOfDofs() { return 6; }
98  virtual int computeNumberOfGlobalDofs() { return 6 + this->numberOfCondensedDofs; }
99  virtual void giveDofManDofIDMask(int inode, IntArray &) const;
100  virtual int giveNumberOfInternalDofManagers() const { return ( ghostNodes [ 0 ] != NULL ) + ( ghostNodes [ 1 ] != NULL ); }
101  virtual DofManager *giveInternalDofManager(int i) const {
102  if ( i == 1 ) {
103  if ( ghostNodes [ 0 ] ) { return ghostNodes [ 0 ]; } else { return ghostNodes [ 1 ]; }
104  } else if ( i == 2 ) { // i==2
105  return ghostNodes [ 1 ];
106  } else {
107  OOFEM_ERROR("No such DOF available on Element %d", number);
108  return NULL;
109  }
110  }
111  virtual void giveInternalDofManDofIDMask(int i, IntArray &answer) const {
112  if ( i == 1 ) {
113  if ( ghostNodes [ 0 ] ) {
114  ghostNodes [ 0 ]->giveCompleteMasterDofIDArray(answer);
115  } else {
116  ghostNodes [ 1 ]->giveCompleteMasterDofIDArray(answer);
117  }
118  } else if ( i == 2 ) { // i==2
119  ghostNodes [ 1 ]->giveCompleteMasterDofIDArray(answer);
120  } else {
121  OOFEM_ERROR("No such DOF available on Element %d", number);
122  }
123  }
124 
125  virtual void giveBoundaryLocationArray(IntArray &locationArray, const IntArray &bNodes, const UnknownNumberingScheme &s, IntArray *dofIds = NULL) {
126  giveLocationArray (locationArray, s, dofIds);
127  }
128 
129  virtual void giveBoundaryLocationArray(IntArray &locationArray, const IntArray &bNodes, const IntArray &dofIDMask, const UnknownNumberingScheme &s, IntArray *dofIds = NULL) {
130  giveLocationArray (locationArray, dofIDMask, s, dofIds);
131  }
132 
133 
134  virtual double computeVolumeAround(GaussPoint *gp);
135  virtual void printOutputAt(FILE *file, TimeStep *tStep);
136 
137  virtual const char *giveClassName() const { return "Beam2d"; }
138  virtual const char *giveInputRecordName() const { return _IFT_Beam2d_Name; }
140 
141 #ifdef __OOFEG
142  virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep);
143  virtual void drawDeformedGeometry(oofegGraphicContext & gc, TimeStep * tStep, UnknownType);
144 #endif
145 
146  virtual void computeStrainVectorInLayer(FloatArray &answer, const FloatArray &masterGpStrain,
147  GaussPoint *masterGp, GaussPoint *slaveGp, TimeStep *tStep);
148 
149  virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep);
150 
151 protected:
152  virtual void computeBoundaryEdgeLoadVector(FloatArray &answer, BoundaryLoad *load, int edge, CharType type, ValueModeType mode, TimeStep *tStep, bool global=true);
153  virtual void computeBmatrixAt(GaussPoint *, FloatMatrix &, int = 1, int = ALL_STRAINS);
154  virtual void computeNmatrixAt(const FloatArray &iLocCoord, FloatMatrix &);
155  virtual bool computeGtoLRotationMatrix(FloatMatrix &answer);
156 
157  virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep);
158  virtual void computeStressVector(FloatArray &answer, const FloatArray &strain, GaussPoint *gp, TimeStep *tStep);
159 
160  void computeBodyLoadVectorAt(FloatArray &answer, Load *load, TimeStep *tStep, ValueModeType mode);
161 
162  double giveKappaCoeff(TimeStep *tStep);
163  virtual double computeLength();
164  double givePitch();
165  virtual void computeGaussPoints();
166  virtual MaterialMode giveMaterialMode() { return _2dBeam; }
167  virtual int giveNumberOfIPForMassMtrxIntegration() { return 4; }
168 
169  bool hasDofs2Condense() { return ( ghostNodes [ 0 ] || ghostNodes [ 1 ] ); }
170 };
171 } // end namespace oofem
172 #endif // beam2d_h
virtual FEInterpolation * giveInterpolation(DofIDItem id) const
Returns the interpolation for the specific dof id.
Definition: beam2d.h:95
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
Definition: beam2d.C:129
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
virtual const char * giveInputRecordName() const
Definition: beam2d.h:138
int number
Component number.
Definition: femcmpnn.h:80
bool hasDofs2Condense()
Definition: beam2d.h:169
virtual void computeNmatrixAt(const FloatArray &iLocCoord, FloatMatrix &)
Computes interpolation matrix for element unknowns.
Definition: beam2d.C:142
virtual Interface * giveInterface(InterfaceType)
Interface requesting service.
Definition: beam2d.C:91
Class and object Domain.
Definition: domain.h:115
double kappa
Definition: beam2d.h:65
virtual int computeNumberOfDofs()
Computes or simply returns total number of element's local DOFs.
Definition: beam2d.h:97
virtual bool computeGtoLRotationMatrix(FloatMatrix &answer)
Returns transformation matrix from global c.s.
Definition: beam2d.C:207
virtual void computeBoundaryEdgeLoadVector(FloatArray &answer, BoundaryLoad *load, int edge, CharType type, ValueModeType mode, TimeStep *tStep, bool global=true)
Computes the contribution of the given load at the given boundary edge.
Definition: beam2d.C:444
virtual void computeStrainVectorInLayer(FloatArray &answer, const FloatArray &masterGpStrain, GaussPoint *masterGp, GaussPoint *slaveGp, TimeStep *tStep)
Computes full 3D strain vector in element layer.
Definition: beam2d.C:271
ValueModeType
Type representing the mode of UnknownType or CharType, or similar types.
Definition: valuemodetype.h:78
virtual void computeStressVector(FloatArray &answer, const FloatArray &strain, GaussPoint *gp, TimeStep *tStep)
Computes the stress vector of receiver at given integration point, at time step tStep.
Definition: beam2d.C:200
Class representing a 2d line with Hermitian interpolation.
virtual void giveInternalDofManDofIDMask(int i, IntArray &answer) const
Returns internal dofmanager dof mask for node.
Definition: beam2d.h:111
virtual void giveDofManDofIDMask(int inode, IntArray &) const
Returns dofmanager dof mask for node.
Definition: beam2d.C:292
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
#define _IFT_Beam2d_Name
Definition: beam2d.h:44
virtual int computeNumberOfGlobalDofs()
Computes the total number of element's global dofs.
Definition: beam2d.h:98
This class implements a 2-dimensional beam element with cubic lateral displacement, quadratic rotations, and linear longitudinal displacements and geometry.
Definition: beam2d.h:62
double pitch
Definition: beam2d.h:65
Base class for dof managers.
Definition: dofmanager.h:113
virtual double computeLength()
Computes the length (zero for all but 1D geometries)
Definition: beam2d.C:301
int numberOfCondensedDofs
number of condensed DOFs
Definition: beam2d.h:74
MaterialMode
Type representing material mode of integration point.
Definition: materialmode.h:89
Class implementing an array of integers.
Definition: intarray.h:61
MatResponseMode
Describes the character of characteristic material matrix.
Beam2d(int n, Domain *aDomain)
Definition: beam2d.C:66
double givePitch()
Definition: beam2d.C:319
virtual FEInterpolation * giveInterpolation() const
Definition: beam2d.C:87
Class representing a general abstraction for finite element interpolation class.
Definition: feinterpol.h:132
virtual void giveBoundaryLocationArray(IntArray &locationArray, const IntArray &bNodes, const UnknownNumberingScheme &s, IntArray *dofIds=NULL)
Returns the location array for the boundary of the element.
Definition: beam2d.h:125
virtual int testElementExtension(ElementExtension ext)
Tests if the element implements required extension.
Definition: beam2d.h:90
virtual int giveNumberOfInternalDofManagers() const
Definition: beam2d.h:100
Abstract base class representing a boundary load (force, momentum, ...) that acts directly on a bound...
Definition: boundaryload.h:110
virtual void computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep)
Computes numerically stiffness matrix of receiver.
Definition: beam2d.C:176
virtual MaterialMode giveMaterialMode()
Returns material mode for receiver integration points.
Definition: beam2d.h:166
virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
Computes constitutive matrix of receiver.
Definition: beam2d.C:193
DofManager * ghostNodes[2]
Ghost nodes are used to introduce additional DOFs at element.
Definition: beam2d.h:72
#define OOFEM_ERROR(...)
Definition: error.h:61
ElementExtension
Type representing element extension.
DofIDItem
Type representing particular dof type.
Definition: dofiditem.h:86
virtual const char * giveClassName() const
Definition: beam2d.h:137
UnknownType
Type representing particular unknown (its physical meaning).
Definition: unknowntype.h:55
void giveLocationArray(IntArray &locationArray, const UnknownNumberingScheme &s, IntArray *dofIds=NULL) const
Returns the location array (array of code numbers) of receiver for given numbering scheme...
Definition: element.C:390
Abstract base class allowing to control the way, how equations are assigned to individual DOFs...
void giveCompleteMasterDofIDArray(IntArray &dofIDArray) const
Returns the full dof ID array of receiver.
Definition: dofmanager.C:249
Class representing a 2d line with linear interpolation.
Definition: fei2dlinelin.h:46
#define ALL_STRAINS
virtual void giveBoundaryLocationArray(IntArray &locationArray, const IntArray &bNodes, const IntArray &dofIDMask, const UnknownNumberingScheme &s, IntArray *dofIds=NULL)
Definition: beam2d.h:129
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Prints output of receiver to stream, for given time step.
Definition: beam2d.C:541
Class representing vector of real numbers.
Definition: floatarray.h:82
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
virtual void giveEndForcesVector(FloatArray &answer, TimeStep *tStep)
Definition: beam2d.C:427
IRResultType
Type defining the return values of InputRecord reading operations.
Definition: irresulttype.h:47
virtual void giveInternalForcesVector(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord=0)
Returns equivalent nodal forces vectors.
Definition: beam2d.C:420
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
Definition: beam2d.C:264
CharType
Definition: chartype.h:87
Class representing the general Input Record.
Definition: inputrecord.h:101
Class Interface.
Definition: interface.h:82
virtual ~Beam2d()
Definition: beam2d.C:80
double length
Definition: beam2d.h:65
static FEI2dLineLin interp_geom
Definition: beam2d.h:76
virtual void computeConsistentMassMatrix(FloatMatrix &answer, TimeStep *tStep, double &mass, const double *ipDensity=NULL)
Computes consistent mass matrix of receiver using numerical integration over element volume...
Definition: beam2d.C:572
InterfaceType
Enumerative type, used to identify interface type.
Definition: interfacetype.h:43
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Definition: beam2d.C:381
virtual void computeInitialStressMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes initial stress matrix for linear stability problem.
Definition: beam2d.C:622
Load is base abstract class for all loads.
Definition: load.h:61
The element interface required by LayeredCrossSection.
the oofem namespace is to define a context or scope in which all oofem names are defined.
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
Definition: beam2d.C:499
virtual DofManager * giveInternalDofManager(int i) const
Returns i-th internal element dof manager of the receiver.
Definition: beam2d.h:101
static FEI2dLineHermite interp_beam
Definition: beam2d.h:77
virtual void computeBmatrixAt(GaussPoint *, FloatMatrix &, int=1, int=ALL_STRAINS)
Computes the geometrical matrix of receiver in given integration point.
Definition: beam2d.C:102
virtual int giveLocalCoordinateSystem(FloatMatrix &answer)
Returns local coordinate system of receiver Required by material models with ortho- and anisotrophy...
Definition: beam2d.C:356
Class representing integration point in finite element program.
Definition: gausspoint.h:93
Class representing solution step.
Definition: timestep.h:80
virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep)
Definition: beam2d.C:665
This class implements a base beam intented to be a base class for beams based on lagrangian interpola...
virtual void drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType)
Definition: beam2d.C:691
void computeBodyLoadVectorAt(FloatArray &answer, Load *load, TimeStep *tStep, ValueModeType mode)
Computes the load vector due to body load acting on receiver, at given time step. ...
Definition: beam2d.C:490
Element extension for edge loads.
double giveKappaCoeff(TimeStep *tStep)
Definition: beam2d.C:339
virtual int giveNumberOfIPForMassMtrxIntegration()
Return desired number of integration points for consistent mass matrix computation, if required.
Definition: beam2d.h:167

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