stfc / PSyclone

Domain-specific compiler and code transformation system for Finite Difference/Volume/Element Earth-system models in Fortran
BSD 3-Clause "New" or "Revised" License
104 stars 27 forks source link

[InlineTrans] An error occured during the call to Reference2ArrayRangeTrans apply while inlining a subroutine. #2726

Open hbrunie opened 1 day ago

hbrunie commented 1 day ago

Here is the code snippet to highlight the bug

def test_psyclone_limits_on_ref2arrayrangetrans():
    from psyclone.psyir import nodes as psyir_nodes
    from psyclone.psyir.frontend.fortran import FortranReader
    from psyclone.psyir.transformations import Reference2ArrayRangeTrans

    ref2arraytrans = Reference2ArrayRangeTrans()
    code = """subroutine sub_caller()
      integer, dimension(3), intent(in) :: logc
      real, dimension(logc(3):,:,logc(1):,logc(2):), intent(out) :: arg1_caller
      call sub_callee(arg1_caller, logc)
    end subroutine sub_caller

    subroutine sub_callee(arg1, logc)
      integer, dimension(3), intent(in) :: logc
      real, dimension(logc(3):,:,logc(1):,logc(2):), intent(inout) :: arg1
      arg1(1,1,1,1) = logc(1)
    end subroutine sub_callee
    """
    reader = FortranReader()
    psyir_tree: psyir_nodes.Node = reader.psyir_from_source(code)
    node: psyir_nodes.Call = psyir_tree.walk(psyir_nodes.Call)[0]
    assert not node.routine is None
    routine: psyir_nodes.Reference = node.routine
    assert routine.name == "sub_callee"
    for child in node.arguments:
        ref2arraytrans.apply(child)

Here is the message of error I got:

../venv/lib/python3.10/site-packages/psyclone/psyir/transformations/reference2arrayrange_trans.py:189: in apply
    Reference2ArrayRangeTrans._get_array_bound(symbol, idx)
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 

symbol = <psyclone.psyir.symbols.datasymbol.DataSymbol object at 0x7d75a3ca7970>, index = 0

    @staticmethod
    def _get_array_bound(symbol, index):
        '''A utility function that returns the appropriate loop bounds (lower,
        upper and step) for an array dimension.  If the array
        dimension is declared with known bounds (an integer or a
        symbol) then these bound values are used. If the size is
        unknown (a deferred or attribute type) then the LBOUND and
        UBOUND PSyIR nodes are used.

        :param symbol: the symbol that we are interested in.
        :type symbol: :py:class:`psyir.symbols.DataSymbol`
        :param int index: the (array) reference index that we are \
            interested in.

        :returns: the loop bounds for this array index.
        :rtype: Tuple(:py:class:`psyclone.psyir.nodes.Literal`, \
                      :py:class:`psyclone.psyir.nodes.Literal`, \
                      :py:class:`psyclone.psyir.nodes.Literal`) or \
                Tuple(:py:class:`psyclone.psyir.nodes.BinaryOperation`, \
                      :py:class:`psyclone.psyir.nodes.BinaryOperation`, \
                      :py:class:`psyclone.psyir.nodes.Literal`)

        '''
        # Look for explicit bounds in the array declaration.
        my_dim = symbol.shape[index]
        if isinstance(my_dim, ArrayType.ArrayBounds):
            lower_bound = my_dim.lower.copy()
>           upper_bound = my_dim.upper.copy()
E           AttributeError: 'Extent' object has no attribute 'copy'

../venv/lib/python3.10/site-packages/psyclone/psyir/transformations/reference2arrayrange_trans.py:112: AttributeError

Here is the code snippet I wrote in reference2arrayrange_trans.py line 108 to avoid the bug:

        # Look for explicit bounds in the array declaration.
        my_dim = symbol.shape[index]
        if isinstance(my_dim, ArrayType.ArrayBounds):
            if isinstance(my_dim.lower, Reference):
                lower_bound = my_dim.lower.copy()
            else:
                lower_bound = IntrinsicCall.create(
                    IntrinsicCall.Intrinsic.LBOUND,
                    [Reference(symbol), ("dim", Literal(str(index+1), INTEGER_TYPE))])
            if isinstance(my_dim.upper, Reference):
                upper_bound = my_dim.upper.copy()
            else:
                upper_bound = IntrinsicCall.create(
                    IntrinsicCall.Intrinsic.UBOUND,
                    [Reference(symbol), ("dim", Literal(str(index+1), INTEGER_TYPE))])
            step = Literal("1", INTEGER_TYPE)
            return (lower_bound, upper_bound, step)

        # No explicit array bound information could be found so use the
        # LBOUND and UBOUND intrinsics.
        lower_bound = IntrinsicCall.create(
            IntrinsicCall.Intrinsic.LBOUND,
            [Reference(symbol), ("dim", Literal(str(index+1), INTEGER_TYPE))])
        upper_bound = IntrinsicCall.create(
            IntrinsicCall.Intrinsic.UBOUND,
            [Reference(symbol), ("dim", Literal(str(index+1), INTEGER_TYPE))])
        step = Literal("1", INTEGER_TYPE)
        return (lower_bound, upper_bound, step)

The corresponding git patch I made:

From 368e37a78071b926f38ae51344824728a8ccf4cd Mon Sep 17 00:00:00 2001
From: Hugo Brunie <hugo.brunie@univ-grenoble-alpes.fr>
Date: Tue, 1 Oct 2024 16:03:13 +0200
Subject: [PATCH] [Ref2ArrayRangeTrans] fix _get_array_bound

AttributeError: 'Extent' object has no attribute 'copy'
real, dimension(logc(3):,:,logc(1):,logc(2):), intent(inout) :: hy_starstate.
An array is used as of lower bound.
Exlicit the bug.
---
 .../transformations/reference2arrayrange_trans.py  | 14 ++++++++++++--
 1 file changed, 12 insertions(+), 2 deletions(-)

diff --git a/src/psyclone/psyir/transformations/reference2arrayrange_trans.py b/src/psyclone/psyir/transformations/reference2arrayrange_trans.py
index 4a32c93f9..206ddb1b5 100644
--- a/src/psyclone/psyir/transformations/reference2arrayrange_trans.py
+++ b/src/psyclone/psyir/transformations/reference2arrayrange_trans.py
@@ -108,8 +108,18 @@ class Reference2ArrayRangeTrans(Transformation):
         # Look for explicit bounds in the array declaration.
         my_dim = symbol.shape[index]
         if isinstance(my_dim, ArrayType.ArrayBounds):
-            lower_bound = my_dim.lower.copy()
-            upper_bound = my_dim.upper.copy()
+            if isinstance(my_dim.lower, Reference):
+                lower_bound = my_dim.lower.copy()
+            else:
+                lower_bound = IntrinsicCall.create(
+                    IntrinsicCall.Intrinsic.LBOUND,
+                    [Reference(symbol), ("dim", Literal(str(index+1), INTEGER_TYPE))])
+            if isinstance(my_dim.upper, Reference):
+                upper_bound = my_dim.upper.copy()
+            else:
+                upper_bound = IntrinsicCall.create(
+                    IntrinsicCall.Intrinsic.UBOUND,
+                    [Reference(symbol), ("dim", Literal(str(index+1), INTEGER_TYPE))])
             step = Literal("1", INTEGER_TYPE)
             return (lower_bound, upper_bound, step)

-- 
2.34.1