Context Navigation

network.c @ c981e5b

Visit:

main

Last change on this file since c981e5b was 8048d405, checked in by Olly Betts <olly@…>, 5 weeks ago

Separate linkcommon reverse and flag bits

We had at least 2 spare bytes of padding, so we can store the reverse
direction separately which seems cleaner (and slightly more efficient
but the difference is so small it's only measureable using cachegrind).

Property mode set to 100644

File size: 21.8 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.bits \|= FLAG_ARTICULATION;
99	rev_leg->l.bits \|= 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;
156	newleg->l.bits = FLAG_DATAHERE \| FLAG_REPLACEMENTLEG;
157
158	trav->join[1] = stn4->leg[dirn4];
159	newleg2->l.to = stn3;
160	newleg2->l.reverse = dirn3;
161	newleg2->l.bits = FLAG_REPLACEMENTLEG;
162
163	stn3->leg[dirn3] = newleg;
164	stn4->leg[dirn4] = newleg2;
165
166	trav->next = reduction_stack;
167	#if PRINT_NETBITS
168	printf("remove lollipop\n");
169	#endif
170	reduction_stack = trav;
171	fMore = true;
172	}
173	}
174	}
175
176	if (optimize & BITA('p')) {
177	#if PRINT_NETBITS
178	printf("replacing parallel legs\n");
179	#endif
180	FOR_EACH_STN(stn, stnlist) {
181	/*
182	* : :
183	* * stn3 * stn3
184	* \| \|
185	* * stn \|
186	* ( ) --> \|
187	* * stn2 \|
188	* \| \|
189	* * stn4 * stn4
190	* : :
191	*/
192	if (three_node(stn)) {
193	stn2 = stn->leg[0]->l.to;
194	if (stn2 == stn->leg[1]->l.to) {
195	dirn = 2;
196	} else if (stn2 == stn->leg[2]->l.to) {
197	dirn = 1;
198	} else {
199	if (stn->leg[1]->l.to != stn->leg[2]->l.to) continue;
200	stn2 = stn->leg[1]->l.to;
201	dirn = 0;
202	}
203
204	/* stn == stn2 => lollipop */
205	if (stn == stn2 \|\| fixed(stn2)) continue;
206
207	SVX_ASSERT(three_node(stn2));
208
209	stn3 = stn->leg[dirn]->l.to;
210	/* 3 parallel legs (=> nothing else) so leave */
211	if (stn3 == stn2) continue;
212
213	dirn3 = reverse_leg_dirn(stn->leg[dirn]);
214	dirn2 = (0 + 1 + 2 - reverse_leg_dirn(stn->leg[(dirn + 1) % 3])
215	- reverse_leg_dirn(stn->leg[(dirn + 2) % 3]));
216
217	stn4 = stn2->leg[dirn2]->l.to;
218	dirn4 = reverse_leg_dirn(stn2->leg[dirn2]);
219
220	trav = allocate_reduction(2);
221	trav->type = TYPE_PARALLEL;
222
223	newleg = copy_link(stn->leg[(dirn + 1) % 3]);
224	/* use newleg2 for scratch */
225	newleg2 = copy_link(stn->leg[(dirn + 2) % 3]);
226	{
227	#ifdef NO_COVARIANCES
228	vars sum;
229	var prod;
230	delta temp, temp2;
231	addss(&sum, &newleg->v, &newleg2->v);
232	SVX_ASSERT2(!fZeros(&sum), "loop of zero variance found");
233	mulss(&prod, &newleg->v, &newleg2->v);
234	mulsd(&temp, &newleg2->v, &newleg->d);
235	mulsd(&temp2, &newleg->v, &newleg2->d);
236	adddd(&temp, &temp, &temp2);
237	divds(&newleg->d, &temp, &sum);
238	sdivvs(&newleg->v, &prod, &sum);
239	#else
240	svar inv1, inv2, sum;
241	delta temp, temp2;
242	/* if leg one is an equate, we can just ignore leg two
243	* whatever it is */
244	if (invert_svar(&inv1, &newleg->v)) {
245	if (invert_svar(&inv2, &newleg2->v)) {
246	addss(&sum, &inv1, &inv2);
247	if (!invert_svar(&newleg->v, &sum)) {
248	BUG("matrix singular in parallel legs replacement");
249	}
250
251	mulsd(&temp, &inv1, &newleg->d);
252	mulsd(&temp2, &inv2, &newleg2->d);
253	adddd(&temp, &temp, &temp2);
254	mulsd(&newleg->d, &newleg->v, &temp);
255	} else {
256	/* leg two is an equate, so just ignore leg 1 */
257	linkfor *tmpleg;
258	tmpleg = newleg;
259	newleg = newleg2;
260	newleg2 = tmpleg;
261	}
262	}
263	#endif
264	}
265	free(newleg2);
266	newleg2 = (linkfor*)osnew(linkcommon);
267
268	addto_link(newleg, stn2->leg[dirn2]);
269	addto_link(newleg, stn3->leg[dirn3]);
270
271	#if 0
272	printf("Parallel found with stn...stn4 = \n");
273	(dump_node)(stn); (dump_node)(stn2); (dump_node)(stn3); (dump_node)(stn4);
274	printf("dirns = %d %d %d %d\n", dirn, dirn2, dirn3, dirn4);
275	#endif
276	SVX_ASSERT2(stn3->leg[dirn3]->l.to == stn, "stn3 end of \|\| doesn't recip");
277	SVX_ASSERT2(stn4->leg[dirn4]->l.to == stn2, "stn4 end of \|\| doesn't recip");
278	SVX_ASSERT2(stn->leg[(dirn+1)%3]->l.to == stn2 && stn->leg[(dirn + 2) % 3]->l.to == stn2, "\|\| legs aren't");
279
280	/* remove stn and stn2 (already discarded triple parallel) */
281	/* so stn!=stn4 <=> stn2!=stn3 */
282	remove_stn_from_list(&stnlist, stn);
283	remove_stn_from_list(&stnlist, stn2);
284
285	/* stack parallel and replace with a leg between stn3 and stn4 */
286	trav->join[0] = stn3->leg[dirn3];
287	newleg->l.to = stn4;
288	newleg->l.reverse = dirn4;
289	newleg->l.bits = FLAG_DATAHERE \| FLAG_REPLACEMENTLEG;
290
291	trav->join[1] = stn4->leg[dirn4];
292	newleg2->l.to = stn3;
293	newleg2->l.reverse = dirn3;
294	newleg2->l.bits = FLAG_REPLACEMENTLEG;
295
296	stn3->leg[dirn3] = newleg;
297	stn4->leg[dirn4] = newleg2;
298
299	trav->next = reduction_stack;
300	#if PRINT_NETBITS
301	printf("remove parallel\n");
302	#endif
303	reduction_stack = trav;
304	fMore = true;
305	}
306	}
307	}
308
309	if (optimize & BITA('d')) {
310	node stn5, stn6;
311	int dirn5, dirn6;
312	linkfor legAB, legBC, *legCA;
313	#if PRINT_NETBITS
314	printf("replacing deltas with stars\n");
315	#endif
316	FOR_EACH_STN(stn, stnlist) {
317	/* printf("");/
318	/*
319	* : :
320	* * stn5 * stn5
321	* \| \|
322	* * stn2 \|
323	* / \ --> O stnZ
324	* \| \| / \
325	* stn --- stn3 / \
326	* / \ / \
327	* stn4 * * stn6 stn4 * * stn6
328	* : : : :
329	*/
330	if (three_node(stn)) {
331	for (int dirn12 = 0; dirn12 <= 2; dirn12++) {
332	stn2 = stn->leg[dirn12]->l.to;
333	if (stn2 == stn \|\| fixed(stn2)) continue;
334	SVX_ASSERT(three_node(stn2));
335	int dirn13 = (dirn12 + 1) % 3;
336	stn3 = stn->leg[dirn13]->l.to;
337	if (stn3 == stn \|\| stn3 == stn2 \|\| fixed(stn3)) continue;
338	SVX_ASSERT(three_node(stn3));
339	int dirn23 = reverse_leg_dirn(stn->leg[dirn12]);
340	dirn23 = (dirn23 + 1) % 3;
341	if (stn2->leg[dirn23]->l.to != stn3) {
342	dirn23 = (dirn23 + 1) % 3;
343	if (stn2->leg[dirn23]->l.to != stn3) {
344	continue;
345	}
346	}
347	legAB = copy_link(stn->leg[dirn12]);
348	legBC = copy_link(stn2->leg[dirn23]);
349	legCA = copy_link(stn3->leg[reverse_leg_dirn(stn->leg[dirn13])]);
350	dirn = (0 + 1 + 2) - dirn12 - dirn13;
351	dirn2 = (0 + 1 + 2) - dirn23 - reverse_leg_dirn(stn->leg[dirn12]);
352	dirn3 = (0 + 1 + 2) - reverse_leg_dirn(stn->leg[dirn13]) - reverse_leg_dirn(stn2->leg[dirn23]);
353	stn4 = stn->leg[dirn]->l.to;
354	stn5 = stn2->leg[dirn2]->l.to;
355	stn6 = stn3->leg[dirn3]->l.to;
356	if (stn4 == stn2 \|\| stn4 == stn3 \|\| stn5 == stn3) continue;
357	dirn4 = reverse_leg_dirn(stn->leg[dirn]);
358	dirn5 = reverse_leg_dirn(stn2->leg[dirn2]);
359	dirn6 = reverse_leg_dirn(stn3->leg[dirn3]);
360	#if 0
361	printf("delta-star, stn ... stn6 are:\n");
362	(dump_node)(stn);
363	(dump_node)(stn2);
364	(dump_node)(stn3);
365	(dump_node)(stn4);
366	(dump_node)(stn5);
367	(dump_node)(stn6);
368	#endif
369	SVX_ASSERT(stn4->leg[dirn4]->l.to == stn);
370	SVX_ASSERT(stn5->leg[dirn5]->l.to == stn2);
371	SVX_ASSERT(stn6->leg[dirn6]->l.to == stn3);
372
373	trav = allocate_reduction(3);
374	trav->type = TYPE_DELTASTAR;
375	{
376	linkfor legAZ, legBZ, *legCZ;
377	node *stnZ;
378	prefix *nameZ;
379	svar invAB, invBC, invCA, tmp, sum, inv;
380	var vtmp;
381	svar sumAZBZ, sumBZCZ, sumCZAZ;
382	delta temp, temp2;
383
384	/* FIXME: ought to handle cases when some legs are
385	* equates, but handle as a special case maybe? */
386	if (!invert_svar(&invAB, &legAB->v)) break;
387	if (!invert_svar(&invBC, &legBC->v)) break;
388	if (!invert_svar(&invCA, &legCA->v)) break;
389
390	addss(&sum, &legBC->v, &legCA->v);
391	addss(&tmp, &sum, &legAB->v);
392	if (!invert_svar(&inv, &tmp)) {
393	/* impossible - loop of zero variance */
394	BUG("loop of zero variance found");
395	}
396
397	legAZ = osnew(linkfor);
398	legBZ = osnew(linkfor);
399	legCZ = osnew(linkfor);
400
401	/* AZBZ */
402	/* done above: addvv(&sum, &legBC->v, &legCA->v); */
403	mulss(&vtmp, &sum, &inv);
404	smulvs(&sumAZBZ, &vtmp, &legAB->v);
405
406	adddd(&temp, &legBC->d, &legCA->d);
407	divds(&temp2, &temp, &sum);
408	mulsd(&temp, &invAB, &legAB->d);
409	subdd(&temp, &temp2, &temp);
410	mulsd(&legBZ->d, &sumAZBZ, &temp);
411
412	/* leg vectors after transform are determined up to
413	* a constant addition, so arbitrarily fix AZ = 0 */
414	legAZ->d[2] = legAZ->d[1] = legAZ->d[0] = 0;
415
416	/* BZCZ */
417	addss(&sum, &legCA->v, &legAB->v);
418	mulss(&vtmp, &sum, &inv);
419	smulvs(&sumBZCZ, &vtmp, &legBC->v);
420
421	/* CZAZ */
422	addss(&sum, &legAB->v, &legBC->v);
423	mulss(&vtmp, &sum, &inv);
424	smulvs(&sumCZAZ, &vtmp, &legCA->v);
425
426	adddd(&temp, &legAB->d, &legBC->d);
427	divds(&temp2, &temp, &sum);
428	mulsd(&temp, &invCA, &legCA->d);
429	/* NB: swapped arguments to negate answer for legCZ->d */
430	subdd(&temp, &temp, &temp2);
431	mulsd(&legCZ->d, &sumCZAZ, &temp);
432
433	free(legAB);
434	free(legBC);
435	free(legCA);
436
437	/* Now add two, subtract third, and scale by 0.5 */
438	addss(&sum, &sumAZBZ, &sumCZAZ);
439	subss(&sum, &sum, &sumBZCZ);
440	mulsc(&legAZ->v, &sum, 0.5);
441
442	addss(&sum, &sumBZCZ, &sumAZBZ);
443	subss(&sum, &sum, &sumCZAZ);
444	mulsc(&legBZ->v, &sum, 0.5);
445
446	addss(&sum, &sumCZAZ, &sumBZCZ);
447	subss(&sum, &sum, &sumAZBZ);
448	mulsc(&legCZ->v, &sum, 0.5);
449
450	nameZ = osnew(prefix);
451	nameZ->pos = osnew(pos);
452	nameZ->ident.p = NULL;
453	stnZ = osnew(node);
454	stnZ->name = nameZ;
455	nameZ->stn = stnZ;
456	nameZ->up = NULL;
457	nameZ->filename = NULL;
458	nameZ->min_export = nameZ->max_export = 0;
459	nameZ->sflags = 0;
460	unfix(stnZ);
461	add_stn_to_list(&stnlist, stnZ);
462	legAZ->l.to = stnZ;
463	legAZ->l.reverse = 0;
464	legAZ->l.bits = FLAG_DATAHERE \| FLAG_REPLACEMENTLEG;
465	legBZ->l.to = stnZ;
466	legBZ->l.reverse = 1;
467	legBZ->l.bits = FLAG_DATAHERE \| FLAG_REPLACEMENTLEG;
468	legCZ->l.to = stnZ;
469	legCZ->l.reverse = 2;
470	legCZ->l.bits = FLAG_DATAHERE \| FLAG_REPLACEMENTLEG;
471	stnZ->leg[0] = (linkfor*)osnew(linkcommon);
472	stnZ->leg[1] = (linkfor*)osnew(linkcommon);
473	stnZ->leg[2] = (linkfor*)osnew(linkcommon);
474	stnZ->leg[0]->l.to = stn4;
475	stnZ->leg[0]->l.reverse = dirn4;
476	stnZ->leg[0]->l.bits = 0;
477	stnZ->leg[1]->l.to = stn5;
478	stnZ->leg[1]->l.reverse = dirn5;
479	stnZ->leg[1]->l.bits = 0;
480	stnZ->leg[2]->l.to = stn6;
481	stnZ->leg[2]->l.reverse = dirn6;
482	stnZ->leg[2]->l.bits = 0;
483	addto_link(legAZ, stn4->leg[dirn4]);
484	addto_link(legBZ, stn5->leg[dirn5]);
485	addto_link(legCZ, stn6->leg[dirn6]);
486	/* stack stuff */
487	trav->join[0] = stn4->leg[dirn4];
488	trav->join[1] = stn5->leg[dirn5];
489	trav->join[2] = stn6->leg[dirn6];
490	trav->next = reduction_stack;
491	#if PRINT_NETBITS
492	printf("remove delta*\n");
493	#endif
494	reduction_stack = trav;
495	fMore = true;
496
497	remove_stn_from_list(&stnlist, stn);
498	remove_stn_from_list(&stnlist, stn2);
499	remove_stn_from_list(&stnlist, stn3);
500	stn4->leg[dirn4] = legAZ;
501	stn5->leg[dirn5] = legBZ;
502	stn6->leg[dirn6] = legCZ;
503	}
504	break;
505	}
506	}
507	}
508	}
509
510	}
511	}
512
513	extern void
514	replace_subnets(void)
515	{
516	node stn2, stn3, *stn4;
517	int dirn2, dirn3, dirn4;
518
519	/* help to catch bad accesses */
520	stn2 = stn3 = stn4 = NULL;
521	dirn2 = dirn3 = dirn4 = 0;
522
523	out_current_action(msg(/Calculating network/130));
524
525	while (reduction_stack != NULL) {
526	/* printf("replace_subnets() type %d\n", reduction_stack->type);*/
527
528	#if PRINT_NETBITS
529	printf("replace_subnets\n");
530	if (reduction_stack->type == TYPE_LOLLIPOP) printf("islollipop\n");
531	if (reduction_stack->type == TYPE_PARALLEL) printf("isparallel\n");
532	if (reduction_stack->type == TYPE_DELTASTAR) printf("isdelta*\n");
533	#endif
534
535	if (reduction_stack->type != TYPE_DELTASTAR) {
536	linkfor *leg;
537	leg = reduction_stack->join[0]; leg = reverse_leg(leg);
538	stn3 = leg->l.to; dirn3 = reverse_leg_dirn(leg);
539	leg = reduction_stack->join[1]; leg = reverse_leg(leg);
540	stn4 = leg->l.to; dirn4 = reverse_leg_dirn(leg);
541
542	if (!fixed(stn3) \|\| !fixed(stn4)) {
543	SVX_ASSERT(!fixed(stn3) && !fixed(stn4));
544	goto skip;
545	}
546	SVX_ASSERT(data_here(stn3->leg[dirn3]));
547	}
548
549	if (reduction_stack->type == TYPE_LOLLIPOP) {
550	node *stn;
551	delta e;
552	linkfor *leg;
553	int zero;
554
555	leg = stn3->leg[dirn3];
556	stn2 = reduction_stack->join[0]->l.to;
557	dirn2 = reverse_leg_dirn(reduction_stack->join[0]);
558
559	zero = fZeros(&leg->v);
560	if (!zero) {
561	delta tmp;
562	subdd(&e, &POSD(stn4), &POSD(stn3));
563	subdd(&tmp, &e, &leg->d);
564	divds(&e, &tmp, &leg->v);
565	}
566	if (data_here(reduction_stack->join[0])) {
567	adddd(&POSD(stn2), &POSD(stn3), &reduction_stack->join[0]->d);
568	if (!zero) {
569	delta tmp;
570	mulsd(&tmp, &reduction_stack->join[0]->v, &e);
571	adddd(&POSD(stn2), &POSD(stn2), &tmp);
572	}
573	} else {
574	subdd(&POSD(stn2), &POSD(stn3), &stn2->leg[dirn2]->d);
575	if (!zero) {
576	delta tmp;
577	mulsd(&tmp, &stn2->leg[dirn2]->v, &e);
578	adddd(&POSD(stn2), &POSD(stn2), &tmp);
579	}
580	}
581	dirn2 = (dirn2 + 2) % 3; /* point back at stn again */
582	stn = stn2->leg[dirn2]->l.to;
583	#if 0
584	printf("Replacing lollipop with stn...stn4 = \n");
585	print_prefix(stn->name); putnl();
586	print_prefix(stn2->name); putnl();
587	print_prefix(stn3->name); putnl();
588	print_prefix(stn4->name); putnl();
589	#endif
590	if (data_here(stn2->leg[dirn2]))
591	adddd(&POSD(stn), &POSD(stn2), &stn2->leg[dirn2]->d);
592	else
593	subdd(&POSD(stn), &POSD(stn2), &reverse_leg(stn2->leg[dirn2])->d);
594
595	/* The "stick" of the lollipop is a new articulation. */
596	stn2->leg[dirn2]->l.bits \|= FLAG_ARTICULATION;
597	reverse_leg(stn2->leg[dirn2])->l.bits \|= FLAG_ARTICULATION;
598
599	add_stn_to_list(&fixedlist, stn);
600	add_stn_to_list(&fixedlist, stn2);
601
602	free(stn3->leg[dirn3]);
603	stn3->leg[dirn3] = reduction_stack->join[0];
604	free(stn4->leg[dirn4]);
605	stn4->leg[dirn4] = reduction_stack->join[1];
606	} else if (reduction_stack->type == TYPE_PARALLEL) {
607	/* parallel legs */
608	node *stn;
609	delta e, e2;
610	linkfor *leg;
611	int dirn;
612
613	stn = reduction_stack->join[0]->l.to;
614	stn2 = reduction_stack->join[1]->l.to;
615
616	dirn = reverse_leg_dirn(reduction_stack->join[0]);
617	dirn2 = reverse_leg_dirn(reduction_stack->join[1]);
618
619	leg = stn3->leg[dirn3];
620
621	if (leg->l.bits & FLAG_ARTICULATION) {
622	reduction_stack->join[0]->l.bits \|= FLAG_ARTICULATION;
623	stn->leg[dirn]->l.bits \|= FLAG_ARTICULATION;
624	reduction_stack->join[1]->l.bits \|= FLAG_ARTICULATION;
625	stn2->leg[dirn2]->l.bits \|= FLAG_ARTICULATION;
626	}
627
628	if (fZeros(&leg->v))
629	e[0] = e[1] = e[2] = 0.0;
630	else {
631	delta tmp;
632	subdd(&e, &POSD(stn4), &POSD(stn3));
633	subdd(&tmp, &e, &leg->d);
634	divds(&e, &tmp, &leg->v);
635	}
636
637	if (data_here(reduction_stack->join[0])) {
638	leg = reduction_stack->join[0];
639	adddd(&POSD(stn), &POSD(stn3), &leg->d);
640	} else {
641	leg = stn->leg[dirn];
642	subdd(&POSD(stn), &POSD(stn3), &leg->d);
643	}
644	mulsd(&e2, &leg->v, &e);
645	adddd(&POSD(stn), &POSD(stn), &e2);
646
647	if (data_here(reduction_stack->join[1])) {
648	leg = reduction_stack->join[1];
649	adddd(&POSD(stn2), &POSD(stn4), &leg->d);
650	} else {
651	leg = stn2->leg[dirn2];
652	subdd(&POSD(stn2), &POSD(stn4), &leg->d);
653	}
654	mulsd(&e2, &leg->v, &e);
655	subdd(&POSD(stn2), &POSD(stn2), &e2);
656	#if 0
657	printf("Replacing parallel with stn...stn4 = \n");
658	print_prefix(stn->name); putnl();
659	print_prefix(stn2->name); putnl();
660	print_prefix(stn3->name); putnl();
661	print_prefix(stn4->name); putnl();
662	#endif
663
664	add_stn_to_list(&fixedlist, stn);
665	add_stn_to_list(&fixedlist, stn2);
666
667	free(stn3->leg[dirn3]);
668	stn3->leg[dirn3] = reduction_stack->join[0];
669	free(stn4->leg[dirn4]);
670	stn4->leg[dirn4] = reduction_stack->join[1];
671	} else if (reduction_stack->type == TYPE_DELTASTAR) {
672	node *stnZ;
673	node *stn[3];
674	int dirn[3];
675	int i;
676	linkfor *leg;
677
678	/* work out ends as we don't bother stacking them */
679	leg = reverse_leg(reduction_stack->join[0]);
680	stn[0] = leg->l.to;
681	dirn[0] = reverse_leg_dirn(leg);
682	stnZ = stn[0]->leg[dirn[0]]->l.to;
683	if (!fixed(stnZ)) {
684	SVX_ASSERT(!fixed(stn[0]));
685	goto skip;
686	}
687	stn[1] = stnZ->leg[1]->l.to;
688	dirn[1] = reverse_leg_dirn(stnZ->leg[1]);
689	stn[2] = stnZ->leg[2]->l.to;
690	dirn[2] = reverse_leg_dirn(stnZ->leg[2]);
691	/print_prefix(stnZ->name);printf(" %p\n",(void)stnZ);*/
692
693	for (i = 0; i < 3; i++) {
694	SVX_ASSERT2(fixed(stn[i]), "stn not fixed for D*");
695
696	leg = stn[i]->leg[dirn[i]];
697
698	SVX_ASSERT2(data_here(leg), "data not on leg for D*");
699	SVX_ASSERT2(leg->l.to == stnZ, "bad sub-network for D*");
700
701	stn2 = reduction_stack->join[i]->l.to;
702
703	if (data_here(reduction_stack->join[i])) {
704	adddd(&POSD(stn2), &POSD(stn[i]), &reduction_stack->join[i]->d);
705	} else {
706	subdd(&POSD(stn2), &POSD(stn[i]), &reverse_leg(reduction_stack->join[i])->d);
707	}
708
709	if (!fZeros(&leg->v)) {
710	delta e, tmp;
711	subdd(&e, &POSD(stnZ), &POSD(stn[i]));
712	subdd(&e, &e, &leg->d);
713	divds(&tmp, &e, &leg->v);
714	if (data_here(reduction_stack->join[i])) {
715	mulsd(&e, &reduction_stack->join[i]->v, &tmp);
716	} else {
717	mulsd(&e, &reverse_leg(reduction_stack->join[i])->v, &tmp);
718	}
719	adddd(&POSD(stn2), &POSD(stn2), &e);
720	}
721	add_stn_to_list(&fixedlist, stn2);
722	free(leg);
723	stn[i]->leg[dirn[i]] = reduction_stack->join[i];
724	/* transfer the articulation status of the radial legs */
725	if (stnZ->leg[i]->l.bits & FLAG_ARTICULATION) {
726	reduction_stack->join[i]->l.bits \|= FLAG_ARTICULATION;
727	reverse_leg(reduction_stack->join[i])->l.bits \|= FLAG_ARTICULATION;
728	}
729	free(stnZ->leg[i]);
730	stnZ->leg[i] = NULL;
731	}
732	/printf("---%f %f %f\n",POS(stnZ, 0), POS(stnZ, 1), POS(stnZ, 2));/
733	remove_stn_from_list(&fixedlist, stnZ);
734	free(stnZ->name);
735	free(stnZ);
736	} else {
737	BUG("reduction_stack has unknown type");
738	}
739
740	skip:;
741	reduction *ptrOld = reduction_stack;
742	reduction_stack = reduction_stack->next;
743	free(ptrOld);
744	}
745	}

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

Context Navigation

source: git/src/network.c @ c981e5b

Download in other formats: