My Project
Functions | Variables
testsvdunit Namespace Reference

Functions

template<unsigned int Precision>
bool testsvd (bool silent)
 
template<unsigned int Precision>
void fillsparsea (ap::template_2d_array< amp::ampf< Precision > > &a, int m, int n, amp::ampf< Precision > sparcity)
 
template<unsigned int Precision>
void getsvderror (const ap::template_2d_array< amp::ampf< Precision > > &a, int m, int n, const ap::template_2d_array< amp::ampf< Precision > > &u, const ap::template_1d_array< amp::ampf< Precision > > &w, const ap::template_2d_array< amp::ampf< Precision > > &vt, amp::ampf< Precision > &materr, amp::ampf< Precision > &orterr, bool &wsorted)
 
template<unsigned int Precision>
void testsvdproblem (const ap::template_2d_array< amp::ampf< Precision > > &a, int m, int n, amp::ampf< Precision > &materr, amp::ampf< Precision > &orterr, amp::ampf< Precision > &othererr, bool &wsorted, bool &wfailed)
 
template<unsigned int Precision>
bool testsvdunit_test_silent ()
 
template<unsigned int Precision>
bool testsvdunit_test ()
 

Variables

STATIC_VAR int failcount
 
STATIC_VAR int succcount
 

Function Documentation

◆ fillsparsea()

template<unsigned int Precision>
void testsvdunit::fillsparsea ( ap::template_2d_array< amp::ampf< Precision > > &  a,
int  m,
int  n,
amp::ampf< Precision >  sparcity 
)

Definition at line 239 of file testsvdunit.h.

243 {
244 int i;
245 int j;
246
247
248 for(i=0; i<=m-1; i++)
249 {
250 for(j=0; j<=n-1; j++)
251 {
252 if( amp::ampf<Precision>::getRandom()>=sparcity )
253 {
254 a(i,j) = 2*amp::ampf<Precision>::getRandom()-1;
255 }
256 else
257 {
258 a(i,j) = 0;
259 }
260 }
261 }
262 }
m
int m
Definition: cfEzgcd.cc:128
i
int i
Definition: cfEzgcd.cc:132
amp::ampf
Definition: amp.h:82
amp::ampf::getRandom
static const ampf getRandom()
Definition: amp.h:593
j
int j
Definition: facHensel.cc:110

◆ getsvderror()

template<unsigned int Precision>
void testsvdunit::getsvderror ( const ap::template_2d_array< amp::ampf< Precision > > &  a,
int  m,
int  n,
const ap::template_2d_array< amp::ampf< Precision > > &  u,
const ap::template_1d_array< amp::ampf< Precision > > &  w,
const ap::template_2d_array< amp::ampf< Precision > > &  vt,
amp::ampf< Precision > &  materr,
amp::ampf< Precision > &  orterr,
bool &  wsorted 
)

Definition at line 266 of file testsvdunit.h.

275 {
276 int i;
277 int j;
278 int k;
279 int minmn;
280 amp::ampf<Precision> locerr;
281 amp::ampf<Precision> sm;
282
283
284 minmn = ap::minint(m, n);
285
286 //
287 // decomposition error
288 //
289 locerr = 0;
290 for(i=0; i<=m-1; i++)
291 {
292 for(j=0; j<=n-1; j++)
293 {
294 sm = 0;
295 for(k=0; k<=minmn-1; k++)
296 {
297 sm = sm+w(k)*u(i,k)*vt(k,j);
298 }
299 locerr = amp::maximum<Precision>(locerr, amp::abs<Precision>(a(i,j)-sm));
300 }
301 }
302 materr = amp::maximum<Precision>(materr, locerr);
303
304 //
305 // orthogonality error
306 //
307 locerr = 0;
308 for(i=0; i<=minmn-1; i++)
309 {
310 for(j=i; j<=minmn-1; j++)
311 {
312 sm = ap::vdotproduct(u.getcolumn(i, 0, m-1), u.getcolumn(j, 0, m-1));
313 if( i!=j )
314 {
315 locerr = amp::maximum<Precision>(locerr, amp::abs<Precision>(sm));
316 }
317 else
318 {
319 locerr = amp::maximum<Precision>(locerr, amp::abs<Precision>(sm-1));
320 }
321 sm = ap::vdotproduct(vt.getrow(i, 0, n-1), vt.getrow(j, 0, n-1));
322 if( i!=j )
323 {
324 locerr = amp::maximum<Precision>(locerr, amp::abs<Precision>(sm));
325 }
326 else
327 {
328 locerr = amp::maximum<Precision>(locerr, amp::abs<Precision>(sm-1));
329 }
330 }
331 }
332 orterr = amp::maximum<Precision>(orterr, locerr);
333
334 //
335 // values order error
336 //
337 for(i=1; i<=minmn-1; i++)
338 {
339 if( w(i)>w(i-1) )
340 {
341 wsorted = false;
342 }
343 }
344 }
k
int k
Definition: cfEzgcd.cc:99
ap::template_2d_array::getrow
raw_vector< T > getrow(int iRow, int iColumnStart, int iColumnEnd)
Definition: ap.h:939
ap::template_2d_array::getcolumn
raw_vector< T > getcolumn(int iColumn, int iRowStart, int iRowEnd)
Definition: ap.h:931
w
const CanonicalForm & w
Definition: facAbsFact.cc:51
ap::vdotproduct
T vdotproduct(const_raw_vector< T > v1, const_raw_vector< T > v2)
Definition: ap.h:181
ap::minint
int minint(int m1, int m2)
Definition: ap.cpp:167

◆ testsvd()

template<unsigned int Precision>
bool testsvdunit::testsvd ( bool  silent)

Definition at line 58 of file testsvdunit.h.

