From: Michael Orlitzky Date: Fri, 2 Oct 2020 22:41:30 +0000 (-0400) Subject: mjo/ldlt.py: fix two bugs in the imperative block_ldlt() function. X-Git-Url: http://gitweb.michael.orlitzky.com/?a=commitdiff_plain;h=78a924d7f6b518daa3886b3356ae9a4af86ebee1;p=sage.d.git mjo/ldlt.py: fix two bugs in the imperative block_ldlt() function. --- diff --git a/mjo/ldlt.py b/mjo/ldlt.py index 6b99083..c8483c7 100644 --- a/mjo/ldlt.py +++ b/mjo/ldlt.py @@ -409,28 +409,28 @@ def block_ldlt(A): p[k] = p[s] p[s] = p_k - # Now the pivot is in the (k,k)th position. - d.append( matrix(ring, 1, [[A[k,k]]]) ) - - # Compute the Schur complement that we'll work on during - # the following iteration, and store it back in the lower- - # right-hand corner of "A". - for i in range(n-k-1): - for j in range(i+1): - A[k+1+j,k+1+i] = ( A[k+1+j,k+1+i] - - A[k,k+1+j]*A[k,k+1+i]/alpha ) - A[k+1+i,k+1+j] = A[k+1+j,k+1+i] # keep it symmetric! - - for i in range(n-k-1): - # Store the new (kth) column of "L" within the lower- - # left-hand corner of "A", being sure to set the lower- - # left entries from the upper-right ones to avoid - #collisions. - A[k+i+1,k] = A[k,k+1+i]/alpha - - # No return value, only the desired side effects of updating - # p, d, and A. - return + # Now the pivot is in the (k,k)th position. + d.append( matrix(ring, 1, [[A[k,k]]]) ) + + # Compute the Schur complement that we'll work on during + # the following iteration, and store it back in the lower- + # right-hand corner of "A". + for i in range(n-k-1): + for j in range(i+1): + A[k+1+j,k+1+i] = ( A[k+1+j,k+1+i] - + A[k,k+1+j]*A[k,k+1+i]/A[k,k] ) + A[k+1+i,k+1+j] = A[k+1+j,k+1+i] # keep it symmetric! + + for i in range(n-k-1): + # Store the new (kth) column of "L" within the lower- + # left-hand corner of "A", being sure to set the lower- + # left entries from the upper-right ones to avoid + # collisions. + A[k+i+1,k] = A[k,k+1+i]/A[k,k] + + # No return value, only the desired side effects of updating + # p, d, and A. + return k = 0 while k < n: