Test Array in mat_rclist_*

inst/@sym/private/mat_rclist_access.m

@@ -1,4 +1,5 @@
 %% Copyright (C) 2014, 2016, 2019, 2022 Colin B. Macdonald
+%% Copyright (C) 2022 Alex Vong
 %%
 %% This file is part of OctSymPy.
 %%
@@ -32,15 +33,12 @@
     error('this routine is for a list of rows and cols');
   end
 
-  cmd = { '(A, rr, cc) = _ins'
-          'if A is None or not A.is_Matrix:'
-          '    A = sp.Matrix([A])'
-          'n = len(rr)'
-          'M = [[0] for i in range(n)]'
-          'for i in range(0, n):'
-          '    M[i][0] = A[rr[i],cc[i]]'
-          'M = sp.Matrix(M)'
-          'return M,' };
+  cmd = {'dbg_no_array = True'
+         '(A, rr, cc) = _ins'
+         'AA = A.tolist() if isinstance(A, (MatrixBase, NDimArray)) else [[A]]'
+         'MM = [[AA[i][j]] for i, j in zip(rr, cc)]'
+         'M = make_matrix_or_array(MM)'
+         'return M,'};
 
   rr = num2cell(int32(r-1));
   cc = num2cell(int32(c-1));

inst/@sym/private/mat_rclist_asgn.m

@@ -1,6 +1,7 @@
 %% Copyright (C) 2014, 2016-2017, 2019, 2022 Colin B. Macdonald
 %% Copyright (C) 2020 Mike Miller
 %% Copyright (C) 2020 Fernando Alvarruiz
+%% Copyright (C) 2022 Alex Vong
 %%
 %% This file is part of OctSymPy.
 %%
@@ -59,40 +60,38 @@
   %    AA[0, 0] = A
   % Also usefil: .copyin_matrix
 
-  cmd = { '(A, r, c, B) = _ins'
-          '# B linear access fix, transpose for sympy row-based'
-          'if B is None or not B.is_Matrix:'
-          '    B = sp.Matrix([[B]])'
-          'BT = B.T'
-          '# make a resized copy of A, and copy existing stuff in'
-          'if isinstance(A, list):'
-          '    assert len(A) == 0, "unexpectedly non-empty list: report bug!"'
-          '    n = max(max(r) + 1, 1)'
-          '    m = max(max(c) + 1, 1)'
-          '    AA = [[0]*m for i in range(n)]'
-          'elif A is None or not isinstance(A, MatrixBase):'
-          '    # we have non-matrix, put in top-left'
-          '    n = max(max(r) + 1, 1)'
-          '    m = max(max(c) + 1, 1)'
-          '    AA = [[0]*m for i in range(n)]'
-          '    AA[0][0] = A'
-          'else:'
-          '    # build bigger matrix'
-          '    n = max(max(r) + 1, A.rows)'
-          '    m = max(max(c) + 1, A.cols)'
-          '    AA = [[0]*m for i in range(n)]'
-          '    # copy current matrix in'
-          '    for i in range(A.rows):'
-          '        for j in range(A.cols):'
-          '            AA[i][j] = A[i, j]'
-          '# now insert the new bits from B'
-          'for i, (r, c) in enumerate(zip(r, c)):'
-          '    AA[r][c] = BT[i]'
-          'return sp.Matrix(AA),' };
+  cmd = {'dbg_no_array = True'
+         '(A, rr, cc, b) = _ins'
+         'assert A == [] or not isinstance(A, list), "unexpectedly non-empty list: report bug!"'
+         'if A == []:'
+         '    AA = []'
+         '    (nrows_A, ncols_A) = (0, 0)'
+         'elif isinstance(A, (MatrixBase, NDimArray)):'
+         '    AA = A.tolist()'
+         '    (nrows_A, ncols_A) = A.shape'
+         'else:'
+         '    AA = [[A]]'
+         '    (nrows_A, ncols_A) = (1, 1)'
+         'bb = b.tolist() if isinstance(b, (MatrixBase, NDimArray)) else [[b]]'
+         'bb_elts = itertools.chain.from_iterable(bb)'
+         'entries = dict(zip(zip(rr, cc), bb_elts))'
+         'def entry(i, j):'
+         '    if (i, j) in entries:'
+         '        return entries[i, j]'
+         '    elif i < nrows_A and j < ncols_A:'
+         '        return AA[i][j]'
+         '    else:'
+         '        return 0'
+         'n = max(max(rr) + 1, nrows_A)'
+         'm = max(max(cc) + 1, ncols_A)'
+         'MM = [[entry(i, j) for j in range(m)] for i in range(n)]'
+         'M = make_matrix_or_array(MM)'
+         'return M,'};
 
   rr = num2cell(int32(r-1));
   cc = num2cell(int32(c-1));
-  z = pycall_sympy__ (cmd, A, rr, cc, B);
+  b = vec (B); # B is accessed with linear indexing, as a column vector
+  z = pycall_sympy__ (cmd, A, rr, cc, b);
 
   % a simpler earlier version, but only for scalar r,c
   %cmd = { '(A, r, c, b) = _ins'

inst/private/python_header.py

@@ -235,6 +235,30 @@ def octoutput(x, et):
 except:
     echo_exception_stdout("in python_header defining fcns block 4")
     raise
+
+
+try:
+    def make_matrix_or_array(ls_of_ls, dbg_no_array=True):
+        # should be kept in sync with the same function
+        # defined in inst/private/python_ipc_native.m
+        # FIXME: dbg_no_array is currently used for debugging,
+        # remove it once sympy drops non-Expr supp in Matrix
+        """
+        Given a list of list of syms LS_OF_LS
+        If all elements of LS_OF_LS are Expr,
+        construct the corresponding Matrix.
+        Otherwise, construct the corresponding 2D array.
+        """
+        elts = itertools.chain.from_iterable(ls_of_ls)
+        if (dbg_no_array
+            or all(isinstance(elt, Expr) for elt in elts)):
+            return Matrix(ls_of_ls)
+        else:
+            dbout(f"make_matrix_or_array: making 2D Array...")
+            return Array(ls_of_ls)
+except:
+    echo_exception_stdout("in python_header defining fcns block 5")
+    raise
 # end of python header, now couple blank lines
 
 

inst/private/python_ipc_native.m

@@ -111,6 +111,24 @@
                     '    # should be kept in sync with the same function'
                     '    # defined in inst/private/python_header.py'
                     '    sys.stderr.write("pydebug: " + str(l) + "\n")'
+                    'def make_matrix_or_array(ls_of_ls, dbg_no_array=True):'
+                    '    # should be kept in sync with the same function'
+                    '    # defined in inst/private/python_header.py'
+                    '    # FIXME: dbg_no_array is currently used for debugging,'
+                    '    # remove it once sympy drops non-Expr supp in Matrix'
+                    '    """'
+                    '    Given a list of list of syms LS_OF_LS'
+                    '    If all elements of LS_OF_LS are Expr,'
+                    '    construct the corresponding Matrix.'
+                    '    Otherwise, construct the corresponding 2D array.'
+                    '    """'
+                    '    elts = itertools.chain.from_iterable(ls_of_ls)'
+                    '    if (dbg_no_array'
+                    '        or all(isinstance(elt, Expr) for elt in elts)):'
+                    '        return Matrix(ls_of_ls)'
+                    '    else:'
+                    '        dbout(f"make_matrix_or_array: making 2D Array...")'
+                    '        return Array(ls_of_ls)'
                   }, newl))
     have_headers = true;
   end