diff --git a/python/dune/perftool/pdelab/driver.py b/python/dune/perftool/pdelab/driver.py
index 524ce57967e4c4c42f94d08b30c433ba2aab40ed..26de3b67de909ee8156d33aefafcaff3727c718d 100644
--- a/python/dune/perftool/pdelab/driver.py
+++ b/python/dune/perftool/pdelab/driver.py
@@ -1471,6 +1471,46 @@ def evaluate_residual_timer():
     return evaluation
 
 
+@preamble
+def apply_jacobian_timer():
+    # Set the matrix_free option to True!
+    from dune.perftool.options import set_option
+    set_option("matrix_free", True)
+
+    formdata = _driver_data['formdata']
+    n_go = name_gridoperator(formdata)
+    v = name_vector(formdata)
+    t_v = type_vector(formdata)
+    setup_timer()
+
+    if get_option('instrumentation_level') >= 2:
+        # Write back times
+        from dune.perftool.generation import post_include
+        post_include("HP_DECLARE_TIMER(apply_jacobian);", filetag="driver")
+        timestream = name_timing_stream()
+        print_times = []
+
+    from dune.perftool.generation import get_global_context_value
+    formdatas = get_global_context_value("formdatas")
+    for formdata in formdatas:
+        lop_name = name_localoperator(formdata)
+        if get_option('instrumentation_level') >= 3:
+            print_times.append("{}.dump_timers({}, argv[0], true);".format(lop_name, timestream))
+
+    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(timestream)]
+        evaluation.extend(print_times)
+    else:
+        evaluation = ["{}.jacobian_apply({}, j);".format(n_go, v)]
+
+    evaluation = ["{} j({});".format(t_v, v), "j=0.0;"] + evaluation
+
+    return evaluation
+
+
 @preamble
 def assemble_matrix_timer():
     formdata = _driver_data['formdata']
@@ -1743,6 +1783,7 @@ def generate_driver(formdatas, data):
         # In case of operator conunting we only assemble the matrix and evaluate the residual
         #assemble_matrix_timer()
         evaluate_residual_timer()
+        apply_jacobian_timer()
     elif is_stationary():
         # We could also use solve if we are not interested in visualization
         vtkoutput()
diff --git a/python/dune/perftool/pdelab/signatures.py b/python/dune/perftool/pdelab/signatures.py
index d2c48409d746278ea35d9abb22d99a9bad6675e8..97560a248a4aa28f3db489859c64dcb48d54c58c 100644
--- a/python/dune/perftool/pdelab/signatures.py
+++ b/python/dune/perftool/pdelab/signatures.py
@@ -274,7 +274,7 @@ def jacobian_apply_skeleton_templates():
     lfsvt = type_testfunctionspace()
     cct = type_coefficientcontainer()
     avt = type_accumulation_variable()
-    return (geot, lfsut, cct, lfsvt, lfsut, cct, lfsvt, avt)
+    return (geot, lfsut, cct, lfsvt, lfsut, cct, lfsvt, avt, avt)
 
 
 def jacobian_apply_skeleton_args():
diff --git a/python/dune/perftool/sumfact/sumfact.py b/python/dune/perftool/sumfact/sumfact.py
index 216ce95c0f7a06ab99cedd7210f01b740b69999b..2420747f3d7ceb323ce71ac0dfab61ca21e80c36 100644
--- a/python/dune/perftool/sumfact/sumfact.py
+++ b/python/dune/perftool/sumfact/sumfact.py
@@ -269,7 +269,7 @@ def generate_accumulation_instruction(visitor, accterm, measure, subdomain_id):
         # variable.
         if get_option('fastdg'):
             ft = get_global_context_value("form_type")
-            if ft == 'residual':
+            if ft == 'residual' or ft == 'jacobian_apply':
                 accum = accum + ".data()"
                 size = basis_functions_per_direction() ** world_dimension()
                 globalarg(accum, dtype=np.float64, shape=(size,), managed=False)