Skip to content

Commit fb60fe3

Browse files
committed
also update docs
1 parent afa64bc commit fb60fe3

File tree

1 file changed

+174
-7
lines changed

1 file changed

+174
-7
lines changed

docs/src/user_interface/decompositions.md

Lines changed: 174 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -169,38 +169,205 @@ PolarNewton
169169

170170
## Orthogonal Subspaces
171171

172-
Often it is useful to compute orthogonal bases for a particular subspace defined by a matrix.
173-
Given a matrix `A` we can compute an orthonormal basis for its image or coimage, and factorize the matrix accordingly.
172+
Often it is useful to compute orthogonal bases for particular subspaces defined by a matrix.
173+
Given a matrix `A`, we can compute an orthonormal basis for its image or coimage, and factorize the matrix accordingly.
174174
These bases are accessible through [`left_orth`](@ref) and [`right_orth`](@ref) respectively.
175-
This is implemented through a combination of the decompositions mentioned above, and serves as a convenient interface to these operations.
175+
176+
### Overview
177+
178+
The [`left_orth`](@ref) function computes an orthonormal basis `V` for the image (column space) of `A`, along with a corestriction matrix `C` such that `A = V * C`.
179+
The resulting `V` has orthonormal columns (`V' * V ≈ I` or `isisometric(V)`).
180+
181+
Similarly, [`right_orth`](@ref) computes an orthonormal basis for the coimage (row space) of `A`, i.e., the image of `A'`.
182+
It returns matrices `C` and `Vᴴ` such that `A = C * Vᴴ`, where `V = (Vᴴ)'` has orthonormal columns (`isisometric(Vᴴ; side = :right)`).
183+
184+
These functions serve as high-level interfaces that automatically select the most appropriate decomposition based on the specified options, making them convenient for users who want orthonormalization without worrying about the underlying implementation details.
176185

177186
```@docs; canonical=false
178187
left_orth
179188
right_orth
180189
```
181190

182-
In order to dispatch to the underlying factorizations, the following wrapper algorithms are used:
191+
### Algorithm Selection
192+
193+
Both functions support multiple decomposition drivers, which can be selected through the `alg` keyword argument:
194+
195+
**For `left_orth`:**
196+
- `alg = :qr` (default without truncation): Uses QR decomposition via [`qr_compact`](@ref)
197+
- `alg = :polar`: Uses polar decomposition via [`left_polar`](@ref)
198+
- `alg = :svd` (default with truncation): Uses SVD via [`svd_compact`](@ref) or [`svd_trunc`](@ref)
199+
200+
**For `right_orth`:**
201+
- `alg = :lq` (default without truncation): Uses LQ decomposition via [`lq_compact`](@ref)
202+
- `alg = :polar`: Uses polar decomposition via [`right_polar`](@ref)
203+
- `alg = :svd` (default with truncation): Uses SVD via [`svd_compact`](@ref) or [`svd_trunc`](@ref)
204+
205+
When `alg` is not specified, the function automatically selects `:qr`/`:lq` for exact orthogonalization, or `:svd` when a truncation strategy is provided.
206+
207+
### Extending with Custom Algorithms
208+
209+
To register a custom algorithm type for use with these functions, you need to define the appropriate conversion function, for example:
210+
211+
```julia
212+
# For left_orth
213+
MatrixAlgebraKit.left_orth_alg(alg::MyCustomAlgorithm) = LeftOrthAlgorithm{:qr}(alg)
214+
215+
# For right_orth
216+
MatrixAlgebraKit.right_orth_alg(alg::MyCustomAlgorithm) = RightOrthAlgorithm{:lq}(alg)
217+
```
218+
219+
The type parameter (`:qr`, `:lq`, `:polar`, or `:svd`) indicates which factorization backend will be used.
220+
The wrapper algorithm types handle the dispatch to the appropriate implementation:
183221

184222
```@docs; canonical=false
223+
left_orth_alg
224+
right_orth_alg
185225
LeftOrthAlgorithm
186226
RightOrthAlgorithm
187227
```
188228

229+
### Examples
230+
231+
Basic orthogonalization:
232+
233+
```jldoctest orthnull; output=false
234+
using MatrixAlgebraKit
235+
using LinearAlgebra
236+
237+
A = [1.0 2.0; 3.0 4.0; 5.0 6.0]
238+
V, C = left_orth(A)
239+
(V' * V) ≈ I && A ≈ V * C
240+
241+
# output
242+
true
243+
```
244+
245+
Using different algorithms:
246+
247+
```jldoctest orthnull; output=false
248+
A = randn(4, 3)
249+
V1, C1 = left_orth(A; alg = :qr)
250+
V2, C2 = left_orth(A; alg = :polar)
251+
V3, C3 = left_orth(A; alg = :svd)
252+
A ≈ V1 * C1 ≈ V2 * C2 ≈ V3 * C3
253+
254+
# output
255+
true
256+
```
257+
258+
With truncation:
259+
260+
```jldoctest orthnull; output=false
261+
A = [1.0 0.0; 0.0 1e-10; 0.0 0.0]
262+
V, C = left_orth(A; trunc = (atol = 1e-8,))
263+
size(V, 2) == 1 # Only one column retained
264+
265+
# output
266+
true
267+
```
268+
189269

