Repository: stfc/PSyclone (python)

Issue:
[PSyAD] Add support for logical scalar arguments in test harness generation.
`tl_poly1d_vert_adv_kernel_type` has a scalar argument of `logical` type. As it's a scalar it is read only but the test harness still needs to be extended to support this case.

Pull request:
(Closes #2013) Support logical, scalar arguments in psyad test-harness generation

Patch:
diff --git a/changelog b/changelog
index d76dc1944..62509bfa3 100644
--- a/changelog
+++ b/changelog
@@ -392,6 +392,9 @@
 	standalone driver for an extracted kernel. All required modules are
 	now inlined into the driver source file.
 
+	132) PR #2086 for #2013. Adds support for logical scalar arguments
+	in PSyAD test harness generation.
+
 release 2.3.1 17th of June 2022
 
 	1) PR #1747 for #1720. Adds support for If blocks to PSyAD.
diff --git a/doc/psyad/user_guide/implementation.rst b/doc/psyad/user_guide/implementation.rst
index 647511cee..544620175 100644
--- a/doc/psyad/user_guide/implementation.rst
+++ b/doc/psyad/user_guide/implementation.rst
@@ -556,14 +556,16 @@ Steps 1, 3, 5 and 6 are described in more detail below.
 Initialisation
 --------------
 
-All arguments to the TL kernel are initialised with pseudo-random numbers
+All real arguments to the TL kernel are initialised with pseudo-random numbers
 in the interval :math:`[0.0,1.0]` using the Fortran `random_number` intrinsic
-function. If the LFRic API is selected then only scalar and field arguments
+function. If the LFRic API is selected then only real scalar and field arguments
 are initialised in this way since arguments such as dof-maps contain
 essential information derived from the model configuration. In addition,
 those arguments flagged by the user (see :ref:`geom_kernel_args`) as
 containing geometric information (i.e. mesh coordinates or panel IDs) are
 passed through to the kernel from the Algorithm layer without modification.
+Integer and logical scalar arguments are currently just set to `1` and `.False`,
+respectively. This limitation is the subject of issue #2087.
 
 Inner Products
 --------------
diff --git a/doc/psyad/user_guide/psyad_command.rst b/doc/psyad/user_guide/psyad_command.rst
index 69cfb6753..e057fc933 100644
--- a/doc/psyad/user_guide/psyad_command.rst
+++ b/doc/psyad/user_guide/psyad_command.rst
@@ -1,7 +1,7 @@
 .. -----------------------------------------------------------------------------
 .. BSD 3-Clause License
 ..
-.. Copyright (c) 2021-2022, Science and Technology Facilities Council.
+.. Copyright (c) 2021-2023, Science and Technology Facilities Council.
 .. All rights reserved.
 ..
 .. Redistribution and use in source and binary forms, with or without
@@ -182,8 +182,10 @@ that must be called from within a full LFRic mini-app. Please see the
 Kernel Arguments Containing Geometric Information
 +++++++++++++++++++++++++++++++++++++++++++++++++
 
