OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
transportelement.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 transportelement_h
36 #define transportelement_h
37 
38 #include "element.h"
39 #include "floatmatrix.h"
40 #include "primaryfield.h"
41 #include "matresponsemode.h"
42 
43 namespace oofem {
44 class TransportCrossSection;
45 
52 {
53 public:
55 
56 protected:
59  static const double stefanBoltzmann;
60 
61 public:
63  virtual ~TransportElement();
64 
65  virtual void giveCharacteristicMatrix(FloatMatrix &answer, CharType type, TimeStep *tStep);
66  virtual void giveCharacteristicVector(FloatArray &answer, CharType type, ValueModeType mode, TimeStep *tStep);
67 
68  virtual void giveDofManDofIDMask(int inode, IntArray &answer) const;
69 
70  virtual void computeInternalForcesVector(FloatArray &answer, TimeStep *tStep);
71  virtual void computeExternalForcesVector(FloatArray &answer, TimeStep *tStep, ValueModeType mode);
72  virtual void computeInertiaForcesVector(FloatArray &answer, TimeStep *tStep);
73  virtual void computeLumpedCapacityVector(FloatArray &answer, TimeStep *tStep);
74 
75  //Compute volumetric load from element
76  virtual void computeLoadVector(FloatArray &answer, BodyLoad *load, CharType type, ValueModeType mode, TimeStep *tStep);
77  //Boundary load by prescribed flux, convection, or radiation over surface
78  virtual void computeBoundarySurfaceLoadVector(FloatArray &answer, BoundaryLoad *load, int boundary, CharType type, ValueModeType mode, TimeStep *tStep, bool global=true);
79  //Contribution to conductivity matrix from convection
80  virtual void computeTangentFromSurfaceLoad(FloatMatrix &answer, SurfaceLoad *load, int boundary, MatResponseMode rmode, TimeStep *tStep);
81  virtual void computeTangentFromEdgeLoad(FloatMatrix &answer, EdgeLoad *load, int boundary, MatResponseMode rmode, TimeStep *tStep);
82  //Boundary load by prescribed flux, convection, or radiation over length
83  virtual void computeBoundaryEdgeLoadVector(FloatArray &answer, BoundaryLoad *load, int edge, CharType type, ValueModeType mode, TimeStep *tStep, bool global=true);
84 
85  virtual void computeField(ValueModeType mode, TimeStep *tStep, const FloatArray &lcoords, FloatArray &answer);
86 
88  virtual Material * giveMaterial();
89 
97  virtual double giveThicknessAt(const FloatArray &gcoords) { return 1.0; }
98 
100  virtual void computeCapacityMatrix(FloatMatrix &answer, TimeStep *tStep);
102  virtual void computeConductivityMatrix(FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep);
104  virtual void computeBCVectorAt(FloatArray &answer, TimeStep *tStep, ValueModeType mode);
106  virtual void computeBCMtrxAt(FloatMatrix &answer, TimeStep *tStep, ValueModeType mode);
108  virtual void computeInternalSourceRhsVectorAt(FloatArray &answer, TimeStep *tStep, ValueModeType mode);
110  virtual void computeIntSourceLHSMatrix(FloatMatrix &answer, TimeStep *tStep);
112  virtual void computeIntSourceLHSSubMatrix(FloatMatrix &answer, MatResponseMode rmode, int iri, TimeStep *tStep);
119  virtual void computeFlow(FloatArray &answer, GaussPoint *gp, TimeStep *tStep);
120 
121  // time step termination
122  virtual void updateInternalState(TimeStep *tStep);
123  virtual int checkConsistency();
124 
126  const FloatArray &coords, IntArray &dofId, ValueModeType mode,
127  TimeStep *tStep);
128 
129 #ifdef __OOFEG
131  int node, TimeStep *tStep);
132  // Graphics output
133  //virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep) {}
134  //virtual void drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType) {}
135 #endif
136 
144  virtual void computeConstitutiveMatrixAt(FloatMatrix &answer,
145  MatResponseMode rMode, GaussPoint *gp,
146  TimeStep *tStep);
147 
148 protected:
157  virtual void computeCapacitySubMatrix(FloatMatrix &answer, MatResponseMode rmode, int iri, TimeStep *tStep);
167  virtual void computeConductivitySubMatrix(FloatMatrix &answer, int iri, MatResponseMode rmode, TimeStep *tStep);
176  void computeBCSubVectorAt(FloatArray &answer, TimeStep *tStep, ValueModeType mode, int indx);
185  void computeBCSubMtrxAt(FloatMatrix &answer, TimeStep *tStep, ValueModeType mode, int indx);
186  void computeBodyBCSubVectorAt(FloatArray &answer, Load *load, TimeStep *tStep, ValueModeType mode, int indx);
198  void computeEdgeBCSubVectorAt(FloatArray &answer, Load *load, int iEdge,
199  TimeStep *tStep, ValueModeType mode, int indx);
200 
212  void computeSurfaceBCSubVectorAt(FloatArray &answer, Load *load, int iSurf,
213  TimeStep *tStep, ValueModeType mode, int indx);
214 
220  virtual void computeNAt(FloatArray &answer, const FloatArray &lcoord);
225  virtual void computeNmatrixAt(FloatMatrix &answer, const FloatArray &lcoords);
226  virtual void computeBmatrixAt(FloatMatrix &answer, const FloatArray &lcoords);
230  virtual void computeGradientMatrixAt(FloatMatrix &answer, const FloatArray &lcoords);
234  virtual void computeInternalSourceRhsSubVectorAt(FloatArray &answer, TimeStep *, ValueModeType mode, int indx);
235 
239  virtual void computeEgdeNAt(FloatArray &answer, int iEdge, const FloatArray &lcoord);
243  virtual void giveEdgeDofMapping(IntArray &mask, int iEdge);
247  virtual double computeEdgeVolumeAround(GaussPoint *gp, int iEdge) = 0;
248 
249  virtual IntegrationRule *GetSurfaceIntegrationRule(int approxOrder) { return NULL; }
250  virtual void computeSurfaceNAt(FloatArray &answer, int iSurf, const FloatArray &lcoord);
251  virtual double computeSurfaceVolumeAround(GaussPoint *gp, int iSurf) { return 0.; }
252  virtual void giveSurfaceDofMapping(IntArray &mask, int iSurf);
253 
254  virtual int giveApproxOrder(int unknownIndx);
255 
269  void assembleLocalContribution(FloatMatrix &answer, FloatMatrix &src, int ndofs, int rdof, int cdof);
282  void assembleLocalContribution(FloatArray &answer, FloatArray &src, int ndofs, int rdof);
283 };
284 } // end namespace oofem
285 #endif // transportelement_h
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
virtual void computeTangentFromSurfaceLoad(FloatMatrix &answer, SurfaceLoad *load, int boundary, MatResponseMode rmode, TimeStep *tStep)
Computes the tangent contribution of the given load at the given boundary.
int giveInternalStateAtNode(FloatArray &answer, InternalStateType type, InternalStateMode mode, int node, TimeStep *tStep)
Returns internal state variable (like stress,strain) at node of element in Reduced form...
void assembleLocalContribution(FloatMatrix &answer, FloatMatrix &src, int ndofs, int rdof, int cdof)
Assembles the given source matrix of size (ndofs, ndofs) into target matrix answer.
Class and object Domain.
Definition: domain.h:115
virtual void computeBCVectorAt(FloatArray &answer, TimeStep *tStep, ValueModeType mode)
Computes the RHS contribution to balance equation(s) due to boundary conditions.
virtual void computeCapacityMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes the capacity matrix of the receiver.
Abstract class representing field of primary variables (those, which are unknown and are typically as...
Definition: primaryfield.h:104
virtual int checkConsistency()
Performs consistency check.
Transort cross-section.
ValueModeType
Type representing the mode of UnknownType or CharType, or similar types.
Definition: valuemodetype.h:78
Class implementing element body load, acting over whole element volume (e.g., the dead weight)...
Definition: bodyload.h:49
virtual void giveCharacteristicMatrix(FloatMatrix &answer, CharType type, TimeStep *tStep)
Computes characteristic matrix of receiver of requested type in given time step.
virtual void computeInternalSourceRhsVectorAt(FloatArray &answer, TimeStep *tStep, ValueModeType mode)
Computes the contribution to balance equation(s) due to internal sources.
virtual void computeExternalForcesVector(FloatArray &answer, TimeStep *tStep, ValueModeType mode)
virtual int giveApproxOrder(int unknownIndx)
Abstract base class for all finite elements.
Definition: element.h:145
virtual void computeBmatrixAt(FloatMatrix &answer, const FloatArray &lcoords)
virtual void computeIntSourceLHSMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes the LHS contribution to balance equation(s) due to material internal source.
virtual void computeNAt(FloatArray &answer, const FloatArray &lcoord)
Computes the basis functions.
Class implementing an array of integers.
Definition: intarray.h:61
MatResponseMode
Describes the character of characteristic material matrix.
virtual void computeTangentFromEdgeLoad(FloatMatrix &answer, EdgeLoad *load, int boundary, MatResponseMode rmode, TimeStep *tStep)
Computes the tangent contribution of the given load at the given boundary.
virtual void computeConductivityMatrix(FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep)
Computes the conductivity matrix of the receiver.
virtual void computeNmatrixAt(FloatMatrix &answer, const FloatArray &lcoords)
Computes the interpolation matrix corresponding to all unknowns.
virtual void computeGradientMatrixAt(FloatMatrix &answer, const FloatArray &lcoords)
Computes the gradient matrix corresponding to one unknown.
virtual void computeSurfaceNAt(FloatArray &answer, int iSurf, const FloatArray &lcoord)
virtual void computeInertiaForcesVector(FloatArray &answer, TimeStep *tStep)
Abstract base class representing integration rule.
virtual void computeCapacitySubMatrix(FloatMatrix &answer, MatResponseMode rmode, int iri, TimeStep *tStep)
Computes the matrix .
void computeBodyBCSubVectorAt(FloatArray &answer, Load *load, TimeStep *tStep, ValueModeType mode, int indx)
Element interface class.
Definition: primaryfield.h:58
virtual void computeInternalForcesVector(FloatArray &answer, TimeStep *tStep)
void computeSurfaceBCSubVectorAt(FloatArray &answer, Load *load, int iSurf, TimeStep *tStep, ValueModeType mode, int indx)
Computes the part of RHS due to applied BCs on particular surface.
virtual double computeEdgeVolumeAround(GaussPoint *gp, int iEdge)=0
Computes the length around a integration point on a edge.
virtual void computeInternalSourceRhsSubVectorAt(FloatArray &answer, TimeStep *, ValueModeType mode, int indx)
Computes the contribution to balance equation(s) due to internal sources.
Abstract base class representing a boundary load (force, momentum, ...) that acts directly on a bound...
Definition: boundaryload.h:110
virtual void giveSurfaceDofMapping(IntArray &mask, int iSurf)
This abstract class represent a general base element class for transport problems.
virtual void computeFlow(FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
Computes a flow vector in an integration point.
virtual int EIPrimaryFieldI_evaluateFieldVectorAt(FloatArray &answer, PrimaryField &pf, const FloatArray &coords, IntArray &dofId, ValueModeType mode, TimeStep *tStep)
Evaluates the value of field at given point of interest (should be located inside receiver's volume) ...
Abstract base class representing an edge load (force, momentum, ...) that acts directly on a edge bou...
Definition: boundaryload.h:200
virtual void computeField(ValueModeType mode, TimeStep *tStep, const FloatArray &lcoords, FloatArray &answer)
Computes the unknown vector interpolated at the specified local coordinates.
Abstract base class for all material models.
Definition: material.h:95
virtual void computeLumpedCapacityVector(FloatArray &answer, TimeStep *tStep)
void computeBCSubVectorAt(FloatArray &answer, TimeStep *tStep, ValueModeType mode, int indx)
Computes the part of RHS due to applied BCs.
Class representing vector of real numbers.
Definition: floatarray.h:82
virtual double computeSurfaceVolumeAround(GaussPoint *gp, int iSurf)
Abstract base class representing a surface load (force, momentum, ...) that acts directly on a surfac...
Definition: boundaryload.h:218
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
virtual void computeLoadVector(FloatArray &answer, BodyLoad *load, CharType type, ValueModeType mode, TimeStep *tStep)
Computes the contribution of the given body load (volumetric).
virtual void computeConductivitySubMatrix(FloatMatrix &answer, int iri, MatResponseMode rmode, TimeStep *tStep)
Computes the matrix .
virtual double giveThicknessAt(const FloatArray &gcoords)
Gives the thickness at some global coordinate.
CharType
Definition: chartype.h:87
virtual void computeIntSourceLHSSubMatrix(FloatMatrix &answer, MatResponseMode rmode, int iri, TimeStep *tStep)
Computes the part of internal source LHS contribution corresponding to unknown identified by rmode pa...
virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
Computes constitutive matrix of receiver.
static const double stefanBoltzmann
Stefan–Boltzmann constant W/m2/K4.
virtual void computeBCMtrxAt(FloatMatrix &answer, TimeStep *tStep, ValueModeType mode)
Computes the LHS contribution to balance equation(s) due to boundary conditions.
void computeEdgeBCSubVectorAt(FloatArray &answer, Load *load, int iEdge, TimeStep *tStep, ValueModeType mode, int indx)
Computes the part of RHS due to applied BCs on particular edge.
virtual IntegrationRule * GetSurfaceIntegrationRule(int approxOrder)
Load is base abstract class for all loads.
Definition: load.h:61
virtual void giveEdgeDofMapping(IntArray &mask, int iEdge)
Gives the node indexes for given edge.
TransportCrossSection * giveTransportCrossSection()
the oofem namespace is to define a context or scope in which all oofem names are defined.
virtual void giveCharacteristicVector(FloatArray &answer, CharType type, ValueModeType mode, TimeStep *tStep)
Computes characteristic vector of receiver of requested type in given time step.
virtual void computeBoundarySurfaceLoadVector(FloatArray &answer, BoundaryLoad *load, int boundary, CharType type, ValueModeType mode, TimeStep *tStep, bool global=true)
Computes the contribution of the given load at the given boundary surface in global coordinate system...
TransportElement(int n, Domain *d, ElementMode em=HeatTransferEM)
virtual Material * giveMaterial()
virtual void updateInternalState(TimeStep *tStep)
Updates element state after equilibrium in time step has been reached.
Class representing integration point in finite element program.
Definition: gausspoint.h:93
Class representing solution step.
Definition: timestep.h:80
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.
InternalStateMode
Determines the mode of internal variable.
virtual void giveDofManDofIDMask(int inode, IntArray &answer) const
Returns dofmanager dof mask for node.
void computeBCSubMtrxAt(FloatMatrix &answer, TimeStep *tStep, ValueModeType mode, int indx)
Computes the part of LHS due to applied BCs.
virtual void computeEgdeNAt(FloatArray &answer, int iEdge, const FloatArray &lcoord)
Computes the basis functions at the edge for one unknown.

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