Context Navigation

source: git/kernel/ncSAMult.cc @ cdec33

spielwiese

Last change on this file since cdec33 was 599326, checked in by Kai Krüger <krueger@…>, 14 years ago
Anne, Kai, Frank: - changes to #include "..." statements to allow cleaner build structure - affected directories: omalloc, kernel, Singular - not yet done: IntergerProgramming git-svn-id: file:///usr/local/Singular/svn/trunk@13032 2c84dea3-7e68-4137-9b89-c4e89433aadc
Property mode set to `100644`
File size: 24.3 KB

Line
1	/****************************************
2	* Computer Algebra System SINGULAR *
3	****************************************/
4	/***************************************************************
5	* File: ncSAMult.cc
6	* Purpose: implementation of multiplication in simple NC subalgebras
7	* Author: motsak
8	* Created:
9	* Version: $Id$
10	*******************************************************************/
11
12	#define MYTEST 0
13
14	#if MYTEST
15	#define OM_CHECK 4
16	#define OM_TRACK 5
17	// these settings must be before "mod2.h" in order to work!!!
18	#endif
19
20
21	#include <kernel/mod2.h>
22
23	#ifdef HAVE_PLURAL
24
25	#ifndef NDEBUG
26	#define OUTPUT 1
27	#else
28	#define OUTPUT 0
29	#endif
30
31
32	#include <kernel/options.h>
33	#include <kernel/p_Procs.h>
34	#include <kernel/ncSAMult.h> // for CMultiplier etc classes
35	#include <kernel/sca.h> // for SCA
36	#include <kernel/polys.h> // for p_One
37
38
39
40	// poly functions defined in p_Procs: ;
41	static poly gnc_pp_Mult_mm(const poly p, const poly m, const ring r, poly& last)
42	{
43	if( (p == NULL) \|\| (m == NULL) )
44	return NULL;
45
46	assume( (p != NULL) && (m != NULL) && (r != NULL) );
47
48	#if OUTPUT
49	PrintS("VVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVV gnc_pp_Mult_mm(p, m) VVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVV ");
50	PrintLn();
51	PrintS("p: "); p_Write(p, r);
52	PrintS("m: "); p_Write(m, r);
53	#endif
54	poly pResult;
55
56	if (p_IsConstant(m, r))
57	pResult = pp_Mult_nn(p, p_GetCoeff(m,r),r);
58	else
59	{
60	CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier();
61	assume( pMultiplier != NULL );
62
63	poly pMonom = pMultiplier->LM(m, r);
64	pResult = pMultiplier->MultiplyPE(p, pMonom);
65	p_Delete(&pMonom, r);
66	p_Test(pResult, r);
67	pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r);
68	}
69
70	#if OUTPUT
71	p_Test(pResult, r);
72
73	PrintS("gnc_pp_Mult_mm(p, m) => "); p_Write(pResult, r);
74	PrintS("p: "); p_Write(p, r);
75	PrintS("m: "); p_Write(m, r);
76	PrintS("^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ ");
77	PrintLn();
78	#endif
79
80	return pResult;
81
82	}
83
84	static poly gnc_p_Mult_mm(poly p, const poly m, const ring r)
85	{
86	if( (p == NULL) \|\| (m == NULL) )
87	{
88	p_Delete(&p, r);
89	return NULL;
90	}
91
92	assume( (p != NULL) && (m != NULL) && (r != NULL) );
93
94	#if OUTPUT
95	PrintS("VVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVV gnc_p_Mult_mm(p, m) VVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVV ");
96	PrintLn();
97	PrintS("p: ");
98	p_Write(p, r);
99	PrintS("m: ");
100	p_Write(m, r);
101	#endif
102
103	poly pResult;
104
105	if (p_IsConstant(m, r))
106	pResult = p_Mult_nn(p, p_GetCoeff(m,r),r);
107	else
108	{
109	CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier();
110	assume( pMultiplier != NULL );
111
112	poly pMonom = pMultiplier->LM(m, r);
113	pResult = pMultiplier->MultiplyPEDestroy(p, pMonom);
114	p_Delete(&pMonom, r);
115	p_Test(pResult, r);
116	pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r);
117	}
118
119	#if OUTPUT
120	p_Test(pResult, r);
121
122	PrintS("gnc_p_Mult_mm(p, m) => "); p_Write(pResult, r);
123	// PrintS("p: "); p_Write(p, r);
124	PrintS("m: "); p_Write(m, r);
125	PrintS("^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ ");
126	PrintLn();
127	#endif
128
129	return pResult;
130
131	}
132
133	static poly gnc_mm_Mult_p(const poly m, poly p, const ring r)
134	{
135	if( (p == NULL) \|\| (m == NULL) )
136	{
137	p_Delete(&p, r);
138	return NULL;
139	}
140
141	assume( (p != NULL) && (m != NULL) && (r != NULL) );
142
143	p_Test(m, r);
144	p_Test(p, r);
145
146	#if OUTPUT
147	PrintS("VVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVV gnc_mm_Mult_p(m, p) VVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVV ");
148	PrintLn();
149	PrintS("m: "); p_Write(m, r);
150	PrintS("p: "); p_Write(p, r);
151	#endif
152
153	poly pResult;
154
155	if (p_IsConstant(m, r))
156	pResult = p_Mult_nn(p, p_GetCoeff(m,r),r);
157	else
158	{
159	CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier();
160	assume( pMultiplier != NULL );
161
162	poly pMonom = pMultiplier->LM(m, r);
163	pResult = pMultiplier->MultiplyEPDestroy(pMonom, p);
164	p_Delete(&pMonom, r);
165	p_Test(pResult, r);
166	pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r);
167	}
168
169	#if OUTPUT
170	p_Test(pResult, r);
171
172	PrintS("gnc_mm_Mult_p(m, p) => "); p_Write(pResult, r);
173	// PrintS("p: "); p_Write(p, r);
174	PrintS("m: "); p_Write(m, r);
175	PrintS("^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ ");
176	PrintLn();
177	#endif
178
179	return pResult;
180	}
181
182	static poly gnc_mm_Mult_pp(const poly m, const poly p, const ring r)
183	{
184	if( (p == NULL) \|\| (m == NULL) )
185	{
186	return NULL;
187	}
188
189	assume( (p != NULL) && (m != NULL) && (r != NULL) );
190
191	p_Test(m, r);
192	p_Test(p, r);
193
194	#if OUTPUT
195	PrintS("VVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVV gnc_mm_Mult_pp(m, p) VVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVVV ");
196	PrintLn();
197	PrintS("m: "); p_Write(m, r);
198	PrintS("p: "); p_Write(p, r);
199	#endif
200
201	poly pResult;
202
203	if (p_IsConstant(m, r))
204	pResult = pp_Mult_nn(p, p_GetCoeff(m,r),r);
205	else
206	{
207	CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier();
208	assume( pMultiplier != NULL );
209
210	poly pMonom = pMultiplier->LM(m, r);
211	pResult = pMultiplier->MultiplyEP(pMonom, p);
212	p_Delete(&pMonom, r);
213	p_Test(pResult, r);
214	pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r);
215	}
216
217	#if OUTPUT
218	p_Test(pResult, r);
219
220	PrintS("gnc_mm_Mult_pp(m, p) => "); p_Write(pResult, r);
221	PrintS("p: "); p_Write(p, r);
222	PrintS("m: "); p_Write(m, r);
223	PrintS("^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ ");
224	PrintLn();
225	#endif
226
227	return pResult;
228	}
229
230	static void gnc_p_ProcsSet(ring rGR, p_Procs_s* p_Procs = NULL)
231	{
232	#if OUTPUT
233	PrintS("\|gnc_p_ProcsSet()");
234	PrintLn();
235	#endif
236
237	if( p_Procs == NULL )
238	p_Procs = rGR->p_Procs;
239
240	// "commutative"
241	p_Procs->p_Mult_mm = rGR->p_Procs->p_Mult_mm = gnc_p_Mult_mm;
242	p_Procs->pp_Mult_mm = rGR->p_Procs->pp_Mult_mm = gnc_pp_Mult_mm;
243
244	p_Procs->p_Minus_mm_Mult_qq = rGR->p_Procs->p_Minus_mm_Mult_qq = NULL;
245
246	// non-commutaitve multiplication by monomial from the left
247	rGR->GetNC()->p_Procs.mm_Mult_p = gnc_mm_Mult_p;
248	rGR->GetNC()->p_Procs.mm_Mult_pp = gnc_mm_Mult_pp;
249
250	}
251
252	bool ncInitSpecialPairMultiplication(ring r)
253	{
254	#if OUTPUT
255	PrintS("ncInitSpecialPairMultiplication(ring), ring: \n");
256	rWrite(r);
257	PrintLn();
258	#endif
259
260	assume(rIsPluralRing(r));
261	assume(!rIsSCA(r));
262
263	if( r->GetNC()->GetGlobalMultiplier() != NULL )
264	{
265	WarnS("Already defined!");
266	return false;
267	}
268
269	r->GetNC()->GetGlobalMultiplier() = new CGlobalMultiplier(r);
270
271	gnc_p_ProcsSet(r);
272	return true;
273	}
274
275
276	CGlobalMultiplier::CGlobalMultiplier(ring r):
277	CMultiplier<poly>(r), m_RingFormulaMultiplier(GetFormulaPowerMultiplier(r))
278	{
279	#if OUTPUT
280	PrintS("CGlobalMultiplier::CGlobalMultiplier(ring)!");
281	PrintLn();
282	#endif
283
284	m_cache = new CGlobalCacheHash(r);
285	m_powers = new CPowerMultiplier(r);
286	}
287
288
289	CGlobalMultiplier::~CGlobalMultiplier()
290	{
291	#if OUTPUT
292	PrintS("CGlobalMultiplier::~CGlobalMultiplier()!");
293	PrintLn();
294	#endif
295
296	delete m_cache;
297	delete m_powers;
298
299	// we cannot delete m_RingFormulaMultiplier as it belongs to the ring!
300	}
301
302
303
304	// Exponent * Exponent
305	// TODO: handle components!!!
306	poly CGlobalMultiplier::MultiplyEE(const CGlobalMultiplier::CExponent expLeft, const CGlobalMultiplier::CExponent expRight)
307	{
308
309	const ring r = GetBasering();
310
311	#if OUTPUT
312	PrintS("CGlobalMultiplier::MultiplyEE(expLeft, expRight)!");
313	PrintLn();
314	PrintS("expL: "); p_Write(expLeft, GetBasering());
315	PrintS("expR: "); p_Write(expRight, GetBasering());
316	#endif
317
318	CCacheHash<poly>::CCacheItem* pLookup;
319
320	int b = m_cache->LookupEE(expLeft, expRight, pLookup);
321	// TODO!!!
322
323	// up to now:
324	assume( b == -1 );
325
326	// TODO: use PowerMultiplier!!!!
327
328	poly product = NULL;
329
330	const int N = NVars();
331	int j = N;
332	int i = 1;
333
334	int ej = p_GetExp(expLeft, j, r);
335	int ei = p_GetExp(expRight, i, r);
336
337	while( (i < j) && !((ej != 0) && (ei != 0)) )
338	{
339	if( ei == 0 )
340	ei = p_GetExp(expRight, ++i, r);
341
342	if( ej == 0 )
343	ej = p_GetExp(expLeft, --j, r);
344	}
345
346
347	#if OUTPUT
348	PrintS("<CGlobalMultiplier::MultiplyEE>");
349	PrintLn();
350	Print("i: %d, j: %d", i, j);
351	PrintLn();
352	Print("ei: %d, ej: %d", ei, ej);
353	PrintLn();
354	#endif
355
356
357	// \| expLeft \| * \| expRight \|
358	// \|<<<< ej 0..0\| , \|0..0 ei >>>>\|
359	// \|<<<< j <<<N\| , \|1>>> i >>>>\|
360
361	if( i >= j ) // BUG here!!!???
362	{
363	// either i == j or i = j + 1 => commutative multiple!
364	// TODO: it can be done more efficiently! ()
365	product = p_Head(expRight, r);
366
367	// \| expLeft \| * \| expRight \|
368	// \|<<<< ej 0....0\| , \|0..00 ei >>>>\|
369	// \|<<<< j i <<<N\| , \|1>>>j i >>>>\|
370
371	if(i > j)
372	{
373	--i;
374	ei = 0;
375	}
376
377	if( i == j )
378	{
379	if( ej != 0 )
380	p_SetExp(product, i, ei + ej, r);
381	}
382
383	--i;
384
385	for(; i > 0; --i)
386	{
387	const int e = p_GetExp(expLeft, i, r);
388
389	if( e > 0 )
390	p_SetExp(product, i, e, r);
391	}
392
393	p_Setm(product, r);
394
395	} else
396	{ // i < j, ei != 0, ej != 0
397
398	Enum_ncSAType PairType = _ncSA_notImplemented;
399
400	if( m_RingFormulaMultiplier != NULL )
401	PairType = m_RingFormulaMultiplier->GetPair(i, j);
402
403
404	if( PairType == _ncSA_notImplemented )
405	product = m_powers->MultiplyEE( CPower(j, ej), CPower(i, ei) );
406	// return gnc_uu_Mult_ww_vert(i, a, j, b, r);
407	else
408	// return m_RingFormulaMultiplier->Multiply(j, i, b, a);
409	product = CFormulaPowerMultiplier::Multiply( PairType, i, j, ei, ej, GetBasering());
410
411
412	#if OUTPUT
413	PrintS("<CGlobalMultiplier::MultiplyEE> ==> ");
414	PrintLn();
415	Print("i: %d, j: %d", i, j);
416	PrintLn();
417	Print("ei: %d, ej: %d", ei, ej);
418	PrintLn();
419	PrintS("<product>: "); p_Write(product, GetBasering());
420	#endif
421
422
423	// TODO: Choose some multiplication strategy!!!
424
425	while( (product != NULL) && !((i == NVars()) && (j == 1)) )
426	{
427
428	// make some choice here!:
429
430	if( i < NVars() )
431	{
432	ei = p_GetExp(expRight, ++i, r);
433
434	while( (ei == 0) && (i < NVars()) )
435	ei = p_GetExp(expRight, ++i, r);
436
437	if( ei != 0 )
438	product = m_powers->MultiplyPEDestroy(product, CPower(i, ei));
439	}
440
441	if( j > 1 )
442	{
443	ej = p_GetExp(expLeft, --j, r);
444
445	while( (ej == 0) && (1 < j) )
446	ej = p_GetExp(expLeft, --j, r);
447
448	if( ej != 0 )
449	product = m_powers->MultiplyEPDestroy(CPower(j, ej), product);
450	}
451
452
453	#if OUTPUT
454	PrintS("<CGlobalMultiplier::MultiplyEE> ==> ");
455	PrintLn();
456	Print("i: %d, j: %d", i, j);
457	PrintLn();
458	Print("ei: %d, ej: %d", ei, ej);
459	PrintLn();
460	PrintS("<product>: "); p_Write(product, GetBasering());
461	#endif
462
463	}
464
465	}
466
467	// TODO!
468
469
470	m_cache->StoreEE( expLeft, expRight, product);
471	// up to now:
472	return product;
473	}
474
475	// Monom * Exponent
476	poly CGlobalMultiplier::MultiplyME(const poly pMonom, const CGlobalMultiplier::CExponent expRight)
477	{
478	#if OUTPUT
479	PrintS("CGlobalMultiplier::MultiplyME(monom, expR)!");
480	PrintLn();
481	PrintS("Monom: "); p_Write(pMonom, GetBasering());
482	PrintS("expR: "); p_Write(expRight, GetBasering());
483	#endif
484
485	return MultiplyEE(pMonom, expRight);
486	}
487
488	// Exponent * Monom
489	poly CGlobalMultiplier::MultiplyEM(const CGlobalMultiplier::CExponent expLeft, const poly pMonom)
490	{
491	#if OUTPUT
492	PrintS("CGlobalMultiplier::MultiplyEM(expL, monom)!");
493	PrintLn();
494	PrintS("expL: "); p_Write(expLeft, GetBasering());
495	PrintS("Monom: "); p_Write(pMonom, GetBasering());
496	#endif
497
498	return MultiplyEE(expLeft, pMonom);
499	}
500
501
502
503
504
505	///////////////////////////////////////////////////////////////////////////////////////////
506	CCommutativeSpecialPairMultiplier::CCommutativeSpecialPairMultiplier(ring r, int i, int j):
507	CSpecialPairMultiplier(r, i, j)
508	{
509	#if OUTPUT
510	Print("CCommutativeSpecialPairMultiplier::CCommutativeSpecialPairMultiplier(ring, i: %d, j: %d)!", i, j);
511	PrintLn();
512	#endif
513	};
514
515
516	CCommutativeSpecialPairMultiplier::~CCommutativeSpecialPairMultiplier()
517	{
518	#if OUTPUT
519	PrintS("CCommutativeSpecialPairMultiplier::~CCommutativeSpecialPairMultiplier()");
520	PrintLn();
521	#endif
522	}
523
524	// Exponent * Exponent
525	poly CCommutativeSpecialPairMultiplier::MultiplyEE(const int expLeft, const int expRight)
526	{
527	#if OUTPUT
528	Print("CCommutativeSpecialPairMultiplier::MultiplyEE(var(%d)^{%d}, var(%d)^{%d})!", GetJ(), expLeft, GetI(), expRight);
529	PrintLn();
530	#endif
531
532	const ring r = GetBasering();
533
534	return CFormulaPowerMultiplier::ncSA_1xy0x0y0(GetI(), GetJ(), expRight, expLeft, r);
535	}
536
537	///////////////////////////////////////////////////////////////////////////////////////////
538	CAntiCommutativeSpecialPairMultiplier::CAntiCommutativeSpecialPairMultiplier(ring r, int i, int j):
539	CSpecialPairMultiplier(r, i, j)
540	{
541	#if OUTPUT
542	Print("CAntiCommutativeSpecialPairMultiplier::CAntiCommutativeSpecialPairMultiplier(ring, i: %d, j: %d)!", i, j);
543	PrintLn();
544	#endif
545	};
546
547
548	CAntiCommutativeSpecialPairMultiplier::~CAntiCommutativeSpecialPairMultiplier()
549	{
550	#if OUTPUT
551	PrintS("CAntiCommutativeSpecialPairMultiplier::~CAntiCommutativeSpecialPairMultiplier()");
552	PrintLn();
553	#endif
554	}
555
556	// Exponent * Exponent
557	poly CAntiCommutativeSpecialPairMultiplier::MultiplyEE(const int expLeft, const int expRight)
558	{
559	#if OUTPUT
560	Print("CAntiCommutativeSpecialPairMultiplier::MultiplyEE(var(%d)^{%d}, var(%d)^{%d})!", GetJ(), expLeft, GetI(), expRight);
561	PrintLn();
562	#endif
563
564	const ring r = GetBasering();
565
566	return CFormulaPowerMultiplier::ncSA_Mxy0x0y0(GetI(), GetJ(), expRight, expLeft, r);
567	}
568
569	///////////////////////////////////////////////////////////////////////////////////////////
570	CQuasiCommutativeSpecialPairMultiplier::CQuasiCommutativeSpecialPairMultiplier(ring r, int i, int j, number q):
571	CSpecialPairMultiplier(r, i, j), m_q(q)
572	{
573	#if OUTPUT
574	Print("CQuasiCommutativeSpecialPairMultiplier::CQuasiCommutativeSpecialPairMultiplier(ring, i: %d, j: %d, q)!", i, j);
575	PrintLn();
576	PrintS("Parameter q: ");
577	n_Write(q, r);
578	#endif
579	};
580
581
582	CQuasiCommutativeSpecialPairMultiplier::~CQuasiCommutativeSpecialPairMultiplier()
583	{
584	#if OUTPUT
585	PrintS("CQuasiCommutativeSpecialPairMultiplier::~CQuasiCommutativeSpecialPairMultiplier()");
586	PrintLn();
587	#endif
588	}
589
590	// Exponent * Exponent
591	poly CQuasiCommutativeSpecialPairMultiplier::MultiplyEE(const int expLeft, const int expRight)
592	{
593	#if OUTPUT
594	Print("CQuasiCommutativeSpecialPairMultiplier::MultiplyEE(var(%d)^{%d}, var(%d)^{%d})!", GetJ(), expLeft, GetI(), expRight);
595	PrintLn();
596	#endif
597
598	const ring r = GetBasering();
599
600	return CFormulaPowerMultiplier::ncSA_Qxy0x0y0(GetI(), GetJ(), expRight, expLeft, m_q, r);
601	}
602
603
604	///////////////////////////////////////////////////////////////////////////////////////////
605	CWeylSpecialPairMultiplier::CWeylSpecialPairMultiplier(ring r, int i, int j, number g):
606	CSpecialPairMultiplier(r, i, j), m_g(g)
607	{
608	#if OUTPUT
609	Print("CWeylSpecialPairMultiplier::CWeylSpecialPairMultiplier(ring, i: %d, j: %d, g)!", i, j);
610	PrintLn();
611	PrintS("Parameter g: ");
612	n_Write(g, r);
613	#endif
614	};
615
616
617	CWeylSpecialPairMultiplier::~CWeylSpecialPairMultiplier()
618	{
619	#if OUTPUT
620	PrintS("CWeylSpecialPairMultiplier::~CWeylSpecialPairMultiplier()");
621	PrintLn();
622	#endif
623	}
624
625	// Exponent * Exponent
626	poly CWeylSpecialPairMultiplier::MultiplyEE(const int expLeft, const int expRight)
627	{
628	#if OUTPUT
629	Print("CWeylSpecialPairMultiplier::MultiplyEE(var(%d)^{%d}, var(%d)^{%d})!", GetJ(), expLeft, GetI(), expRight);
630	PrintLn();
631	#endif
632	// Char == 0, otherwise - problem!
633
634
635	const ring r = GetBasering();
636
637	assume( expLeft*expRight > 0 );
638
639	return CFormulaPowerMultiplier::ncSA_1xy0x0yG(GetI(), GetJ(), expRight, expLeft, m_g, r);
640	}
641
642
643	///////////////////////////////////////////////////////////////////////////////////////////
644	CShiftSpecialPairMultiplier::CShiftSpecialPairMultiplier(ring r, int i, int j, int s, number c):
645	CSpecialPairMultiplier(r, i, j), m_shiftCoef(c), m_shiftVar(s)
646	{
647	#if OUTPUT
648	Print("CShiftSpecialPairMultiplier::CShiftSpecialPairMultiplier(ring, i: %d, j: %d, s: %d, c)!", i, j, s);
649	PrintLn();
650	PrintS("Parameter c: "); n_Write(c, r);
651	#endif
652	};
653
654
655	CShiftSpecialPairMultiplier::~CShiftSpecialPairMultiplier()
656	{
657	#if OUTPUT
658	PrintS("CShiftSpecialPairMultiplier::~CShiftSpecialPairMultiplier()");
659	PrintLn();
660	#endif
661	}
662
663	// Exponent * Exponent
664	poly CShiftSpecialPairMultiplier::MultiplyEE(const int expLeft, const int expRight)
665	{
666	#if OUTPUT
667	Print("CShiftSpecialPairMultiplier::MultiplyEE(var(%d)^{%d}, var(%d)^{%d})!", GetJ(), expLeft, GetI(), expRight);
668	PrintLn();
669	#endif
670	// Char == 0, otherwise - problem!
671
672	assume( expLeft*expRight > 0 );
673
674	const ring r = GetBasering();
675
676	if( m_shiftVar != GetI() ) // YX = XY + b*Y?
677	return CFormulaPowerMultiplier::ncSA_1xy0xBy0(GetI(), GetJ(), expRight, expLeft, m_shiftCoef, r); // case: (1, 0, beta, 0, 0)
678	else
679	return CFormulaPowerMultiplier::ncSA_1xyAx0y0(GetI(), GetJ(), expRight, expLeft, m_shiftCoef, r); // case: (1, alpha, 0, 0)
680
681	}
682
683
684
685	///////////////////////////////////////////////////////////////////////////////////////////
686	CExternalSpecialPairMultiplier::CExternalSpecialPairMultiplier(ring r, int i, int j, Enum_ncSAType type):
687	CSpecialPairMultiplier(r, i, j), m_ncSAtype(type)
688	{
689	#if OUTPUT
690	Print("CExternalSpecialPairMultiplier::CExternalSpecialPairMultiplier(ring, i: %d, j: %d, type: %d, c)!", i, j, (int)type);
691	PrintLn();
692	#endif
693	};
694
695
696	CExternalSpecialPairMultiplier::~CExternalSpecialPairMultiplier()
697	{
698	#if OUTPUT
699	PrintS("CExternalSpecialPairMultiplier::~CExternalSpecialPairMultiplier()");
700	PrintLn();
701	#endif
702	}
703
704	// Exponent * Exponent
705	poly CExternalSpecialPairMultiplier::MultiplyEE(const int expLeft, const int expRight)
706	{
707	#if OUTPUT
708	Print("CExternalSpecialPairMultiplier::MultiplyEE(var(%d)^{%d}, var(%d)^{%d})!", GetJ(), expLeft, GetI(), expRight);
709	PrintLn();
710	#endif
711	// Char == 0, otherwise - problem!
712
713	assume( expLeft*expRight > 0 );
714
715	const ring r = GetBasering();
716
717	return CFormulaPowerMultiplier::Multiply(m_ncSAtype, GetI(), GetJ(), expRight, expLeft, r);
718
719	}
720
721
722
723	///////////////////////////////////////////////////////////////////////////////////////////
724
725	// factory method!
726	CSpecialPairMultiplier* AnalyzePair(const ring r, int i, int j)
727	{
728	#if OUTPUT
729	Print("AnalyzePair(ring, i: %d, j: %d)!", i, j);
730	PrintLn();
731	#endif
732
733	Enum_ncSAType type = CFormulaPowerMultiplier::AnalyzePair(r, i, j);
734
735	if( type == _ncSA_notImplemented ) return NULL;
736
737
738	// last possibility:
739	return new CExternalSpecialPairMultiplier(r, i, j, type); // For tests!
740
741
742	if( type == _ncSA_1xy0x0y0 )
743	return new CCommutativeSpecialPairMultiplier(r, i, j);
744
745	if( type == _ncSA_Mxy0x0y0 )
746	return new CAntiCommutativeSpecialPairMultiplier(r, i, j);
747
748	if( type == _ncSA_Qxy0x0y0 )
749	{
750	const number q = p_GetCoeff(GetC(r, i, j), r);
751	return new CQuasiCommutativeSpecialPairMultiplier(r, i, j, q);
752	}
753
754	const poly d = GetD(r, i, j);
755
756	assume( d != NULL );
757	assume( pNext(d) == NULL );
758
759	const number g = p_GetCoeff(d, r);
760
761	if( type == _ncSA_1xy0x0yG ) // Weyl
762	return new CWeylSpecialPairMultiplier(r, i, j, g);
763
764	if( type == _ncSA_1xyAx0y0 ) // Shift 1
765	return new CShiftSpecialPairMultiplier(r, i, j, i, g);
766
767	if( type == _ncSA_1xy0xBy0 ) // Shift 2
768	return new CShiftSpecialPairMultiplier(r, i, j, j, g);
769
770	}
771
772
773
774
775
776
777	CPowerMultiplier::CPowerMultiplier(ring r): CMultiplier<CPower>(r)
778	{
779	#if OUTPUT
780	PrintS("CPowerMultiplier::CPowerMultiplier(ring)!");
781	PrintLn();
782	#endif
783
784	m_specialpairs = (CSpecialPairMultiplier*)omAlloc0( ((NVars() (NVars()-1)) / 2) * sizeof(CSpecialPairMultiplier*) );
785
786	for( int i = 1; i < NVars(); i++ )
787	for( int j = i + 1; j <= NVars(); j++ )
788	GetPair(i, j) = AnalyzePair(GetBasering(), i, j); // factory method!
789	}
790
791
792	CPowerMultiplier::~CPowerMultiplier()
793	{
794	#if OUTPUT
795	PrintS("CPowerMultiplier::~CPowerMultiplier()!");
796	PrintLn();
797	#endif
798
799	omFreeSize((ADDRESS)m_specialpairs, ((NVars() * (NVars()-1)) / 2) * sizeof(CSpecialPairMultiplier*) );
800	}
801
802
803	// Monom * Exponent
804	// pMonom may NOT be of the form: var(j)^{n}!
805	poly CPowerMultiplier::MultiplyME(const poly pMonom, const CExponent expRight)
806	{
807	const int j = expRight.Var;
808	const int n = expRight.Power;
809
810	const ring r = GetBasering();
811
812	#if OUTPUT
813	Print("CPowerMultiplier::MultiplyME(monom * var(%d)^{%d})!", j, n);
814	PrintLn();
815	PrintS("Monom: "); p_Write(pMonom, r);
816	#endif
817
818	assume( (j > 0) && (j <= NVars()));
819
820	if( n == 0 )
821	return p_Head(pMonom, r); // Copy?!?
822
823
824	int v = NVars();
825	int e = p_GetExp(pMonom, v, r);
826
827	while((v > j) && (e == 0))
828	e = p_GetExp(pMonom, --v, r);
829
830	// TODO: review this!
831	if( (v == j) )
832	{
833	poly p = p_Head(pMonom, r);
834	p_SetExp(p, v, e + n, r);
835	p_Setm(p, r);
836
837	return p;
838	}
839
840	assume( v > j );
841	assume( e > 0 );
842
843	// And now the General Case: v > j!
844
845	poly p = MultiplyEE( CPower(v, e), expRight ); // Easy way!
846
847	--v;
848
849	while(v > 0)
850	{
851	e = p_GetExp(pMonom, v, GetBasering());
852
853	if( e > 0 )
854	p = MultiplyEPDestroy(CPower(v, e), p);
855
856	--v;
857	}
858
859	#if OUTPUT
860	PrintS("CPowerMultiplier::MultiplyME() ===> ");
861	p_Write(p, GetBasering());
862	#endif
863
864	return p;
865	}
866
867	// Exponent * Monom
868	// pMonom may NOT be of the form: var(i)^{m}!
869	poly CPowerMultiplier::MultiplyEM(const CExponent expLeft, const poly pMonom)
870	{
871	const ring r = GetBasering();
872
873	// TODO: as above! (difference due to Left/Right semmantics!)
874	const int j = expLeft.Var;
875	const int n = expLeft.Power;
876
877	#if OUTPUT
878	Print("CPowerMultiplier::MultiplyEM(var(%d)^{%d} * monom)!", j, n);
879	PrintLn();
880	PrintS("Monom: "); p_Write(pMonom, r);
881	#endif
882
883	assume( (j > 0) && (j <= NVars()));
884
885	if( n == 0 )
886	return p_Head(pMonom, r); // Copy?!?
887
888
889	int v = 1; // NVars();
890	int e = p_GetExp(pMonom, v, r);
891
892	while((v < j) && (e == 0))
893	e = p_GetExp(pMonom, ++v, r);
894
895	if( v == j )
896	{
897	poly p = p_Head(pMonom, r);
898	p_SetExp(p, j, e + n, r);
899	p_Setm(p, r);
900
901	return p;
902	}
903
904	assume( v < j );
905	assume( e > 0 );
906
907
908	// And now the General Case: v > j!
909
910	poly p = MultiplyEE( expLeft, CPower(v, e) ); // Easy way!
911
912	++v;
913
914	while(v <= NVars())
915	{
916	e = p_GetExp(pMonom, v, r);
917
918	if( e > 0 )
919	p = MultiplyPEDestroy(p, CPower(v, e));
920
921	++v;
922	}
923
924	#if OUTPUT
925	PrintS("CPowerMultiplier::MultiplyEM() ===> ");
926	p_Write(p, r);
927	#endif
928
929	return p;
930
931	}
932
933
934	// Exponent * Exponent
935	// Computes: var(j)^{expLeft} * var(i)^{expRight}
936	poly CPowerMultiplier::MultiplyEE(const CExponent expLeft, const CExponent expRight)
937	{
938	#if OUTPUT
939	PrintS("CPowerMultiplier::MultiplyEE)!");
940	PrintLn();
941	#endif
942
943	const int i = expRight.Var, j = expLeft.Var;
944	const int ei = expRight.Power, ej = expLeft.Power;
945
946	#if OUTPUT
947	Print("Input: var(%d)^{%d} * var(%d)^{%d}", j, ej, i, ei);
948	PrintLn();
949	#endif
950
951	assume(1 <= i);
952	assume(j <= NVars());
953	assume(1 <= j);
954	assume(i <= NVars());
955	assume(ei > 0);
956	assume(ej > 0);
957
958	if( i >= j )
959	{
960	const ring r = GetBasering();
961
962	poly product = p_One(r);
963	p_SetExp(product, j, ej, r);
964	p_SetExp(product, i, ei, r);
965	p_Setm(product, r);
966
967	return product;
968
969	} else
970	{
971	assume(i < j);
972
973	// No Cache Lookup!? :(
974
975	CSpecialPairMultiplier* pSpecialMultiplier = GetPair(i, j);
976
977	poly product = NULL;
978
979	// Special case?
980	if( pSpecialMultiplier != NULL )
981	{
982	assume( pSpecialMultiplier->GetI() == i );
983	assume( pSpecialMultiplier->GetJ() == j );
984	assume( pSpecialMultiplier->GetBasering() == GetBasering() );
985
986	return pSpecialMultiplier->MultiplyEE(ej, ei);
987	} else
988	{
989	// Perform general NC Multiplication:
990	// TODO
991
992	// product = NULL;
993	}
994	}
995
996	return NULL;
997	}
998
999
1000
1001
1002
1003
1004	CSpecialPairMultiplier::CSpecialPairMultiplier(ring r, int i, int j):
1005	CMultiplier<int>(r), m_i(i), m_j(j)
1006	{
1007	#if OUTPUT
1008	Print("CSpecialPairMultiplier::CSpecialPairMultiplier(ring, i: %d, j: %d)!", i, j);
1009	PrintLn();
1010	#endif
1011
1012	assume(i < j);
1013	assume(i > 0);
1014	assume(j <= NVars());
1015	}
1016
1017
1018	CSpecialPairMultiplier::~CSpecialPairMultiplier()
1019	{
1020	#if OUTPUT
1021	PrintS("CSpecialPairMultiplier::~CSpecialPairMultiplier()!");
1022	PrintLn();
1023	#endif
1024	}
1025
1026
1027
1028	// Monom * Exponent
1029	poly CSpecialPairMultiplier::MultiplyME(const poly pMonom, const CExponent expRight)
1030	{
1031	#if OUTPUT
1032	Print("CSpecialPairMultiplier::MultiplyME(monom, var(%d)^{%d})!", GetI(), expRight);
1033	PrintLn();
1034	PrintS("Monom: "); p_Write(pMonom, GetBasering());
1035	#endif
1036
1037	return MultiplyEE(p_GetExp(pMonom, GetJ(), GetBasering()), expRight);
1038	}
1039
1040	// Exponent * Monom
1041	poly CSpecialPairMultiplier::MultiplyEM(const CExponent expLeft, const poly pMonom)
1042	{
1043	#if OUTPUT
1044	Print("CSpecialPairMultiplier::MultiplyEM(var(%d)^{%d}, monom)!", GetJ(), expLeft);
1045	PrintLn();
1046	PrintS("Monom: "); p_Write(pMonom, GetBasering());
1047	#endif
1048
1049	return MultiplyEE(expLeft, p_GetExp(pMonom, GetI(), GetBasering()));
1050	}
1051	#endif

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

Download in other formats: