wolffd@0
|
1 function pot = direct_combine_pots(pot1, pot2)
|
wolffd@0
|
2 % DIRECTED_COMBINE_POTS The combination operation corresponds to ordinary composition of conditional distributions.
|
wolffd@0
|
3 % In some sense is similar to that of forming disjoint union of set.
|
wolffd@0
|
4 % pot = direct_combine_pots(pot1, pot2)
|
wolffd@0
|
5
|
wolffd@0
|
6 % directed combine can be performed under the conditon that the head node set of pot1 is disjoint from the domain of
|
wolffd@0
|
7 % pot2 or vice versa. if the last conditon was satisfied we exchange the pot1 and pot2 firstly then perform the operation.
|
wolffd@0
|
8 % If neither of them was satified the directed combine is undifined.
|
wolffd@0
|
9
|
wolffd@0
|
10
|
wolffd@0
|
11 if isempty( myintersect(pot1.domain, pot2.cheaddom) )
|
wolffd@0
|
12 pot1 = pot1;
|
wolffd@0
|
13 pot2 = pot2;
|
wolffd@0
|
14 elseif isempty( myintersect(pot2.domain, pot1.cheaddom))
|
wolffd@0
|
15 temppot = pot1;
|
wolffd@0
|
16 pot1 = pot2;
|
wolffd@0
|
17 pot2 = temppot;
|
wolffd@0
|
18 else
|
wolffd@0
|
19 assert(0);
|
wolffd@0
|
20 return;
|
wolffd@0
|
21 end
|
wolffd@0
|
22
|
wolffd@0
|
23 domain = myunion(pot1.domain, pot2.domain);
|
wolffd@0
|
24 nodesizes = zeros(1,max(domain));
|
wolffd@0
|
25 nodesizes(pot2.ctaildom) = pot2.ctailsizes;
|
wolffd@0
|
26 nodesizes(pot2.cheaddom) = pot2.cheadsizes;
|
wolffd@0
|
27 nodesizes(pot2.ddom) = pot2.dsizes;
|
wolffd@0
|
28 nodesizes(pot1.ctaildom) = pot1.ctailsizes;
|
wolffd@0
|
29 nodesizes(pot1.cheaddom) = pot1.cheadsizes;
|
wolffd@0
|
30 nodesizes(pot1.ddom) = pot1.dsizes;
|
wolffd@0
|
31
|
wolffd@0
|
32 dom_u = mysetdiff(pot2.ctaildom, pot1.cheaddom);
|
wolffd@0
|
33 if ~isempty(dom_u) & ~mysubset(dom_u, pot1.ctaildom)
|
wolffd@0
|
34 pot1 = extension_pot(pot1, [], [], dom_u, nodesizes(dom_u));
|
wolffd@0
|
35 end
|
wolffd@0
|
36
|
wolffd@0
|
37 dom_u = myunion(pot1.cheaddom, pot1.ctaildom);
|
wolffd@0
|
38 if ~isempty(dom_u) & ~mysubset(dom_u, pot2.ctaildom)
|
wolffd@0
|
39 pot2 = extension_pot(pot2, [], [], dom_u, nodesizes(dom_u));
|
wolffd@0
|
40 end
|
wolffd@0
|
41
|
wolffd@0
|
42
|
wolffd@0
|
43 cheaddom = myunion(pot1.cheaddom, pot2.cheaddom);
|
wolffd@0
|
44 ctaildom = mysetdiff(myunion(pot1.ctaildom, pot2.ctaildom), cheaddom);
|
wolffd@0
|
45 cdom = myunion(cheaddom, ctaildom);
|
wolffd@0
|
46 ddom = mysetdiff(domain, cdom);
|
wolffd@0
|
47 dsizes = nodesizes(ddom);
|
wolffd@0
|
48 dsize = prod(nodesizes(ddom));
|
wolffd@0
|
49 cheadsizes = nodesizes(cheaddom);
|
wolffd@0
|
50 cheadsize = sum(nodesizes(cheaddom));
|
wolffd@0
|
51 ctailsizes = nodesizes(ctaildom);
|
wolffd@0
|
52 ctailsize = sum(nodesizes(ctaildom));
|
wolffd@0
|
53
|
wolffd@0
|
54 r1 = pot1.cheadsize;
|
wolffd@0
|
55 s1 = pot1.ctailsize;
|
wolffd@0
|
56 scpot = cell(1, dsize);
|
wolffd@0
|
57 mask1 = [];
|
wolffd@0
|
58 mask2 = [];
|
wolffd@0
|
59 if ~isempty(pot1.ddom)
|
wolffd@0
|
60 mask1 = find_equiv_posns(pot1.ddom, ddom);
|
wolffd@0
|
61 end
|
wolffd@0
|
62 if ~isempty(pot2.ddom)
|
wolffd@0
|
63 mask2 = find_equiv_posns(pot2.ddom, ddom);
|
wolffd@0
|
64 end
|
wolffd@0
|
65 cmask1 = [];
|
wolffd@0
|
66 cmask2 = [];
|
wolffd@0
|
67 if ~isempty(pot1.cheaddom)
|
wolffd@0
|
68 cmask1 = find_equiv_posns(pot1.cheaddom, cheaddom);
|
wolffd@0
|
69 end
|
wolffd@0
|
70 if ~isempty(pot2.cheaddom)
|
wolffd@0
|
71 cmask2 = find_equiv_posns(pot2.cheaddom, cheaddom);
|
wolffd@0
|
72 end
|
wolffd@0
|
73
|
wolffd@0
|
74 u1 = block(cmask1, cheadsizes);
|
wolffd@0
|
75 u2 = block(cmask2, cheadsizes);
|
wolffd@0
|
76
|
wolffd@0
|
77 fmaskh = find_equiv_posns(pot1.cheaddom, pot2.ctaildom);
|
wolffd@0
|
78 fmaskt = find_equiv_posns(pot1.ctaildom, pot2.ctaildom);
|
wolffd@0
|
79
|
wolffd@0
|
80 fh = block(fmaskh, pot2.ctailsizes);
|
wolffd@0
|
81 ft = block(fmaskt, pot2.ctailsizes);
|
wolffd@0
|
82
|
wolffd@0
|
83 for i=1:dsize
|
wolffd@0
|
84 sub = ind2subv(dsizes, i);
|
wolffd@0
|
85 sub1 = sub(mask1);
|
wolffd@0
|
86 sub2 = sub(mask2);
|
wolffd@0
|
87 ind1 = subv2ind(pot1.dsizes, sub1);
|
wolffd@0
|
88 ind2 = subv2ind(pot2.dsizes, sub2);
|
wolffd@0
|
89
|
wolffd@0
|
90 if isempty(ind1)
|
wolffd@0
|
91 ind1 = 1;
|
wolffd@0
|
92 end
|
wolffd@0
|
93 if isempty(ind2)
|
wolffd@0
|
94 ind2 = 1;
|
wolffd@0
|
95 end
|
wolffd@0
|
96 potc1 = struct(pot1.scgpotc{ind1});
|
wolffd@0
|
97 potc2 = struct(pot2.scgpotc{ind2});
|
wolffd@0
|
98 p = potc1.p;
|
wolffd@0
|
99 q = potc2.p;
|
wolffd@0
|
100 ro = p*q;
|
wolffd@0
|
101
|
wolffd@0
|
102 A = potc1.A;
|
wolffd@0
|
103 B = potc1.B;
|
wolffd@0
|
104 C = potc1.C;
|
wolffd@0
|
105
|
wolffd@0
|
106 E = potc2.A;
|
wolffd@0
|
107 F = potc2.B;
|
wolffd@0
|
108 G = potc2.C;
|
wolffd@0
|
109
|
wolffd@0
|
110 F1 = F(:, fh);
|
wolffd@0
|
111 F2 = F(:, ft);
|
wolffd@0
|
112
|
wolffd@0
|
113 if ~isempty(F1)
|
wolffd@0
|
114 K1 = F1*A;
|
wolffd@0
|
115 K2 = F1*B;
|
wolffd@0
|
116 FCF = F1*C*F1';
|
wolffd@0
|
117 FC = F1*C;
|
wolffd@0
|
118 CFT = C*F1';
|
wolffd@0
|
119 else
|
wolffd@0
|
120 K1 = zeros(size(E));
|
wolffd@0
|
121 K2 = zeros(size(F2));
|
wolffd@0
|
122 FCF = zeros(size(G));
|
wolffd@0
|
123 FC = zeros(size(C, 1), size(G, 2));
|
wolffd@0
|
124 CFT = zeros(size(G, 2), size(C, 1));
|
wolffd@0
|
125 end
|
wolffd@0
|
126
|
wolffd@0
|
127
|
wolffd@0
|
128 U = zeros(cheadsize,1);
|
wolffd@0
|
129 W = zeros(cheadsize,cheadsize);
|
wolffd@0
|
130 V = zeros(cheadsize,ctailsize);
|
wolffd@0
|
131
|
wolffd@0
|
132 if cheadsize > 0
|
wolffd@0
|
133 U(u1) = A;
|
wolffd@0
|
134 U(u2) = E + K1;
|
wolffd@0
|
135 W(u1, u1) = C;
|
wolffd@0
|
136 W(u2, u2) = G + FCF;
|
wolffd@0
|
137 W(u1, u2) = CFT;
|
wolffd@0
|
138 W(u2, u1) = FC;
|
wolffd@0
|
139 else
|
wolffd@0
|
140 U = zeros(cheadsize,1);
|
wolffd@0
|
141 W = zeros(cheadsize,cheadsize);
|
wolffd@0
|
142 end
|
wolffd@0
|
143 if cheadsize > 0 | ctailsize > 0
|
wolffd@0
|
144 if ~isempty(u1)
|
wolffd@0
|
145 V(u1, :) = B;
|
wolffd@0
|
146 else
|
wolffd@0
|
147 V(u1, :) = zeros(potc1.cheadsize, ctailsize);
|
wolffd@0
|
148 end
|
wolffd@0
|
149 if ~isempty(u2)
|
wolffd@0
|
150 V(u2, :) = F2 + K2;
|
wolffd@0
|
151 else
|
wolffd@0
|
152 V(u2, :) = zeros(potc2.cheadsize, ctailsize);
|
wolffd@0
|
153 end
|
wolffd@0
|
154 else
|
wolffd@0
|
155 V = zeros(cheadsize,ctailsize);
|
wolffd@0
|
156 end
|
wolffd@0
|
157
|
wolffd@0
|
158 scpot{i} = scgcpot(cheadsize, ctailsize, ro, U, V, W);
|
wolffd@0
|
159 end
|
wolffd@0
|
160
|
wolffd@0
|
161 pot = scgpot(ddom, cheaddom, ctaildom, nodesizes, scpot);
|