1 /****************************************
2 * Computer Algebra System SINGULAR *
3 ****************************************/
4 /***************************************************************
5 * File: pDebug.h
6 * Purpose: implementation of debug related poly routines
7 * Author: obachman (Olaf Bachmann)
8 * Created: 8/00
9 *******************************************************************/
10
11 #ifndef PDEBUG_CC
12 #define PDEBUG_CC
13
14 #include <stdarg.h>
15 #include <stdio.h>
16
17
18
19
20
21 #include "misc/auxiliary.h"
22
23
24 #ifdef PDEBUG
25
26 // do the following to always enforce checking of pSetm
27 // #undef PDEBUG
28 // #define PDEBUG 2
29
30 #include "polys/monomials/ring.h"
31 #include "polys/monomials/p_polys.h"
32
33 #include "coeffs/coeffs.h"
34
35 /***************************************************************
36 *
37 * Error reporting
38 *
39 ***************************************************************/
40 // avoid recursive calls
41 STATIC_VAR BOOLEAN d_poly_error_reporting = FALSE;
dPolyReportError(poly p,ring r,const char * fmt,...)42 BOOLEAN dPolyReportError(poly p, ring r, const char* fmt, ...)
43 {
44 if (d_poly_error_reporting) return FALSE;
45 d_poly_error_reporting = TRUE;
46 va_list ap;
47 va_start(ap, fmt);
48
49 fprintf(stderr, "\n// ***dPolyReportError: ");
50 vfprintf(stderr, fmt, ap);
51 fprintf(stderr, "\n occurred at\n");
52 omPrintCurrentBackTraceMax(stderr, 8);
53 if (p != NULL)
54 {
55 fprintf(stderr, " occurred for poly: ");
56 p_wrp(p, r);
57 omPrintAddrInfo(stderr, p, " ");
58 }
59 dErrorBreak();
60 d_poly_error_reporting = FALSE;
61 return FALSE;
62 }
63
64 /***************************************************************
65 *
66 * checking for ring stuff
67 *
68 ***************************************************************/
p_LmCheckIsFromRing(poly p,ring r)69 BOOLEAN p_LmCheckIsFromRing(poly p, ring r)
70 {
71 if (p != NULL)
72 {
73 #if (OM_TRACK > 0) && defined(OM_TRACK_CUSTOM)
74 void* custom = omGetCustomOfAddr(p);
75 if (custom != NULL)
76 {
77 pPolyAssumeReturnMsg(custom == r ||
78 // be more sloppy for qrings
79 (r->qideal != NULL &&
80 omIsBinPageAddr(p) &&
81 omSizeWOfAddr(p)==omSizeWOfBin(r->PolyBin)) ||
82 rSamePolyRep((ring) custom, r),
83 "monomial not from specified ring",p,r);
84 return TRUE;
85 }
86 else
87 #endif
88 #ifndef X_OMALLOC
89 {
90 _pPolyAssumeReturn(omIsBinPageAddr(p),p,r);
91 _pPolyAssumeReturn(omSizeWOfAddr(p)==omSizeWOfBin(r->PolyBin),p,r);
92 return TRUE;
93 }
94 return FALSE;
95 #endif
96 }
97 return TRUE;
98 }
99
p_CheckIsFromRing(poly p,ring r)100 BOOLEAN p_CheckIsFromRing(poly p, ring r)
101 {
102 while (p!=NULL)
103 {
104 pFalseReturn(p_LmCheckIsFromRing(p, r));
105 pIter(p);
106 }
107 return TRUE;
108 }
109
p_CheckPolyRing(poly p,ring r)110 BOOLEAN p_CheckPolyRing(poly p, ring r)
111 {
112 #ifndef X_OMALLOC
113 pAssumeReturn(r != NULL && r->PolyBin != NULL);
114 #endif
115 return p_CheckIsFromRing(p, r);
116 }
117
p_LmCheckPolyRing(poly p,ring r)118 BOOLEAN p_LmCheckPolyRing(poly p, ring r)
119 {
120 #ifndef X_OMALLOC
121 pAssumeReturn(r != NULL && r->PolyBin != NULL);
122 #endif
123 pAssumeReturn(p != NULL);
124 return p_LmCheckIsFromRing(p, r);
125 }
p_CheckRing(ring r)126 BOOLEAN p_CheckRing(ring r)
127 {
128 #ifndef X_OMALLOC
129 pAssumeReturn(r != NULL && r->PolyBin != NULL);
130 #endif
131 return TRUE;
132 }
133
134 /***************************************************************
135 *
136 * Debugging/statistics of pDivisibleBy
137 *
138 ***************************************************************/
p_DebugLmDivisibleByNoComp(poly a,poly b,ring r)139 BOOLEAN p_DebugLmDivisibleByNoComp(poly a, poly b, ring r)
140 {
141 int i=r->N;
142
143 do
144 {
145 if (p_GetExp(a,i,r) > p_GetExp(b,i,r))
146 return FALSE;
147 i--;
148 }
149 while (i);
150 #ifdef HAVE_RINGS
151 return n_DivBy(pGetCoeff(b), pGetCoeff(a), r->cf);
152 #else
153 return TRUE;
154 #endif
155 }
156
157
158 /***************************************************************
159 *
160 * Misc things helpful for debugging
161 *
162 ***************************************************************/
pIsMonomOf(poly p,poly m)163 BOOLEAN pIsMonomOf(poly p, poly m)
164 {
165 if (m == NULL) return TRUE;
166 while (p != NULL)
167 {
168 if (p == m) return TRUE;
169 pIter(p);
170 }
171 return FALSE;
172 }
pHaveCommonMonoms(poly p,poly q)173 BOOLEAN pHaveCommonMonoms(poly p, poly q)
174 {
175 while (p != NULL)
176 {
177 if (pIsMonomOf(q, p))
178 {
179 return TRUE;
180 }
181 pIter(p);
182 }
183 return FALSE;
184 }
185
186 /***************************************************************
187 *
188 * Testing of polys
189 *
190 ***************************************************************/
191 extern void p_Setm_General(poly p, ring r);
192
p_DebugInit(poly p,ring src_ring,ring dest_ring)193 static poly p_DebugInit(poly p, ring src_ring, ring dest_ring)
194 {
195 poly d_p = p_Init(dest_ring);
196 int i;
197 assume(dest_ring->N == src_ring->N);
198
199 for (i=1; i<= src_ring->N; i++)
200 {
201 p_SetExp(d_p, i, p_GetExp(p, i, src_ring), dest_ring);
202 }
203 if (rRing_has_Comp(dest_ring))
204 p_SetComp(d_p, p_GetComp(p, src_ring), dest_ring);
205
206 p_Setm_General(d_p, dest_ring);
207 return d_p;
208 }
209
_p_Test(poly p,ring r,int level)210 BOOLEAN _p_Test(poly p, ring r, int level)
211 {
212 assume(r->cf !=NULL);
213
214 if (PDEBUG > level) level = PDEBUG;
215 if (level < 0 || p == NULL) return TRUE;
216
217 poly p_prev = NULL;
218
219 #ifndef OM_NDEBUG
220 #ifndef X_OMALLOC
221 // check addr with level+1 so as to check bin/page of addr
222 _pPolyAssumeReturnMsg(omTestBinAddrSize(p, (omSizeWOfBin(r->PolyBin))*SIZEOF_LONG, level+1)
223 == omError_NoError, "memory error",p,r);
224 #endif
225 #endif
226
227 pFalseReturn(p_CheckRing(r));
228
229 // this checks that p does not contain a loop: rather expensive O(length^2)
230 #ifndef OM_NDEBUG
231 if (level > 1)
232 pFalseReturn(omTestList(p, level) == omError_NoError);
233 #endif
234
235 int ismod = p_GetComp(p, r) != 0;
236
237 while (p != NULL)
238 {
239 // ring check
240 pFalseReturn(p_LmCheckIsFromRing(p, r));
241 #ifndef OM_NDEBUG
242 #ifndef X_OMALLOC
243 // omAddr check
244 _pPolyAssumeReturnMsg(omTestBinAddrSize(p, (omSizeWOfBin(r->PolyBin))*SIZEOF_LONG, 1)
245 == omError_NoError, "memory error",p,r);
246 #endif
247 #endif
248 // number/coef check
249 _pPolyAssumeReturnMsg(p->coef != NULL || (n_GetChar(r->cf) >= 2), "NULL coef",p,r);
250
251 #ifdef LDEBUG
252 _pPolyAssumeReturnMsg(n_Test(p->coef,r->cf),"coeff err",p,r);
253 #endif
254 _pPolyAssumeReturnMsg(!n_IsZero(p->coef, r->cf), "Zero coef",p,r);
255
256 // check for valid comp
257 _pPolyAssumeReturnMsg(p_GetComp(p, r) >= 0 && (p_GetComp(p, r)<65000), "component out of range ?",p,r);
258 // check for mix poly/vec representation
259 _pPolyAssumeReturnMsg(ismod == (p_GetComp(p, r) != 0), "mixed poly/vector",p,r);
260
261 // special check for ringorder_s/S
262 if ((r->typ!=NULL) && (r->typ[0].ord_typ == ro_syzcomp))
263 {
264 long c1, cc1, ccc1, ec1;
265 sro_ord* o = &(r->typ[0]);
266
267 c1 = p_GetComp(p, r);
268 if (o->data.syzcomp.Components!=NULL)
269 {
270 cc1 = o->data.syzcomp.Components[c1];
271 ccc1 = o->data.syzcomp.ShiftedComponents[cc1];
272 }
273 else { cc1=0; ccc1=0; }
274 _pPolyAssumeReturnMsg(c1 == 0 || cc1 != 0, "Component <-> TrueComponent zero mismatch",p,r);
275 _pPolyAssumeReturnMsg(c1 == 0 || ccc1 != 0,"Component <-> ShiftedComponent zero mismatch",p,r);
276 ec1 = p->exp[o->data.syzcomp.place];
277 //pPolyAssumeReturnMsg(ec1 == ccc1, "Shifted comp out of sync. should %d, is %d");
278 if (ec1 != ccc1)
279 {
280 dPolyReportError(p,r,"Shifted comp out of sync. should %d, is %d",ccc1,ec1);
281 return FALSE;
282 }
283 }
284
285 // check that p_Setm works ok
286 if (level > 0)
287 {
288 poly p_should_equal = p_DebugInit(p, r, r);
289 _pPolyAssumeReturnMsg(p_ExpVectorEqual(p, p_should_equal, r), "p_Setm field(s) out of sync",p,r);
290 p_LmFree(p_should_equal, r);
291 }
292
293 // check order
294 if (p_prev != NULL)
295 {
296 int cmp = p_LmCmp(p_prev, p, r);
297 if (cmp == 0)
298 {
299 _pPolyAssumeReturnMsg(0, "monoms p and p->next are equal", p_prev, r);
300 }
301 else
302 _pPolyAssumeReturnMsg(p_LmCmp(p_prev, p, r) == 1, "wrong order", p_prev, r);
303
304 // check that compare worked sensibly
305 if (level > 1 && p_GetComp(p_prev, r) == p_GetComp(p, r))
306 {
307 int i;
308 for (i=r->N; i>0; i--)
309 {
310 if (p_GetExp(p_prev, i, r) != p_GetExp(p, i, r)) break;
311 }
312 _pPolyAssumeReturnMsg(i > 0, "Exponents equal but compare different", p_prev, r);
313 }
314 }
315 p_prev = p;
316 pIter(p);
317 }
318 return TRUE;
319 }
320
_p_LmTest(poly p,ring r,int level)321 BOOLEAN _p_LmTest(poly p, ring r, int level)
322 {
323 if (level < 0 || p == NULL) return TRUE;
324 poly pnext = pNext(p);
325 pNext(p) = NULL;
326 BOOLEAN test_res = _p_Test(p, r, level);
327 pNext(p) = pnext;
328 return test_res;
329 }
330
_pp_Test(poly p,ring lmRing,ring tailRing,int level)331 BOOLEAN _pp_Test(poly p, ring lmRing, ring tailRing, int level)
332 {
333 if (PDEBUG > level) level = PDEBUG;
334 if (level < 0 || p == NULL) return TRUE;
335 if (pNext(p) == NULL || lmRing == tailRing) return _p_Test(p, lmRing, level);
336
337 pFalseReturn(_p_LmTest(p, lmRing, level));
338 pFalseReturn(_p_Test(pNext(p), tailRing, level));
339
340 // check that lm > Lm(tail)
341 if (level > 1)
342 {
343 poly lm = p;
344 poly tail = p_DebugInit(pNext(p), tailRing, lmRing);
345 poly pnext = pNext(lm);
346 pNext(lm) = tail;
347 BOOLEAN cmp = p_LmCmp(lm, tail, lmRing);
348 if (cmp != 1)
349 dPolyReportError(lm, lmRing, "wrong order: lm <= Lm(tail)");
350 p_LmFree(tail, lmRing);
351 pNext(lm) = pnext;
352 return (cmp == 1);
353 }
354 return TRUE;
355 }
356
357 #endif // PDEBUG
358
359 #if defined(PDEBUG) || defined(PDIV_DEBUG)
360 STATIC_VAR unsigned long pDivisibleBy_number = 1;
361 STATIC_VAR unsigned long pDivisibleBy_FALSE = 1;
362 STATIC_VAR unsigned long pDivisibleBy_ShortFalse = 1;
363
pDebugLmShortDivisibleBy(poly p1,unsigned long sev_1,ring r_1,poly p2,unsigned long not_sev_2,ring r_2)364 BOOLEAN pDebugLmShortDivisibleBy(poly p1, unsigned long sev_1, ring r_1,
365 poly p2, unsigned long not_sev_2, ring r_2)
366 {
367 _pPolyAssume(p_GetShortExpVector(p1, r_1) == sev_1, p1, r_1);
368 _pPolyAssume(p_GetShortExpVector(p2, r_2) == ~ not_sev_2, p2, r_2);
369
370 pDivisibleBy_number++;
371 BOOLEAN ret;
372 if (r_1 == r_2)
373 ret = p_LmDivisibleBy(p1, p2, r_1);
374 else
375 ret = p_LmDivisibleBy(p1, r_1, p2, r_2);
376
377 if (! ret) pDivisibleBy_FALSE++;
378 if (sev_1 & not_sev_2)
379 {
380 pDivisibleBy_ShortFalse++;
381 if (ret)
382 dReportError("p1 divides p2, but sev's are wrong");
383 }
384 return ret;
385 }
386
pDebugLmShortDivisibleByNoComp(poly p1,unsigned long sev_1,ring r_1,poly p2,unsigned long not_sev_2,ring r_2)387 BOOLEAN pDebugLmShortDivisibleByNoComp(poly p1, unsigned long sev_1, ring r_1,
388 poly p2, unsigned long not_sev_2, ring r_2)
389 {
390 // _pPolyAssume((p_GetComp(p1, r_1) == p_GetComp(p2, r_2)) || (p_GetComp(p1, r_1) == 0));
391 _pPolyAssume(p_GetShortExpVector(p1, r_1) == sev_1, p1, r_1);
392 _pPolyAssume(p_GetShortExpVector(p2, r_2) == ~ not_sev_2, p2, r_2);
393
394 pDivisibleBy_number++;
395 BOOLEAN ret;
396 if (r_1 == r_2)
397 ret = p_LmDivisibleByNoComp(p1, p2, r_1);
398 else
399 ret = p_LmDivisibleByNoComp(p1, r_1, p2, r_2);
400
401 if (! ret) pDivisibleBy_FALSE++;
402 if (sev_1 & not_sev_2)
403 {
404 pDivisibleBy_ShortFalse++;
405 if (ret)
406 dReportError("p1 divides p2, but sev's are wrong");
407 }
408 return ret;
409 }
410
pPrintDivisbleByStat()411 void pPrintDivisbleByStat()
412 {
413 Print("#Tests: %ld; #FALSE %ld(%ld); #SHORT %ld(%ld)\n",
414 pDivisibleBy_number,
415 pDivisibleBy_FALSE, (unsigned long) ((double)pDivisibleBy_FALSE*((double) 100)/(double)pDivisibleBy_number),
416 pDivisibleBy_ShortFalse, (unsigned long) ((double)pDivisibleBy_ShortFalse*((double)100)/(double)pDivisibleBy_FALSE));
417 }
418
419 #endif // PDEBUG
420
421 #endif // PDEBUG_CC
422