matrix manipulation

matrix manipulation

[hnptz]

What is an elegant way to program exchange of columns in a matrix. 
Application: 
C is a matrix. Cj be the j-th column of C.
Sort C1,C2, ..., Cn in non increasing order with respect to their first component.
I tried and ended up in a very messy program, as I had to sort the first components first, remember their new index and then reconstruct the matrix.

Re: matrix manipulation

[Jacob Sparre Andersen]

hnptz@yahoo.de writes:

> What is an elegant way to program exchange of columns in a matrix.

One way could be to declare a matrix type as an ordered set of columns,
and insert the columns one at a time.

Greetings,

Jacob
-- 
No trees were killed in the sending of this message.
However a large number of electrons were terribly inconvenienced.

Re: matrix manipulation

[Stephen Leake]

On Saturday, November 12, 2016 at 12:16:47 PM UTC-6, hn...@yahoo.de wrote:
> What is an elegant way to program exchange of columns in a matrix. 
> Application: 
> C is a matrix. Cj be the j-th column of C.
> Sort C1,C2, ..., Cn in non increasing order with respect to their first component.
> I tried and ended up in a very messy program, as I had to sort the first components first, remember their new index and then reconstruct the matrix.

declare the matrix as an array of columns:

type Column is array (1 .. 10) of float;
type Matrix is array (1 .. 10) of Column;

A : Matrix := ...;
Tmp : Column;

Tmp := A (2);
A (2) := A (1);
A (1) := Tmp;

Re: matrix manipulation

[hnptz]

On Sunday, November 13, 2016 at 12:30:16 AM UTC+1, Stephen Leake wrote:
> On Saturday, November 12, 2016 at 12:16:47 PM UTC-6, hn...@yahoo.de wrote:
> > What is an elegant way to program exchange of columns in a matrix. 
> > Application: 
> > C is a matrix. Cj be the j-th column of C.
> > Sort C1,C2, ..., Cn in non increasing order with respect to their first component.
> > I tried and ended up in a very messy program, as I had to sort the first components first, remember their new index and then reconstruct the matrix.
> 
> declare the matrix as an array of columns:
> 
> type Column is array (1 .. 10) of float;
> type Matrix is array (1 .. 10) of Column;
> 
> A : Matrix := ...;
> Tmp : Column;
> 
> Tmp := A (2);
> A (2) := A (1);
> A (1) := Tmp;

How to access a single element of matrix A ? Is it A(1)(1) ? I think, I cannot use it as a Matrix anymore applying build in Ada-procedures, like for solving  equations or transposing the matrix. I would have to transform it into a 2-dimensional array, right ?





On Sunday, November 13, 2016 at 12:30:16 AM UTC+1, Stephen Leake wrote:
> On Saturday, November 12, 2016 at 12:16:47 PM UTC-6, hn...@yahoo.de wrote:
> > What is an elegant way to program exchange of columns in a matrix. 
> > Application: 
> > C is a matrix. Cj be the j-th column of C.
> > Sort C1,C2, ..., Cn in non increasing order with respect to their first component.
> > I tried and ended up in a very messy program, as I had to sort the first components first, remember their new index and then reconstruct the matrix.
> 
> declare the matrix as an array of columns:
> 
> type Column is array (1 .. 10) of float;
> type Matrix is array (1 .. 10) of Column;
> 
> A : Matrix := ...;
> Tmp : Column;
> 
> Tmp := A (2);
> A (2) := A (1);
> A (1) := Tmp;

Re: matrix manipulation

[Stephen Leake]

On Monday, November 14, 2016 at 11:39:56 AM UTC-6, hn...@yahoo.de wrote:
> On Sunday, November 13, 2016 at 12:30:16 AM UTC+1, Stephen Leake wrote:
> > On Saturday, November 12, 2016 at 12:16:47 PM UTC-6, hn...@yahoo.de wrote:
> > > What is an elegant way to program exchange of columns in a matrix. 
> > > Application: 
> > > C is a matrix. Cj be the j-th column of C.
> > > Sort C1,C2, ..., Cn in non increasing order with respect to their first component.
> > > I tried and ended up in a very messy program, as I had to sort the first components first, remember their new index and then reconstruct the matrix.
> > 
> > declare the matrix as an array of columns:
> > 
> > type Column is array (1 .. 10) of float;
> > type Matrix is array (1 .. 10) of Column;
> > 
> > A : Matrix := ...;
> > Tmp : Column;
> > 
> > Tmp := A (2);
> > A (2) := A (1);
> > A (1) := Tmp;
> 
> How to access a single element of matrix A ? Is it A(1)(1) ? 

Yes.

> I think, I cannot use it as a Matrix anymore applying build in 
> Ada-procedures, like for solving  equations or transposing the matrix. I 
> would have to transform it into a 2-dimensional array, right ?

I guess you mean the type Ada.Numerics.Generic_Complex_Arrays.Complex_Matrix, or something similar.

It helps to have _all_ of the requirements at hand when designing a solution :)

You have two choices here; write code two swap two rows of a 2D array, or copy the 2D array to an array of arrays, swap, and copy back. You'd have to measure to be sure which is more efficient, and the result will depend on matrix size.

For an initial design, it is best to go with code that is easiest to understand and implement; get it right first, then measure it, then improve only if necessary.

Re: matrix manipulation

[Robert Eachus]

On Saturday, November 12, 2016 at 1:16:47 PM UTC-5, hn...@yahoo.de wrote:
> What is an elegant way to program exchange of columns in a matrix. 
> Application: 
> C is a matrix. Cj be the j-th column of C.
> Sort C1,C2, ..., Cn in non increasing order with respect to their first component.
> I tried and ended up in a very messy program, as I had to sort the first components first, remember their new index and then reconstruct the matrix.

I don't consider the first part messy.  I have lots of code lying around that creates an index of a set of records based on a greater than function which is passed in.  Reusing the generic would require aliasing the array as a matrix, easier just to change the spec and recompile--none of the generic code depends on the contents of the record.

The second part is more interesting.  If you have a large array and even if you don't it is more efficient to move columns to the correct location. You are going to need to copy one or more columns to a temporary location to do the copying, but you can manage to copy most columns once.

After you have the index array walk through it.  If the index matches the current location go on.  If it doesn't copy the row at the new location to the temporary vector (and keep track of where it goes), copy the current location to the emptied location, change its index to 0 or the right index (your choice) and start an inner loop.

The inner loop recursively walks the indexes until it finds itself a zero (or identity if you took that approach). Now unwind the recursion moving columns as you go.  Are you finished?  You can keep a count (remember to include columns not moved along with those moved to the temporary vector) or just continue through the indexes looking for another column that needs moving.

matrix manipulation

[hnptz]

Hi,

I would like to solve the following Ada-Programming demand:

bk :=3; Dv := 4;

H : Real_Matrix(1..bk*Dv,1..bk*Dv);

type Block_Matrix
   is array (Positive range <>, Positive range <>) of Real_Matrix(1..Dv,1..Dv);

HB : Block_Matrix(1..bk,1..bk);

Now, there are 12 Block_Matrices which make up H. After calculation of all HB's I want to place them into H and ended up in coding of 8 different for loops. Is there a more elegant solution? I noticed that slicing would only work with a definition of linear arrays of a linear array, which does not fit for the rest of my program.

Re: matrix manipulation

[mockturtle]

On Saturday, April 27, 2019 at 4:50:13 PM UTC+2, hn...@yahoo.de wrote:
> Hi,
> 
> I would like to solve the following Ada-Programming demand:
> 
> bk :=3; Dv := 4;
> 
> H : Real_Matrix(1..bk*Dv,1..bk*Dv);
> 
> type Block_Matrix
>    is array (Positive range <>, Positive range <>) of Real_Matrix(1..Dv,1..Dv);
> 
> HB : Block_Matrix(1..bk,1..bk);
> 
> Now, there are 12 Block_Matrices which make up H. After calculation of all HB's I want to place them into H and ended up in coding of 8 different for loops. Is there a more elegant solution? I noticed that slicing would only work with a definition of linear arrays of a linear array, which does not fit for the rest of my program.

