Merge branch 'feature/measure-nonlinear-jacobian-apply' into 'master'

Feature/measure nonlinear jacobian apply

Closes #111 and #119

See merge request dominic/dune-perftool!244

applications/knl/poisson_dg/knl_poisson_dg.mini

@@ -3,6 +3,7 @@ __exec_suffix = deg{formcompiler.ufl_variants.degree}_{opcount_suffix}_level{for
 
 opcount_suffix = opcount, nonopcount | expand opcount
 {opcount_suffix} == opcount and {formcompiler.instrumentation_level} != 4 | exclude
+dune-opcounter_FOUND, 1 | expand opcount | cmake_guard
 
 # Calculate the size of the grid to equlibritate it to 100 MB/rank
 # Input parameters

applications/knl/poisson_dg_tensor/knl_poisson_dg_tensor.mini

@@ -3,6 +3,7 @@ __exec_suffix = deg{formcompiler.ufl_variants.degree}_{opcount_suffix}_level{for
 
 opcount_suffix = opcount, nonopcount | expand opcount
 {opcount_suffix} == opcount and {formcompiler.instrumentation_level} != 4 | exclude
+dune-opcounter_FOUND, 1 | expand opcount | cmake_guard
 
 # Calculate the size of the grid to equlibritate it to 100 MB/rank
 # Input parameters

applications/poisson_dg/poisson_dg.mini

@@ -3,6 +3,7 @@ __exec_suffix = deg{formcompiler.ufl_variants.degree}_{opcount_suffix}_level{for
 
 opcount_suffix = opcount, nonopcount | expand opcount
 {opcount_suffix} == opcount and {formcompiler.instrumentation_level} != 4 | exclude
+dune-opcounter_FOUND, 1 | expand opcount | cmake_guard
 
 # Calculate the size of the grid to equlibritate it to 100 MB/rank
 # Input parameters

applications/poisson_dg_tensor/poisson_dg_tensor.mini

@@ -3,6 +3,7 @@ __exec_suffix = deg{formcompiler.ufl_variants.degree}_{opcount_suffix}_level{for
 
 opcount_suffix = opcount, nonopcount | expand opcount
 {opcount_suffix} == opcount and {formcompiler.instrumentation_level} != 4 | exclude
+dune-opcounter_FOUND, 1 | expand opcount | cmake_guard
 
 # Calculate the size of the grid to equlibritate it to 100 MB/rank
 # Input parameters

applications/stokes_dg/stokes_dg.mini

@@ -3,6 +3,7 @@ __exec_suffix = deg{formcompiler.ufl_variants.v_degree}_{opcount_suffix}_level{f
 
 opcount_suffix = opcount, nonopcount | expand opcount
 {opcount_suffix} == opcount and {formcompiler.instrumentation_level} != 4 | exclude
+dune-opcounter_FOUND, 1 | expand opcount | cmake_guard
 
 # Calculate the size of the grid to equlibritate it to 100 MB/rank
 # Input parameters

cmake/modules/GeneratedSystemtests.cmake

@@ -96,7 +96,7 @@ function(dune_add_formcompiler_system_test)
                )
 
       set_tests_properties(${tname} PROPERTIES SKIP_RETURN_CODE 77)
-      set_tests_properties(${tname} PROPERTIES TIMEOUT 60)
+      set_tests_properties(${tname} PROPERTIES TIMEOUT 120)
     endif()
   endforeach()
 endfunction()

dune.module

@@ -8,3 +8,4 @@ Version: 0.0
 Maintainer: dominic.kempf@iwr.uni-heidelberg.de
 #depending on
 Depends: dune-testtools dune-pdelab dune-alugrid
+Suggests: dune-opcounter

dune/perftool/common/CMakeLists.txt

 install(FILES muladd_workarounds.hh
-              opcounter.hh
               timer.hh
+              timer_tsc.hh
+              timer_chrono.hh
               tsc.hh
               vectorclass.hh
         DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/perftool/common

dune/perftool/common/muladd_workarounds.hh

@@ -3,36 +3,14 @@
 
 /* We are currently having some issues with FMA nodes not being
  * eliminated correctly upon code generation. We "solve" the problem
- * for now with overloads of the mul_add function for scalars.
+ * for now with a generic implementation of  the mul_add function.
  */
 
-#include<dune/perftool/common/opcounter.hh>
 
-
-inline double mul_add(double op1, double& op2, double op3)
-{
-	return op1 * op2 + op3;
-}
-
-
-inline float mul_add(float op1, float& op2, float op3)
+template<typename T>
+inline T mul_add(T op1, T op2, T op3)
 {
-	return op1 * op2 + op3;
+  return op1 * op2 + op3;
 }
 
-
-#ifdef ENABLE_COUNTER
-
-oc::OpCounter<double> mul_add(oc::OpCounter<double> op1, oc::OpCounter<double>& op2, oc::OpCounter<double> op3)
-{
-	return op1 * op2 + op3;
-}
-
-oc::OpCounter<float> mul_add(oc::OpCounter<float> op1, oc::OpCounter<float>& op2, oc::OpCounter<float> op3)
-{
-	return op1 * op2 + op3;
-}
-
-#endif
-
 #endif

dune/perftool/common/timer_chrono.hh

@@ -7,7 +7,7 @@
 
 #include <chrono>
 
-#include <dune/perftool/common/opcounter.hh>
+#include <dune/opcounter/opcounter.hh>
 
 #define HP_TIMER_OPCOUNTER oc::OpCounter<double>
 

dune/perftool/common/timer_tsc.hh

@@ -6,7 +6,7 @@
 #endif
 
 #include <dune/perftool/common/tsc.hh>
-#include <dune/perftool/common/opcounter.hh>
+#include <dune/opcounter/opcounter.hh>
 
 #define HP_TIMER_DURATION(name) __hp_timer_##name##_duration
 #define HP_TIMER_STARTTIME(name) __hp_timer_##name##_start
@@ -76,16 +76,22 @@
 #ifdef ENABLE_COUNTER
 
 #define DUMP_TIMER(level,name,os,reset)\
-  if (HP_TIMER_DURATION(name) > 1e-12) \
-    os << #level << " " << ident << " " << #name << " time " << Dune::PDELab::TSC::seconds(HP_TIMER_DURATION(name)) << std::endl; \
-  HP_TIMER_OPCOUNTERS(name).reportOperations(os,#level,ident,#name,reset);
+  { \
+    std::string prefix = std::string(#level) + " " + ident + " " + std::string(#name); \
+    if (HP_TIMER_DURATION(name) > 1e-12) \
+      os << prefix << " time " << Dune::PDELab::TSC::seconds(HP_TIMER_DURATION(name)) << std::endl; \
+    HP_TIMER_OPCOUNTERS(name).reportOperations(os,prefix,reset); \
+  }
 
 #define DUMP_AND_ACCUMULATE_TIMER(level,name,os,reset,time,ops)  \
-  if (HP_TIMER_DURATION(name) > 1e-12) \
-    os << #level << " " << ident << " " << #name << " time " << Dune::PDELab::TSC::seconds(HP_TIMER_DURATION(name)) << std::endl; \
-  time += HP_TIMER_DURATION(name); \
-  ops += HP_TIMER_OPCOUNTERS(name); \
-  HP_TIMER_OPCOUNTERS(name).reportOperations(os,#level,ident,#name,reset);
+  { \
+    std::string prefix = std::string(#level) + " " + ident + " " + std::string(#name); \
+    if (HP_TIMER_DURATION(name) > 1e-12) \
+      os << prefix << " time " << Dune::PDELab::TSC::seconds(HP_TIMER_DURATION(name)) << std::endl; \
+    time += HP_TIMER_DURATION(name); \
+    ops += HP_TIMER_OPCOUNTERS(name); \
+    HP_TIMER_OPCOUNTERS(name).reportOperations(os,prefix,reset); \
+  }
 
 #elif defined ENABLE_HP_TIMERS
 

python/dune/perftool/loopy/target.py

@@ -30,7 +30,8 @@ import cgen
 
 
 def _type_to_op_counter_type(name):
-    return "oc::OpCounter<{}>".format(name)
+    include_file("dune/opcounter/opcounter.hh")
+    return "OpCounter::OpCounter<{}>".format(name)
 
 
 def dtype_floatingpoint():
@@ -63,6 +64,10 @@ def type_floatingpoint():
     return numpy_to_cpp_dtype(NumpyType(dtype).dtype.name)
 
 
+def type_context_floatingpoint():
+    return {np.float32: 'f', np.float64: 'd'}.get(dtype_floatingpoint())
+
+
 class DuneExpressionToCExpressionMapper(ExpressionToCExpressionMapper):
     def map_subscript(self, expr, type_context):
         arr = self.find_array(expr)
@@ -88,6 +93,10 @@ class DuneExpressionToCExpressionMapper(ExpressionToCExpressionMapper):
         return expr
 
     def map_constant(self, expr, type_context):
+        # We correct the type context to force all floating point literals to be of
+        # the type that we use throughout the computation.
+        if type_context in ("f", "d"):
+            type_context = type_context_floatingpoint()
         ret = ExpressionToCExpressionMapper.map_constant(self, expr, type_context)
         if get_option('opcounter'):
             if type_context in ("f", "d"):

python/dune/perftool/loopy/transformations/instrumentation.py

@@ -22,7 +22,7 @@ def _union(a):
     return frozenset.union(*a)
 
 
-def add_instrumentation(knl, match, identifier, level, filetag='operatorfile', operator=False):
+def add_instrumentation(knl, match, identifier, level, filetag='operatorfile', operator=False, depends_on=frozenset()):
     """ Transform loopy kernel to contain instrumentation code
 
     Arguments:
@@ -32,6 +32,9 @@ def add_instrumentation(knl, match, identifier, level, filetag='operatorfile', o
     identifier : The name of the counter to start and stop
     level : The instrumentation level this measurement is defined at
     filetag : The tag of the file that should contain the counter definitions
+    depends_on: Additional dependencies to add to the start instruction. This is used to correct
+                currently wrong behaviour of the transformation in cases where a lot of structure
+                of the instrumentation is known a priori.
     """
     # If the instrumentation level is not high enough, this is a no-op
     if level > get_option("instrumentation_level"):
@@ -53,6 +56,7 @@ def add_instrumentation(knl, match, identifier, level, filetag='operatorfile', o
     insn_inames = _intersect(tuple(i.within_inames for i in insns))
     other_inames = _union(tuple(i.within_inames for i in lp.find_instructions(knl, lp.match.Not(match))))
     within = _intersect((insn_inames, other_inames))
+    uniontags = _intersect(tuple(i.tags for i in insns))
 
     # Get a unique identifer - note that the same timer could be started and stopped several times
     # within one kernel...
@@ -67,8 +71,9 @@ def add_instrumentation(knl, match, identifier, level, filetag='operatorfile', o
                                  "HP_TIMER_START({});".format(identifier),
                                  id=start_id,
                                  within_inames=within,
-                                 depends_on=start_depends,
+                                 depends_on=depends_on.union(start_depends),
                                  boostable_into=frozenset(),
+                                 tags=uniontags,
                                  )
 
     # Add dependencies on the timing instructions
@@ -82,6 +87,7 @@ def add_instrumentation(knl, match, identifier, level, filetag='operatorfile', o
                                 within_inames=within,
                                 depends_on=frozenset(i.id for i in insns),
                                 boostable_into=frozenset(),
+                                tags=uniontags,
                                 )
 
     # Find all the instructions that should depend on stop
@@ -98,4 +104,7 @@ def add_instrumentation(knl, match, identifier, level, filetag='operatorfile', o
     other_insns = list(filter(lambda i: i.id not in [j.id for j in rewritten_insns], knl.instructions))
 
     # Add all the modified instructions into the kernel object
-    return knl.copy(instructions=rewritten_insns + other_insns + [start_insn, stop_insn])
+    knl = knl.copy(instructions=rewritten_insns + other_insns + [start_insn, stop_insn])
+
+    from loopy.kernel.creation import resolve_dependencies
+    return resolve_dependencies(knl)

python/dune/perftool/pdelab/driver/error.py

@@ -161,7 +161,7 @@ def accumulate_L2_squared():
 @preamble(section="error")
 def define_accumulated_L2_error(name):
     t = type_range()
-    return "{} {}(0.0);".format(t, name)
+    return "Dune::FieldVector<{}, 1> {}(0.0);".format(t, name)
 
 
 def name_accumulated_L2_error():
@@ -182,7 +182,7 @@ def compare_L2_squared():
             "if ({}.comm().rank() == 0){{".format(gv),
             "  std::cout << \"\\nl2errorsquared: \" << {} << std::endl << std::endl;".format(accum_error),
             "}",
-            "if (isnan({0}) or abs({0})>{1})".format(accum_error, get_option("compare_l2errorsquared")),
+            "if (isnan({0}[0]) or abs({0}[0])>{1})".format(accum_error, get_option("compare_l2errorsquared")),
             "  {} = true;".format(fail)]
 
 

python/dune/perftool/pdelab/driver/timings.py

@@ -8,6 +8,7 @@ from dune.perftool.generation import (cached,
                                       preamble,
                                       )
 from dune.perftool.pdelab.driver import (get_form_ident,
+                                         is_linear,
                                          name_initree,
                                          name_mpihelper,
                                          )
@@ -131,18 +132,19 @@ def apply_jacobian_timer():
     if get_option('instrumentation_level') >= 3:
         print_times.append("{}.dump_timers({}, {}, true);".format(lop_name, timestream, name_timing_identifier()))
 
-    if get_option('instrumentation_level') >= 2:
-        evaluation = ["HP_TIMER_START(apply_jacobian);",
-                      "{}.jacobian_apply({}, j);".format(n_go, v),
-                      "HP_TIMER_STOP(apply_jacobian);",
-                      "DUMP_TIMER({}, apply_jacobian, {}, true);".format(get_option("instrumentation_level"), timestream)]
-        evaluation.extend(print_times)
-    else:
+    if is_linear():
+        declaration = ["{} j({});".format(t_v, v), "j=0.0;"]
         evaluation = ["{}.jacobian_apply({}, j);".format(n_go, v)]
+    else:
+        declaration = ["{} j0({});".format(t_v, v), "j0=0.0;",
+                       "{} j1({});".format(t_v, v), "j1=0.0;"]
+        evaluation = ["{}.nonlinear_jacobian_apply({}, j0, j1);".format(n_go, v)]
 
-    evaluation = ["{} j({});".format(t_v, v), "j=0.0;"] + evaluation
+    if get_option('instrumentation_level') >= 2:
+        evaluation = ["HP_TIMER_START(apply_jacobian);"] + evaluation + ["HP_TIMER_STOP(apply_jacobian);", "DUMP_TIMER({}, apply_jacobian, {}, true);".format(get_option("instrumentation_level"), timestream)]
+        evaluation.extend(print_times)
 
-    return evaluation
+    return declaration + evaluation
 
 
 @preamble(section="timings")

python/dune/perftool/pdelab/driver/visitor.py

@@ -4,11 +4,19 @@ from dune.perftool.ufl.visitor import UFL2LoopyVisitor
 import pymbolic.primitives as prim
 
 
-@preamble
+@preamble(section="init")
 def driver_using_statement(what):
     return "using {};".format(what)
 
 
+@preamble(section="gridoperator")
+def set_lop_to_starting_time():
+    from dune.perftool.pdelab.driver import get_form_ident
+    from dune.perftool.pdelab.driver.gridoperator import name_localoperator
+    lop = name_localoperator(get_form_ident())
+    return "{}.setTime(0.0);".format(lop)
+
+
 class DriverUFL2PymbolicVisitor(UFL2LoopyVisitor):
     def __init__(self):
         from dune.perftool.pdelab import PDELabInterface
@@ -33,6 +41,7 @@ class DriverUFL2PymbolicVisitor(UFL2LoopyVisitor):
             from dune.perftool.pdelab.driver import get_form_ident
             from dune.perftool.pdelab.driver.gridoperator import name_localoperator
             lop = name_localoperator(get_form_ident())
+            set_lop_to_starting_time()
             return prim.Call(prim.Variable("{}.getTime".format(lop)), ())
         else:
             return UFL2LoopyVisitor.coefficient(self, o)

python/dune/perftool/pdelab/localoperator.py

@@ -572,8 +572,8 @@ def extract_kernel_from_cache(tag, name, signature, wrap_in_cgen=True, add_timin
     if add_timings and get_form_option("sumfact"):
         from dune.perftool.pdelab.signatures import assembler_routine_name
         kernel = add_instrumentation(kernel, lp.match.Tagged("sumfact_stage1"), "{}_kernel_stage1".format(assembler_routine_name()), 4)
-        kernel = add_instrumentation(kernel, lp.match.Tagged("sumfact_stage2"), "{}_kernel_quadratureloop".format(assembler_routine_name()), 4)
-        kernel = add_instrumentation(kernel, lp.match.Tagged("sumfact_stage3"), "{}_kernel_stage3".format(assembler_routine_name()), 4)
+        kernel = add_instrumentation(kernel, lp.match.Tagged("sumfact_stage2"), "{}_kernel_quadratureloop".format(assembler_routine_name()), 4, depends_on=frozenset({lp.match.Tagged("sumfact_stage1")}))
+        kernel = add_instrumentation(kernel, lp.match.Tagged("sumfact_stage3"), "{}_kernel_stage3".format(assembler_routine_name()), 4, depends_on=frozenset({lp.match.Tagged("sumfact_stage2")}))
 
     if wrap_in_cgen:
         # Wrap the kernel in something which can generate code

test/poisson/CMakeLists.txt

@@ -44,10 +44,12 @@ dune_add_formcompiler_system_test(UFLFILE poisson_dg.ufl
                                   )
 
 # 8. Poisson with operator counting
-dune_add_formcompiler_system_test(UFLFILE opcount_poisson_dg.ufl
-                                  BASENAME opcount_poisson_dg_symdiff
-                                  INIFILE opcount_poisson_dg_symdiff.mini
-                                  )
+if(dune-opcounter_FOUND)
+  dune_add_formcompiler_system_test(UFLFILE opcount_poisson_dg.ufl
+                                    BASENAME opcount_poisson_dg_symdiff
+                                    INIFILE opcount_poisson_dg_symdiff.mini
+                                    )
+endif()
 
 # 9. Poisson Test Case: DG quadrilaterals
 dune_add_formcompiler_system_test(UFLFILE poisson_dg_quadrilateral.ufl