OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
nodalspringelement.C
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 #include "../sm/Elements/nodalspringelement.h"
36 #include "floatmatrix.h"
37 #include "intarray.h"
38 #include "floatarray.h"
39 #include "classfactory.h"
40 
41 #ifdef __OOFEG
42  #include "oofeggraphiccontext.h"
43 #endif
44 
45 namespace oofem {
46 REGISTER_Element(NodalSpringElement);
47 
49 {
50  numberOfDofMans = 1;
51 }
52 
53 void
55 {
56  answer.beDiagonal(this->springConstants);
57 }
58 
59 
60 void
61 NodalSpringElement :: giveInternalForcesVector(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord)
62 {
63  int ndofs = computeNumberOfDofs();
64  FloatArray u;
65  this->computeVectorOf(VM_Total, tStep, u);
66 
67  answer.resize(ndofs);
68 
69  for (int i = 1; i<=ndofs; i++) {
70  answer.at(i) = this->springConstants.at(i)*u.at(i);
71  }
72 }
73 
74 
75 bool
77 {
78  return false;
79 }
80 
81 
82 void
84 {
85  answer= this->dofMask;
86 }
87 
88 void
90 {
91  if (this->masses.isNotEmpty()) {
92  answer.beDiagonal(masses);
93  } else {
94  int ndofs = this->computeNumberOfDofs();
95  answer.resize(ndofs, ndofs);
96  answer.zero();
97  }
98 }
99 
100 int
102 {
103  return this->computeNumberOfDofs();
104 }
105 
106 
109 {
110  IRResultType result; // Required by IR_GIVE_FIELD macro
111 
114 
115  this->masses.resize(dofMask.giveSize());
116  this->masses.zero();
118 
119  int ndofs = dofMask.giveSize();
120  if (ndofs != springConstants.giveSize()) {
121  OOFEM_ERROR ("Spring constants size not equal to number of DOFs");
122  }
123 
124  if (ndofs != masses.giveSize()) {
125  OOFEM_ERROR ("Masses array size not equal to number of DOFs");
126  }
127 
128  // from element
132 
133  return IRRT_OK;
134 }
135 
137 {
138  FloatArray F;
139  this->giveInternalForcesVector(F, tStep);
140 
141  fprintf(File, "NodalSpring element %d (%8d) :\n", this->giveLabel(), this->giveNumber() );
142  fprintf(File, " Force/Moment ");
143  for (int i = 1; i<= F.giveSize(); i++) {
144  fprintf(File, "%.4e ", F.at(i));
145  }
146  fprintf(File, "\n");
147 }
148 } // end namespace oofem
IntArray dofManArray
Array containing dofmanager numbers.
Definition: element.h:151
IntArray dofMask
Dof mask.
Class and object Domain.
Definition: domain.h:115
void computeVectorOf(ValueModeType u, TimeStep *tStep, FloatArray &answer)
Returns local vector of unknowns.
Definition: element.C:86
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
virtual bool computeGtoLRotationMatrix(FloatMatrix &answer)
Returns transformation matrix from global c.s.
FloatArray springConstants
Spring constants.
void beDiagonal(const FloatArray &diag)
Modifies receiver to be a diagonal matrix with the components specified in diag.
Definition: floatmatrix.C:1433
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Prints output of receiver to stream, for given time step.
Class implementing an array of integers.
Definition: intarray.h:61
MatResponseMode
Describes the character of characteristic material matrix.
#define _IFT_NodalSpringElement_dofmask
virtual int computeNumberOfDofs()
Computes or simply returns total number of element&#39;s local DOFs.
virtual void computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep)
Computes numerically stiffness matrix of receiver.
int activityTimeFunction
Element activity time function. If defined, nonzero value indicates active receiver, zero value inactive element.
Definition: element.h:176
Abstract base class for all "structural" finite elements.
virtual void computeLumpedMassMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes lumped mass matrix of receiver.
#define OOFEM_ERROR(...)
Definition: error.h:61
#define _IFT_NodalSpringElement_springConstants
REGISTER_Element(LSpace)
#define _IFT_Element_nodes
Definition: element.h:65
virtual void giveDofManDofIDMask(int inode, IntArray &answer) const
Returns dofmanager dof mask for node.
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
NodalSpringElement(int n, Domain *d)
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
void resize(int rows, int cols)
Checks size of receiver towards requested bounds.
Definition: floatmatrix.C:1358
Class representing the general Input Record.
Definition: inputrecord.h:101
void zero()
Zeroes all coefficients of receiver.
Definition: floatarray.C:658
FloatArray masses
total mass of the spring; to be distributed to nodes
#define _IFT_Element_activityTimeFunction
Definition: element.h:71
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
#define _IFT_NodalSpringElement_masses
#define IR_GIVE_OPTIONAL_FIELD(__ir, __value, __id)
Macro facilitating the use of input record reading methods.
Definition: inputrecord.h:78
int giveSize() const
Definition: intarray.h:203
int giveSize() const
Returns the size of receiver.
Definition: floatarray.h:218
int giveLabel() const
Definition: element.h:1055
the oofem namespace is to define a context or scope in which all oofem names are defined.
bool isNotEmpty() const
Returns true if receiver is not empty.
Definition: floatarray.h:220
#define IR_GIVE_FIELD(__ir, __value, __id)
Macro facilitating the use of input record reading methods.
Definition: inputrecord.h:69
int giveNumber() const
Definition: femcmpnn.h:107
virtual void giveInternalForcesVector(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord=0)
Returns equivalent nodal forces vectors.
Class representing solution step.
Definition: timestep.h:80
int numberOfDofMans
Number of dofmanagers.
Definition: element.h:149
void resize(int s)
Resizes receiver towards requested size.
Definition: floatarray.C:631
virtual int computeNumberOfGlobalDofs()
Computes the total number of element&#39;s global dofs.

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