Changeset aad4ca4 in git for Singular/svd

Timestamp:

Dec 6, 2017, 7:05:25 PM (6 years ago)

Author:

Hans Schoenemann <hannes@…>

Branches:

(u'fieker-DuVal', '117eb8c30fc9e991c4decca4832b1d19036c4c65')(u'spielwiese', '38dfc5131670d387a89455159ed1e071997eec94')

Children:

06b3e84d21e073249ada6d3b863bb3b8ce91d3be

Parents:

3803c6a39fad57f9a8e12942b2928c215504e027

Message:

format: Warn -> WarnS, trailing spaces

Location:

Singular/svd

Files:

• bdsvd.h (modified) (58 diffs)
• bidiagonal.h (modified) (49 diffs)
• blas.h (modified) (13 diffs)
• libs/amp.cpp (modified) (4 diffs)
• libs/amp.english.html (modified) (6 diffs)
• libs/amp.h (modified) (68 diffs)
• libs/amp.russian.html (modified) (7 diffs)
• libs/ap.cpp (modified) (1 diff)
• libs/ap.english.html (modified) (6 diffs)
• libs/ap.h (modified) (12 diffs)
• libs/ap.russian.html (modified) (1 diff)
• lq.h (modified) (13 diffs)
• makeheader (modified) (2 diffs)
• mkinstalldirs (modified) (1 diff)
• qr.h (modified) (14 diffs)
• reflections.h (modified) (8 diffs)
• rotations.h (modified) (10 diffs)
• tests/testsvdunit.h (modified) (12 diffs)

Singular/svd/bdsvd.h

@@ -346 +346 @@
             return result;
         }
-        
+
         //
         // init
@@ -363 +363 @@
         rightside = true;
         fwddir = true;
-        
+
         //
         // resize E from N-1 to N
@@ -379 +379 @@
         e(n) = 0;
         idir = 0;
-        
+
         //
         // Get machine constants
@@ -385 +385 @@
         eps = amp::ampf<Precision>::getAlgoPascalEpsilon();
         unfl = amp::ampf<Precision>::getAlgoPascalMinNumber();
-        
+
         //
         // If matrix lower bidiagonal, rotate to be upper bidiagonal
@@ -401 +401 @@
                 work1(i) = sn;
             }
-            
+
             //
             // Update singular vectors if desired
@@ -414 +414 @@
             }
         }
-        
+
         //
         // Compute singular values to relative accuracy TOL
@@ -426 +426 @@
             tol = -tol;
         }
-        
+
         //
         // Compute approximate maximum, minimum singular values
@@ -442 +442 @@
         if( tol>=0 )
         {
-            
+
             //
             // Relative accuracy desired
@@ -465 +465 @@
         else
         {
-            
+
             //
             // Absolute accuracy desired
@@ -471 +471 @@
             thresh = amp::maximum<Precision>(amp::abs<Precision>(tol)*smax, maxitr*n*n*unfl);
         }
-        
+
         //
         // Prepare for main iteration loop for the singular values
@@ -481 +481 @@
         oldll = -1;
         oldm = -1;
-        
+
         //
         // M points to last element of unconverged part of matrix
         //
         m = n;
-        
+
         //
         // Begin main iteration loop
@@ -492 +492 @@
         while( true )
         {
-            
+
             //
             // Check for convergence or exceeding iteration count
@@ -505 +505 @@
                 return result;
             }
-            
+
             //
             // Find diagonal block of matrix to work on
@@ -539 +539 @@
             else
             {
-                
+
                 //
                 // Matrix splits since E(LL) = 0
@@ -546 +546 @@
                 if( ll==m-1 )
                 {
-                    
+
                     //
                     // Convergence of bottom singular value, return to top of loop
@@ -555 +555 @@
             }
             ll = ll+1;
-            
+
             //
             // E(LL) through E(M-1) are nonzero, E(LL-1) is zero
@@ -561 +561 @@
             if( ll==m-1 )
             {
-                
+
                 //
                 // 2 by 2 block, handle separately
@@ -569 +569 @@
                 e(m-1) = 0;
                 d(m) = sigmn;
-                
+
                 //
                 // Compute singular vectors, if desired
@@ -606 +606 @@
                 continue;
             }
-            
+
             //
             // If working on new submatrix, choose shift direction
@@ -629 +629 @@
                 if( amp::abs<Precision>(d(ll))>=amp::abs<Precision>(d(m)) )
                 {
-                    
+
                     //
                     // Chase bulge from top (big end) to bottom (small end)
@@ -637 +637 @@
                 else
                 {
-                    
+
                     //
                     // Chase bulge from bottom (big end) to top (small end)
@@ -644 +644 @@
                 }
             }
-            
+
             //
             // Apply convergence tests
@@ -650 +650 @@
             if( idir==1 )
             {
-                
+
                 //
                 // Run convergence test in forward direction
@@ -662 +662 @@
                 if( tol>=0 )
                 {
-                    
+
                     //
                     // If relative accuracy desired,
@@ -690 +690 @@
             else
             {
-                
+
                 //
                 // Run convergence test in backward direction
@@ -702 +702 @@
                 if( tol>=0 )
                 {
-                    
+
                     //
                     // If relative accuracy desired,
@@ -730 +730 @@
             oldll = ll;
             oldm = m;
-            
+
             //
             // Compute shift.  First, test if shifting would ruin relative
@@ -737 +737 @@
             if( tol>=0 && n*tol*(sminl/smax)<=amp::maximum<Precision>(eps, amp::ampf<Precision>("0.01")*tol) )
             {
-                
+
                 //
                 // Use a zero shift to avoid loss of relative accuracy
@@ -745 +745 @@
             else
             {
-                
+
                 //
                 // Compute the shift from 2-by-2 block at end of matrix
@@ -759 +759 @@
                     svd2x2<Precision>(d(ll), e(ll), d(ll+1), shift, r);
                 }
-                
+
                 //
                 // Test if shift negligible, and if so set to zero
@@ -771 +771 @@
                 }
             }
-            
+
             //
             // Increment iteration count
             //
             iter = iter+m-ll;
-            
+
             //
             // If SHIFT = 0, do simplified QR iteration
@@ -784 +784 @@
                 if( idir==1 )
                 {
-                    
+
                     //
                     // Chase bulge from top to bottom
@@ -808 +808 @@
                     d(m) = h*oldcs;
                     e(m-1) = h*oldsn;
-                    
+
                     //
                     // Update singular vectors
@@ -824 +824 @@
                         rotations::applyrotationsfromtheleft<Precision>(fwddir, ll+cstart-1, m+cstart-1, cstart, cend, work2, work3, c, ctemp);
                     }
-                    
+
                     //
                     // Test convergence
@@ -835 +835 @@
                 else
                 {
-                    
+
                     //
                     // Chase bulge from bottom to top
@@ -859 +859 @@
                     d(ll) = h*oldcs;
                     e(ll) = h*oldsn;
-                    
+
                     //
                     // Update singular vectors
@@ -875 +875 @@
                         rotations::applyrotationsfromtheleft<Precision>(!fwddir, ll+cstart-1, m+cstart-1, cstart, cend, work0, work1, c, ctemp);
                     }
-                    
+
                     //
                     // Test convergence
@@ -887 +887 @@
             else
             {
-                
+
                 //
                 // Use nonzero shift
@@ -893 +893 @@
                 if( idir==1 )
                 {
-                    
+
                     //
                     // Chase bulge from top to bottom
@@ -926 +926 @@
                     }
                     e(m-1) = f;
-                    
+
                     //
                     // Update singular vectors
@@ -942 +942 @@
                         rotations::applyrotationsfromtheleft<Precision>(fwddir, ll+cstart-1, m+cstart-1, cstart, cend, work2, work3, c, ctemp);
                     }
-                    
+
                     //
                     // Test convergence
@@ -953 +953 @@
                 else
                 {
-                    
+
                     //
                     // Chase bulge from bottom to top
@@ -986 +986 @@
                     }
                     e(ll) = f;
-                    
+
                     //
                     // Test convergence
@@ -994 +994 @@
                         e(ll) = 0;
                     }
-                    
+
                     //
                     // Update singular vectors if desired
@@ -1012 +1012 @@
                 }
             }
-            
+
             //
             // QR iteration finished, go back and check convergence
@@ -1018 +1018 @@
             continue;
         }
-        
+
         //
         // All singular values converged, so make them positive
@@ -1027 +1027 @@
             {
                 d(i) = -d(i);
-                
+
                 //
                 // Change sign of singular vectors, if desired
@@ -1037 +1037 @@
             }
         }
