OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
libeam3d.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 libeam3d_h
36 #define libeam3d_h
37 
38 #include "../sm/Elements/structuralelement.h"
39 #include "../sm/CrossSections/fiberedcs.h"
40 
42 
43 #define _IFT_LIBeam3d_Name "libeam3d"
44 #define _IFT_LIBeam3d_refnode "refnode"
45 
46 
47 namespace oofem {
53 {
54 private:
55  double length;
57 
58 public:
59  LIBeam3d(int n, Domain * d);
60  virtual ~LIBeam3d() { }
61 
63 
64  virtual void computeLumpedMassMatrix(FloatMatrix &answer, TimeStep *tStep);
65  virtual void computeMassMatrix(FloatMatrix &answer, TimeStep *tStep)
66  { computeLumpedMassMatrix(answer, tStep); }
67  virtual void computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep);
68  virtual bool computeGtoLRotationMatrix(FloatMatrix &answer);
69 
70  virtual int testElementExtension(ElementExtension ext);
71 
72  virtual int computeNumberOfDofs() { return 12; }
73  virtual void giveDofManDofIDMask(int inode, IntArray &answer) const;
74  virtual double computeVolumeAround(GaussPoint *gp);
75  virtual int giveLocalCoordinateSystem(FloatMatrix &answer);
76  virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep);
77  virtual int computeGlobalCoordinates(FloatArray &answer, const FloatArray &lcoords);
78 
79  // Fibered cross section support functions
80  virtual void FiberedCrossSectionInterface_computeStrainVectorInFiber(FloatArray &answer, const FloatArray &masterGpStrain,
81  GaussPoint *slaveGp, TimeStep *tStep);
82 
84 
85  // definition & identification
86  virtual const char *giveInputRecordName() const { return _IFT_LIBeam3d_Name; }
87  virtual const char *giveClassName() const { return "LIBeam3d"; }
88 
89 #ifdef __OOFEG
90  virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep);
92 #endif
93 
94  virtual integrationDomain giveIntegrationDomain() const { return _Line; }
95  virtual MaterialMode giveMaterialMode() { return _3dBeam; }
96 
97 protected:
98  // edge load support
99  virtual void giveEdgeDofMapping(IntArray &answer, int iEdge) const;
100  virtual double computeEdgeVolumeAround(GaussPoint *gp, int iEdge);
101  virtual int computeLoadLEToLRotationMatrix(FloatMatrix &answer, int, GaussPoint *gp);
102  virtual int computeLoadGToLRotationMtrx(FloatMatrix &answer);
103  virtual void computeBodyLoadVectorAt(FloatArray &answer, Load *load, TimeStep *tStep, ValueModeType mode);
104 
105  virtual void computeBmatrixAt(GaussPoint *gp, FloatMatrix &answer, int = 1, int = ALL_STRAINS);
106  virtual void computeNmatrixAt(const FloatArray &iLocCoord, FloatMatrix &answer);
107  virtual void computeGaussPoints();
108  virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep);
109  virtual void computeStressVector(FloatArray &answer, const FloatArray &strain, GaussPoint *gp, TimeStep *tStep);
110 
111  virtual double computeLength();
112 };
113 } // end namespace oofem
114 #endif // libeam3d_h
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
integrationDomain
Used by integrator class to supply integration points for proper domain to be integrated (Area...
Class and object Domain.
Definition: domain.h:115
virtual void computeLumpedMassMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes lumped mass matrix of receiver.
Definition: libeam3d.C:133
virtual const char * giveInputRecordName() const
Definition: libeam3d.h:86
int referenceNode
Definition: libeam3d.h:56
virtual int computeLoadLEToLRotationMatrix(FloatMatrix &answer, int, GaussPoint *gp)
Returns transformation matrix from local edge c.s to element local coordinate system of load vector c...
Definition: libeam3d.C:364
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
Definition: libeam3d.C:107
ValueModeType
Type representing the mode of UnknownType or CharType, or similar types.
Definition: valuemodetype.h:78
virtual Interface * giveInterface(InterfaceType it)
Interface requesting service.
Definition: libeam3d.C:438
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
Definition: libeam3d.C:217
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
virtual int computeNumberOfDofs()
Computes or simply returns total number of element's local DOFs.
Definition: libeam3d.h:72
virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
Computes constitutive matrix of receiver.
Definition: libeam3d.C:119
virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep)
Definition: libeam3d.C:450
MaterialMode
Type representing material mode of integration point.
Definition: materialmode.h:89
LIBeam3d(int n, Domain *d)
Definition: libeam3d.C:55
Class implementing an array of integers.
Definition: intarray.h:61
MatResponseMode
Describes the character of characteristic material matrix.
virtual void computeNmatrixAt(const FloatArray &iLocCoord, FloatMatrix &answer)
Computes interpolation matrix for element unknowns.
Definition: libeam3d.C:148
virtual void drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType)
Definition: libeam3d.C:477
Abstract base class for all "structural" finite elements.
virtual double computeLength()
Computes the length (zero for all but 1D geometries)
Definition: libeam3d.C:266
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Definition: libeam3d.C:286
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
Definition: libeam3d.C:207
virtual void computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep)
Computes numerically stiffness matrix of receiver.
Definition: libeam3d.C:177
ElementExtension
Type representing element extension.
#define _IFT_LIBeam3d_Name
Definition: libeam3d.h:43
UnknownType
Type representing particular unknown (its physical meaning).
Definition: unknowntype.h:55
virtual void giveEdgeDofMapping(IntArray &answer, int iEdge) const
Assembles edge dof mapping mask, which provides mapping between edge local DOFs and "global" element ...
Definition: libeam3d.C:308
double length
Definition: libeam3d.h:55
virtual int computeGlobalCoordinates(FloatArray &answer, const FloatArray &lcoords)
Computes the global coordinates from given element's local coordinates.
Definition: libeam3d.C:241
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: libeam3d.C:126
virtual int computeLoadGToLRotationMtrx(FloatMatrix &answer)
Returns transformation matrix from global coordinate system to local element coordinate system for el...
Definition: libeam3d.C:337
#define ALL_STRAINS
virtual integrationDomain giveIntegrationDomain() const
Returns integration domain for receiver, used to initialize integration point over receiver volume...
Definition: libeam3d.h:94
virtual 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: libeam3d.C:378
virtual void computeBmatrixAt(GaussPoint *gp, FloatMatrix &answer, int=1, int=ALL_STRAINS)
Computes the geometrical matrix of receiver in given integration point.
Definition: libeam3d.C:65
virtual ~LIBeam3d()
Definition: libeam3d.h:60
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 int giveLocalCoordinateSystem(FloatMatrix &answer)
Returns local coordinate system of receiver Required by material models with ortho- and anisotrophy...
Definition: libeam3d.C:386
Class representing the general Input Record.
Definition: inputrecord.h:101
Class Interface.
Definition: interface.h:82
virtual void FiberedCrossSectionInterface_computeStrainVectorInFiber(FloatArray &answer, const FloatArray &masterGpStrain, GaussPoint *slaveGp, TimeStep *tStep)
Computes full 3d strain vector in element fiber.
Definition: libeam3d.C:422
This class implements a 3-dimensional mindlin theory Linear Isoparametric beam element, with reduced integration.
Definition: libeam3d.h:52
The element interface required by FiberedCrossSection.
Definition: fiberedcs.h:220
InterfaceType
Enumerative type, used to identify interface type.
Definition: interfacetype.h:43
virtual const char * giveClassName() const
Definition: libeam3d.h:87
Load is base abstract class for all loads.
Definition: load.h:61
virtual double computeEdgeVolumeAround(GaussPoint *gp, int iEdge)
Computes volume related to integration point on local edge.
Definition: libeam3d.C:326
the oofem namespace is to define a context or scope in which all oofem names are defined.
virtual int testElementExtension(ElementExtension ext)
Tests if the element implements required extension.
Definition: libeam3d.C:259
virtual void computeMassMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes mass matrix of receiver.
Definition: libeam3d.h:65
Class representing integration point in finite element program.
Definition: gausspoint.h:93
Class representing solution step.
Definition: timestep.h:80
virtual MaterialMode giveMaterialMode()
Returns material mode for receiver integration points.
Definition: libeam3d.h:95
virtual bool computeGtoLRotationMatrix(FloatMatrix &answer)
Returns transformation matrix from global c.s.
Definition: libeam3d.C:186
virtual void giveDofManDofIDMask(int inode, IntArray &answer) const
Returns dofmanager dof mask for node.
Definition: libeam3d.C:234

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