59 {
60 bool result;
61 ap::template_2d_array< amp::ampf<Precision> > a;
62 int m;
63 int n;
64 int maxmn;
65 int i;
66 int j;
67 int gpass;
68 int pass;
69 bool waserrors;
70 bool wsorted;
71 bool wfailed;
72 amp::ampf<Precision> materr;
73 amp::ampf<Precision> orterr;
74 amp::ampf<Precision> othererr;
75 amp::ampf<Precision> threshold;
76 amp::ampf<Precision> failthreshold;
77 amp::ampf<Precision> failr;
78
79
80 failcount = 0;
81 succcount = 0;
82 materr = 0;
83 orterr = 0;
84 othererr = 0;
85 wsorted = true;
86 wfailed = false;
87 waserrors = false;
88 maxmn = 30;
89 threshold = 5*100*amp::ampf<Precision>::getAlgoPascalEpsilon();
90 failthreshold = amp::ampf<Precision>("5.0E-3");
91 a.setbounds(0, maxmn-1, 0, maxmn-1);
92
93 //
94 // TODO: div by zero fail, convergence fail
95 //
96 for(gpass=1; gpass<=1; gpass++)
97 {
98
99 //
100 // zero matrix, several cases
101 //
102 for(i=0; i<=maxmn-1; i++)
103 {
104 for(j=0; j<=maxmn-1; j++)
105 {
106 a(i,j) = 0;
107 }
108 }
109 for(i=1; i<=ap::minint(5, maxmn); i++)
110 {
111 for(j=1; j<=ap::minint(5, maxmn); j++)
112 {
113 testsvdproblem<Precision>(a, i, j, materr, orterr, othererr, wsorted, wfailed);
114 }
115 }
116
117 //
118 // Long dense matrix
119 //
120 for(i=0; i<=maxmn-1; i++)
121 {
122 for(j=0; j<=ap::minint(5, maxmn)-1; j++)
123 {
124 a(i,j) = 2*amp::ampf<Precision>::getRandom()-1;
125 }
126 }
127 for(i=1; i<=maxmn; i++)
128 {
129 for(j=1; j<=ap::minint(5, maxmn); j++)
130 {
131 testsvdproblem<Precision>(a, i, j, materr, orterr, othererr, wsorted, wfailed);
132 }
133 }
134 for(i=0; i<=ap::minint(5, maxmn)-1; i++)
135 {
136 for(j=0; j<=maxmn-1; j++)
137 {
138 a(i,j) = 2*amp::ampf<Precision>::getRandom()-1;
139 }
140 }
141 for(i=1; i<=ap::minint(5, maxmn); i++)
142 {
143 for(j=1; j<=maxmn; j++)
144 {
145 testsvdproblem<Precision>(a, i, j, materr, orterr, othererr, wsorted, wfailed);
146 }
147 }
148
149 //
150 // Dense matrices
151 //
152 for(m=1; m<=ap::minint(10, maxmn); m++)
153 {
154 for(n=1; n<=ap::minint(10, maxmn); n++)
155 {
156 for(i=0; i<=m-1; i++)
157 {
158 for(j=0; j<=n-1; j++)
159 {
160 a(i,j) = 2*amp::ampf<Precision>::getRandom()-1;
161 }
162 }
163 testsvdproblem<Precision>(a, m, n, materr, orterr, othererr, wsorted, wfailed);
164 }
165 }
166
167 //
168 // Sparse matrices, very sparse matrices, incredible sparse matrices
169 //
170 for(m=1; m<=10; m++)
171 {
172 for(n=1; n<=10; n++)
173 {
174 for(pass=1; pass<=2; pass++)
175 {
176 fillsparsea<Precision>(a, m, n, amp::ampf<Precision>("0.8"));
177 testsvdproblem<Precision>(a, m, n, materr, orterr, othererr, wsorted, wfailed);
178 fillsparsea<Precision>(a, m, n, amp::ampf<Precision>("0.9"));
179 testsvdproblem<Precision>(a, m, n, materr, orterr, othererr, wsorted, wfailed);
180 fillsparsea<Precision>(a, m, n, amp::ampf<Precision>("0.95"));
181 testsvdproblem<Precision>(a, m, n, materr, orterr, othererr, wsorted, wfailed);
182 }
183 }
184 }
185 }
186
187 //
188 // report
189 //
190 failr = amp::ampf<Precision>(failcount)/(amp::ampf<Precision>(succcount+failcount));
191 waserrors = materr>threshold || orterr>threshold || othererr>threshold || !wsorted || failr>failthreshold;
192 if( !silent )
193 {
194 printf("TESTING SVD DECOMPOSITION\n");
195 printf("SVD decomposition error: %5.3le\n",
196 double(amp::ampf<Precision>(materr).toDouble()));
197 printf("SVD orthogonality error: %5.3le\n",
198 double(amp::ampf<Precision>(orterr).toDouble()));
199 printf("SVD with different parameters error: %5.3le\n",
200 double(amp::ampf<Precision>(othererr).toDouble()));
201 printf("Singular values order: ");
202 if( wsorted )
203 {
204 printf("OK\n");
205 }
206 else
207 {
208 printf("FAILED\n");
209 }
210 printf("Always converged: ");
211 if( !wfailed )
212 {
213 printf("YES\n");
214 }
215 else
216 {
217 printf("NO\n");
218 printf("Fail ratio: %5.3lf\n",
219 double(amp::ampf<Precision>(failr).toDouble()));
220 }
221 printf("Threshold: %5.3le\n",
222 double(amp::ampf<Precision>(threshold).toDouble()));
223 if( waserrors )
224 {
225 printf("TEST FAILED\n");
226 }
227 else
228 {
229 printf("TEST PASSED\n");
230 }
231 printf("\n\n");
232 }
233 result = !waserrors;
234 return result;
235 }
amp::ampf::getAlgoPascalEpsilon
static const ampf getAlgoPascalEpsilon()
Definition: amp.h:565
ap::template_2d_array
Definition: ap.h:807
ap::template_2d_array::setbounds
void setbounds(int iLow1, int iHigh1, int iLow2, int iHigh2)
Definition: ap.h:890
result
return result
Definition: facAbsBiFact.cc:75
testsvdunit::failcount
STATIC_VAR int failcount
Definition: testsvdunit.h:50
testsvdunit::succcount
STATIC_VAR int succcount
Definition: testsvdunit.h:51

◆ testsvdproblem()

template<unsigned int Precision>
void testsvdunit::testsvdproblem ( const ap::template_2d_array< amp::ampf< Precision > > &  a,
int  m,
int  n,
amp::ampf< Precision > &  materr,
amp::ampf< Precision > &  orterr,
amp::ampf< Precision > &  othererr,
bool &  wsorted,
bool &  wfailed 
)

Definition at line 348 of file testsvdunit.h.

356 {
357 ap::template_2d_array< amp::ampf<Precision> > u;
358 ap::template_2d_array< amp::ampf<Precision> > vt;
359 ap::template_2d_array< amp::ampf<Precision> > u2;
360 ap::template_2d_array< amp::ampf<Precision> > vt2;
361 ap::template_1d_array< amp::ampf<Precision> > w;
362 ap::template_1d_array< amp::ampf<Precision> > w2;
363 int i;
364 int j;
365 int k;
366 int ujob;
367 int vtjob;
368 int memjob;
369 int ucheck;
370 int vtcheck;
371 amp::ampf<Precision> v;
372 amp::ampf<Precision> mx;
373
374
375
376 //
377 // Main SVD test
378 //
379 if( !svd::rmatrixsvd<Precision>(a, m, n, 2, 2, 2, w, u, vt) )
380 {
381 failcount = failcount+1;
382 wfailed = true;
383 return;
384 }
385 getsvderror<Precision>(a, m, n, u, w, vt, materr, orterr, wsorted);
386
387 //
388 // Additional SVD tests
389 //
390 for(ujob=0; ujob<=2; ujob++)
391 {
392 for(vtjob=0; vtjob<=2; vtjob++)
393 {
394 for(memjob=0; memjob<=2; memjob++)
395 {
396 if( !svd::rmatrixsvd<Precision>(a, m, n, ujob, vtjob, memjob, w2, u2, vt2) )
397 {
398 failcount = failcount+1;
399 wfailed = true;
400 return;
401 }
402 ucheck = 0;
403 if( ujob==1 )
404 {
405 ucheck = ap::minint(m, n);
406 }
407 if( ujob==2 )
408 {
409 ucheck = m;
410 }
411 vtcheck = 0;
412 if( vtjob==1 )
413 {
414 vtcheck = ap::minint(m, n);
415 }
416 if( vtjob==2 )
417 {
418 vtcheck = n;
419 }
420 for(i=0; i<=m-1; i++)
421 {
422 for(j=0; j<=ucheck-1; j++)
423 {
424 othererr = amp::maximum<Precision>(othererr, amp::abs<Precision>(u(i,j)-u2(i,j)));
425 }
426 }
427 for(i=0; i<=vtcheck-1; i++)
428 {
429 for(j=0; j<=n-1; j++)
430 {
431 othererr = amp::maximum<Precision>(othererr, amp::abs<Precision>(vt(i,j)-vt2(i,j)));
432 }
433 }
434 for(i=0; i<=ap::minint(m, n)-1; i++)
435 {
436 othererr = amp::maximum<Precision>(othererr, amp::abs<Precision>(w(i)-w2(i)));
437 }
438 }
439 }
440 }
441
442 //
443 // update counter
444 //
445 succcount = succcount+1;
446 }
ap::template_1d_array
Definition: ap.h:657
v
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:39

◆ testsvdunit_test()

template<unsigned int Precision>
bool testsvdunit::testsvdunit_test ( )

Definition at line 467 of file testsvdunit.h.

468 {
469 bool result;
470
471
472 result = testsvd<Precision>(false);
473 return result;
474 }

◆ testsvdunit_test_silent()

template<unsigned int Precision>
bool testsvdunit::testsvdunit_test_silent ( )

Definition at line 453 of file testsvdunit.h.

454 {
455 bool result;
456
457
458 result = testsvd<Precision>(true);
459 return result;
460 }

Variable Documentation

◆ failcount

STATIC_VAR int testsvdunit::failcount

Definition at line 50 of file testsvdunit.h.

◆ succcount

STATIC_VAR int testsvdunit::succcount

Definition at line 51 of file testsvdunit.h.

Generated on Thu Apr 27 2023 08:16:02 for My Project by   doxygen 1.9.4 for Singular