-        
+
         //
         // Sort the singular values into decreasing order (insertion sort on
@@ -1044 +1044 @@
         for(i=1; i<=n-1; i++)
         {
-            
+
             //
             // Scan for smallest D(I)
@@ -1060 +1060 @@
             if( isub!=n+1-i )
             {
-                
+
                 //
                 // Swap singular values and vectors
@@ -1163 +1163 @@
                 if( au==0 )
                 {
-                    
+
                     //
                     // Avoid possible harmful underflow if exponent range
@@ -1228 +1228 @@
         ht = h;
         ha = amp::abs<Precision>(h);
-        
+
         //
         // PMAX points to the maximum absolute element of matrix
@@ -1239 +1239 @@
         if( swp )
         {
-            
+
             //
             // Now FA .ge. HA
@@ -1255 +1255 @@
         if( ga==0 )
         {
-            
+
             //
             // Diagonal matrix
@@ -1274 +1274 @@
                 if( fa/ga<amp::ampf<Precision>::getAlgoPascalEpsilon() )
                 {
-                    
+
                     //
                     // Case of very large GA
@@ -1298 +1298 @@
             if( gasmal )
             {
-                
+
                 //
                 // Normal case
@@ -1329 +1329 @@
                 if( mm==0 )
                 {
-                    
+
                     //
                     // Note that M is very tiny
@@ -1368 +1368 @@
             snr = srt;
         }
-        
+
         //
         // Correct signs of SSMAX and SSMIN

Singular/svd/bidiagonal.h

@@ -207 +207 @@
 
 
-        
+
         //
         // Prepare
@@ -231 +231 @@
         if( m>=n )
         {
-            
+
             //
             // Reduce to upper bidiagonal form
@@ -237 +237 @@
             for(i=0; i<=n-1; i++)
             {
-                
+
                 //
                 // Generate elementary reflector H(i) to annihilate A(i+1:m-1,i)
@@ -246 +246 @@
                 ap::vmove(a.getcolumn(i, i, m-1), t.getvector(1, m-i));
                 t(1) = 1;
-                
+
                 //
                 // Apply H(i) to A(i:m-1,i+1:n-1) from the left
@@ -253 +253 @@
                 if( i<n-1 )
                 {
-                    
+
                     //
                     // Generate elementary reflector G(i) to annihilate
@@ -263 +263 @@
                     ap::vmove(a.getrow(i, i+1, n-1), t.getvector(1, n-1-i));
                     t(1) = 1;
-                    
+
                     //
                     // Apply G(i) to A(i+1:m-1,i+1:n-1) from the right
@@ -277 +277 @@
         else
         {
-            
+
             //
             // Reduce to lower bidiagonal form
@@ -283 +283 @@
             for(i=0; i<=m-1; i++)
             {
-                
+
                 //
                 // Generate elementary reflector G(i) to annihilate A(i,i+1:n-1)
@@ -292 +292 @@
                 ap::vmove(a.getrow(i, i, n-1), t.getvector(1, n-i));
                 t(1) = 1;
-                
+
                 //
                 // Apply G(i) to A(i+1:m-1,i:n-1) from the right
@@ -299 +299 @@
                 if( i<m-1 )
                 {
-                    
+
                     //
                     // Generate elementary reflector H(i) to annihilate
@@ -309 +309 @@
                     ap::vmove(a.getcolumn(i, i+1, m-1), t.getvector(1, m-1-i));
                     t(1) = 1;
-                    
+
                     //
                     // Apply H(i) to A(i+1:m-1,i+1:n-1) from the left
@@ -363 +363 @@
             return;
         }
-        
+
         //
         // prepare Q
@@ -382 +382 @@
             }
         }
-        
+
         //
         // Calculate
@@ -444 +444 @@
         }
         ap::ap_error::make_assertion(fromtheright && zcolumns==m || !fromtheright && zrows==m);
-        
+
         //
         // init
@@ -455 +455 @@
         if( m>=n )
         {
-            
+
             //
             // setup
@@ -478 +478 @@
                 istep = -istep;
             }
-            
+
             //
             // Process
@@ -501 +501 @@
         else
         {
-            
+
             //
             // setup
@@ -524 +524 @@
                 istep = -istep;
             }
-            
+
             //
             // Process
@@ -590 +590 @@
             return;
         }
-        
+
         //
         // prepare PT
@@ -609 +609 @@
             }
         }
-        
+
         //
         // Calculate
@@ -671 +671 @@
         }
         ap::ap_error::make_assertion(fromtheright && zcolumns==n || !fromtheright && zrows==n);
-        
+
         //
         // init
@@ -684 +684 @@
         if( m>=n )
         {
-            
+
             //
             // setup
@@ -707 +707 @@
                 istep = -istep;
             }
-            
+
             //
             // Process
@@ -733 +733 @@
         else
         {
-            
+
             //
             // setup
@@ -756 +756 @@
                 istep = -istep;
             }
-            
+
             //
             // Process
@@ -875 +875 @@
         if( m>=n )
         {
-            
+
             //
             // Reduce to upper bidiagonal form
@@ -881 +881 @@
             for(i=1; i<=n; i++)
             {
-                
+
                 //
                 // Generate elementary reflector H(i) to annihilate A(i+1:m,i)
@@ -891 +891 @@
                 ap::vmove(a.getcolumn(i, i, m), t.getvector(1, mmip1));
                 t(1) = 1;
-                
+
                 //
                 // Apply H(i) to A(i:m,i+1:n) from the left
@@ -898 +898 @@
                 if( i<n )
                 {
-                    
+
                     //
                     // Generate elementary reflector G(i) to annihilate
@@ -910 +910 @@
                     ap::vmove(a.getrow(i, ip1, n), t.getvector(1, nmi));
                     t(1) = 1;
-                    
+
                     //
                     // Apply G(i) to A(i+1:m,i+1:n) from the right
@@ -924 +924 @@
         else
         {
-            
+
             //
             // Reduce to lower bidiagonal form
@@ -930 +930 @@
             for(i=1; i<=m; i++)
             {
-                
+
                 //
                 // Generate elementary reflector G(i) to annihilate A(i,i+1:n)
@@ -940 +940 @@
                 ap::vmove(a.getrow(i, i, n), t.getvector(1, nmip1));
                 t(1) = 1;
-                
+
                 //
                 // Apply G(i) to A(i+1:m,i:n) from the right
@@ -947 +947 @@
                 if( i<m )
                 {
-                    
+
                     //
                     // Generate elementary reflector H(i) to annihilate
@@ -959 +959 @@
                     ap::vmove(a.getcolumn(i, ip1, m), t.getvector(1, mmi));
                     t(1) = 1;
-                    
+
                     //
                     // Apply H(i) to A(i+1:m,i+1:n) from the left
@@ -999 +999 @@
             return;
         }
-        
+
         //
         // init
@@ -1006 +1006 @@
         v.setbounds(1, m);
         work.setbounds(1, qcolumns);
-        
+
         //
         // prepare Q
@@ -1079 +1079 @@
         }
         ap::ap_error::make_assertion(fromtheright && zcolumns==m || !fromtheright && zrows==m);
-        
+
         //
         // init
@@ -1090 +1090 @@
         if( m>=n )
         {
-            
+
             //
             // setup
@@ -1113 +1113 @@
                 istep = -istep;
             }
-            
+
             //
             // Process
@@ -1137 +1137 @@
         else
         {
-            
+
             //
             // setup
@@ -1160 +1160 @@
                 istep = -istep;
             }
-            
+
             //
             // Process
@@ -1214 +1214 @@
             return;
         }
-        
+
         //
         // init
@@ -1221 +1221 @@
         v.setbounds(1, n);
         work.setbounds(1, ptrows);
-        
+
         //
         // prepare PT
@@ -1294 +1294 @@
         }
         ap::ap_error::make_assertion(fromtheright && zcolumns==n || !fromtheright && zrows==n);
-        
+
         //
         // init
@@ -1307 +1307 @@
         if( m>=n )
         {
-            
+
             //
             // setup
@@ -1330 +1330 @@
                 istep = -istep;
             }
-            
+
             //
             // Process
@@ -1358 +1358 @@
         else
         {
-            
+
             //
             // setup
@@ -1381 +1381 @@
                 istep = -istep;
             }
-            
+
             //
             // Process

Singular/svd/blas.h

@@ -400 +400 @@
         if( !trans )
         {
-            
+
             //
             // y := alpha*A*x + beta*y;
@@ -410 +410 @@
             ap::ap_error::make_assertion(j2-j1==ix2-ix1);
             ap::ap_error::make_assertion(i2-i1==iy2-iy1);
-            
+
             //
             // beta*y
@@ -425 +425 @@
                 ap::vmul(y.getvector(iy1, iy2), beta);
             }
-            
+
             //
             // alpha*A*x
@@ -437 +437 @@
         else
         {
-            
+
             //
             // y := alpha*A'*x + beta*y;
@@ -447 +447 @@
             ap::ap_error::make_assertion(i2-i1==ix2-ix1);
             ap::ap_error::make_assertion(j2-j1==iy2-iy1);
-            
+
             //
             // beta*y
@@ -462 +462 @@
                 ap::vmul(y.getvector(iy1, iy2), beta);
             }
-            
+
             //
             // alpha*A'*x
@@ -538 +538 @@
 
 
-        
+
         //
         // Setup
@@ -569 +569 @@
         crows = arows;
         ccols = bcols;
-        
+
         //
         // Test WORK
@@ -578 +578 @@
         work(1) = 0;
         work(i) = 0;
-        
+
         //
         // Prepare C
@@ -599 +599 @@
             }
         }
-        
+
         //
         // A*B
@@ -616 +616 @@
             return;
         }
-        
+
         //
         // A*B'
@@ -647 +647 @@
             }
         }
-        
+
         //
         // A'*B
@@ -664 +664 @@
             return;
         }
-        
+
         //
         // A'*B'

Singular/svd/libs/amp.cpp

@@ -20 +20 @@
         lst = rec;
     }
-    
+
     amp::mpfr_record *p = lst;
     p->refCount = 1;
@@ -111 +111 @@
         free();
 }
-        
+
 void amp::mpfr_reference::initialize(int Precision)
 {
@@ -130 +130 @@
     ref = NULL;
 }
-        
+
 mpfr_srcptr amp::mpfr_reference::getReadPtr() const
 {
@@ -149 +149 @@
     amp::mpfr_record *newref = amp::mpfr_storage::newMpfr(ref->Precision);
     mpfr_set(newref->value, ref->value, GMP_RNDN);
-    
+
     free();
     ref = newref;

Singular/svd/libs/amp.english.html

@@ -128 +128 @@
         ~ampf();
 
-        ampf ();        
+        ampf ();
         ampf (<b>long double</b> v);
         ampf (<b>double</b> v);
@@ -142 +142 @@
         ampf (<b>const</b> std::string &amp;s);
         ampf (<b>const char</b> *s);
-        
+
         ampf(<b>const</b> ampf&amp; r);
         <b>template</b>&lt;<b>unsigned int</b> Precision2&gt;
@@ -163 +163 @@
         <b>template</b>&lt;<b>unsigned int</b> Precision2&gt;
         ampf&amp; <b>operator</b>= (<b>const</b> ampf&lt;Precision2&gt;&amp; r);
-        
+
         <b>template</b>&lt;<b>class</b> T&gt;
         ampf&amp; <b>operator</b>+=(<b>const</b> T&amp; v);
@@ -172 +172 @@
         <b>template</b>&lt;class T&gt;
         ampf&amp; <b>operator</b>/=(<b>const</b> T&amp; v);
-        
+
         mpfr_srcptr getReadPtr() <b>const</b>;
         mpfr_ptr getWritePtr();
-        
+
         <b>bool</b> isFiniteNumber() <b>const</b>;
         <b>bool</b> isPositiveNumber() <b>const</b>;
@@ -185 +185 @@
         std::string toHex() <b>const</b>;
         std::string toDec() <b>const</b>;
-        
+
         <b>static const</b> ampf getUlpOf(<b>const</b> ampf &amp;x);
         <b>static const</b> ampf getUlp();
@@ -499 +499 @@
 Calculates the arc tangent of y/x. It produces correct results even when the resulting angle is near pi/2 or -pi/2 (x near 0).
 </p>
-    
+
 <p>
 <span class=func>

Singular/svd/libs/amp.h

@@ -19 +19 @@
 namespace amp
 {
-    class exception {}; 
+    class exception {};
     class incorrectPrecision    : public exception {};
     class overflow              : public exception {};
@@ -28 +28 @@
     class internalError         : public exception {};
     class domainError           : public exception {};
-    
+
     typedef unsigned long unsigned32;
     typedef signed long   signed32;
-    
+
     struct mpfr_record
     {
@@ -39 +39 @@
         mpfr_record *next;
     };
-    
+
     typedef mpfr_record* mpfr_record_ptr;
-    
+
     //
     // storage for mpfr_t instances
@@ -47 +47 @@
     class mpfr_storage
     {
-    public:        
+    public:
         static mpfr_record* newMpfr(unsigned int Precision);
         static void deleteMpfr(mpfr_record* ref);
@@ -55 +55 @@
         static mpfr_record_ptr& getList(unsigned int Precision);
     };
-    
+
     //
     // mpfr_t reference
@@ -66 +66 @@
         mpfr_reference& operator= (const mpfr_reference &r);
         ~mpfr_reference();
-        
+
         void initialize(int Precision);
         void free();
-        
+
         mpfr_srcptr getReadPtr() const;
         mpfr_ptr getWritePtr();
@@ -75 +75 @@
         mpfr_record *ref;
     };
-        
+
     //
     // ampf template
@@ -92 +92 @@
                 mpfr_storage::deleteMpfr(rval);
         }
-        
+
         //
         // Initializing
