1function C = bitset (A, B, arg3, arg4)
2%BITSET set bit.
3% C = bitset (A,B) sets a bit in A to 1, where the bit position is
4% determined by B.  A is an integer array.  If B(i,j) is an integer in the
5% range 1 (the least significant bit) to the number of bits in the data
6% type of A, then C(i,j) is equal to the value of A(i,j) after setting the
7% bit to 1.  If B(i,j) is outside this range, C(i,j) is set to A(i,j),
8% unmodified; note that this behavior is an extension of the built-in
9% MATLAB bigset, which results in an error for this case.  This modified
10% rule allows the inputs A and B to be sparse.
11%
12% If A and B are matrices, the pattern of C is the set union of A
13% and B.  If one of A or B is a nonzero scalar, the scalar is expanded
14% into a sparse matrix with the same pattern as the other matrix, and the
15% result is a sparse matrix.
16%
17% If the last input argument is a string, C = bigset (A,B,assumedtype)
18% provides a data type to convert A to if it has a floating-point type.
19% If A already has an integer type, then it is not modified.  Otherwise, A
20% is converted to assumedtype, which can be 'int8', 'int16', 'int32',
21% 'int64', 'uint8', 'uint16', 'uint32' or 'uint64'.  The default is
22% 'uint64'.
23%
24% C = bitset (A,B,V) sets the bit in A(i,j) at position B(i,j) to 0 if
25% V(i,j) is zero, or to 1 if V(i,j) is nonzero.  If V is a scalar, it
26% is implicitly expanded to V * spones (B).
27%
28% All four arguments may be used, as C = bitset (A,B,V,assumedtype).
29%
30% Example:
31%
32%   A = GrB (magic (4), 'uint8')
33%   B = reshape ([1:8 1:8], 4, 4)
34%   C = bitset (A, B)
35%   fprintf ('\nA: ') ; fprintf ('%3x ', A) ; fprintf ('\n') ;
36%   fprintf ('\nB: ') ; fprintf ('%3x ', B) ; fprintf ('\n') ;
37%   fprintf ('\nC: ') ; fprintf ('%3x ', C) ; fprintf ('\n') ;
38%   % in MATLAB:
39%   C2 = bitset (uint8 (A), B)
40%   isequal (C2, C)
41%
42% See also GrB/bitor, GrB/bitand, GrB/bitxor, GrB/bitcmp, GrB/bitshift,
43% GrB/bitset, GrB/bitclr.
44
45% SuiteSparse:GraphBLAS, Timothy A. Davis, (c) 2017-2021, All Rights Reserved.
46% SPDX-License-Identifier: GPL-3.0-or-later
47
48% FUTURE: bitset(A,B,V) for two matrices A and B is slower than it could be.
49% See comments in gb_union_op.
50
51if (isobject (A))
52    A = A.opaque ;
53end
54
55if (isobject (B))
56    B = B.opaque ;
57end
58
59[am, an, atype] = gbsize (A) ;
60[bm, bn, btype] = gbsize (B) ;
61
62if (contains (atype, 'complex') || contains (btype, 'complex'))
63    error ('inputs must be real') ;
64end
65
66if (isequal (atype, 'logical') || isequal (btype, 'logical'))
67    error ('inputs must not be logical') ;
68end
69
70a_is_scalar = (am == 1) && (an == 1) ;
71b_is_scalar = (bm == 1) && (bn == 1) ;
72
73% get the optional input arguments
74if (nargin == 4)
75    V = arg3 ;
76    assumedtype = arg4 ;
77elseif (nargin == 3)
78    if (ischar (arg3))
79        V = 1 ;
80        assumedtype = arg3 ;
81    else
82        V = arg3 ;
83        assumedtype = 'uint64' ;
84    end
85else
86    V = 1 ;
87    assumedtype = 'uint64' ;
88end
89
90if (~contains (assumedtype, 'int'))
91    error ('assumedtype must be an integer type') ;
92end
93
94% C will have the same type as A on input
95ctype = atype ;
96
97% determine the type of A
98if (isequal (atype, 'double') || isequal (atype, 'single'))
99    A = gbnew (A, assumedtype) ;
100    atype = assumedtype ;
101end
102
103% ensure B has the same type as A
104if (~isequal (btype, atype))
105    B = gbnew (B, atype) ;
106end
107
108% get the matrix or scalar V
109if (isobject (V))
110    V = V.opaque ;
111end
112[m, n] = gbsize (V) ;
113V_is_scalar = (m == 1) && (n == 1) ;
114
115if (V_is_scalar)
116
117    % V is a scalar:  all bits in A indexed by B are either cleared or set.
118    if (gb_scalar (V) == 0)
119        % any bit reference by B(i,j) is set to 0 in A
120        op = ['bitclr.' atype] ;
121    else
122        % any bit reference by B(i,j) is set to 1 in A
123        op = ['bitset.' atype] ;
124    end
125
126    if (a_is_scalar)
127        % A is a scalar
128        if (b_is_scalar)
129            % both A and B are scalars
130            C = gb_union_op (op, A, B) ;
131        else
132            % A is a scalar, B is a matrix
133            C = gbapply2 (op, A, B) ;
134        end
135    else
136        % A is a matrix
137        if (b_is_scalar)
138            % A is a matrix, B is scalar
139            C = gbapply2 (op, A, B) ;
140        else
141            % both A and B are matrices
142            C = gb_union_op (op, A, B) ;
143        end
144    end
145
146else
147
148    % V is a matrix: A and B can be scalars or matrices, but if they
149    % are matrices, they must have the same size as V.
150
151    % if B(i,j) is nonzero and V(i,j)=1, then:
152    % C(i,j) = bitset (A (i,j), B (i,j)).
153
154    % if B(i,j) is nonzero and V(i,j)=0 (implicit or explicit), then:
155    % C(i,j) = bitclr (A (i,j), B (i,j)).
156
157    if (a_is_scalar)
158        % expand A to a full matrix the same size as V.
159        A = gb_scalar_to_full (m, n, atype, gb_fmt (V), A) ;
160    end
161    if (b_is_scalar)
162        % expand B to a full matrix the same size as V.
163        B = gb_scalar_to_full (m, n, atype, gb_fmt (V), B) ;
164    end
165
166    % Set all bits referenced by B(i,j) to 1, even those that need to be
167    % set to 0, without considering V(i,j).
168    C = gb_union_op (['bitset.', atype], A, B) ;
169
170    % The pattern of C is now the set intersection of A and B, but
171    % bits referenced by B(i,j) have been set to 1, not 0.  Construct B0
172    % as the bits in B(i,j) that must be set to 0; B0<~V>=B defines the
173    % pattern of bit positions B0 to set to 0 in A.
174    d.mask = 'complement' ;
175    B0 = gbassign (gbnew (m, n, atype), V, B, d) ;
176
177    % Clear the bits in C, referenced by B0(i,j), where V(i,j) is zero.
178    C = gbeadd (['bitclr.', atype], C, B0) ;
179
180end
181
182% return result
183if (isequal (gbtype (C), ctype))
184    C = GrB (C) ;
185else
186    C = GrB (gbnew (C, ctype)) ;
187end
188
189
190