Why 8?  I would say 4: two to select the block and two to select the element inside the block.  You could use just 2: iterate over the rows and columns of H and compute (using mod and division) the block where the element of H goes.

Re: matrix manipulation

[hnptz]

On Saturday, April 27, 2019 at 4:50:13 PM UTC+2, hn...@yahoo.de wrote:
> Hi,
> 
> I would like to solve the following Ada-Programming demand:
> 
> bk :=3; Dv := 4;
> 
> H : Real_Matrix(1..bk*Dv,1..bk*Dv);
> 
> type Block_Matrix
>    is array (Positive range <>, Positive range <>) of Real_Matrix(1..Dv,1..Dv);
> 
> HB : Block_Matrix(1..bk,1..bk);
> 
> Now, there are 12 Block_Matrices which make up H. After calculation of all HB's I want to place them into H and ended up in coding of 8 different for loops. Is there a more elegant solution? I noticed that slicing would only work with a definition of linear arrays of a linear array, which does not fit for the rest of my program.

Sorry, I meant six. Two for H, two for the block and two inside the block;

H(i,j) := HB(i1,j1)(i2,j2);

The mod solution sounds interesting. I then might switch to indices starting with zero, which may be easier to handle in connection with the mod function.

Re: matrix manipulation

[J-P. Rosen]

Le 27/04/2019 à 17:21, hnptz@yahoo.de a écrit :
> The mod solution sounds interesting. I then might switch to indices
> starting with zero, which may be easier to handle in connection with
> the mod function.

It just means that you'll end up doing by hand the index calculation
that the compiler would do for you with nested loops.

Replacing multiple loops with a single loop buys you nothing as long as
the total number of iterations is the same...
-- 
J-P. Rosen
Adalog
2 rue du Docteur Lombard, 92441 Issy-les-Moulineaux CEDEX
Tel: +33 1 45 29 21 52, Fax: +33 1 45 29 25 00
http://www.adalog.fr

Re: matrix manipulation

[mockturtle]

On Saturday, April 27, 2019 at 6:37:45 PM UTC+2, J-P. Rosen wrote:
> Le 27/04/2019 à 17:21, hnptz@yahoo.de a écrit :
> > The mod solution sounds interesting. I then might switch to indices
> > starting with zero, which may be easier to handle in connection with
> > the mod function.
> 
> It just means that you'll end up doing by hand the index calculation
> that the compiler would do for you with nested loops.

Correct me if I am wrong, but it seems to me that you would end to do some calculation anyway

* If you iterate over rows and column of the full matrix (that is, the non-block matrix) you need to compute block position and position inside the block with mod/division.  This requires two loops

* If you iterate over blocks and row/column inside the block you need to compute the corresponding position inside the full matrix.  This requires 4 loops

* To be honest, I cannot imagine how you would use 6 loops: the 2 indexes  used by the full matrix are in bijective correspondence with the 4 indexes used by the block matrix, therefore you loop over one matrix and compute the position in the other.

I guess it is more or less a matter of taste and readability. Personally, I would _slightly_ prefer the 2-loop solution because I do not like too many nested loops.


> 
> Replacing multiple loops with a single loop buys you nothing as long as
> the total number of iterations is the same...
> -- 
> J-P. Rosen
> Adalog
> 2 rue du Docteur Lombard, 92441 Issy-les-Moulineaux CEDEX
> Tel: +33 1 45 29 21 52, Fax: +33 1 45 29 25 00
> http://www.adalog.fr

Re: matrix manipulation

[Simon Wright]

hnptz@yahoo.de writes:

> Hi,
>
> I would like to solve the following Ada-Programming demand:
>
> bk :=3; Dv := 4;
>
> H : Real_Matrix(1..bk*Dv,1..bk*Dv);
>
> type Block_Matrix
>    is array (Positive range <>, Positive range <>) of Real_Matrix(1..Dv,1..Dv);
>
> HB : Block_Matrix(1..bk,1..bk);
>
> Now, there are 12 Block_Matrices which make up H. After calculation of
> all HB's I want to place them into H and ended up in coding of 8
> different for loops. Is there a more elegant solution? I noticed that
> slicing would only work with a definition of linear arrays of a linear
> array, which does not fit for the rest of my program.

You might find this package helpful. I haven't been able to prove the
postconditions, due to lack of SPARK smarts, but see the demo.

--

package Block_Matrices
  with SPARK_Mode
is

   pragma Assertion_Policy (Check);

   type Matrix is array (Positive range <>, Positive range <>) of Integer;

   function Submatrix (Of_Matrix : Matrix;
                       At_Row    : Positive;
                       At_Column : Positive;
                       Rows      : Natural;
                       Columns   : Natural) return Matrix
   with
     Pre => At_Row in Of_Matrix'Range (1)
            and At_Column in Of_Matrix'Range (2)
            and At_Row + (Rows - 1) in Of_Matrix'Range (1)
            and At_Column + (Columns - 1) in Of_Matrix'Range (2),
     Post => Submatrix'Result'First (1) = At_Row
             and Submatrix'Result'Length (1) = Rows
             and Submatrix'Result'First (2) = At_Column
             and Submatrix'Result'Length (2) = Columns
             and (for all J in Submatrix'Result'Range (1) =>
                    (for all K in Submatrix'Result'Range (2) =>
                       Submatrix'Result (J, K) = Of_Matrix (J, K)));

   procedure Assign (Submatrix : Matrix; To : in out Matrix)
   with
     Pre => Submatrix'First (1) in To'Range (1)
            and Submatrix'Last (1) in To'Range (1)
            and Submatrix'First (2) in To'Range (2)
            and Submatrix'Last (2) in To'Range (2),
     Post => (for all J in Submatrix'Range (1) =>
                (for all K in Submatrix'Range (2) =>
                   To (J, K) = Submatrix (J, K)));

end Block_Matrices;

package body Block_Matrices
with SPARK_Mode
is

   function Submatrix (Of_Matrix : Matrix;
                       At_Row    : Positive;
                       At_Column : Positive;
                       Rows      : Natural;
                       Columns   : Natural) return Matrix
   is
   begin
      return Result : Matrix
        (At_Row .. At_Row - 1 + Rows, At_Column .. At_Column - 1 + Columns) do
         for J in Result'Range (1) loop
            for K in Result'Range (2) loop
               Result (J, K) := Of_Matrix (J, K);
            end loop;
         end loop;
      end return;
   end Submatrix;

   procedure Assign (Submatrix : Matrix; To : in out Matrix)
   is
   begin
      for J in Submatrix'Range (1) loop
         for K in Submatrix'Range (2) loop
            To (J, K) := Submatrix (J, K);
         end loop;
      end loop;
   end Assign;

end Block_Matrices;

with Ada.Text_IO; use Ada.Text_IO;
procedure Block_Matrices.Test
is

   procedure Print (M : Matrix) is
   begin
      for R in M'Range (1) loop
         for C in M'Range (2) loop
            Put (" " & M (R, C)'Image);
         end loop;
         New_Line;
      end loop;
   end Print;

   M : Matrix (1 .. 8, 1 .. 8);

begin

   for J in M'Range (1) loop
      for K in M'Range (2) loop
         M (J, K) := J * K;
      end loop;
   end loop;

   declare
      S : Matrix := Submatrix (M, 3, 5, 2, 3);
   begin
      Print (S);
      New_Line;
      for El of S loop
         El := - El;
      end loop;
      Assign (S, M);
      Print (M);
   end;

end Block_Matrices.Test;