@@ -98 +98 @@
         ampf ()                 { InitializeAsZero(); }
         ampf(mpfr_record *v)    { rval = v; }
-        
+
         ampf (long double v)    { InitializeAsDouble(v); }
         ampf (double v)         { InitializeAsDouble(v); }
@@ -110 +110 @@
         ampf (signed char v)    { InitializeAsSLong(v); }
         ampf (unsigned char v)  { InitializeAsULong(v); }
-        
+
         //
         // initializing from string
@@ -117 +117 @@
         ampf (const std::string &s) { InitializeAsString(s.c_str()); }
         ampf (const char *s)        { InitializeAsString(s); }
-        
+
         //
         // copy constructors
@@ -177 +177 @@
         }
 #endif
-        
+
         //
         // in-place operators
@@ -186 +186 @@
         template<class T> ampf& operator*=(const T& v){ *this = *this * v; return *this; };
         template<class T> ampf& operator/=(const T& v){ *this = *this / v; return *this; };
-        
+
         //
         // MPFR access
@@ -192 +192 @@
         mpfr_srcptr getReadPtr() const;
         mpfr_ptr getWritePtr();
-        
+
         //
         // properties and information
@@ -209 +209 @@
         std::string toDec() const;
         char * toString() const;
-        
-        
+
+
         //
         // static methods
-        //        
+        //
         static const ampf getUlpOf(const ampf &x);
         static const ampf getUlp();
@@ -231 +231 @@
         void InitializeAsDouble(long double v);
         void InitializeAsString(const char *s);
-        
+
         //mpfr_reference  ref;
         mpfr_record *rval;
@@ -264 +264 @@
         mpfr_set_si(getWritePtr(), sv, GMP_RNDN);
     }
-    
+
     template<unsigned int Precision>
     void ampf<Precision>::InitializeAsULong(unsigned long v)
@@ -272 +272 @@
         mpfr_set_ui(getWritePtr(), v, GMP_RNDN);
     }
-                           
+
     template<unsigned int Precision>
     void ampf<Precision>::InitializeAsDouble(long double v)
@@ -287 +287 @@
         rval = mpfr_storage::newMpfr(Precision);
         mpfr_strtofr(getWritePtr(), s, NULL, 0, GMP_RNDN);
-    } 
+    }
 
     template<unsigned int Precision>
     mpfr_srcptr ampf<Precision>::getReadPtr() const
     {
-        // TODO: ïîäóìàòü, íóæíî ëè ñäåëàòü, ÷òîáû è ïðè getRead, è ïðè 
+        // TODO: ïîäóìàòü, íóæíî ëè ñäåëàòü, ÷òîáû è ïðè getRead, è ïðè
         //       getWrite ñîçäàâàëàñü íîâàÿ instance mpfr_t.
         //       ýòî ìîæåò áûòü íóæíî äëÿ êîððåêòíîé îáðàáîòêè ñèòóàöèé âèäà
@@ -299 +299 @@
         return rval->value;
     }
-    
+
     template<unsigned int Precision>
     mpfr_ptr ampf<Precision>::getWritePtr()
@@ -311 +311 @@
         return rval->value;
     }
-        
+
     template<unsigned int Precision>
     bool ampf<Precision>::isFiniteNumber() const
@@ -325 +325 @@
         return mpfr_sgn(getReadPtr())>0;
     }
-    
+
     template<unsigned int Precision>
     bool ampf<Precision>::isZero() const
@@ -331 +331 @@
         return mpfr_zero_p(getReadPtr())!=0;
     }
-    
+
     template<unsigned int Precision>
     bool ampf<Precision>::isNegativeNumber() const
@@ -345 +345 @@
         return getUlpOf(*this);
     }
-    
+
     template<unsigned int Precision>
     double ampf<Precision>::toDouble() const
@@ -351 +351 @@
         return mpfr_get_d(getReadPtr(), GMP_RNDN);
     }
-    
+
     template<unsigned int Precision>
     std::string ampf<Precision>::toHex() const
@@ -368 +368 @@
             return r;
         }
-        
+
         //
         // general case
@@ -396 +396 @@
         return r;
     }
-    
+
     template<unsigned int Precision>
     std::string ampf<Precision>::toDec() const
     {
         // TODO: advanced output formatting (zero, integers)
-        
+
         //
         // some special cases
@@ -415 +415 @@
             return r;
         }
-        
+
         //
         // general case
@@ -459 +459 @@
             return toString_Block;
         }
-        
+
         //
         // general case
@@ -508 +508 @@
         return r;
     }
-        
+
     template<unsigned int Precision>
     const ampf<Precision> ampf<Precision>::getUlp()
@@ -517 +517 @@
         return r;
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> ampf<Precision>::getUlp256()
@@ -531 +531 @@
         return r;
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> ampf<Precision>::getUlp512()
@@ -544 +544 @@
             GMP_RNDN);
         return r;
-    } 
+    }
 
     template<unsigned int Precision>
@@ -554 +554 @@
         return r;
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> ampf<Precision>::getMinNumber()
@@ -568 +568 @@
         return getUlp256();
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> ampf<Precision>::getAlgoPascalMaxNumber()
@@ -579 +579 @@
         return r;
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> ampf<Precision>::getAlgoPascalMinNumber()
@@ -598 +598 @@
         return r;
     }
-    
+
     //
     // comparison operators
@@ -613 +613 @@
         return mpfr_cmp(op1.getReadPtr(), op2.getReadPtr())!=0;
     }
-    
+
     template<unsigned int Precision>
     const bool operator<(const ampf<Precision>& op1, const ampf<Precision>& op2)
@@ -619 +619 @@
         return mpfr_cmp(op1.getReadPtr(), op2.getReadPtr())<0;
     }
-    
+
     template<unsigned int Precision>
     const bool operator>(const ampf<Precision>& op1, const ampf<Precision>& op2)
@@ -625 +625 @@
         return mpfr_cmp(op1.getReadPtr(), op2.getReadPtr())>0;
     }
-    
+
     template<unsigned int Precision>
     const bool operator<=(const ampf<Precision>& op1, const ampf<Precision>& op2)
@@ -631 +631 @@
         return mpfr_cmp(op1.getReadPtr(), op2.getReadPtr())<=0;
     }
-    
+
     template<unsigned int Precision>
     const bool operator>=(const ampf<Precision>& op1, const ampf<Precision>& op2)
@@ -637 +637 @@
         return mpfr_cmp(op1.getReadPtr(), op2.getReadPtr())>=0;
     }
-    
+
     //
     // arithmetic operators
@@ -646 +646 @@
         return op1;
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> operator-(const ampf<Precision>& op1)
@@ -654 +654 @@
         return v;
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> operator+(const ampf<Precision>& op1, const ampf<Precision>& op2)
@@ -662 +662 @@
         return v;
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> operator-(const ampf<Precision>& op1, const ampf<Precision>& op2)
@@ -670 +670 @@
         return v;
     }
-        
-    
+
+
     template<unsigned int Precision>
     const ampf<Precision> operator*(const ampf<Precision>& op1, const ampf<Precision>& op2)
@@ -679 +679 @@
         return v;
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> operator/(const ampf<Precision>& op1, const ampf<Precision>& op2)
@@ -687 +687 @@
         return v;
     }
- 
+
     //
     // basic functions
@@ -843 +843 @@
         return r;
     }
- 
+
     //
     // different types of arguments
@@ -918 +918 @@
     __AMP_BINARY_OPF(long double)
     #undef __AMP_BINARY_OPF
-    
+
     //
     // transcendent functions
@@ -1003 +1003 @@
         return v;
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> log(const ampf<Precision> &x)
@@ -1077 +1077 @@
         campf():x(0),y(0){};
         campf(long double v)    { x=v; y=0; }
-        campf(double v)         { x=v; y=0; } 
+        campf(double v)         { x=v; y=0; }
         campf(float v)          { x=v; y=0; }
         campf(signed long v)    { x=v; y=0; }
@@ -1126 +1126 @@
         ampf<Precision> x, y;
     };
-    
+
     //
     // complex operations
@@ -1133 +1133 @@
     const bool operator==(const campf<Precision>& lhs, const campf<Precision>& rhs)
     { return lhs.x==rhs.x && lhs.y==rhs.y; }
-    
+
     template<unsigned int Precision>
     const bool operator!=(const campf<Precision>& lhs, const campf<Precision>& rhs)
     { return lhs.x!=rhs.x || lhs.y!=rhs.y; }
-    
+
     template<unsigned int Precision>
     const campf<Precision> operator+(const campf<Precision>& lhs)
     { return lhs; }
-    
+
     template<unsigned int Precision>
     campf<Precision>& operator+=(campf<Precision>& lhs, const campf<Precision>& rhs)
     { lhs.x += rhs.x; lhs.y += rhs.y; return lhs; }
-    
+
     template<unsigned int Precision>
     const campf<Precision> operator+(const campf<Precision>& lhs, const campf<Precision>& rhs)
     { campf<Precision> r = lhs; r += rhs; return r; }
-    
+
     template<unsigned int Precision>
     const campf<Precision> operator-(const campf<Precision>& lhs)
     { return campf<Precision>(-lhs.x, -lhs.y); }
-    
+
     template<unsigned int Precision>
     campf<Precision>& operator-=(campf<Precision>& lhs, const campf<Precision>& rhs)
     { lhs.x -= rhs.x; lhs.y -= rhs.y; return lhs; }
-    
+
     template<unsigned int Precision>
     const campf<Precision> operator-(const campf<Precision>& lhs, const campf<Precision>& rhs)
     { campf<Precision> r = lhs; r -= rhs; return r; }
-    
+
     template<unsigned int Precision>
     campf<Precision>& operator*=(campf<Precision>& lhs, const campf<Precision>& rhs)
@@ -1170 +1170 @@
         return lhs;
     }
-    
+
     template<unsigned int Precision>
     const campf<Precision> operator*(const campf<Precision>& lhs, const campf<Precision>& rhs)
     { campf<Precision> r = lhs; r *= rhs; return r; }
-    
+
     template<unsigned int Precision>
     const campf<Precision> operator/(const campf<Precision>& lhs, const campf<Precision>& rhs)
@@ -1197 +1197 @@
         return result;
     }
-    
+
     template<unsigned int Precision>
     campf<Precision>& operator/=(campf<Precision>& lhs, const campf<Precision>& rhs)
@@ -1204 +1204 @@
         return lhs;
     }
-    
+
     template<unsigned int Precision>
     const ampf<Precision> abscomplex(const campf<Precision> &z)
     {
         ampf<Precision> w, xabs, yabs, v;
-    
+
         xabs = abs(z.x);
         yabs = abs(z.y);
         w = xabs>yabs ? xabs : yabs;
-        v = xabs<yabs ? xabs : yabs; 
+        v = xabs<yabs ? xabs : yabs;
         if( v==0 )
             return w;
@@ -1222 +1222 @@
         }
     }
-    
+
     template<unsigned int Precision>
     const campf<Precision> conj(const campf<Precision> &z)
     {
-        return campf<Precision>(z.x, -z.y); 
-    }
-    
+        return campf<Precision>(z.x, -z.y);
+    }
+
     template<unsigned int Precision>
     const campf<Precision> csqr(const campf<Precision> &z)
     {
-        ampf<Precision> t = z.x*z.y; 
-        return campf<Precision>(sqr(z.x)-sqr(z.y), t+t); 
-    }
-    
+        ampf<Precision> t = z.x*z.y;
+        return campf<Precision>(sqr(z.x)-sqr(z.y), t+t);
+    }
+
     //
     // different types of arguments
@@ -1281 +1281 @@
         template<unsigned int Precision>                   bool operator==(const campf<Precision>& op1, const type& op2)             { return op1.x==op2 && op1.y==0; }   \
         template<unsigned int Precision>                   bool operator!=(const type& op1,             const campf<Precision>& op2) { return op1!=op2.x || op2.y!=0; }   \
-        template<unsigned int Precision>                   bool operator!=(const campf<Precision>& op1, const type& op2)             { return op1.x!=op2 || op1.y!=0; }   
+        template<unsigned int Precision>                   bool operator!=(const campf<Precision>& op1, const type& op2)             { return op1.x!=op2 || op1.y!=0; }
     __AMP_BINARY_OPF(float)
     __AMP_BINARY_OPF(double)
@@ -1287 +1287 @@
     __AMP_BINARY_OPF(ampf<Precision>)
     #undef __AMP_BINARY_OPF
-    
+
     //
     // Real linear algebra
@@ -1341 +1341 @@
         }
     }
-    
+
     template<unsigned int Precision>
     void vMoveNeg(ap::raw_vector< ampf<Precision> > vDst, ap::const_raw_vector< ampf<Precision> > vSrc)
@@ -1358 +1358 @@
         }
     }
-    
+
     template<unsigned int Precision, class T2>
     void vMove(ap::raw_vector< ampf<Precision> > vDst, ap::const_raw_vector< ampf<Precision> > vSrc, T2 alpha)
@@ -1376 +1376 @@
         }
     }
-    
+
     template<unsigned int Precision>
     void vAdd(ap::raw_vector< ampf<Precision> > vDst, ap::const_raw_vector< ampf<Precision> > vSrc)
@@ -1393 +1393 @@
         }
     }
-    
+
     template<unsigned int Precision, class T2>
     void vAdd(ap::raw_vector< ampf<Precision> > vDst, ap::const_raw_vector< ampf<Precision> > vSrc, T2 alpha)
@@ -1412 +1412 @@
         }
     }
-    
+
     template<unsigned int Precision>
     void vSub(ap::raw_vector< ampf<Precision> > vDst, ap::const_raw_vector< ampf<Precision> > vSrc)
@@ -1429 +1429 @@
         }
     }
-    
+
     template<unsigned int Precision, class T2>
     void vSub(ap::raw_vector< ampf<Precision> > vDst, ap::const_raw_vector< ampf<Precision> > vSrc, T2 alpha)
@@ -1435 +1435 @@
         vAdd(vDst, vSrc, -alpha);
     }
-    
+
     template<unsigned int Precision, class T2>
     void vMul(ap::raw_vector< ampf<Precision> > vDst, T2 alpha)
@@ -1450 +1450 @@
     }
 }
- 
+
 #endif

Singular/svd/libs/amp.russian.html

@@ -77 +77 @@
 
 <p>
