<< Chapter < Page | Chapter >> Page > |
and
When this factored form of or any of the variations alluded to above, is used, the number of additions incurred is given by
where .
Although the use of operations of the form is simple, it does not exactly separate the circular convolution into smaller disjoint convolutions.In other words, its use does not give rise in [link] and [link] to block diagonal matrices whose diagonalblocks are the form . However, by reorganizing the arrangement of theoperations we can obtain the block diagonal form we seek.
First, suppose , and are matrices of sizes , and respectively. If
where and are matrices of sizes and , then
where
Here denotes the first rows and all the columns of and similarly for . Note that
That these two expressions are not equal explains why the arrangement of operations must be reorganizedin order to obtain the desired block diagonal form. The appropriate reorganization is described by theexpression in [link] . Therefore, we must modify the transformation of [link] appropriately. It should be noted that this reorganization ofoperations does not change their computational cost. It is still given by [link] .
For example, we can use this observation and the expression in [link] to arrive at the following similarity transformation:
where
is a column vector of 1's
and is the matrix:
In general we have
where is given by
with , and
and are as given in [link] and [link] .
Notice the distinction between this example and example "Reduction Operations" . In example "Reduction Operations" we obtained from a Kronecker product of two block diagonal matrices, but here we obtaineda block diagonal matrix whose diagonal blocks are the Kronecker product of cyclotomic companionmatrices. Each block in [link] represents a multi-dimensional cyclotomic convolution.
A Matlab program that carries out the operation
in
[link] is
KRED()
.
function x = KRED(P,E,K,x)
% x = KRED(P,E,K,x);% P : P = [P(1),...,P(K)];% E : E = [E(K),...,E(K)];for i = 1:K
a = prod(P(1:i-1).^E(1:i-1));c = prod(P(i+1:K).^E(i+1:K));
p = P(i);e = E(i);
for j = e-1:-1:0x(1:a*c*(p^(j+1))) = RED(p,a,c*(p^j),x(1:a*c*(p^(j+1))));
endend
It calls the Matlab program
RED()
.
function y = RED(p,a,c,x)
% y = RED(p,a,c,x);y = zeros(a*c*p,1);
for i = 0:c:(a-1)*cfor j = 0:c-1
y(i+j+1) = x(i*p+j+1);for k = 0:c:c*(p-2)
y(i+j+1) = y(i+j+1) + x(i*p+j+k+c+1);y(i*(p-1)+j+k+a*c+1) = x(i*p+j+k+1) - x(i*p+j+c*(p-1)+1);
endend
end
These two Matlab programs are not written to run as fast as they could be.They are a `naive' coding of and are meant to serve as a basis for more efficient programs. In particular, the indexing and the loop counters canbe modified to improve the efficiency. However, the modifications that minimize the overheadincurred by indexing operations depends on the programming language, the compiler and the computer used.These two programs are written with simple loop counters and complicated indexing operationsso that appropriate modifications can be easily made.
The inverse of has the form
and
Because the inverse of in [link] is given by
the inverse of the matrix described by eqs [link] , [link] and [link] is given by
with , and
where
denotes the matrix in
[link] without the first column.
A Matlab program for
is
tKRED()
,
it calls the Matlab program
tRED()
.
A Matlab program for
is
itKRED()
,
it calls the Matlab program
itRED()
.
These programs all appear in one of the appendices.
Notification Switch
Would you like to follow the 'Automatic generation of prime length fft programs' conversation and receive update notifications?