New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

add support to sparse linear algebra. #85

Merged

cnpetra merged 19 commits into develop from nywrk

Nov 15, 2020

Collaborator

nychiang commented Sep 23, 2020 •

edited by cnpetra

Adds an NLP general sparse problem formulation and interface, new sparse linear algebra, and backend linear solvers (MA57 and StrumPack). The PR focuses on software and some math aspects and does not fully support device computations. These will be addressed by a future PR.

pelesh reviewed

View reviewed changes

src/LinAlg/hiopMatrixSparseTripletStorage.hpp Outdated

                   Tidx *irow, *jcol;
-                  Tval *values;
+                  Tval *values;

Collaborator

pelesh Sep 23, 2020

Suggested change

      
                Tval *values;
          
                Tval* values_;

Perhaps it would be good to follow HiOp style guidelines.

pelesh reviewed

View reviewed changes

src/LinAlg/hiopMatrixSparseTriplet.cpp Outdated

+               * @pre 'this' must have exactly, or more than 'n_rows' rows
+               * @pre 'this' must have exactly, or more cols than 'src'
+               */
+              void hiopMatrixSparseTriplet::copyRowsFromSrcToDest(const hiopMatrix& src_gen,

Collaborator

pelesh Sep 23, 2020

Please add a unit test for this kernel.

pelesh reviewed

View reviewed changes

src/LinAlg/hiopMatrixSparseTriplet.cpp Outdated

                 assert(false && "not needed");
               }
+              void hiopMatrixSparseTriplet::copySubDiagonalEleFromVec(const long long& start_on_dest_diag, const long long& num_elems,

Collaborator

pelesh Sep 23, 2020

I think this kernel needs a unit test.

Collaborator Author

nychiang Oct 2, 2020

add a unit test into test/testMatrixSparse.cpp.
The corresponding code is added into test/LinAlg/matrixTestsSparseTriplet.cpp
However, the corresponding function verifyAnser(...) hasn't been implemented. Will set a follow-up meeting with Asher @ashermancinelli to discuss this issue.

Contributor

ashermancinelli Oct 5, 2020

Update: @nychiang and I discussed our unit test conventions, so he should now be able to update accordingly.

Collaborator

pelesh commented Sep 23, 2020

It would be good to take a look at RAJA snapshot branch to get a heads up on possible merge conflicts. RAJA PR will be submitted soon, so it would be good to coordinate.

Collaborator

cnpetra commented Sep 23, 2020

@pelesh RAJA PR will be merged first in master

cnpetra reviewed

View reviewed changes

src/Drivers/nlpMDS_ex5.hpp Outdated

               		 const double* x, bool new_x, double* cons)
                 {
                   const double* s = x+ns_;
                   const double* y = x+2*ns_;

Collaborator

cnpetra Sep 23, 2020

@nychiang : end of line is changed in very many places. This will make merging extremely difficult and laborious. I know why it happened ;) I suggest you somehow fix this.

Collaborator Author

nychiang Sep 24, 2020

I tried to run command dos2unit, but it seems that I need to fix this issue manually. :-)

Collaborator

cnpetra Sep 24, 2020

ouch... maybe core.autocrlf settings of git will help?

cnpetra reviewed

View reviewed changes

src/Drivers/nlpMDS_ex5.hpp Outdated

                 //sum 0.5 {x_i*(x_{i}-1) : i=1,...,ns} + 0.5 y'*Qd*y + 0.5 s^T s
                 bool eval_grad_f(const long long& n, const double* x, bool new_x, double* gradf)
                 {
                   //! assert(ns>=4); assert(Q->n()==ns/4); assert(Q->m()==ns/4);
-                  //x_i - 0.5

Collaborator

cnpetra Sep 23, 2020

@nychiang : also here, I guess \n changed with \r\n

Collaborator Author

nychiang commented Sep 26, 2020

Just push another commit, in which I fix the bug and white-space issue.
Ex6 now is working fine. Its usage is:
nlpSparse_ex6.exe n
where n>=3 is a the number of variables provided by the user

Collaborator Author

nychiang commented Sep 28, 2020

NOTE: need to fix the hard-coded path for Metis/Ma57 in CMake file before merging into the master branch

nychiang commented

View reviewed changes

tests/LinAlg/matrixTestsSparse.hpp

                 virtual const local_ordinal_type* getRowIndices(const hiop::hiopMatrixSparse* a) = 0;
                 virtual const local_ordinal_type* getColumnIndices(const hiop::hiopMatrixSparse* a) = 0;
                 virtual local_ordinal_type getLocalSize(const hiop::hiopVector* x) = 0;
                 virtual int verifyAnswer(hiop::hiopMatrix* A, real_type answer) = 0;
+                virtual int verifyAnswer(hiop::hiopMatrix* A, local_ordinal_type nnz_st, local_ordinal_type nnz_ed, const double answer) = 0;

Collaborator Author

nychiang Oct 7, 2020

@ashermancinelli Please verify. I add one function to verify the answers for a sparse matrix, since the other "verifyAnswer" requires a dense input matrix "A". See file matrixTestsSparseTriplet.cpp for examples. NOTE that since this branch is developed based on 'master', it doesn't has the RAJA implementation.

Contributor

ashermancinelli Oct 7, 2020

I think if you're just checking a single element, probably better to use getLocalElement and check result

ashermancinelli reviewed

View reviewed changes

tests/LinAlg/matrixTestsSparse.hpp

+                    }
+                  }
+                  assert(A.n() >= B.n());

