Kronecker product: Difference between revisions

</pre>

=={{header|Fortran}}==

The plan is to pass the two arrays to a subroutine, which will return their Kronecker product as a third parameter. This relies on the expanded array-handling facilities introduced with F90, especially the ability of a subroutine to allocate an array of a size of its choosing, this array being a parameter to the subroutine. Some compilers offering the "allocate" statement do not handle this! Further features of the MODULE protocol of F90 allow arrays passed to a subroutine to have their sizes ascertained in the subroutine (via function UBOUND, ''etc.'') rather than this information being supplied via additional parameters. This is not all to the good: multi-dimensional arrays must therefore be the actual size of their usage rather than say A(100,100) but only using the first fifty elements (in one place) and the first thirty in another. Thus, for such usage the array must be re-allocated the correct size each time, and, the speed of access to such arrays is reduced - see [[Sequence_of_primorial_primes#Fixed-size_data_aggregates]] for an example. A further detail of the MODULE protocol when passing arrays is that if the parameter's declaration does not specify the lower bound, it will be treated as if it were one even if the actual array is declared otherwise - see [[Array_length#Fortran]] for example.

In older-style Fortran, the arrays would be of some "surely-big-enough" size, fixed at compile time, and there would be additional parameters describing the bounds in use for each invocation. Since no array-assignment statements were available, there would be additional DO-loops to copy each block of values. In all versions of Fortran, the ordering of array elements in storage is "column-major" so that the DATA statement appears to initialise the arrays with their transpose - see [[Matrix_transposition#Fortran]] for example. As a result, the default output order for an array, if written as just <code>WRITE (6,*) A</code> will be that of the transposed order, just as with the default order of the DATA statement's data. To show the desired order, the array must be written with explicit specification of the order of elements, as done by subroutine SHOW: columns across the page, rows running down the page. <lang Fortran> MODULE ARRAYMUSH !A rather small collection.

CONTAINS !For the specific problem only.

SUBROUTINE KPRODUCT(A,B,AB) !AB = Kronecker product of A and B, both two-dimensional arrays.

Considers the arrays to be addressed as A(row,column), despite any storage order arrangements. .

Creating array AB to fit here, adjusting the caller's array AB, may not work on some compilers.

INTEGER A(:,:),B(:,:) !Two-dimensional arrays, lower bound one.

INTEGER, ALLOCATABLE:: AB(:,:) !To be created to fit.

INTEGER R,RA,RB,C,CA,CB,I !Assistants.

RA = UBOUND(A,DIM = 1) !Ascertain the upper bounds of the incoming arrays.

CA = UBOUND(A,DIM = 2) !Their lower bounds will be deemed one,

RB = UBOUND(B,DIM = 1) !And the upper bound as reported will correspond.

CB = UBOUND(B,DIM = 2) !UBOUND(A) would give an array of two values, RA and CA, more for higher dimensionality.

WRITE (6,1) "A",RA,CA,"B",RB,CB,"A.k.B",RA*RB,CA*CB !Announce.

1 FORMAT (3(A," is ",I0,"x",I0,1X)) !Three sets of sizes.

IF (ALLOCATED(AB)) DEALLOCATE(AB) !Discard any lingering storage.

ALLOCATE (AB(RA*RB,CA*CB)) !Obtain the exact desired size.

R = 0 !Syncopation: start the row offset.

DO I = 1,RA !Step down the rows of A.

C = 0 !For each row, start the column offset.

DO J = 1,CA !Step along the columns of A.

AB(R + 1:R + RB,C + 1:C + CB) = A(I,J)*B !Place a block of B values.

C = C + CB !Advance a block of columns.

END DO !On to the next column of A.

R = R + RB !Advance a block of rows.

END DO !On to the next row of A.

END SUBROUTINE KPRODUCT !No tests for bad parametes, or lack of staorage...

SUBROUTINE SHOW(F,A) !Write array A in row,column order.

INTEGER F !Output file unit number.

INTEGER A(:,:) !The 2-D array, lower bound one.

INTEGER R !The row stepper.

DO R = 1,UBOUND(A,DIM = 1) !Each row gets its own line.

WRITE (F,1) A(R,:) !Write all the columns of that row.

1 FORMAT (666I3) !This suffices for the example.

END DO !On to the next row.

END SUBROUTINE SHOW !WRITE (F,*) A or similar would show A as if transposed.

END MODULE ARRAYMUSH !That was simple enough.

PROGRAM POKE

USE ARRAYMUSH

INTEGER A(2,2),B(2,2) !First test: square arrays.

INTEGER, ALLOCATABLE:: AB(:,:) !To be created for each result.

INTEGER C(3,3),D(3,4) !Second test: some rectilinearity.

DATA A/1,3, 2,4/,B/0,6, 5,7/ !Furrytran uses "column-major" order for successive storage elements.

DATA C/0,1,0, 1,1,1, 0,1,0/ !So, the first three values go down the rows of the first column.

DATA D/1,1,1, 1,0,1, 1,0,1, 1,1,1/!And then follow the values for the next column, etc.

WRITE (6,*) "First test..."

CALL KPRODUCT(A,B,AB)

CALL SHOW (6,AB)

WRITE (6,*)

WRITE (6,*) "Second test..."

CALL KPRODUCT(C,D,AB)

CALL SHOW (6,AB)

END</lang>

Output:

<pre>

First test...

A is 2x2 B is 2x2 A.k.B is 4x4

0 5 0 10

6 7 12 14

0 15 0 20

18 21 24 28

Second test...

A is 3x3 B is 3x4 A.k.B is 9x12

0 0 0 0 1 1 1 1 0 0 0 0

0 0 0 0 1 0 0 1 0 0 0 0

0 0 0 0 1 1 1 1 0 0 0 0

1 1 1 1 1 1 1 1 1 1 1 1

1 0 0 1 1 0 0 1 1 0 0 1

1 1 1 1 1 1 1 1 1 1 1 1

0 0 0 0 1 1 1 1 0 0 0 0

0 0 0 0 1 0 0 1 0 0 0 0

0 0 0 0 1 1 1 1 0 0 0 0