Re: Dynamic mutidimensional arrays

[Michael Riedl]

Andreas,

sorry to say that but I am still under the impression that you do not distguish between a code example and and a full library implementation.

In the samples I could have called e.g. Lapack or Blas or part of my own huger libraies - but that does not make the undelaying principle clear.

Again, a forgotten feature in ISO Modula-2 versus a bit spoiled and unreadbly way to code simple things.

What I see here as samples is even a way back behing Algol60 published more than half a century ago.

If what I can read in this discussen is the solution I would ask to get my problem back :-)

Gruß

Michael

Am 05.04.23 um 08:25 schrieb Fischlin Andreas:

Dear all,

Yes, you can implement this in PIM and ISO and need no extra feature. Michael, as I tried to express earlier, I fear you are hung up in doing the client side of things and the implementation in the same go, i.e. all in a program module (e.g. your M2 and Fortran pseudo code samples). Instead I advocate for splitting the task into the servicing part, which goes into a library module, and the client part, the program module where you use that functionality. Yes, this requires a bit more of coding and a design effort, but the overall result is more robust, offers reusability and is safe to use. The way to do it in M2.

Benjamin, I prefer our approach (cf. LgMatrices) over what you sketched as being less limited (MaxElems being maximally large, with LMatrix = POINTER TO LMat and can allocate the maximum storage for any given compiler that is addressable via CARDINAL*). And more importantly, it avoids the messy address arithmetic – I really advise against – as the notation when accessing elements or vectors of a matrix is quite familiar: E.g. vec[i]-th element is denoted by vec^[i] or mat(i,j)-th element is denoted by mat^[i]^[j] ). This notation helps when it comes to implement more complex calculations, e.g. when computing eigen values etc. I would not want to implement that with the address arithmetic your approach would require.

Andreas

* I must admit that there are some aspects of LgMatrices specific coding I do not like. E.g. type of procedure arguments nRows or nCols should be of course at least CARDINAL and not INTEGER as given in procedures such as AllocMatrix or AllocVector. Even better would be to use TYPE Index = 1..MaxElems;

ETH Zurich

Prof. em. Dr. Andreas Fischlin

IPCC Vice-Chair WGII

Systems Ecology - Institute of Biogeochemistry and Pollutant Dynamics

CHN E 24

Universitaetstrasse 16

8092 Zurich

SWITZERLAND

andreas.fischlin@env.ethz.ch

www.sysecol.ethz.ch/people/andreas.fischlin.html

+41 44 633-6090 phone

+41 44 633-1136 fax

+41 79 595-4050 mobile

             Make it as simple as possible, but distrust it!

________________________________________________________________________

On Wed, 05.04.23, at 01:52, Benjamin Kowarsch <trijezdci@gmail.com> wrote:

On Wed, 5 Apr 2023 at 04:43, Michael Riedl <udo-michael.riedl@t-online.de> wrote:

actually I do not see your approach to be "easy" at all.

So please can you line out how the short sample I lined out should look like ?

I'll give you the public interface to the entire library instead.

With that you can allocate, initialise, introspect and operate on arbitrarily sized matrices.

DEFINITION MODULE Matrix;

IMPORT Vector;

FROM Vector IMPORT VectorT;

TYPE Matrix; (* OPAQUE *)

TYPE MatrixT = Matrix; (* for unqualified use *)

(* allocator, passes back null matrix *)

PROCEDURE New ( VAR m : Matrix; rows, columns : CARDINAL );

(* introspection *)

PROCEDURE rowCount ( m : Matrix ) : CARDINAL;

PROCEDURE columnCount ( m : Matrix ) : CARDINAL;

PROCEDURE isInvertible ( m : Matrix ) : BOOLEAN;

(* constructors *)

PROCEDURE NewWithRows ( VAR m : Matrix; rowValues : ARRAY OF VectorT );

PROCEDURE NewWithColumns ( VAR m : Matrix; colValues : ARRAY OF VectorT );

(* mutators *)

PROCEDURE StoreValue ( m : Matrix; row, col : CARDINAL; value : REAL );

PROCEDURE StoreRow ( m : Matrix; rowIndex : CARDINAL; rowValues : ARRAY OF VectorT );

PROCEDURE StoreColumn ( m : Matrix; colIndex : CARDINAL; colValues : ARRAY OF VectorT );

(* accessors *)

PROCEDURE value ( m : Matrix; row, col : CARDINAL ) : REAL;

PROCEDURE row ( m : Matrix; rowIndex : CARDINAL ) : VectorT;

PROCEDURE column ( m : Matrix; colIndex : CARDINAL ) : Vector T;

(* operations *)

PROCEDURE eq ( m1, m2 : Matrix ) : BOOLEAN;

PROCEDURE sum ( m1, m2 : Matrix ) : Matrix;

PROCEDURE diff ( m1, m2 : Matrix ) : Matrix;

PROCEDURE prod ( m1, m2 : Matrix ) : Matrix;

PROCEDURE scalarProd ( m : Matrix; scalar : REAL ) : Matrix;

PROCEDURE vectorProd ( m : Matrix; vector : VectorT ) : VectorT;

PROCEDURE inverse ( m : Matrix ) : Matrix;

PROCEDURE transpose ( m : Matrix ) : Matrix;

(* destructor, deallocates and passes back NIL *)

PROCEDURE Release ( VAR m : Matrix );

END Matrix.

The Vector type module follows the same pattern, so I omitted it.

I believe you can do the implementation for this yourself, so I omitted that, too.

Note that you can implement this in PIM and ISO without the need for any additional features. However, since neither PIM nor ISO provide any facility to declare a record field of indeterminate size, you need to use an ADDRESS field and then use address arithmetic to access the individual values.

TYPE Matrix = RECORD

  rows, columns : CARDINAL;

  capacity : LONGCARD;

  data : ADDRESS

END; (* Matrix *)

Doing the address arithmetic will look very messy and represents plenty of opportunity for error. In this way it is even less safe than doing it in C. But it can be done and the unholy mess is hidden in the implementation module. Once you've got it working, you will only ever use the public interface.

By contrast, with the facility of records with indeterminate fields that I described does not rely on address arithmetic, it allows the use of array subscripts to address the values of the array field, clean and safe, as it should be in Modula-2.

(* type declaration in implementation module *)

TYPE Matrix = POINTER TO RECORD

  rows, columns : CARDINAL

+ data : ARRAY capacity OF REAL

END; (* Matrix *)

I hope this clarifies.

regards

benjamin