Context Navigation

network.c @ b826dcf

Visit:

main

Last change on this file since b826dcf was 0b99107, checked in by Olly Betts <olly@…>, 8 weeks ago

Eliminate old FSF addresses

Update GPL/LGPL licence files and boilerplate to direct people who
didn't receive the licence text to the FSF website, as the current
versions of the FSF licence texts now do, rather than giving a postal
address.

Property mode set to 100644

File size: 21.6 KB

Line
1	/* network.c
2	* Survex network reduction - find patterns and apply network reductions
3	* Copyright (C) 1991-2002,2005,2024 Olly Betts
4	*
5	* This program is free software; you can redistribute it and/or modify
6	* it under the terms of the GNU General Public License as published by
7	* the Free Software Foundation; either version 2 of the License, or
8	* (at your option) any later version.
9	*
10	* This program is distributed in the hope that it will be useful,
11	* but WITHOUT ANY WARRANTY; without even the implied warranty of
12	* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13	* GNU General Public License for more details.
14	*
15	* You should have received a copy of the GNU General Public License
16	* along with this program; if not, see
17	* <https://www.gnu.org/licenses/>.
18	*/
19
20	#if 0
21	#define DEBUG_INVALID 1
22	#define VALIDATE 1
23	#define DUMP_NETWORK 1
24	#endif
25
26	#include <config.h>
27
28	#include "validate.h"
29	#include "debug.h"
30	#include "cavern.h"
31	#include "message.h"
32	#include "netbits.h"
33	#include "network.h"
34	#include "osalloc.h"
35	#include "out.h"
36
37	typedef struct reduction {
38	struct reduction *next;
39	enum {
40	TYPE_PARALLEL,
41	TYPE_LOLLIPOP,
42	TYPE_DELTASTAR,
43	} type;
44	linkfor* join[];
45	} reduction;
46
47	#define allocate_reduction(N) osmalloc(sizeof(reduction) + (N) * sizeof(linkfor*))
48
49	// Head of linked list of reductions.
50	static reduction *reduction_stack;
51
52	/* can be altered by -z<letters> on command line */
53	unsigned long optimize = BITA('l') \| BITA('p') \| BITA('d');
54	/* Lollipops, Parallel legs, Iterate mx, Delta* */
55
56	extern void
57	remove_subnets(void)
58	{
59	node stn, stn2, stn3, stn4;
60	int dirn, dirn2, dirn3, dirn4;
61	reduction *trav;
62	linkfor newleg, newleg2;
63	bool fMore = true;
64
65	reduction_stack = NULL;
66
67	out_current_action(msg(/Simplifying network/129));
68
69	while (fMore) {
70	fMore = false;
71	if (optimize & BITA('l')) {
72	#if PRINT_NETBITS
73	printf("replacing lollipops\n");
74	#endif
75	/* NB can have non-fixed 0 nodes */
76	FOR_EACH_STN(stn, stnlist) {
77	if (three_node(stn)) {
78	dirn = -1;
79	if (stn->leg[1]->l.to == stn) dirn++;
80	if (stn->leg[0]->l.to == stn) dirn += 2;
81	if (dirn < 0) continue;
82
83	stn2 = stn->leg[dirn]->l.to;
84	if (fixed(stn2)) {
85	/* _
86	* ( )
87	* * stn
88	* \|
89	* * stn2 (fixed)
90	* : (may have other connections)
91	*
92	* The leg forming the "stick" of the lollipop is
93	* articulating so we can just fix stn with coordinates
94	* calculated by adding or subtracting the leg's vector.
95	*/
96	linkfor *leg = stn->leg[dirn];
97	linkfor *rev_leg = reverse_leg(leg);
98	leg->l.reverse \|= FLAG_ARTICULATION;
99	rev_leg->l.reverse \|= FLAG_ARTICULATION;
100	if (data_here(leg)) {
101	subdd(&POSD(stn), &POSD(stn2), &leg->d);
102	} else {
103	adddd(&POSD(stn), &POSD(stn2), &rev_leg->d);
104	}
105	remove_stn_from_list(&stnlist, stn);
106	add_stn_to_list(&fixedlist, stn);
107	continue;
108	}
109
110	SVX_ASSERT(three_node(stn2));
111
112	/* _
113	* ( )
114	* * stn
115	* \|
116	* * stn2
117	* / \
118	* stn4 * * stn3 --> stn4 --- stn3
119	* : : : :
120	*/
121	dirn2 = reverse_leg_dirn(stn->leg[dirn]);
122	dirn2 = (dirn2 + 1) % 3;
123	stn3 = stn2->leg[dirn2]->l.to;
124	if (stn2 == stn3) continue; /* dumb-bell - leave alone */
125
126	dirn3 = reverse_leg_dirn(stn2->leg[dirn2]);
127
128	trav = allocate_reduction(2);
129	trav->type = TYPE_LOLLIPOP;
130
131	newleg2 = (linkfor*)osnew(linkcommon);
132
133	newleg = copy_link(stn3->leg[dirn3]);
134
135	dirn2 = (dirn2 + 1) % 3;
136	stn4 = stn2->leg[dirn2]->l.to;
137	dirn4 = reverse_leg_dirn(stn2->leg[dirn2]);
138	#if 0
139	printf("Lollipop found with stn...stn4 = \n");
140	print_prefix(stn->name); putnl();
141	print_prefix(stn2->name); putnl();
142	print_prefix(stn3->name); putnl();
143	print_prefix(stn4->name); putnl();
144	#endif
145
146	addto_link(newleg, stn2->leg[dirn2]);
147
148	/* remove stn and stn2 */
149	remove_stn_from_list(&stnlist, stn);
150	remove_stn_from_list(&stnlist, stn2);
151
152	/* stack lollipop and replace with a leg between stn3 and stn4 */
153	trav->join[0] = stn3->leg[dirn3];
154	newleg->l.to = stn4;
155	newleg->l.reverse = dirn4 \| FLAG_DATAHERE \| FLAG_REPLACEMENTLEG;
156
157	trav->join[1] = stn4->leg[dirn4];
158	newleg2->l.to = stn3;
159	newleg2->l.reverse = dirn3 \| FLAG_REPLACEMENTLEG;
160
161	stn3->leg[dirn3] = newleg;
162	stn4->leg[dirn4] = newleg2;
163
164	trav->next = reduction_stack;
165	#if PRINT_NETBITS
166	printf("remove lollipop\n");
167	#endif
168	reduction_stack = trav;
169	fMore = true;
170	}
171	}
172	}
173
174	if (optimize & BITA('p')) {
175	#if PRINT_NETBITS
176	printf("replacing parallel legs\n");
177	#endif
178	FOR_EACH_STN(stn, stnlist) {
179	/*
180	* : :
181	* * stn3 * stn3
182	* \| \|
183	* * stn \|
184	* ( ) --> \|
185	* * stn2 \|
186	* \| \|
187	* * stn4 * stn4
188	* : :
189	*/
190	if (three_node(stn)) {
191	stn2 = stn->leg[0]->l.to;
192	if (stn2 == stn->leg[1]->l.to) {
193	dirn = 2;
194	} else if (stn2 == stn->leg[2]->l.to) {
195	dirn = 1;
196	} else {
197	if (stn->leg[1]->l.to != stn->leg[2]->l.to) continue;
198	stn2 = stn->leg[1]->l.to;
199	dirn = 0;
200	}
201
202	/* stn == stn2 => lollipop */
203	if (stn == stn2 \|\| fixed(stn2)) continue;
204
205	SVX_ASSERT(three_node(stn2));
206
207	stn3 = stn->leg[dirn]->l.to;
208	/* 3 parallel legs (=> nothing else) so leave */
209	if (stn3 == stn2) continue;
210
211	dirn3 = reverse_leg_dirn(stn->leg[dirn]);
212	dirn2 = (0 + 1 + 2 - reverse_leg_dirn(stn->leg[(dirn + 1) % 3])
213	- reverse_leg_dirn(stn->leg[(dirn + 2) % 3]));
214
215	stn4 = stn2->leg[dirn2]->l.to;
216	dirn4 = reverse_leg_dirn(stn2->leg[dirn2]);
217
218	trav = allocate_reduction(2);
219	trav->type = TYPE_PARALLEL;
220
221	newleg = copy_link(stn->leg[(dirn + 1) % 3]);
222	/* use newleg2 for scratch */
223	newleg2 = copy_link(stn->leg[(dirn + 2) % 3]);
224	{
225	#ifdef NO_COVARIANCES
226	vars sum;
227	var prod;
228	delta temp, temp2;
229	addss(&sum, &newleg->v, &newleg2->v);
230	SVX_ASSERT2(!fZeros(&sum), "loop of zero variance found");
231	mulss(&prod, &newleg->v, &newleg2->v);
232	mulsd(&temp, &newleg2->v, &newleg->d);
233	mulsd(&temp2, &newleg->v, &newleg2->d);
234	adddd(&temp, &temp, &temp2);
235	divds(&newleg->d, &temp, &sum);
236	sdivvs(&newleg->v, &prod, &sum);
237	#else
238	svar inv1, inv2, sum;
239	delta temp, temp2;
240	/* if leg one is an equate, we can just ignore leg two
241	* whatever it is */
242	if (invert_svar(&inv1, &newleg->v)) {
243	if (invert_svar(&inv2, &newleg2->v)) {
244	addss(&sum, &inv1, &inv2);
245	if (!invert_svar(&newleg->v, &sum)) {
246	BUG("matrix singular in parallel legs replacement");
247	}
248
249	mulsd(&temp, &inv1, &newleg->d);
250	mulsd(&temp2, &inv2, &newleg2->d);
251	adddd(&temp, &temp, &temp2);
252	mulsd(&newleg->d, &newleg->v, &temp);
253	} else {
254	/* leg two is an equate, so just ignore leg 1 */
255	linkfor *tmpleg;
256	tmpleg = newleg;
257	newleg = newleg2;
258	newleg2 = tmpleg;
259	}
260	}
261	#endif
262	}
263	free(newleg2);
264	newleg2 = (linkfor*)osnew(linkcommon);
265
266	addto_link(newleg, stn2->leg[dirn2]);
267	addto_link(newleg, stn3->leg[dirn3]);
268
269	#if 0
270	printf("Parallel found with stn...stn4 = \n");
271	(dump_node)(stn); (dump_node)(stn2); (dump_node)(stn3); (dump_node)(stn4);
272	printf("dirns = %d %d %d %d\n", dirn, dirn2, dirn3, dirn4);
273	#endif
274	SVX_ASSERT2(stn3->leg[dirn3]->l.to == stn, "stn3 end of \|\| doesn't recip");
275	SVX_ASSERT2(stn4->leg[dirn4]->l.to == stn2, "stn4 end of \|\| doesn't recip");
276	SVX_ASSERT2(stn->leg[(dirn+1)%3]->l.to == stn2 && stn->leg[(dirn + 2) % 3]->l.to == stn2, "\|\| legs aren't");
277
278	/* remove stn and stn2 (already discarded triple parallel) */
279	/* so stn!=stn4 <=> stn2!=stn3 */
280	remove_stn_from_list(&stnlist, stn);
281	remove_stn_from_list(&stnlist, stn2);
282
283	/* stack parallel and replace with a leg between stn3 and stn4 */
284	trav->join[0] = stn3->leg[dirn3];
285	newleg->l.to = stn4;
286	newleg->l.reverse = dirn4 \| FLAG_DATAHERE \| FLAG_REPLACEMENTLEG;
287
288	trav->join[1] = stn4->leg[dirn4];
289	newleg2->l.to = stn3;
290	newleg2->l.reverse = dirn3 \| FLAG_REPLACEMENTLEG;
291
292	stn3->leg[dirn3] = newleg;
293	stn4->leg[dirn4] = newleg2;
294
295	trav->next = reduction_stack;
296	#if PRINT_NETBITS
297	printf("remove parallel\n");
298	#endif
299	reduction_stack = trav;
300	fMore = true;
301	}
302	}
303	}
304
305	if (optimize & BITA('d')) {
306	node stn5, stn6;
307	int dirn5, dirn6;
308	linkfor legAB, legBC, *legCA;
309	#if PRINT_NETBITS
310	printf("replacing deltas with stars\n");
311	#endif
312	FOR_EACH_STN(stn, stnlist) {
313	/* printf("");/
314	/*
315	* : :
316	* * stn5 * stn5
317	* \| \|
318	* * stn2 \|
319	* / \ --> O stnZ
320	* \| \| / \
321	* stn --- stn3 / \
322	* / \ / \
323	* stn4 * * stn6 stn4 * * stn6
324	* : : : :
325	*/
326	if (three_node(stn)) {
327	for (int dirn12 = 0; dirn12 <= 2; dirn12++) {
328	stn2 = stn->leg[dirn12]->l.to;
329	if (stn2 == stn \|\| fixed(stn2)) continue;
330	SVX_ASSERT(three_node(stn2));
331	int dirn13 = (dirn12 + 1) % 3;
332	stn3 = stn->leg[dirn13]->l.to;
333	if (stn3 == stn \|\| stn3 == stn2 \|\| fixed(stn3)) continue;
334	SVX_ASSERT(three_node(stn3));
335	int dirn23 = reverse_leg_dirn(stn->leg[dirn12]);
336	dirn23 = (dirn23 + 1) % 3;
337	if (stn2->leg[dirn23]->l.to != stn3) {
338	dirn23 = (dirn23 + 1) % 3;
339	if (stn2->leg[dirn23]->l.to != stn3) {
340	continue;
341	}
342	}
343	legAB = copy_link(stn->leg[dirn12]);
344	legBC = copy_link(stn2->leg[dirn23]);
345	legCA = copy_link(stn3->leg[reverse_leg_dirn(stn->leg[dirn13])]);
346	dirn = (0 + 1 + 2) - dirn12 - dirn13;
347	dirn2 = (0 + 1 + 2) - dirn23 - reverse_leg_dirn(stn->leg[dirn12]);
348	dirn3 = (0 + 1 + 2) - reverse_leg_dirn(stn->leg[dirn13]) - reverse_leg_dirn(stn2->leg[dirn23]);
349	stn4 = stn->leg[dirn]->l.to;
350	stn5 = stn2->leg[dirn2]->l.to;
351	stn6 = stn3->leg[dirn3]->l.to;
352	if (stn4 == stn2 \|\| stn4 == stn3 \|\| stn5 == stn3) continue;
353	dirn4 = reverse_leg_dirn(stn->leg[dirn]);
354	dirn5 = reverse_leg_dirn(stn2->leg[dirn2]);
355	dirn6 = reverse_leg_dirn(stn3->leg[dirn3]);
356	#if 0
357	printf("delta-star, stn ... stn6 are:\n");
358	(dump_node)(stn);
359	(dump_node)(stn2);
360	(dump_node)(stn3);
361	(dump_node)(stn4);
362	(dump_node)(stn5);
363	(dump_node)(stn6);
364	#endif
365	SVX_ASSERT(stn4->leg[dirn4]->l.to == stn);
366	SVX_ASSERT(stn5->leg[dirn5]->l.to == stn2);
367	SVX_ASSERT(stn6->leg[dirn6]->l.to == stn3);
368
369	trav = allocate_reduction(3);
370	trav->type = TYPE_DELTASTAR;
371	{
372	linkfor legAZ, legBZ, *legCZ;
373	node *stnZ;
374	prefix *nameZ;
375	svar invAB, invBC, invCA, tmp, sum, inv;
376	var vtmp;
377	svar sumAZBZ, sumBZCZ, sumCZAZ;
378	delta temp, temp2;
379
380	/* FIXME: ought to handle cases when some legs are
381	* equates, but handle as a special case maybe? */
382	if (!invert_svar(&invAB, &legAB->v)) break;
383	if (!invert_svar(&invBC, &legBC->v)) break;
384	if (!invert_svar(&invCA, &legCA->v)) break;
385
386	addss(&sum, &legBC->v, &legCA->v);
387	addss(&tmp, &sum, &legAB->v);
388	if (!invert_svar(&inv, &tmp)) {
389	/* impossible - loop of zero variance */
390	BUG("loop of zero variance found");
391	}
392
393	legAZ = osnew(linkfor);
394	legBZ = osnew(linkfor);
395	legCZ = osnew(linkfor);
396
397	/* AZBZ */
398	/* done above: addvv(&sum, &legBC->v, &legCA->v); */
399	mulss(&vtmp, &sum, &inv);
400	smulvs(&sumAZBZ, &vtmp, &legAB->v);
401
402	adddd(&temp, &legBC->d, &legCA->d);
403	divds(&temp2, &temp, &sum);
404	mulsd(&temp, &invAB, &legAB->d);
405	subdd(&temp, &temp2, &temp);
406	mulsd(&legBZ->d, &sumAZBZ, &temp);
407
408	/* leg vectors after transform are determined up to
409	* a constant addition, so arbitrarily fix AZ = 0 */
410	legAZ->d[2] = legAZ->d[1] = legAZ->d[0] = 0;
411
412	/* BZCZ */
413	addss(&sum, &legCA->v, &legAB->v);
414	mulss(&vtmp, &sum, &inv);
415	smulvs(&sumBZCZ, &vtmp, &legBC->v);
416
417	/* CZAZ */
418	addss(&sum, &legAB->v, &legBC->v);
419	mulss(&vtmp, &sum, &inv);
420	smulvs(&sumCZAZ, &vtmp, &legCA->v);
421
422	adddd(&temp, &legAB->d, &legBC->d);
423	divds(&temp2, &temp, &sum);
424	mulsd(&temp, &invCA, &legCA->d);
425	/* NB: swapped arguments to negate answer for legCZ->d */
426	subdd(&temp, &temp, &temp2);
427	mulsd(&legCZ->d, &sumCZAZ, &temp);
428
429	free(legAB);
430	free(legBC);
431	free(legCA);
432
433	/* Now add two, subtract third, and scale by 0.5 */
434	addss(&sum, &sumAZBZ, &sumCZAZ);
435	subss(&sum, &sum, &sumBZCZ);
436	mulsc(&legAZ->v, &sum, 0.5);
437
438	addss(&sum, &sumBZCZ, &sumAZBZ);
439	subss(&sum, &sum, &sumCZAZ);
440	mulsc(&legBZ->v, &sum, 0.5);
441
442	addss(&sum, &sumCZAZ, &sumBZCZ);
443	subss(&sum, &sum, &sumAZBZ);
444	mulsc(&legCZ->v, &sum, 0.5);
445
446	nameZ = osnew(prefix);
447	nameZ->pos = osnew(pos);
448	nameZ->ident.p = NULL;
449	stnZ = osnew(node);
450	stnZ->name = nameZ;
451	nameZ->stn = stnZ;
452	nameZ->up = NULL;
453	nameZ->min_export = nameZ->max_export = 0;
454	nameZ->sflags = 0;
455	unfix(stnZ);
456	add_stn_to_list(&stnlist, stnZ);
457	legAZ->l.to = stnZ;
458	legAZ->l.reverse = 0 \| FLAG_DATAHERE \| FLAG_REPLACEMENTLEG;
459	legBZ->l.to = stnZ;
460	legBZ->l.reverse = 1 \| FLAG_DATAHERE \| FLAG_REPLACEMENTLEG;
461	legCZ->l.to = stnZ;
462	legCZ->l.reverse = 2 \| FLAG_DATAHERE \| FLAG_REPLACEMENTLEG;
463	stnZ->leg[0] = (linkfor*)osnew(linkcommon);
464	stnZ->leg[1] = (linkfor*)osnew(linkcommon);
465	stnZ->leg[2] = (linkfor*)osnew(linkcommon);
466	stnZ->leg[0]->l.to = stn4;
467	stnZ->leg[0]->l.reverse = dirn4;
468	stnZ->leg[1]->l.to = stn5;
469	stnZ->leg[1]->l.reverse = dirn5;
470	stnZ->leg[2]->l.to = stn6;
471	stnZ->leg[2]->l.reverse = dirn6;
472	addto_link(legAZ, stn4->leg[dirn4]);
473	addto_link(legBZ, stn5->leg[dirn5]);
474	addto_link(legCZ, stn6->leg[dirn6]);
475	/* stack stuff */
476	trav->join[0] = stn4->leg[dirn4];
477	trav->join[1] = stn5->leg[dirn5];
478	trav->join[2] = stn6->leg[dirn6];
479	trav->next = reduction_stack;
480	#if PRINT_NETBITS
481	printf("remove delta*\n");
482	#endif
483	reduction_stack = trav;
484	fMore = true;
485
486	remove_stn_from_list(&stnlist, stn);
487	remove_stn_from_list(&stnlist, stn2);
488	remove_stn_from_list(&stnlist, stn3);
489	stn4->leg[dirn4] = legAZ;
490	stn5->leg[dirn5] = legBZ;
491	stn6->leg[dirn6] = legCZ;
492	}
493	break;
494	}
495	}
496	}
497	}
498
499	}
500	}
501
502	extern void
503	replace_subnets(void)
504	{
505	node stn2, stn3, *stn4;
506	int dirn2, dirn3, dirn4;
507
508	/* help to catch bad accesses */
509	stn2 = stn3 = stn4 = NULL;
510	dirn2 = dirn3 = dirn4 = 0;
511
512	out_current_action(msg(/Calculating network/130));
513
514	while (reduction_stack != NULL) {
515	/* printf("replace_subnets() type %d\n", reduction_stack->type);*/
516
517	#if PRINT_NETBITS
518	printf("replace_subnets\n");
519	if (reduction_stack->type == TYPE_LOLLIPOP) printf("islollipop\n");
520	if (reduction_stack->type == TYPE_PARALLEL) printf("isparallel\n");
521	if (reduction_stack->type == TYPE_DELTASTAR) printf("isdelta*\n");
522	#endif
523
524	if (reduction_stack->type != TYPE_DELTASTAR) {
525	linkfor *leg;
526	leg = reduction_stack->join[0]; leg = reverse_leg(leg);
527	stn3 = leg->l.to; dirn3 = reverse_leg_dirn(leg);
528	leg = reduction_stack->join[1]; leg = reverse_leg(leg);
529	stn4 = leg->l.to; dirn4 = reverse_leg_dirn(leg);
530
531	if (!fixed(stn3) \|\| !fixed(stn4)) {
532	SVX_ASSERT(!fixed(stn3) && !fixed(stn4));
533	goto skip;
534	}
535	SVX_ASSERT(data_here(stn3->leg[dirn3]));
536	}
537
538	if (reduction_stack->type == TYPE_LOLLIPOP) {
539	node *stn;
540	delta e;
541	linkfor *leg;
542	int zero;
543
544	leg = stn3->leg[dirn3];
545	stn2 = reduction_stack->join[0]->l.to;
546	dirn2 = reverse_leg_dirn(reduction_stack->join[0]);
547
548	zero = fZeros(&leg->v);
549	if (!zero) {
550	delta tmp;
551	subdd(&e, &POSD(stn4), &POSD(stn3));
552	subdd(&tmp, &e, &leg->d);
553	divds(&e, &tmp, &leg->v);
554	}
555	if (data_here(reduction_stack->join[0])) {
556	adddd(&POSD(stn2), &POSD(stn3), &reduction_stack->join[0]->d);
557	if (!zero) {
558	delta tmp;
559	mulsd(&tmp, &reduction_stack->join[0]->v, &e);
560	adddd(&POSD(stn2), &POSD(stn2), &tmp);
561	}
562	} else {
563	subdd(&POSD(stn2), &POSD(stn3), &stn2->leg[dirn2]->d);
564	if (!zero) {
565	delta tmp;
566	mulsd(&tmp, &stn2->leg[dirn2]->v, &e);
567	adddd(&POSD(stn2), &POSD(stn2), &tmp);
568	}
569	}
570	dirn2 = (dirn2 + 2) % 3; /* point back at stn again */
571	stn = stn2->leg[dirn2]->l.to;
572	#if 0
573	printf("Replacing lollipop with stn...stn4 = \n");
574	print_prefix(stn->name); putnl();
575	print_prefix(stn2->name); putnl();
576	print_prefix(stn3->name); putnl();
577	print_prefix(stn4->name); putnl();
578	#endif
579	if (data_here(stn2->leg[dirn2]))
580	adddd(&POSD(stn), &POSD(stn2), &stn2->leg[dirn2]->d);
581	else
582	subdd(&POSD(stn), &POSD(stn2), &reverse_leg(stn2->leg[dirn2])->d);
583
584	/* The "stick" of the lollipop is a new articulation. */
585	stn2->leg[dirn2]->l.reverse \|= FLAG_ARTICULATION;
586	reverse_leg(stn2->leg[dirn2])->l.reverse \|= FLAG_ARTICULATION;
587
588	add_stn_to_list(&fixedlist, stn);
589	add_stn_to_list(&fixedlist, stn2);
590
591	free(stn3->leg[dirn3]);
592	stn3->leg[dirn3] = reduction_stack->join[0];
593	free(stn4->leg[dirn4]);
594	stn4->leg[dirn4] = reduction_stack->join[1];
595	} else if (reduction_stack->type == TYPE_PARALLEL) {
596	/* parallel legs */
597	node *stn;
598	delta e, e2;
599	linkfor *leg;
600	int dirn;
601
602	stn = reduction_stack->join[0]->l.to;
603	stn2 = reduction_stack->join[1]->l.to;
604
605	dirn = reverse_leg_dirn(reduction_stack->join[0]);
606	dirn2 = reverse_leg_dirn(reduction_stack->join[1]);
607
608	leg = stn3->leg[dirn3];
609
610	if (leg->l.reverse & FLAG_ARTICULATION) {
611	reduction_stack->join[0]->l.reverse \|= FLAG_ARTICULATION;
612	stn->leg[dirn]->l.reverse \|= FLAG_ARTICULATION;
613	reduction_stack->join[1]->l.reverse \|= FLAG_ARTICULATION;
614	stn2->leg[dirn2]->l.reverse \|= FLAG_ARTICULATION;
615	}
616
617	if (fZeros(&leg->v))
618	e[0] = e[1] = e[2] = 0.0;
619	else {
620	delta tmp;
621	subdd(&e, &POSD(stn4), &POSD(stn3));
622	subdd(&tmp, &e, &leg->d);
623	divds(&e, &tmp, &leg->v);
624	}
625
626	if (data_here(reduction_stack->join[0])) {
627	leg = reduction_stack->join[0];
628	adddd(&POSD(stn), &POSD(stn3), &leg->d);
629	} else {
630	leg = stn->leg[dirn];
631	subdd(&POSD(stn), &POSD(stn3), &leg->d);
632	}
633	mulsd(&e2, &leg->v, &e);
634	adddd(&POSD(stn), &POSD(stn), &e2);
635
636	if (data_here(reduction_stack->join[1])) {
637	leg = reduction_stack->join[1];
638	adddd(&POSD(stn2), &POSD(stn4), &leg->d);
639	} else {
640	leg = stn2->leg[dirn2];
641	subdd(&POSD(stn2), &POSD(stn4), &leg->d);
642	}
643	mulsd(&e2, &leg->v, &e);
644	subdd(&POSD(stn2), &POSD(stn2), &e2);
645	#if 0
646	printf("Replacing parallel with stn...stn4 = \n");
647	print_prefix(stn->name); putnl();
648	print_prefix(stn2->name); putnl();
649	print_prefix(stn3->name); putnl();
650	print_prefix(stn4->name); putnl();
651	#endif
652
653	add_stn_to_list(&fixedlist, stn);
654	add_stn_to_list(&fixedlist, stn2);
655
656	free(stn3->leg[dirn3]);
657	stn3->leg[dirn3] = reduction_stack->join[0];
658	free(stn4->leg[dirn4]);
659	stn4->leg[dirn4] = reduction_stack->join[1];
660	} else if (reduction_stack->type == TYPE_DELTASTAR) {
661	node *stnZ;
662	node *stn[3];
663	int dirn[3];
664	int i;
665	linkfor *leg;
666
667	/* work out ends as we don't bother stacking them */
668	leg = reverse_leg(reduction_stack->join[0]);
669	stn[0] = leg->l.to;
670	dirn[0] = reverse_leg_dirn(leg);
671	stnZ = stn[0]->leg[dirn[0]]->l.to;
672	if (!fixed(stnZ)) {
673	SVX_ASSERT(!fixed(stn[0]));
674	goto skip;
675	}
676	stn[1] = stnZ->leg[1]->l.to;
677	dirn[1] = reverse_leg_dirn(stnZ->leg[1]);
678	stn[2] = stnZ->leg[2]->l.to;
679	dirn[2] = reverse_leg_dirn(stnZ->leg[2]);
680	/print_prefix(stnZ->name);printf(" %p\n",(void)stnZ);*/
681
682	for (i = 0; i < 3; i++) {
683	SVX_ASSERT2(fixed(stn[i]), "stn not fixed for D*");
684
685	leg = stn[i]->leg[dirn[i]];
686
687	SVX_ASSERT2(data_here(leg), "data not on leg for D*");
688	SVX_ASSERT2(leg->l.to == stnZ, "bad sub-network for D*");
689
690	stn2 = reduction_stack->join[i]->l.to;
691
692	if (data_here(reduction_stack->join[i])) {
693	adddd(&POSD(stn2), &POSD(stn[i]), &reduction_stack->join[i]->d);
694	} else {
695	subdd(&POSD(stn2), &POSD(stn[i]), &reverse_leg(reduction_stack->join[i])->d);
696	}
697
698	if (!fZeros(&leg->v)) {
699	delta e, tmp;
700	subdd(&e, &POSD(stnZ), &POSD(stn[i]));
701	subdd(&e, &e, &leg->d);
702	divds(&tmp, &e, &leg->v);
703	if (data_here(reduction_stack->join[i])) {
704	mulsd(&e, &reduction_stack->join[i]->v, &tmp);
705	} else {
706	mulsd(&e, &reverse_leg(reduction_stack->join[i])->v, &tmp);
707	}
708	adddd(&POSD(stn2), &POSD(stn2), &e);
709	}
710	add_stn_to_list(&fixedlist, stn2);
711	free(leg);
712	stn[i]->leg[dirn[i]] = reduction_stack->join[i];
713	/* transfer the articulation status of the radial legs */
714	if (stnZ->leg[i]->l.reverse & FLAG_ARTICULATION) {
715	reduction_stack->join[i]->l.reverse \|= FLAG_ARTICULATION;
716	reverse_leg(reduction_stack->join[i])->l.reverse \|= FLAG_ARTICULATION;
717	}
718	free(stnZ->leg[i]);
719	stnZ->leg[i] = NULL;
720	}
721	/printf("---%f %f %f\n",POS(stnZ, 0), POS(stnZ, 1), POS(stnZ, 2));/
722	remove_stn_from_list(&fixedlist, stnZ);
723	free(stnZ->name);
724	free(stnZ);
725	} else {
726	BUG("reduction_stack has unknown type");
727	}
728
729	skip:;
730	reduction *ptrOld = reduction_stack;
731	reduction_stack = reduction_stack->next;
732	free(ptrOld);
733	}
734	}

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: git/src/network.c @ b826dcf

Download in other formats: