~ubuntu-branches/ubuntu/wily/julia/wily
» Revision
16
: .pc/suitesparse-3.4.patch/base/linalg/cholmod.jl
« back to all changes in this revision
Viewing changes to .pc/suitesparse-3.4.patch/base/linalg/cholmod.jl
- browse files at revision 16
- compare with another revision
- download diff
- download tarball
- view history from revision 16
- Committer: Package Import Robot
- Author(s): Sébastien Villemot
- Date: 2013-11-17 19:32:52 UTC
- mfrom: (1.1.12)
- Revision ID: package-import@ubuntu.com-20131117193252-tkrpclguqqebqa35
Tags: 0.2.0+dfsg-3
testsuite-i386.patch: loosen the numerical precision for yet another test.
- files added:
- .pc/sysconfdir-install.patch
- files removed:
- .pc/fhs.patch
- .pc/fhs.patch/base
- .pc/suitesparse-3.4.patch
- .pc/suitesparse-3.4.patch/base
- .pc/suitesparse-3.4.patch/base/linalg
- files modified:
- .build_info
added
removed
.pc/suitesparse-3.4.patch/base/linalg/cholmod.jl
1
module CHOLMOD
2
3
export # types
4
CholmodDense,
5
CholmodFactor,
6
CholmodSparse,
7
CholmodTriplet,
8
9
CholmodSparse!, # destructive constructors
10
CholmodDense!,
11
12
etree
13
14
using Base.LinAlg.UMFPACK # for decrement, increment, etc.
15
16
import Base: (*), convert, copy, ctranspose, eltype, findnz, getindex, hcat,
17
isvalid, nnz, show, size, sort!, transpose, vcat
18
19
import ..LinAlg: (\), A_mul_Bc, A_mul_Bt, Ac_ldiv_B, Ac_mul_B, At_ldiv_B, At_mul_B,
20
cholfact, cholfact!, copy, dense, det, diag,
21
full, logdet, norm, scale, scale!, solve, sparse
22
23
include("cholmod_h.jl")
24
25
const chm_ver = Array(Cint, 3)
26
if dlsym(dlopen("libcholmod"), :cholmod_version) != C_NULL
27
ccall((:cholmod_version, "libcholmod.so.1.7.1"), Cint, (Ptr{Cint},), chm_ver)
28
else
29
ccall((:jl_cholmod_version, :libsuitesparse_wrapper), Cint, (Ptr{Cint},), chm_ver)
30
end
31
const chm_com_sz = ccall((:jl_cholmod_common_size,:libsuitesparse_wrapper),Int,())
32
const chm_com = fill(0xff, chm_com_sz)
33
const chm_l_com = fill(0xff, chm_com_sz)
34
## chm_com and chm_l_com must be initialized at runtime because they contain pointers
35
## to functions in libc.so, whose addresses can change
36
function cmn(::Type{Int32})
37
if isnan(reinterpret(Float64,chm_com[1:8])[1])
38
status = ccall((:cholmod_start, "libcholmod.so.1.7.1"), Cint, (Ptr{Uint8},), chm_com)
39
if status != CHOLMOD_TRUE throw(CholmodException) end
40
end
41
chm_com
42
end
43
function cmn(::Type{Int64})
44
if isnan(reinterpret(Float64,chm_l_com[1:8])[1])
45
status = ccall((:cholmod_l_start, "libcholmod.so.1.7.1"), Cint, (Ptr{Uint8},), chm_l_com)
46
if status != CHOLMOD_TRUE throw(CholmodException) end
47
end
48
chm_l_com
49
end
50
51
typealias CHMITypes Union(Int32,Int64)
52
typealias CHMVTypes Union(Complex64, Complex128, Float32, Float64)
53
54
type CholmodException <: Exception end
55
56
## macro to generate the name of the C function according to the integer type
57
macro chm_nm(nm,typ) string("cholmod_", eval(typ) == Int64 ? "l_" : "", nm) end
58
59
### A way of examining some of the fields in chm_com
60
### Probably better to make this a Dict{ASCIIString,Tuple} and
61
### save the offsets and the lengths and the types. Then the names can be checked.
62
type ChmCommon
63
dbound::Float64
64
maxrank::Int
65
supernodal_switch::Float64
66
supernodal::Int32
67
final_asis::Int32
68
final_super::Int32
69
final_ll::Int32
70
final_pack::Int32
71
final_monotonic::Int32
72
final_resymbol::Int32
73
prefer_zomplex::Int32 # should always be false
74
prefer_upper::Int32
75
print::Int32 # print level. Default: 3
76
precise::Int32 # print 16 digits, otherwise 5
77
nmethods::Int32 # number of ordering methods
78
selected::Int32
79
postorder::Int32
80
itype::Int32
81
dtype::Int32
82
end
83
84
### These offsets should be reconfigured to be less error-prone in matches
85
const chm_com_offsets = Array(Int, length(ChmCommon.types))
86
ccall((:jl_cholmod_common_offsets, :libsuitesparse_wrapper),
87
Void, (Ptr{Uint8},), chm_com_offsets)
88
const chm_final_ll_inds = (1:4) + chm_com_offsets[7]
89
const chm_prt_inds = (1:4) + chm_com_offsets[13]
90
const chm_ityp_inds = (1:4) + chm_com_offsets[18]
91
92
### there must be an easier way but at least this works.
93
function ChmCommon(aa::Array{Uint8,1})
94
typs = ChmCommon.types
95
sz = map(sizeof, typs)
96
args = map(i->reinterpret(typs[i], aa[chm_com_offsets[i] + (1:sz[i])])[1], 1:length(sz))
97
eval(Expr(:call, unshift!(args, :ChmCommon), Any))
98
end
99
100
function set_chm_prt_lev(cm::Array{Uint8}, lev::Integer) # can probably be removed
101
cm[(1:4) + chm_com_offsets[13]] = reinterpret(Uint8, [int32(lev)])
102
end
103
104
## cholmod_dense pointers passed to or returned from C functions are of Julia type
105
## Ptr{c_CholmodDense}. The CholmodDense type contains a c_CholmodDense object and other
106
## fields then ensure the memory pointed to is freed when it should be and not before.
107
type c_CholmodDense{T<:CHMVTypes}
108
m::Int
109
n::Int
110
nzmax::Int
111
lda::Int
112
xpt::Ptr{T}
113
zpt::Ptr{Void}
114
xtype::Cint
115
dtype::Cint
116
end
117
118
type CholmodDense{T<:CHMVTypes}
119
c::c_CholmodDense
120
mat::Matrix{T}
121
end
122
123
if (1000chm_ver[1]+chm_ver[2]) >= 2001 # CHOLMOD version 2.1.0 or later
124
type c_CholmodFactor{Tv<:CHMVTypes,Ti<:CHMITypes}
125
n::Int
126
minor::Int
127
Perm::Ptr{Ti} # this pointer was added in verison 2.1.0
128
ColCount::Ptr{Ti}
129
IPerm::Ptr{Ti}
130
nzmax::Int
131
p::Ptr{Ti}
132
i::Ptr{Ti}
133
x::Ptr{Tv}
134
z::Ptr{Void}
135
nz::Ptr{Ti}
136
next::Ptr{Ti}
137
prev::Ptr{Ti}
138
nsuper::Int
139
ssize::Int
140
xsize::Int
141
maxcsize::Int
142
maxesize::Int
143
super::Ptr{Ti}
144
pi::Ptr{Ti}
145
px::Ptr{Ti}
146
s::Ptr{Ti}
147
ordering::Cint
148
is_ll::Cint
149
is_super::Cint
150
is_monotonic::Cint
151
itype::Cint
152
xtype::Cint
153
dtype::Cint
154
end
155
156
type CholmodFactor{Tv<:CHMVTypes,Ti<:CHMITypes}
157
c::c_CholmodFactor{Tv,Ti}
158
Perm::Vector{Ti}
159
ColCount::Vector{Ti}
160
IPerm::Vector{Ti}
161
p::Vector{Ti}
162
i::Vector{Ti}
163
x::Vector{Tv}
164
nz::Vector{Ti}
165
next::Vector{Ti}
166
prev::Vector{Ti}
167
super::Vector{Ti}
168
pi::Vector{Ti}
169
px::Vector{Ti}
170
s::Vector{Ti}
171
end
172
173
function CholmodFactor{Tv<:CHMVTypes,Ti<:CHMITypes}(cp::Ptr{c_CholmodFactor{Tv,Ti}})
174
cfp = unsafe_load(cp)
175
Perm = pointer_to_array(cfp.Perm, (cfp.n,), true)
176
ColCount = pointer_to_array(cfp.ColCount, (cfp.n,), true)
177
IPerm = pointer_to_array(cfp.IPerm, (cfp.IPerm == C_NULL ? 0 : cfp.n + 1,), true)
178
p = pointer_to_array(cfp.p, (cfp.p == C_NULL ? 0 : cfp.n + 1,), true)
179
i = pointer_to_array(cfp.i, (cfp.i == C_NULL ? 0 : cfp.nzmax,), true)
180
x = pointer_to_array(cfp.x, (cfp.x == C_NULL ? 0 : max(cfp.nzmax,cfp.xsize),), true)
181
nz = pointer_to_array(cfp.nz, (cfp.nz == C_NULL ? 0 : cfp.n,), true)
182
next = pointer_to_array(cfp.next, (cfp.next == C_NULL ? 0 : cfp.n + 2,), true)
183
prev = pointer_to_array(cfp.prev, (cfp.prev == C_NULL ? 0 : cfp.n + 2,), true)
184
super = pointer_to_array(cfp.super, (cfp.super == C_NULL ? 0 : cfp.nsuper + 1,), true)
185
pi = pointer_to_array(cfp.pi, (cfp.pi == C_NULL ? 0 : cfp.nsuper + 1,), true)
186
px = pointer_to_array(cfp.px, (cfp.px == C_NULL ? 0 : cfp.nsuper + 1,), true)
187
s = pointer_to_array(cfp.s, (cfp.s == C_NULL ? 0 : cfp.ssize + 1,), true)
188
cf = CholmodFactor{Tv,Ti}(cfp, Perm, ColCount, IPerm, p, i, x, nz, next, prev,
189
super, pi, px, s)
190
c_free(cp)
191
cf
192
end
193
else
194
type c_CholmodFactor{Tv<:CHMVTypes,Ti<:CHMITypes}
195
n::Int
196
minor::Int
197
ColCount::Ptr{Ti}
198
IPerm::Ptr{Ti}
199
nzmax::Int
200
p::Ptr{Ti}
201
i::Ptr{Ti}
202
x::Ptr{Tv}
203
z::Ptr{Void}
204
nz::Ptr{Ti}
205
next::Ptr{Ti}
206
prev::Ptr{Ti}
207
nsuper::Int
208
ssize::Int
209
xsize::Int
210
maxcsize::Int
211
maxesize::Int
212
super::Ptr{Ti}
213
pi::Ptr{Ti}
214
px::Ptr{Ti}
215
s::Ptr{Ti}
216
ordering::Cint
217
is_ll::Cint
218
is_super::Cint
219
is_monotonic::Cint
220
itype::Cint
221
xtype::Cint
222
dtype::Cint
223
end
224
225
type CholmodFactor{Tv<:CHMVTypes,Ti<:CHMITypes}
226
c::c_CholmodFactor{Tv,Ti}
227
Perm::Vector{Ti}
228
ColCount::Vector{Ti}
229
p::Vector{Ti}
230
i::Vector{Ti}
231
x::Vector{Tv}
232
nz::Vector{Ti}
233
next::Vector{Ti}
234
prev::Vector{Ti}
235
super::Vector{Ti}
236
pi::Vector{Ti}
237
px::Vector{Ti}
238
s::Vector{Ti}
239
end
240
241
function CholmodFactor{Tv<:CHMVTypes,Ti<:CHMITypes}(cp::Ptr{c_CholmodFactor{Tv,Ti}})
242
cfp = unsafe_load(cp)
243
Perm = pointer_to_array(cfp.Perm, (cfp.n,), true)
244
ColCount = pointer_to_array(cfp.ColCount, (cfp.n,), true)
245
p = pointer_to_array(cfp.p, (cfp.p == C_NULL ? 0 : cfp.n + 1,), true)
246
i = pointer_to_array(cfp.i, (cfp.i == C_NULL ? 0 : cfp.nzmax,), true)
247
x = pointer_to_array(cfp.x, (cfp.x == C_NULL ? 0 : max(cfp.nzmax,cfp.xsize),), true)
248
nz = pointer_to_array(cfp.nz, (cfp.nz == C_NULL ? 0 : cfp.n,), true)
249
next = pointer_to_array(cfp.next, (cfp.next == C_NULL ? 0 : cfp.n + 2,), true)
250
prev = pointer_to_array(cfp.prev, (cfp.prev == C_NULL ? 0 : cfp.n + 2,), true)
251
super = pointer_to_array(cfp.super, (cfp.super == C_NULL ? 0 : cfp.nsuper + 1,), true)
252
pi = pointer_to_array(cfp.pi, (cfp.pi == C_NULL ? 0 : cfp.nsuper + 1,), true)
253
px = pointer_to_array(cfp.px, (cfp.px == C_NULL ? 0 : cfp.nsuper + 1,), true)
254
s = pointer_to_array(cfp.s, (cfp.s == C_NULL ? 0 : cfp.ssize + 1,), true)
255
cf = CholmodFactor{Tv,Ti}(cfp, Perm, ColCount, IPerm, p, i, x, nz, next, prev,
256
super, pi, px, s)
257
c_free(cp)
258
cf
259
end
260
end
261
262
type c_CholmodSparse{Tv<:CHMVTypes,Ti<:CHMITypes}
263
m::Int
264
n::Int
265
nzmax::Int
266
ppt::Ptr{Ti}
267
ipt::Ptr{Ti}
268
nzpt::Ptr{Void}
269
xpt::Ptr{Tv}
270
zpt::Ptr{Void}
271
stype::Cint
272
itype::Cint
273
xtype::Cint
274
dtype::Cint
275
sorted::Cint
276
packed::Cint
277
end
278
279
type CholmodSparse{Tv<:CHMVTypes,Ti<:CHMITypes}
280
c::c_CholmodSparse{Tv,Ti}
281
colptr0::Vector{Ti}
282
rowval0::Vector{Ti}
283
nzval::Vector{Tv}
284
end
285
286
type c_CholmodTriplet{Tv<:CHMVTypes,Ti<:CHMITypes}
287
m::Int
288
n::Int
289
nzmax::Int
290
nnz::Int
291
i::Ptr{Ti}
292
j::Ptr{Ti}
293
x::Ptr{Tv}
294
z::Ptr{Void}
295
stype:Cint
296
itype::Cint
297
xtype::Cint
298
dtype::Cint
299
end
300
301
type CholmodTriplet{Tv<:CHMVTypes,Ti<:CHMITypes}
302
c::c_CholmodTriplet{Tv,Ti}
303
i::Vector{Ti}
304
j::Vector{Ti}
305
x::Vector{Tv}
306
end
307
308
eltype{T<:CHMVTypes}(A::CholmodDense{T}) = T
309
eltype{T<:CHMVTypes}(A::CholmodFactor{T}) = T
310
eltype{T<:CHMVTypes}(A::CholmodSparse{T}) = T
311
eltype{T<:CHMVTypes}(A::CholmodTriplet{T}) = T
312
313
## The CholmodDense! constructor does not copy the contents, which is generally what you
314
## want as most uses of CholmodDense objects are read-only.
315
function CholmodDense!{T<:CHMVTypes}(aa::VecOrMat{T}) # uses the memory from Julia
316
m = size(aa,1); n = size(aa,2)
317
CholmodDense(c_CholmodDense{T}(m, n, m*n, stride(aa,2), convert(Ptr{T}, aa),
318
C_NULL, xtyp(T), dtyp(T)),
319
length(size(aa)) == 2 ? aa : reshape(aa, (m,n)))
320
end
321
322
## The CholmodDense constructor copies the contents
323
function CholmodDense{T<:CHMVTypes}(aa::VecOrMat{T})
324
m = size(aa,1); n = size(aa,2)
325
acp = length(size(aa)) == 2 ? copy(aa) : reshape(copy(aa), (m,n))
326
CholmodDense(c_CholmodDense{T}(m, n, m*n, stride(aa,2), convert(Ptr{T}, acp),
327
C_NULL, xtyp(T), dtyp(T)), acp)
328
end
329
330
function CholmodDense{T<:CHMVTypes}(c::Ptr{c_CholmodDense{T}})
331
cp = unsafe_load(c)
332
if cp.lda != cp.m || cp.nzmax != cp.m * cp.n
333
error("overallocated cholmod_dense returned object of size $(cp.m) by $(cp.n) with leading dim $(cp.lda) and nzmax $(cp.nzmax)")
334
end
335
## the true in the call to pointer_to_array means Julia will free the memory
336
val = CholmodDense(cp, pointer_to_array(cp.xpt, (cp.m,cp.n), true))
337
c_free(c)
338
val
339
end
340
CholmodDense!{T<:CHMVTypes}(c::Ptr{c_CholmodDense{T}}) = CholmodDense(c) # no distinction
341
342
function isvalid{T<:CHMVTypes}(cd::CholmodDense{T})
343
bool(ccall((:cholmod_check_dense, "libcholmod.so.1.7.1"), Cint,
344
(Ptr{c_CholmodDense{T}}, Ptr{Uint8}), &cd.c, cmn(Int32)))
345
end
346
347
function chm_eye{T<:Union(Float64,Complex128)}(m::Integer, n::Integer, t::T)
348
CholmodDense(ccall((:cholmod_eye, "libcholmod.so.1.7.1"), Ptr{c_CholmodDense{T}},
349
(Int, Int, Cint, Ptr{Uint8}),
350
m, n,xtyp(T),cmn(Int32)))
351
end
352
chm_eye(m::Integer, n::Integer) = chm_eye(m, n, 1.)
353
chm_eye(n::Integer) = chm_eye(n, n, 1.)
354
355
function chm_ones{T<:Union(Float64,Complex128)}(m::Integer, n::Integer, t::T)
356
CholmodDense(ccall((:cholmod_ones, "libcholmod.so.1.7.1"), Ptr{c_CholmodDense{T}},
357
(Int, Int, Cint, Ptr{Uint8}),
358
m, n, xtyp(T), cmn(Int32)))
359
end
360
chm_ones(m::Integer, n::Integer) = chm_ones(m, n, 1.)
361
362
function chm_zeros{T<:Union(Float64,Complex128)}(m::Integer, n::Integer, t::T)
363
CholmodDense(ccall((:cholmod_zeros, "libcholmod.so.1.7.1"), Ptr{c_CholmodDense{T}},
364
(Int, Int, Cint, Ptr{Uint8}),
365
m, n, xtyp(T), cmn(Int32)))
366
end
367
chm_zeros(m::Integer, n::Integer) = chm_zeros(m, n, 1.)
368
369
function chm_print{T<:CHMVTypes}(cd::CholmodDense{T}, lev::Integer, nm::ASCIIString)
370
cm = cmn(Int32)
371
orig = cm[chm_prt_inds]
372
cm[chm_prt_inds] = reinterpret(Uint8, [int32(lev)])
373
status = ccall((:cholmod_print_dense, "libcholmod.so.1.7.1"), Cint,
374
(Ptr{c_CholmodDense{T}}, Ptr{Uint8}, Ptr{Uint8}),
375
&cd.c, nm, cm)
376
cm[chm_prt_inds] = orig
377
if status != CHOLMOD_TRUE throw(CholmodException) end
378
end
379
chm_print(cd::CholmodDense, lev::Integer) = chm_print(cd, lev, "")
380
chm_print(cd::CholmodDense) = chm_print(cd, int32(4), "")
381
show(io::IO,cd::CholmodDense) = chm_print(cd, int32(4), "")
382
383
function copy{Tv<:CHMVTypes}(B::CholmodDense{Tv})
384
CholmodDense(ccall((:cholmod_copy_dense,"libcholmod.so.1.7.1"), Ptr{c_CholmodDense{Tv}},
385
(Ptr{c_CholmodDense{Tv}},Ptr{Uint8}), &B.c, cmn(Int32)))
386
end
387
388
function norm{Tv<:CHMVTypes}(D::CholmodDense{Tv},p::Number)
389
ccall((:cholmod_norm_dense, "libcholmod.so.1.7.1"), Float64,
390
(Ptr{c_CholmodDense{Tv}}, Cint, Ptr{Uint8}),
391
&D.c, p == 1 ? 1 :(p == Inf ? 1 : error("p must be 1 or Inf")),cmn(Int32))
392
end
393
norm{Tv<:CHMVTypes}(D::CholmodDense{Tv}) = norm(D,1)
394
395
function CholmodSparse!{Tv<:CHMVTypes,Ti<:CHMITypes}(colpt::Vector{Ti},
396
rowval::Vector{Ti},
397
nzval::Vector{Tv},
398
m::Integer,
399
n::Integer,
400
stype::Signed)
401
bb = colpt[1]
402
if bb != 0 && bb != 1 error("colpt[1] is $bb, must be 0 or 1") end
403
if any(diff(colpt) .< 0) error("elements of colpt must be non-decreasing") end
404
if length(colpt) != n + 1 error("length(colptr) = $(length(colpt)), should be $(n+1)") end
405
if bool(bb) # one-based
406
decrement!(colpt)
407
decrement!(rowval)
408
end
409
nz = colpt[end]
410
if length(rowval) != nz || length(nzval) != nz
411
error("length(rowval) = $(length(rowval)) and length(nzval) = $(length(nzval)) should be $nz")
412
end
413
if any(rowval .< 0) || any(rowval .>= m)
414
error("all elements of rowval must be in the range [0,$(m-1)]")
415
end
416
it = ityp(Ti)
417
sort!(CholmodSparse(c_CholmodSparse{Tv,Ti}(m,n,int(nz),convert(Ptr{Ti},colpt),
418
convert(Ptr{Ti},rowval), C_NULL,
419
convert(Ptr{Tv},nzval), C_NULL,
420
int32(stype), ityp(Ti),
421
xtyp(Tv),dtyp(Tv),
422
CHOLMOD_FALSE,CHOLMOD_TRUE),
423
colpt,rowval,nzval))
424
end
425
function CholmodSparse{Tv<:CHMVTypes,Ti<:CHMITypes}(colpt::Vector{Ti},
426
rowval::Vector{Ti},
427
nzval::Vector{Tv},
428
m::Integer,
429
n::Integer,
430
stype::Signed)
431
CholmodSparse!(copy(colpt),copy(rowval),copy(nzval),m,n,stype)
432
end
433
function CholmodSparse!{Tv<:CHMVTypes,Ti<:CHMITypes}(A::SparseMatrixCSC{Tv,Ti}, stype::Signed)
434
CholmodSparse!(A.colptr,A.rowval,A.nzval,size(A,1),size(A,2),stype)
435
end
436
function CholmodSparse{Tv<:CHMVTypes,Ti<:CHMITypes}(A::SparseMatrixCSC{Tv,Ti}, stype::Signed)
437
CholmodSparse!(copy(A.colptr),copy(A.rowval),copy(A.nzval),size(A,1),size(A,2),stype)
438
end
439
function CholmodSparse(A::SparseMatrixCSC)
440
stype = ishermitian(A) ? 1 : 0
441
CholmodSparse(stype > 0 ? triu(A) : A, stype)
442
end
443
function CholmodSparse{Tv<:CHMVTypes,Ti<:CHMITypes}(cp::Ptr{c_CholmodSparse{Tv,Ti}})
444
csp = unsafe_load(cp)
445
colptr0 = pointer_to_array(csp.ppt, (csp.n + 1,), true)
446
nnz = int(colptr0[end])
447
cms = CholmodSparse{Tv,Ti}(csp, colptr0,
448
pointer_to_array(csp.ipt, (nnz,), true),
449
pointer_to_array(csp.xpt, (nnz,), true))
450
c_free(cp)
451
cms
452
end
453
CholmodSparse!{Tv<:CHMVTypes,Ti<:CHMITypes}(cp::Ptr{c_CholmodSparse{Tv,Ti}}) = CholmodSparse(cp)
454
CholmodSparse{Tv<:CHMVTypes}(D::CholmodDense{Tv}) = CholmodSparse(D,1) # default Ti is Int
455
456
function CholmodTriplet{Tv<:CHMVTypes,Ti<:CHMITypes}(tp::Ptr{c_CholmodTriplet{Tv,Ti}})
457
ctp = unsafe_load(tp)
458
i = pointer_to_array(ctp.i, (ctp.nnz,), true)
459
j = pointer_to_array(ctp.j, (ctp.nnz,), true)
460
x = pointer_to_array(ctp.x, (ctp.x == C_NULL ? 0 : ctp.nnz), true)
461
ct = CholmodTriplet{Tv,Ti}(ctp, i, j, x)
462
c_free(tp)
463
ct
464
end
465
466
function chm_rdsp(fnm::String)
467
fd = ccall(:fopen, Ptr{Void}, (Ptr{Uint8},Ptr{Uint8}), fnm, "r")
468
res = ccall((:cholmod_read_sparse,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Float64,Cint}},
469
(Ptr{Void},Ptr{Uint8}),fd,cmn(Cint))
470
ccall(:fclose, Cint, (Ptr{Void},), fd) # should do this in try/finally/end
471
CholmodSparse(res)
472
end
473
474
for Ti in (:Int32,:Int64)
475
@eval begin
476
function (*){Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti},
477
B::CholmodSparse{Tv,$Ti})
478
CholmodSparse(ccall((@chm_nm "ssmult" $Ti
479
, "libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
480
(Ptr{c_CholmodSparse{Tv,$Ti}},Ptr{c_CholmodSparse{Tv,$Ti}},
481
Cint,Cint,Cint,Ptr{Uint8}), &A.c,&B.c,0,true,true,cmn($Ti)))
482
end
483
function A_mul_Bc{Tv<:Union(Float32,Float64)}(A::CholmodSparse{Tv,$Ti},
484
B::CholmodSparse{Tv,$Ti})
485
cm = cmn($Ti)
486
aa = Array(Ptr{c_CholmodSparse{Tv,$Ti}}, 2)
487
if !is(A,B)
488
aa[1] = ccall((@chm_nm "transpose" $Ti
489
,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
490
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}),
491
&B.c,cm)
492
aa[2] = ccall((@chm_nm "ssmult" $Ti
493
,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
494
(Ptr{c_CholmodSparse{Tv,$Ti}},
495
Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}),
496
&A.c, aa[1], cmn($Ti))
497
status = ccall((@chm_nm "free_sparse" $Ti
498
,"libcholmod.so.1.7.1"), Cint,
499
(Ptr{Ptr{c_CholmodSparse{Tv,$Ti}}}, Ptr{Uint8}), aa, cm)
500
if status != CHOLMOD_TRUE throw(CholmodException) end
501
return CholmodSparse(aa[2])
502
end
503
## The A*A' case is handled by cholmod_aat. Strangely the matrix returned by
504
## cholmod_aat is not marked as symmetric. The code following the call to
505
## cholmod_aat is to create the symmetric-storage version of the result then
506
## transpose it to provide sorted columns. The result is stored in the upper
507
## triangle
508
aa[1] = ccall((@chm_nm "aat" $Ti
509
, "libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
510
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Void}, Int, Cint, Ptr{Uint8}),
511
&A.c, C_NULL, 0, 1, cm)
512
## Create the lower triangle unsorted
513
aa[2] = ccall((@chm_nm "copy" $Ti
514
, "libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
515
(Ptr{c_CholmodSparse{Tv,$Ti}}, Cint, Cint, Ptr{Uint8}),
516
aa[1], -1, 1, cm)
517
status = ccall((@chm_nm "free_sparse" $Ti
518
, "libcholmod.so.1.7.1"), Cint,
519
(Ptr{Ptr{c_CholmodSparse{Tv,$Ti}}}, Ptr{Uint8}), aa, cm)
520
if status != CHOLMOD_TRUE throw(CholmodException) end
521
aa[1] = aa[2]
522
r = unsafe_load(aa[1])
523
## Now transpose the lower triangle to the upper triangle to do the sorting
524
rpt = ccall((@chm_nm "allocate_sparse" $Ti
525
,"libcholmod.so.1.7.1"),Ptr{c_CholmodSparse{Tv,$Ti}},
526
(Csize_t,Csize_t,Csize_t,Cint,Cint,Cint,Cint,Ptr{Cuchar}),
527
r.m,r.n,r.nzmax,r.sorted,r.packed,-r.stype,r.xtype,cm)
528
status = ccall((@chm_nm "transpose_sym" $Ti
529
,"libcholmod.so.1.7.1"),Cint,
530
(Ptr{c_CholmodSparse{Tv,$Ti}}, Cint, Ptr{$Ti},
531
Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}),
532
aa[1],1,C_NULL,rpt,cm)
533
if status != CHOLMOD_TRUE throw(CholmodException) end
534
status = ccall((@chm_nm "free_sparse" $Ti
535
, "libcholmod.so.1.7.1"), Cint,
536
(Ptr{Ptr{c_CholmodSparse{Tv,$Ti}}}, Ptr{Uint8}), aa, cm)
537
if status != CHOLMOD_TRUE throw(CholmodException) end
538
CholmodSparse(rpt)
539
end
540
function Ac_mul_B{Tv<:Union(Float32,Float64)}(A::CholmodSparse{Tv,$Ti},
541
B::CholmodSparse{Tv,$Ti})
542
cm = cmn($Ti)
543
aa = Array(Ptr{c_CholmodSparse{Tv,$Ti}}, 2)
544
aa[1] = ccall((@chm_nm "transpose" $Ti
545
,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
546
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}),
547
&B.c,cm)
548
if is(A,B)
549
Ac = CholmodSparse(aa[1])
550
return A_mul_Bc(Ac,Ac)
551
end
552
aa[2] = ccall((@chm_nm "ssmult" $Ti
553
,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
554
(Ptr{c_CholmodSparse{Tv,$Ti}},
555
Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}),
556
aa[1],&B.c,cm)
557
status = ccall((@chm_nm "free_sparse" $Ti
558
, "libcholmod.so.1.7.1"), Cint,
559
(Ptr{Ptr{c_CholmodSparse{Tv,$Ti}}}, Ptr{Uint8}), aa, cm)
560
if status != CHOLMOD_TRUE throw(CholmodException) end
561
CholmodSparse(aa[2])
562
end
563
function CholmodDense{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti})
564
CholmodDense(ccall((@chm_nm "sparse_to_dense" $Ti
565
,"libcholmod.so.1.7.1"), Ptr{c_CholmodDense{Tv,$Ti}},
566
(Ptr{c_CholmodSparse{Tv,Ti}},Ptr{Uint8}),
567
&A.c,cmn($Ti)))
568
end
569
function CholmodSparse{Tv<:CHMVTypes}(D::CholmodDense{Tv},i::$Ti)
570
CholmodSparse(ccall((@chm_nm "dense_to_sparse" $Ti
571
,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,Ti}},
572
(Ptr{c_CholmodDense{Tv,Ti}},Ptr{Uint8}),
573
&D.c,cmn($Ti)))
574
end
575
function CholmodSparse{Tv<:CHMVTypes}(L::CholmodFactor{Tv,$Ti})
576
if bool(L.c.is_ll)
577
return CholmodSparse(ccall((@chm_nm "factor_to_sparse" $Ti
578
,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
579
(Ptr{c_CholmodFactor{Tv,$Ti}},Ptr{Uint8}),
580
&L.c,cmn($Ti)))
581
end
582
cm = cmn($Ti)
583
Lcll = ccall((@chm_nm "copy_factor" $Ti
584
,"libcholmod.so.1.7.1"), Ptr{c_CholmodFactor{Tv,$Ti}},
585
(Ptr{c_CholmodFactor{Tv,$Ti}},Ptr{Uint8}),
586
&L.c,cm)
587
status = ccall((@chm_nm "change_factor" $Ti
588
,"libcholmod.so.1.7.1"), Cint,
589
(Cint,Cint,Cint,Cint,Cint,Ptr{c_CholmodFactor{Tv,$Ti}},Ptr{Uint8}),
590
L.c.xtype,true,L.c.is_super,true,true,Lcll,cm)
591
if status != CHOLMOD_TRUE throw(CholmodException) end
592
val = CholmodSparse(ccall((@chm_nm "factor_to_sparse" $Ti
593
,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
594
(Ptr{c_CholmodFactor{Tv,$Ti}},Ptr{Uint8}),
595
Lcll,cmn($Ti)))
596
status = ccall((@chm_nm "free_factor" $Ti
597
,"libcholmod.so.1.7.1"), Cint,
598
(Ptr{Ptr{c_CholmodFactor{Tv,$Ti}}},Ptr{Uint8}),
599
[Lcll],cm)
600
if status != CHOLMOD_TRUE throw(CholmodException) end
601
val
602
end
603
function CholmodSparse{Tv<:CHMVTypes,Ti<:$Ti}(T::CholmodTriplet{Tv,Ti})
604
CholmodSparse(ccall((@chm_nm "triplet_to_sparse" $Ti
605
,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,Ti}},
606
(Ptr{c_CholmodTriplet{Tv,Ti}},Ptr{Uint8}),
607
&T.c,cmn($Ti)))
608
end
609
function CholmodTriplet{Tv<:CHMVTypes,Ti<:$Ti}(A::CholmodSparse{Tv,Ti})
610
CholmodTriplet(ccall((@chm_nm "sparse_to_triplet" $Ti
611
,"libcholmod.so.1.7.1"), Ptr{c_CholmodTriplet{Tv,Ti}},
612
(Ptr{c_CholmodSparse{Tv,Ti}},Ptr{Uint8}),
613
&A.c,cmn($Ti)))
614
end
615
function isvalid{Tv<:CHMVTypes}(L::CholmodFactor{Tv,$Ti})
616
bool(ccall((@chm_nm "check_factor" $Ti
617
,"libcholmod.so.1.7.1"), Cint,
618
(Ptr{c_CholmodFactor{Tv,$Ti}}, Ptr{Uint8}),
619
&L.c, cmn($Ti)))
620
end
621
function isvalid{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti})
622
bool(ccall((@chm_nm "check_sparse" $Ti
623
,"libcholmod.so.1.7.1"), Cint,
624
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}),
625
&A.c, cmn($Ti)))
626
end
627
function isvalid{Tv<:CHMVTypes}(T::CholmodTriplet{Tv,$Ti})
628
bool(ccall((@chm_nm "check_triplet" $Ti
629
,"libcholmod.so.1.7.1"), Cint,
630
(Ptr{c_CholmodTriplet{Tv,$Ti}}, Ptr{Uint8}),
631
&T.c, cmn($Ti)))
632
end
633
function cholfact{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti}, ll::Bool)
634
cm = cmn($Ti)
635
## may need to change final_asis as well as final_ll
636
if ll cm[chm_final_ll_inds] = reinterpret(Uint8, [one(Cint)]) end
637
Lpt = ccall((@chm_nm "analyze" $Ti
638
,"libcholmod.so.1.7.1"), Ptr{c_CholmodFactor{Tv,$Ti}},
639
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}), &A.c, cm)
640
status = ccall((@chm_nm "factorize" $Ti
641
,"libcholmod.so.1.7.1"), Cint,
642
(Ptr{c_CholmodSparse{Tv,$Ti}},
643
Ptr{c_CholmodFactor{Tv,$Ti}}, Ptr{Uint8}),
644
&A.c, Lpt, cm)
645
if status != CHOLMOD_TRUE throw(CholmodException) end
646
CholmodFactor(Lpt)
647
end
648
function cholfact{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti},beta::Tv,ll::Bool)
649
cm = cmn($Ti)
650
## may need to change final_asis as well as final_ll
651
if ll cm[chm_final_ll_inds] = reinterpret(Uint8, [one(Cint)]) end
652
Lpt = ccall((@chm_nm "analyze" $Ti
653
,"libcholmod.so.1.7.1"), Ptr{c_CholmodFactor{Tv,$Ti}},
654
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}), &A.c, cm)
655
status = ccall((@chm_nm "factorize_p" $Ti
656
,"libcholmod.so.1.7.1"), Cint,
657
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Tv}, Ptr{Cint}, Csize_t,
658
Ptr{c_CholmodFactor{Tv,$Ti}}, Ptr{Uint8}),
659
&A.c, &beta, C_NULL, zero(Csize_t), Lpt, cmn($Ti))
660
if status != CHOLMOD_TRUE throw(CholmodException) end
661
CholmodFactor(Lpt)
662
end
663
function chm_analyze{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti})
664
CholmodFactor(ccall((@chm_nm "analyze" $Ti
665
,"libcholmod.so.1.7.1"), Ptr{c_CholmodFactor{Tv,$Ti}},
666
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}), &A.c, cmn($Ti)))
667
end
668
function cholfact!{Tv<:CHMVTypes}(L::CholmodFactor{Tv,$Ti},A::CholmodSparse{Tv,$Ti})
669
status = ccall((@chm_nm "factorize" $Ti
670
,"libcholmod.so.1.7.1"), Cint,
671
(Ptr{c_CholmodSparse{Tv,$Ti}},
672
Ptr{c_CholmodFactor{Tv,$Ti}}, Ptr{Uint8}),
673
&A.c, &L.c, cmn($Ti))
674
if status != CHOLMOD_TRUE throw(CholmodException) end
675
L
676
end
677
function cholfact!{Tv<:Float64}(L::CholmodFactor{Tv,$Ti},A::CholmodSparse{Tv,$Ti},
678
beta::Tv)
679
status = ccall((@chm_nm "factorize_p" $Ti
680
,"libcholmod.so.1.7.1"), Cint,
681
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Tv}, Ptr{Cint}, Csize_t,
682
Ptr{c_CholmodFactor{Tv,$Ti}}, Ptr{Uint8}),
683
&A.c, &beta, C_NULL, zero(Csize_t), &L.c, cmn($Ti))
684
if status != CHOLMOD_TRUE throw(CholmodException) end
685
L
686
end
687
function chm_pack!{Tv<:CHMVTypes}(L::CholmodFactor{Tv,$Ti})
688
status = ccall((@chm_nm "pack_factor" $Ti
689
,"libcholmod.so.1.7.1"), Cint,
690
(Ptr{c_CholmodFactor{Tv,$Ti}}, Ptr{Uint8}),
691
&L.c,cmn($Ti))
692
if status != CHOLMOD_TRUE throw(CholmodException) end
693
L
694
end
695
function chm_print{Tv<:CHMVTypes}(L::CholmodFactor{Tv,$Ti},lev,nm)
696
cm = cmn($Ti)
697
orig = cm[chm_prt_inds]
698
cm[chm_prt_inds] = reinterpret(Uint8, [int32(lev)])
699
status = ccall((@chm_nm "print_factor" $Ti
700
,"libcholmod.so.1.7.1"), Cint,
701
(Ptr{c_CholmodFactor{Tv,$Ti}}, Ptr{Uint8}, Ptr{Uint8}),
702
&L.c, nm, cm)
703
cm[chm_prt_inds] = orig
704
if status != CHOLMOD_TRUE throw(CholmodException) end
705
end
706
function chm_print{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti},lev,nm)
707
cm = cmn($Ti)
708
orig = cm[chm_prt_inds]
709
cm[chm_prt_inds] = reinterpret(Uint8, [int32(lev)])
710
status = ccall((@chm_nm "print_sparse" $Ti
711
,"libcholmod.so.1.7.1"), Cint,
712
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}, Ptr{Uint8}),
713
&A.c, nm, cm)
714
cm[chm_prt_inds] = orig
715
if status != CHOLMOD_TRUE throw(CholmodException) end
716
end
717
function chm_scale!{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti},
718
S::CholmodDense{Tv},
719
typ::Integer)
720
status = ccall((@chm_nm "scale" $Ti
721
,"libcholmod.so.1.7.1"), Cint,
722
(Ptr{c_CholmodDense{Tv}},Cint,Ptr{c_CholmodSparse{Tv,$Ti}},
723
Ptr{Uint8}), &S.c, typ, &A.c, cmn($Ti))
724
if status != CHOLMOD_TRUE throw(CholmodException) end
725
A
726
end
727
function chm_sdmult{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti},
728
trans::Bool,
729
alpha::Real,
730
beta::Real,
731
X::CholmodDense{Tv})
732
m,n = size(A)
733
nc = trans ? m : n
734
nr = trans ? n : m
735
if nc != size(X,1)
736
error("Incompatible dimensions, $nc and $(size(X,1)), in chm_sdmult")
737
end
738
aa = float64([alpha, 0.])
739
bb = float64([beta, 0.])
740
Y = CholmodDense(zeros(Tv,nr,size(X,2)))
741
status = ccall((@chm_nm "sdmult" $Ti
742
,"libcholmod.so.1.7.1"), Cint,
743
(Ptr{c_CholmodSparse{Tv,$Ti}},Cint,Ptr{Cdouble},Ptr{Cdouble},
744
Ptr{c_CholmodDense{Tv}}, Ptr{c_CholmodDense{Tv}}, Ptr{Uint8}),
745
&A.c,trans,aa,bb,&X.c,&Y.c,cmn($Ti))
746
if status != CHOLMOD_TRUE throw(CholmodException) end
747
Y
748
end
749
function chm_speye{Tv<:CHMVTypes,Ti<:$Ti}(m::Ti, n::Ti, x::Tv)
750
CholmodSparse(ccall((@chm_nm "speye" $Ti
751
, "libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
752
(Int, Int, Cint, Ptr{Uint8}),
753
m, n, xtyp{Tv}, cmn($Ti)))
754
end
755
function chm_spzeros{Tv<:Union(Float64,Complex128)}(m::$Ti, n::$Ti, nzmax::$Ti, x::Tv)
756
CholmodSparse(ccall((@chm_nm "spzeros" $Ti
757
, "libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
758
(Int, Int, Int, Ptr{Uint8}),
759
m, n, nzmax, xtyp{Tv}, cmn($Ti)))
760
end
761
## add chm_xtype and chm_pack
762
function copy{Tv<:CHMVTypes}(L::CholmodFactor{Tv,$Ti})
763
CholmodFactor(ccall((@chm_nm "copy_factor" $Ti
764
,"libcholmod.so.1.7.1"), Ptr{c_CholmodFactor{Tv,$Ti}},
765
(Ptr{c_CholmodFactor{Tv,$Ti}},Ptr{Uint8}), &L.c, cmn($Ti)))
766
end
767
function copy{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti})
768
CholmodSparse(ccall((@chm_nm "copy_sparse" $Ti
769
,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
770
(Ptr{c_CholmodSparse{Tv,$Ti}},Ptr{Uint8}), &A.c, cmn($Ti)))
771
end
772
function copy{Tv<:CHMVTypes}(T::CholmodTriplet{Tv,$Ti})
773
CholmodTriplet(ccall((@chm_nm "copy_triplet" $Ti
774
,"libcholmod.so.1.7.1"), Ptr{c_CholmodTriplet{Tv,$Ti}},
775
(Ptr{c_CholmodTriplet{Tv,$Ti}},Ptr{Uint8}), &T.c, cmn($Ti)))
776
end
777
function ctranspose{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti})
778
CholmodSparse(ccall((@chm_nm "transpose" $Ti
779
,"libcholmod.so.1.7.1"),Ptr{c_CholmodSparse{Tv,$Ti}},
780
(Ptr{c_CholmodSparse{Tv,$Ti}}, Cint, Ptr{Uint8}),
781
&A.c, 2, cmn($Ti)))
782
end
783
function etree{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti})
784
tr = Array($Ti,size(A,2))
785
status = ccall((@chm_nm "etree" $Ti
786
,"libcholmod.so.1.7.1"), Cint,
787
(Ptr{c_CholmodSparse{Tv,$Ti}},Ptr{$Ti},Ptr{Uint8}),
788
&A.c,tr,cmn($Ti))
789
if status != CHOLMOD_TRUE throw(CholmodException) end
790
tr
791
end
792
function hcat{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti},B::CholmodSparse{Tv,$Ti})
793
ccall((@chm_nm "horzcat" $Ti
794
, "libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
795
(Ptr{c_CholmodSparse{Tv,$Ti}},Ptr{c_CholmodSparse{Tv,$Ti}},Cint,Ptr{Uint8}),
796
&A.c,&B.c,true,cmn($Ti))
797
end
798
function nnz{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti})
799
ccall((@chm_nm "nnz" $Ti
800
,"libcholmod.so.1.7.1"), Int, (Ptr{c_CholmodSparse{Tv,$Ti}},Ptr{Uint8}),&A.c,cmn($Ti))
801
end
802
function norm{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti},p::Number)
803
ccall((@chm_nm "norm_sparse" $Ti
804
, "libcholmod.so.1.7.1"), Float64,
805
(Ptr{c_CholmodSparse{Tv,$Ti}}, Cint, Ptr{Uint8}),
806
&A.c,p == 1 ? 1 :(p == Inf ? 1 : error("p must be 1 or Inf")),cmn($Ti))
807
end
808
function solve{Tv<:CHMVTypes}(L::CholmodFactor{Tv,$Ti},
809
B::CholmodDense{Tv}, typ::Integer)
810
CholmodDense(ccall((@chm_nm "solve" $Ti
811
,"libcholmod.so.1.7.1"), Ptr{c_CholmodDense{Tv}},
812
(Cint, Ptr{c_CholmodFactor{Tv,$Ti}},
813
Ptr{c_CholmodDense{Tv}}, Ptr{Uint8}),
814
typ, &L.c, &B.c, cmn($Ti)))
815
end
816
function solve{Tv<:CHMVTypes}(L::CholmodFactor{Tv,$Ti},
817
B::CholmodSparse{Tv,$Ti},
818
typ::Integer)
819
CholmodSparse(ccall((@chm_nm "spsolve" $Ti
820
,"libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
821
(Cint, Ptr{c_CholmodFactor{Tv,$Ti}},
822
Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}),
823
typ, &L.c, &B.c, cmn($Ti)))
824
end
825
function sort!{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti})
826
status = ccall((@chm_nm "sort" $Ti
827
,"libcholmod.so.1.7.1"), Cint,
828
(Ptr{c_CholmodSparse{Tv,$Ti}}, Ptr{Uint8}),
829
&A.c, cmn($Ti))
830
if status != CHOLMOD_TRUE throw(CholmodException) end
831
A
832
end
833
function copysym{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti})
834
CholmodSparse(ccall((@chm_nm "copy" $Ti
835
,"libcholmod.so.1.7.1"),Ptr{c_CholmodSparse{Tv,$Ti}},
836
(Ptr{c_CholmodSparse{Tv,$Ti}},Cint,Cint,Ptr{Uint8}),
837
&A.c,0,1,cmn($Ti)))
838
end
839
function transpose{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti})
840
CholmodSparse(ccall((@chm_nm "transpose" $Ti
841
,"libcholmod.so.1.7.1"),Ptr{c_CholmodSparse{Tv,$Ti}},
842
(Ptr{c_CholmodSparse{Tv,$Ti}}, Cint, Ptr{Uint8}),
843
&A.c, 1, cmn($Ti)))
844
end
845
function vcat{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti},B::CholmodSparse{Tv,$Ti})
846
ccall((@chm_nm "vertcat" $Ti
847
, "libcholmod.so.1.7.1"), Ptr{c_CholmodSparse{Tv,$Ti}},
848
(Ptr{c_CholmodSparse{Tv,$Ti}},Ptr{c_CholmodSparse{Tv,$Ti}},Cint,Ptr{Uint8}),
849
&A.c,&B.c,true,cmn($Ti))
850
end
851
end
852
end
853
(*){Tv<:CHMVTypes}(A::CholmodSparse{Tv},B::CholmodDense{Tv}) = chm_sdmult(A,false,1.,0.,B)
854
(*){Tv<:CHMVTypes}(A::CholmodSparse{Tv},B::VecOrMat{Tv}) = chm_sdmult(A,false,1.,0.,CholmodDense(B))
855
856
(\){T<:CHMVTypes}(L::CholmodFactor{T},B::CholmodDense{T}) = solve(L,B,CHOLMOD_A)
857
(\){T<:CHMVTypes}(L::CholmodFactor{T},B::VecOrMat{T}) = solve(L,CholmodDense!(B),CHOLMOD_A).mat
858
function (\){Tv<:CHMVTypes,Ti<:CHMITypes}(L::CholmodFactor{Tv,Ti},B::CholmodSparse{Tv,Ti})
859
solve(L,B,CHOLMOD_A)
860
end
861
function (\){Tv<:CHMVTypes,Ti<:CHMITypes}(L::CholmodFactor{Tv,Ti},B::SparseMatrixCSC{Tv,Ti})
862
sparse!(solve(L,CholmodSparse(B),CHOLMOD_A))
863
end
864
865
function A_mul_Bt{Tv<:Union(Float32,Float64),Ti<:CHMITypes}(A::CholmodSparse{Tv,Ti},
866
B::CholmodSparse{Tv,Ti})
867
A_mul_Bc(A,B) # in the unlikely event of writing A*B.' instead of A*B'
868
end
869
870
Ac_ldiv_B{T<:CHMVTypes}(L::CholmodFactor{T},B::CholmodDense{T}) = solve(L,B,CHOLMOD_A)
871
Ac_ldiv_B{T<:CHMVTypes}(L::CholmodFactor{T},B::VecOrMat{T}) = solve(L,CholmodDense!(B),CHOLMOD_A).mat
872
function Ac_ldiv_B{Tv<:CHMVTypes,Ti<:CHMITypes}(L::CholmodFactor{Tv,Ti},B::CholmodSparse{Tv,Ti})
873
solve(L,B,CHOLMOD_A)
874
end
875
function Ac_ldiv_B{Tv<:CHMVTypes,Ti<:CHMITypes}(L::CholmodFactor{Tv,Ti},B::SparseMatrixCSC{Tv,Ti})
876
sparse!(solve(L,CholmodSparse(B),CHOLMOD_A))
877
end
878
879
function Ac_mul_B{Tv<:CHMVTypes}(A::CholmodSparse{Tv},B::CholmodDense{Tv})
880
chm_sdmult(A,true,1.,0.,B)
881
end
882
function Ac_mul_B{Tv<:CHMVTypes}(A::CholmodSparse{Tv},B::VecOrMat{Tv})
883
chm_sdmult(A,true,1.,0.,CholmodDense(B))
884
end
885
886
function At_mul_B{Tv<:Union(Float32,Float64),Ti<:CHMITypes}(A::CholmodSparse{Tv,Ti},
887
B::CholmodSparse{Tv,Ti})
888
Ac_mul_B(A,B) # in the unlikely event of writing A.'*B instead of A'*B
889
end
890
891
cholfact{T<:CHMVTypes}(A::CholmodSparse{T},beta::T) = cholfact(A,beta,false)
892
cholfact(A::CholmodSparse) = cholfact(A,false)
893
cholfact(A::SparseMatrixCSC,ll::Bool) = cholfact(CholmodSparse(A),ll)
894
cholfact(A::SparseMatrixCSC) = cholfact(CholmodSparse(A),false)
895
function cholfact!{T<:CHMVTypes}(L::CholmodFactor{T},A::CholmodSparse{T},beta::Number)
896
cholfact!(L,A,convert(T,beta))
897
end
898
function cholfact!{T<:CHMVTypes}(L::CholmodFactor{T},A::SparseMatrixCSC{T},beta::Number)
899
cholfact!(L,CholmodSparse(A),convert(T,beta))
900
end
901
cholfact!{T<:CHMVTypes}(L::CholmodFactor{T},A::SparseMatrixCSC{T}) = cholfact!(L,CholmodSparse(A))
902
903
chm_analyze(A::SparseMatrixCSC) = chm_analyze(CholmodSparse(A))
904
905
chm_print(A::CholmodSparse, lev::Integer) = chm_print(A, lev, "")
906
chm_print(A::CholmodFactor, lev::Integer) = chm_print(L, lev, "")
907
908
function chm_scale!{T<:CHMVTypes}(A::CholmodSparse{T},b::Vector{T},typ::Integer)
909
chm_scale!(A,CholmodDense(b),typ)
910
A
911
end
912
chm_scale{T<:CHMVTypes}(A::CholmodSparse{T},b::Vector{T},typ::Integer) = chm_scale!(copy(A),b,typ)
913
914
chm_speye(m::Integer, n::Integer) = chm_speye(m, n, 1., 1) # default element type is Float32
915
chm_speye(n::Integer) = chm_speye(n, n, 1.) # default shape is square
916
917
chm_spzeros(m::Integer,n::Integer,nzmax::Integer) = chm_spzeros(m,n,nzmax,1.)
918
919
function scale!{T<:CHMVTypes}(b::Vector{T}, A::CholmodSparse{T})
920
chm_scale!(A,CholmodDense(b),CHOLMOD_ROW)
921
A
922
end
923
scale{T<:CHMVTypes}(b::Vector{T}, A::CholmodSparse{T}) = scale!(b,copy(A))
924
function scale!{T<:CHMVTypes}(A::CholmodSparse{T},b::Vector{T})
925
chm_scale!(A,CholmodDense(b),CHOLMOD_COL)
926
A
927
end
928
scale{T<:CHMVTypes}(A::CholmodSparse{T},b::Vector{T}) = scale!(copy(A), b)
929
930
norm(A::CholmodSparse) = norm(A,1)
931
932
show(io::IO,L::CholmodFactor) = chm_print(L,int32(4),"")
933
show(io::IO,A::CholmodSparse) = chm_print(A,int32(4),"")
934
935
size(B::CholmodDense) = size(B.mat)
936
size(B::CholmodDense,d) = size(B.mat,d)
937
size(A::CholmodSparse) = (int(A.c.m), int(A.c.n))
938
function size(A::CholmodSparse, d::Integer)
939
d == 1 ? A.c.m : (d == 2 ? A.c.n : 1)
940
end
941
size(L::CholmodFactor) = (n = int(L.c.n); (n,n))
942
size(L::CholmodFactor,d::Integer) = d < 1 ? error("dimension out of range") : (d <= 2 ? int(L.c.n) : 1)
943
944
function solve{Tv<:CHMVTypes,Ti<:CHMITypes}(L::CholmodFactor{Tv,Ti},
945
B::SparseMatrixCSC{Tv,Ti},typ::Integer)
946
sparse!(solve(L,CholmodSparse(B),typ))
947
end
948
function solve{Tv<:CHMVTypes,Ti<:CHMITypes}(L::CholmodFactor{Tv,Ti},B::SparseMatrixCSC{Tv,Ti})
949
sparse!(solve(L,CholmodSparse(B),CHOLMOD_A))
950
end
951
function solve{Tv<:CHMVTypes,Ti<:CHMITypes}(L::CholmodFactor{Tv,Ti},B::CholmodSparse{Tv,Ti})
952
solve(L,B,CHOLMOD_A)
953
end
954
function (\){Tv<:CHMVTypes,Ti<:CHMITypes}(L::CholmodFactor{Tv,Ti},B::CholmodSparse{Tv,Ti})
955
solve(L,B,CHOLMOD_A)
956
end
957
function solve{T<:CHMVTypes}(L::CholmodFactor{T},B::Vector{T},typ::Integer)
958
solve(L,CholmodDense!(B),typ).mat[:]
959
end
960
function solve{T<:CHMVTypes}(L::CholmodFactor{T},B::Matrix{T},typ::Integer)
961
solve(L,CholmodDense!(B),typ).mat
962
end
963
solve{T<:CHMVTypes}(L::CholmodFactor{T},B::CholmodDense{T}) = solve(L,B,CHOLMOD_A)
964
965
function findnz{Tv,Ti}(A::CholmodSparse{Tv,Ti})
966
jj = similar(A.rowval0) # expand A.colptr0 to a vector of indices
967
for j in 1:A.c.n, k in (A.colptr0[j]+1):A.colptr0[j+1]
968
jj[k] = j
969
end
970
971
ind = similar(A.rowval0)
972
ipos = 1
973
count = 0
974
for k in 1:length(A.nzval)
975
if A.nzval[k] != 0
976
ind[ipos] = k
977
ipos += 1
978
count += 1
979
else
980
println("Warning: sparse matrix contains explicitly stored zeros.")
981
end
982
end
983
ind = ind[1:count] # ind is the indices of nonzeros in A.nzval
984
(increment!(A.rowval0[ind]), jj[ind], A.nzval[ind])
985
end
986
987
findnz(L::CholmodFactor) = findnz(CholmodSparse(L))
988
989
function diag{Tv}(A::CholmodSparse{Tv})
990
minmn = minimum(size(A))
991
res = zeros(Tv,minmn)
992
cp0 = A.colptr0
993
rv0 = A.rowval0
994
anz = A.nzval
995
for j in 1:minmn, k in (cp0[j]+1):cp0[j+1]
996
if rv0[k] == j-1
997
res[j] += anz[k]
998
end
999
end
1000
res
1001
end
1002
1003
function diag{Tv}(L::CholmodFactor{Tv})
1004
res = zeros(Tv,L.c.n)
1005
xv = L.x
1006
if bool(L.c.is_super)
1007
nec = decrement!(diff(L.super)) # number of excess columns per supernode
1008
dstride = increment!(diff(L.pi)) # stride of diagonal elements (nrow + 1)
1009
px = L.px
1010
pos = 1
1011
for i in 1:length(nec)
1012
base = px[i] + 1
1013
res[pos] = xv[base]
1014
pos += 1
1015
for j in 1:nec[i]
1016
res[pos] = xv[base + j*dstride[i]]
1017
pos += 1
1018
end
1019
end
1020
else
1021
c0 = L.p
1022
r0 = L.i
1023
xv = L.x
1024
for j in 1:L.c.n
1025
jj = c0[j]+1
1026
assert(r0[jj] == j-1)
1027
res[j] = xv[jj]
1028
end
1029
end
1030
res
1031
end
1032
1033
function logdet{Tv,Ti}(L::CholmodFactor{Tv,Ti},sub)
1034
res = zero(Tv)
1035
for d in diag(L)[sub] res += log(abs(d)) end
1036
bool(L.c.is_ll) ? 2res : res
1037
end
1038
logdet(L::CholmodFactor) = logdet(L, 1:L.c.n)
1039
det(L::CholmodFactor) = exp(logdet(L))
1040
1041
dense(A::CholmodSparse) = CholmodDense(A).mat
1042
dense(A::CholmodDense) = A.mat
1043
full(A::CholmodSparse) = dense(A)
1044
function sparse(A::CholmodSparse)
1045
if bool(A.c.stype) return sparse!(copysym(A)) end
1046
SparseMatrixCSC(A.c.m, A.c.n, increment(A.colptr0), increment(A.rowval0), copy(A.nzval))
1047
end
1048
function sparse!(A::CholmodSparse)
1049
SparseMatrixCSC(A.c.m, A.c.n, increment!(A.colptr0), increment!(A.rowval0), A.nzval)
1050
end
1051
sparse(L::CholmodFactor) = sparse!(CholmodSparse(L))
1052
sparse(D::CholmodDense) = sparse!(CholmodSparse(D))
1053
sparse(T::CholmodTriplet) = sparse!(CholmodSparse(T))
1054
end #module
Loggerhead is a web-based interface for Breezy
Version: 2.0.1