[Trilinos-Users] Question on use of Epetra_BlockMap with Epetra_Vector

Heroux, Michael A maherou at sandia.gov
Tue May 21 12:49:35 MDT 2013


Robert,

I think the code you sent is accurate, but I don't know that you will want to use a blockmap.  Blockmap only makes sense if you are going to use a VBR matrix with a compatible block partitioning.

If you are going to use a CRS matrix, you will want to create a standard Epetra map.  You will need unique GIDs for the map.  There is a utility function BlockMap2PointMap in Epetra_VbrMatrix that shows how to do this:

http://trilinos.sandia.gov/packages/docs/r11.2/packages/epetra/browser/doc/html/classEpetra__VbrMatrix.html#ae135a223d42bc85db28b23141861d11c

Let me know if you have any issues.

Mike

On May 21, 2013, at 1:18 PM, "Robert Crockett" <rkcrockett at txcorp.com<mailto:rkcrockett at txcorp.com>> wrote:

Hello,

The context:
I would like to be able to instantiate an Epetra_Vector from data in
existing memory, but only using a subset of that data. If anyone can
provide answers to the following questions, I would be most thankful.

The questions:
(1) Is this possible and, if so, is my below proposal a good way
    to do it?
(2) Is there existing example code that does this?

The details:
We have data in an existing (parallel) data structure, which holds it in
a contiguous swath of memory, eg.

int m_dataSize;  // size of data
double* m_data;  // pointer to data

I would like to pass this memory into Epetra_Vector for performing
solves elsewhere in Trilinos. The catch is that only parts of m_data are
valid degrees of freedom for the solves. These parts are contiguous, so
I believe using an Epetra_BlockMap fits the bill. Here is my proposed
pseudo-code for how I would do this, which is the basis for my questions
above:

<code>
Epetra_Comm* m_epCommPtr;  // initialized elsewhere
Epetra_BlockMap* m_mapPtr;
Epetra_Vector* m_vecPtr;
vector<int> myGlobalElements;
vector<int> elementSizeList;
// iterate only over valid DoF
DataIterator it(m_data, m_dataSize, VALID_CELLS);
/// add all valid cells to global element list for this proc,
//+ and add a one for element size (will optimize for contiguous
//+ blocks -- ie. element size > 1 -- later)
for (it.begin(); it.ok(); ++it) {
  myGlobalElements.push_back(it.globalIndex());
  elementSizeList.push_back(1);
}

// compute all element counts
int numMyElements = myGlobalElements.size();
int numGlobalElements;
m_comm.allReduceSum(&numMyElements, &numGlobalElements);
int indexBase = 0;
m_mapPtr = new Epetra_BlockMap(numGlobalElements, numMyElements,
                               myGlobalElements, elementSizeList,
                               indexBase, *m_epCommPtr);
m_vecPtr = new Epetra_Vector(View, *m_mapPtr, m_dataPtr);
</code>

That's it! Thanks,
Robert

--
"The mind is not a vessel to be filled but a fire to be kindled."
-Plutarch
------------------------------------------------------------------------------------------------
Robert Crockett, Ph.D.            Research Scientist, Tech-X Corporation
rkcrockett at txcorp.com<mailto:rkcrockett at txcorp.com>             5621 Arapahoe Avenue, Suite A
303/996-2037                      Boulder, CO 80303
------------------------------------------------------------------------------------------------

_______________________________________________
Trilinos-Users mailing list
Trilinos-Users at software.sandia.gov<mailto:Trilinos-Users at software.sandia.gov>
http://software.sandia.gov/mailman/listinfo/trilinos-users

-------------- next part --------------
An HTML attachment was scrubbed...
URL: https://software.sandia.gov/pipermail/trilinos-users/attachments/20130521/5b02ff36/attachment.html 


More information about the Trilinos-Users mailing list