-By default, the test harness code generated by PSyAD populates all scalar
-and field arguments to a kernel with pseudo-random data. However, certain
+By default, the test harness code generated by PSyAD populates all real,
+scalar and field arguments to a kernel with pseudo-random data. (Integer and
+logical scalar arguments are currently set to `1` and `.False.`, respectively
+- see issue #2087.) However, certain
 LFRic kernels have arguments which carry information on the geometry of the
 simulation mesh (either or both of panel IDs and mesh coordinates) and these
 must be preserved if the kernel is to execute
diff --git a/src/psyclone/domain/lfric/kern_call_invoke_arg_list.py b/src/psyclone/domain/lfric/kern_call_invoke_arg_list.py
index 9e9530a07..a11713b2b 100644
--- a/src/psyclone/domain/lfric/kern_call_invoke_arg_list.py
+++ b/src/psyclone/domain/lfric/kern_call_invoke_arg_list.py
@@ -153,6 +153,8 @@ class KernCallInvokeArgList(ArgOrdering):
             datatype = LFRicTypes("LFRicRealScalarDataType")()
         elif scalar_arg.intrinsic_type == "integer":
             datatype = LFRicTypes("LFRicIntegerScalarDataType")()
+        elif scalar_arg.intrinsic_type == "logical":
+            datatype = LFRicTypes("LFRicLogicalScalarDataType")()
         else:
             raise NotImplementedError(
                 f"Scalar of type '{scalar_arg.intrinsic_type}' not supported.")
diff --git a/src/psyclone/psyad/domain/lfric/lfric_adjoint_harness.py b/src/psyclone/psyad/domain/lfric/lfric_adjoint_harness.py
index 2b17e7c2f..08bcc65e4 100644
--- a/src/psyclone/psyad/domain/lfric/lfric_adjoint_harness.py
+++ b/src/psyclone/psyad/domain/lfric/lfric_adjoint_harness.py
@@ -38,7 +38,8 @@
 
 from fparser import api as fpapi
 from psyclone.core import AccessType
-from psyclone.domain.lfric import LFRicConstants
+from psyclone.domain.lfric import LFRicConstants, LFRicTypes
+
 from psyclone.domain.lfric.algorithm.psyir import (
     LFRicAlgorithmInvokeCall, LFRicBuiltinFunctorFactory, LFRicKernelFunctor)
 from psyclone.domain.lfric.algorithm.lfric_alg import LFRicAlg
@@ -46,8 +47,9 @@ from psyclone.errors import InternalError, GenerationError
 from psyclone.psyad.domain.common.adjoint_utils import (create_adjoint_name,
                                                         create_real_comparison,
                                                         find_container)
-from psyclone.psyir.nodes import (Call, Reference, ArrayReference, Assignment,
-                                  Literal, BinaryOperation, Routine)
+from psyclone.psyir.nodes import (
+    IntrinsicCall, Reference, ArrayReference, Assignment,
+    Literal, BinaryOperation, Routine)
 from psyclone.psyir.symbols import (ImportInterface, ContainerSymbol,
                                     ScalarType, ArrayType, RoutineSymbol,
                                     DataTypeSymbol, DataSymbol, DeferredType)
@@ -79,6 +81,9 @@ def _compute_lfric_inner_products(prog, scalars, field_sums, sum_sym):
                                     Literal("0.0", sum_sym.datatype)))
     for scalar in scalars:
         # Compute the product of the pair of scalars: scalar[0]*scalar[1]
+        # (unless they are boolean).
+        if scalar[0].datatype.intrinsic == ScalarType.Intrinsic.BOOLEAN:
+            continue
         prod = BinaryOperation.create(BinaryOperation.Operator.MUL,
                                       Reference(scalar[0]),
                                       Reference(scalar[1]))
@@ -324,6 +329,54 @@ def _init_operators_random(operators, table):
     return kernel_list
 
 
