1 /*
2 
3     Copyright (C) 2014, The University of Texas at Austin
4 
5     This file is part of libflame and is available under the 3-Clause
6     BSD license, which can be found in the LICENSE file at the top-level
7     directory, or at http://opensource.org/licenses/BSD-3-Clause
8 
9 */
10 
11 #include "FLAME.h"
12 
FLA_Obj_datatype(FLA_Obj obj)13 FLA_Datatype FLA_Obj_datatype( FLA_Obj obj )
14 {
15   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
16     FLA_Obj_datatype_check( obj );
17 
18   return obj.base->datatype;
19 }
20 
21 
22 
FLA_Obj_datatype_proj_to_real(FLA_Obj A)23 FLA_Datatype FLA_Obj_datatype_proj_to_real( FLA_Obj A )
24 {
25 	FLA_Datatype datatype;
26 
27 	if ( FLA_Obj_is_single_precision( A ) )
28 		datatype = FLA_FLOAT;
29 	else
30 		datatype = FLA_DOUBLE;
31 
32 	return datatype;
33 }
34 
35 
36 
FLA_Obj_datatype_proj_to_complex(FLA_Obj A)37 FLA_Datatype FLA_Obj_datatype_proj_to_complex( FLA_Obj A )
38 {
39 	FLA_Datatype datatype;
40 
41 	if ( FLA_Obj_is_single_precision( A ) )
42 		datatype = FLA_COMPLEX;
43 	else
44 		datatype = FLA_DOUBLE_COMPLEX;
45 
46 	return datatype;
47 }
48 
49 
50 
FLA_Obj_elemtype(FLA_Obj obj)51 FLA_Elemtype FLA_Obj_elemtype( FLA_Obj obj )
52 {
53   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
54     FLA_Obj_elemtype_check( obj );
55 
56   return obj.base->elemtype;
57 }
58 
59 
60 
FLA_Obj_datatype_size(FLA_Datatype datatype)61 dim_t FLA_Obj_datatype_size( FLA_Datatype datatype )
62 {
63   dim_t datatype_size = 0;
64 
65   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
66     FLA_Obj_datatype_size_check( datatype );
67 
68   switch( datatype )
69   {
70     case FLA_INT:
71       datatype_size = sizeof( int );
72       break;
73     case FLA_FLOAT:
74       datatype_size = sizeof( float );
75       break;
76     case FLA_DOUBLE:
77       datatype_size = sizeof( double );
78       break;
79     case FLA_COMPLEX:
80       datatype_size = sizeof( scomplex );
81       break;
82     case FLA_DOUBLE_COMPLEX:
83       datatype_size = sizeof( dcomplex );
84       break;
85     case FLA_CONSTANT:
86       datatype_size = FLA_CONSTANT_SIZE;
87       break;
88   }
89 
90   return datatype_size;
91 }
92 
93 
94 
FLA_Obj_elem_size(FLA_Obj obj)95 dim_t FLA_Obj_elem_size( FLA_Obj obj )
96 {
97   dim_t elem_size = 0;
98 
99   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
100     FLA_Obj_elem_size_check( obj );
101 
102   if ( FLA_Obj_elemtype( obj ) == FLA_MATRIX )
103   {
104     elem_size = sizeof( FLA_Obj );
105   }
106   else // if ( FLA_Obj_elemtype( obj ) == FLA_SCALAR )
107   {
108     elem_size = FLA_Obj_datatype_size( FLA_Obj_datatype( obj ) );
109   }
110 
111   return elem_size;
112 }
113 
114 
115 
FLA_Obj_length(FLA_Obj obj)116 dim_t FLA_Obj_length( FLA_Obj obj )
117 {
118   return obj.m;
119 }
120 
121 
122 
FLA_Obj_width(FLA_Obj obj)123 dim_t FLA_Obj_width( FLA_Obj obj )
124 {
125   return obj.n;
126 }
127 
128 
129 
FLA_Obj_structure(FLA_Obj obj)130 FLA_Uplo FLA_Obj_structure( FLA_Obj obj )
131 {
132   return obj.base->uplo;
133 }
134 
135 
136 
FLA_Obj_vector_dim(FLA_Obj obj)137 dim_t FLA_Obj_vector_dim( FLA_Obj obj )
138 {
139   return ( obj.m == 1 ? obj.n
140                       : obj.m );
141 }
142 
143 
144 
FLA_Obj_vector_inc(FLA_Obj obj)145 dim_t FLA_Obj_vector_inc( FLA_Obj obj )
146 {
147   return ( obj.m == 1 ? (obj.base)->cs
148                       : (obj.base)->rs );
149 }
150 
151 
152 
FLA_Obj_min_dim(FLA_Obj obj)153 dim_t FLA_Obj_min_dim( FLA_Obj obj )
154 {
155   return min( obj.m, obj.n );
156 }
157 
158 
159 
FLA_Obj_max_dim(FLA_Obj obj)160 dim_t FLA_Obj_max_dim( FLA_Obj obj )
161 {
162   return max( obj.m, obj.n );
163 }
164 
165 
166 
FLA_Obj_row_stride(FLA_Obj obj)167 dim_t FLA_Obj_row_stride( FLA_Obj obj )
168 {
169   return (obj.base)->rs;
170 }
171 
172 
173 
FLA_Obj_col_stride(FLA_Obj obj)174 dim_t FLA_Obj_col_stride( FLA_Obj obj )
175 {
176   return (obj.base)->cs;
177 }
178 
179 
FLA_Obj_row_offset(FLA_Obj obj)180 dim_t FLA_Obj_row_offset( FLA_Obj obj )
181 {
182 	return obj.offm;
183 }
184 
185 
FLA_Obj_col_offset(FLA_Obj obj)186 dim_t FLA_Obj_col_offset( FLA_Obj obj )
187 {
188 	return obj.offn;
189 }
190 
191 
FLA_Obj_base_length(FLA_Obj obj)192 dim_t FLA_Obj_base_length( FLA_Obj obj )
193 {
194 	return (obj.base)->m;
195 }
196 
197 
FLA_Obj_base_width(FLA_Obj obj)198 dim_t FLA_Obj_base_width( FLA_Obj obj )
199 {
200 	return (obj.base)->n;
201 }
202 
203 
FLA_Obj_num_elem_alloc(FLA_Obj obj)204 dim_t FLA_Obj_num_elem_alloc( FLA_Obj obj )
205 {
206   return (obj.base)->n_elem_alloc;
207 }
208 
209 
FLA_Obj_base_buffer(FLA_Obj obj)210 void* FLA_Obj_base_buffer( FLA_Obj obj )
211 {
212   return (obj.base)->buffer;
213 }
214 
FLA_Obj_buffer_at_view(FLA_Obj obj)215 void* FLA_Obj_buffer_at_view( FLA_Obj obj )
216 {
217   char*  buffer;
218   size_t elem_size, offm, offn, rs, cs;
219   size_t byte_offset;
220 
221   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
222     FLA_Obj_buffer_at_view_check( obj );
223 
224   elem_size   = ( size_t ) FLA_Obj_elem_size( obj );
225   rs          = ( size_t ) FLA_Obj_row_stride( obj );
226   cs          = ( size_t ) FLA_Obj_col_stride( obj );
227   offm        = ( size_t ) obj.offm;
228   offn        = ( size_t ) obj.offn;
229 
230   byte_offset = elem_size * ( offn * cs + offm * rs );
231 
232   buffer      = ( char * ) (obj.base)->buffer;
233 
234   return ( void* ) ( buffer + byte_offset );
235 }
236 
237 
238 
FLA_Obj_buffer_is_null(FLA_Obj obj)239 FLA_Bool FLA_Obj_buffer_is_null( FLA_Obj obj )
240 {
241   FLA_Bool r_val;
242 
243   if      ( obj.base == NULL )
244     r_val = TRUE;
245   else if ( (obj.base)->buffer == NULL )
246     r_val = TRUE;
247   else
248     r_val = FALSE;
249 
250   return r_val;
251 }
252 
253 
254 
FLA_Obj_is_int(FLA_Obj A)255 FLA_Bool FLA_Obj_is_int( FLA_Obj A )
256 {
257   FLA_Datatype datatype;
258   FLA_Bool     r_val;
259 
260   datatype = FLA_Obj_datatype( A );
261 
262   if ( datatype == FLA_INT )
263     r_val = TRUE;
264   else
265     r_val = FALSE;
266 
267   return r_val;
268 }
269 
270 
271 
FLA_Obj_is_floating_point(FLA_Obj A)272 FLA_Bool FLA_Obj_is_floating_point( FLA_Obj A )
273 {
274   FLA_Datatype datatype;
275   FLA_Bool     r_val;
276 
277   datatype = FLA_Obj_datatype( A );
278 
279   if ( datatype == FLA_FLOAT || datatype == FLA_COMPLEX ||
280        datatype == FLA_DOUBLE || datatype == FLA_DOUBLE_COMPLEX )
281     r_val = TRUE;
282   else
283     r_val = FALSE;
284 
285   return r_val;
286 }
287 
288 
289 
FLA_Obj_is_constant(FLA_Obj A)290 FLA_Bool FLA_Obj_is_constant( FLA_Obj A )
291 {
292   FLA_Datatype datatype;
293   FLA_Bool     r_val;
294 
295   datatype = FLA_Obj_datatype( A );
296 
297   if ( datatype == FLA_CONSTANT )
298     r_val = TRUE;
299   else
300     r_val = FALSE;
301 
302   return r_val;
303 }
304 
305 
306 
FLA_Obj_is_real(FLA_Obj A)307 FLA_Bool FLA_Obj_is_real( FLA_Obj A )
308 {
309   FLA_Datatype datatype;
310   FLA_Bool     r_val;
311 
312   datatype = FLA_Obj_datatype( A );
313 
314   if ( datatype == FLA_CONSTANT || datatype == FLA_FLOAT || datatype == FLA_DOUBLE )
315     r_val = TRUE;
316   else
317     r_val = FALSE;
318 
319   return r_val;
320 }
321 
322 
323 
FLA_Obj_is_complex(FLA_Obj A)324 FLA_Bool FLA_Obj_is_complex( FLA_Obj A )
325 {
326   FLA_Datatype datatype;
327   FLA_Bool     r_val;
328 
329   datatype = FLA_Obj_datatype( A );
330 
331   if ( datatype == FLA_CONSTANT || datatype == FLA_COMPLEX || datatype == FLA_DOUBLE_COMPLEX )
332     r_val = TRUE;
333   else
334     r_val = FALSE;
335 
336   return r_val;
337 }
338 
339 
340 
FLA_Obj_is_single_precision(FLA_Obj A)341 FLA_Bool FLA_Obj_is_single_precision( FLA_Obj A )
342 {
343   FLA_Datatype datatype;
344   FLA_Bool     r_val;
345 
346   datatype = FLA_Obj_datatype( A );
347 
348   if ( datatype == FLA_CONSTANT || datatype == FLA_FLOAT || datatype == FLA_COMPLEX )
349     r_val = TRUE;
350   else
351     r_val = FALSE;
352 
353   return r_val;
354 }
355 
356 
357 
FLA_Obj_is_double_precision(FLA_Obj A)358 FLA_Bool FLA_Obj_is_double_precision( FLA_Obj A )
359 {
360   FLA_Datatype datatype;
361   FLA_Bool     r_val;
362 
363   datatype = FLA_Obj_datatype( A );
364 
365   if ( datatype == FLA_CONSTANT || datatype == FLA_DOUBLE || datatype == FLA_DOUBLE_COMPLEX )
366     r_val = TRUE;
367   else
368     r_val = FALSE;
369 
370   return r_val;
371 }
372 
373 
374 
FLA_Obj_is_scalar(FLA_Obj A)375 FLA_Bool FLA_Obj_is_scalar( FLA_Obj A )
376 {
377   FLA_Bool r_val = FALSE;
378 
379   if ( FLA_Obj_length( A ) == 1 &&
380        FLA_Obj_width( A ) == 1 )
381     r_val = TRUE;
382 
383   return r_val;
384 }
385 
386 
387 
FLA_Obj_is_vector(FLA_Obj A)388 FLA_Bool FLA_Obj_is_vector( FLA_Obj A )
389 {
390   FLA_Bool r_val = FALSE;
391 
392   if ( FLA_Obj_length( A ) == 1 || FLA_Obj_width( A ) == 1 )
393     r_val = TRUE;
394 
395   return r_val;
396 }
397 
398 
399 
FLA_Obj_has_zero_dim(FLA_Obj A)400 FLA_Bool FLA_Obj_has_zero_dim( FLA_Obj A )
401 {
402   FLA_Bool r_val = FALSE;
403 
404   if ( FLA_Obj_length( A ) == 0 || FLA_Obj_width( A ) == 0 )
405     r_val = TRUE;
406 
407   return r_val;
408 }
409 
410 
411 
FLA_Obj_is_col_major(FLA_Obj A)412 FLA_Bool FLA_Obj_is_col_major( FLA_Obj A )
413 {
414   FLA_Bool r_val = FALSE;
415 
416   // A row stride of 1 indicates column-major storage.
417   if ( FLA_Obj_row_stride( A ) == 1 )
418     r_val = TRUE;
419 
420   return r_val;
421 }
422 
423 
424 
FLA_Obj_is_row_major(FLA_Obj A)425 FLA_Bool FLA_Obj_is_row_major( FLA_Obj A )
426 {
427   FLA_Bool r_val = FALSE;
428 
429   // A column stride of 1 indicates row-major storage.
430   if ( FLA_Obj_col_stride( A ) == 1 )
431     r_val = TRUE;
432 
433   return r_val;
434 }
435 
436 
437 
FLA_Obj_is_conformal_to(FLA_Trans trans,FLA_Obj A,FLA_Obj B)438 FLA_Bool FLA_Obj_is_conformal_to( FLA_Trans trans, FLA_Obj A, FLA_Obj B )
439 {
440   FLA_Bool r_val = TRUE;
441 
442   if ( trans == FLA_NO_TRANSPOSE || trans == FLA_CONJ_NO_TRANSPOSE )
443   {
444     if ( FLA_Obj_length( A ) != FLA_Obj_length( B ) ||
445          FLA_Obj_width( A )  != FLA_Obj_width( B ) )
446       r_val = FALSE;
447   }
448   else
449   {
450     if ( FLA_Obj_width( A )  != FLA_Obj_length( B ) ||
451          FLA_Obj_length( A ) != FLA_Obj_width( B ) )
452       r_val = FALSE;
453   }
454 
455   return r_val;
456 }
457 
458 
459 
FLA_Obj_is(FLA_Obj A,FLA_Obj B)460 FLA_Bool FLA_Obj_is( FLA_Obj A, FLA_Obj B )
461 {
462   FLA_Bool r_val = FALSE;
463 
464   if ( A.base == B.base )
465     r_val = TRUE;
466 
467   return r_val;
468 }
469 
FLA_Obj_is_identical(FLA_Obj A,FLA_Obj B)470 FLA_Bool FLA_Obj_is_identical( FLA_Obj A, FLA_Obj B )
471 {
472   FLA_Bool r_val = FALSE;
473 
474   // For LU_piv, if A and B are identical, we do not need copy.
475   // Elemtype should be checked as they can have the same buffer pointer
476   // but elemtype can be either scalar or matrix.
477   if ( A.base != NULL && A.base != NULL )
478     if ( ( A.base == B.base ) ||
479          ( A.base->elemtype == B.base->elemtype  &&
480            A.base->datatype == B.base->datatype ) )
481       if ( FLA_Obj_buffer_at_view( A ) == FLA_Obj_buffer_at_view( B ) )
482         if ( A.m == B.m && A.n == B.n )
483           r_val = TRUE;
484 
485   return r_val;
486 }
487 
FLA_Obj_is_overlapped(FLA_Obj A,FLA_Obj B)488 FLA_Bool FLA_Obj_is_overlapped( FLA_Obj A, FLA_Obj B )
489 {
490   FLA_Bool r_val = FALSE;
491 
492   // For form_Q, if A and B are not overlapped, we do not use in-place forming Q.
493   if ( A.base != NULL && A.base != NULL )
494     if ( ( A.base == B.base ) ||
495          ( A.base->elemtype == B.base->elemtype &&
496            A.base->datatype == B.base->datatype ) )
497       if ( FLA_Obj_buffer_at_view( A ) == FLA_Obj_buffer_at_view( B ) )
498         if ( ( ( A.offm <= B.offm && B.offm < ( A.offm + A.m ) ) &&
499                ( A.offn <= B.offn && B.offn < ( A.offn + A.n ) ) ) ||
500              ( ( B.offm <= A.offm && A.offm < ( B.offm + B.m ) ) &&
501                ( B.offn <= A.offn && A.offn < ( B.offn + B.n ) ) ) )
502           r_val = TRUE;
503 
504   return r_val;
505 }
506 
FLA_Obj_equals(FLA_Obj A,FLA_Obj B)507 FLA_Bool FLA_Obj_equals( FLA_Obj A, FLA_Obj B )
508 {
509   FLA_Datatype datatype_A;
510   FLA_Datatype datatype_B;
511   FLA_Datatype datatype;
512   dim_t        m, n;
513   dim_t        rs_A, cs_A;
514   dim_t        rs_B, cs_B;
515   dim_t        i, j;
516 
517   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
518     FLA_Obj_equals_check( A, B );
519 
520   m      = FLA_Obj_length( A );
521   n      = FLA_Obj_width( A );
522   rs_A   = FLA_Obj_row_stride( A );
523   cs_A   = FLA_Obj_col_stride( A );
524   rs_B   = FLA_Obj_row_stride( B );
525   cs_B   = FLA_Obj_col_stride( B );
526 
527   datatype_A = FLA_Obj_datatype( A );
528   datatype_B = FLA_Obj_datatype( B );
529 
530   // If A is a non-FLA_CONSTANT object, then we should proceed based on the
531   // value of datatype_A. In such a situation, either datatype_B is an exact
532   // match and we're fine, or datatype_B is FLA_CONSTANT, in which case we're
533   // also covered since FLA_CONSTANT encompassas all numerical types.
534   // If A is an FLA_CONSTANT object, then we should proceed based on the value
535   // of datatype_B. In this case, datatype_B is either a non-FLA_CONSTANT type,
536   // which mirrors the second sub-case above, or datatype_B is FLA_CONSTANT,
537   // in which case both types are FLA_CONSTANT and therefore we have to handle
538   // that case. Only if both are FLA_CONSTANTs does the FLA_CONSTANT case
539   // statement below execute.
540   if ( datatype_A != FLA_CONSTANT )
541     datatype = datatype_A;
542   else
543     datatype = datatype_B;
544 
545   switch ( datatype )
546   {
547     case FLA_CONSTANT:
548     {
549       // We require ALL floating-point fields to be the same.
550       float*    buffs_A = ( float    * ) FLA_FLOAT_PTR( A );
551       float*    buffs_B = ( float    * ) FLA_FLOAT_PTR( B );
552       double*   buffd_A = ( double   * ) FLA_DOUBLE_PTR( A );
553       double*   buffd_B = ( double   * ) FLA_DOUBLE_PTR( B );
554       scomplex* buffc_A = ( scomplex * ) FLA_COMPLEX_PTR( A );
555       scomplex* buffc_B = ( scomplex * ) FLA_COMPLEX_PTR( B );
556       dcomplex* buffz_A = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( A );
557       dcomplex* buffz_B = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( B );
558 
559       if ( *buffs_A != *buffs_B ||
560            *buffd_A != *buffd_B ||
561            buffc_A->real != buffc_B->real ||
562            buffc_A->imag != buffc_B->imag ||
563            buffz_A->real != buffz_B->real ||
564            buffz_A->imag != buffz_B->imag )
565       {
566         return FALSE;
567       }
568 
569       break;
570     }
571 
572     case FLA_INT:
573     {
574       int *buff_A = ( int * ) FLA_INT_PTR( A );
575       int *buff_B = ( int * ) FLA_INT_PTR( B );
576 
577       for ( j = 0; j < n; j++ )
578         for ( i = 0; i < m; i++ )
579           if ( buff_A[ j * cs_A + i * rs_A ] !=
580                buff_B[ j * cs_B + i * rs_B ] )
581           {
582             return FALSE;
583           }
584 
585       break;
586     }
587 
588     case FLA_FLOAT:
589     {
590       float *buff_A = ( float * ) FLA_FLOAT_PTR( A );
591       float *buff_B = ( float * ) FLA_FLOAT_PTR( B );
592 
593       for ( j = 0; j < n; j++ )
594         for ( i = 0; i < m; i++ )
595           if ( buff_A[ j * cs_A + i * rs_A ] !=
596                buff_B[ j * cs_B + i * rs_B ] )
597           {
598             return FALSE;
599           }
600 
601       break;
602     }
603 
604     case FLA_DOUBLE:
605     {
606       double *buff_A = ( double * ) FLA_DOUBLE_PTR( A );
607       double *buff_B = ( double * ) FLA_DOUBLE_PTR( B );
608 
609       for ( j = 0; j < n; j++ )
610         for ( i = 0; i < m; i++ )
611           if ( buff_A[ j * cs_A + i * rs_A ] !=
612                buff_B[ j * cs_B + i * rs_B ] )
613           {
614             return FALSE;
615           }
616 
617       break;
618     }
619 
620     case FLA_COMPLEX:
621     {
622       scomplex *buff_A = ( scomplex * ) FLA_COMPLEX_PTR( A );
623       scomplex *buff_B = ( scomplex * ) FLA_COMPLEX_PTR( B );
624 
625       for ( j = 0; j < n; j++ )
626         for ( i = 0; i < m; i++ )
627           if ( buff_A[ j * cs_A + i * rs_A ].real != buff_B[ j * cs_B + i * rs_B ].real ||
628                buff_A[ j * cs_A + i * rs_A ].imag != buff_B[ j * cs_B + i * rs_B ].imag )
629           {
630             return FALSE;
631           }
632 
633       break;
634     }
635 
636     case FLA_DOUBLE_COMPLEX:
637     {
638       dcomplex *buff_A = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( A );
639       dcomplex *buff_B = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( B );
640 
641       for ( j = 0; j < n; j++ )
642         for ( i = 0; i < m; i++ )
643           if ( buff_A[ j * cs_A + i * rs_A ].real != buff_B[ j * cs_B + i * rs_B ].real ||
644                buff_A[ j * cs_A + i * rs_A ].imag != buff_B[ j * cs_B + i * rs_B ].imag )
645           {
646             return FALSE;
647           }
648 
649       break;
650     }
651   }
652 
653   return TRUE;
654 }
655 
656 
657 
FLA_Obj_gt(FLA_Obj A,FLA_Obj B)658 FLA_Bool FLA_Obj_gt( FLA_Obj A, FLA_Obj B )
659 {
660   FLA_Datatype datatype_A;
661   FLA_Datatype datatype_B;
662   FLA_Datatype datatype;
663 
664   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
665     FLA_Obj_gt_check( A, B );
666 
667   datatype_A = FLA_Obj_datatype( A );
668   datatype_B = FLA_Obj_datatype( B );
669 
670   if ( datatype_A != FLA_CONSTANT ) datatype = datatype_A;
671   else                              datatype = datatype_B;
672 
673   switch ( datatype )
674   {
675     case FLA_CONSTANT:
676     {
677       // We require ALL floating-point fields to be the same.
678       float*    buffs_A = ( float    * ) FLA_FLOAT_PTR( A );
679       float*    buffs_B = ( float    * ) FLA_FLOAT_PTR( B );
680       double*   buffd_A = ( double   * ) FLA_DOUBLE_PTR( A );
681       double*   buffd_B = ( double   * ) FLA_DOUBLE_PTR( B );
682       scomplex* buffc_A = ( scomplex * ) FLA_COMPLEX_PTR( A );
683       scomplex* buffc_B = ( scomplex * ) FLA_COMPLEX_PTR( B );
684       dcomplex* buffz_A = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( A );
685       dcomplex* buffz_B = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( B );
686 
687       if ( !( *buffs_A      > *buffs_B      &&
688               *buffd_A      > *buffd_B      &&
689               buffc_A->real > buffc_B->real &&
690               buffc_A->imag > buffc_B->imag &&
691               buffz_A->real > buffz_B->real &&
692               buffz_A->imag > buffz_B->imag ) )
693       {
694         return FALSE;
695       }
696 
697       break;
698     }
699 
700     case FLA_INT:
701     {
702       int *buff_A = ( int * ) FLA_INT_PTR( A );
703       int *buff_B = ( int * ) FLA_INT_PTR( B );
704 
705       if ( !( *buff_A > *buff_B ) ) return FALSE;
706 
707       break;
708     }
709 
710     case FLA_FLOAT:
711     {
712       float *buff_A = ( float * ) FLA_FLOAT_PTR( A );
713       float *buff_B = ( float * ) FLA_FLOAT_PTR( B );
714 
715       if ( !( *buff_A > *buff_B ) ) return FALSE;
716 
717       break;
718     }
719 
720     case FLA_DOUBLE:
721     {
722       double *buff_A = ( double * ) FLA_DOUBLE_PTR( A );
723       double *buff_B = ( double * ) FLA_DOUBLE_PTR( B );
724 
725       if ( !( *buff_A > *buff_B ) ) return FALSE;
726 
727       break;
728     }
729 
730   }
731 
732   return TRUE;
733 }
734 
735 
FLA_Obj_ge(FLA_Obj A,FLA_Obj B)736 FLA_Bool FLA_Obj_ge( FLA_Obj A, FLA_Obj B )
737 {
738   FLA_Datatype datatype_A;
739   FLA_Datatype datatype_B;
740   FLA_Datatype datatype;
741 
742   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
743     FLA_Obj_ge_check( A, B );
744 
745   datatype_A = FLA_Obj_datatype( A );
746   datatype_B = FLA_Obj_datatype( B );
747 
748   if ( datatype_A != FLA_CONSTANT ) datatype = datatype_A;
749   else                              datatype = datatype_B;
750 
751   switch ( datatype )
752   {
753     case FLA_CONSTANT:
754     {
755       // We require ALL floating-point fields to be the same.
756       float*    buffs_A = ( float    * ) FLA_FLOAT_PTR( A );
757       float*    buffs_B = ( float    * ) FLA_FLOAT_PTR( B );
758       double*   buffd_A = ( double   * ) FLA_DOUBLE_PTR( A );
759       double*   buffd_B = ( double   * ) FLA_DOUBLE_PTR( B );
760       scomplex* buffc_A = ( scomplex * ) FLA_COMPLEX_PTR( A );
761       scomplex* buffc_B = ( scomplex * ) FLA_COMPLEX_PTR( B );
762       dcomplex* buffz_A = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( A );
763       dcomplex* buffz_B = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( B );
764 
765       if ( !( *buffs_A      >= *buffs_B      &&
766               *buffd_A      >= *buffd_B      &&
767               buffc_A->real >= buffc_B->real &&
768               buffc_A->imag >= buffc_B->imag &&
769               buffz_A->real >= buffz_B->real &&
770               buffz_A->imag >= buffz_B->imag ) )
771       {
772         return FALSE;
773       }
774 
775       break;
776     }
777 
778     case FLA_INT:
779     {
780       int *buff_A = ( int * ) FLA_INT_PTR( A );
781       int *buff_B = ( int * ) FLA_INT_PTR( B );
782 
783       if ( !( *buff_A >= *buff_B ) ) return FALSE;
784 
785       break;
786     }
787 
788     case FLA_FLOAT:
789     {
790       float *buff_A = ( float * ) FLA_FLOAT_PTR( A );
791       float *buff_B = ( float * ) FLA_FLOAT_PTR( B );
792 
793       if ( !( *buff_A >= *buff_B ) ) return FALSE;
794 
795       break;
796     }
797 
798     case FLA_DOUBLE:
799     {
800       double *buff_A = ( double * ) FLA_DOUBLE_PTR( A );
801       double *buff_B = ( double * ) FLA_DOUBLE_PTR( B );
802 
803       if ( !( *buff_A >= *buff_B ) ) return FALSE;
804 
805       break;
806     }
807 
808   }
809 
810   return TRUE;
811 }
812 
FLA_Obj_lt(FLA_Obj A,FLA_Obj B)813 FLA_Bool FLA_Obj_lt( FLA_Obj A, FLA_Obj B )
814 {
815   FLA_Datatype datatype_A;
816   FLA_Datatype datatype_B;
817   FLA_Datatype datatype;
818 
819   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
820     FLA_Obj_lt_check( A, B );
821 
822   datatype_A = FLA_Obj_datatype( A );
823   datatype_B = FLA_Obj_datatype( B );
824 
825   if ( datatype_A != FLA_CONSTANT ) datatype = datatype_A;
826   else                              datatype = datatype_B;
827 
828   switch ( datatype )
829   {
830     case FLA_CONSTANT:
831     {
832       // We require ALL floating-point fields to be the same.
833       float*    buffs_A = ( float    * ) FLA_FLOAT_PTR( A );
834       float*    buffs_B = ( float    * ) FLA_FLOAT_PTR( B );
835       double*   buffd_A = ( double   * ) FLA_DOUBLE_PTR( A );
836       double*   buffd_B = ( double   * ) FLA_DOUBLE_PTR( B );
837       scomplex* buffc_A = ( scomplex * ) FLA_COMPLEX_PTR( A );
838       scomplex* buffc_B = ( scomplex * ) FLA_COMPLEX_PTR( B );
839       dcomplex* buffz_A = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( A );
840       dcomplex* buffz_B = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( B );
841 
842       if ( !( *buffs_A      < *buffs_B      &&
843               *buffd_A      < *buffd_B      &&
844               buffc_A->real < buffc_B->real &&
845               buffc_A->imag < buffc_B->imag &&
846               buffz_A->real < buffz_B->real &&
847               buffz_A->imag < buffz_B->imag ) )
848       {
849         return FALSE;
850       }
851 
852       break;
853     }
854 
855     case FLA_INT:
856     {
857       int *buff_A = ( int * ) FLA_INT_PTR( A );
858       int *buff_B = ( int * ) FLA_INT_PTR( B );
859 
860       if ( !( *buff_A < *buff_B ) ) return FALSE;
861 
862       break;
863     }
864 
865     case FLA_FLOAT:
866     {
867       float *buff_A = ( float * ) FLA_FLOAT_PTR( A );
868       float *buff_B = ( float * ) FLA_FLOAT_PTR( B );
869 
870       if ( !( *buff_A < *buff_B ) ) return FALSE;
871 
872       break;
873     }
874 
875     case FLA_DOUBLE:
876     {
877       double *buff_A = ( double * ) FLA_DOUBLE_PTR( A );
878       double *buff_B = ( double * ) FLA_DOUBLE_PTR( B );
879 
880       if ( !( *buff_A < *buff_B ) ) return FALSE;
881 
882       break;
883     }
884 
885   }
886 
887   return TRUE;
888 }
889 
FLA_Obj_le(FLA_Obj A,FLA_Obj B)890 FLA_Bool FLA_Obj_le( FLA_Obj A, FLA_Obj B )
891 {
892   FLA_Datatype datatype_A;
893   FLA_Datatype datatype_B;
894   FLA_Datatype datatype;
895 
896   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
897     FLA_Obj_le_check( A, B );
898 
899   datatype_A = FLA_Obj_datatype( A );
900   datatype_B = FLA_Obj_datatype( B );
901 
902   if ( datatype_A != FLA_CONSTANT ) datatype = datatype_A;
903   else                              datatype = datatype_B;
904 
905   switch ( datatype )
906   {
907     case FLA_CONSTANT:
908     {
909       // We require ALL floating-point fields to be the same.
910       float*    buffs_A = ( float    * ) FLA_FLOAT_PTR( A );
911       float*    buffs_B = ( float    * ) FLA_FLOAT_PTR( B );
912       double*   buffd_A = ( double   * ) FLA_DOUBLE_PTR( A );
913       double*   buffd_B = ( double   * ) FLA_DOUBLE_PTR( B );
914       scomplex* buffc_A = ( scomplex * ) FLA_COMPLEX_PTR( A );
915       scomplex* buffc_B = ( scomplex * ) FLA_COMPLEX_PTR( B );
916       dcomplex* buffz_A = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( A );
917       dcomplex* buffz_B = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( B );
918 
919       if ( !( *buffs_A      <= *buffs_B      &&
920               *buffd_A      <= *buffd_B      &&
921               buffc_A->real <= buffc_B->real &&
922               buffc_A->imag <= buffc_B->imag &&
923               buffz_A->real <= buffz_B->real &&
924               buffz_A->imag <= buffz_B->imag ) )
925       {
926         return FALSE;
927       }
928 
929       break;
930     }
931 
932     case FLA_INT:
933     {
934       int *buff_A = ( int * ) FLA_INT_PTR( A );
935       int *buff_B = ( int * ) FLA_INT_PTR( B );
936 
937       if ( !( *buff_A <= *buff_B ) ) return FALSE;
938 
939       break;
940     }
941 
942     case FLA_FLOAT:
943     {
944       float *buff_A = ( float * ) FLA_FLOAT_PTR( A );
945       float *buff_B = ( float * ) FLA_FLOAT_PTR( B );
946 
947       if ( !( *buff_A <= *buff_B ) ) return FALSE;
948 
949       break;
950     }
951 
952     case FLA_DOUBLE:
953     {
954       double *buff_A = ( double * ) FLA_DOUBLE_PTR( A );
955       double *buff_B = ( double * ) FLA_DOUBLE_PTR( B );
956 
957       if ( !( *buff_A <= *buff_B ) ) return FALSE;
958 
959       break;
960     }
961 
962   }
963 
964   return TRUE;
965 }
966 
967 
968 
FLA_Submatrix_at(FLA_Datatype datatype,void * buffer,dim_t i,dim_t j,dim_t rs,dim_t cs)969 void* FLA_Submatrix_at( FLA_Datatype datatype, void* buffer, dim_t i, dim_t j, dim_t rs, dim_t cs )
970 {
971   void* r_val = buffer;
972 
973   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
974     FLA_Submatrix_at_check( datatype, buffer, i, j, rs, cs );
975 
976   switch( datatype )
977   {
978     case FLA_INT:
979       r_val = ( void* ) ( (      ( int* ) buffer ) + i * rs + j * cs );
980       break;
981 
982     case FLA_FLOAT:
983       r_val = ( void* ) ( (    ( float* ) buffer ) + i * rs + j * cs );
984       break;
985 
986     case FLA_DOUBLE:
987       r_val = ( void* ) ( (   ( double* ) buffer ) + i * rs + j * cs );
988       break;
989 
990     case FLA_COMPLEX:
991       r_val = ( void* ) ( ( ( scomplex* ) buffer ) + i * rs + j * cs );
992       break;
993 
994     case FLA_DOUBLE_COMPLEX:
995       r_val = ( void* ) ( ( ( dcomplex* ) buffer ) + i * rs + j * cs );
996       break;
997   }
998 
999   return r_val;
1000 }
1001 
FLA_Obj_has_nan(FLA_Obj A)1002 FLA_Bool FLA_Obj_has_nan( FLA_Obj A )
1003 {
1004   FLA_Datatype datatype;
1005   dim_t        i, j, m, n, cs, rs;
1006 
1007   if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
1008     FLA_Obj_has_nan_check( A );
1009 
1010   datatype = FLA_Obj_datatype( A );
1011   m        = FLA_Obj_length( A );
1012   n        = FLA_Obj_width( A );
1013   cs       = FLA_Obj_col_stride( A );
1014   rs       = FLA_Obj_row_stride( A );
1015 
1016   switch ( datatype )
1017   {
1018     case FLA_FLOAT:
1019     {
1020       float *buff = ( float * ) FLA_FLOAT_PTR( A );
1021 
1022       for ( j=0; j<n; ++j )
1023         for ( i=0; i<m; ++i )
1024         {
1025           float val = buff[i*cs + j*rs];
1026           if ( val != val ) return TRUE;
1027         }
1028       break;
1029     }
1030     case FLA_DOUBLE:
1031     {
1032       double *buff = ( double * ) FLA_DOUBLE_PTR( A );
1033 
1034       for ( j=0; j<n; ++j )
1035         for ( i=0; i<m; ++i )
1036         {
1037           double val = buff[i*cs + j*rs];
1038           if ( val != val ) return TRUE;
1039         }
1040       break;
1041     }
1042     case FLA_COMPLEX:
1043     {
1044       scomplex *buff = ( scomplex * ) FLA_COMPLEX_PTR( A );
1045 
1046       for ( j=0; j<n; ++j )
1047         for ( i=0; i<m; ++i )
1048         {
1049           scomplex val = buff[i*cs + j*rs];
1050           if ( val.real != val.real || val.imag != val.imag ) return TRUE;
1051         }
1052       break;
1053     }
1054     case FLA_DOUBLE_COMPLEX:
1055     {
1056       dcomplex *buff = ( dcomplex * ) FLA_DOUBLE_COMPLEX_PTR( A );
1057 
1058       for ( j=0; j<n; ++j )
1059         for ( i=0; i<m; ++i )
1060         {
1061           dcomplex val = buff[i*cs + j*rs];
1062           if ( val.real != val.real || val.imag != val.imag ) return TRUE;
1063         }
1064       break;
1065     }
1066   }
1067 
1068   return FALSE;
1069 }
1070 
1071