source: git/libpolys/polys/pDebug.cc @ 85bcd6

jengelh-datetimespielwiese
Last change on this file since 85bcd6 was 85bcd6, checked in by Oleksandr Motsak <motsak@…>, 11 years ago
added p_LmShortDivisibleByNoComp by analogy add: also needed pDebugLmShortDivisibleByNoComp and p_LmDivisibleByNoComp (with 2 rings) for completeness
  • Property mode set to 100644
File size: 10.6 KB
Line 
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#include "config.h"
19#include <misc/auxiliary.h>
20
21
22#ifdef PDEBUG
23
24// do the following to always enforce checking of pSetm
25// #undef PDEBUG
26// #define PDEBUG 2
27
28#include <omalloc/omalloc.h>
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
41static BOOLEAN d_poly_error_reporting = FALSE;
42BOOLEAN 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// ***dPolyError: ");
50  vfprintf(stderr, fmt, ap);
51  fprintf(stderr, "\n occured at\n");
52  omPrintCurrentBackTraceMax(stderr, 8);
53  if (p != NULL)
54  {
55    fprintf(stderr, " occured 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 ***************************************************************/
69BOOLEAN 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    {
89      _pPolyAssumeReturn(omIsBinPageAddr(p) && omSizeWOfAddr(p)==omSizeWOfBin(r->PolyBin),p,r);
90      return TRUE;
91    }
92    return FALSE;
93  }
94  return TRUE;
95}
96
97BOOLEAN p_CheckIsFromRing(poly p, ring r)
98{
99  while (p!=NULL)
100  {
101    pFalseReturn(p_LmCheckIsFromRing(p, r));
102    pIter(p);
103  }
104  return TRUE;
105}
106
107BOOLEAN p_CheckPolyRing(poly p, ring r)
108{
109  pAssumeReturn(r != NULL && r->PolyBin != NULL);
110  return p_CheckIsFromRing(p, r);
111}
112
113BOOLEAN p_LmCheckPolyRing(poly p, ring r)
114{
115  pAssumeReturn(r != NULL && r->PolyBin != NULL);
116  pAssumeReturn(p != NULL);
117  return p_LmCheckIsFromRing(p, r);
118}
119BOOLEAN p_CheckRing(ring r)
120{
121  pAssumeReturn(r != NULL && r->PolyBin != NULL);
122  return TRUE;
123}
124
125/***************************************************************
126 *
127 * Debugging/statistics of pDivisibleBy
128 *
129 ***************************************************************/
130BOOLEAN p_DebugLmDivisibleByNoComp(poly a, poly b, ring r)
131{
132  int i=r->N;
133
134  do
135  {
136    if (p_GetExp(a,i,r) > p_GetExp(b,i,r))
137      return FALSE;
138    i--;
139  }
140  while (i);
141#ifdef HAVE_RINGS
142  return n_DivBy(pGetCoeff(b), pGetCoeff(a), r->cf);
143#else
144  return TRUE;
145#endif
146     }
147
148
149/***************************************************************
150 *
151 * Misc things helpful for debugging
152 *
153 ***************************************************************/
154BOOLEAN pIsMonomOf(poly p, poly m)
155{
156  if (m == NULL) return TRUE;
157  while (p != NULL)
158  {
159    if (p == m) return TRUE;
160    pIter(p);
161  }
162  return FALSE;
163}
164BOOLEAN pHaveCommonMonoms(poly p, poly q)
165{
166  while (p != NULL)
167  {
168    if (pIsMonomOf(q, p))
169    {
170      return TRUE;
171    }
172    pIter(p);
173  }
174  return FALSE;
175}
176
177/***************************************************************
178 *
179 * Testing of polys
180 *
181 ***************************************************************/
182extern void p_Setm_General(poly p, ring r);
183
184static poly p_DebugInit(poly p, ring src_ring, ring dest_ring)
185{
186  poly d_p = p_Init(dest_ring);
187  int i;
188  assume(dest_ring->N == src_ring->N);
189
190  for (i=1; i<= src_ring->N; i++)
191  {
192    p_SetExp(d_p, i, p_GetExp(p, i, src_ring), dest_ring);
193  }
194  if (rRing_has_Comp(dest_ring))
195    p_SetComp(d_p, p_GetComp(p, src_ring), dest_ring);
196
197  p_Setm_General(d_p, dest_ring);
198  return d_p;
199}
200
201BOOLEAN _p_Test(poly p, ring r, int level)
202{
203  assume(r->cf !=NULL);
204
205  if (PDEBUG > level) level = PDEBUG;
206  if (level < 0 || p == NULL) return TRUE;
207
208  poly p_prev = NULL;
209
210  #ifndef OM_NDEBUG
211  // check addr with level+1 so as to check bin/page of addr
212  _pPolyAssumeReturnMsg(omTestBinAddrSize(p, (omSizeWOfBin(r->PolyBin))*SIZEOF_LONG, level+1)
213                        == omError_NoError, "memory error",p,r);
214  #endif
215
216  pFalseReturn(p_CheckRing(r));
217
218  // this checks that p does not contain a loop: rather expensive O(length^2)
219  #ifndef OM_NDEBUG
220  if (level > 1)
221    pFalseReturn(omTestList(p, level) == omError_NoError);
222  #endif
223
224  int ismod = p_GetComp(p, r) != 0;
225
226  while (p != NULL)
227  {
228    // ring check
229    pFalseReturn(p_LmCheckIsFromRing(p, r));
230    #ifndef OM_NDEBUG
231    // omAddr check
232    _pPolyAssumeReturnMsg(omTestBinAddrSize(p, (omSizeWOfBin(r->PolyBin))*SIZEOF_LONG, 1)
233                     == omError_NoError, "memory error",p,r);
234    #endif
235    // number/coef check
236    _pPolyAssumeReturnMsg(p->coef != NULL || (n_GetChar(r->cf) >= 2), "NULL coef",p,r);
237    #ifdef LDEBUG
238    r->cf->cfDBTest(p->coef,__FILE__,__LINE__,r->cf);
239    #endif
240    _pPolyAssumeReturnMsg(!n_IsZero(p->coef, r->cf), "Zero coef",p,r);
241
242    // check for valid comp
243    _pPolyAssumeReturnMsg(p_GetComp(p, r) >= 0 && (p_GetComp(p, r)<65000), "component out of range ?",p,r);
244    // check for mix poly/vec representation
245    _pPolyAssumeReturnMsg(ismod == (p_GetComp(p, r) != 0), "mixed poly/vector",p,r);
246
247    // special check for ringorder_s/S
248    if ((r->typ!=NULL) && (r->typ[0].ord_typ == ro_syzcomp))
249    {
250      long c1, cc1, ccc1, ec1;
251      sro_ord* o = &(r->typ[0]);
252
253      c1 = p_GetComp(p, r);
254      if (o->data.syzcomp.Components!=NULL)
255      {
256        cc1 = o->data.syzcomp.Components[c1];
257        ccc1 = o->data.syzcomp.ShiftedComponents[cc1];
258      }
259      else { cc1=0; ccc1=0; }
260      _pPolyAssumeReturnMsg(c1 == 0 || cc1 != 0, "Component <-> TrueComponent zero mismatch",p,r);
261      _pPolyAssumeReturnMsg(c1 == 0 || ccc1 != 0,"Component <-> ShiftedComponent zero mismatch",p,r);
262      ec1 = p->exp[o->data.syzcomp.place];
263      //pPolyAssumeReturnMsg(ec1 == ccc1, "Shifted comp out of sync. should %d, is %d");
264      if (ec1 != ccc1)
265      {
266        dPolyReportError(p,r,"Shifted comp out of sync. should %d, is %d",ccc1,ec1);
267        return FALSE;
268      }
269    }
270
271    // check that p_Setm works ok
272    if (level > 0)
273    {
274      poly p_should_equal = p_DebugInit(p, r, r);
275      _pPolyAssumeReturnMsg(p_ExpVectorEqual(p, p_should_equal, r), "p_Setm field(s) out of sync",p,r);
276      p_LmFree(p_should_equal, r);
277    }
278
279    // check order
280    if (p_prev != NULL)
281    {
282      int cmp = p_LmCmp(p_prev, p, r);
283      if (cmp == 0)
284      {
285        _pPolyAssumeReturnMsg(0, "monoms p and p->next are equal", p_prev, r);
286      }
287      else
288        _pPolyAssumeReturnMsg(p_LmCmp(p_prev, p, r) == 1, "wrong order", p_prev, r);
289
290      // check that compare worked sensibly
291      if (level > 1 && p_GetComp(p_prev, r) == p_GetComp(p, r))
292      {
293        int i;
294        for (i=r->N; i>0; i--)
295        {
296          if (p_GetExp(p_prev, i, r) != p_GetExp(p, i, r)) break;
297        }
298        _pPolyAssumeReturnMsg(i > 0, "Exponents equal but compare different", p_prev, r);
299      }
300    }
301    p_prev = p;
302    pIter(p);
303  }
304  return TRUE;
305}
306
307BOOLEAN _p_LmTest(poly p, ring r, int level)
308{
309  if (level < 0 || p == NULL) return TRUE;
310  poly pnext = pNext(p);
311  pNext(p) = NULL;
312  BOOLEAN test_res = _p_Test(p, r, level);
313  pNext(p) = pnext;
314  return test_res;
315}
316
317BOOLEAN _pp_Test(poly p, ring lmRing, ring tailRing, int level)
318{
319  if (PDEBUG > level) level = PDEBUG;
320  if (level < 0 || p == NULL) return TRUE;
321  if (pNext(p) == NULL || lmRing == tailRing) return _p_Test(p, lmRing, level);
322
323  pFalseReturn(_p_LmTest(p, lmRing, level));
324  pFalseReturn(_p_Test(pNext(p), tailRing, level));
325
326  // check that lm > Lm(tail)
327  if (level > 1)
328  {
329    poly lm = p;
330    poly tail = p_DebugInit(pNext(p), tailRing, lmRing);
331    poly pnext = pNext(lm);
332    pNext(lm) = tail;
333    BOOLEAN cmp = p_LmCmp(lm, tail, lmRing);
334    if (cmp != 1)
335      dPolyReportError(lm, lmRing, "wrong order: lm <= Lm(tail)");
336    p_LmFree(tail, lmRing);
337    pNext(lm) = pnext;
338    return (cmp == 1);
339  }
340  return TRUE;
341}
342
343#endif // PDEBUG
344
345#if defined(PDEBUG) || defined(PDIV_DEBUG)
346static unsigned long pDivisibleBy_number = 1;
347static unsigned long pDivisibleBy_FALSE = 1;
348static unsigned long pDivisibleBy_ShortFalse = 1;
349
350BOOLEAN pDebugLmShortDivisibleBy(poly p1, unsigned long sev_1, ring r_1,
351                               poly p2, unsigned long not_sev_2, ring r_2)
352{
353  _pPolyAssume(p_GetShortExpVector(p1, r_1) == sev_1, p1, r_1);
354  _pPolyAssume(p_GetShortExpVector(p2, r_2) == ~ not_sev_2, p2, r_2);
355
356  pDivisibleBy_number++;
357  BOOLEAN ret;
358  if (r_1 == r_2)
359    ret = p_LmDivisibleBy(p1, p2, r_1);
360  else
361    ret = p_LmDivisibleBy(p1, r_1, p2, r_2);
362
363  if (! ret) pDivisibleBy_FALSE++;
364  if (sev_1 & not_sev_2)
365  {
366    pDivisibleBy_ShortFalse++;
367    if (ret)
368      dReportError("p1 divides p2, but sev's are wrong");
369  }
370  return ret;
371}
372
373BOOLEAN pDebugLmShortDivisibleByNoComp(poly p1, unsigned long sev_1, ring r_1,
374                                       poly p2, unsigned long not_sev_2, ring r_2)
375{
376//  _pPolyAssume((p_GetComp(p1, r_1) == p_GetComp(p2, r_2)) || (p_GetComp(p1, r_1) == 0));
377  _pPolyAssume(p_GetShortExpVector(p1, r_1) == sev_1, p1, r_1);
378  _pPolyAssume(p_GetShortExpVector(p2, r_2) == ~ not_sev_2, p2, r_2);
379
380  pDivisibleBy_number++;
381  BOOLEAN ret;
382  if (r_1 == r_2)
383    ret = p_LmDivisibleByNoComp(p1, p2, r_1);
384  else
385    ret = p_LmDivisibleByNoComp(p1, r_1, p2, r_2);
386
387  if (! ret) pDivisibleBy_FALSE++;
388  if (sev_1 & not_sev_2)
389  {
390    pDivisibleBy_ShortFalse++;
391    if (ret)
392      dReportError("p1 divides p2, but sev's are wrong");
393  }
394  return ret;
395}
396
397void pPrintDivisbleByStat()
398{
399  Print("#Tests: %ld; #FALSE %ld(%ld); #SHORT %ld(%ld)\n",
400        pDivisibleBy_number,
401        pDivisibleBy_FALSE, (unsigned long) ((double)pDivisibleBy_FALSE*((double) 100)/(double)pDivisibleBy_number),
402        pDivisibleBy_ShortFalse, (unsigned long) ((double)pDivisibleBy_ShortFalse*((double)100)/(double)pDivisibleBy_FALSE));
403}
404
405#endif // PDEBUG
406
407#endif // PDEBUG_CC
Note: See TracBrowser for help on using the repository browser.