[05b386]: src / __gmtimes__.cc Maximize Restore History

Download this file

__gmtimes__.cc    108 lines (91 with data), 3.1 kB

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
/*
Copyright (C) 2011 David Bateman
This program is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2, or (at your option)
any later version.
This program is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
General Public License for more details.
You should have received a copy of the GNU General Public License
along with this software; see the file COPYING. If not, see
<http://www.gnu.org/licenses/>.
In addition to the terms of the GPL, you are permitted to link
this program with any Open Source program, as defined by the
Open Source Initiative (www.opensource.org)
*/
#include <iostream>
#include <iomanip>
#include <sstream>
#include <octave/oct.h>
#include <octave/pager.h>
#include <octave/quit.h>
static inline octave_idx_type
modn(octave_idx_type x, octave_idx_type m, octave_idx_type n)
{
while (x >= n)
{
x -= n;
x = (x >> m) + (x & n);
}
return x;
}
DEFUN_DLD (__gmtimes__, args, nargout,
"-*- texinfo -*-\n"
"@deftypefnx {Loadable Function} {@var{y} =} __mtimes__ (@var{a}, @var{b}, @var{m}, @var{alpha_to}, @var{index_of})\n"
"\n"
"Undocumented internal function\n"
"@end deftypefn")
{
octave_value retval;
int nargin = args.length ();
if (nargin != 5)
{
print_usage ();
return retval;
}
const NDArray a = args(0).array_value();
const NDArray b = args(1).array_value();
octave_idx_type m = args(2).idx_type_value();
octave_idx_type n =
(static_cast<octave_idx_type>(1) << m) - static_cast<octave_idx_type>(1);
const Array<octave_idx_type> alpha_to =
args(3).octave_idx_type_vector_value();
const Array<octave_idx_type> index_of =
args(4).octave_idx_type_vector_value();
octave_idx_type a_nr = a.rows();
octave_idx_type a_nc = a.cols();
octave_idx_type b_nr = b.rows();
octave_idx_type b_nc = b.cols();
if (a_nc != b_nr)
{
gripe_nonconformant ("operator *", a_nr, a_nc, b_nr, b_nc);
}
else
{
NDArray y(dim_vector(a_nr, b_nc), 0.0);
if (a_nr != 0 && a_nc != 0 && b_nc != 0)
{
// This is not optimum for referencing b, but can use vector
// to represent index(a(k,j)). Seems to be the fastest.
Array<octave_idx_type> c(dim_vector(a_nr, 1));
for (octave_idx_type j = 0; j < b_nr; j++) {
for (octave_idx_type k = 0; k < a_nr; k++)
c.xelem(k) = index_of(static_cast<octave_idx_type>(a.xelem(k,j)));
for (octave_idx_type i = 0; i < b_nc; i++)
if (static_cast<octave_idx_type>(b(j,i)) != 0)
{
octave_idx_type tmp = index_of(static_cast<octave_idx_type>(b.xelem(j,i)));
for (octave_idx_type k = 0; k < a_nr; k++)
{
if (static_cast<octave_idx_type>(a.xelem(k,j)) != 0)
y.xelem(k,i) = static_cast<octave_idx_type>(y.xelem(k,i)) ^ alpha_to(modn(tmp + c.xelem(k), m, n));
}
}
}
}
retval = y;
}
return retval;
}