Contributor

ashermancinelli Oct 7, 2020

Should these be assertions or failure conditions for the test?

tests/LinAlg/matrixTestsSparse.hpp Outdated Show resolved Hide resolved

tests/LinAlg/matrixTestsSparse.hpp

@@ @@ -432,10 +500,12 @@ class MatrixTestsSparse : public TestBase @@
                 virtual real_type getLocalElement(const hiop::hiopMatrix* a, local_ordinal_type i, local_ordinal_type j) = 0;
                 virtual real_type getLocalElement(const hiop::hiopVector* x, local_ordinal_type i) = 0;
                 virtual real_type* getMatrixData(hiop::hiopMatrixSparse* a) = 0;
+                virtual real_type getMatrixData(hiop::hiopMatrixSparse* a, local_ordinal_type i, local_ordinal_type j) = 0;

Contributor

ashermancinelli Oct 7, 2020

Why create another function instead of using getLocalElement?

Collaborator Author

nychiang Oct 7, 2020

getLocalElement in class matrixTestsSparseTriplet assumes that the destination matrix "hiop::hiopMatrix* a" is a dense matrix.
See

hiop/tests/LinAlg/matrixTestsSparseTriplet.cpp

Line 86 in fcd29a0

auto mat = dynamic_cast<const hiop::hiopMatrixDense*>(A);

tests/LinAlg/matrixTestsSparseTriplet.cpp Show resolved Hide resolved

tests/LinAlg/matrixTestsSparse.hpp

                 virtual const local_ordinal_type* getRowIndices(const hiop::hiopMatrixSparse* a) = 0;
                 virtual const local_ordinal_type* getColumnIndices(const hiop::hiopMatrixSparse* a) = 0;
                 virtual local_ordinal_type getLocalSize(const hiop::hiopVector* x) = 0;
                 virtual int verifyAnswer(hiop::hiopMatrix* A, real_type answer) = 0;
+                virtual int verifyAnswer(hiop::hiopMatrix* A, local_ordinal_type nnz_st, local_ordinal_type nnz_ed, const double answer) = 0;

Contributor

ashermancinelli Oct 7, 2020

I think if you're just checking a single element, probably better to use getLocalElement and check result

tests/LinAlg/matrixTestsSparse.hpp


		auto val = getMatrixData(&B);

		fail += verifyAnswer(&B,0,B_nnz_st,B_val);

Contributor

ashermancinelli Oct 7, 2020

To be more thorough, we should have another verifyAnswer that takes a lambda for sparse mats, eg

fail += verifyAnswer(&B,
  [=] (local_ordinal_type i, local_ordinal_type j) -> real_type
  {
    if(i==0 && j==B_nnz_st) return B_val;
    else if // other two conditions...
    else if // other two conditions...
    else return 0.;
  });

Collaborator Author

nychiang Oct 7, 2020

I am trying to check multiple values in a sparse matrix, not a single value.
The nonzero indices for these nonzero starts from nnz_st to nnz_ed.
I can change it to a lambda function like you mentioned, but I need to assign it to a different name, or use different declaration.
This is because the current implementation assumes the input matrix is dense

hiop/tests/LinAlg/matrixTestsSparseTriplet.cpp

