diff --git a/loopy/frontend/fortran/__init__.py b/loopy/frontend/fortran/__init__.py index 00dc837e16ad2a414a13c1ceaf4f36f3f3fb3049..10c390ff1e937d5a708e8667bc8598e1c91ef861 100644 --- a/loopy/frontend/fortran/__init__.py +++ b/loopy/frontend/fortran/__init__.py @@ -235,8 +235,13 @@ def parse_transformed_fortran(source, free_form=True, strict=True, def parse_fortran(source, filename="", free_form=True, strict=True, - seq_dependencies=None, auto_dependencies=None, target=None): + seq_dependencies=None, auto_dependencies=None, target=None, + all_names_known=True): """ + :arg all_names_known: if set to *False*, enter an undocumented mode + in which Fortran parsing will try to tolerate unknown names. + If used, ``loopy.frontend.fortran.translator.specialize_fortran_division`` + must be called as soon as all names are known. :returns: a list of :class:`loopy.LoopKernel` objects """ @@ -268,7 +273,8 @@ def parse_fortran(source, filename="", free_form=True, strict=True, "and returned invalid data (Sorry!)") from loopy.frontend.fortran.translator import F2LoopyTranslator - f2loopy = F2LoopyTranslator(filename, target=target) + f2loopy = F2LoopyTranslator( + filename, target=target, all_names_known=all_names_known) f2loopy(tree) return f2loopy.make_kernels(seq_dependencies=seq_dependencies) diff --git a/loopy/frontend/fortran/translator.py b/loopy/frontend/fortran/translator.py index 78eddfb549c4ebabbf933abd7832a449a44dbeec..abd1b888f1781595118cfc1b797b3a9cc0254391 100644 --- a/loopy/frontend/fortran/translator.py +++ b/loopy/frontend/fortran/translator.py @@ -28,11 +28,13 @@ import loopy as lp import numpy as np from warnings import warn from loopy.frontend.fortran.tree import FTreeWalkerBase +from loopy.diagnostic import warn_with_kernel from loopy.frontend.fortran.diagnostic import ( TranslationError, TranslatorWarning) import islpy as isl from islpy import dim_type -from loopy.symbolic import IdentityMapper +from loopy.symbolic import (IdentityMapper, RuleAwareIdentityMapper, + SubstitutionRuleMappingContext) from loopy.diagnostic import LoopyError from loopy.kernel.instruction import LegacyStringInstructionTag from pymbolic.primitives import Wildcard @@ -80,7 +82,9 @@ class SubscriptIndexBaseShifter(IdentityMapper): # {{{ scope class Scope: - def __init__(self, subprogram_name, arg_names=set()): + def __init__(self, subprogram_name, arg_names=None): + if arg_names is None: + arg_names = set() self.subprogram_name = subprogram_name # map name to data @@ -193,13 +197,66 @@ class Scope: # }}} +# {{{ fortran division specializers + +class FortranDivisionToFloorDiv(IdentityMapper): + def map_fortran_division(self, expr, *args): + from warnings import warn + from loopy.diagnostic import LoopyWarning + warn( + "Integer division in Fortran do loop bound. " + "Loopy currently forces this to integers and gets it wrong for " + "negative arguments.", LoopyWarning) + from pymbolic.primitives import FloorDiv + return FloorDiv( + self.rec(expr.numerator, *args), + self.rec(expr.denominator, *args)) + + +class FortranDivisionSpecializer(RuleAwareIdentityMapper): + def __init__(self, rule_mapping_context, kernel): + super().__init__(rule_mapping_context) + from loopy.type_inference import TypeInferenceMapper + self.infer_type = TypeInferenceMapper(kernel) + self.kernel = kernel + + def map_fortran_division(self, expr, *args): + # We remove all these before type inference ever sees them. + num_dtype = self.infer_type(expr.numerator).numpy_dtype + den_dtype = self.infer_type(expr.denominator).numpy_dtype + + from pymbolic.primitives import Quotient, FloorDiv + if num_dtype.kind in "iub" and den_dtype.kind in "iub": + warn_with_kernel(self.kernel, + "fortran_int_div", + "Integer division in Fortran code. Loopy currently gets this " + "wrong for negative arguments.") + return FloorDiv( + self.rec(expr.numerator, *args), + self.rec(expr.denominator, *args)) + + else: + return Quotient( + self.rec(expr.numerator, *args), + self.rec(expr.denominator, *args)) + + +def specialize_fortran_division(knl): + rmc = SubstitutionRuleMappingContext( + knl.substitutions, knl.get_var_name_generator()) + return FortranDivisionSpecializer(rmc, knl).map_kernel(knl) + +# }}} + + # {{{ translator class F2LoopyTranslator(FTreeWalkerBase): - def __init__(self, filename, target=None): + def __init__(self, filename, target=None, all_names_known=True): FTreeWalkerBase.__init__(self, filename) self.target = target + self.all_names_known = all_names_known self.scope_stack = [] @@ -365,7 +422,9 @@ class F2LoopyTranslator(FTreeWalkerBase): for name, data in node.stmts: name, = name assert name not in scope.data - scope.data[name] = [self.parse_expr(node, i) for i in data] + scope.data[name] = [ + scope.process_expression_for_loopy( + self.parse_expr(node, i)) for i in data] return [] @@ -461,7 +520,9 @@ class F2LoopyTranslator(FTreeWalkerBase): cond_var = var(cond_name) self.add_expression_instruction( - cond_var, self.parse_expr(node, node.expr)) + cond_var, + scope.process_expression_for_loopy( + self.parse_expr(node, node.expr))) cond_expr = cond_var if context_cond is not None: @@ -531,9 +592,10 @@ class F2LoopyTranslator(FTreeWalkerBase): % (loop_var, iname_dtype, self.index_dtype)) scope.use_name(loop_var) - loop_bounds = self.parse_expr( - node, - loop_bounds, min_precedence=self.expr_parser._PREC_FUNC_ARGS) + loop_bounds = scope.process_expression_for_loopy( + self.parse_expr( + node, + loop_bounds, min_precedence=self.expr_parser._PREC_FUNC_ARGS)) if len(loop_bounds) == 2: start, stop = loop_bounds @@ -598,7 +660,8 @@ class F2LoopyTranslator(FTreeWalkerBase): isl.Constraint.inequality_from_aff( iname_rel_aff(space, loopy_loop_var, "<=", - aff_from_expr(space, stop-start))))) + aff_from_expr(space, FortranDivisionToFloorDiv()( + stop-start)))))) from pymbolic import var scope.active_iname_aliases[loop_var] = \ @@ -717,6 +780,9 @@ class F2LoopyTranslator(FTreeWalkerBase): lang_version=MOST_RECENT_LANGUAGE_VERSION ) + if self.all_names_known: + knl = specialize_fortran_division(knl) + from loopy.loop import fuse_loop_domains knl = fuse_loop_domains(knl) knl = lp.fold_constants(knl) diff --git a/loopy/frontend/fortran/tree.py b/loopy/frontend/fortran/tree.py index c33578dc844d1a77b084d8cf83cb5009cccc489d..cc41afc6c5bca7d202548490983169374dcf1535 100644 --- a/loopy/frontend/fortran/tree.py +++ b/loopy/frontend/fortran/tree.py @@ -23,6 +23,14 @@ THE SOFTWARE. import re from loopy.diagnostic import LoopyError +from loopy.symbolic import IdentityMapper, FortranDivision + + +class DivisionToFortranDivisionMapper(IdentityMapper): + def map_quotient(self, expr): + return FortranDivision( + self.rec(expr.numerator), + self.rec(expr.denominator)) class FTreeWalkerBase: @@ -84,7 +92,8 @@ class FTreeWalkerBase: def parse_expr(self, node, expr_str, **kwargs): try: - return self.expr_parser(expr_str, **kwargs) + return DivisionToFortranDivisionMapper()( + self.expr_parser(expr_str, **kwargs)) except Exception as e: raise LoopyError( "Error parsing expression '%s' on line %d of '%s': %s" diff --git a/loopy/symbolic.py b/loopy/symbolic.py index 002092eedda1940d8ac5ed9a6393f7f7113fcc28..15570c4607f90fa3fc03ea6d2cb3c4ab11458326 100644 --- a/loopy/symbolic.py +++ b/loopy/symbolic.py @@ -146,6 +146,11 @@ class IdentityMapperMixin: map_rule_argument = map_group_hw_index + def map_fortran_division(self, expr, *args, **kwargs): + return type(expr)( + self.rec(expr.numerator, *args, **kwargs), + self.rec(expr.denominator, *args, **kwargs)) + class IdentityMapper(IdentityMapperBase, IdentityMapperMixin): pass @@ -197,6 +202,8 @@ class WalkMapper(WalkMapperBase): map_rule_argument = map_group_hw_index + map_fortran_division = WalkMapperBase.map_quotient + class CallbackMapper(CallbackMapperBase, IdentityMapper): map_reduction = CallbackMapperBase.map_constant @@ -208,6 +215,8 @@ class CombineMapper(CombineMapperBase): map_linear_subscript = CombineMapperBase.map_subscript + map_fortran_division = CombineMapperBase.map_quotient + class SubstitutionMapper( CSECachingMapperMixin, SubstitutionMapperBase, IdentityMapperMixin): @@ -265,6 +274,11 @@ class StringifyMapper(StringifyMapperBase): return "cast({}, {})".format( repr(expr.type), self.rec(expr.child, PREC_NONE)) + def map_fortran_division(self, expr, enclosing_prec): + from pymbolic.mapper.stringifier import PREC_NONE + result = self.map_quotient(expr, PREC_NONE) + return f"[FORTRANDIV]({result})" + class EqualityPreservingStringifyMapper(StringifyMapperBase): """ @@ -357,6 +371,8 @@ class DependencyMapper(DependencyMapperBase): def map_literal(self, expr): return set() + map_fortran_division = DependencyMapperBase.map_quotient + class SubstitutionRuleExpander(IdentityMapper): def __init__(self, rules): @@ -748,6 +764,20 @@ class RuleArgument(LoopyExpressionBase): mapper_method = intern("map_rule_argument") + +class FortranDivision(p.QuotientBase, LoopyExpressionBase): + """This exists for the benefit of the Fortran frontend, which specializes + to floating point division for floating point inputs and round-to-zero + division for integer inputs. Despite the name, this would also be usable + for C semantics. (:mod:`loopy` division semantics match Python's.) + + .. note:: + + This is not a documented expression node type. It may disappear + at any moment. + """ + mapper_method = "map_fortran_division" + # }}} diff --git a/test/test_fortran.py b/test/test_fortran.py index ff0855fa239218b05fae25b37cbe0f80d4f0362e..c72afa78f08b1b802a608308b8f866d5702d2d38 100644 --- a/test/test_fortran.py +++ b/test/test_fortran.py @@ -507,6 +507,28 @@ def test_precompute_some_exist(ctx_factory): lp.auto_test_vs_ref(ref_knl, ctx, knl, parameters=dict(n=128, m=128, ell=128)) +def test_division_in_shapes(ctx_factory): + fortran_src = """ + subroutine halve(m, a) + implicit none + integer m, i, j + real*8 a(m/2,m/2) + do i = 1,m/2 + do j = 1,m/2 + a(i, j) = 2*a(i, j) + end do + end do + end subroutine + """ + knl, = lp.parse_fortran(fortran_src) + ref_knl = knl + + print(knl) + + ctx = ctx_factory() + lp.auto_test_vs_ref(ref_knl, ctx, knl, parameters=dict(m=128)) + + if __name__ == "__main__": if len(sys.argv) > 1: exec(sys.argv[1]) diff --git a/test/test_numa_diff.py b/test/test_numa_diff.py index 7e708f449c72ab93dc7ee342da055d4b5becb8a5..54fbe066e800f8e953bfefe4a5c2731079a3a63d 100644 --- a/test/test_numa_diff.py +++ b/test/test_numa_diff.py @@ -58,7 +58,9 @@ def test_gnuma_horiz_kernel(ctx_factory, ilp_multiple, Nq, opt_level): # noqa source = source.replace("datafloat", "real*4") hsv_r, hsv_s = [ - knl for knl in lp.parse_fortran(source, filename, seq_dependencies=False) + knl for knl in lp.parse_fortran( + source, filename, seq_dependencies=False, + all_names_known=False) if "KernelR" in knl.name or "KernelS" in knl.name ] hsv_r = lp.tag_instructions(hsv_r, "rknl") @@ -77,6 +79,9 @@ def test_gnuma_horiz_kernel(ctx_factory, ilp_multiple, Nq, opt_level): # noqa hsv = lp.assume(hsv, "elements >= 1") hsv = fix_euler_parameters(hsv, p_p0=1, p_Gamma=1.4, p_R=1) + from loopy.frontend.fortran.translator import specialize_fortran_division + hsv = specialize_fortran_division(hsv) + for name in ["Q", "rhsQ"]: hsv = set_q_storage_format(hsv, name)