source: issm/trunk/src/c/Container/Elements.cpp@ 9874

Last change on this file since 9874 was 9874, checked in by Mathieu Morlighem, 13 years ago

Added support for result_on_vertices
File size: 7.9 KB
Line 
1/*
2 * \file Elements.c
3 * \brief: implementation of the Elements class, derived from DataSet class
4 */
5
6/*Headers: {{{1*/
7#ifdef HAVE_CONFIG_H
8 #include <config.h>
9#else
10#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
11#endif
12
13#include <vector>
14#include <functional>
15#include <algorithm>
16#include <iostream>
17
18#include "./DataSet.h"
19#include "../shared/shared.h"
20#include "../include/include.h"
21#include "../EnumDefinitions/EnumDefinitions.h"
22
23using namespace std;
24/*}}}*/
25
26/*Object constructors and destructor*/
27/*FUNCTION Elements::Elements(){{{1*/
28Elements::Elements(){
29 return;
30}
31/*}}}*/
32/*FUNCTION Elements::Elements(int in_enum){{{1*/
33Elements::Elements(int in_enum): DataSet(in_enum){
34 //do nothing;
35 return;
36}
37/*}}}*/
38/*FUNCTION Elements::~Elements(){{{1*/
39Elements::~Elements(){
40 return;
41}
42/*}}}*/
43
44/*Object management*/
45/*FUNCTION Elements::Configure{{{1*/
46void Elements::Configure(Elements* elements,Loads* loads, Nodes* nodes, Vertices* vertices, Materials* materials,Parameters* parameters){
47
48 vector<Object*>::iterator object;
49 Element* element=NULL;
50
51 for ( object=objects.begin() ; object < objects.end(); object++ ){
52
53 element=(Element*)(*object);
54 element->Configure(elements,loads,nodes,materials,parameters);
55
56 }
57
58}
59/*}}}*/
60/*FUNCTION Elements::ProcessResultsUnits{{{1*/
61void Elements::ProcessResultsUnits(void){
62
63 //Process results to be output in the correct units
64 for(int i=0;i<this->Size();i++){
65 Element* element=(Element*)this->GetObjectByOffset(i);
66 element->ProcessResultsUnits();
67 }
68}
69/*}}}*/
70/*FUNCTION Elements::DeleteResults{{{1*/
71void Elements::DeleteResults(void){
72
73 for (int i=0;i<this->Size();i++){
74 Element* element=(Element*)this->GetObjectByOffset(i);
75 element->DeleteResults();
76 }
77}
78/*}}}*/
79/*FUNCTION Elements::ResultsToPatch{{{1*/
80Patch* Elements::ResultsToPatch(void){
81
82 /*output: */
83 Patch* patch=NULL;
84
85 /*intermediary: */
86 int i;
87 int count;
88 int numrows;
89 int numvertices;
90 int numnodes;
91 int max_numvertices;
92 int max_numnodes;
93 int element_numvertices;
94 int element_numrows;
95 int element_numnodes;
96
97 /*We are going to extract from the results within the elements, the desired results , and create a table
98 * of patch information, that will hold, for each element that computed the result that
99 * we desire, the enum_type of the result, the step and time, the id of the element, the interpolation type, the vertices ids, and the values
100 * at the nodes (could be different from the vertices). This will be used for visualization purposes.
101 * For example, we could build the following patch table, for velocities:
102 *
103 1. on a Beam element, Vx, at step 1, time .5, element id 1, interpolation type P0 (constant), vertices ids 1 and 2, one constant value 4.5
104 VxEnum 1 .5 1 P0 1 2 4.5 NaN NaN
105 2. on a Tria element, Vz, at step 2, time .8, element id 2, interpolation type P1 (linear), vertices ids 1 3 and 4, with values at 3 nodes 4.5, 3.2, 2.5
106 VzEnum 2 .8 2 P1 1 3 4 4.5 3.2 2.5
107 * ... etc ...
108 *
109 * So what do we need to build the table: the maximum number of vertices included in the table,
110 * and the maximum number of nodal values, as well as the number of rows. Once we have that,
111 * we ask the elements to fill their own row in the table, by looping on the elememnts.
112 *
113 * We will use the Patch object, which will store all of the information needed, and will be able
114 * to output itself to disk on its own. See object/Patch.h for format of this object.*/
115
116 /*First, determine maximum number of vertices, nodes, and number of results: */
117 numrows=0;
118 numvertices=0;
119 numnodes=0;
120
121 for(i=0;i<this->Size();i++){
122
123 Element* element=(Element*)this->GetObjectByOffset(i);
124 element->PatchSize(&element_numrows,&element_numvertices,&element_numnodes);
125
126 numrows+=element_numrows;
127 if(element_numvertices>numvertices)numvertices=element_numvertices;
128 if(element_numnodes>numnodes)numnodes=element_numnodes;
129 }
130
131 #ifdef _PARALLEL_
132 /*Synchronize across cluster, so as to not end up with different sizes for each patch on each cpu: */
133 MPI_Reduce (&numvertices,&max_numvertices,1,MPI_INT,MPI_MAX,0,MPI_COMM_WORLD );
134 MPI_Bcast(&max_numvertices,1,MPI_INT,0,MPI_COMM_WORLD);
135 numvertices=max_numvertices;
136
137 MPI_Reduce (&numnodes,&max_numnodes,1,MPI_INT,MPI_MAX,0,MPI_COMM_WORLD );
138 MPI_Bcast(&max_numnodes,1,MPI_INT,0,MPI_COMM_WORLD);
139 numnodes=max_numnodes;
140 #endif
141
142 /*Ok, initialize Patch object: */
143 patch=new Patch(numrows,numvertices,numnodes);
144
145 /*Now, go through elements, and fill the Patch object: */
146 count=0;
147 for(i=0;i<this->Size();i++){
148 Element* element=(Element*)this->GetObjectByOffset(i);
149 element->PatchFill(&count,patch);
150 }
151
152 return patch;
153}
154/*}}}*/
155/*FUNCTION Elements::SetCurrentConfiguration{{{1*/
156void Elements::SetCurrentConfiguration(Elements* elements,Loads* loads, Nodes* nodes, Vertices* vertices, Materials* materials,Parameters* parameters){
157
158 vector<Object*>::iterator object;
159 Element* element=NULL;
160
161 for ( object=objects.begin() ; object < objects.end(); object++ ){
162
163 element=(Element*)(*object);
164 element->SetCurrentConfiguration(elements,loads,nodes,materials,parameters);
165
166 }
167
168}
169/*}}}*/
170/*FUNCTION Elements::ToResults{{{1*/
171void Elements::ToResults(Results* results,Parameters* parameters,int step, double time){
172
173 Patch *patch = NULL;
174 int *resultsenums = NULL;
175 double *vector_serial= NULL;
176 Vec vector = NULL;
177 bool io_gather;
178 bool results_on_vertices;
179 int numberofvertices;
180 int numberofresults;
181
182 /*Recover parameters: */
183 parameters->FindParam(&io_gather,SettingsIoGatherEnum);
184 parameters->FindParam(&results_on_vertices,SettingsResultsOnVerticesEnum);
185 parameters->FindParam(&numberofvertices,MeshNumberofverticesEnum);
186
187 if(results_on_vertices){
188 /*No patch here, we prepare vectors*/
189
190 /*OK, see what the first element of this partition has in stock (this is common to all partitions)*/
191 Element* element=(Element*)this->GetObjectByOffset(0);
192 element->ListResultsEnums(&resultsenums,&numberofresults);
193
194 /*Loop over all results and get nodal vector*/
195 for(int i=0;i<numberofresults;i++){
196
197 /*Get vector for result number i*/
198 vector=NewVec(numberofvertices);
199 for(int j=0;j<this->Size();j++){
200 Element* element=(Element*)this->GetObjectByOffset(j);
201 element->GetVectorFromResults(vector,resultsenums[i]);
202 }
203 VecAssemblyBegin(vector);
204 VecAssemblyEnd(vector);
205
206 /*Serialize and add to results*/
207 VecToMPISerial(&vector_serial,vector);
208 results->AddObject(new DoubleVecExternalResult(results->Size()+1,resultsenums[i],vector_serial,numberofvertices,step,time));
209
210 /*clean up*/
211 VecFree(&vector);
212 xfree((void**)&vector_serial);
213 }
214
215 }
216 else{
217 /*create patch object out of all results in this dataset: */
218 patch=this->ResultsToPatch();
219
220 /*Gather onto master cpu 0, if needed: */
221#ifdef _PARALLEL_
222 if(io_gather)patch->Gather();
223#endif
224
225 /*create result object and add to results dataset:*/
226 if (patch->maxvertices && patch->maxnodes){
227 results->AddObject(new IntExternalResult(results->Size()+1,PatchVerticesEnum,patch->maxvertices,step,time));
228 results->AddObject(new IntExternalResult(results->Size()+1,PatchNodesEnum, patch->maxnodes,step,time));
229 results->AddObject(new DoubleMatExternalResult(results->Size()+1,PatchEnum,patch->values,patch->numrows,patch->numcols,step,time));
230 }
231 }
232
233 /*Free ressources:*/
234 xfree((void**)&resultsenums);
235 delete patch;
236}
237/*}}}*/
238/*FUNCTION Elements::NumberOfElements{{{1*/
239int Elements::NumberOfElements(void){
240
241 int local_nelem=0;
242 int numberofelements;
243
244 #ifdef _PARALLEL_
245 local_nelem=this->Size();
246 MPI_Allreduce ( (void*)&local_nelem,(void*)&numberofelements,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
247 #else
248 numberofelements=this->Size();
249 #endif
250
251 return numberofelements;
252}
253/*}}}*/
254/*FUNCTION Elements::InputCopy{{{1*/
255void Elements::InputDuplicate(int input_enum,int output_enum){
256
257 for(int i=0;i<this->Size();i++){
258 Element* element=(Element*)this->GetObjectByOffset(i);
259 element->InputDuplicate(input_enum,output_enum);
260 }
261}
262/*}}}*/
Note: See TracBrowser for help on using the repository browser.