-
Notifications
You must be signed in to change notification settings - Fork 5
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.
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
Some data in the State variable may be corrupt #625
Comments
The |
Thanks Kyle for the clarification. I looked at these two lines https://github.com/NCAR/micm/blob/main/include/micm/solver/lu_decomposition.inl#L201 and https://github.com/NCAR/micm/blob/main/include/micm/solver/lu_decomposition.inl#L214. Could the |
We set them just before those, but it's in an if block so maybe not all of the values are zeros. Perhaps we should add back a zero of L and U |
Thanks @K20shores for your comments. Hmm, interesting. For the first time when the |
I've learned some things. Thanks to @mattldawson's suggestion, I set the L and U matrices to Inf and found out that sometimes we end up with infs after the LU decomposition. The issue seems to be related to the interplay between block size (which are grid cells) and the number of species. Sizes leading to valid L and U matrices
Sizes leading to L and U matrices with infs
|
Also, the math on geeks for geeks for this algorithm indicates that U should always take a value from the input matrix, but I can see in the debugger that we sometimes don't. I think this means we are missing some micm/include/micm/solver/lu_decomposition.inl Lines 56 to 66 in 94c22da
Missing values mean that when we set the U vector here we are subtracting a value from whatever the value of U is (in this case micm/include/micm/solver/lu_decomposition.inl Lines 197 to 203 in 94c22da
|
Thanks to @K20shores , we track down the issue further. Given an initial Jacobian matrix with size
The expected L matrix should be:
While we are constructing the third element of the fourth row, micm/include/micm/solver/lu_decomposition.inl Line 211 in 94c22da
Since we are initializing the L matrix with inf now, this leads to a inf - some value operation later and corrupts the result.micm/include/micm/solver/lu_decomposition.inl Line 214 in 94c22da
After discussion with Kyle, the lines
|
@jian helped me pin down the problem. The issue is that the sparsity pattern of the LU matrices doesn't necessarily match the sparsity pattern of the jacobian matrix. In the algorithm, the LU matrices need to be filled with the value from the A matrix when the LU matrices are nonzero. The first change we needed to make was to accurately record when to set the LU matrix values based off of when the LU matrices are zero, not the jacobian as was happening before. The second change is to record a sentinel value. This value indicates when the L or U matrix is nonzero and the jacobian is actually zero. In this case, we set the LU matrix to zero rather than copying from the jacobian matrix. |
As identified by this issue https://github.com/NCAR/MUSICA-Performance-Comparison/issues/62, if we reuse the
State
variable between different iterations ofsolve
function, the results are corrupt and we getfail to converge
error message later.Kyle commented that this was also the reason why he had to initialize the LU matrix to zero originally (#587).
This issue indicates that the LU matrix in
State
variable is not overwritten correctly between iterations (or we should not expect it to be overwritten correctly at all?).The text was updated successfully, but these errors were encountered: