1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
|
function [c,phi,flag]=min_qcost_flow(eps,g)
// Copyright INRIA
[lhs,rhs]=argn(0)
if rhs<>2 then error(39), end
// check eps
if prod(size(eps))<>1 then
error('First argument must be a scalar')
end
if eps<=0 then
error('First argument must be strictly positive')
end
p=-log(eps)/log(2)
// check g
check_graph(g)
// check capacities
ma=prod(size(g('tail')))
n=g('node_number')
mincap=g('edge_min_cap')
maxcap=g('edge_max_cap')
if mincap==[] then
mincap=zeros(1,ma)
end
if maxcap==[] then
maxcap=zeros(1,ma)
end
verif=find(mincap>maxcap)
if verif<>[] then
error('Maximum capacities must be greater than minimal capacities')
end
verif=find(mincap<>maxcap)
if verif==[] then
error('Maximum capacities must not be all equal to minimal capacities')
end
// check quadratic costs
qorig=g('edge_q_orig')
qweight=g('edge_q_weight')
if qorig==[] then
qorig=zeros(1,ma)
end
if qweight==[] then
qweight=zeros(1,ma)
end
// check demand
demand=g('node_demand')
if demand==[] then
demand=zeros(1,n)
end
verif=find(demand<>0)
if verif<>[] then
error('Demands must be equal to zero')
end
// compute lp, la and ls
// always consider the graph as undirected!
[lp,la,ls]=m6ta2lpu(g('tail'),g('head'),n+1,n,2*ma)
// compute min quadratic cost flow
[phi,flag]=m6floqua(p,mincap,maxcap,g('head'),g('tail'),la,lp,n,qorig,qweight)
c=sum(qweight.*(phi-qorig).*(phi-qorig))/2
|
