Compile error when composing linear algebra functions

[redhatturtle]

Summary of Problem

Trying to compile the following program in the most recent builds yields the error

linalg-test.chpl:12: error: non-lvalue actual is passed to a 'ref' formal of inv()

which I guess is related to the recent ref-maybe-const deprecation.

There's obviously a simple workaround of creating a temporary constant, which I'm currently using.

I've had a similar experience composing an internal function with reshape that I solved with a temporary variable. But if there's a different solution by changing the function signature it seems like it would be more appropriate.

Steps to Reproduce

Source Code:

use LinearAlgebra;

config const n = 3;

var A, B, AtB_inv : [1..n, 1..n] real;

// This works
var AtB = dot(A,B);
AtB_inv = inv(AtB);

// This doesn't, but it did work for a long time including in v1.31
AtB_inv = inv(dot(A,B));

write(AtB_inv);

Compile command: chpl -lcblas -llapacke linalg-test.chpl;

Execution command: ./linalg-test -nl 1

Configuration Information

• Output of chpl --version:

  chpl version 1.32.0 pre-release (ca8b6f73db)
  built with LLVM version 15.0.7
  available LLVM targets: amdgcn, r600, nvptx64, nvptx, aarch64_32, aarch64_be, aarch64, arm64_32, arm64, x86-64, x86
  Copyright 2020-2023 Hewlett Packard Enterprise Development LP
  Copyright 2004-2019 Cray Inc.
  (See LICENSE file for more details)

• Output of $CHPL_HOME/util/printchplenv --anonymize:

  CHPL_TARGET_PLATFORM: linux64
  CHPL_TARGET_COMPILER: llvm
  CHPL_TARGET_ARCH: x86_64
  CHPL_TARGET_CPU: native +
  CHPL_LOCALE_MODEL: flat
  CHPL_COMM: gasnet +
  CHPL_COMM_SUBSTRATE: smp +
  CHPL_GASNET_SEGMENT: fast
  CHPL_TASKS: qthreads
  CHPL_LAUNCHER: smp
  CHPL_TIMERS: generic
  CHPL_UNWIND: bundled +
  CHPL_MEM: jemalloc
  CHPL_ATOMICS: cstdlib
  CHPL_NETWORK_ATOMICS: none
  CHPL_GMP: bundled
  CHPL_HWLOC: bundled
  CHPL_RE2: bundled +
  CHPL_LLVM: bundled +
  CHPL_AUX_FILESYS: none

• Back-end compiler and version, e.g. gcc --version:

  gcc (SUSE Linux) 13.2.1 20230912 [revision b96e66fd4ef3e36983969fb8cdd1956f551a074b]
  Copyright (C) 2023 Free Software Foundation, Inc.
  This is free software; see the source for copying conditions.  There is NO
  warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

[mppf]

IMO this is an API design issue in the linear algebra module.

Looking at proc inv specifically:

proc inv(ref A: [?Adom] ?eltType, overwrite=false) where usingLAPACK { ... }

Here it seems that the LAPACK calls modify the array in-place. At some point in history, it might have been the case that if overwrite=true, this routine would not return a copy of A. However, today, it copies A either way (with overwrite=false, it does it with var A_clone = A; and with overwrite=true it does it with return A;. So, I don't think the design of having overwrite as an optional argument here is reasonable. It doesn't offer much benefit (in terms of reducing copying in some cases, at least) and it means that the function needs to accept A with ref intent (because it might mutate it).

Instead, I'd recommend that we change proc inv (and any LinearAlgebra routine in a similar situation) into two routines with different names:

• A routine that never modifies the array passed in, always copies it and returns the new array (the result of inverting). E.g. proc inv(const ref A:[]): A.type.
• Another routine that never returns anything and mutates the array in place to compute the inverse. E.g. proc invInPlace(ref A:[]).

[damianmoz]

I agree with Michael that there is an issue with the design of the API. I would go so far as as to say there is an issue with its use.

I am assuming that you are intending these routines for use in real applications, not just toy (or purely pedagogic) applications.

Most of the time, you won't need a copy. In those situations where you really do need a copy of the original matrix as well as the mutated one, or probably more precisely, a factored or transformed one,, you really need to make those copies blatantly obvious in the code, not hidden under the hood so you never know the real overhead!

As a user of LINPACK and now LAPACK for a long time, and having rewritten them in C and C++ for real world use, the default name should refer to modify in place. Having the reverse convention encourages very (very very) bad programming practices.

I regularly get given MATLAB code that gobbles up memory like it is going out of style because there are copies of arrays all over the place. It makes such codes impractical for production work. The bane of my life.

Memory is expensive (from numerous perspectives even if its dollar cost has dropped). Don't waste it on lazy programming practices.

That advice is worth a lot more than my usual 2c.

[vasslitvinov]

@redhatturtle - until we have a better API for inv, your working variant is the same what used to happen in v1.31 under the covers for AtB_inv = inv(dot(A,B));.

I would suggest a minor tweak to your solution: move the declaration of AtB_inv to where it gets the output of inv(). This will ensure that AtB_inv takes over the outcome of inv():

var A, B : [1..n, 1..n] real;

var AtB = dot(A,B);
var AtB_inv = inv(AtB);

Cf. here the array produced by inv() is copied into the array AtB_inv then discarded:

var A, B, AtB_inv : [1..n, 1..n] real;
writeln(AtB_inv);

var AtB = dot(A,B);
AtB_inv = inv(AtB);