-<b>Îïðåäåëåíèå òî÷íîñòè íà ýòàïå êîìïèëÿöèè</b>. Îñîáåííîñòüþ áèáëèîòåêè AMP ÿâëÿåòñÿ òî, ÷òî êëàññ, ïîçâîëÿþùèé îñóùåñòâëÿòü îïåðàöèè ñ âåùåñòâåííûìè ÷èñëàìè âûñîêîé òî÷íîñòè, ÿâëÿåòñÿ êëàññîì-øàáëîíîì, ïàðàìåòðîì êîòîðîãî ÿâëÿåòñÿ öåëîå ÷èñëî, çàäàþùåå òî÷íîñòü õðàíÿùåãîñÿ çíà÷åíèÿ. Òàêèì îáðàçîì, òî÷íîñòü, ñ êîòîðîé îñóùåñòâëÿþòñÿ âû÷èñëåíèÿ, ñòàíîâèòñÿ èçâåñòíà åùž íà ýòàïå êîìïèëÿöèè - ïðè ïàðàìåòðèçàöèè øàáëîíà, è íå ìîæåò áûòü èçìåíåíà â õîäå âûïîëíåíèÿ ïðîãðàììû. Ýòî ñóùåñòâåííîå îãðàíè÷åíèå ôóíêöèîíàëüíîñòè áèáëèîòåêè áûëî ââåäåíî ñîçíàòåëüíî. Ïðè÷èíîé ÿâëÿåòñÿ òî, ÷òî îñíîâíîå íàçíà÷åíèå áèáëèîòåêè AMP - ïðèìåíåíèå â ðàìêàõ ïðîåêòà ALGLIB.  ïðîåêòå ALGLIB áèáëèîòåêà AMP èñïîëüçóåòñÿ ïðîãðàììàìè, ÿâëÿþùèìèñÿ ðåçóëüòàòîì ðàáîòû àâòîìàòè÷åñêîãî òðàíñëÿòîðà. Ãåíåðèðóåìûé òðàíñëÿòîðîì êîä íå ñïîñîáåí ðàáîòàòü ñ âåùåñòâåííûìè ÷èñëàìè, ÷üÿ òî÷íîñòü ìåíÿåòñÿ äèíàìè÷åñêè, ïîñêîëüêó ïðè ýòîì âîçíèêàåò ðÿä ïðîáëåì (íàïðèìåð, êàêóþ òî÷íîñòü ñëåäóåò èñïîëüçîâàòü äëÿ ïðîìåæóòî÷íûõ ðåçóëüòàòîâ, åñëè îäíè ýëåìåíòû ìàòðèöû çàäàíû ñ òî÷íîñòüþ 128 áèò, à äðóãèå - ñ òî÷íîñòüþ 256 áèò). Ïî ýòîé ïðè÷èíå áûëî ðåøåíî ââåñòè âûøåóêàçàííîå îãðàíè÷åíèå. 
+<b>Îïðåäåëåíèå òî÷íîñòè íà ýòàïå êîìïèëÿöèè</b>. Îñîáåííîñòüþ áèáëèîòåêè AMP ÿâëÿåòñÿ òî, ÷òî êëàññ, ïîçâîëÿþùèé îñóùåñòâëÿòü îïåðàöèè ñ âåùåñòâåííûìè ÷èñëàìè âûñîêîé òî÷íîñòè, ÿâëÿåòñÿ êëàññîì-øàáëîíîì, ïàðàìåòðîì êîòîðîãî ÿâëÿåòñÿ öåëîå ÷èñëî, çàäàþùåå òî÷íîñòü õðàíÿùåãîñÿ çíà÷åíèÿ. Òàêèì îáðàçîì, òî÷íîñòü, ñ êîòîðîé îñóùåñòâëÿþòñÿ âû÷èñëåíèÿ, ñòàíîâèòñÿ èçâåñòíà åùž íà ýòàïå êîìïèëÿöèè - ïðè ïàðàìåòðèçàöèè øàáëîíà, è íå ìîæåò áûòü èçìåíåíà â õîäå âûïîëíåíèÿ ïðîãðàììû. Ýòî ñóùåñòâåííîå îãðàíè÷åíèå ôóíêöèîíàëüíîñòè áèáëèîòåêè áûëî ââåäåíî ñîçíàòåëüíî. Ïðè÷èíîé ÿâëÿåòñÿ òî, ÷òî îñíîâíîå íàçíà÷åíèå áèáëèîòåêè AMP - ïðèìåíåíèå â ðàìêàõ ïðîåêòà ALGLIB.  ïðîåêòå ALGLIB áèáëèîòåêà AMP èñïîëüçóåòñÿ ïðîãðàììàìè, ÿâëÿþùèìèñÿ ðåçóëüòàòîì ðàáîòû àâòîìàòè÷åñêîãî òðàíñëÿòîðà. Ãåíåðèðóåìûé òðàíñëÿòîðîì êîä íå ñïîñîáåí ðàáîòàòü ñ âåùåñòâåííûìè ÷èñëàìè, ÷üÿ òî÷íîñòü ìåíÿåòñÿ äèíàìè÷åñêè, ïîñêîëüêó ïðè ýòîì âîçíèêàåò ðÿä ïðîáëåì (íàïðèìåð, êàêóþ òî÷íîñòü ñëåäóåò èñïîëüçîâàòü äëÿ ïðîìåæóòî÷íûõ ðåçóëüòàòîâ, åñëè îäíè ýëåìåíòû ìàòðèöû çàäàíû ñ òî÷íîñòüþ 128 áèò, à äðóãèå - ñ òî÷íîñòüþ 256 áèò). Ïî ýòîé ïðè÷èíå áûëî ðåøåíî ââåñòè âûøåóêàçàííîå îãðàíè÷åíèå.
 </p>
 
@@ -128 +128 @@
         ~ampf();
 
-        ampf ();        
+        ampf ();
         ampf (<b>long double</b> v);
         ampf (<b>double</b> v);
@@ -142 +142 @@
         ampf (<b>const</b> std::string &amp;s);
         ampf (<b>const char</b> *s);
-        
+
         ampf(<b>const</b> ampf&amp; r);
         <b>template</b>&lt;<b>unsigned int</b> Precision2&gt;
@@ -163 +163 @@
         <b>template</b>&lt;<b>unsigned int</b> Precision2&gt;
         ampf&amp; <b>operator</b>= (<b>const</b> ampf&lt;Precision2&gt;&amp; r);
-        
+
         <b>template</b>&lt;<b>class</b> T&gt;
         ampf&amp; <b>operator</b>+=(<b>const</b> T&amp; v);
@@ -172 +172 @@
         <b>template</b>&lt;class T&gt;
         ampf&amp; <b>operator</b>/=(<b>const</b> T&amp; v);
-        
+
         mpfr_srcptr getReadPtr() <b>const</b>;
         mpfr_ptr getWritePtr();
-        
+
         <b>bool</b> isFiniteNumber() <b>const</b>;
         <b>bool</b> isPositiveNumber() <b>const</b>;
@@ -185 +185 @@
         std::string toHex() <b>const</b>;
         std::string toDec() <b>const</b>;
-        
+
         <b>static const</b> ampf getUlpOf(<b>const</b> ampf &amp;x);
         <b>static const</b> ampf getUlp();
@@ -496 +496 @@
 Ýòà ôóíêöèÿ âîçâðàùàåò àðêòàíãåíñ ÷èñëà, ðàâíîãî îòíîøåíèþ àðãóìåíòîâ <code>y/x</code>. Ôóíêöèÿ âîçâðàùàåò êîððåêòíûé ðåçóëüòàò, äàæå åñëè <code>x</code> ðàâíî <code>0</code>.
 </p>
-    
+
 <p>
 <span class=func>

Singular/svd/libs/ap.cpp

@@ -104 +104 @@
     yabs = fabs(z.y);
     w = xabs>yabs ? xabs : yabs;
-    v = xabs<yabs ? xabs : yabs; 
+    v = xabs<yabs ? xabs : yabs;
     if( v==0 )
         return w;

Singular/svd/libs/ap.english.html