+def _init_scalar_value(scalar_arg, routine, input_symbols):
+    '''
+    Extends the supplied Routine with the necessary statements to initialise
+    the supplied scalar argument. If that argument appears in the
+    `input_symbols` dict then its initial value is also assigned to the Symbol
+    in the dict entry.
+
+    :param scalar_arg: the scalar kernel argument to initialise.
+    :type scalar_arg: :py:class:`psyclone.psyir.symbols.DataSymbol`
+    :param routine: the routine to which to add assignments.
+    :type routine: :py:class:`psyclone.psyir.nodes.Routine`
+    :param input_symbols: dict containing those kernel arguments for which we \
+                          need to keep copies of their input values.
+    :type input_symbols: Dict[str, \
+                              :py:class:`psyclone.psyir.symbols.DataSymbol`]
+
+    :raises InternalError: if the type of the scalar argument is not supported.
+
+    '''
+    if scalar_arg.datatype.intrinsic == ScalarType.Intrinsic.REAL:
+        # TODO #2087 - we use the RANDOM_NUMBER intrinsic but make no attempt
+        # to scale the resulting value.
+        routine.addchild(
+            IntrinsicCall.create(IntrinsicCall.Intrinsic.RANDOM_NUMBER,
+                                 [Reference(scalar_arg)]))
+    elif scalar_arg.datatype.intrinsic == ScalarType.Intrinsic.BOOLEAN:
+        # TODO #2087 - just set the variable to False for the moment.
+        routine.addchild(Assignment.create(
+            Reference(scalar_arg),
+            Literal("false", LFRicTypes("LFRicLogicalScalarDataType")())))
+    elif scalar_arg.datatype.intrinsic == ScalarType.Intrinsic.INTEGER:
+        # TODO #2087 - just set the variable to 1 for the moment.
+        routine.addchild(Assignment.create(
+            Reference(scalar_arg),
+            Literal("1", LFRicTypes("LFRicIntegerScalarDataType")())))
+    else:
+        raise InternalError(
+            f"_init_scalar_value: only scalars of REAL, INTEGER or BOOLEAN "
+            f"type are supported but got symbol '{scalar_arg.name}' of type "
+            f"'{scalar_arg.datatype}'.")
+
+    if scalar_arg.name in input_symbols:
+        # We need to keep a copy of the input value of this argument.
+        input_sym = input_symbols[scalar_arg.name]
+        routine.addchild(Assignment.create(Reference(input_sym),
+                                           Reference(scalar_arg)))
+
+
 def _validate_geom_arg(kern, arg_idx, name, valid_spaces, vec_len):
     '''
     Check that the argument at the supplied index is consistent with the
@@ -531,9 +584,8 @@ def generate_lfric_adjoint_harness(tl_psyir, coord_arg_idx=None,
         kernel_input_arg_list.append(sym)
 
     # Initialise argument values and keep copies.
-    # Scalars - we use the Fortran 'random_number' intrinsic directly.
-    # TODO #1345 - this is Fortran specific.
-    random_num = RoutineSymbol("random_number")
+
+    # Scalars.
     for sym in kern_args.scalars:
         idx = kern_args.arglist.index(sym.name)
         if (kern_args.metadata_index_from_actual_index(idx) in
@@ -541,10 +593,8 @@ def generate_lfric_adjoint_harness(tl_psyir, coord_arg_idx=None,
             # This kernel argument is not modified by the test harness
             # because it contains geometry information.
             continue
-        routine.addchild(Call.create(random_num, [Reference(sym)]))
-        input_sym = table.lookup(sym.name+"_input")
-        routine.addchild(Assignment.create(Reference(input_sym),
-                                           Reference(sym)))
+        _init_scalar_value(sym, routine, input_symbols)
+
     # Fields.
     kernel_list = _init_fields_random(kernel_input_arg_list, input_symbols,
                                       table)


Test patch:
diff --git a/src/psyclone/tests/domain/lfric/conftest.py b/src/psyclone/tests/domain/lfric/conftest.py
index b59d41b02..cd0a60ef1 100644
--- a/src/psyclone/tests/domain/lfric/conftest.py
+++ b/src/psyclone/tests/domain/lfric/conftest.py
@@ -58,14 +58,15 @@ def dynkern_fixture():
     mdata_code = '''
 module testkern_field_mod
   type, extends(kernel_type) :: testkern_field_type
-     type(arg_type), meta_args(7) =                               &
+     type(arg_type), meta_args(8) =                               &
           (/ arg_type(gh_scalar, gh_real,    gh_read),            &
              arg_type(gh_field,  gh_real,    gh_readinc, w0),     &
              arg_type(gh_field,  gh_real,    gh_inc,     w1),     &
              arg_type(gh_field*3,gh_integer, gh_read,    w2),     &
              arg_type(gh_field,  gh_integer, gh_write,   wtheta), &
              arg_type(gh_field,  gh_integer, gh_read,    w3),     &
-             arg_type(gh_scalar, gh_integer, gh_read)             &
+             arg_type(gh_scalar, gh_integer, gh_read),            &
+             arg_type(gh_scalar, gh_logical, gh_read)             &
            /)
      type(func_type), dimension(2) :: meta_funcs =  &
           (/ func_type(w1, gh_basis),               &
diff --git a/src/psyclone/tests/domain/lfric/kern_call_invoke_arg_list_test.py b/src/psyclone/tests/domain/lfric/kern_call_invoke_arg_list_test.py
index fa06bc09e..b2377a2c0 100644
--- a/src/psyclone/tests/domain/lfric/kern_call_invoke_arg_list_test.py
+++ b/src/psyclone/tests/domain/lfric/kern_call_invoke_arg_list_test.py
@@ -64,16 +64,16 @@ def test_kcial_generate(dynkern):
     kcial = KernCallInvokeArgList(dynkern, table)
     kcial.generate()
     assert len(kcial.fields) == 5
-    assert len(kcial.scalars) == 2
+    assert len(kcial.scalars) == 3
     # Check that we can call it repeatedly.
     kcial.generate()
     assert len(kcial.fields) == 5
     # Check that an unsupported scalar type gives the expected error.
-    dynkern.arguments.args[0]._intrinsic_type = 'boolean'
+    dynkern.arguments.args[0]._intrinsic_type = 'wrong'
     kcial = KernCallInvokeArgList(dynkern, table)
     with pytest.raises(NotImplementedError) as err:
         kcial.generate()
-    assert "Scalar of type 'boolean' not supported" in str(err.value)
+    assert "Scalar of type 'wrong' not supported" in str(err.value)
 
 
 def test_kcial_generate_operator(dynkern_op):
diff --git a/src/psyclone/tests/psyad/domain/lfric/test_lfric_adjoint_harness.py b/src/psyclone/tests/psyad/domain/lfric/test_lfric_adjoint_harness.py
index 8246675b6..d9ec943d9 100644
--- a/src/psyclone/tests/psyad/domain/lfric/test_lfric_adjoint_harness.py
+++ b/src/psyclone/tests/psyad/domain/lfric/test_lfric_adjoint_harness.py
@@ -39,7 +39,7 @@
 import pytest
 from fparser import api as fpapi
 
-from psyclone.domain.lfric import LFRicSymbolTable, LFRicConstants
+from psyclone.domain.lfric import LFRicSymbolTable, LFRicConstants, LFRicTypes
 from psyclone.domain.lfric.algorithm import (
     LFRicBuiltinFunctor, LFRicAlg, LFRicBuiltinFunctorFactory,
     LFRicKernelFunctor)
@@ -50,10 +50,11 @@ from psyclone.psyad.domain.lfric.lfric_adjoint_harness import (
     _compute_field_inner_products,
     _init_fields_random,
     _init_operators_random,
+    _init_scalar_value,
     _validate_geom_arg,
     generate_lfric_adjoint_harness)
-from psyclone.psyir.nodes import Routine, Literal, Assignment
-from psyclone.psyir.symbols import (DataSymbol, REAL_TYPE,
+from psyclone.psyir import nodes
+from psyclone.psyir.symbols import (DataSymbol, REAL_TYPE, BOOLEAN_TYPE,
                                     ArrayType, DataTypeSymbol, DeferredType,
                                     INTEGER_TYPE, ContainerSymbol,
                                     ImportInterface, ScalarType)
@@ -69,28 +70,32 @@ def lfric_consts_fixture():
 
 def test_compute_inner_products_scalars(fortran_writer):
     '''Test that _compute_lfric_inner_products generates the expected code
-    for scalars.'''
+    for scalars and ignores any of boolean type.'''
     table = LFRicSymbolTable()
-    prog = Routine.create("test_prog", table, [], is_program=True)
+    prog = nodes.Routine.create("test_prog", table, [], is_program=True)
     sum_sym = table.new_symbol(root_name="my_sum",
                                symbol_type=DataSymbol, datatype=REAL_TYPE)
     sym1 = table.new_symbol(root_name="var1", symbol_type=DataSymbol,
                             datatype=REAL_TYPE)
     sym2 = table.new_symbol(root_name="var2", symbol_type=DataSymbol,
                             datatype=REAL_TYPE)
-    _compute_lfric_inner_products(prog, [(sym1, sym1), (sym1, sym2)], [],
-                                  sum_sym)
+    sym3 = table.new_symbol(root_name="var3", symbol_type=DataSymbol,
+                            datatype=BOOLEAN_TYPE)
+    _compute_lfric_inner_products(
+        prog, [(sym1, sym1), (sym1, sym2), (sym3, sym3)], [], sum_sym)
     gen = fortran_writer(prog)
+    # The resulting code should not include var3 since it is boolean.
     assert ("  my_sum = 0.0\n"
             "  my_sum = my_sum + var1 * var1\n"
-            "  my_sum = my_sum + var1 * var2\n" in gen)
+            "  my_sum = my_sum + var1 * var2\n\n"
+            "end program" in gen)
 
 
 def test_compute_inner_products_fields(fortran_writer):
     '''Test that _compute_lfric_inner_products generates the expected code
     when supplied with symbols representing the innerproducts of fields.'''
     table = LFRicSymbolTable()
-    prog = Routine.create("test_prog", table, [], is_program=True)
+    prog = nodes.Routine.create("test_prog", table, [], is_program=True)
     sum_sym = table.new_symbol(root_name="my_sum",
                                symbol_type=DataSymbol, datatype=REAL_TYPE)
     sym1 = table.new_symbol(root_name="ip1", symbol_type=DataSymbol,
@@ -118,7 +123,7 @@ def test_compute_field_inner_products(fortran_writer, type_map):
     assignments and functors for fields.'''
     bin_factory = LFRicBuiltinFunctorFactory.get()
     table = LFRicSymbolTable()
-    prog = Routine.create("test_prog", table, [], is_program=True)
+    prog = nodes.Routine.create("test_prog", table, [], is_program=True)
     csym = table.new_symbol(type_map["field"]["module"],
                             symbol_type=ContainerSymbol)
     fld_type = table.new_symbol(type_map["field"]["type"],
@@ -140,7 +145,7 @@ def test_compute_field_inner_products(fortran_writer, type_map):
                       bin_factory._get_builtin_class("x_innerproduct_x"))
     assert isinstance(functors[1],
                       bin_factory._get_builtin_class("x_innerproduct_y"))
-    assert isinstance(prog.children[0], Assignment)
+    assert isinstance(prog.children[0], nodes.Assignment)
     code = fortran_writer(prog)
     assert "field1_inner_prod = 0.0_r_def" in code
     assert "field1_field2_inner_prod = 0.0_r_def" in code
@@ -151,7 +156,7 @@ def test_compute_field_vector_inner_products(fortran_writer, type_map):
     assignments and functors for field vectors.'''
     bin_factory = LFRicBuiltinFunctorFactory.get()
     table = LFRicSymbolTable()
-    prog = Routine.create("test_prog", table, [], is_program=True)
+    prog = nodes.Routine.create("test_prog", table, [], is_program=True)
     csym = table.new_symbol(type_map["field"]["module"],
                             symbol_type=ContainerSymbol)
     fld_type = table.new_symbol(type_map["field"]["type"],
@@ -185,7 +190,7 @@ def test_compute_field_inner_products_errors(type_map):
     '''Check that _compute_field_inner_products raises the expected errors
     when passed incorrect arguments.'''
     table = LFRicSymbolTable()
-    prog = Routine.create("test_prog", table, [], is_program=True)
+    prog = nodes.Routine.create("test_prog", table, [], is_program=True)
     csym = table.new_symbol(type_map["field"]["module"],
                             symbol_type=ContainerSymbol)
     fld_type = table.new_symbol(type_map["field"]["type"],
@@ -263,7 +268,7 @@ def test_init_fields_random_vector(type_map):
     assert len(kernels) == 6
     for idx in range(3):
         kidx = 2*idx
-        lit = Literal(f"{idx+1}", idef_type)
+        lit = nodes.Literal(f"{idx+1}", idef_type)
         assert isinstance(kernels[kidx], LFRicBuiltinFunctor)
         assert kernels[kidx].symbol.name == "setval_random"
         assert kernels[kidx].children[0].symbol.name == "field1"
@@ -317,6 +322,50 @@ def test_init_operators_random(type_map):
             is csym)
 
 
+# _init_scalar_values
+
+def test_init_scalar_value(monkeypatch):
+    '''Check that _init_scalar_value() adds the expected nodes to the supplied
+    Routine.'''
+    table = LFRicSymbolTable()
+    routine = nodes.Routine("testkern_code", symbol_table=table)
+    sym1 = DataSymbol("my_real1", LFRicTypes("LFRicRealScalarDataType")())
+    sym2 = DataSymbol("my_int2", LFRicTypes("LFRicIntegerScalarDataType")())
+    sym2_input = DataSymbol("my_int2_input",
+                            LFRicTypes("LFRicIntegerScalarDataType")())
+    table.add(sym2_input)
+    _init_scalar_value(sym1, routine, {})
+    assert len(routine.children) == 1
+    # We should get a call to random_number for a real scalar.
+    assert isinstance(routine[0], nodes.Call)
+    _init_scalar_value(sym2, routine, {"my_int2": sym2_input})
+    assert len(routine.children) == 3
+    # An integer should just be assigned the value 1 (TODO #2087)
+    assert isinstance(routine[1], nodes.Assignment)
+    assert routine[1].rhs.value == "1"
+    # and we should store this value as it's listed as an 'input'.
+    assert isinstance(routine[2], nodes.Assignment)
+    assert routine[2].lhs.symbol.name == "my_int2_input"
+    # A logical argument should just be assigned False (TODO #2087)
+    sym3 = DataSymbol("my_bool", LFRicTypes("LFRicLogicalScalarDataType")())
+    _init_scalar_value(sym3, routine, {})
+    assert isinstance(routine[3], nodes.Assignment)
+    assert routine[3].rhs.value == "false"
+    # Unrecognised type of scalar. This is tricky to reproduce so we create
+    # a new class that has a 'name' attribute and monkeypatch the 'intrinsic'
+    # property of the datatype.
+    sym4 = DataSymbol("my_var", LFRicTypes("LFRicRealScalarDataType")())
+
+    class broken_type:
+        def __init__(self):
+            self.name = "wrong"
+    monkeypatch.setattr(sym4.datatype, "intrinsic", broken_type())
+    with pytest.raises(InternalError) as err:
+        _init_scalar_value(sym4, routine, {})
+    assert ("scalars of REAL, INTEGER or BOOLEAN type are supported but got "
+            "symbol 'my_var' of type 'Scalar<wrong" in str(err.value))
+
+
 # _validate_geom_arg
 
 def test_validate_geom_arg():
@@ -445,7 +494,7 @@ def test_generate_lfric_adjoint_harness(fortran_reader, fortran_writer):
     expected test-harness code.'''
     tl_psyir = fortran_reader.psyir_from_source(TL_CODE)
     psyir = generate_lfric_adjoint_harness(tl_psyir)
-    gen = fortran_writer(psyir)
+    gen = fortran_writer(psyir).lower()
     assert "module adjoint_test_mod" in gen
     assert "subroutine adjoint_test(mesh, chi, panel_id)" in gen
     # We should have a field, a copy of that field and an inner-product value
@@ -468,7 +517,7 @@ def test_generate_lfric_adjoint_harness(fortran_reader, fortran_writer):
     # The TL kernel must then be called and the inner-product of the result
     # computed.
     assert "field_2_inner_prod = 0.0_r_def" in gen
-    assert ("    ! Initialise arguments and call the tangent-linear kernel.\n"
+    assert ("    ! initialise arguments and call the tangent-linear kernel.\n"
             "    call invoke(setval_random(field_2), setval_x(field_2_input, "
             "field_2), testkern_type(rscalar_1, field_2), x_innerproduct_x("
             "field_2_inner_prod, field_2))\n" in gen)
@@ -500,7 +549,7 @@ def test_generate_lfric_adj_test_quadrature(fortran_reader):
         "     integer :: gh_shape = gh_quadrature_xyoz\n")
     tl_psyir = fortran_reader.psyir_from_source(new_code)
     psyir = generate_lfric_adjoint_harness(tl_psyir)
-    routine = psyir.walk(Routine)[0]
+    routine = psyir.walk(nodes.Routine)[0]
     for sym in routine.symbol_table.datasymbols:
         # All input variables should be either scalars or fields (operators
         # are never active).