grey/BOSL2
1
0
Fork 0
mirror of https://github.com/BelfrySCAD/BOSL2.git synced 2025-01-01 09:49:45 +00:00
Code Issues Projects Releases Packages Wiki Activity Actions

explore new scope of transpose in linear_solve and _qr_factor

Browse source
This commit is contained in:
RonaldoCMP 2020-08-11 15:21:36 +01:00
parent 0dc8bf6c8f
commit 2e9d87f584
1 changed files with 12 additions and 13 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

25
math.scad
View file

@ -698,9 +698,8 @@ function linear_solve(A,b) =
zeros != [] ? [] : zeros != [] ? [] :
m<n m<n
// avoiding input validation in back_substitute // avoiding input validation in back_substitute
? let( n = len(R), ? let( n = len(R) )
Rt = [for(i=[0:n-1]) [for(j=[0:n-1]) R[n-1-j][n-1-i]]] ) Q*reverse(_back_substitute(transpose(R, reverse=true), reverse(b)))
Q*reverse(_back_substitute(Rt,reverse(b)))
: _back_substitute(R, transpose(Q)*b); : _back_substitute(R, transpose(Q)*b);
// Function: matrix_inverse() // Function: matrix_inverse()
@ -728,13 +727,14 @@ function qr_factor(A) =
n = len(A[0]) n = len(A[0])
) )
let( let(
qr =_qr_factor(A, column=0, m = m, n=n, Q=ident(m)), qr =_qr_factor(A, Q=ident(m), column=0, m = m, n=n),
Rzero = [ Rzero =
for(i=[0:m-1]) [ let( R = qr[1] )
for(j=[0:n-1]) [ for(i=[0:m-1]) [
i>j ? 0 : qr[1][i][j] let( ri = R[i] )
for(j=[0:n-1]) i>j ? 0 : ri[j]
] ]
] ]
) [qr[0],Rzero]; ) [qr[0],Rzero];
function _qr_factor(A,Q, column, m, n) = function _qr_factor(A,Q, column, m, n) =
@ -763,13 +763,12 @@ function back_substitute(R, b, transpose = false) =
let(n=len(R)) let(n=len(R))
assert(is_vector(b,n) || is_matrix(b,n),str("R and b are not compatible in back_substitute ",n, len(b))) assert(is_vector(b,n) || is_matrix(b,n),str("R and b are not compatible in back_substitute ",n, len(b)))
transpose transpose
? reverse(_back_substitute([for(i=[0:n-1]) [for(j=[0:n-1]) R[n-1-j][n-1-i]]], ? reverse(_back_substitute(transpose(R, reverse=true), reverse(b)))
reverse(b)))
: _back_substitute(R,b); : _back_substitute(R,b);
function _back_substitute(R, b, x=[]) = function _back_substitute(R, b, x=[]) =
let(n=len(R)) let(n=len(R))
len(x) == n ? x len(x) == n ? x
: let(ind = n - len(x) - 1) : let(ind = n - len(x) - 1)
R[ind][ind] == 0 ? [] R[ind][ind] == 0 ? []
: let( : let(