Context Navigation

network.c @ 6d6f4e3

Visit:

stereo-2025

Last change on this file since 6d6f4e3 was f78034b, checked in by Olly Betts <olly@…>, 9 months ago

Eliminate EXPLICIT_FIXED_FLAG option

It's been off by default for over 30 years, and the default approach
of setting the x coordinate to HUGE_VAL to indicate a station is not
fixed requires 4 bytes less memory per station while still being
efficient to check.

Property mode set to 100644

File size: 22.2 KB

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

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

Context Navigation

source: git/src/network.c @ 6d6f4e3

Download in other formats: