OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
quasicontinuumvtkxmlexportmodule.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 
36 #include "vtkxmlexportmodule.h"
37 #include "element.h"
38 #include "gausspoint.h"
39 #include "timestep.h"
40 #include "engngm.h"
41 #include "node.h"
42 #include "dof.h"
43 #include "materialinterface.h"
44 #include "mathfem.h"
45 #include "cltypes.h"
46 #include "material.h"
47 #include "classfactory.h"
48 #include "crosssection.h"
49 #include "unknownnumberingscheme.h"
50 
51 
52 namespace oofem {
54 
55 
56 QuasicontinuumVTKXMLExportModule :: QuasicontinuumVTKXMLExportModule(int n, EngngModel *e) : VTKXMLExportModule(n, e), internalVarsToExport(), primaryVarsToExport()
57 {}
58 
59 
61 {}
62 
63 
66 {
67  IRResultType result; // Required by IR_GIVE_FIELD macro
68 
71 
73 }
74 
75 
76 
77 void
79 {
80  // Stores all neccessary data (of a region) in a VTKPiece so it can be exported later.
81 
82  Domain *d = emodel->giveDomain(1);
83  Element *elem;
84  FloatArray *coords;
85 
86  this->giveSmoother(); // make sure smoother is created
87 
88  // output nodes Region By Region
89  int numNodes, numRegionEl;
90  IntArray mapG2L, mapL2G;
91 
92  // Assemble local->global and global->local region map and get number of
93  // single cells to process, the composite cells exported individually.
94  this->initRegionNodeNumbering(mapG2L, mapL2G, numNodes, numRegionEl, d, tStep, region);
95  if ( numNodes > 0 && numRegionEl > 0 ) {
96  // Export nodes as vtk vertices
97  vtkPiece.setNumberOfNodes(numNodes);
98  for ( int inode = 1; inode <= numNodes; inode++ ) {
99  coords = d->giveNode( mapL2G.at(inode) )->giveCoordinates();
100  vtkPiece.setNodeCoords(inode, * coords);
101  }
102 
103 
104  //-------------------------------------------
105  // Export all the cell data for the piece
106  //-------------------------------------------
107  IntArray cellNodes;
108  vtkPiece.setNumberOfCells(numRegionEl);
109 
110  int offset = 0;
111  int cellNum = 0;
112  IntArray elems = this->giveRegionSet(region)->giveElementList();
113  for ( int ei = 1; ei <= elems.giveSize(); ei++ ) {
114  int elNum = elems.at(ei);
115  elem = d->giveElement(elNum);
116 
117 
119  // Skip elements that are deactivated by QC
120  if ( !emodel->isElementActivated(elem) ) {
121  continue;
122  }
123  }
124 
125  // Skip elements that:
126  // are inactivated or of composite type ( these are exported individually later)
127  if ( this->isElementComposite(elem) || !elem->isActivated(tStep) ) {
128  continue;
129  }
130 
131  //skip materials with casting time > current time
132  if ( !elem->isCast(tStep) ) {
133  continue;
134  }
135 
136  if ( elem->giveParallelMode() != Element_local ) {
137  continue;
138  }
139 
140  cellNum++;
141 
142  // Set the connectivity
143  this->giveElementCell(cellNodes, elem); // node numbering of the cell with according to the VTK format
144 
145  // Map from global to local node numbers for the current piece
146  int numElNodes = cellNodes.giveSize();
147  IntArray connectivity(numElNodes);
148  for ( int i = 1; i <= numElNodes; i++ ) {
149  connectivity.at(i) = mapG2L.at( cellNodes.at(i) );
150  }
151 
152  vtkPiece.setConnectivity(cellNum, connectivity);
153 
154  vtkPiece.setCellType( cellNum, this->giveCellType(elem) ); // VTK cell type
155 
156  offset += numElNodes;
157  vtkPiece.setOffset(cellNum, offset);
158  }
159 
160 
161  // Export primary, internal and XFEM variables as nodal quantities
162  this->exportPrimaryVars(vtkPiece, mapG2L, mapL2G, region, tStep);
163  this->exportIntVars(vtkPiece, mapG2L, mapL2G, region, tStep);
164  this->exportExternalForces(vtkPiece, mapG2L, mapL2G, region, tStep);
165 
166  this->exportCellVars(vtkPiece, elems, tStep);
167  } // end of default piece for simple geometry elements
168 }
169 
170 
171 
172 
173 
174 int
176  IntArray &regionL2GNodalNumbers,
177  int &regionDofMans,
178  int &regionSingleCells,
179  Domain *domain, TimeStep *tStep, int reg)
180 {
181  // regionG2LNodalNumbers is array with mapping from global numbering to local region numbering.
182  // The i-th value contains the corresponding local region number (or zero, if global numbar is not in region).
183 
184  // regionL2GNodalNumbers is array with mapping from local to global numbering.
185  // The i-th value contains the corresponding global node number.
186 
187 
188  int nnodes = domain->giveNumberOfDofManagers();
189  int elemNodes;
190  int elementNode, node;
191  int currOffset = 1;
192  Element *element;
193 
194  regionG2LNodalNumbers.resize(nnodes);
195  regionG2LNodalNumbers.zero();
196  regionDofMans = 0;
197  regionSingleCells = 0;
198 
199  IntArray elements = this->giveRegionSet(reg)->giveElementList();
200  for ( int ie = 1; ie <= elements.giveSize(); ie++ ) {
201  int ielem = elements.at(ie);
202  element = domain->giveElement(ielem);
203 
204 
206  if ( !domain->giveEngngModel()->isElementActivated(element) ) {
207  continue; // skip elements deactovated by QC
208  }
209  }
210 
211  if ( this->isElementComposite(element) ) {
212  continue; // composite cells exported individually
213  }
214 
215  if ( !element->isActivated(tStep) ) { //skip inactivated elements
216  continue;
217  }
218 
219  //skip materials with casting time > current time
220  if ( !element->isCast(tStep) ) {
221  continue;
222  }
223 
224  if ( element->giveParallelMode() != Element_local ) {
225  continue;
226  }
227 
228  regionSingleCells++;
229  elemNodes = element->giveNumberOfNodes();
230  // elemSides = element->giveNumberOfSides();
231 
232  // determine local region node numbering
233  for ( elementNode = 1; elementNode <= elemNodes; elementNode++ ) {
234  node = element->giveNode(elementNode)->giveNumber();
235  if ( regionG2LNodalNumbers.at(node) == 0 ) { // assign new number
236  /* mark for assignment. This is done later, as it allows to preserve
237  * natural node numbering.
238  */
239  regionG2LNodalNumbers.at(node) = 1;
240  regionDofMans++;
241  }
242  }
243  }
244 
245  regionL2GNodalNumbers.resize(regionDofMans);
246 
247  for ( int i = 1; i <= nnodes; i++ ) {
248  if ( regionG2LNodalNumbers.at(i) ) {
249  regionG2LNodalNumbers.at(i) = currOffset++;
250  regionL2GNodalNumbers.at( regionG2LNodalNumbers.at(i) ) = i;
251  }
252  }
253 
254  return 1;
255 }
256 } // end namespace oofem
#define _IFT_QuasicontinuumVTKXMLExportModule_ExportDeactivatedElements
virtual bool isActivated(TimeStep *tStep)
Definition: element.C:793
void exportExternalForces(VTKPiece &piece, IntArray &mapG2L, IntArray &mapL2G, int region, TimeStep *tStep)
Export external forces.
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Class and object Domain.
Definition: domain.h:115
bool isElementComposite(Element *elem)
void exportCellVars(VTKPiece &piece, const IntArray &elems, TimeStep *tStep)
virtual int initRegionNodeNumbering(IntArray &mapG2L, IntArray &mapL2G, int &regionDofMans, int &totalcells, Domain *domain, TimeStep *tStep, int reg)
Assembles the region node map.
int giveNumberOfDofManagers() const
Returns number of dof managers in domain.
Definition: domain.h:432
int giveCellType(Element *element)
Returns corresponding element cell_type.
void setNumberOfNodes(int numNodes)
void zero()
Sets all component to zero.
Definition: intarray.C:52
EngngModel * giveEngngModel()
Returns engineering model to which receiver is associated.
Definition: domain.C:433
Abstract base class for all finite elements.
Definition: element.h:145
virtual void setupVTKPiece(VTKPiece &vtkPiece, TimeStep *tStep, int region)
void exportIntVars(VTKPiece &piece, IntArray &mapG2L, IntArray &mapL2G, int ireg, TimeStep *tStep)
Export internal variables by smoothing.
Class implementing an array of integers.
Definition: intarray.h:61
int & at(int i)
Coefficient access function.
Definition: intarray.h:103
void setNumberOfCells(int numCells)
void setNodeCoords(int nodeNum, FloatArray &coords)
virtual int giveNumberOfNodes() const
Returns number of nodes of receiver.
Definition: element.h:662
void giveElementCell(IntArray &answer, Element *elem)
Returns the element cell geometry.
void setCellType(int cellNum, int type)
Element * giveElement(int n)
Service for accessing particular domain fe element.
Definition: domain.C:160
EngngModel * emodel
Problem pointer.
Definition: exportmodule.h:77
Set * giveRegionSet(int i)
Returns element set.
Definition: exportmodule.C:109
void resize(int n)
Checks size of receiver towards requested bounds.
Definition: intarray.C:124
Represents VTK (Visualization Toolkit) export module.
Class representing vector of real numbers.
Definition: floatarray.h:82
elementParallelMode giveParallelMode() const
Return elementParallelMode of receiver.
Definition: element.h:1069
Element is local, there are no contributions from other domains to this element.
Definition: element.h:101
IRResultType
Type defining the return values of InputRecord reading operations.
Definition: irresulttype.h:47
NodalRecoveryModel * giveSmoother()
Returns the internal smoother.
void setOffset(int cellNum, int offset)
virtual bool isElementActivated(int elemNum)
Definition: engngm.h:1118
Class representing the general Input Record.
Definition: inputrecord.h:101
void exportPrimaryVars(VTKPiece &piece, IntArray &mapG2L, IntArray &mapL2G, int region, TimeStep *tStep)
Export primary variables.
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
void setConnectivity(int cellNum, IntArray &nodes)
Abstract base class representing the "problem" under consideration.
Definition: engngm.h:181
#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
Node * giveNode(int n)
Service for accessing particular domain node.
Definition: domain.h:371
the oofem namespace is to define a context or scope in which all oofem names are defined.
Domain * giveDomain(int n)
Service for accessing particular problem domain.
Definition: engngm.C:1720
virtual bool isCast(TimeStep *tStep)
Definition: element.C:809
int giveNumber() const
Definition: femcmpnn.h:107
Node * giveNode(int i) const
Returns reference to the i-th node of element.
Definition: element.h:610
Class representing solution step.
Definition: timestep.h:80
const IntArray & giveElementList()
Returns list of elements within set.
Definition: set.C:138
REGISTER_ExportModule(ErrorCheckingExportModule)

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