OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
tr1bubblestokes.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 "tr1bubblestokes.h"
36 #include "node.h"
37 #include "elementinternaldofman.h"
38 #include "domain.h"
39 #include "gaussintegrationrule.h"
40 #include "gausspoint.h"
41 #include "bcgeomtype.h"
42 #include "load.h"
43 #include "bodyload.h"
44 #include "boundaryload.h"
45 #include "mathfem.h"
46 #include "fluiddynamicmaterial.h"
47 #include "fei2dtrlin.h"
48 #include "masterdof.h"
49 #include "fluidcrosssection.h"
50 #include "classfactory.h"
51 
52 namespace oofem {
53 REGISTER_Element(Tr1BubbleStokes);
54 
55 // Set up interpolation coordinates
56 FEI2dTrLin Tr1BubbleStokes :: interp(1, 2);
57 // Set up ordering vectors (for assembling)
58 IntArray Tr1BubbleStokes :: momentum_ordering = {1, 2, 4, 5, 7, 8, 10, 11};
59 IntArray Tr1BubbleStokes :: conservation_ordering = {3, 6, 9};
60 IntArray Tr1BubbleStokes :: edge_ordering [ 3 ] = {
61  {1, 2, 4, 5},
62  {4, 5, 7, 8},
63  {7, 8, 1, 2}
64 };
65 
67 {
68  this->numberOfDofMans = 3;
69  this->numberOfGaussPoints = 7;
70 
71  this->bubble.reset( new ElementDofManager(1, aDomain, this) );
72  this->bubble->appendDof( new MasterDof(this->bubble.get(), V_u) );
73  this->bubble->appendDof( new MasterDof(this->bubble.get(), V_v) );
74 }
75 
77 {
78 }
79 
81 {
82  if ( integrationRulesArray.size() == 0 ) {
83  integrationRulesArray.resize(1);
84  integrationRulesArray [ 0 ].reset( new GaussIntegrationRule(1, this, 1, 3) );
86  }
87 }
88 
90 {
91  return 11;
92 }
93 
94 void Tr1BubbleStokes :: giveDofManDofIDMask(int inode, IntArray &answer) const
95 {
96  answer = {V_u, V_v, P_f};
97 }
98 
100 {
101  answer = {V_u, V_v};
102 }
103 
105 {
106  return 1;
107 
108 }
109 
111 {
112  return this->bubble.get();
113 }
114 
116 {
117  double detJ = fabs( this->interp.giveTransformationJacobian( gp->giveNaturalCoordinates(), FEIElementGeometryWrapper(this) ) );
118  return detJ *gp->giveWeight();
119 }
120 
122  TimeStep *tStep)
123 {
124  // Compute characteristic vector for this element. I.e the load vector(s)
125  if ( mtrx == ExternalForcesVector ) {
126  this->computeExternalForcesVector(answer, tStep);
127  } else if ( mtrx == InternalForcesVector ) {
128  this->computeInternalForcesVector(answer, tStep);
129  } else {
130  OOFEM_ERROR("Unknown Type of characteristic mtrx.");
131  }
132 }
133 
135  CharType mtrx, TimeStep *tStep)
136 {
137  if ( mtrx == TangentStiffnessMatrix ) {
138  this->computeStiffnessMatrix(answer, TangentStiffness, tStep);
139  } else {
140  OOFEM_ERROR("Unknown Type of characteristic mtrx.");
141  }
142 }
143 
145 {
146  FluidDynamicMaterial *mat = static_cast< FluidCrossSection * >( this->giveCrossSection() )->giveFluidMaterial();
147  FloatArray a_pressure, a_velocity, devStress, epsp, N, dNv(8);
148  double r_vol, pressure;
149  FloatMatrix dN, B(3, 8);
150  B.zero();
151 
152  this->computeVectorOfVelocities(VM_Total, tStep, a_velocity);
153  this->computeVectorOfPressures(VM_Total, tStep, a_pressure);
154 
155  FloatArray momentum, conservation;
156 
157  for ( GaussPoint *gp: *integrationRulesArray [ 0 ] ) {
158  const FloatArray &lcoords = gp->giveNaturalCoordinates();
159 
160  double detJ = fabs( this->interp.evaldNdx( dN, lcoords, FEIElementGeometryWrapper(this) ) );
161  this->interp.evalN( N, lcoords, FEIElementGeometryWrapper(this) );
162  double dA = detJ * gp->giveWeight();
163 
164  for ( int j = 0, k = 0; j < 3; j++, k += 2 ) {
165  dNv(k) = B(0, k) = B(2, k + 1) = dN(j, 0);
166  dNv(k + 1) = B(1, k + 1) = B(2, k) = dN(j, 1);
167  }
168 
169  // Bubble contribution;
170  dNv(6) = B(0, 6) = B(2, 7) = 27. * ( dN(0, 0) * N(1) * N(2) + N(0) * dN(1, 0) * N(2) + N(0) * N(1) * dN(2, 0) );
171  dNv(7) = B(1, 7) = B(2, 6) = 27. * ( dN(0, 1) * N(1) * N(2) + N(0) * dN(1, 1) * N(2) + N(0) * N(1) * dN(2, 1) );
172 
173  pressure = N.dotProduct(a_pressure);
174  epsp.beProductOf(B, a_velocity);
175 
176  mat->computeDeviatoricStressVector(devStress, r_vol, gp, epsp, pressure, tStep);
177 
178  momentum.plusProduct(B, devStress, dA);
179  momentum.add(-pressure * dA, dNv);
180  conservation.add(r_vol * dA, N);
181  }
182 
183  answer.resize(11);
184  answer.zero();
185  answer.assemble(momentum, this->momentum_ordering);
186  answer.assemble(conservation, this->conservation_ordering);
187 }
188 
190 {
191  int load_number, load_id;
192  Load *load;
193  BodyLoad *bLoad;
194  bcGeomType ltype;
195  FloatArray vec;
196 
197  int nLoads = this->boundaryLoadArray.giveSize() / 2;
198  answer.resize(11);
199  answer.zero();
200  for ( int i = 1; i <= nLoads; i++ ) { // For each Neumann boundary condition
201  load_number = this->boundaryLoadArray.at(2 * i - 1);
202  load_id = this->boundaryLoadArray.at(2 * i);
203  load = this->domain->giveLoad(load_number);
204  ltype = load->giveBCGeoType();
205 
206  if ( ltype == EdgeLoadBGT ) {
207  this->computeBoundarySurfaceLoadVector(vec, static_cast< BoundaryLoad * >(load), load_id, ExternalForcesVector, VM_Total, tStep);
208  answer.add(vec);
209  }
210  }
211 
212  nLoads = this->giveBodyLoadArray()->giveSize();
213  for ( int i = 1; i <= nLoads; i++ ) {
214  load = domain->giveLoad( bodyLoadArray.at(i) );
215  if ((bLoad = dynamic_cast<BodyLoad*>(load))) {
216  ltype = load->giveBCGeoType();
217  if ( ltype == BodyLoadBGT && load->giveBCValType() == ForceLoadBVT ) {
218  this->computeLoadVector(vec, bLoad, ExternalForcesVector, VM_Total, tStep);
219  answer.add(vec);
220  }
221  }
222  }
223 }
224 
225 
227 {
228  if ( type != ExternalForcesVector ) {
229  answer.clear();
230  return;
231  }
232 
233  FloatArray N, gVector, temparray(8);
234 
235  load->computeComponentArrayAt(gVector, tStep, VM_Total);
236  temparray.zero();
237  if ( gVector.giveSize() ) {
238  for ( GaussPoint *gp: *integrationRulesArray [ 0 ] ) {
239  const FloatArray &lcoords = gp->giveNaturalCoordinates();
240 
241  double rho = static_cast< FluidCrossSection * >( this->giveCrossSection() )->giveDensity(gp);
242  double detJ = fabs( this->interp.giveTransformationJacobian( lcoords, FEIElementGeometryWrapper(this) ) );
243  double dA = detJ * gp->giveWeight();
244 
245  this->interp.evalN( N, lcoords, FEIElementGeometryWrapper(this) );
246  for ( int j = 0; j < 3; j++ ) {
247  temparray(2 * j) += N(j) * rho * gVector(0) * dA;
248  temparray(2 * j + 1) += N(j) * rho * gVector(1) * dA;
249  }
250 
251  temparray(6) += N(0) * N(1) * N(2) * rho * gVector(0) * dA;
252  temparray(7) += N(0) * N(1) * N(2) * rho * gVector(1) * dA;
253  }
254  }
255 
256  answer.resize(11);
257  answer.zero();
258  answer.assemble(temparray, this->momentum_ordering);
259 }
260 
261  void Tr1BubbleStokes :: computeBoundarySurfaceLoadVector(FloatArray &answer, BoundaryLoad *load, int iEdge, CharType type, ValueModeType mode, TimeStep *tStep, bool global)
262 {
263  if ( type != ExternalForcesVector ) {
264  answer.clear();
265  return;
266  }
267 
268  if ( load->giveType() == TransmissionBC ) { // Neumann boundary conditions (traction)
269 
270  int numberOfEdgeIPs = ( int ) ceil( ( load->giveApproxOrder() + 2. ) / 2. );
271 
272  GaussIntegrationRule iRule(1, this, 1, 1);
273  FloatArray N, t, f(4);
274 
275  f.zero();
276  iRule.SetUpPointsOnLine(numberOfEdgeIPs, _Unknown);
277 
278  for ( GaussPoint *gp: iRule ) {
279  const FloatArray &lcoords = gp->giveNaturalCoordinates();
280 
281  this->interp.edgeEvalN( N, iEdge, lcoords, FEIElementGeometryWrapper(this) );
282  double detJ = fabs( this->interp.boundaryGiveTransformationJacobian( iEdge, lcoords, FEIElementGeometryWrapper(this) ) );
283  double dS = gp->giveWeight() * detJ;
284 
285  if ( load->giveFormulationType() == Load :: FT_Entity ) { // Edge load in xi-eta system
286  load->computeValueAt(t, tStep, lcoords, VM_Total);
287  } else { // Edge load in x-y system
288  FloatArray gcoords;
289  this->interp.boundaryLocal2Global( gcoords, iEdge, lcoords, FEIElementGeometryWrapper(this) );
290  load->computeValueAt(t, tStep, gcoords, VM_Total);
291  }
292 
293  // Reshape the vector
294  for ( int j = 0; j < 2; j++ ) {
295  f(2 * j) += N(j) * t(0) * dS;
296  f(2 * j + 1) += N(j) * t(1) * dS;
297  }
298  }
299 
300  answer.resize(11);
301  answer.zero();
302  answer.assemble(f, this->edge_ordering [ iEdge - 1 ]);
303  } else {
304  OOFEM_ERROR("Strange boundary condition type");
305  }
306 }
307 
309 {
310  // Note: Working with the components; [K, G+Dp; G^T+Dv^T, C] . [v,p]
311  FluidDynamicMaterial *mat = static_cast< FluidCrossSection * >( this->giveCrossSection() )->giveFluidMaterial();
312  FloatMatrix B(3, 8), EdB, K, G, Dp, DvT, C, Ed, dN;
313  FloatArray dNv(8), N, Ep, Cd, tmpA, tmpB;
314  double Cp;
315  B.zero();
316 
317  for ( GaussPoint *gp: *integrationRulesArray [ 0 ] ) {
318  // Compute Gauss point and determinant at current element
319  const FloatArray &lcoords = gp->giveNaturalCoordinates();
320 
321  double detJ = fabs( this->interp.evaldNdx( dN, lcoords, FEIElementGeometryWrapper(this) ) );
322  double dA = detJ * gp->giveWeight();
323  this->interp.evalN( N, lcoords, FEIElementGeometryWrapper(this) );
324  for ( int j = 0, k = 0; j < 3; j++, k += 2 ) {
325  dNv(k) = B(0, k) = B(2, k + 1) = dN(j, 0);
326  dNv(k + 1) = B(1, k + 1) = B(2, k) = dN(j, 1);
327  }
328 
329  // Bubble contribution;
330  dNv(6) = B(0, 6) = B(2, 7) = 27. * ( dN(0, 0) * N(1) * N(2) + N(0) * dN(1, 0) * N(2) + N(0) * N(1) * dN(2, 0) );
331  dNv(7) = B(1, 7) = B(2, 6) = 27. * ( dN(0, 1) * N(1) * N(2) + N(0) * dN(1, 1) * N(2) + N(0) * N(1) * dN(2, 1) );
332 
333  // Computing the internal forces should have been done first.
334  // dsigma_dev/deps_dev dsigma_dev/dp deps_vol/deps_dev deps_vol/dp
335  mat->giveStiffnessMatrices(Ed, Ep, Cd, Cp, mode, gp, tStep);
336 
337  EdB.beProductOf(Ed, B);
338  K.plusProductSymmUpper(B, EdB, dA);
339  G.plusDyadUnsym(dNv, N, -dA);
340  C.plusDyadSymmUpper(N, Cp * dA);
341 
342  tmpA.beTProductOf(B, Ep);
343  Dp.plusDyadUnsym(tmpA, N, dA);
344 
345  tmpB.beTProductOf(B, Cd);
346  DvT.plusDyadUnsym(N, tmpB, dA);
347  }
348 
349  K.symmetrized();
350  C.symmetrized();
351 
352  FloatMatrix GTDvT, GDp;
353  GTDvT.beTranspositionOf(G);
354  GTDvT.add(DvT);
355  GDp = G;
356  GDp.add(Dp);
357 
358  answer.resize(11, 11);
359  answer.zero();
360  answer.assemble(K, this->momentum_ordering);
361  answer.assemble(GTDvT, this->conservation_ordering, this->momentum_ordering);
362  answer.assemble(GDp, this->momentum_ordering, this->conservation_ordering);
363  answer.assemble(C, this->conservation_ordering);
364 }
365 
367 {
368  return & interp;
369 }
370 
372 {
373  return & interp;
374 }
375 
377 {
379 }
380 
381 // Some extension Interfaces to follow:
382 
384 {
385  switch ( it ) {
387  return static_cast< ZZNodalRecoveryModelInterface * >(this);
388 
390  return static_cast< SpatialLocalizerInterface * >(this);
391 
392  default:
393  return FMElement :: giveInterface(it);
394  }
395 }
396 
398 {
399  FloatArray n, n_lin, pressures, velocities;
400  this->interp.evalN( n, lcoords, FEIElementGeometryWrapper(this) );
401  this->interp.evalN( n_lin, lcoords, FEIElementGeometryWrapper(this) );
402  this->computeVectorOf({P_f}, mode, tStep, pressures);
403  this->computeVectorOf({V_u, V_v}, mode, tStep, velocities);
404 
405  answer.resize(3);
406  answer.zero();
407  for ( int i = 1; i <= n.giveSize(); i++ ) {
408  answer(0) += n.at(i) * velocities.at(i*2-1);
409  answer(1) += n.at(i) * velocities.at(i*2);
410  }
411  answer(0) += n.at(1) * n.at(2) * n.at(3) * velocities.at(7);
412  answer(1) += n.at(1) * n.at(2) * n.at(3) * velocities.at(8);
413 
414  for ( int i = 1; i <= n_lin.giveSize(); i++ ) {
415  answer(2) += n_lin.at(i) * pressures.at(i);
416  }
417 }
418 
419 } // end namespace oofem
CrossSection * giveCrossSection()
Definition: element.C:495
virtual void edgeEvalN(FloatArray &answer, int iedge, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the array of edge interpolation functions (shape functions) at given point.
Definition: fei2dtrlin.C:184
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...
virtual double evaldNdx(FloatMatrix &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the matrix of derivatives of interpolation functions (shape functions) at given point...
Definition: fei2dtrlin.C:53
Class and object Domain.
Definition: domain.h:115
virtual void evalN(FloatArray &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the array of interpolation functions (shape functions) at given point.
Definition: fei2dtrlin.C:43
void computeVectorOf(ValueModeType u, TimeStep *tStep, FloatArray &answer)
Returns local vector of unknowns.
Definition: element.C:86
Fluid cross-section.
static FEI2dTrLin interp
Interpolation for pressure.
Abstract base class for all fluid materials.
Domain * domain
Link to domain object, useful for communicating with other FEM components.
Definition: femcmpnn.h:82
std::unique_ptr< ElementDofManager > bubble
The extra dofs from the bubble.
Class implementing internal element dof manager having some DOFs.
The element interface required by ZZNodalRecoveryModel.
bcGeomType
Type representing the geometric character of loading.
Definition: bcgeomtype.h:40
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
virtual void computeLoadVector(FloatArray &answer, BodyLoad *load, CharType type, ValueModeType mode, TimeStep *tStep)
Computes the contribution of the given body load (volumetric).
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
void computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode mode, TimeStep *tStep)
virtual FEInterpolation * giveInterpolation() const
void clear()
Clears receiver (zero size).
Definition: floatarray.h:206
virtual void updateYourself(TimeStep *tStep)
Updates element state after equilibrium in time step has been reached.
Definition: element.C:779
Base class for dof managers.
Definition: dofmanager.h:113
virtual void computeComponentArrayAt(FloatArray &answer, TimeStep *tStep, ValueModeType mode)
Computes boundary condition value - its components values at given time.
Definition: load.C:82
void plusDyadUnsym(const FloatArray &a, const FloatArray &b, double dV)
Adds to the receiver the product .
Definition: floatmatrix.C:811
virtual void giveStiffnessMatrices(FloatMatrix &dsdd, FloatArray &dsdp, FloatArray &dedd, double &dedp, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
Computes the 4 tangents of the compressible material response in 3D.
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 void updateYourself(TimeStep *tStep)
Updates element state after equilibrium in time step has been reached.
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
Distributed body load.
Definition: bcgeomtype.h:43
Class representing a general abstraction for finite element interpolation class.
Definition: feinterpol.h:132
Class representing "master" degree of freedom.
Definition: masterdof.h:92
virtual void computeVectorOfPressures(ValueModeType mode, TimeStep *tStep, FloatArray &pressures)
Definition: fmelement.C:55
void computeInternalForcesVector(FloatArray &answer, TimeStep *tStep)
static IntArray conservation_ordering
static IntArray edge_ordering[3]
Ordering of dofs on edges. Used to assemble edge loads.
virtual void computeDeviatoricStressVector(FloatArray &stress_dev, double &epsp_vol, GaussPoint *gp, const FloatArray &eps, double pressure, TimeStep *tStep)
Computes the deviatoric stress vector and volumetric strain rate from given deviatoric strain and pre...
Distributed edge load.
Definition: bcgeomtype.h:44
Abstract base class representing a boundary load (force, momentum, ...) that acts directly on a bound...
Definition: boundaryload.h:110
virtual Interface * giveInterface(InterfaceType it)
Interface requesting service.
virtual void computeVectorOfVelocities(ValueModeType mode, TimeStep *tStep, FloatArray &velocities)
Definition: fmelement.C:48
double dotProduct(const FloatArray &x) const
Computes the dot product (or inner product) of receiver and argument.
Definition: floatarray.C:463
#define OOFEM_ERROR(...)
Definition: error.h:61
virtual void giveDofManDofIDMask(int inode, IntArray &answer) const
Returns dofmanager dof mask for node.
REGISTER_Element(LSpace)
void computeExternalForcesVector(FloatArray &answer, TimeStep *tStep)
DofIDItem
Type representing particular dof type.
Definition: dofiditem.h:86
virtual double giveWeight()
Returns integration weight of receiver.
Definition: gausspoint.h:181
virtual FormulationType giveFormulationType()
Specifies is load should take local or global coordinates.
Definition: load.h:155
Neumann type (prescribed flux).
Definition: bctype.h:43
virtual double boundaryGiveTransformationJacobian(int boundary, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the determinant of the transformation Jacobian on the requested boundary.
Definition: feinterpol2d.C:154
Wrapper around element definition to provide FEICellGeometry interface.
Definition: feinterpol.h:95
IntArray * giveBodyLoadArray()
Returns array containing load numbers of loads acting on element.
Definition: element.C:372
void beProductOf(const FloatMatrix &aMatrix, const FloatArray &anArray)
Receiver becomes the result of the product of aMatrix and anArray.
Definition: floatarray.C:676
IntArray bodyLoadArray
Array containing indexes of loads (body loads and boundary loads are kept separately), that apply on receiver.
Definition: element.h:160
#define N(p, q)
Definition: mdm.C:367
void plusProductSymmUpper(const FloatMatrix &a, const FloatMatrix &b, double dV)
Adds to the receiver the product .
Definition: floatmatrix.C:698
virtual int giveApproxOrder()=0
void beTProductOf(const FloatMatrix &aMatrix, const FloatArray &anArray)
Receiver becomes the result of the product of aMatrix^T and anArray.
Definition: floatarray.C:708
virtual bcType giveType() const
Returns receiver load type.
Definition: boundaryload.h:166
virtual int setupIntegrationPoints(IntegrationRule &irule, int npoints, Element *element)
Sets up integration rule for the given element.
Definition: crosssection.C:54
int numberOfGaussPoints
Number of integration points as specified by nip.
Definition: element.h:188
Class representing vector of real numbers.
Definition: floatarray.h:82
void plusDyadSymmUpper(const FloatArray &a, double dV)
Adds to the receiver the dyadic product .
Definition: floatmatrix.C:756
virtual void giveCharacteristicMatrix(FloatMatrix &answer, CharType type, TimeStep *tStep)
Computes characteristic matrix of receiver of requested type in given time step.
This abstract class represent a general base element class for fluid dynamic problems.
Definition: fmelement.h:54
virtual void computeField(ValueModeType u, TimeStep *tStep, const FloatArray &coords, FloatArray &answer)
Computes the unknown vector interpolated at the specified local coordinates.
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
virtual double giveTransformationJacobian(const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the determinant of the transformation.
Definition: fei2dtrlin.C:147
virtual int computeNumberOfDofs()
Computes or simply returns total number of element&#39;s local DOFs.
void assemble(const FloatArray &fe, const IntArray &loc)
Assembles the array fe (typically, the load vector of a finite element) into the receiver, using loc as location array.
Definition: floatarray.C:551
static IntArray momentum_ordering
Ordering of dofs in element. Used to assemble the element stiffness.
CharType
Definition: chartype.h:87
void resize(int rows, int cols)
Checks size of receiver towards requested bounds.
Definition: floatmatrix.C:1358
void add(const FloatMatrix &a)
Adds matrix to the receiver.
Definition: floatmatrix.C:1023
Class Interface.
Definition: interface.h:82
void zero()
Zeroes all coefficients of receiver.
Definition: floatarray.C:658
virtual void giveInternalDofManDofIDMask(int i, IntArray &answer) const
Returns internal dofmanager dof mask for node.
virtual bcGeomType giveBCGeoType() const
Returns geometry character of boundary condition.
virtual int giveNumberOfInternalDofManagers() const
virtual Interface * giveInterface(InterfaceType t)
Interface requesting service.
Definition: femcmpnn.h:179
void beTranspositionOf(const FloatMatrix &src)
Assigns to the receiver the transposition of parameter.
Definition: floatmatrix.C:323
The spatial localizer element interface associated to spatial localizer.
std::vector< std::unique_ptr< IntegrationRule > > integrationRulesArray
List of integration rules of receiver (each integration rule contains associated integration points a...
Definition: element.h:170
virtual void giveCharacteristicVector(FloatArray &answer, CharType type, ValueModeType mode, TimeStep *tStep)
Computes characteristic vector of receiver of requested type in given time step.
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
InterfaceType
Enumerative type, used to identify interface type.
Definition: interfacetype.h:43
Tr1BubbleStokes(int n, Domain *d)
virtual DofManager * giveInternalDofManager(int i) const
Returns i-th internal element dof manager of the receiver.
Load is base abstract class for all loads.
Definition: load.h:61
void beProductOf(const FloatMatrix &a, const FloatMatrix &b)
Assigns to the receiver product of .
Definition: floatmatrix.C:337
int giveSize() const
Definition: intarray.h:203
int giveSize() const
Returns the size of receiver.
Definition: floatarray.h:218
the oofem namespace is to define a context or scope in which all oofem names are defined.
void assemble(const FloatMatrix &src, const IntArray &loc)
Assembles the contribution using localization array into receiver.
Definition: floatmatrix.C:251
virtual bcValType giveBCValType() const
Returns receiver load type.
virtual void boundaryLocal2Global(FloatArray &answer, int boundary, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Maps the local boundary coordinates to global.
Definition: feinterpol2d.C:159
virtual void computeValueAt(FloatArray &answer, TimeStep *tStep, const FloatArray &coords, ValueModeType mode)
Computes components values of load at given point - global coordinates (coordinates given)...
Definition: boundaryload.C:58
Load * giveLoad(int n)
Service for accessing particular domain load.
Definition: domain.C:222
virtual int SetUpPointsOnLine(int nPoints, MaterialMode mode)
Sets up receiver&#39;s integration points on unit line integration domain.
void symmetrized()
Initializes the lower half of the receiver according to the upper half.
Definition: floatmatrix.C:1576
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
Class representing integration point in finite element program.
Definition: gausspoint.h:93
IntArray boundaryLoadArray
Definition: element.h:160
Class representing solution step.
Definition: timestep.h:80
int numberOfDofMans
Number of dofmanagers.
Definition: element.h:149
void add(const FloatArray &src)
Adds array src to receiver.
Definition: floatarray.C:156
const FloatArray & giveNaturalCoordinates()
Returns coordinate array of receiver.
Definition: gausspoint.h:138
Class representing Gaussian-quadrature integration rule.
void resize(int s)
Resizes receiver towards requested size.
Definition: floatarray.C:631
void plusProduct(const FloatMatrix &b, const FloatArray &s, double dV)
Adds the product .
Definition: floatarray.C:226

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