[Trilinos-Users] Epetra - Dense Matrix-Vector multiply basics
Gyorgy Matyasfalvi
matyasfalvi at gmail.com
Sat Apr 11 19:44:43 MDT 2015
Dear Erik and Michael:
thanks for your responses.
I've tried out some test cases and I have some further questions regarding
reading in data and rebalancing with Isorropia.
1)
For EpetraExt::MatrixMarketFileToCrsMatrix(datafile, Comm, A, false, true)
the more processors I use the slower the "File Read time (secs)" gets. E.g.
for a 193MB file I get 3.47222 sec for 1 proc and 5.84088 sec for 2 proc
etc. I think this makes sense because only processor 0 is reading the data
file but I just wanted to double check with you, if this is really makes
sense.
2)
I've created a sparse matrix A using:
EpetraExt::MatrixMarketFileToCrsMatrix(datafile, comm, A, false, true);
This is matrix has 19996 rows, 1355191 columns and 9097916 nonzeros.
And then used Isorropia to rebalance it with:
Teuchos::ParameterList paramlist;
paramlist.set("BALANCE OBJECTIVE", "NONZEROS");
balA = Isorropia::Epetra::createBalancedCopy(*A, paramlist);
I got the following result for 4 cores:
Proc 1 nonzeros = 2501926
Proc 2 nonzeros = 2365455
Proc 3 nonzeros = 2501919
Proc 4 nonzeros = 1728616
Which means that proc 4 has about 30% less nonzeros than the others.
So I tried to improve this by forcing a split along the columns when
creating the matrix by using:
EpetraExt::MatrixMarketFileToCrsMatrix(datafile, rowLocalMap, colMap, A,
false, true);
where rowLocalMap and colMap are the following:
Epetra_LocalMap rowLocalMap(numRows, 0, comm);
Epetra_Map colMap(numCols, 0, comm);
When I tried to use Isorropia on this I got the following error messages:
**************************************************************************************
Epetra ERROR -2,
/home/gyorgy/Trilinos/Trilinos-11.4.3-Source/packages/epetra/src/Epetra_CrsMatrix.cpp,
line 751
Epetra ERROR -2,
/home/gyorgy/Trilinos/Trilinos-11.4.3-Source/packages/epetra/src/Epetra_CrsMatrix.cpp,
line 574
Epetra ERROR -2,
/home/gyorgy/Trilinos/Trilinos-11.4.3-Source/packages/epetra/src/Epetra_CrsMatrix.cpp,
line 2333
Epetra ERROR -2,
/home/gyorgy/Trilinos/Trilinos-11.4.3-Source/packages/epetra/src/Epetra_CrsMatrix.cpp,
line 2175
Epetra ERROR -2,
/home/gyorgy/Trilinos/Trilinos-11.4.3-Source/packages/epetra/src/Epetra_DistObject.cpp,
line 242
Epetra ERROR -2,
/home/gyorgy/Trilinos/Trilinos-11.4.3-Source/packages/epetra/src/Epetra_DistObject.cpp,
line 137
**************************************************************************************
I've looked a the Errors they have to do with inserting global values,
somehow the indices are messed up. But I don't know why. I think the way I
have created matrix A should be legitimate, or?
Thanks a lot! Best,
Gyorgy
On Mon, Apr 6, 2015 at 11:30 AM, Erik Boman <egboman at sandia.gov> wrote:
> Gyorgy,
>
> Epetra does not do load balancing, except create maps in the trivial way.
> For sparse matrices, you should use Zoltan or Isorropia. For your case.
> I'd recommend the PHG hypergraph partitioner in Zoltan (also available in
> Isorropia with Epetra interfaces).
> You can balance with respect to nonzeros instead of rows, which should
> take care of the dense row/column problem in most cases. (Parameter BALANCE
> OBJECTIVE = nonzeros in Isorropia):
>
> http://trilinos.org/docs/r11.14/packages/isorropia/doc/html/index.html
>
> It is possible to split rows in Epetra (overlapping row maps) but often
> this is not worth the extra complexity.
>
> Hope this helps,
> Erik
>
> Gyorgy Matyasfalvi wrote:
>
>> Dear Michael:
>>
>> thanks for your response. I thought that maybe reading in overlapping
>> chuncks of the text file using MPI I/O might work. But I never implemented
>> it, so I don't know how it would actually perform. I will look into HDF5
>> more thoroughly.
>> I have one more question in terms of the data distribution. Say a sparse
>> matrix has some very dense columns or rows. In that case a processor could
>> be assigned substantially more work then everyone else creating a
>> bottleneck for the algorithm. How does Epetra deal with this situation?
>> I've encountered the following issues:
>> i) Assigning columns and rows randomly may not work too well, say if
>> there is only one completely dense column than the processor that gets it
>> has more work then anyone else. ii) Coming up with an optimal load balance
>> is NP hard. iii) Splitting data according to number of non-zeros will
>> require columns or rows to be split between processors, which would violate
>> the unique global ID (GID) requirement of Epetra.
>>
>> Thanks a lot!
>> Best,
>> Gyorgy
>>
>> On Mon, Apr 6, 2015 at 12:08 AM, Heroux, Mike <MHeroux at csbsju.edu
>> <mailto:MHeroux at csbsju.edu>> wrote:
>>
>> Gyorgy,
>>
>> The Matrix Market I/O functions are parallel, but only process 0
>> reads from the file. The file contents are used to construction
>> an Epetra object on process 0, which is then exported to a
>> distributed object that is split across all processes according to
>> an Epetra map. The default map splits the data by a uniform
>> distribution. You can prescribe your own distribution too.
>>
>> This is the only portable approach to I/O that I know of for
>> reading and writing a plain text file.
>> In general, text-based file I/O from a single file is not very
>> scalable. It’s primarily meant as a way to enable problem
>> studies. You would probably be better off using the HDF5
>> functionality from EpetraExt.
>>
>> Mike
>>
>> From: Gyorgy Matyasfalvi <matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com>>
>> Date: Sunday, April 5, 2015 at 10:52 AM
>> To: "maherou at sandia.gov <mailto:maherou at sandia.gov>"
>> <maherou at sandia.gov <mailto:maherou at sandia.gov>>
>> Cc: "trilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov>"
>> <trilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov>>
>> Subject: Re: [Trilinos-Users] Epetra - Dense Matrix-Vector
>> multiply basics
>>
>> Dear Michael:
>>
>> thanks for your response. Ok that makes sense. One more question
>> I have in terms of the MatrixMarket I/O routines. To me it seems that
>> the MatrixMarketFileToCrsMatrix in EpetraExt
>> is actually a serial operation since only processor 0 does the
>> data read-in and then broadcasts the data to the other processors.
>> Is that correct?
>>
>> Thanks! Best,
>> Gyorgy
>>
>> On Sat, Apr 4, 2015 at 8:23 PM, Heroux, Michael A
>> <maherou at sandia.gov <mailto:maherou at sandia.gov>> wrote:
>>
>> Gyorgy,
>>
>> I believe the operations you need are supported. In the case
>> of a wide matrix, for either a sparse or dense matrix, you
>> want to think about storing A^T, and then apply the transpose
>> operation to get the action of A. This is basically what I
>> described in the previous message.
>>
>> You are also right about the square case, that doing a 2D
>> decomposition would be appropriate, and Epetra (and similarly
>> Tpetra) is unique in its support for this kind of
>> decomposition, compared to other popular libraries. Even so,
>> the infrastructure to support construction and management of
>> 2D decompositions can be substantially improved.
>>
>> If you need more details on what I describe in the first
>> paragraph, let me know.
>>
>> Mike
>>
>> On Apr 3, 2015, at 1:33 PM, Gyorgy Matyasfalvi
>> <matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com>>> wrote:
>>
>> Dear Michael:
>>
>> We're solving the following optimization problem:
>>
>> min ||A(x^+ - x^-) - b||^2_2 + \nu * 1^T (x^+ + x^-)
>> x^+, x^- \in R^n
>>
>> s.t. x^+, x^- >= 0
>>
>> aka (LASSO) using a first order method called the Nonmonotone
>> Spectral Projected Gradient (SPG).
>> ( You can find more details in this paper:
>> http://www.optimization-online.org/DB_FILE/2015/01/4748.pdf )
>>
>> The matrix A in the objective function comes in many different
>> forms: tall, wide, square. In one of our test problems A is a
>> wide matrix, with dimensions: 19,996 rows and 1.355,191 columns.
>> Our most "expensive" operations are mtx-vec multiplies when
>> computing the objective function value or the gradient.
>>
>> What I'd like to achieve is that our code can handle matrices
>> of arbitrary form "appropriately". For example if A is wide
>> i.e. few rows but many columns then we partition A along it's
>> columns (domain decomposition), if A is tall i.e. many rows
>> and few columns then we partition A along it's rows (range
>> decomposition). If A is square and huge than maybe do both
>> range and domain decomposition.
>> It seems to me this is actually doable for sparse matrices
>> using the Epetra_CrsMatrix class. My question is can I do this
>> in the dense case with the Epetra_MultiVector class?
>>
>> Thanks a lot for your help!
>> Best,
>> Gyorgy
>>
>> On Fri, Apr 3, 2015 at 1:24 AM, Heroux, Michael A
>> <maherou at sandia.gov
>> <mailto:maherou at sandia.gov><mailto:maherou at sandia.gov
>> <mailto:maherou at sandia.gov>>> wrote:
>> Gyorgy,
>>
>> Please give me some details about your problem: Algorithm,
>> size of matrices, type of data distribution. Just the basics
>> to start.
>>
>> Thanks.
>>
>> Mike
>>
>> On Apr 2, 2015, at 11:03 PM, Gyorgy Matyasfalvi
>> <matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com>><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com>>>> wrote:
>>
>> Dear Mike:
>>
>> Thanks for your help. Indeed I got an error code of -2 and I
>> suspect that, as you pointed out, the issue is that *this is
>> not a row vector.
>>
>> My other question then would be: Is it possible to define a
>> "Domain Map" for Multivectors? i.e I want processors to own
>> columns of the matrix instead of rows. So basically spread out
>> the columns among different processors. I know think this is
>> possible for sparse matrices.
>>
>> Thanks! Best,
>> Gyorgy
>>
>>
>>
>> On Thu, Apr 2, 2015 at 9:02 AM, Heroux, Michael A
>> <maherou at sandia.gov
>> <mailto:maherou at sandia.gov><mailto:maherou at sandia.gov
>> <mailto:maherou at sandia.gov>><mailto:maherou at sandia.gov
>> <mailto:maherou at sandia.gov><mailto:maherou at sandia.gov
>> <mailto:maherou at sandia.gov>>>> wrote:
>> Gyorgy,
>>
>> Assuming your x, y and A are dense objects, you are using the
>> correct function, Multiply(). Chances are that you have some
>> detail incorrect. First you should check for error codes that
>> are being returned. Multiply() will return integer error
>> codes if it runs into issues:
>>
>> Int ierr = y.Multiply(...)
>>
>> If (ierr!=0) ...
>>
>> Regarding your transpose operation, if I understand the
>> situation correctly, x should also be transposed to be a row
>> vector. Multiply() is probably complaining about this. There
>> is no way in the Epetra framework to specify that the "this"
>> argument is transposed.
>>
>> I think it is best to perform the transpose operation as
>>
>> x = A^T y
>>
>> This expression keeps the dimensions correct.
>>
>> Try this and see how it goes.
>>
>> Mike
>>
>> From: Gyorgy Matyasfalvi <matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com>><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com>>><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com>><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com><mailto:matyasfalvi at gmail.com
>> <mailto:matyasfalvi at gmail.com>>>>>
>> Date: Wednesday, April 1, 2015 at 10:41 PM
>> To: "trilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov><mailto:tr
>> ilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov>><mailto:t
>> rilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov><mailto:tr
>> ilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov>>><mailto:
>> trilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov><mailto:tr
>> ilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov>><mailto:t
>> rilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov><mailto:tr
>> ilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov>>>>"
>> <trilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov><mailto:tr
>> ilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov>><mailto:t
>> rilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov><mailto:tr
>> ilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov>>><mailto:
>> trilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov><mailto:tr
>> ilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov>><mailto:t
>> rilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov><mailto:tr
>> ilinos-users at software.sandia.gov
>> <mailto:trilinos-users at software.sandia.gov>>>>>
>> Subject: [EXTERNAL] [Trilinos-Users] Epetra - Dense
>> Matrix-Vector multiply basics
>>
>> Dear Developers/Users:
>>
>> I have a question concerning matrix-vector multiplies:
>>
>> 1) Is there a document or example file that explains the
>> basics of matrix-vector multiplies? Preferably isolated
>> examples of mtx-vec mults. instead of being embedded in some
>> algorithm.
>>
>> If this doesn't exist. Then I'd be interested in finding out:
>>
>> 2) How to do dense matrix, dense vector multiplies?
>> I assume it's done via the Epetra_MultiVector class. However,
>> the only operation that seems to have an effect is the following:
>>
>> y = A x
>> using
>> y.Multiply('N', 'N', 1.0, A, x, 0.0)
>> (x was constructed using an Epetra_LocalMap
>> A, and y were constructed using the same Epetra_Map)
>>
>> The problem is that if I want to compute for example:
>>
>> x = y^T A
>> using
>> x.Multiply('T', 'N', 1.0, y, A, 0.0)
>>
>> it seems the Multiply() operation doesn't do anything x
>> remains unchanged.
>>
>> What am I doing incorrectly here?
>> Thanks a lot! Best,
>> Gyorgy
>>
>>
>>
>>
>>
>>
>>
>>
>> ------------------------------------------------------------------------
>>
>> _______________________________________________
>> Trilinos-Users mailing list
>> Trilinos-Users at software.sandia.gov
>> https://software.sandia.gov/mailman/listinfo/trilinos-users
>>
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://software.sandia.gov/pipermail/trilinos-users/attachments/20150411/3cd00253/attachment.html>
More information about the Trilinos-Users
mailing list