Lines 211 to 224 in fcd29a0

    
           int MatrixTestsSparseTriplet::verifyAnswer( 
        
               hiop::hiopMatrix* Amat, 
        
               std::function<real_type(local_ordinal_type, local_ordinal_type)> expect) 
        
           { 
        
             auto A = dynamic_cast<hiop::hiopMatrixDense*>(Amat); 
        
             assert(A->get_local_size_n() == A->n() && "Matrix should not be distributed"); 
        
             const local_ordinal_type M = A->get_local_size_m(); 
        
             const local_ordinal_type N = A->get_local_size_n(); 
        
             int fail = 0; 
        
             for (local_ordinal_type i=0; i<M; i++) 
        
             { 
        
               for (local_ordinal_type j=0; j<N; j++) 
        
               { 
        
                 if (!isEqual(getLocalElement(A, i, j), expect(i, j)))

Contributor

ashermancinelli Oct 7, 2020

I'll submit a pr against your branch to fix this

Collaborator Author

nychiang Oct 7, 2020

sure. I tried to not touch your existing code and only add my code into the repository.
I believer it will be really helpful if we can fix this issue.
About verifying the answer for a sparse matrix, I do have some different idea. Currently we only verify the nonzero values given by the nonzero index, but I think we may also want to verify the nonzero pattern, For example, verify the kth nonzero is located on some particular row and col. This may be not urgent, since now we assume the sparsity pattern won't change. If we use sparse matrices correctly, we can skip check the pattern. However, if there is an error in the sparsity pattern when copy to/from another matrix/vector, it will be hard to detect. Please let me know your opinion about this. Cheers,

Collaborator

pelesh Oct 8, 2020

The way I see it there is little use for a getElement method for a sparse matrix simply because such method is cumbersome to implement and use. Method verifyAnswer is always implementation specific, so it should access matrix data directly, no need to use getElement or similar.

Collaborator

pelesh Oct 8, 2020

Verifying sparsity pattern should not be difficult either. Again, I would implement that by accessing matrix data directly and not try to implement getElement method for that. Method getElement is design to retrieve and check only one single element in a vector or matrix. If you are checking multiple elements, you should be accessing matrix data directly. This is also true for dense matrices.

Collaborator Author

nychiang commented Oct 7, 2020 via email •

edited

I thought about it. The problem is that in this routine, there is no way to get the ‘correct’ number of entries in each row. See line 117 in the main function testMatrixSparse.cpp: initializeSparseTriplet(mxn_sparse, entries_per_row); which sets the number of entries in each row of a given matrix. If we can pass this information to the LinAlg layer, I don’t need to do this loop to know how many entries will be copied from A to B. Instead, I can simply set auto nnz_A_need_to_copy = n_rows * entries_per_row.

Contributor

ashermancinelli commented Oct 7, 2020

Could you respond via github in the future? Email responses create a lot of garbage. It would be appropriate to pass num entries per row to the test from the test driver.

Collaborator Author

nychiang commented Oct 7, 2020

I noticed it. :-)
just deleted all the email garbage in my reply.

nychiang mentioned this pull request

addToSymDenseMatrixUpperTriangle and transAddToSymDenseMatrixUpperTriangle for sparse matrix have never been used #98

Closed

nychiang force-pushed the nywrk branch 3 times, most recently from d2ffa0d to 50cf63b Compare

October 20, 2020 22:05

Collaborator

pelesh commented Nov 4, 2020

I believe this PR should target develop rather than master branch.

cnpetra marked this pull request as ready for review

November 10, 2020 23:52

cnpetra changed the base branch from master to develop

November 10, 2020 23:57

Nai-Yuan Chiang and others added 9 commits

November 10, 2020 19:12


          add support to sparse linear algebra. Code can be compiled. Without t…

…est for functionalities


          Add example for sparse case; Fix bugs and white space issue.

b250b8f


          remove hard-coded path from CMake file

6fae57d


          correct IpoptAdapter for sparse matrix. Add unit test

3fb3391


          add example 7

01b1d6e


          add unit test; add opt obj for selfcheck

4d70343


          add unit test; build sparse solver by default

ee480cd


          fix conflics due to git rebase

469243a


          add support for STRUMPACK

56b4fec

nychiang added 5 commits

November 10, 2020 19:14


          add STRUMPACK for sparse indefinite sys. fix CMAKE files

2de4e42


          fix bugs for full sparse KKT system

fee555e


          fix CMake file and bugs in formulating CSR matrix

d84033d


          change function names

8319d18


          correct index when print out the sparse matrix

nychiang force-pushed the nywrk branch from a850fd0 to 4898619 Compare

November 11, 2020 03:25

nychiang and others added 5 commits

November 10, 2020 19:40


          fix index when print out dense matrix

ad0a6dd


          build system updates

acaf4c5

- builds without STRUMPACK (SparseKKT class was updated)
- add "dl" when coinhsl is linked with
- left some further todo items


          fix bugs in system xycyd, add make test, find STRUMPACK module

268b1e6


          small touches related to STRUMPACK

9a41e03


          fix conditional jump or move depends on uninitialised value in MA57

d9fd2e3

Collaborator

cnpetra commented Nov 14, 2020

Does this also address issue #79 ?

Collaborator Author

nychiang commented Nov 14, 2020

NO. The patch is ready. Do you want to address issue #79 in this PR?

Collaborator

cnpetra commented Nov 14, 2020 •

edited

NO. The patch is ready. Do you want to address issue #79 in this PR?

let's address it after the merge as this PR blocks other efforts as well and it is better to close it asap

cnpetra merged commit ec5bf1d into develop

cnpetra deleted the nywrk branch

April 7, 2021 21:34

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment