• Real_Arrays on heap with overloaded operators and clean syntax

    From Jim Paloander@21:1/5 to All on Sun Jan 22 13:34:18 2023
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in runtime)
    and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Joakim Strandberg@21:1/5 to All on Sun Jan 22 13:56:27 2023
    söndag 22 januari 2023 kl. 22:34:19 UTC+1 skrev dhmos.a...@gmail.com:
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in runtime)
    and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!

    Easiest solution is probably to declare a new task and specify the stack size using the Storage_Size aspect. Allocate as much stack space as you need to be able to do the calculations and do all the allocations on the declared task, not on the
    environment task. You will avoid the unnecessary heap allocations and have nice clean syntax.

    Best regards,
    Joakim

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to All on Sun Jan 22 14:07:51 2023
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in runtime)
    and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!
    Easiest solution is probably to declare a new task and specify the stack size using the Storage_Size aspect. Allocate as much stack space as you need to be able to do the calculations and do all the allocations on the declared task, not on the
    environment task. You will avoid the unnecessary heap allocations and have nice clean syntax.

    Best regards,
    Joakim

    Thank you for your reply,
    since I am a newbie I was under the impression that tasks are used only when you want to write a parallel code that takes advantage of multicore architectures. You suggest I have a single task and single thread something like this? I see, but there
    should be a way to do this also for the main program. But thanks anyway. Are you aware of any libraries similar to Real_Arrays, but who allocated memory internally using heap? This is the natural way to do such things. Similarly to the Containers.Vector.
    But Vector has such an awful syntax. There should be something like an indexer [i] similarly to the C++ std::vector to make things simpler and overloaded operators similarly to Real_Arrays. It is a no brainer. Most programs need to allocate on the heap,
    why did they restrict Real_Arrays on the stack?

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Dmitry A. Kazakov@21:1/5 to Jim Paloander on Sun Jan 22 23:13:14 2023
    On 2023-01-22 22:34, Jim Paloander wrote:
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in runtime)
    and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?

    You can define "+" on the access type, which should probably be an arena pointer for performance reasons:

    Arena : Mark_And_Release_Pool;
    type Real_Vector_Ptr is access Real_Vector;
    for Real_Vector_Ptr'Storage_Pool use Arena;

    function "+" (Left, Right : Real_Vector_Ptr)
    return Real_Vector_Ptr is
    begin
    if Left'Length /= Right'Length then
    raise Constraint_Error;
    end if;
    return Result : Real_Vector_Ptr := new Real_Vector (Left'Range) do
    for I in Result'Range loop
    Result (I) :=
    Left (I) + Right (I - Left'First + Right'First);
    end loop;
    end return;
    end "+";

    You can overload that with

    function "+" (Left : Real_Vector_Ptr; Right : Real_Vector)
    return Real_Vector_Ptr is
    begin
    if Left'Length /= Right'Length then
    raise Constraint_Error;
    end if;
    return Result : Real_Vector_Ptr := new Real_Vector (Left'Range) do
    for I in Result'Range loop
    Result (I) :=
    Left (I) + Right (I - Left'First + Right'First);
    end loop;
    end return;
    end "+";

    and with

    function "+" (Left : Real_Vector; Right : Real_Vector_Ptr)
    return Real_Vector_Ptr;

    Then you will be able to write:

    X := A + B + C + C + A + B;
    -- Use X
    Free (X); -- Pop all arena garbage

    But of course, the optimal way to work large linear algebra problems is
    by using in-place operations, e.g.

    procedure Add (Left : in out Real_Vector; Right : Real_Vector);

    etc.

    --
    Regards,
    Dmitry A. Kazakov
    http://www.dmitry-kazakov.de

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Joakim Strandberg@21:1/5 to All on Sun Jan 22 14:42:50 2023
    söndag 22 januari 2023 kl. 23:07:53 UTC+1 skrev dhmos.a...@gmail.com:
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in
    runtime) and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!
    Easiest solution is probably to declare a new task and specify the stack size using the Storage_Size aspect. Allocate as much stack space as you need to be able to do the calculations and do all the allocations on the declared task, not on the
    environment task. You will avoid the unnecessary heap allocations and have nice clean syntax.

    Best regards,
    Joakim
    Thank you for your reply,
    since I am a newbie I was under the impression that tasks are used only when you want to write a parallel code that takes advantage of multicore architectures. You suggest I have a single task and single thread something like this? I see, but there
    should be a way to do this also for the main program. But thanks anyway. Are you aware of any libraries similar to Real_Arrays, but who allocated memory internally using heap? This is the natural way to do such things. Similarly to the Containers.Vector.
    But Vector has such an awful syntax. There should be something like an indexer [i] similarly to the C++ std::vector to make things simpler and overloaded operators similarly to Real_Arrays. It is a no brainer. Most programs need to allocate on the heap,
    why did they restrict Real_Arrays on the stack?

    It my impression that in the Ada community the preferred way of working is in general stack only. Heap allocations are avoided for a number of reasons for example performance, the application needs to ask the operating system for memory which one doesn't
    know how much time that will take nor if it always will succeed. In addition, applications that use the heap may be susceptible to heap memory fragmentation. In Ada, it is easy to specify stack sizes when declaring tasks. It is not part of the Ada
    standard to specify stack size of the environment task. The Ada way is to declare new tasks and do work on them. Prefer the bounded containers over the unbounded containers. If you really need to allocate objects I recommend using storage pools with pre-
    allocated memory also known as arena pools.

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to All on Sun Jan 22 14:36:28 2023
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in runtime)
    and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    You can define "+" on the access type, which should probably be an arena pointer for performance reasons:

    Arena : Mark_And_Release_Pool;
    type Real_Vector_Ptr is access Real_Vector;
    for Real_Vector_Ptr'Storage_Pool use Arena;

    function "+" (Left, Right : Real_Vector_Ptr)
    return Real_Vector_Ptr is
    begin
    if Left'Length /= Right'Length then
    raise Constraint_Error;
    end if;
    return Result : Real_Vector_Ptr := new Real_Vector (Left'Range) do
    for I in Result'Range loop
    Result (I) :=
    Left (I) + Right (I - Left'First + Right'First);
    end loop;
    end return;
    end "+";

    You can overload that with

    function "+" (Left : Real_Vector_Ptr; Right : Real_Vector)
    return Real_Vector_Ptr is
    begin
    if Left'Length /= Right'Length then
    raise Constraint_Error;
    end if;
    return Result : Real_Vector_Ptr := new Real_Vector (Left'Range) do
    for I in Result'Range loop
    Result (I) :=
    Left (I) + Right (I - Left'First + Right'First);
    end loop;
    end return;
    end "+";

    and with

    function "+" (Left : Real_Vector; Right : Real_Vector_Ptr)
    return Real_Vector_Ptr;

    Then you will be able to write:

    X := A + B + C + C + A + B;
    -- Use X
    Free (X); -- Pop all arena garbage

    But of course, the optimal way to work large linear algebra problems is
    by using in-place operations, e.g.

    procedure Add (Left : in out Real_Vector; Right : Real_Vector);

    etc.

    --
    Regards,
    Dmitry A. Kazakov
    http://www.dmitry-kazakov.de
    Privet Dmitry,
    thank you for your answer, I was thinking about that, but I was not sure whether or not it can be avoided with Implicit_Dereference,

    type Accessor (Data: not null access Element) is limited private
    with Implicit_Dereference => Data;

    somehow... Otherwise what you described for operator+ one has to do for every operator overloaded inside Real_Arrays package. The optimal way to work large linear algebra problem is what you describe because unfortunately ADA does not allow what Fortran
    does since 30 years ago or more. But in C++ you can reproduce the same functionality as Fortran using Expression Templates and Template Metaprogramming. Perhaps ADA should allow something like that. Because for maintainability reasons the best would be
    to write the mathematical expressions as close as possible to the mathematical formulas.

    Spacebo Bolshoe,
    Jim.

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Joakim Strandberg@21:1/5 to All on Sun Jan 22 15:11:00 2023
    söndag 22 januari 2023 kl. 23:49:11 UTC+1 skrev dhmos.a...@gmail.com:
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in
    runtime) and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!
    Easiest solution is probably to declare a new task and specify the stack size using the Storage_Size aspect. Allocate as much stack space as you need to be able to do the calculations and do all the allocations on the declared task, not on the
    environment task. You will avoid the unnecessary heap allocations and have nice clean syntax.

    Best regards,
    Joakim
    Thank you for your reply,
    since I am a newbie I was under the impression that tasks are used only when you want to write a parallel code that takes advantage of multicore architectures. You suggest I have a single task and single thread something like this? I see, but there
    should be a way to do this also for the main program. But thanks anyway. Are you aware of any libraries similar to Real_Arrays, but who allocated memory internally using heap? This is the natural way to do such things. Similarly to the Containers.Vector.
    But Vector has such an awful syntax. There should be something like an indexer [i] similarly to the C++ std::vector to make things simpler and overloaded operators similarly to Real_Arrays. It is a no brainer. Most programs need to allocate on the heap,
    why did they restrict Real_Arrays on the stack?
    It my impression that in the Ada community the preferred way of working is in general stack only. Heap allocations are avoided for a number of reasons for example performance, the application needs to ask the operating system for memory which one
    doesn't know how much time that will take nor if it always will succeed. In addition, applications that use the heap may be susceptible to heap memory fragmentation. In Ada, it is easy to specify stack sizes when declaring tasks. It is not part of the
    Ada standard to specify stack size of the environment task. The Ada way is to declare new tasks and do work on them. Prefer the bounded containers over the unbounded containers. If you really need to allocate objects I recommend using storage pools with
    pre-allocated memory also known as arena pools.
    With great depression I realized that the preferred way is of stack only. This is very restrictive excluding all scientific modelling involving solvers for partial differential equations, linear algebra kernels, etc. It is insane. Completely insane. 3D
    simulations of physical phenomena may involve billions of grid-cells and at each grid-cell several unknowns are defined (velocity, pressure, temperature, energy, density, etc). That is why they are using Fortran or C++, but ADA has really cool stuff for
    so many things, why not vectors and matrices and heap allocation? Would you please give me an example, I googled and I cannot find a single example demonstrating how to use a task with the declaration of stack size. Why is there so little information
    online about so important things such as allocation?

    For what you described above being able to write
    X := A + B + C + C + A + B;
    working with a bigger stack is the easiest solution. Don't understand why you think it would restrict all scientific modelling? Maybe I haven't described it well enough.

    From the top of my head:
    task T with Storage_Size => 10_000_000 is

    end T;

    task body T is
    begin
    null;
    end T;

    Tasks need to be defined inside Ada packages, they cannot be stand-alone.

    An optional type is for example in Ada:
    type Optional_Integer (Exists : Boolean := False) is record
    case Exists is
    when True => Value : Integer;
    when False => null;
    end Exists;
    end record;

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to All on Sun Jan 22 14:49:09 2023
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in
    runtime) and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!
    Easiest solution is probably to declare a new task and specify the stack size using the Storage_Size aspect. Allocate as much stack space as you need to be able to do the calculations and do all the allocations on the declared task, not on the
    environment task. You will avoid the unnecessary heap allocations and have nice clean syntax.

    Best regards,
    Joakim
    Thank you for your reply,
    since I am a newbie I was under the impression that tasks are used only when you want to write a parallel code that takes advantage of multicore architectures. You suggest I have a single task and single thread something like this? I see, but there
    should be a way to do this also for the main program. But thanks anyway. Are you aware of any libraries similar to Real_Arrays, but who allocated memory internally using heap? This is the natural way to do such things. Similarly to the Containers.Vector.
    But Vector has such an awful syntax. There should be something like an indexer [i] similarly to the C++ std::vector to make things simpler and overloaded operators similarly to Real_Arrays. It is a no brainer. Most programs need to allocate on the heap,
    why did they restrict Real_Arrays on the stack?
    It my impression that in the Ada community the preferred way of working is in general stack only. Heap allocations are avoided for a number of reasons for example performance, the application needs to ask the operating system for memory which one doesn'
    t know how much time that will take nor if it always will succeed. In addition, applications that use the heap may be susceptible to heap memory fragmentation. In Ada, it is easy to specify stack sizes when declaring tasks. It is not part of the Ada
    standard to specify stack size of the environment task. The Ada way is to declare new tasks and do work on them. Prefer the bounded containers over the unbounded containers. If you really need to allocate objects I recommend using storage pools with pre-
    allocated memory also known as arena pools.

    With great depression I realized that the preferred way is of stack only. This is very restrictive excluding all scientific modelling involving solvers for partial differential equations, linear algebra kernels, etc. It is insane. Completely insane. 3D
    simulations of physical phenomena may involve billions of grid-cells and at each grid-cell several unknowns are defined (velocity, pressure, temperature, energy, density, etc). That is why they are using Fortran or C++, but ADA has really cool stuff for
    so many things, why not vectors and matrices and heap allocation? Would you please give me an example, I googled and I cannot find a single example demonstrating how to use a task with the declaration of stack size. Why is there so little information
    online about so important things such as allocation?

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Gautier write-only address@21:1/5 to All on Sun Jan 22 15:14:03 2023
    Note that Real_Arrays does not specify where things are allocated (heap or stack).
    Only when you define "x : Real_Vector (1 .. n)", it is on stack. You can always write something like the snippet below.
    Anyway, after a certain size, you may have to find compromises, like avoiding operators (they do too many allocations & deallocations in the background, even assuming elegant heap-allocated objects) and also give up plain matrices, against sparse
    matrices or band-stored matrices, typically for solving Partial Differential Equations.

    with Ada.Numerics.Generic_Real_Arrays;

    procedure Test_Large is

    type Float_15 is digits 15;

    package F15_R_A is new Ada.Numerics.Generic_Real_Arrays (Float_15);

    use F15_R_A;

    procedure Solve_it
    (x : in Real_Vector;
    y : out Real_Vector;
    A : in Real_Matrix) is
    begin
    null; -- Here, the big number-crunching
    end;

    n : constant := 10_000;

    type Vector_Access is access Real_Vector;
    type Matrix_Access is access Real_Matrix;

    x, y : Vector_Access := new Real_Vector (1 .. n);
    A : Matrix_Access := new Real_Matrix (1 .. n, 1 .. n);

    begin
    Solve_it (x.all, y.all, A.all);
    -- !! Deallocation here
    end;

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Rod Kay@21:1/5 to Jim Paloander on Mon Jan 23 10:34:47 2023
    On 23/1/23 10:20, Jim Paloander wrote:
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in runtime)
    and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!
    If you are on linux, then you could set the stack size with

    $ ulimit -s unlimited
    $ launch_my_app



    Regards.
    On Windows 10 with mingw64?

    Not sure. I don't have a windows machine.

    What happens when try ?

    $ ulimit -a

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Rod Kay@21:1/5 to Jim Paloander on Mon Jan 23 10:18:52 2023
    On 23/1/23 08:34, Jim Paloander wrote:
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in runtime)
    and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!



    If you are on linux, then you could set the stack size with

    $ ulimit -s unlimited
    $ launch_my_app



    Regards.

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to All on Sun Jan 22 15:20:12 2023
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in runtime)
    and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!
    If you are on linux, then you could set the stack size with

    $ ulimit -s unlimited
    $ launch_my_app



    Regards.
    On Windows 10 with mingw64?

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Joakim Strandberg@21:1/5 to All on Sun Jan 22 15:53:13 2023
    måndag 23 januari 2023 kl. 00:34:31 UTC+1 skrev roda...@gmail.com:
    On 23/1/23 10:20, Jim Paloander wrote:
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in
    runtime) and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!
    If you are on linux, then you could set the stack size with

    $ ulimit -s unlimited
    $ launch_my_app



    Regards.
    On Windows 10 with mingw64?
    Not sure. I don't have a windows machine.

    What happens when try ?

    $ ulimit -a

    Something came up and I had to send my previous reply/e-mail as is. I wanted to find the video where Jean Pierre Rosen talks about how memory is handled in the Ada language from FOSDEM perhaps 2018-2019. Unfortunately I have been unable to find it.

    Secondly, example of Ada code where Storage pool usage is demonstrated: https://github.com/joakim-strandberg/advent_of_code

    Example of task with rendez-vous mechanism:
    task T with Storage_Size => 10_000_000 is
    entry Run (I : Integer);
    end T;

    task body T is
    begin
    loop
    begin
    select
    accept Run (I : Integer) do
    null; -- Here save the value of the integer I for later processing in the task
    end Run;
    null; -- Whenever another task calls Run on this task, the work to be done should be put here.
    -- Put the mathematical calcuations here.
    or
    terminate; -- To end the select-statement with "or terminate;" means the task will terminate when the environment task has finished execution of the "Main" procedure of the application. No need to tell the task T that now it is time to
    shutdown.
    end select;
    exception
    when Error : others =>
    -- Print error to standard out here using subprograms from Ada.Exceptions and Ada.Text_IO.
    end;
    end loop;
    end T;

    Best regards,
    Joakim

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Leo Brewin@21:1/5 to All on Mon Jan 23 12:14:47 2023
    Here is a slight variation on the solution suggested by Gautier. It uses
    Aad's "rename" syntax so that you can avoid all the .all stuff. I use
    this construction extensively in my large scale scientific computations.

    with Ada.Numerics.Generic_Real_Arrays;
    with Ada.Unchecked_Deallocation;

    procedure Test_Large is

    type Float_15 is digits 15;

    package F15_R_A is new Ada.Numerics.Generic_Real_Arrays (Float_15);

    use F15_R_A;

    procedure Solve_it
    (x : in Real_Vector;
    y : out Real_Vector;
    A : in Real_Matrix) is
    begin
    null; -- Here, the big number-crunching
    end;

    n : constant := 10_000;

    type Vector_Access is access Real_Vector;
    type Matrix_Access is access Real_Matrix;

    x_ptr, y_ptr : Vector_Access := new Real_Vector (1 .. n);
    A_ptr : Matrix_Access := new Real_Matrix (1 .. n, 1 .. n);

    x : Real_Vector renames x_ptr.all;
    y : Real_Vector renames y_ptr.all;
    A : Real_Matrix renames A_ptr.all;

    procedure FreeVector is new
    Ada.Unchecked_Deallocation (Real_Vector,Vector_Access);
    procedure FreeMatrix is new
    Ada.Unchecked_Deallocation (Real_Matrix,Matrix_Access);

    begin
    Solve_it (x, y, A);
    -- Deallocation here
    FreeVector (x_ptr);
    FreeVector (y_ptr);
    FreeMatrix (A_ptr);
    end;

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to Leo Brewin on Sun Jan 22 22:01:58 2023
    On Monday, January 23, 2023 at 2:14:58 AM UTC+1, Leo Brewin wrote:
    Here is a slight variation on the solution suggested by Gautier. It uses Aad's "rename" syntax so that you can avoid all the .all stuff. I use
    this construction extensively in my large scale scientific computations.

    with Ada.Numerics.Generic_Real_Arrays;
    with Ada.Unchecked_Deallocation;
    procedure Test_Large is

    type Float_15 is digits 15;

    package F15_R_A is new Ada.Numerics.Generic_Real_Arrays (Float_15);

    use F15_R_A;

    procedure Solve_it
    (x : in Real_Vector;
    y : out Real_Vector;
    A : in Real_Matrix) is
    begin
    null; -- Here, the big number-crunching
    end;

    n : constant := 10_000;

    type Vector_Access is access Real_Vector;
    type Matrix_Access is access Real_Matrix;
    x_ptr, y_ptr : Vector_Access := new Real_Vector (1 .. n);
    A_ptr : Matrix_Access := new Real_Matrix (1 .. n, 1 .. n);

    x : Real_Vector renames x_ptr.all;
    y : Real_Vector renames y_ptr.all;
    A : Real_Matrix renames A_ptr.all;

    procedure FreeVector is new
    Ada.Unchecked_Deallocation (Real_Vector,Vector_Access);
    procedure FreeMatrix is new
    Ada.Unchecked_Deallocation (Real_Matrix,Matrix_Access);

    begin
    Solve_it (x, y, A);
    -- Deallocation here
    FreeVector (x_ptr);
    FreeVector (y_ptr);
    FreeMatrix (A_ptr);
    end;

    Thank you very much, would a for Real_Vector_Access'Storage_Pool use localPool; save you from the need to free the vectors and matrix yourself?

    On the other hand, is there any way of avoiding temporaries? Possibly a modern version of the Real_Array using expression generic syntax or something similar? Since you are using scientific computationas extensively, you must be aware of Fortran. Have
    you compared Fortran's complex numbers with ADA's for inner products or similar computations to see who is faster? You see, I like a lot of things about ADA, but the syntax is really difficult to follow. Sometimes it gives me the impression that it is
    more difficult than really needed to be. For example there should be a way for Real_Arrays to allocate memory internally and not on the stack directly like containers. And for containers to provide an indexer Vector(i) and overloaded operators similarly
    to Real_Vectors. But the fact that they do not give me the impression that this Language although being designed by the army for mission critical applications never realized that modern mission critical need to simplify mathematical calculations
    providing an easy syntax. I am surprised that after so many years and so many updates to the Standard the designers of the Language did not realized that such mathematical syntax should be simplified to attract more people from scientific computing, who
    are tired with Fortran 10000 ways of declaring something a variable.

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to All on Sun Jan 22 22:56:38 2023
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I >>>> get STACK_OVERFLOW ERROR while trying to check how fast operator
    overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I >>>> lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using
    Real_Arrays for implementing some linear algebra method who relies
    heavilly on matrix vector products and vector updates, you do need
    to allocate on the heap (sizes are determined in runtime) and you do >>>> need a clean syntax. So, is there any way to simplify my life
    without using the .all or even without declaring A,B,C,X as access
    Real_Vector?
    Thanks for your time!
    If you are on linux, then you could set the stack size with

    $ ulimit -s unlimited
    $ launch_my_app



    Regards.
    On Windows 10 with mingw64?

    Not sure. I don't have a windows machine.

    What happens when try ?

    $ ulimit -a
    ulimit is available on cygwin.

    It is not available on mingw64 then ?
    It is, but I am not sure if it works since -D400m -d400m together with ulimit were not able to solve the problem. So I thought that under windows 10 and mingw things related to stack size settings are not equivalent to Linux. I checked your repository.
    Looks good. We need some elegant solution, because sometimes it is not easy to predict the stack size that will be needed. Imagine you write a commercial application that is running from a GUI and the user loads a particular case from a file. The memory
    that will be needed may even exceed the available memory installed on the computer. There must be some elegant way to set the stack size globally. On the other hand storage pool seem that does not need deallocation. Isnt there any storage pool for stack
    based allocation so that you enjoy best of both worlds?

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Rod Kay@21:1/5 to Rod Kay on Mon Jan 23 17:34:06 2023
    On 23/1/23 10:34, Rod Kay wrote:
    On 23/1/23 10:20, Jim Paloander wrote:
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I
    get STACK_OVERFLOW ERROR while trying to check how fast operator
    overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I
    lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using
    Real_Arrays for implementing some linear algebra method who relies
    heavilly on matrix vector products and vector updates, you do need
    to allocate on the heap (sizes are determined in runtime) and you do
    need a clean syntax. So, is there any way to simplify my life
    without using the .all or even without declaring A,B,C,X as access
    Real_Vector?
    Thanks for your time!
    If you are on linux, then you could set the stack size with

    $ ulimit -s unlimited
    $ launch_my_app



    Regards.
    On Windows 10 with mingw64?

       Not sure. I don't have a windows machine.

       What happens when try ?

          $ ulimit -a


    ulimit is available on cygwin.

    It is not available on mingw64 then ?

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Rod Kay@21:1/5 to All on Mon Jan 23 18:31:25 2023
    If you are on linux, then you could set the stack size with

    $ ulimit -s unlimited
    $ launch_my_app



    Regards.
    On Windows 10 with mingw64?

    Not sure. I don't have a windows machine.

    What happens when try ?

    $ ulimit -a
    ulimit is available on cygwin.

    It is not available on mingw64 then ?
    It is, but I am not sure if it works since -D400m -d400m together with ulimit were not able to solve the problem. So I thought that under windows 10 and mingw things related to stack size settings are not equivalent to Linux. I checked your repository.
    Looks good. We need some elegant solution, because sometimes it is not easy to predict the stack size that will be needed. Imagine you write a commercial application that is running from a GUI and the user loads a particular case from a file. The memory
    that will be needed may even exceed the available memory installed on the computer. There must be some elegant way to set the stack size globally. On the other hand storage pool seem that does not need deallocation. Isnt there any storage pool for stack
    based allocation so that you enjoy best of both worlds?

    Afraid I'm not familiar with '-D400m -d400m'. I did see a few examples
    on the net which seemed to set stack size as a linker parameter.

    My repository ? Which one ?


    Regards.

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Egil H H@21:1/5 to joak...@kth.se on Sun Jan 22 23:50:11 2023
    On Monday, January 23, 2023 at 12:53:15 AM UTC+1, joak...@kth.se wrote:
    Something came up and I had to send my previous reply/e-mail as is. I wanted to find the video where Jean Pierre Rosen talks about how memory is handled in the Ada language from FOSDEM perhaps 2018-2019. Unfortunately I have been unable to find it.

    It was in 2016:

    https://archive.fosdem.org/2016/schedule/event/ada_memory/

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From G.B.@21:1/5 to Jim Paloander on Mon Jan 23 09:39:39 2023
    On 22.01.23 23:07, Jim Paloander wrote:

    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    X := A + B + C + C + A + B, where
    A,B,C,X are all Real_Vector ( 1 .. N ).

    So my only option was to allocate on the heap using new. But then I lost the clean syntax

    X := A + B + C + C + A + B

    and I had to write instead:

    X.all := A.all + B.all + C.all + C.all + A.all + B.all.

    This is really ugly and annoying because when you are using Real_Arrays for implementing some linear algebra method who relies heavilly on matrix vector products and vector updates, you do need to allocate on the heap (sizes are determined in runtime)
    and you do need a clean syntax. So, is there any way to simplify my life without using the .all or even without declaring A,B,C,X as access Real_Vector?
    Thanks for your time!
    Easiest solution is probably to declare a new task and specify the stack size using the Storage_Size aspect. Allocate as much stack space as you need to be able to do the calculations and do all the allocations on the declared task, not on the
    environment task. You will avoid the unnecessary heap allocations and have nice clean syntax.

    Best regards,
    Joakim

    Thank you for your reply,
    since I am a newbie I was under the impression that tasks are used only when you want to write a parallel code that takes advantage of multicore architectures.

    Tasks---their name is an indication that they just
    perform something that is specified---are not restricted to
    parallelism, or to multiple cores. The can do the job on
    a single processor, or move to some other, they are more
    abstract than contemporary hardware, or hardware for that
    matter.

    You suggest I have a single task and single thread something like this? I see, but there should be a way to do this also for the main program.

    It is unfortunate that GCC had to introduce a tiny scale
    stack(frame) model some time ago.

    But thanks anyway. Are you aware of any libraries similar to Real_Arrays, but who allocated memory internally using heap? This is the natural way to do such things.

    The most natural way to work with an array of FPT numbers is
    for the programmer to declare an array indexed by some index type.
    Done. If GNAT gets in the way there, it might be worth a note sent to
    its maintainers. Whenever a programmer is tasked with considering
    memory allocation, then depending on one's propensity towards working
    on memory allocation it is inconvenient and distracting.
    Math programs don't make you do this, I think.

    Also, std::vector an its relatives shield the programmer from
    the absurdly clever, yet unreadable memory allocation that
    needs to be stuffed behind the scenes.
    More importantly, though, C++ introduced std::move semantics
    after a few decades of its existence, to address copying when
    using chains of +. It might be interesting to see Ada's in-situ
    construction of return values in comparison.


    Similarly to the Containers.Vector. But Vector has such an awful syntax. There should be something like an indexer [i] similarly to the C++ std::vector to make things simpler

    .at() does some of what Ada does.
    Is

    v.at(k) = 4;

    less awful than

    v(k) := 4;

    ?

    Another thing: Mathematical notation has ellipsis, thus

    A + B + ... + Y + Z;

    Most general purpose languages don't have ellipsis for this kind
    of expression. However, even mathematical formulas use what
    programmers can usually achieve, too. The usual

    \sum_k A_k.

    No "+" at all, and an array of vectors, not single ones.
    Going further, some like to write

    reduce("+", A);

    In Ada, you could have a generic function for this, or use a
    function pointer.

    The .all thing vanishes automatically whenever you want to refer
    to a particular component of the pointed-at object, as opposed
    to all of them. So, A.all(K) is the same as A(K).
    Likewise, .all can be dropped if want to invoke the pointed-at
    subprogram if it has parameters.

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From J-P. Rosen@21:1/5 to All on Mon Jan 23 09:51:55 2023
    Le 23/01/2023 à 08:50, Egil H H a écrit :
    On Monday, January 23, 2023 at 12:53:15 AM UTC+1, joak...@kth.se wrote:
    Something came up and I had to send my previous reply/e-mail as is. I wanted to find the video where Jean Pierre Rosen talks about how memory is handled in the Ada language from FOSDEM perhaps 2018-2019. Unfortunately I have been unable to find it.

    It was in 2016:

    https://archive.fosdem.org/2016/schedule/event/ada_memory/

    Thanks Egil, you were faster than me...
    I have also a full tutorial at several Ada-Europe conferences. No video,
    but I can send the slides to those interested.

    --
    J-P. Rosen
    Adalog
    2 rue du Docteur Lombard, 92441 Issy-les-Moulineaux CEDEX
    https://www.adalog.fr https://www.adacontrol.fr

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Dmitry A. Kazakov@21:1/5 to Jim Paloander on Mon Jan 23 09:28:46 2023
    On 2023-01-22 23:36, Jim Paloander wrote:

    I was not sure whether or not it can be avoided with Implicit_Dereference,

    type Accessor (Data: not null access Element) is limited private
    with Implicit_Dereference => Data;

    If you create a new wrapper type, anyway, then it is easier to define operations directly on that new type.

    Otherwise what you described for operator+ one has to do for every operator overloaded inside Real_Arrays package.

    You should not use the standard library anyway. It is not intended for
    large problems, which require specific approaches and methods, like
    sparse matrices, concurrent processing and so on.

    The optimal way to work large linear algebra problem is what you describe because unfortunately ADA does not allow what Fortran does since 30 years ago or more.

    I am not sure what you mean. It is quite possible to design a wrapper
    datatype allocating vectors/matrices in the pool. E.g. Ada's
    Unbounded_String is such a thing. Real_Arrays were not designed this way because see above.

    But in C++ you can reproduce the same functionality as Fortran using Expression Templates and Template Metaprogramming.

    Nothing prevents you from wrapping Real_Array in a generic way:

    generic
    with package Real_Arrays is new Numerics.Generic_Real_Arrays (<>);
    package Generic_Pool_Real_Arrays is
    ...
    end Generic_Pool_Real_Arrays;

    Perhaps ADA should allow something like that. Because for maintainability reasons the best would be to write the mathematical expressions as close as possible to the mathematical formulas.

    There is no problem with that as you can define operations on pointers.

    --
    Regards,
    Dmitry A. Kazakov
    http://www.dmitry-kazakov.de

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to All on Mon Jan 23 17:04:39 2023
    Dmitry,
    you mentioned in a previous email

    You should not use the standard library anyway. It is not intended for
    large problems, which require specific approaches and methods, like
    sparse matrices, concurrent processing and so on.
    The optimal way to work large linear algebra problem is what you describe because unfortunately ADA does not allow what Fortran does since 30 years ago or more.
    I am not sure what you mean. It is quite possible to design a wrapper

    that the optimal way is to work with specific approaches and methods. Nevertheless, at least for vector operations, and/or dense matrices operations, ADA should provide native support of vector-math or Vector.Numerics, or AdvancedNumerics or Numerics.
    LinearAlgebra and/or Numerics.SparseLinearAlgebra. Honestly with all the mission critical applications ADA is used for, I would expect something like that to enable people working in scientific computing and in general wherever Linear Algebra kernels are
    essential tool of the core algorithms development.

    So it seems that as you correctly implied even dense linear vector math are not enough. But it is at least what I hoped for. Sparse matrix times vector could be stored in a temporary vector. I would not expect that an expression such as

    y := alpha*a + beta*(z-w*gamma) + alpha*beta*SparseMat*x;

    is computed in a way similar to Fortran, which is

    for I in 1 .. N loop
    y[I] := alpha*a[I] + beta*(z[I] - w[I]*gamma) + alpha*beta*Sum_j ( SparseMat[I,j]*x[j])
    end

    where Sum_j is implemented in another loop. I would instead compute v := SparseMat*x, and then replace SparseMat*x with v in the above expression. So I would be happy if ADA would give me a library that allows Vector Math with performance similar to
    Fortran without temporaries for both Complex and Real numbers. If they could provide that also for Sparse Linear Algebra even better even for CSR format, because there are plenty Sparse Matrix formats for Sparse linear Algebra.

    Thanks.

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From J-P. Rosen@21:1/5 to All on Tue Jan 24 11:42:11 2023
    Le 24/01/2023 à 02:04, Jim Paloander a écrit :
    ADA should provide native support of vector-math or Vector.Numerics, or AdvancedNumerics or Numerics.LinearAlgebra and/or Numerics.SparseLinearAlgebra. Honestly with all the mission critical applications ADA is used for, I would expect something like
    that to enable people working in scientific computing and in general wherever Linear Algebra kernels are essential tool of the core algorithms development.

    Did you consider the facilities from Annex G? Moreover, there are
    bindings to BLAS, LAPACK, and more. Try "Ada mathematical libraries" on Google...

    --
    J-P. Rosen
    Adalog
    2 rue du Docteur Lombard, 92441 Issy-les-Moulineaux CEDEX
    https://www.adalog.fr https://www.adacontrol.fr

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Gautier write-only address@21:1/5 to All on Tue Jan 24 11:47:13 2023
    There are plenty of legitimate complaints and good ideas in this thread.

    - need for a package for large matrices (issue with stack allocation)
    - should support approximations (at least, floating-point) real and complex numbers (and perhaps others)
    - should support various matrix storages: sparse storages, band storages, ...

    Fortunately Ada (not "ADA" = "American Dentist Association" and some others...) is powerful enough to have such components written in Ada itself (with GNAT, it is done exactly like this).
    So _we_ don't need to wait that _they_ do anything about it :-)

    So, if I summarize the ideas discussed and combine with stuff grabbed from a toolbox of mine ( https://github.com/zertovitch/mathpaqs
    ), I obtain the following specification.
    The choice between different kind of matrix storages would be trivial.
    Comments are welcome!

    ---8<--------8<--------8<--------8<--------8<--------8<-----
    -- Draft of a specification for an universal matrix package.

    -- The elements of the vectors and matrices can be of any algebraic ring
    -- and their implementation, eventually approximate:
    --
    -- - integers (stored as Integer, Integer_n or Big_Integer)
    -- - modular integers
    -- - real numbers (approximated as floating-point, fixed-point or Big_Real) -- - complex numbers
    -- - rational numbers
    -- - polynomials (of real, complex, rational, ...)
    -- - other matrices
    -- - ...

    with Ada.Containers.Vectors;

    generic

    type Ring_Element is private; -- Element of an algebraic ring

    zero, one : Ring_Element; -- 0 and 1 elements

    with function "-" (a : Ring_Element) return Ring_Element; -- Unary operators

    with function "+" (a, b : Ring_Element) return Ring_Element; -- Binary operators
    with function "-" (a, b : Ring_Element) return Ring_Element;
    with function "*" (a, b : Ring_Element) return Ring_Element;

    package Universal_Matrices is

    package UM_Vectors is new Ada.Containers.Vectors (Positive, Ring_Element);

    subtype Vector is UM_Vectors.Vector;

    -------------------------
    -- Vector operations --
    -------------------------

    function "*" (v : Vector; factor : Ring_Element) return Vector;
    function "*" (factor : Ring_Element; v : Vector) return Vector;

    -- v := factor * v :
    procedure Scale_Left (v : in out Vector; factor : Ring_Element);

    function "-" (v : Vector) return Vector;

    function "+" (v, w : Vector) return Vector;
    -- v := v + w :
    procedure Add (v : in out Vector; w : Vector);
    -- v := v + factor * w :
    procedure Add_Left_Scaled (v : in out Vector; factor : Ring_Element; w : Vector);

    function "-" (v, w : Vector) return Vector;
    -- v := v - w :
    procedure Subtract (v : in out Vector; w : Vector);

    function "*" (v, w : Vector) return Ring_Element;

    -- Euclidean norm and distance:
    function Square_L2_Norm (v : Vector) return Ring_Element;
    function Square_L2_Distance (v, w : Vector) return Ring_Element;

    -------------------------------------------------------------------
    -- Root matrix type. --
    -- Possible derivations: dense, sparse, band storage matrices. --
    -------------------------------------------------------------------

    type Matrix is interface;

    -------------------------
    -- Matrix operations --
    -------------------------

    function Transpose (A : Matrix) return Matrix is abstract;
    function Identity (order : Positive) return Matrix is abstract;
    function "*" (factor : Ring_Element; A : Matrix) return Matrix is abstract;
    function "*" (A : Matrix; factor : Ring_Element) return Matrix is abstract;
    function "*" (A, B : Matrix) return Matrix is abstract;
    function "+" (A, B : Matrix) return Matrix is abstract;
    function "-" (A, B : Matrix) return Matrix is abstract;

    -- Matrix-Vector operations

    function "*" (A : Matrix; x : Vector) return Vector is abstract;

    end Universal_Matrices;

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Gautier write-only address@21:1/5 to All on Tue Jan 24 15:02:21 2023
    Current version of the draft is located here: https://github.com/zertovitch/mathpaqs/blob/master/lin_alg/universal_matrices.ads

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to All on Wed Jan 25 01:50:56 2023
    Dear Gautier,

    thank you for the draft. Indeed something like that would be great possibly without polynomials or rational numbers etc. You just need doubles and/or floats and complex numbers. Sparse matrices are another animal, and once adopted they are another
    dimension of complexity. This is because apart from the triplet format which is easily expandable, one may also need to convert the triplet to CSR and there multiple entries with same row and column index need first to be summed so we have a single entry
    (i, j, value) for each (i, j) before converting to CSR. So I mean, one would like to build a library then and I am not so sure how easy would be to hardcode this library so that operations on matrices/vectors are as efficient as in Fortran. For Dense
    matrices it may still be possible to support mathematical expressions without introducing temporary objects. But if you also bring sparse into the game then?
    (alpha,beta,gamma scalars, while x,y,z,w vectors, A, B_sparse matrices)

    does this have a chance to be executed without temporaries?
    y := alpha*x + beta*A*x + gamma*B_sparse*(z-w)

    this has a chance to be executed without temporaries Fortran already does it.
    y := alpha*x + beta*A*x

    I would create a SparseMatrix library and break the statement above in 2
    y1 := B_sparse*(z-w)
    y := alpha*x + beta*A*x + gamma*y1;

    I could survive with that. The problem dear Gautier is the temporary objects in such libraries. Very correctly J.P.Rosen suggested BLAS. The problem is that for simple vector operations (BLAS1) there is no speedup at all and the operations of BLAS1
    cannot be executed for the whole expression, just in an `axpy (a*x+y)` manner one by one. This would require transversal of the arrays many times in a long expression involving vectors as (y := x1 + x2*alpha + beta*(x3-gamma*x4) + alpha*beta*(x1-x2)).
    For BLAS2 operations (matrix times vector) you may see a performance gain using BLAS but for large matrices, very large. Where you definitely see such a performance gain is for BLAS3 and there we definitely need something like that.

    This discussion started about techniques to avoid generation of temporaries. Even when you manipulate complex numbers you introduce temporaries. I need to come back to you guys with some number to understand the problem better. I will be back with these
    numbers. But until then, please let me know, do you feel that there is a need of a performance uplifting and mathematics enabling ADA library which is not Real_Arrays or Containers, but a real Math library allowing Real and Complex math manipulation and
    integrating possibly some Lapack routines as Real_Arrays do at the moment?

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to All on Wed Jan 25 01:52:57 2023
    ADA should provide native support of vector-math or Vector.Numerics, or AdvancedNumerics or Numerics.LinearAlgebra and/or Numerics.SparseLinearAlgebra. Honestly with all the mission critical applications ADA is used for, I would expect something like
    that to enable people working in scientific computing and in general wherever Linear Algebra kernels are essential tool of the core algorithms development.
    Did you consider the facilities from Annex G? Moreover, there are
    bindings to BLAS, LAPACK, and more. Try "Ada mathematical libraries" on Google...

    Thank you, very nice and positive suggestion. Please read my answer below to Gautier. Thank you again for the nice feedback.By the way, is there any other forum/discussion for ADA in github or elsewhere where we can also post code with colours etc?

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From J-P. Rosen@21:1/5 to All on Wed Jan 25 13:21:15 2023
    Le 25/01/2023 à 10:52, Jim Paloander a écrit :

    ADA should provide native support of vector-math or Vector.Numerics, or AdvancedNumerics or Numerics.LinearAlgebra and/or Numerics.SparseLinearAlgebra. Honestly with all the mission critical applications ADA is used for, I would expect something like
    that to enable people working in scientific computing and in general wherever Linear Algebra kernels are essential tool of the core algorithms development.
    Did you consider the facilities from Annex G? Moreover, there are
    bindings to BLAS, LAPACK, and more. Try "Ada mathematical libraries" on
    Google...

    Thank you, very nice and positive suggestion. Please read my answer below to Gautier. Thank you again for the nice feedback.By the way, is there any other forum/discussion for ADA in github or elsewhere where we can also post code with colours etc?

    Gitter: https://gitter.im/ada-lang
    LinkedIn: https://www.linkedin.com/groups/114211
    Reddit: https://www.reddit.com/r/ada
    StackOverflow: https://www.stackoverflow.com/questions/tagged/ada
    Twitter: https://twitter.com/search?q=adaprogramming

    --
    J-P. Rosen
    Adalog
    2 rue du Docteur Lombard, 92441 Issy-les-Moulineaux CEDEX
    https://www.adalog.fr https://www.adacontrol.fr

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Gautier write-only address@21:1/5 to All on Wed Jan 25 14:41:59 2023
    ...and: https://forum.ada-lang.io/

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to All on Thu Jan 26 11:08:46 2023
    ...and: https://forum.ada-lang.io/

    Which one do you think is better for general discussions on both ADA and ADA libraries?
    Better in the sense that it is preferred by experienced open-minded people with engineering background or
    people who may not necessarily hold an engineering degree or PhDs and Masters but who
    have worked for many years on real world engineering problems ADA is supposed to assist.

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jerry@21:1/5 to ...@gmail.com on Thu Jan 26 12:39:29 2023
    On Sunday, January 22, 2023 at 2:34:19 PM UTC-7, ...@gmail.com wrote:
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression

    Leo's answer is responsive to the OP but seems to have gotten buried in an amazingly long discussion which I haven't read. This answer appeared here years ago (was it you, Leo?) and I have used it ever since. If there needs to be a long discussion maybe
    it could be about why this workaround is necessary. I will pull out Leo's answer from his post and put my version here:

    type Real_Vector_Access is access Real_Vector;
    then
    x_Ptr : Real_Vector_Access := new Real_Vector(0 .. N - 1);
    x : Real_Vector renames x_Ptr.all;

    That's all. All the overloaded operators work as though the vector was declared from the stack. If you have code which formerly used stack allocation (of x, in this example), it will work without modification with this trick.

    Jerry

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jim Paloander@21:1/5 to All on Thu Jan 26 13:52:58 2023
    Dear ADA lovers,
    with stack allocation of Real_Vector ( 1 .. N ) when N >= 100,000 I get STACK_OVERFLOW ERROR while trying to check how fast operator overloading is working for an expression
    Leo's answer is responsive to the OP but seems to have gotten buried in an amazingly long discussion which I haven't read. This answer appeared here years ago (was it you, Leo?) and I have used it ever since. If there needs to be a long discussion
    maybe it could be about why this workaround is necessary. I will pull out Leo's answer from his post and put my version here:

    type Real_Vector_Access is access Real_Vector;
    then
    x_Ptr : Real_Vector_Access := new Real_Vector(0 .. N - 1);
    x : Real_Vector renames x_Ptr.all;

    That's all. All the overloaded operators work as though the vector was declared from the stack. If you have code which formerly used stack allocation (of x, in this example), it will work without modification with this trick.

    Jerry

    Thank you Jerry, that's true, this trick does the job indeed. Thanks. But to be honest among us, it seems to me that the Real_Arrays library has been implemented probably (note the probably) by people not targeting real world engineering / industrial
    applications. Examples are Finite Volume solutions of semiconductor device equations, or weather prediction problems. You cannot expect to solve these on the stack and of course ADA containers are impractical for such computations. Possibly a library as
    discussed above that provides Real_Arrays functionality but memory allocation is done internally on the heap or even on the stack if the user chooses to do so to take advantage of fast stack access. But then the mechanism for changing the stack
    allocation size transparently and easily should also be provided as simply as setting the memory pool in an access type.

    Do you agree?

    In the meantime, I am wondering whether or not a hybrid approach (stack/heap) would work provided such a mechanism to set the stack size is easy accessible as simply as a memory pool. The functionality would be during operator overloading we create fast
    stack allocated arrays to transfer the intermediate results stored in temporary objects. So all the temporaries described above would reside in the stack if this would translate to performance improvement until a smarter approach is introduced similar to
    C++ expression templates.

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From Jerry@21:1/5 to All on Thu Feb 2 13:59:42 2023
    Thank you Jerry, that's true, this trick does the job indeed. Thanks. But to be honest among us, it seems to me that the Real_Arrays library has been implemented probably (note the probably) by people not targeting real world engineering / industrial
    applications. Examples are Finite Volume solutions of semiconductor device equations, or weather prediction problems. You cannot expect to solve these on the stack and of course ADA containers are impractical for such computations. Possibly a library as
    discussed above that provides Real_Arrays functionality but memory allocation is done internally on the heap or even on the stack if the user chooses to do so to take advantage of fast stack access. But then the mechanism for changing the stack
    allocation size transparently and easily should also be provided as simply as setting the memory pool in an access type.

    Do you agree?

    Yes, with the minor caveat that I have no idea what I'm talking about.

    When I ran into this limitation a few years ago, I learned how to increase the stack size before running my program. This helped but only a little, as large arrays still were not accommodated. (Define "large" as anything too big for stack expanded to OS
    limits.) I have never compared speed with stack versus heap allocation because, for me, it would serve no purpose. I can't just use smaller arrays in my work because large ones are too slow, if you follow my drift here.

    The applications you mention are great but the list must be nearly endless. My applications run towards audio signals and radar (simulation) images.

    Jerry

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)