Unexported (Private) functions.
This is a list of all unexported functions and types from Laplacians.
Laplacians.ApproxCholPQ
Laplacians.IJV
Laplacians.LDLinv
Laplacians.LLmatp
Laplacians.LLp
Laplacians.addToGPrime
Laplacians.approxCholPQDec!
Laplacians.approxCholPQInc!
Laplacians.approxchol_lapChol
Laplacians.augmentTreeOpt
Laplacians.blockSolver
Laplacians.complete_bipartite_graph_ijv
Laplacians.complete_graph_ijv
Laplacians.empty_graph_ijv
Laplacians.extendMatrix
Laplacians.firstn
Laplacians.forceLap
Laplacians.generalizedNecklace
Laplacians.getCutSet
Laplacians.initDictCol!
Laplacians.initGPrime
Laplacians.lapWrapComponents
Laplacians.lapWrapConnected
Laplacians.ldli2Chol
Laplacians.localBlockFlow
Laplacians.localFlow
Laplacians.path_graph_ijv
Laplacians.print_ll_col
Laplacians.print_ll_col
Laplacians.pure_random_ijv
Laplacians.pushSpeedResult!
Laplacians.rand_regular_bipartite
Laplacians.ring_graph_ijv
Laplacians.sampleByWeight
Laplacians.sddmWrapLap
Laplacians.sortSet
Laplacians.star_graph_ijv
Laplacians.testZeroDiag
Laplacians.treeDepthDFS
Laplacians.uniformWeight_ver
Laplacians.wrapCapture
Laplacians.wrapCaptureRhs
Laplacians.wrapInterface
Laplacians.ApproxCholPQ
— Type.An approximate priority queue. Items are bundled together into doubly-linked lists with all approximately the same key. minlist is the min list we know to be non-empty. It should always be a lower bound. keyMap maps keys to lists
Laplacians.IJV
— Method.ijv = IJV(A::SparseMatrixCSC)
Convert a sparse matrix to an IJV.
Laplacians.LDLinv
— Type.LDLinv contains the information needed to solve the Laplacian systems. It does it by applying Linv, then Dinv, then Linv (transpose). But, it is specially constructed for this particular solver. It does not explicitly make the matrix triangular. Rather, col[i] is the name of the ith col to be eliminated
Laplacians.LLmatp
— Type.LLmatp is the data structure used to maintain the matrix during elimination. It stores the elements in each column in a singly linked list (only next ptrs) Each element is an LLp (linked list pointer). The head of each column is pointed to by cols.
We probably can get rid of degs - as it is only used to store initial degrees.
Laplacians.LLp
— Type.LLp elements are all in the same column. row tells us the row, and val is the entry. val is set to zero for some edges that we should remove. next gives the next in the column. It points to itself to terminate. reverse is the index into lles of the other copy of this edge, since every edge is stored twice as we do not know the order of elimination in advance.
Laplacians.addToGPrime
— Method.Add a new vertex to GPrime
Laplacians.approxCholPQDec!
— Method.Decrement the key of element i
This could crash if i exceeds the maxkey
Laplacians.approxCholPQInc!
— Method.Increment the key of element i
This could crash if i exceeds the maxkey
Laplacians.approxchol_lapChol
— Method.This variation of approxChol creates a cholesky factor to do the elimination. It has not yet been optimized, and does not yet make the cholesky factor lower triangular
Laplacians.augmentTreeOpt
— Method.B = augmentTreeOpt{Tv,Ti}(tree, A, params)
Takes as input a tree and an adjacency matrix of a graph. It then computes the stretch of every edge of the graph wrt the tree. It uses cholmod to decide how many edge to add back, shooting for nnzLfac times n entries in the factored augmented tree, with a number of flops to factor equal to nnz(a)*flopsfac. The edges to add back are then choen at random.
Laplacians.blockSolver
— Method.Apply the ith solver on the ith component
Laplacians.complete_bipartite_graph_ijv
— Method.ijv = complete_bipartite_graph_ijv(n)
Laplacians.complete_graph_ijv
— Method.ijv = complete_graph_ijv(n)
Laplacians.empty_graph_ijv
— Method.ijv = empty_graph_ijv(n)
Laplacians.extendMatrix
— Method.Add a new vertex to a with weights to the other vertices corresponding to diagonal surplus weight.
This is an efficient way of writing [a d; d' 0]
Laplacians.firstn
— Method.b = firstn(a::IJV, n::Integer)
Only keep the first n vertices of a.
Laplacians.forceLap
— Method.la = forceLap(a)
Create a Laplacian matrix from an adjacency matrix. If the input looks like a Laplacian, throw a warning and convert it.
Laplacians.generalizedNecklace
— Method.graph = generalizedNecklace(A, H, k::Int64)
Constructs a generalized necklace graph starting with two graphs A and H. The resulting new graph will be constructed by expanding each vertex in H to an instance of A. k random edges will be generated between components. Thus, the resulting graph may have weighted edges.
Laplacians.getCutSet
— Method.Get the min cut from the source - return all vertices in the cut besides the source
Laplacians.initDictCol!
— Method.initDictCol!(dic, name, typ)
For a dictionary in which each key indexes an array. If dic does not contain an entry of name
, create with set to Array(typ,0)
.
Laplacians.initGPrime
— Method.Initialize GPrime with the set A and edges of type s->u
Laplacians.lapWrapComponents
— Method.f = lapWrapComponents(solver, a::AbstractArray; tol::Real=1e-6, maxits=Inf, maxtime=Inf, verbose=false, pcgIts=Int[], params...)
Applies a Laplacian solver
that satisfies our interface to each connected component of the graph with adjacency matrix a
. Passes kwargs on the solver.
Laplacians.lapWrapConnected
— Method.f = lapWrapConnected(sddmSolver, a::AbstractMatrix; kwargs...)
Applies a sddmSolver
to the Laplacian of the adjacency matrix a
of a connected graph. Passes on kwargs to the solver. sddmSolver
should be a solver that obeys the interface.
Laplacians.ldli2Chol
— Method.L = ldli2Chol(ldli)
This produces a matrix L so that L L^T approximate the original Laplacians. It is not quite a Cholesky factor, because it is off by a perm (and the all-1s vector orthogonality.
Laplacians.localBlockFlow
— Method.Compute block flow between s and t
Laplacians.localFlow
— Method.The LocalFlow function, from the Orecchia-Zhu paper
Laplacians.path_graph_ijv
— Method.ijv = path_graph_ijv(n::Int64)
Laplacians.print_ll_col
— Method.Print a column in an LLMatOrd matrix. This is here for diagnostics.
Laplacians.print_ll_col
— Method.Print a column in an LLmatp matrix. This is here for diagnostics.
Laplacians.pure_random_ijv
— Method.a = pure_random_ijv(n::Integer; verbose=false, prefix="", ver=Vcur)
Chooses among pathgraph, ringgraph, gridgraph, completebinarytree, randgenring, growngraph and ErdosRenyiClusterFix. It can produce a disconnected graph. For code that always produces a connected graph (and is the same as with Julia v0.6, use purerandomijv_v6)
Laplacians.pushSpeedResult!
— Method.ret
is the answer returned by a speed test. This pushed it into the dictionary on which we are storing the tests.
Laplacians.rand_regular_bipartite
— Method.a = rand_regular_bipartite(n,k)
Random k-regular bipartite graph between two sets of n vertices. No repeat edges, so can take a long time to build of k is close to n.
Returns a (possibly) asymmetric matrix that contains the upper-right block.
Laplacians.ring_graph_ijv
— Method.ijv = ring_graph_ijv(n)
Laplacians.sampleByWeight
— Method.ind = sampleByWeight(wt; ver=Vcur)
sample an index with probability proportional to its weight given here
Laplacians.sddmWrapLap
— Method.f = sddmWrapLap(lapSolver, sddm::AbstractArray; tol::Real=1e-6, maxits=Inf, maxtime=Inf, verbose=false, pcgIts=Int[], params...)
Uses a lapSolver
to solve systems of linear equations in sddm matrices.
Laplacians.sortSet
— Method.Given a set of integers, set
between 1 and n, return a sorted version of them
Laplacians.star_graph_ijv
— Method.ijv = star_graph_ijv(n::Int64)
Laplacians.testZeroDiag
— Method.testZeroDiag(a)
Returns true if a
has zero diagonal, false otherwise
Laplacians.treeDepthDFS
— Method.Compute the vector of depths in a tree that is in DFS order, with the root at the first position, and the leaves at the end
Laplacians.uniformWeight_ver
— Method.wted = uniformWeight(unwted)
Put a uniform [0,1] weight on every edge. This is an example of how to use mapweight.
Laplacians.wrapCapture
— Method.f = wrapCapture(solver::Function, mats, rhss; tol::Real=1e-6, maxits=Inf, maxtime=Inf, verbose=false, pcgIts=Int[], params...)
This wraps a solver so that we can capture all the matrices that it solves and all the right-hand-sides. Those are pushed into the arrays mats
and rhss
. For example
julia> mats = []
julia> rhss = []
julia> solver = wrapCapture(approxchol_lap, mats, rhss)
julia> a = chimera(10)
julia> f = solver(a);
julia> size(mats[1])
(10,10)
julia> b = randn(10)
julia> x = f(b);
julia> rhss
1-element Array{Any,1}:
[0.404962,-0.827718,0.704616,-0.403223,0.204891,-0.505589,0.907015,1.90266,-0.438115,0.0464351]
Laplacians.wrapCaptureRhs
— Method.f = wrapCaptureRhs(sola::Function, rhss; tol::Real=1e-6, maxits=Inf, maxtime=Inf, verbose=false, pcgIts=Int[], params...)
Captures all the right-hand-sides that are passed to the solver sola
. It pushes them into an array called rhhs. For example
julia> rhss = []
julia> a = wted_chimera(100)
julia> sola = approxchol_lap(a)
julia> wrappedSolver = wrapCaptureRhs(sola,rhss)
julia> b = randn(100)
julia> x = wrappedSolver(b,verbose=true)
PCG BLAS stopped after: 0.0 seconds and 11 iterations with relative error 3.160275810360986e-7.
julia> length(rhss[1])
100
Laplacians.wrapInterface
— Method.solveA = wrapInterface(solver::Function, A::AbstractMatrix; tol, maxits, maxtime, verbose, pcgIts=Int[],params...)
solverConstructor = wrapInterface(A::AbstractMatrix; tol, maxits, maxtime, verbose, pcgIts=Int[],params...)
Returns a function that discards tol
, maxits
, maxtime
and verbose
, sets pcgIts
to 0 (because it might not be using pcg), and passes whatever params
are left to the solver.
Examples
julia> a = randn(5,5);
julia> a = a * a';
julia> solvea = wrapInterface(X->cholesky(X,Val(true)), a, maxits=100, verbose=true);
julia> b = randn(5,1);
julia> norm(a*solvea(b, verbose=false)-b)
1.575705319704736e-14
julia> f = wrapInterface(X->cholesky(X,Val(true)))
julia> solvea = f(a, maxits=1000, maxtime = 1)
julia> norm(a*solvea(b, verbose=false, maxtime = 10)-b)
1.575705319704736e-14