void Exporter<MeshType>::exportRegionMarkerID( boost::shared_ptr<MeshType> mesh, boost::shared_ptr<Epetra_Comm> comm ) { // TODO: use FESpace M_spacemap for generality const ReferenceFE* refFEPtr; // Need a factory!!!! // @todo Need a factory! switch ( MeshType::S_geoDimensions ) { case 3: refFEPtr = &feTetraP0; break; case 2: refFEPtr = &feTriaP0; break; case 1: refFEPtr = &feSegP0; break; default: ERROR_MSG ( "Dimension not supported " ); } // Useless quadrature rule const QuadratureRule & dummyQR = quadRuleDummy; const feSpacePtr_Type regionMarkerID_FESpacePtr( new feSpace_Type( mesh, *refFEPtr, dummyQR, dummyQR, 1, comm ) ); vectorPtr_Type regionMarkerIDData ( new vector_Type ( regionMarkerID_FESpacePtr->map() ) ); for ( UInt iElem( 0 ); iElem < mesh->numElements(); ++iElem ) { const ID globalElem = mesh->element(iElem).id(); (*regionMarkerIDData)[ globalElem ] = mesh->element(iElem).markerID(); } addVariable( exporterData_Type::ScalarField, "regionMarkerID", regionMarkerID_FESpacePtr, regionMarkerIDData, 0, exporterData_Type::SteadyRegime, exporterData_Type::Cell ); } // exportRegionMarkerID
UInt ExporterVTK<MeshType>::whichCellType ( const feSpacePtr_Type& _feSpacePtr ) { ASSERT ( _feSpacePtr.get(), "\nA pointer to a valid FE object is required!"); UInt vtkCellType (0); switch ( _feSpacePtr->fe().refFE().type() ) { case FE_P1_2D: // case FE_P1bubble_2D: vtkCellType = VTK_TRIANGLE; break; case FE_P2_2D: vtkCellType = VTK_QUADRATIC_TRIANGLE; break; case FE_P0_3D: case FE_P1_3D: vtkCellType = VTK_TETRA; break; case FE_P2_3D: case FE_P2tilde_3D: vtkCellType = VTK_QUADRATIC_TETRA; break; case FE_Q1_2D: vtkCellType = VTK_QUAD; break; case FE_Q2_2D: vtkCellType = VTK_QUADRATIC_QUAD; break; case FE_Q0_3D: case FE_Q1_3D: vtkCellType = VTK_HEXAHEDRON; break; case FE_Q2_3D: vtkCellType = VTK_QUADRATIC_HEXAHEDRON; break; default: ERROR_MSG ( "WARNING: the element is not yet implemented in ExporterVTK\n" ) break; } return vtkCellType; }
void Exporter<MeshType>::addVariable(const FieldTypeEnum& type, const std::string& variableName, const feSpacePtr_Type& feSpacePtr, const vectorPtr_Type& vectorPtr, const UInt& start, const FieldRegimeEnum& regime, const WhereEnum& where ) { M_dataVector.push_back( exporterData_Type(type, variableName, feSpacePtr, vectorPtr, start, regime, where) ); M_whereToDataIdMap.insert( std::pair<WhereEnum,UInt >(where, M_dataVector.size()-1 ) ); M_feTypeToDataIdMap.insert( std::pair<FE_TYPE,UInt >(feSpacePtr->fe().refFE().type(), M_dataVector.size()-1 ) ); }
ExporterData<MeshType>::ExporterData( const FieldTypeEnum& type, const std::string& variableName, const feSpacePtr_Type& feSpacePtr, const vectorPtr_Type& vectorPtr, const UInt& start, const FieldRegimeEnum& regime, const WhereEnum& where ): M_variableName ( variableName ), M_feSpacePtr ( feSpacePtr ), M_storedArrayPtr ( vectorPtr ), M_numDOF ( feSpacePtr->dim() ), M_start ( start ), M_fieldType ( type ), M_regime ( regime ), M_where ( where ) {}
void computeP0pressure (const feSpacePtr_Type& pFESpacePtr, const feSpacePtr_Type& p0FESpacePtr, const feSpacePtr_Type& uFESpacePtr, const vector_Type& velAndPressure, vector_Type& P0pres, Real /*time*/) { int MyPID; MPI_Comm_rank (MPI_COMM_WORLD, &MyPID); UInt offset = 3 * uFESpacePtr->dof().numTotalDof(); std::vector<int> gid0Vec (0); gid0Vec.reserve (p0FESpacePtr->mesh()->numVolumes() ); std::vector<Real> val0Vec (0); val0Vec.reserve (p0FESpacePtr->mesh()->numVolumes() ); for (UInt ivol = 0; ivol < pFESpacePtr->mesh()->numVolumes(); ++ivol) { pFESpacePtr->fe().update ( pFESpacePtr->mesh()->volumeList ( ivol ), UPDATE_DPHI ); p0FESpacePtr->fe().update ( p0FESpacePtr->mesh()->volumeList ( ivol ) ); UInt eleID = pFESpacePtr->fe().currentLocalId(); double tmpsum = 0.; for (UInt iNode = 0; iNode < (UInt) pFESpacePtr->fe().nbFEDof(); iNode++) { int ig = pFESpacePtr->dof().localToGlobalMap ( eleID, iNode ); tmpsum += velAndPressure (ig + offset); gid0Vec.push_back ( p0FESpacePtr->fe().currentId() ); val0Vec.push_back ( tmpsum / (double) pFESpacePtr->fe().nbFEDof() ); } } P0pres.setCoefficients (gid0Vec, val0Vec); P0pres.globalAssemble(); MPI_Barrier (MPI_COMM_WORLD); }