@@ -142 +142 @@
 
 <p align=justify>
-Classes of the standard library allow operations with matrixes and vectors (one-dimensional and two-dimensional arrays) of variable size and with variable numeration of elements, that is, the array numeration can start at any number, end at any number and change dynamically. Because the array classes are templates, the arrays of the same dimension have the same set of member functions. And as the member functions of arrays with different dimensions differ only by the number of arguments, there is little difference between two-dimensional and one-dimenstional arrays. 
+Classes of the standard library allow operations with matrixes and vectors (one-dimensional and two-dimensional arrays) of variable size and with variable numeration of elements, that is, the array numeration can start at any number, end at any number and change dynamically. Because the array classes are templates, the arrays of the same dimension have the same set of member functions. And as the member functions of arrays with different dimensions differ only by the number of arguments, there is little difference between two-dimensional and one-dimenstional arrays.
 </p>
 
@@ -154 +154 @@
 
 <p align=justify>
-To address the array elements, an overloaded <code:>operator()</code:> is used. That is, the code addressing the element of array <code:>a</code:> with indexes <code:>a(i,j,k)</code:> will look like <code:>a(i,j,k)</code:>. Below is given an example of factorial array calculation, illustrating the work with arrays. 
+To address the array elements, an overloaded <code:>operator()</code:> is used. That is, the code addressing the element of array <code:>a</code:> with indexes <code:>a(i,j,k)</code:> will look like <code:>a(i,j,k)</code:>. Below is given an example of factorial array calculation, illustrating the work with arrays.
 </p>
 
@@ -201 +201 @@
 
 <p align=justify>
-<span class=func>T&amp; operator()(<b>int</b> i)</span><br> Addressing the array element number i 
+<span class=func>T&amp; operator()(<b>int</b> i)</span><br> Addressing the array element number i
 </p>
 
@@ -213 +213 @@
 
 <p align=justify>
-<span class=func>T* getcontent()</span><br> Gets pointer to the array. The data pointed by the returned pointer can be changed, and the array content will be changed as well. 
+<span class=func>T* getcontent()</span><br> Gets pointer to the array. The data pointed by the returned pointer can be changed, and the array content will be changed as well.
 </p>
 
@@ -225 +225 @@
 
 <p align=justify>
-<span class=func>const_raw_vector&lt;T&gt; getvector(<b>int</b> iStart, <b>int</b> iEnd) <b>const</b></span><br> The method is used by the basic subroutines of linear algebra to get access to the internal memory of the array. The method returns an object, holding the pointer to a vector part (starting from the element with iStart index value and finishing with iEnd index value). If iEnd&lt;iStart, then an empty vector is considered to be set. The returned object is for read only. 
+<span class=func>const_raw_vector&lt;T&gt; getvector(<b>int</b> iStart, <b>int</b> iEnd) <b>const</b></span><br> The method is used by the basic subroutines of linear algebra to get access to the internal memory of the array. The method returns an object, holding the pointer to a vector part (starting from the element with iStart index value and finishing with iEnd index value). If iEnd&lt;iStart, then an empty vector is considered to be set. The returned object is for read only.
 </p>
 
@@ -273 +273 @@
 
 <p align=justify>
-<span class=func>T* getcontent()</span><br> Gets pointer to the array. The data pointed by the returned pointer can be changed, and the array content will be changed as well. 
+<span class=func>T* getcontent()</span><br> Gets pointer to the array. The data pointed by the returned pointer can be changed, and the array content will be changed as well.
 </p>
 

Singular/svd/libs/ap.h

@@ -28 +28 @@
 
 /********************************************************************
-This symbol is used for debugging. Do not define it and do not remove 
-comments. 
+This symbol is used for debugging. Do not define it and do not remove
+comments.
 ********************************************************************/
 //#define UNSAFE_MEM_COPY
@@ -124 +124 @@
 
 /********************************************************************
-Template defining vector in memory. It is used by the basic 
+Template defining vector in memory. It is used by the basic
 subroutines of linear algebra.
 
-Vector consists of Length elements of type T, starting from an element, 
-which Data is pointed to. Interval between adjacent elements equals 
+Vector consists of Length elements of type T, starting from an element,
+which Data is pointed to. Interval between adjacent elements equals
 the value of Step.
 
@@ -158 +158 @@
 It is used by the basic subroutines of linear algebra.
 
-Vector consists of Length elements of type T, starting from an element, 
-which Data is pointed to. Interval between adjacent elements equals 
+Vector consists of Length elements of type T, starting from an element,
+which Data is pointed to. Interval between adjacent elements equals
 the value of Step.
 
@@ -688 +688 @@
     };
 
-    
+
     const template_1d_array& operator=(const template_1d_array &rhs)
     {
@@ -714 +714 @@
     };
 
-    
+
     const T& operator()(int i) const
     {
@@ -723 +723 @@
     };
 
-    
+
     T& operator()(int i)
     {
@@ -732 +732 @@
     };
 
-    
+
     void setbounds( int iLow, int iHigh )
     {
@@ -743 +743 @@
     };
 
-    
+
     void setcontent( int iLow, int iHigh, const T *pContent )
     {
@@ -751 +751 @@
     };
 
-    
+
     T* getcontent()
     {
@@ -762 +762 @@
     };
 
-    
+
     int getlowbound(int iBoundNum = 0) const
     {
@@ -768 +768 @@
     };
 
-    
+
     int gethighbound(int iBoundNum = 0) const
     {
@@ -782 +782 @@
     };
 
-    
+
     const_raw_vector<T> getvector(int iStart, int iEnd) const
     {

Singular/svd/libs/ap.russian.html

@@ -184 +184 @@
 
 <p align=justify>
-Äëÿ îáðàùåíèÿ ê ýëåìåíòàì ìàññèâà èñïîëüçóåòñÿ ïåðåãðóæåííûé <code:>operator()</code:>. Ò.å. êîä, îáðàùàþùèéñÿ ê ýëåìåíòó ìàññèâà <code:>a</code:> ñ èíäåêñàìè <code:>i, j, k</code:> áóäåò âûãëÿäåòü êàê <code:>a(i,j,k)</code:>. Íèæå ïðèâåäåí ïðèìåð âû÷èñëåíèÿ ìàññèâà ôàêòîðèàëîâ, èëëþñòðèðóþùèé ðàáîòó ñ ìàññèâàìè. 
+Äëÿ îáðàùåíèÿ ê ýëåìåíòàì ìàññèâà èñïîëüçóåòñÿ ïåðåãðóæåííûé <code:>operator()</code:>. Ò.å. êîä, îáðàùàþùèéñÿ ê ýëåìåíòó ìàññèâà <code:>a</code:> ñ èíäåêñàìè <code:>i, j, k</code:> áóäåò âûãëÿäåòü êàê <code:>a(i,j,k)</code:>. Íèæå ïðèâåäåí ïðèìåð âû÷èñëåíèÿ ìàññèâà ôàêòîðèàëîâ, èëëþñòðèðóþùèé ðàáîòó ñ ìàññèâàìè.
 </p>
 

Singular/svd/lq.h

@@ -133 +133 @@
         for(i=0; i<=k-1; i++)
         {
-            
+
             //
             // Generate elementary reflector H(i) to annihilate A(i,i+1:n-1)
@@ -144 +144 @@
             if( i<n )
             {
-                
+
                 //
                 // Apply H(i) to A(i+1:m,i:n) from the right
@@ -195 +195 @@
             return;
         }
-        
+
         //
         // init
@@ -218 +218 @@
             }
         }
-        
+
         //
         // unpack Q
@@ -224 +224 @@
         for(i=k-1; i>=0; i--)
         {
-            
+
             //
             // Apply H(i)
@@ -306 +306 @@
         t.setbounds(1, n);
         tau.setbounds(1, minmn);
-        
+
         //
         // Test the input arguments
@@ -313 +313 @@
         for(i=1; i<=k; i++)
         {
-            
+
             //
             // Generate elementary reflector H(i) to annihilate A(i,i+1:n)
@@ -325 +325 @@
             if( i<n )
             {
-                
+
                 //
                 // Apply H(i) to A(i+1:m,i:n) from the right
@@ -361 +361 @@
             return;
         }
-        
+
         //
         // init
@@ -384 +384 @@
             }
         }
-        
+
         //
         // unpack Q
@@ -390 +390 @@
         for(i=k; i>=1; i--)
         {
-            
+
             //
             // Apply H(i)
@@ -423 +423 @@
         q.setbounds(1, n, 1, n);
         l.setbounds(1, m, 1, n);
-        
+
         //
         // LQDecomposition
         //
         lqdecomposition<Precision>(a, m, n, tau);
-        
+
         //
         // L
@@ -446 +446 @@
             }
         }
-        
+
         //
         // Q

Singular/svd/makeheader

@@ -18 +18 @@
 
 To include a header file into the generated header file
-(= <outfile>) use the two lines 
+(= <outfile>) use the two lines
 
 /*MAKEHEADER*/
