My Project  debian-1:4.1.1-p2+ds-4build2
ap.h
Go to the documentation of this file.
1 /********************************************************************
2 AP library.
3 See www.alglib.net or alglib.sources.ru for details.
4 ********************************************************************/
5 
6 #ifndef AP_H
7 #define AP_H
8 
9 #include <stdlib.h>
10 #include <math.h>
11 
12 /********************************************************************
13 Checking of the array boundaries mode.
14 ********************************************************************/
15 //#define NO_AP_ASSERT
16 #define AP_ASSERT
17 
18 #ifndef AP_ASSERT
19 #define NO_AP_ASSERT
20 #endif
21 
22 #ifdef NO_AP_ASSERT
23 #ifdef AP_ASSERT
24 #undef NO_AP_ASSERT
25 #endif
26 #endif
27 
28 
29 /********************************************************************
30 This symbol is used for debugging. Do not define it and do not remove
31 comments.
32 ********************************************************************/
33 //#define UNSAFE_MEM_COPY
34 
35 
36 /********************************************************************
37 Namespace of a standard library AlgoPascal.
38 ********************************************************************/
39 namespace ap
40 {
41 
42 
43 /********************************************************************
44 Exception class.
45 ********************************************************************/
46 class ap_error
47 {
48 public:
49  static void make_assertion(bool bClause)
50  { if(!bClause) throw ap_error(); };
51 private:
52 };
53 
54 /********************************************************************
55 Class defining a complex number with double precision.
56 ********************************************************************/
57 class complex;
58 
59 class complex
60 {
61 public:
62  complex():x(0.0),y(0.0){};
63  complex(const double &_x):x(_x),y(0.0){};
64  complex(const double &_x, const double &_y):x(_x),y(_y){};
65  complex(const complex &z):x(z.x),y(z.y){};
66 
67  complex& operator= (const double& v){ x = v; y = 0.0; return *this; };
68  complex& operator+=(const double& v){ x += v; return *this; };
69  complex& operator-=(const double& v){ x -= v; return *this; };
70  complex& operator*=(const double& v){ x *= v; y *= v; return *this; };
71  complex& operator/=(const double& v){ x /= v; y /= v; return *this; };
72 
73  complex& operator= (const complex& z){ x = z.x; y = z.y; return *this; };
74  complex& operator+=(const complex& z){ x += z.x; y += z.y; return *this; };
75  complex& operator-=(const complex& z){ x -= z.x; y -= z.y; return *this; };
76  complex& operator*=(const complex& z){ double t = x*z.x-y*z.y; y = x*z.y+y*z.x; x = t; return *this; };
77  complex& operator/=(const complex& z)
78  {
80  double e;
81  double f;
82  if( fabs(z.y)<fabs(z.x) )
83  {
84  e = z.y/z.x;
85  f = z.x+z.y*e;
86  result.x = (z.x+z.y*e)/f;
87  result.y = (z.y-z.x*e)/f;
88  }
89  else
90  {
91  e = z.x/z.y;
92  f = z.y+z.x*e;
93  result.x = (z.y+z.x*e)/f;
94  result.y = (-z.x+z.y*e)/f;
95  }
96  *this = result;
97  return *this;
98  };
99 
100  double x, y;
101 };
102 
103 const complex operator/(const complex& lhs, const complex& rhs);
104 const bool operator==(const complex& lhs, const complex& rhs);
105 const bool operator!=(const complex& lhs, const complex& rhs);
106 const complex operator+(const complex& lhs);
107 const complex operator-(const complex& lhs);
108 const complex operator+(const complex& lhs, const complex& rhs);
109 const complex operator+(const complex& lhs, const double& rhs);
110 const complex operator+(const double& lhs, const complex& rhs);
111 const complex operator-(const complex& lhs, const complex& rhs);
112 const complex operator-(const complex& lhs, const double& rhs);
113 const complex operator-(const double& lhs, const complex& rhs);
114 const complex operator*(const complex& lhs, const complex& rhs);
115 const complex operator*(const complex& lhs, const double& rhs);
116 const complex operator*(const double& lhs, const complex& rhs);
117 const complex operator/(const complex& lhs, const complex& rhs);
118 const complex operator/(const double& lhs, const complex& rhs);
119 const complex operator/(const complex& lhs, const double& rhs);
120 const double abscomplex(const complex &z);
121 const complex conj(const complex &z);
122 const complex csqr(const complex &z);
123 
124 
125 /********************************************************************
126 Template defining vector in memory. It is used by the basic
127 subroutines of linear algebra.
128 
129 Vector consists of Length elements of type T, starting from an element,
130 which Data is pointed to. Interval between adjacent elements equals
131 the value of Step.
132 
133 The class provides an access for reading only.
134 ********************************************************************/
135 template<class T>
136 class const_raw_vector
137 {
138 public:
139  const_raw_vector(const T *Data, int Length, int Step):
140  pData(const_cast<T*>(Data)),iLength(Length),iStep(Step){};
141 
142  const T* GetData() const
143  { return pData; };
144 
145  int GetLength() const
146  { return iLength; };
147 
148  int GetStep() const
149  { return iStep; };
150 protected:
151  T *pData;
152  int iLength, iStep;
153 };
154 
155 
156 /********************************************************************
157 Template defining vector in memory, derived from const_raw_vector.
158 It is used by the basic subroutines of linear algebra.
159 
160 Vector consists of Length elements of type T, starting from an element,
161 which Data is pointed to. Interval between adjacent elements equals
162 the value of Step.
163 
164 The class provides an access both for reading and writing.
165 ********************************************************************/
166 template<class T>
167 class raw_vector : public const_raw_vector<T>
168 {
169 public:
170  raw_vector(T *Data, int Length, int Step):const_raw_vector<T>(Data, Length, Step){};
171 
173  { return const_raw_vector<T>::pData; };
174 };
175 
176 
177 /********************************************************************
178 Scalar product
179 ********************************************************************/
180 template<class T>
181 T vdotproduct(const_raw_vector<T> v1, const_raw_vector<T> v2)
182 {
183  ap_error::make_assertion(v1.GetLength()==v2.GetLength());
184  if( v1.GetStep()==1 && v2.GetStep()==1 )
185  {
186  //
187  // fast
188  //
189  T r = 0;
190  const T *p1 = v1.GetData();
191  const T *p2 = v2.GetData();
192  int imax = v1.GetLength()/4;
193  int i;
194  for(i=imax; i!=0; i--)
195  {
196  r += (*p1)*(*p2) + p1[1]*p2[1] + p1[2]*p2[2] + p1[3]*p2[3];
197  p1+=4;
198  p2+=4;
199  }
200  for(i=0; i<v1.GetLength()%4; i++)
201  r += (*(p1++))*(*(p2++));
202  return r;
203  }
204  else
205  {
206  //
207  // general
208  //
209  int offset11 = v1.GetStep(), offset12 = 2*offset11, offset13 = 3*offset11, offset14 = 4*offset11;
210  int offset21 = v2.GetStep(), offset22 = 2*offset21, offset23 = 3*offset21, offset24 = 4*offset21;
211  T r = 0;
212  const T *p1 = v1.GetData();
213  const T *p2 = v2.GetData();
214  int imax = v1.GetLength()/4;
215  int i;
216  for(i=0; i<imax; i++)
217  {
218  r += (*p1)*(*p2) + p1[offset11]*p2[offset21] + p1[offset12]*p2[offset22] + p1[offset13]*p2[offset23];
219  p1+=offset14;
220  p2+=offset24;
221  }
222  for(i=0; i<v1.GetLength()%4; i++)
223  {
224  r += (*p1)*(*p2);
225  p1+=offset11;
226  p2+=offset21;
227  }
228  return r;
229  }
230 }
231 
232 
233 /********************************************************************
234 Copy one vector into another
235 ********************************************************************/
236 template<class T>
237 void vmove(raw_vector<T> vdst, const_raw_vector<T> vsrc)
238 {
239  ap_error::make_assertion(vdst.GetLength()==vsrc.GetLength());
240  if( vdst.GetStep()==1 && vsrc.GetStep()==1 )
241  {
242  //
243  // fast
244  //
245  T *p1 = vdst.GetData();
246  const T *p2 = vsrc.GetData();
247  int imax = vdst.GetLength()/2;
248  int i;
249  for(i=imax; i!=0; i--)
250  {
251  *p1 = *p2;
252  p1[1] = p2[1];
253  p1 += 2;
254  p2 += 2;
255  }
256  if(vdst.GetLength()%2 != 0)
257  *p1 = *p2;
258  return;
259  }
260  else
261  {
262  //
263  // general
264  //
265  int offset11 = vdst.GetStep(), offset12 = 2*offset11, offset13 = 3*offset11, offset14 = 4*offset11;
266  int offset21 = vsrc.GetStep(), offset22 = 2*offset21, offset23 = 3*offset21, offset24 = 4*offset21;
267  T *p1 = vdst.GetData();
268  const T *p2 = vsrc.GetData();
269  int imax = vdst.GetLength()/4;
270  int i;
271  for(i=0; i<imax; i++)
272  {
273  *p1 = *p2;
274  p1[offset11] = p2[offset21];
275  p1[offset12] = p2[offset22];
276  p1[offset13] = p2[offset23];
277  p1 += offset14;
278  p2 += offset24;
279  }
280  for(i=0; i<vdst.GetLength()%4; i++)
281  {
282  *p1 = *p2;
283  p1 += vdst.GetStep();
284  p2 += vsrc.GetStep();
285  }
286  return;
287  }
288 }
289 
290 
291 /********************************************************************
292 Copy one vector multiplied by -1 into another.
293 ********************************************************************/
294 template<class T>
295 void vmoveneg(raw_vector<T> vdst, const_raw_vector<T> vsrc)
296 {
297  ap_error::make_assertion(vdst.GetLength()==vsrc.GetLength());
298  if( vdst.GetStep()==1 && vsrc.GetStep()==1 )
299  {
300  //
301  // fast
302  //
303  T *p1 = vdst.GetData();
304  const T *p2 = vsrc.GetData();
305  int imax = vdst.GetLength()/2;
306  int i;
307  for(i=0; i<imax; i++)
308  {
309  *p1 = -*p2;
310  p1[1] = -p2[1];
311  p1 += 2;
312  p2 += 2;
313  }
314  if(vdst.GetLength()%2 != 0)
315  *p1 = -*p2;
316  return;
317  }
318  else
319  {
320  //
321  // general
322  //
323  int offset11 = vdst.GetStep(), offset12 = 2*offset11, offset13 = 3*offset11, offset14 = 4*offset11;
324  int offset21 = vsrc.GetStep(), offset22 = 2*offset21, offset23 = 3*offset21, offset24 = 4*offset21;
325  T *p1 = vdst.GetData();
326  const T *p2 = vsrc.GetData();
327  int imax = vdst.GetLength()/4;
328  int i;
329  for(i=imax; i!=0; i--)
330  {
331  *p1 = -*p2;
332  p1[offset11] = -p2[offset21];
333  p1[offset12] = -p2[offset22];
334  p1[offset13] = -p2[offset23];
335  p1 += offset14;
336  p2 += offset24;
337  }
338  for(i=0; i<vdst.GetLength()%4; i++)
339  {
340  *p1 = -*p2;
341  p1 += vdst.GetStep();
342  p2 += vsrc.GetStep();
343  }
344  return;
345  }
346 }
347 
348 
349 /********************************************************************
350 Copy one vector multiplied by a number into another vector.
351 ********************************************************************/
352 template<class T, class T2>
353 void vmove(raw_vector<T> vdst, const_raw_vector<T> vsrc, T2 alpha)
354 {
355  ap_error::make_assertion(vdst.GetLength()==vsrc.GetLength());
356  if( vdst.GetStep()==1 && vsrc.GetStep()==1 )
357  {
358  //
359  // fast
360  //
361  T *p1 = vdst.GetData();
362  const T *p2 = vsrc.GetData();
363  int imax = vdst.GetLength()/4;
364  int i;
365  for(i=imax; i!=0; i--)
366  {
367  *p1 = alpha*(*p2);
368  p1[1] = alpha*p2[1];
369  p1[2] = alpha*p2[2];
370  p1[3] = alpha*p2[3];
371  p1 += 4;
372  p2 += 4;
373  }
374  for(i=0; i<vdst.GetLength()%4; i++)
375  *(p1++) = alpha*(*(p2++));
376  return;
377  }
378  else
379  {
380  //
381  // general
382  //
383  int offset11 = vdst.GetStep(), offset12 = 2*offset11, offset13 = 3*offset11, offset14 = 4*offset11;
384  int offset21 = vsrc.GetStep(), offset22 = 2*offset21, offset23 = 3*offset21, offset24 = 4*offset21;
385  T *p1 = vdst.GetData();
386  const T *p2 = vsrc.GetData();
387  int imax = vdst.GetLength()/4;
388  int i;
389  for(i=0; i<imax; i++)
390  {
391  *p1 = alpha*(*p2);
392  p1[offset11] = alpha*p2[offset21];
393  p1[offset12] = alpha*p2[offset22];
394  p1[offset13] = alpha*p2[offset23];
395  p1 += offset14;
396  p2 += offset24;
397  }
398  for(i=0; i<vdst.GetLength()%4; i++)
399  {
400  *p1 = alpha*(*p2);
401  p1 += vdst.GetStep();
402  p2 += vsrc.GetStep();
403  }
404  return;
405  }
406 }
407 
408 
409 /********************************************************************
410 Vector addition
411 ********************************************************************/
412 template<class T>
413 void vadd(raw_vector<T> vdst, const_raw_vector<T> vsrc)
414 {
415  ap_error::make_assertion(vdst.GetLength()==vsrc.GetLength());
416  if( vdst.GetStep()==1 && vsrc.GetStep()==1 )
417  {
418  //
419  // fast
420  //
421  T *p1 = vdst.GetData();
422  const T *p2 = vsrc.GetData();
423  int imax = vdst.GetLength()/4;
424  int i;
425  for(i=imax; i!=0; i--)
426  {
427  *p1 += *p2;
428  p1[1] += p2[1];
429  p1[2] += p2[2];
430  p1[3] += p2[3];
431  p1 += 4;
432  p2 += 4;
433  }
434  for(i=0; i<vdst.GetLength()%4; i++)
435  *(p1++) += *(p2++);
436  return;
437  }
438  else
439  {
440  //
441  // general
442  //
443  int offset11 = vdst.GetStep(), offset12 = 2*offset11, offset13 = 3*offset11, offset14 = 4*offset11;
444  int offset21 = vsrc.GetStep(), offset22 = 2*offset21, offset23 = 3*offset21, offset24 = 4*offset21;
445  T *p1 = vdst.GetData();
446  const T *p2 = vsrc.GetData();
447  int imax = vdst.GetLength()/4;
448  int i;
449  for(i=0; i<imax; i++)
450  {
451  *p1 += *p2;
452  p1[offset11] += p2[offset21];
453  p1[offset12] += p2[offset22];
454  p1[offset13] += p2[offset23];
455  p1 += offset14;
456  p2 += offset24;
457  }
458  for(i=0; i<vdst.GetLength()%4; i++)
459  {
460  *p1 += *p2;
461  p1 += vdst.GetStep();
462  p2 += vsrc.GetStep();
463  }
464  return;
465  }
466 }
467 
468 
469 /********************************************************************
470 Add one vector multiplied by a number to another vector.
471 ********************************************************************/
472 template<class T, class T2>
473 void vadd(raw_vector<T> vdst, const_raw_vector<T> vsrc, T2 alpha)
474 {
475  ap_error::make_assertion(vdst.GetLength()==vsrc.GetLength());
476  if( vdst.GetStep()==1 && vsrc.GetStep()==1 )
477  {
478  //
479  // fast
480  //
481  T *p1 = vdst.GetData();
482  const T *p2 = vsrc.GetData();
483  int imax = vdst.GetLength()/4;
484  int i;
485  for(i=imax; i!=0; i--)
486  {
487  *p1 += alpha*(*p2);
488  p1[1] += alpha*p2[1];
489  p1[2] += alpha*p2[2];
490  p1[3] += alpha*p2[3];
491  p1 += 4;
492  p2 += 4;
493  }
494  for(i=0; i<vdst.GetLength()%4; i++)
495  *(p1++) += alpha*(*(p2++));
496  return;
497  }
498  else
499  {
500  //
501  // general
502  //
503  int offset11 = vdst.GetStep(), offset12 = 2*offset11, offset13 = 3*offset11, offset14 = 4*offset11;
504  int offset21 = vsrc.GetStep(), offset22 = 2*offset21, offset23 = 3*offset21, offset24 = 4*offset21;
505  T *p1 = vdst.GetData();
506  const T *p2 = vsrc.GetData();
507  int imax = vdst.GetLength()/4;
508  int i;
509  for(i=0; i<imax; i++)
510  {
511  *p1 += alpha*(*p2);
512  p1[offset11] += alpha*p2[offset21];
513  p1[offset12] += alpha*p2[offset22];
514  p1[offset13] += alpha*p2[offset23];
515  p1 += offset14;
516  p2 += offset24;
517  }
518  for(i=0; i<vdst.GetLength()%4; i++)
519  {
520  *p1 += alpha*(*p2);
521  p1 += vdst.GetStep();
522  p2 += vsrc.GetStep();
523  }
524  return;
525  }
526 }
527 
528 
529 /********************************************************************
530 Vector subtraction
531 ********************************************************************/
532 template<class T>
533 void vsub(raw_vector<T> vdst, const_raw_vector<T> vsrc)
534 {
535  ap_error::make_assertion(vdst.GetLength()==vsrc.GetLength());
536  if( vdst.GetStep()==1 && vsrc.GetStep()==1 )
537  {
538  //
539  // fast
540  //
541  T *p1 = vdst.GetData();
542  const T *p2 = vsrc.GetData();
543  int imax = vdst.GetLength()/4;
544  int i;
545  for(i=imax; i!=0; i--)
546  {
547  *p1 -= *p2;
548  p1[1] -= p2[1];
549  p1[2] -= p2[2];
550  p1[3] -= p2[3];
551  p1 += 4;
552  p2 += 4;
553  }
554  for(i=0; i<vdst.GetLength()%4; i++)
555  *(p1++) -= *(p2++);
556  return;
557  }
558  else
559  {
560  //
561  // general
562  //
563  int offset11 = vdst.GetStep(), offset12 = 2*offset11, offset13 = 3*offset11, offset14 = 4*offset11;
564  int offset21 = vsrc.GetStep(), offset22 = 2*offset21, offset23 = 3*offset21, offset24 = 4*offset21;
565  T *p1 = vdst.GetData();
566  const T *p2 = vsrc.GetData();
567  int imax = vdst.GetLength()/4;
568  int i;
569  for(i=0; i<imax; i++)
570  {
571  *p1 -= *p2;
572  p1[offset11] -= p2[offset21];
573  p1[offset12] -= p2[offset22];
574  p1[offset13] -= p2[offset23];
575  p1 += offset14;
576  p2 += offset24;
577  }
578  for(i=0; i<vdst.GetLength()%4; i++)
579  {
580  *p1 -= *p2;
581  p1 += vdst.GetStep();
582  p2 += vsrc.GetStep();
583  }
584  return;
585  }
586 }
587 
588 
589 /********************************************************************
590 Subtract one vector multiplied by a number from another vector.
591 ********************************************************************/
592 template<class T, class T2>
593 void vsub(raw_vector<T> vdst, const_raw_vector<T> vsrc, T2 alpha)
594 {
595  vadd(vdst, vsrc, -alpha);
596 }
597 
598 
599 /********************************************************************
600 In-place vector multiplication
601 ********************************************************************/
602 template<class T, class T2>
603 void vmul(raw_vector<T> vdst, T2 alpha)
604 {
605  if( vdst.GetStep()==1 )
606  {
607  //
608  // fast
609  //
610  T *p1 = vdst.GetData();
611  int imax = vdst.GetLength()/4;
612  int i;
613  for(i=imax; i!=0; i--)
614  {
615  *p1 *= alpha;
616  p1[1] *= alpha;
617  p1[2] *= alpha;
618  p1[3] *= alpha;
619  p1 += 4;
620  }
621  for(i=0; i<vdst.GetLength()%4; i++)
622  *(p1++) *= alpha;
623  return;
624  }
625  else
626  {
627  //
628  // general
629  //
630  int offset11 = vdst.GetStep(), offset12 = 2*offset11, offset13 = 3*offset11, offset14 = 4*offset11;
631  T *p1 = vdst.GetData();
632  int imax = vdst.GetLength()/4;
633  int i;
634  for(i=0; i<imax; i++)
635  {
636  *p1 *= alpha;
637  p1[offset11] *= alpha;
638  p1[offset12] *= alpha;
639  p1[offset13] *= alpha;
640  p1 += offset14;
641  }
642  for(i=0; i<vdst.GetLength()%4; i++)
643  {
644  *p1 *= alpha;
645  p1 += vdst.GetStep();
646  }
647  return;
648  }
649 }
650 
651 
652 /********************************************************************
653 Template of a dynamical one-dimensional array
654 ********************************************************************/
655 template<class T>
656 class template_1d_array
657 {
658 public:
660  {
661  m_Vec=0;
662  m_iVecSize = 0;
663  };
664 
666  {
667  if(m_Vec)
668  delete[] m_Vec;
669  };
670 
672  {
673  m_iVecSize = rhs.m_iVecSize;
674  m_iLow = rhs.m_iLow;
675  m_iHigh = rhs.m_iHigh;
676  if(rhs.m_Vec)
677  {
679  #ifndef UNSAFE_MEM_COPY
680  for(int i=0; i<m_iVecSize; i++)
681  m_Vec[i] = rhs.m_Vec[i];
682  #else
683  memcpy(m_Vec, rhs.m_Vec, m_iVecSize*sizeof(T));
684  #endif
685  }
686  else
687  m_Vec=0;
688  };
689 
690 
692  {
693  if( this==&rhs )
694  return *this;
695 
696  m_iLow = rhs.m_iLow;
697  m_iHigh = rhs.m_iHigh;
698  m_iVecSize = rhs.m_iVecSize;
699  if(m_Vec)
700  delete[] m_Vec;
701  if(rhs.m_Vec)
702  {
703  m_Vec = new T[m_iVecSize];
704  #ifndef UNSAFE_MEM_COPY
705  for(int i=0; i<m_iVecSize; i++)
706  m_Vec[i] = rhs.m_Vec[i];
707  #else
708  memcpy(m_Vec, rhs.m_Vec, m_iVecSize*sizeof(T));
709  #endif
710  }
711  else
712  m_Vec=0;
713  return *this;
714  };
715 
716 
717  const T& operator()(int i) const
718  {
719  #ifndef NO_AP_ASSERT
721  #endif
722  return m_Vec[ i-m_iLow ];
723  };
724 
725 
726  T& operator()(int i)
727  {
728  #ifndef NO_AP_ASSERT
730  #endif
731  return m_Vec[ i-m_iLow ];
732  };
733 
734 
735  void setbounds( int iLow, int iHigh )
736  {
737  if(m_Vec)
738  delete[] m_Vec;
739  m_iLow = iLow;
740  m_iHigh = iHigh;
741  m_iVecSize = iHigh-iLow+1;
743  };
744 
745 
746  void setcontent( int iLow, int iHigh, const T *pContent )
747  {
748  setbounds(iLow, iHigh);
749  for(int i=iLow; i<=iHigh; i++)
750  (*this)(i) = pContent[i-iLow];
751  };
752 
753 
754  T* getcontent()
755  {
756  return m_Vec;
757  };
758 
759  const T* getcontent() const
760  {
761  return m_Vec;
762  };
763 
764 
765  int getlowbound(int iBoundNum = 0) const
766  {
767  return m_iLow;
768  };
769 
770 
771  int gethighbound(int iBoundNum = 0) const
772  {
773  return m_iHigh;
774  };
775 
776  raw_vector<T> getvector(int iStart, int iEnd)
777  {
778  if( iStart>iEnd || wrongIdx(iStart) || wrongIdx(iEnd) )
779  return raw_vector<T>(0, 0, 1);
780  else
781  return raw_vector<T>(m_Vec+iStart-m_iLow, iEnd-iStart+1, 1);
782  };
783 
784 
785  const_raw_vector<T> getvector(int iStart, int iEnd) const
786  {
787  if( iStart>iEnd || wrongIdx(iStart) || wrongIdx(iEnd) )
788  return const_raw_vector<T>(0, 0, 1);
789  else
790  return const_raw_vector<T>(m_Vec+iStart-m_iLow, iEnd-iStart+1, 1);
791  };
792 private:
793  bool wrongIdx(int i) const { return i<m_iLow || i>m_iHigh; };
794 
795  T *m_Vec;
796  long m_iVecSize;
797  long m_iLow, m_iHigh;
798 };
799 
800 
801 
802 /********************************************************************
803 Template of a dynamical two-dimensional array
804 ********************************************************************/
805 template<class T>
806 class template_2d_array
807 {
808 public:
810  {
811  m_Vec=0;
812  m_iVecSize=0;
813  };
814 
816  {
817  if(m_Vec)
818  delete[] m_Vec;
819  };
820 
822  {
823  m_iVecSize = rhs.m_iVecSize;
824  m_iLow1 = rhs.m_iLow1;
825  m_iLow2 = rhs.m_iLow2;
826  m_iHigh1 = rhs.m_iHigh1;
827  m_iHigh2 = rhs.m_iHigh2;
830  if(rhs.m_Vec)
831  {
832  m_Vec = new T[m_iVecSize];
833  #ifndef UNSAFE_MEM_COPY
834  for(int i=0; i<m_iVecSize; i++)
835  m_Vec[i] = rhs.m_Vec[i];
836  #else
837  memcpy(m_Vec, rhs.m_Vec, m_iVecSize*sizeof(T));
838  #endif
839  }
840  else
841  m_Vec=0;
842  };
844  {
845  if( this==&rhs )
846  return *this;
847 
848  m_iLow1 = rhs.m_iLow1;
849  m_iLow2 = rhs.m_iLow2;
850  m_iHigh1 = rhs.m_iHigh1;
851  m_iHigh2 = rhs.m_iHigh2;
852  m_iConstOffset = rhs.m_iConstOffset;
853  m_iLinearMember = rhs.m_iLinearMember;
854  m_iVecSize = rhs.m_iVecSize;
855  if(m_Vec)
856  delete[] m_Vec;
857  if(rhs.m_Vec)
858  {
859  m_Vec = new T[m_iVecSize];
860  #ifndef UNSAFE_MEM_COPY
861  for(int i=0; i<m_iVecSize; i++)
862  m_Vec[i] = rhs.m_Vec[i];
863  #else
864  memcpy(m_Vec, rhs.m_Vec, m_iVecSize*sizeof(T));
865  #endif
866  }
867  else
868  m_Vec=0;
869  return *this;
870  };
871 
872  const T& operator()(int i1, int i2) const
873  {
874  #ifndef NO_AP_ASSERT
877  #endif
878  return m_Vec[ m_iConstOffset + i2 +i1*m_iLinearMember];
879  };
880 
881  T& operator()(int i1, int i2)
882  {
883  #ifndef NO_AP_ASSERT
886  #endif
887  return m_Vec[ m_iConstOffset + i2 +i1*m_iLinearMember];
888  };
889 
890  void setbounds( int iLow1, int iHigh1, int iLow2, int iHigh2 )
891  {
892  if(m_Vec)
893  delete[] m_Vec;
894  m_iVecSize = (iHigh1-iLow1+1)*(iHigh2-iLow2+1);
895  m_Vec = new T[m_iVecSize];
896  m_iLow1 = iLow1;
897  m_iHigh1 = iHigh1;
898  m_iLow2 = iLow2;
899  m_iHigh2 = iHigh2;
902  };
903 
904  void setcontent( int iLow1, int iHigh1, int iLow2, int iHigh2, const T *pContent )
905  {
906  setbounds(iLow1, iHigh1, iLow2, iHigh2);
907  for(int i=0; i<m_iVecSize; i++)
908  m_Vec[i]=pContent[i];
909  };
910 
912  {
913  return m_Vec;
914  };
915 
916  const T* getcontent() const
917  {
918  return m_Vec;
919  };
920 
921  int getlowbound(int iBoundNum) const
922  {
923  return iBoundNum==1 ? m_iLow1 : m_iLow2;
924  };
925 
926  int gethighbound(int iBoundNum) const
927  {
928  return iBoundNum==1 ? m_iHigh1 : m_iHigh2;
929  };
930 
931  raw_vector<T> getcolumn(int iColumn, int iRowStart, int iRowEnd)
932  {
933  if( (iRowStart>iRowEnd) || wrongColumn(iColumn) || wrongRow(iRowStart) ||wrongRow(iRowEnd) )
934  return raw_vector<T>(0, 0, 1);
935  else
936  return raw_vector<T>(&((*this)(iRowStart, iColumn)), iRowEnd-iRowStart+1, m_iLinearMember);
937  };
938 
939  raw_vector<T> getrow(int iRow, int iColumnStart, int iColumnEnd)
940  {
941  if( (iColumnStart>iColumnEnd) || wrongRow(iRow) || wrongColumn(iColumnStart) || wrongColumn(iColumnEnd))
942  return raw_vector<T>(0, 0, 1);
943  else
944  return raw_vector<T>(&((*this)(iRow, iColumnStart)), iColumnEnd-iColumnStart+1, 1);
945  };
946 
947  const_raw_vector<T> getcolumn(int iColumn, int iRowStart, int iRowEnd) const
948  {
949  if( (iRowStart>iRowEnd) || wrongColumn(iColumn) || wrongRow(iRowStart) ||wrongRow(iRowEnd) )
950  return const_raw_vector<T>(0, 0, 1);
951  else
952  return const_raw_vector<T>(&((*this)(iRowStart, iColumn)), iRowEnd-iRowStart+1, m_iLinearMember);
953  };
954 
955  const_raw_vector<T> getrow(int iRow, int iColumnStart, int iColumnEnd) const
956  {
957  if( (iColumnStart>iColumnEnd) || wrongRow(iRow) || wrongColumn(iColumnStart) || wrongColumn(iColumnEnd))
958  return const_raw_vector<T>(0, 0, 1);
959  else
960  return const_raw_vector<T>(&((*this)(iRow, iColumnStart)), iColumnEnd-iColumnStart+1, 1);
961  };
962 private:
963  bool wrongRow(int i) const { return i<m_iLow1 || i>m_iHigh1; };
964  bool wrongColumn(int j) const { return j<m_iLow2 || j>m_iHigh2; };
965 
966  T *m_Vec;
967  long m_iVecSize;
970 };
971 
972 
981 
982 
983 /********************************************************************
984 Constants and functions introduced for compatibility with AlgoPascal
985 ********************************************************************/
986 extern const double machineepsilon;
987 extern const double maxrealnumber;
988 extern const double minrealnumber;
989 
990 int sign(double x);
991 double randomreal();
992 int randominteger(int maxv);
993 int round(double x);
994 int trunc(double x);
995 int ifloor(double x);
996 int iceil(double x);
997 double pi();
998 double sqr(double x);
999 int maxint(int m1, int m2);
1000 int minint(int m1, int m2);
1001 double maxreal(double m1, double m2);
1002 double minreal(double m1, double m2);
1003 
1004 };//namespace ap
1005 
1006 
1007 #endif
ap::conj
const complex conj(const complex &z)
Definition: ap.cpp:117
ap::template_2d_array::getlowbound
int getlowbound(int iBoundNum) const
Definition: ap.h:928
ap::template_2d_array::getcontent
T * getcontent()
Definition: ap.h:918
ap::machineepsilon
const double machineepsilon
Definition: svd_si.h:999
ap::operator+
const complex operator+(const complex &lhs)
Definition: ap.cpp:17
ap::operator/
const complex operator/(const complex &lhs, const complex &rhs)
Definition: ap.cpp:50
ap::template_2d_array::template_2d_array
template_2d_array()
Definition: ap.h:816
ap::vmove
void vmove(raw_vector< T > vdst, const_raw_vector< T > vsrc)
Definition: ap.h:242
j
int j
Definition: facHensel.cc:105
f
FILE * f
Definition: checklibs.c:9
ap::operator==
const bool operator==(const complex &lhs, const complex &rhs)
Definition: ap.cpp:11
ap::round
int round(double x)
Definition: ap.cpp:144
x
Variable x
Definition: cfModGcd.cc:4023
ap::vadd
void vadd(raw_vector< T > vdst, const_raw_vector< T > vsrc)
Definition: ap.h:418
ap::template_1d_array::getvector
raw_vector< T > getvector(int iStart, int iEnd)
Definition: ap.h:783
result
return result
Definition: facAbsBiFact.cc:76
ap::template_1d_array::getcontent
T * getcontent()
Definition: ap.h:761
ap::template_2d_array::m_Vec
T * m_Vec
Definition: ap.h:971
ap::template_1d_array
Definition: ap.h:661
ap::template_2d_array
Definition: ap.h:811
ap::template_2d_array::m_iLinearMember
long m_iLinearMember
Definition: ap.h:976
ap::template_2d_array::gethighbound
int gethighbound(int iBoundNum) const
Definition: ap.h:933
ap::const_raw_vector::const_raw_vector
const_raw_vector(const T *Data, int Length, int Step)
Definition: ap.h:153
ap::abscomplex
const double abscomplex(const complex &z)
Definition: ap.cpp:97
ap::template_2d_array::wrongColumn
bool wrongColumn(int j) const
Definition: ap.h:971
ap::template_1d_array::getlowbound
int getlowbound(int iBoundNum=0) const
Definition: ap.h:772
ap::template_1d_array::gethighbound
int gethighbound(int iBoundNum=0) const
Definition: ap.h:778
ap::template_2d_array::m_iConstOffset
long m_iConstOffset
Definition: ap.h:976
ap::maxreal
double maxreal(double m1, double m2)
Definition: ap.cpp:172
ap::template_2d_array::m_iLow1
long m_iLow1
Definition: ap.h:975
ap::template_1d_array::operator=
const template_1d_array & operator=(const template_1d_array &rhs)
Definition: ap.h:698
ap::complex::operator=
complex & operator=(const double &v)
Definition: ap.h:72
ap::template_2d_array::setbounds
void setbounds(int iLow1, int iHigh1, int iLow2, int iHigh2)
Definition: ap.h:897
ap::raw_vector
Definition: ap.h:172
ap::complex::y
double y
Definition: ap.h:105
ap::randomreal
double randomreal()
Definition: ap.cpp:133
ap::operator-
const complex operator-(const complex &lhs)
Definition: ap.cpp:20
ap::template_1d_array::m_iVecSize
long m_iVecSize
Definition: ap.h:803
ap::template_2d_array::m_iHigh2
long m_iHigh2
Definition: ap.h:975
ap::complex::operator-=
complex & operator-=(const double &v)
Definition: ap.h:74
ap::minrealnumber
const double minrealnumber
Definition: svd_si.h:1001
ap::sqr
double sqr(double x)
Definition: ap.cpp:159
ap
Definition: ap.h:39
ap::vmul
void vmul(raw_vector< T > vdst, T2 alpha)
Definition: ap.h:608
ap::const_raw_vector::iStep
int iStep
Definition: ap.h:166
ap::operator*
const complex operator*(const complex &lhs, const complex &rhs)
Definition: ap.cpp:41
ap::operator!=
const bool operator!=(const complex &lhs, const complex &rhs)
Definition: ap.cpp:14
ap::iceil
int iceil(double x)
Definition: ap.cpp:153
ap::template_2d_array::getcolumn
raw_vector< T > getcolumn(int iColumn, int iRowStart, int iRowEnd)
Definition: ap.h:938
ap::complex::operator*=
complex & operator*=(const double &v)
Definition: ap.h:75
ap::real_2d_array
template_2d_array< double > real_2d_array
Definition: ap.h:983
ap::minreal
double minreal(double m1, double m2)
Definition: ap.cpp:177
i
int i
Definition: cfEzgcd.cc:125
ap::complex::operator+=
complex & operator+=(const double &v)
Definition: ap.h:73
ap::real_1d_array
template_1d_array< double > real_1d_array
Definition: ap.h:979
alpha
Variable alpha
Definition: facAbsBiFact.cc:52
ap::template_2d_array::wrongRow
bool wrongRow(int i) const
Definition: ap.h:970
T
static jList * T
Definition: janet.cc:31
ap::template_2d_array::m_iLow2
long m_iLow2
Definition: ap.h:975
ap::boolean_1d_array
template_1d_array< bool > boolean_1d_array
Definition: ap.h:981
ap::complex
Definition: ap.h:64
ap::template_1d_array::m_iHigh
long m_iHigh
Definition: ap.h:804
ap::template_2d_array::m_iVecSize
long m_iVecSize
Definition: ap.h:974
ap::raw_vector::raw_vector
raw_vector(T *Data, int Length, int Step)
Definition: ap.h:184
ap::boolean_2d_array
template_2d_array< bool > boolean_2d_array
Definition: ap.h:985
T2
void T2(ideal h)
Definition: cohomo.cc:2754
ap::const_raw_vector::GetLength
int GetLength() const
Definition: ap.h:159
ap::template_2d_array::getrow
raw_vector< T > getrow(int iRow, int iColumnStart, int iColumnEnd)
Definition: ap.h:946
ap::sign
int sign(double x)
Definition: ap.cpp:126
ap::vsub
void vsub(raw_vector< T > vdst, const_raw_vector< T > vsrc)
Definition: ap.h:538
ap::const_raw_vector::iLength
int iLength
Definition: ap.h:166
ap::ifloor
int ifloor(double x)
Definition: ap.cpp:150
ap::template_2d_array::~template_2d_array
~template_2d_array()
Definition: ap.h:822
ap::complex::x
double x
Definition: ap.h:103
ap::integer_2d_array
template_2d_array< int > integer_2d_array
Definition: ap.h:982
ap::minint
int minint(int m1, int m2)
Definition: ap.cpp:167
ap::ap_error::make_assertion
static void make_assertion(bool bClause)
Definition: ap.h:61
ap::pi
double pi()
Definition: ap.cpp:156
ap::complex::complex
complex()
Definition: ap.h:67
ap::maxint
int maxint(int m1, int m2)
Definition: ap.cpp:162
ap::maxrealnumber
const double maxrealnumber
Definition: svd_si.h:1000
ap::csqr
const complex csqr(const complex &z)
Definition: ap.cpp:120
ap::complex_2d_array
template_2d_array< complex > complex_2d_array
Definition: ap.h:984
ap::trunc
int trunc(double x)
Definition: ap.cpp:147
ap::template_1d_array::wrongIdx
bool wrongIdx(int i) const
Definition: ap.h:800
ap::template_2d_array::operator()
const T & operator()(int i1, int i2) const
Definition: ap.h:879
ap::template_1d_array::setbounds
void setbounds(int iLow, int iHigh)
Definition: ap.h:742
ap::const_raw_vector::GetStep
int GetStep() const
Definition: ap.h:162
ap::vmoveneg
void vmoveneg(raw_vector< T > vdst, const_raw_vector< T > vsrc)
Definition: ap.h:300
ap::vdotproduct
T vdotproduct(const_raw_vector< T > v1, const_raw_vector< T > v2)
Definition: ap.h:186
ap::template_1d_array::m_iLow
long m_iLow
Definition: ap.h:804
ap::ap_error
Definition: ap.h:51
ap::randominteger
int randominteger(int maxv)
Definition: ap.cpp:141
v
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:37
ap::template_2d_array::operator=
const template_2d_array & operator=(const template_2d_array &rhs)
Definition: ap.h:850
ap::template_2d_array::m_iHigh1
long m_iHigh1
Definition: ap.h:975
ap::complex::operator/=
complex & operator/=(const double &v)
Definition: ap.h:76
ap::template_1d_array::operator()
const T & operator()(int i) const
Definition: ap.h:724
ap::integer_1d_array
template_1d_array< int > integer_1d_array
Definition: ap.h:978
ap::template_2d_array::setcontent
void setcontent(int iLow1, int iHigh1, int iLow2, int iHigh2, const T *pContent)
Definition: ap.h:911
ap::complex_1d_array
template_1d_array< complex > complex_1d_array
Definition: ap.h:980
ap::const_raw_vector
Definition: ap.h:141
ap::template_1d_array::template_1d_array
template_1d_array()
Definition: ap.h:666
ap::const_raw_vector::GetData
const T * GetData() const
Definition: ap.h:156
ap::const_raw_vector::pData
T * pData
Definition: ap.h:163
ap::template_1d_array::m_Vec
T * m_Vec
Definition: ap.h:800
ap::template_1d_array::setcontent
void setcontent(int iLow, int iHigh, const T *pContent)
Definition: ap.h:753
ap::template_1d_array::~template_1d_array
~template_1d_array()
Definition: ap.h:672
ap::raw_vector::GetData
T * GetData()
Definition: ap.h:186