190270
## Null Spaces
191271

192272
Similarly, it can be convenient to obtain an orthogonal basis for the kernel or cokernel of a matrix.
193-
These are the compliments of the coimage and image, respectively, and can be computed using the [`left_null`](@ref) and [`right_null`](@ref) functions.
194-
Again, this is typically implemented through a combination of the decompositions mentioned above, and serves as a convenient interface to these operations.
273+
These are the complements of the coimage and image, respectively, and can be computed using the [`left_null`](@ref) and [`right_null`](@ref) functions.
274+
275+
### Overview
276+
277+
The [`left_null`](@ref) function computes an orthonormal basis `N` for the cokernel (left nullspace) of `A`, which is the nullspace of `A'`.
278+
This means `A' * N ≈ 0` and `N' * N ≈ I`.
279+
280+
Similarly, [`right_null`](@ref) computes an orthonormal basis for the kernel (right nullspace) of `A`.
281+
It returns `Nᴴ` such that `A * Nᴴ' ≈ 0` and `Nᴴ * Nᴴ' ≈ I`, where `N = (Nᴴ)'` has orthonormal columns.
282+
283+
These functions automatically handle rank determination and provide convenient access to nullspace computation without requiring detailed knowledge of the underlying decomposition methods.
195284

196285
```@docs; canonical=false
197286
left_null
198287
right_null
199288
```
200289

201-
Again, dispatching happens through the following wrapper algorithm types:
290+
### Algorithm Selection
291+
292+
Both functions support multiple decomposition drivers, which can be selected through the `alg` keyword argument:
293+
294+
**For `left_null`:**
295+
- `alg = :qr` (default without truncation): Uses QR-based nullspace computation via [`qr_null`](@ref)
296+
- `alg = :svd` (default with truncation): Uses SVD via [`svd_full`](@ref) with appropriate truncation
297+
298+
**For `right_null`:**
299+
- `alg = :lq` (default without truncation): Uses LQ-based nullspace computation via [`lq_null`](@ref)
300+
- `alg = :svd` (default with truncation): Uses SVD via [`svd_full`](@ref) with appropriate truncation
301+
302+
When `alg` is not specified, the function automatically selects `:qr`/`:lq` for exact nullspace computation, or `:svd` when a truncation strategy is provided to handle numerical rank determination.
303+
304+
!!! note
305+
For nullspace functions, [`notrunc`](@ref) has special meaning when used with the default QR/LQ algorithms.
306+
It indicates that the nullspace should be computed from the exact zeros determined by the additional rows/columns of the extended matrix, without any tolerance-based truncation.
307+
308+
### Extending with Custom Algorithms
309+
310+
To register a custom algorithm type for use with these functions, you need to define the appropriate conversion function:
311+
312+
```julia
313+
# For left_null
314+
MatrixAlgebraKit.left_null_alg(alg::MyCustomAlgorithm) = LeftNullAlgorithm{:qr}(alg)
315+
316+
# For right_null
317+
MatrixAlgebraKit.right_null_alg(alg::MyCustomAlgorithm) = RightNullAlgorithm{:lq}(alg)
318+
```
319+
320+
The type parameter (`:qr`, `:lq`, or `:svd`) indicates which factorization backend will be used.
321+
The wrapper algorithm types handle the dispatch to the appropriate implementation:
202322

203323
```@docs; canonical=false
204324
LeftNullAlgorithm
205325
RightNullAlgorithm
326+
left_null_alg
327+
right_null_alg
328+
```
329+
330+
### Examples
331+
332+
Basic nullspace computation:
333+
334+
```jldoctest orthnull; output=false
335+
A = [1.0 2.0 3.0; 4.0 5.0 6.0] # Rank 2 matrix
336+
N = left_null(A)
337+
size(N) == (2, 0)
338+
339+
# output
340+
true
341+
```
342+
343+
```jldoctest orthnull; output=false
344+
Nᴴ = right_null(A)
345+
size(Nᴴ) == (1, 3) && norm(A * Nᴴ') < 1e-14 && isisometric(Nᴴ; side = :right)
346+
347+
# output
348+
true
349+
```
350+
351+
Computing nullspace with rank detection:
352+
353+
```jldoctest orthnull; output=false
354+
A = [1.0 2.0; 2.0 4.0; 3.0 6.0] # Rank 1 matrix (second column = 2*first)
355+
N = left_null(A; alg = :svd, trunc = (atol = 1e-10,))
356+
size(N) == (3, 2) && norm(A' * N) < 1e-12 && isisometric(N)
357+
358+
# output
359+
true
360+
```
361+
362+
Using different algorithms:
363+
364+
```jldoctest orthnull; output=false
365+
A = [1.0 0.0 0.0; 0.0 1.0 0.0]
366+
N1 = right_null(A; alg = :lq)
367+
N2 = right_null(A; alg = :svd)
368+
norm(A * N1') < 1e-14 && norm(A * N2') < 1e-14 &&
369+
isisometric(N1; side = :right) && isisometric(N2; side = :right)
370+
371+
# output
372+
true
206373
```

0 commit comments

Comments
 (0)