Fix HyKKT solver reuse and add export support - #464
Conversation
andrewxu319
left a comment
There was a problem hiding this comment.
I can't run it right now, but I left comments.
| << ", allocated = " << allocated_ | ||
| << "\n"; | ||
| std::cout << "ERROR: Changing J_d between empty and nonempty is not " | ||
| "supported when reusing HyKKT.\n"; |
There was a problem hiding this comment.
Why not put this check directly inside setMatrixBlocks()?
Also, I'm pretty sure that if J_d is empty but J_d_flag_ is true, everything would still work. So if on the first solve J_d is nonempty, then on the second solve J_d is empty, you can keep J_d_flag_ = true and maybe throw a warning about algorithm inefficiency.
There was a problem hiding this comment.
setMatrixBlocks() already performs the compatibility check and records the result in status_. Since it returns void, solve() is where the error is returned.
I tried keeping J_d_flag_ = true after replacing a nonempty J_d with an empty one. The next solve segfaults in Sparse::copyValues(). Reading that path, J_d_scaled_ still has the original nonempty structure while the source is now empty.
I've left the empty/nonempty J_d transition unsupported here. The proposed setup() refactor would be the right place to rebuild the cached state that would make it work.
There was a problem hiding this comment.
solve() returns real_type, so it seems strange to return 1 as an error code. I was suggesting putting the error message directly inside setMatrixBlocks(), and if you want to return something, it makes more sense to return it as an int in setMatrixBlocks().
There was a problem hiding this comment.
I had kept the existing status_ flow because setMatrixBlocks() returned void, so the compatibility result was carried into solve(). I understand now that you were suggesting changing the setter API itself.
setMatrixBlocks() now returns an int and performs the check before updating the stored pointers, so a rejected call leaves the solver in its previous valid state. The error is reported through out::error() rather than std::cout. I removed status_ and the guard from solve(), and the test now checks the setter result directly.
There was a problem hiding this comment.
Remove these two lines because they're redundant with 413-414
There was a problem hiding this comment.
Removed them. I also moved loadResultMatrix() into computeSpGEMMHgamma(), since the HIP implementation needs the product and sum dimensions before it initializes the result descriptor.
| matrixHandler_->matvec(J_perm_, omega_perm_, schur_, &ONE, &MINUS_ONE, memspace_); | ||
|
|
||
| sccg_ = new SchurComplementConjugateGradient(J_->getNumRows(), J_->getNumColumns(), cholesky_, matrixHandler_, vectorHandler_, memspace_); | ||
| if (!allocated_) |
There was a problem hiding this comment.
Move if (!allocated_) to inside hykkt::HyKKTSolver::solve() (line 207) for consistency.
| } | ||
|
|
||
| y_->setToZero(memspace_); | ||
| z_->setToZero(memspace_); |
There was a problem hiding this comment.
I think you can move lines 106-115 to SchurComplementConjugateGradient::solve() and remove if (!y_). Nobody is going to call this multiple times on the same RHS. Probably doesn't matter that much, but I feel like that's more logical because otherwise the user needs to call setup() before every solve().
|
|
||
| void setAlpha(real_type alpha) | ||
| { | ||
| impl_->setAlpha(alpha); |
There was a problem hiding this comment.
Move this to SpGEMM.cpp
There was a problem hiding this comment.
Same for SpGEMMCpu.hpp, etc
There was a problem hiding this comment.
Also, might as well rename it setCoefficients or setConstants and let it set both alpha and beta
There was a problem hiding this comment.
Renamed it to setCoefficients(alpha, beta) and updated all three backends to set both values. I added a test that changes beta and checks that the result changes on CPU, CUDA, and HIP. I also removed beta_ = beta from the CPU backend to make sure the test catches it, and it failed as expected.
| status *= validateResult(error, tol); | ||
|
|
||
| // Replace D_s and restore data modified by the first solve. | ||
| std::ifstream D_s_reuse_file(D_s_file_name); |
There was a problem hiding this comment.
Why make these reuse matrices/vectors instead of just modifying the data in J, r_x, r_s, etc?
There was a problem hiding this comment.
J, r_x, r_s, r_y, and r_yd now update in place instead of creating separate reuse copies. I kept D_s as a replacement because the fix refreshes its cached values pointer. Updating the same object would not test that path.
| real_type second_error = hykktSolver.solve(); | ||
| status *= validateResult(second_error, tol); | ||
|
|
||
| // Verify an exact zero RHS is handled without producing NaNs. |
There was a problem hiding this comment.
| // Verify an exact zero RHS is handled without producing NaNs. | |
| // Check that a zero RHS doesn't result in NaNs. |
| real_type zero_rhs_error = hykktSolver.solve(); | ||
| status *= validateResult(zero_rhs_error, tol); | ||
|
|
||
| // Changing J_d between nonempty and empty invalidates cached solver data. |
There was a problem hiding this comment.
| // Changing J_d between nonempty and empty invalidates cached solver data. | |
| // Check that the solver raises an error when trying to change J_d from nonempty to empty. |
There was a problem hiding this comment.
| * @pre J_d_flag_ determines if variables used for Spgemm H_tilde |
There was a problem hiding this comment.
Applied. I reworded it slightly for clarity.
| status_ = true; // when using API, we can't check if sparsity pattern changed | ||
| J_d_flag_ = J_d_flag; | ||
| // Arbitrary sparsity changes remain the caller's responsibility, but | ||
| // switching between empty and nonempty J_d invalidates cached HyKKT data. |
There was a problem hiding this comment.
Probably for another PR, but I think we should make a separate HyKKTSolver::setup() function and remove the allocated_ flag. This matches other parts of the codebase better. This also allows switching J_d_ from nonempty to empty as long as the user calls setup() first.
There was a problem hiding this comment.
This also means you don't need to make a new solver at line 224 of HyKKTSolverTests.hpp, for example.
There was a problem hiding this comment.
Agreed. I think that makes sense as a follow-up. For now, I kept the separate solver in the test so the empty J_d case still covers initialization and reuse.
There was a problem hiding this comment.
Never mind, not possible. Some of the "setup" steps depend on the results of earlier "compute" steps.
|
Thanks for the review @andrewxu319 I appreciate it! I made the SCCG lifecycle, SpGEMM, and reuse-test changes. All tests are passing on my end. |
shakedregev
left a comment
There was a problem hiding this comment.
Looks good, testing.
The branch has been updated since this comment to move the matrix compatibility error handling into setMatrixBlocks(). Please use the latest version for testing. |
shakedregev
left a comment
There was a problem hiding this comment.
Looks good and works, but fix these warnings.
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:21:12: warning: 'loadProductMatrices' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
21 | void loadProductMatrices(matrix::Csr* A, matrix::Csr* B);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:28:20: note: overridden virtual function is here
28 | virtual void loadProductMatrices(matrix::Csr* A, matrix::Csr* B) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:22:12: warning: 'loadSumMatrix' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
22 | void loadSumMatrix(matrix::Csr* D);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:29:20: note: overridden virtual function is here
29 | virtual void loadSumMatrix(matrix::Csr* D) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:23:12: warning: 'loadResultMatrix' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
23 | void loadResultMatrix(matrix::Csr** E_ptr);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:30:20: note: overridden virtual function is here
30 | virtual void loadResultMatrix(matrix::Csr** E_ptr) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:25:12: warning: 'compute' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
25 | void compute();
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:32:20: note: overridden virtual function is here
32 | virtual void compute() = 0;
| ^
4 warnings generated when compiling for gfx90a.
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMM.cpp:8:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:21:12: warning: 'loadProductMatrices' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
21 | void loadProductMatrices(matrix::Csr* A, matrix::Csr* B);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:28:20: note: overridden virtual function is here
28 | virtual void loadProductMatrices(matrix::Csr* A, matrix::Csr* B) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMM.cpp:8:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:22:12: warning: 'loadSumMatrix' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
22 | void loadSumMatrix(matrix::Csr* D);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:29:20: note: overridden virtual function is here
29 | virtual void loadSumMatrix(matrix::Csr* D) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMM.cpp:8:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:23:12: warning: 'loadResultMatrix' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
23 | void loadResultMatrix(matrix::Csr** E_ptr);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:30:20: note: overridden virtual function is here
30 | virtual void loadResultMatrix(matrix::Csr** E_ptr) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMM.cpp:8:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:25:12: warning: 'compute' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
25 | void compute();
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:32:20: note: overridden virtual function is here
32 | virtual void compute() = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMM.cpp:12:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.hpp:27:12: warning: 'loadProductMatrices' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
27 | void loadProductMatrices(matrix::Csr* A, matrix::Csr* B);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:28:20: note: overridden virtual function is here
28 | virtual void loadProductMatrices(matrix::Csr* A, matrix::Csr* B) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMM.cpp:12:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.hpp:28:12: warning: 'loadSumMatrix' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
28 | void loadSumMatrix(matrix::Csr* D);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:29:20: note: overridden virtual function is here
29 | virtual void loadSumMatrix(matrix::Csr* D) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMM.cpp:12:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.hpp:29:12: warning: 'loadResultMatrix' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
29 | void loadResultMatrix(matrix::Csr** E_ptr);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:30:20: note: overridden virtual function is here
30 | virtual void loadResultMatrix(matrix::Csr** E_ptr) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMM.cpp:12:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.hpp:31:12: warning: 'compute' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
31 | void compute();
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:32:20: note: overridden virtual function is here
32 | virtual void compute() = 0;
| ^
8 warnings generated when compiling for gfx90a.
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.hpp:27:12: warning: 'loadProductMatrices' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
27 | void loadProductMatrices(matrix::Csr* A, matrix::Csr* B);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:28:20: note: overridden virtual function is here
28 | virtual void loadProductMatrices(matrix::Csr* A, matrix::Csr* B) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.hpp:28:12: warning: 'loadSumMatrix' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
28 | void loadSumMatrix(matrix::Csr* D);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:29:20: note: overridden virtual function is here
29 | virtual void loadSumMatrix(matrix::Csr* D) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.hpp:29:12: warning: 'loadResultMatrix' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
29 | void loadResultMatrix(matrix::Csr** E_ptr);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:30:20: note: overridden virtual function is here
30 | virtual void loadResultMatrix(matrix::Csr** E_ptr) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMHip.hpp:31:12: warning: 'compute' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
31 | void compute();
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:32:20: note: overridden virtual function is here
32 | virtual void compute() = 0;
| ^
4 warnings generated when compiling for gfx90a.
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:21:12: warning: 'loadProductMatrices' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
21 | void loadProductMatrices(matrix::Csr* A, matrix::Csr* B);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:28:20: note: overridden virtual function is here
28 | virtual void loadProductMatrices(matrix::Csr* A, matrix::Csr* B) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:22:12: warning: 'loadSumMatrix' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
22 | void loadSumMatrix(matrix::Csr* D);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:29:20: note: overridden virtual function is here
29 | virtual void loadSumMatrix(matrix::Csr* D) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:23:12: warning: 'loadResultMatrix' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
23 | void loadResultMatrix(matrix::Csr** E_ptr);
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:30:20: note: overridden virtual function is here
30 | virtual void loadResultMatrix(matrix::Csr** E_ptr) = 0;
| ^
In file included from /ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.cpp:7:
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMCpu.hpp:25:12: warning: 'compute' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
25 | void compute();
| ^
/ccs/home/regevs/ReSolve_dir/ReSolve/resolve/hykkt/spgemm/SpGEMMImpl.hpp:32:20: note: overridden virtual function is here
32 | virtual void compute() = 0;
Fixed the missing |
Description
This PR exports the HyKKT libraries and fixes issues found while reusing
HyKKTSolverthrough the HiOp integration.The existing HyKKT tests primarily covered a single solve. Repeated solves with the same sparsity pattern and updated numerical values exposed stale matrix data, cached solver state, repeated allocations, undersized GPU transpose workspaces, and zero-residual handling issues.
@shakedregev
Proposed changes
HyKKTSolveris reused:D_svalues pointer before each solve.gammachanges.C_nnz_variable shadowed the class member, so the stored output nonzero count was not updated.H_tilde_allocation when solving withoutJ_d.J_dstate and reject reuse whenJ_dchanges between empty and nonempty. Reuse still requires unchanged sparsity patterns.Jand RHS data, replacing and updatingD_s, and changinggamma.J_d.J_dchanges from nonempty to empty.Checklist
make testandmake test_installper testing instructions). Code tested on./examples/<your_example>.exe -hto get instructions how to run examples). Code tested on:-Wall -Wpedantic -Wconversion -Wextra.Further comments
The transpose workspace change keeps the current workspace design and grows the CUDA and HIP transpose buffers when needed rather than implementing the broader redesign proposed in #345.