/* * \file Elements.c * \brief: implementation of the Elements class, derived from DataSet class */ /*Headers: {{{1*/ #ifdef HAVE_CONFIG_H #include "config.h" #else #error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!" #endif #include #include #include #include #include "./DataSet.h" #include "../shared/shared.h" #include "../include/include.h" #include "../EnumDefinitions/EnumDefinitions.h" using namespace std; /*}}}*/ /*Object constructors and destructor*/ /*FUNCTION Elements::Elements(){{{1*/ Elements::Elements(){ return; } /*}}}*/ /*FUNCTION Elements::Elements(int in_enum){{{1*/ Elements::Elements(int in_enum): DataSet(in_enum){ //do nothing; return; } /*}}}*/ /*FUNCTION Elements::~Elements(){{{1*/ Elements::~Elements(){ return; } /*}}}*/ /*Object management*/ /*FUNCTION Elements::ProcessResultsUnits{{{1*/ void Elements::ProcessResultsUnits(void){ int i; //Process results to be output in the correct units for(i=0;iSize();i++){ Element* element=(Element*)this->GetObjectByOffset(i); element->ProcessResultsUnits(); } } /*}}}*/ /*FUNCTION Elements::DeleteResults{{{1*/ void Elements::DeleteResults(void){ int i; for (i=0;iSize();i++){ Element* element=(Element*)this->GetObjectByOffset(i); element->DeleteResults(); } } /*}}}*/ /*FUNCTION Elements::ToResults{{{1*/ void Elements::ToResults(Results* results,Parameters* parameters,int step, double time){ /*output: */ Patch* patch=NULL; /*I/O strategy: */ bool io_gather=true; //the default /*Recover parameters: */ parameters->FindParam(&io_gather,IoGatherEnum); /*create patch object out of all results in this dataset: */ patch=this->ResultsToPatch(); /*Gather onto master cpu 0, if needed: */ #ifdef _PARALLEL_ if(io_gather)patch->MPI_Gather(); #endif /*create result object and add to results dataset:*/ if (patch->maxvertices && patch->maxnodes){ results->AddObject(new IntExternalResult(results->Size()+1,PatchVerticesEnum,patch->maxvertices,step,time)); results->AddObject(new IntExternalResult(results->Size()+1,PatchNodesEnum, patch->maxnodes,step,time)); results->AddObject(new DoubleMatExternalResult(results->Size()+1,PatchEnum,patch->values,patch->numrows,patch->numcols,step,time)); } /*Free ressources:*/ delete patch; } /*}}}*/ /*FUNCTION Elements::ResultsToPatch{{{1*/ Patch* Elements::ResultsToPatch(void){ /*output: */ Patch* patch=NULL; /*intermediary: */ int i; int count; int numrows; int numvertices; int numnodes; int max_numvertices; int max_numnodes; int element_numvertices; int element_numrows; int element_numnodes; /*We are going to extract from the results within the elements, the desired results , and create a table * of patch information, that will hold, for each element that computed the result that * 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 * at the nodes (could be different from the vertices). This will be used for visualization purposes. * For example, we could build the following patch table, for velocities: * 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 VxEnum 1 .5 1 P0 1 2 4.5 NaN NaN 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 VzEnum 2 .8 2 P1 1 3 4 4.5 3.2 2.5 * ... etc ... * * So what do we need to build the table: the maximum number of vertices included in the table, * and the maximum number of nodal values, as well as the number of rows. Once we have that, * we ask the elements to fill their own row in the table, by looping on the elememnts. * * We will use the Patch object, which will store all of the information needed, and will be able * to output itself to disk on its own. See object/Patch.h for format of this object.*/ /*First, determine maximum number of vertices, nodes, and number of results: */ numrows=0; numvertices=0; numnodes=0; for(i=0;iSize();i++){ Element* element=(Element*)this->GetObjectByOffset(i); element->PatchSize(&element_numrows,&element_numvertices,&element_numnodes); numrows+=element_numrows; if(element_numvertices>numvertices)numvertices=element_numvertices; if(element_numnodes>numnodes)numnodes=element_numnodes; } #ifdef _PARALLEL_ /*Synchronize across cluster, so as to not end up with different sizes for each patch on each cpu: */ MPI_Reduce (&numvertices,&max_numvertices,1,MPI_INT,MPI_MAX,0,MPI_COMM_WORLD ); MPI_Bcast(&max_numvertices,1,MPI_INT,0,MPI_COMM_WORLD); numvertices=max_numvertices; MPI_Reduce (&numnodes,&max_numnodes,1,MPI_INT,MPI_MAX,0,MPI_COMM_WORLD ); MPI_Bcast(&max_numnodes,1,MPI_INT,0,MPI_COMM_WORLD); numnodes=max_numnodes; #endif /*Ok, initialize Patch object: */ patch=new Patch(numrows,numvertices,numnodes); /*Now, go through elements, and fill the Patch object: */ count=0; for(i=0;iSize();i++){ Element* element=(Element*)this->GetObjectByOffset(i); element->PatchFill(&count,patch); } return patch; } /*}}}*/