@@ -73 +73 @@
     exit 1
 fi
-    
+
 # try to create $outfile if it does not exist
 if test -w "$outfile" || cp /dev/null "$outfile" > /dev/null 2> /dev/null; then

Singular/svd/mkinstalldirs

@@ -8 +8 @@
 errstatus=0
 
-for file in ${1+"$@"} ; do 
+for file in ${1+"$@"} ; do
    set fnord `echo ":$file" | sed -ne 's/^:\//#/;s/^://;s/\// /g;s/^#/\//;p'`
    shift

Singular/svd/qr.h

@@ -142 +142 @@
         t.setbounds(1, m);
         tau.setbounds(0, minmn-1);
-        
+
         //
         // Test the input arguments
@@ -149 +149 @@
         for(i=0; i<=k-1; i++)
         {
-            
+
             //
             // Generate elementary reflector H(i) to annihilate A(i+1:m,i)
@@ -160 +160 @@
             if( i<n )
             {
-                
+
                 //
                 // Apply H(i) to A(i:m-1,i+1:n-1) from the left
@@ -211 +211 @@
             return;
         }
-        
+
         //
         // init
@@ -234 +234 @@
             }
         }
-        
+
         //
         // unpack Q
@@ -240 +240 @@
         for(i=k-1; i>=0; i--)
         {
-            
+
             //
             // Apply H(i)
@@ -319 +319 @@
         t.setbounds(1, m);
         tau.setbounds(1, minmn);
-        
+
         //
         // Test the input arguments
@@ -326 +326 @@
         for(i=1; i<=k; i++)
         {
-            
+
             //
             // Generate elementary reflector H(i) to annihilate A(i+1:m,i)
@@ -338 +338 @@
             if( i<n )
             {
-                
+
                 //
                 // Apply H(i) to A(i:m,i+1:n) from the left
@@ -373 +373 @@
             return;
         }
-        
+
         //
         // init
@@ -396 +396 @@
             }
         }
-        
+
         //
         // unpack Q
@@ -402 +402 @@
         for(i=k; i>=1; i--)
         {
-            
+
             //
             // Apply H(i)
@@ -440 +440 @@
         q.setbounds(1, m, 1, m);
         r.setbounds(1, m, 1, n);
-        
+
         //
         // QRDecomposition
         //
         qrdecomposition<Precision>(a, m, n, tau);
-        
+
         //
         // R
@@ -461 +461 @@
             ap::vmove(r.getrow(i, i, n), a.getrow(i, i, n));
         }
-        
+
         //
         // Q

Singular/svd/reflections.h

@@ -122 +122 @@
 
 
-        
+
         //
         // Executable Statements ..
@@ -131 +131 @@
             return;
         }
-        
+
         //
         // XNORM = DNRM2( N-1, X, INCX )
@@ -152 +152 @@
         if( xnorm==0 )
         {
-            
+
             //
             // H  =  I
@@ -159 +159 @@
             return;
         }
-        
+
         //
         // general case
@@ -223 +223 @@
             return;
         }
-        
+
         //
         // w := C' * v
@@ -237 +237 @@
             ap::vadd(work.getvector(n1, n2), c.getrow(i, n1, n2), t);
         }
-        
+
         //
         // C := C - tau * v * w'
@@ -296 +296 @@
             return;
         }
-        
+
         //
         // w := C * v
@@ -306 +306 @@
             work(i) = t;
         }
-        
+
         //
         // C := C - w * v'

Singular/svd/rotations.h

@@ -119 +119 @@
             return;
         }
-        
+
         //
         // Form  P * A
@@ -127 +127 @@
             if( n1!=n2 )
             {
-                
+
                 //
                 // Common case: N1<>N2
@@ -148 +148 @@
             else
             {
-                
+
                 //
                 // Special case: N1=N2
@@ -169 +169 @@
             if( n1!=n2 )
             {
-                
+
                 //
                 // Common case: N1<>N2
@@ -190 +190 @@
             else
             {
-                
+
                 //
                 // Special case: N1=N2
@@ -253 +253 @@
 
 
-        
+
         //
         // Form A * P'
@@ -261 +261 @@
             if( m1!=m2 )
             {
-                
+
                 //
                 // Common case: M1<>M2
@@ -282 +282 @@
             else
             {
-                
+
                 //
                 // Special case: M1=M2
@@ -303 +303 @@
             if( m1!=m2 )
             {
-                
+
                 //
                 // Common case: M1<>M2
@@ -324 +324 @@
             else
             {
-                
+
                 //
                 // Special case: M1=M2

Singular/svd/tests/testsvdunit.h

@@ -90 +90 @@
         failthreshold = amp::ampf<Precision>("5.0E-3");
         a.setbounds(0, maxmn-1, 0, maxmn-1);
-        
+
         //
         // TODO: div by zero fail, convergence fail
@@ -96 +96 @@
         for(gpass=1; gpass<=1; gpass++)
         {
-            
+
             //
             // zero matrix, several cases
@@ -114 +114 @@
                 }
             }
-            
+
             //
             // Long dense matrix
@@ -146 +146 @@
                 }
             }
-            
+
             //
             // Dense matrices
@@ -164 +164 @@
                 }
             }
-            
+
             //
             // Sparse matrices, very sparse matrices, incredible sparse matrices
@@ -184 +184 @@
             }
         }
-        
+
         //
         // report
@@ -283 +283 @@
 
         minmn = ap::minint(m, n);
-        
+
         //
         // decomposition error
@@ -301 +301 @@
         }
         materr = amp::maximum<Precision>(materr, locerr);
-        
+
         //
         // orthogonality error
@@ -331 +331 @@
         }
         orterr = amp::maximum<Precision>(orterr, locerr);
-        
+
         //
         // values order error
@@ -373 +373 @@
 
 
-        
+
         //
         // Main SVD test
@@ -384 +384 @@
         }
         getsvderror<Precision>(a, m, n, u, w, vt, materr, orterr, wsorted);
-        
+
         //
         // Additional SVD tests
@@ -439 +439 @@
             }
         }
-        
+
         //
         // update counter