Context Navigation

source: git/libpolys/polys/nc/ncSAMult.cc @ 76fead

spielwiese

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

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

Download in other formats: