Context Navigation

source: git/kernel/ncSAMult.cc @ 690e21e

spielwiese

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

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

Download in other formats: