diff --git a/cmake/modules/DunePerftoolMacros.cmake b/cmake/modules/DunePerftoolMacros.cmake
index c2bda2cfcfdf3ae0f7cb77de06a6f6074e1e23c9..c9b67421f9b111e76ed82a25044ada1f7e8a82d6 100644
--- a/cmake/modules/DunePerftoolMacros.cmake
+++ b/cmake/modules/DunePerftoolMacros.cmake
@@ -8,17 +8,18 @@
 #
 #       The UFL file to create the executable from.
 #
-#    .. cmake_param:: TARGET
+#    .. cmake_param:: INIFILE
 #       :single:
 #       :required:
 #
-#       The name given to the added executable target.
+#       The ini file that controls the form compilation process.
+#       It is expected to contain a [formcompiler] section
 #
-#    .. cmake_param:: OPERATOR
+#    .. cmake_param:: TARGET
 #       :single:
+#       :required:
 #
-#       The local operator file name to generate. Defaults
-#       to a suitably mangled, but not easily readable name.
+#       The name given to the added executable target.
 #
 #    .. cmake_param:: SOURCE
 #
@@ -50,6 +51,7 @@ add_custom_target(generation)
 # to have correct retriggers of generated executables
 if(CMAKE_PROJECT_NAME STREQUAL dune-perftool)
   set(UFL2PDELAB_GLOB_PATTERN "${CMAKE_SOURCE_DIR}/python/*.py")
+  set(perftool_path ${CMAKE_SOURCE_DIR}/cmake/modules)
 else()
   dune_module_path(MODULE dune-perftool
                    RESULT perftool_path
@@ -64,7 +66,7 @@ file(GLOB_RECURSE UFL2PDELAB_SOURCES ${UFL2PDELAB_GLOB_PATTERN})
 
 function(add_generated_executable)
   set(OPTIONS)
-  set(SINGLE TARGET OPERATOR SOURCE UFLFILE)
+  set(SINGLE TARGET SOURCE UFLFILE INIFILE)
   set(MULTI FORM_COMPILER_ARGS DEPENDS)
   include(CMakeParseArguments)
   cmake_parse_arguments(GEN "${OPTIONS}" "${SINGLE}" "${MULTI}" ${ARGN})
@@ -83,9 +85,8 @@ function(add_generated_executable)
   if(NOT IS_ABSOLUTE GEN_UFLFILE)
     set(GEN_UFLFILE ${CMAKE_CURRENT_SOURCE_DIR}/${GEN_UFLFILE})
   endif()
-  if(NOT GEN_OPERATOR)
-    set(GEN_OPERATOR ${GEN_TARGET}_operator.hh)
-    set(GEN_OPERATOR ${CMAKE_CURRENT_BINARY_DIR}/${GEN_OPERATOR})
+  if(EXISTS ${CMAKE_CURRENT_SOURCE_DIR}/${GEN_INIFILE})
+    set(GEN_INIFILE ${CMAKE_CURRENT_SOURCE_DIR}/${GEN_INIFILE})
   endif()
   if(NOT GEN_SOURCE)
     set(GEN_DRIVER ${GEN_TARGET}_driver.hh)
@@ -97,17 +98,23 @@ function(add_generated_executable)
     configure_file(${perftool_path}/StandardMain.cmake ${GEN_SOURCE})
   endif()
 
-  add_custom_command(OUTPUT ${GEN_OPERATOR} ${GEN_DRIVER}
+  # Generate a list of generated headers to teach CMake about
+  dune_execute_process(COMMAND ${CMAKE_BINARY_DIR}/run-in-dune-env python ${perftool_path}/deplist.py ${GEN_INIFILE} ${GEN_TARGET}
+                       OUTPUT_VARIABLE header_deps
+                       )
+
+  add_custom_command(OUTPUT ${GEN_DRIVER} ${header_deps}
                      COMMAND ${CMAKE_BINARY_DIR}/run-in-dune-env ufl2pdelab
                              --project-basedir ${CMAKE_BINARY_DIR}
-                             --operator-file ${GEN_OPERATOR}
                              ${GEN_FORM_COMPILER_ARGS}
                              --uflfile ${GEN_UFLFILE}
+                             --ini-file ${GEN_INIFILE}
+                             --target-name ${GEN_TARGET}
                      DEPENDS ${GEN_UFLFILE} ${UFL2PDELAB_SOURCES} ${GEN_DEPENDS} ${DUNE_PERFTOOL_ADDITIONAL_PYTHON_SOURCES}
                      COMMENT "Running ufl2pdelab for the target ${GEN_TARGET}"
                     )
 
-  add_executable(${GEN_TARGET} ${GEN_SOURCE} ${GEN_OPERATOR})
+  add_executable(${GEN_TARGET} ${GEN_SOURCE} ${GEN_DRIVER} ${header_deps})
   target_include_directories(${GEN_TARGET} PUBLIC ${CMAKE_CURRENT_BINARY_DIR})
   add_dependencies(generation ${GEN_TARGET})
 endfunction()
diff --git a/cmake/modules/GeneratedSystemtests.cmake b/cmake/modules/GeneratedSystemtests.cmake
index d42b8aebf18ec27890e5bbdaf995a642e642d81c..4badbae9eff808b718cfd9b7a35b6202652c6416 100644
--- a/cmake/modules/GeneratedSystemtests.cmake
+++ b/cmake/modules/GeneratedSystemtests.cmake
@@ -49,7 +49,7 @@ function(dune_add_formcompiler_system_test)
 
     add_generated_executable(TARGET ${tname}
                              UFLFILE ${SYSTEMTEST_UFLFILE}
-                             FORM_COMPILER_ARGS --ini-file ${inifile}
+                             INIFILE "${CMAKE_CURRENT_BINARY_DIR}/${inifile}"
                              DEPENDS ${SYSTEMTEST_INIFILE}
                              )
 
diff --git a/python/dune/perftool/blockstructured/geometry.py b/python/dune/perftool/blockstructured/geometry.py
index 951ada695f3c8d29d78cbf8992d2fe7172470083..b933c06ce6e79e717ab3116f39f30d4a5971394f 100644
--- a/python/dune/perftool/blockstructured/geometry.py
+++ b/python/dune/perftool/blockstructured/geometry.py
@@ -2,7 +2,7 @@ from dune.perftool.generation import (get_backend,
                                       temporary_variable,
                                       instruction)
 from dune.perftool.tools import get_pymbolic_basename
-from dune.perftool.options import (get_option,
+from dune.perftool.options import (get_form_option,
                                    option_switch)
 from dune.perftool.pdelab.geometry import (name_jacobian_determinant,
                                            local_dimension,
@@ -28,7 +28,7 @@ def pymbolic_jacobian_inverse_transposed(i, j, restriction):
 # scale determinant according to the order of the blockstructure
 def pymbolic_facet_jacobian_determinant():
     return prim.Quotient(prim.Variable(name_facet_jacobian_determinant()),
-                         prim.Power(get_option("number_of_blocks"), local_dimension()))
+                         prim.Power(get_form_option("number_of_blocks"), local_dimension()))
 
 
 # translate a point in the micro element into macro coordinates
@@ -45,7 +45,7 @@ def define_point_in_macro(name, point_in_micro):
         else:
             expr = prim.Subscript(point_in_micro, (i,))
         expr = prim.Sum((expr, prim.Variable(subelem_inames[i]),))
-        expr = prim.Quotient(expr, get_option('number_of_blocks'))
+        expr = prim.Quotient(expr, get_form_option('number_of_blocks'))
         instruction(assignee=prim.Subscript(prim.Variable(name), (i,)),
                     expression=expr,
                     within_inames=frozenset(subelem_inames),
diff --git a/python/dune/perftool/blockstructured/tools.py b/python/dune/perftool/blockstructured/tools.py
index 5fa0d32d2369de9c133c9022434f7f1a9d55c070..31f57e70f7bfbdafc061a0b248e6614a7ee6480c 100644
--- a/python/dune/perftool/blockstructured/tools.py
+++ b/python/dune/perftool/blockstructured/tools.py
@@ -11,7 +11,7 @@ from dune.perftool.pdelab.geometry import (local_dimension,
 
 from dune.perftool.pdelab.quadrature import quadrature_inames
 from dune.perftool.generation.counter import get_counted_variable
-from dune.perftool.options import get_option
+from dune.perftool.options import get_form_option
 import pymbolic.primitives as prim
 
 
@@ -25,7 +25,7 @@ def sub_element_inames():
     for i in range(dim):
         name_counted = get_counted_variable(name)
         inames = inames + (name_counted,)
-        domain(name_counted, get_option("number_of_blocks"))
+        domain(name_counted, get_form_option("number_of_blocks"))
     return inames
 
 
@@ -59,7 +59,7 @@ def sub_facet_inames():
                     predicates=frozenset([prim.LogicalNot(predicate(index))])
                     )
 
-    k = get_option("number_of_blocks")
+    k = get_form_option("number_of_blocks")
 
     inames = ("x",)
     temporary_variable(inames[0])
@@ -113,7 +113,7 @@ def micro_index_to_macro_index(element, iname):
     elif it == "exterior_facet" or it == "interior_facet":
         subelem_inames = sub_facet_inames()
 
-    k = get_option("number_of_blocks")
+    k = get_form_option("number_of_blocks")
     p = element.degree()
     modified_index = prim.Sum((tensor_index_to_sequential_index(to_tensor_index(iname, p + 1), p * k + 1),
                                tensor_index_to_sequential_index(tuple(prim.Variable(iname) * p for iname in subelem_inames), p * k + 1)))
diff --git a/python/dune/perftool/compile.py b/python/dune/perftool/compile.py
index 8a95f17a3e8f2a22fd99aadb2b1c28d1e4e6bec0..8886a4d85b1f590f739082244cfcbb55b355654e 100644
--- a/python/dune/perftool/compile.py
+++ b/python/dune/perftool/compile.py
@@ -16,16 +16,14 @@ from ufl.algorithms.formfiles import interpret_ufl_namespace
 from dune.perftool.generation import (delete_cache_items,
                                       global_context,
                                       )
-from dune.perftool.options import (get_option,
+from dune.perftool.options import (get_form_option,
+                                   get_option,
                                    initialize_options,
                                    )
-from dune.perftool.pdelab.driver import (generate_driver,
-                                         set_driver_data,
-                                         )
-from dune.perftool.pdelab.localoperator import (generate_localoperator_basefile,
-                                                generate_localoperator_file,
+from dune.perftool.pdelab.driver import generate_driver
+from dune.perftool.pdelab.localoperator import (generate_localoperator_file,
                                                 generate_localoperator_kernels,
-                                                name_localoperator_file)
+                                                )
 from dune.perftool.ufl.preprocess import preprocess_form
 
 from os.path import splitext, basename, join, dirname, abspath
@@ -102,46 +100,27 @@ def read_ufl(uflfile):
     for name in magic_names:
         data.object_by_name[name] = namespace.get(name, None)
 
-    formdatas = []
-    forms = data.forms
-    for index, form in enumerate(forms):
-        formdatas.append(preprocess_form(form))
-        forms[index] = formdatas[index].preprocessed_form
-
-    # We expect at least one form
-    assert len(data.forms) >= 1
-
-    return formdatas, data
+    return data
 
 
 # This function is the entrypoint of the ufl2pdelab executable
 def compile_form():
     initialize_options()
-    formdatas, data = read_ufl(get_option("uflfile"))
-
-    with global_context(data=data, formdatas=formdatas):
-        # The driver module uses a global dictionary for storing necessary data
-        set_driver_data(formdatas, data)
+    data = read_ufl(get_option("uflfile"))
 
+    with global_context(data=data):
         # Generate driver file
         if get_option("driver_file"):
-            generate_driver(formdatas, data)
-
-    # In case of multiple forms: Genarate one file that includes all localoperator files
-    if len(formdatas) > 1:
-        generate_localoperator_basefile(formdatas, data)
-
-    # Generate local operator files
-    for formdata in formdatas:
-        with global_context(data=data, formdata=formdata):
-            # Make sure cache is empty
-            delete_cache_items()
-
-            # Create localoperator kernels
-            if get_option("operator_file"):
-                kernels = generate_localoperator_kernels(formdata, data)
-
-            # Create c++ file from kernels
-            if get_option("operator_file"):
-                filename = name_localoperator_file(formdata, data)
-                generate_localoperator_file(formdata, kernels, filename)
+            generate_driver()
+
+        for operator in get_option("operators").split(","):
+            with global_context(form_identifier=operator):
+                # Make sure cache is empty
+                delete_cache_items()
+
+                # Choose the form from the UFL input
+                kernels = generate_localoperator_kernels(operator)
+
+                # Write the result to a file
+                filename = get_form_option("filename")
+                generate_localoperator_file(kernels, filename)
diff --git a/python/dune/perftool/loopy/transformations/disjointgroups.py b/python/dune/perftool/loopy/transformations/disjointgroups.py
index 78ea4c9d6574f294416c5bb0f6c440083ddd3603..791a7b35ea26cf96af224cd5ee0623cc2d7f673a 100644
--- a/python/dune/perftool/loopy/transformations/disjointgroups.py
+++ b/python/dune/perftool/loopy/transformations/disjointgroups.py
@@ -1,12 +1,12 @@
 """ A helper transformation that makes all groups conflicting """
 
-from dune.perftool.options import get_option
+from dune.perftool.options import get_form_option
 
 
 def make_groups_conflicting(knl):
     # As this transformation introduces a performance bug that basically
     # kills our CI, we only apply it if really needed - meaning in production.
-    if get_option("assure_statement_ordering"):
+    if get_form_option("assure_statement_ordering"):
         groups = frozenset().union(*tuple(i.groups for i in knl.instructions))
         return knl.copy(instructions=[i.copy(conflicts_with_groups=groups - i.groups) for i in knl.instructions])
     else:
diff --git a/python/dune/perftool/options.py b/python/dune/perftool/options.py
index 8ed41f795bc33ceec6c0a89a9cc96772416202bd..aeb6baef9fcaf292c48326993ba85006f72f1626 100644
--- a/python/dune/perftool/options.py
+++ b/python/dune/perftool/options.py
@@ -50,6 +50,7 @@ class PerftoolGlobalOptionsArray(ImmutableRecord):
     precision_bits = PerftoolOption(default=64, helpstr="The number of bits for the floating point type")
     overlapping = PerftoolOption(default=False, helpstr="Use an overlapping solver and constraints. You still need to make sure to construct a grid with overlap! The parallel option will be set automatically.")
     operators = PerftoolOption(default="operator", helpstr="A comma separated list of operators, each name will be interpreted as a subsection name within the formcompiler section")
+    target_name = PerftoolOption(default=None, helpstr="The target name from CMake")
 
     # Arguments that are mainly to be set by logic depending on other options
     max_vector_width = PerftoolOption(default=256, helpstr=None)
@@ -67,6 +68,9 @@ class PerftoolFormOptionsArray(ImmutableRecord):
         ImmutableRecord.__init__(self, **opts)
 
     # Form specific options
+    form = PerftoolOption(default="r", helpstr="The name of the UFL object representing the form in the UFL file")
+    filename = PerftoolOption(default=None, helpstr="The filename to use for this LocalOperator")
+    classname = PerftoolOption(default=None, helpstr="The name of the C++ class to generate")
     numerical_jacobian = PerftoolOption(default=False, helpstr="use numerical jacobians (only makes sense, if uflpdelab for some reason fails to generate analytic jacobians)")
     matrix_free = PerftoolOption(default=False, helpstr="Use iterative solver with matrix free jacobian application")
     print_transformations = PerftoolOption(default=False, helpstr="print out dot files after ufl tree transformations")
@@ -160,6 +164,12 @@ def process_form_options(opt):
     if opt.numerical_jacobian:
         opt = opt.copy(generate_jacobians=False)
 
+    if opt.classname is None:
+        opt = opt.copy(classname="{}Operator".format(opt.form))
+
+    if opt.filename is None:
+        opt = opt.copy(filename="{}_{}_file.hh".format(get_option("target_name"), opt.classname))
+
     return opt
 
 
@@ -187,17 +197,29 @@ def set_option(key, value):
 
 def get_option(key):
     processed_global_opts = process_global_options(_global_options)
-    if hasattr(processed_global_opts, key):
-        return getattr(processed_global_opts, key)
-    else:
-        processed_form_opts = process_form_options(_form_options["operator"])
-        return getattr(processed_form_opts, key)
+    return getattr(processed_global_opts, key)
+
+
+def get_form_option(key, form=None):
+    if form is None:
+        from dune.perftool.generation import get_global_context_value
+        form = get_global_context_value("form_identifier", 0)
+    if isinstance(form, int):
+        form = get_option("operators").split(",")[form]
+    processed_form_opts = process_form_options(_form_options[form])
+    return getattr(processed_form_opts, key)
 
 
 def option_switch(opt):
     def _switch():
-        if get_option(opt):
-            return opt
-        else:
-            return "default"
+        try:
+            if get_option(opt):
+                return opt
+            else:
+                return "default"
+        except AttributeError:
+            if get_form_option(opt):
+                return opt
+            else:
+                return "default"
     return _switch
diff --git a/python/dune/perftool/pdelab/argument.py b/python/dune/perftool/pdelab/argument.py
index 848d9d9ad275b4ff910da9e9f6514fe883f9b436..e4dfb972a48b34428b54bef25f142eafd2c86866 100644
--- a/python/dune/perftool/pdelab/argument.py
+++ b/python/dune/perftool/pdelab/argument.py
@@ -5,7 +5,6 @@ Namely:
 * accumulation object (r, jac...)
 """
 
-from dune.perftool.options import get_option
 from dune.perftool.generation import (domain,
                                       function_mangler,
                                       iname,
diff --git a/python/dune/perftool/pdelab/basis.py b/python/dune/perftool/pdelab/basis.py
index 0c7ab18a63f12b9e870992aad88143d65229ef67..53917205e71c983e7036cffaf2d9c98aec30dc5a 100644
--- a/python/dune/perftool/pdelab/basis.py
+++ b/python/dune/perftool/pdelab/basis.py
@@ -9,7 +9,7 @@ from dune.perftool.generation import (backend,
                                       temporary_variable,
                                       )
 from dune.perftool.options import (option_switch,
-                                   get_option
+                                   get_form_option,
                                    )
 from dune.perftool.pdelab.spaces import (lfs_iname,
                                          lfs_inames,
@@ -172,7 +172,7 @@ def evaluate_coefficient(visitor, element, name, container, restriction, index):
     basis = visitor.interface.pymbolic_basis(sub_element, restriction, 0, context='trial')
     basisindex, _ = get_pymbolic_indices(basis)
 
-    if get_option("blockstructured"):
+    if get_form_option("blockstructured"):
         from dune.perftool.blockstructured.argument import pymbolic_coefficient
         coeff = pymbolic_coefficient(container, lfs, sub_element, basisindex)
     else:
@@ -211,7 +211,7 @@ def evaluate_coefficient_gradient(visitor, element, name, container, restriction
     from dune.perftool.tools import maybe_wrap_subscript
     basis = maybe_wrap_subscript(basis, Variable(dimindex))
 
-    if get_option("blockstructured"):
+    if get_form_option("blockstructured"):
         from dune.perftool.blockstructured.argument import pymbolic_coefficient
         coeff = pymbolic_coefficient(container, lfs, sub_element, basisindex)
     else:
diff --git a/python/dune/perftool/pdelab/driver/__init__.py b/python/dune/perftool/pdelab/driver/__init__.py
index e47446537d95a57b428cc832a5abddef6aac36f6..f11875d043302a9f94ded59ec06514ba4b61a294 100644
--- a/python/dune/perftool/pdelab/driver/__init__.py
+++ b/python/dune/perftool/pdelab/driver/__init__.py
@@ -19,84 +19,61 @@ from dune.perftool.generation import (generator_factory,
                                       pre_include,
                                       preamble,
                                       )
-from dune.perftool.options import get_option
+from dune.perftool.options import (get_form_option,
+                                   get_option,
+                                   )
 
-# Have a global variable with the entire form data. This allows functions that depend
-# deterministically on the entire data set to directly access it instead of passing it
-# through the entire generator chain.
-_driver_data = {}
+#
+# The following functions are not doing anything useful, but providing easy access
+# to quantities that are needed throughout the process of generating the driver! 
+#
 
+def get_form_ident(index=0):
+    return get_option("operators").split(",")[index]
 
-# Have a function access this global data structure
-def set_driver_data(formdatas, data):
-    assert (len(formdatas) <= 2)
-    if len(formdatas) == 1:
-        _driver_data['form'] = formdatas[0].preprocessed_form
-        _driver_data['formdata'] = formdatas[0]
-    else:
-        mass_index = mass_form_index(formdatas, data)
-        if mass_index is None:
-            raise NotImplementedError("Form for mass matrix needs to have name 'mass' in ufl file.")
-        _driver_data['mass_form'] = formdatas[mass_index].preprocessed_form
-        _driver_data['mass_formdata'] = formdatas[mass_index]
-        _driver_data['form'] = formdatas[1 - mass_index].preprocessed_form
-        _driver_data['formdata'] = formdatas[1 - mass_index]
 
-    _driver_data['data'] = data
+def get_form(what=None):
+    """ Return the ith form specified """
+    if what is None:
+        what = get_global_context_value("form_identifier", 0)
+    if isinstance(what, int):
+        what = get_form_ident(what)
+    data = get_global_context_value("data")
+    return data.object_by_name[get_form_option("form", what)]
+
+
+def get_preprocessed_form(what=None):
+    from dune.perftool.ufl.preprocess import preprocess_form
+    return preprocess_form(get_form(what)).preprocessed_form
 
 
 def get_dimension():
-    return _driver_data['form'].ufl_cell().geometric_dimension()
+    return get_form().ufl_cell().geometric_dimension()
 
 
 def get_cell():
-    return _driver_data['form'].ufl_cell().cellname()
+    return get_form().ufl_cell().cellname()
 
 
 def get_test_element():
-    return _driver_data['form'].arguments()[0].ufl_element()
+    return get_form().arguments()[0].ufl_element()
 
 
 def get_trial_element():
-    return _driver_data['form'].coefficients()[0].ufl_element()
-
-
-def get_formdata():
-    return _driver_data['formdata']
-
-
-def get_mass_formdata():
-    return _driver_data["mass_formdata"]
+    return get_form().coefficients()[0].ufl_element()
 
 
 def is_stationary():
-    return 'mass_form' not in _driver_data
-
-
-def form_name_suffix(name, formdata):
-    from dune.perftool.pdelab.localoperator import name_form
-    data = get_global_context_value('data')
-    form_name = name_form(formdata, data)
-    return name + '_' + form_name
-
-
-def get_object(name):
-    return _driver_data['data'].object_by_name.get(name, None)
-
-
-def mass_form_index(formdatas, data):
-    for index, formdata in enumerate(formdatas):
-        try:
-            if data.object_names[id(formdata.original_form)] == 'mass':
-                return index
-        except KeyError:
-            continue
+    # TODO I am completely unsure how this should work in the future
+    # This only fixes instationary stuff, it will break Renes adjoint stuff
+    return len(get_option("operators").split(",")) == 1
+#     return 'mass_form' not in _driver_data
 
 
 def is_linear(form=None):
     '''Test if form is linear in trial function'''
     if form is None:
-        form = get_formdata().original_form
+        form = get_form()
     from ufl import derivative
     from ufl.algorithms import expand_derivatives
     jacform = expand_derivatives(derivative(form, form.coefficients()[0]))
@@ -194,7 +171,7 @@ def unroll_list_tensors(data):
 
 
 def preprocess_leaf_data(element, data, applyZeroDefault=True):
-    data = get_object(data)
+    data = get_global_context_value("data").object_by_name.get(data, None)
     if data is None and not applyZeroDefault:
         return None
 
@@ -262,7 +239,7 @@ def check_parallel_execution():
             "}"]
 
 
-def generate_driver(formdatas, data):
+def generate_driver():
     # Add check to c++ file if this program should only be used in parallel mode
     if get_option("parallel"):
         check_parallel_execution()
@@ -273,7 +250,7 @@ def generate_driver(formdatas, data):
             assert(not get_option("opcounter"))
         assert(any(_driver_data['form'].ufl_cell().cellname() in x for x in
                    ["vertex", "interval", "quadrilateral", "hexahedron"]))
-        # In case of operator conunting we only assemble the matrix and evaluate the residual
+        # In case of operator counting we only assemble the matrix and evaluate the residual
         # assemble_matrix_timer()
         from dune.perftool.pdelab.driver.timings import apply_jacobian_timer, evaluate_residual_timer
         from dune.perftool.loopy.target import type_floatingpoint
diff --git a/python/dune/perftool/pdelab/driver/constraints.py b/python/dune/perftool/pdelab/driver/constraints.py
index 901133c853819c45f7cae6af21ba212346d4adf2..6e9d7e7a9fa42af8db9177302f0418fb37ef6a39 100644
--- a/python/dune/perftool/pdelab/driver/constraints.py
+++ b/python/dune/perftool/pdelab/driver/constraints.py
@@ -3,7 +3,6 @@ from dune.perftool.generation import (global_context,
                                       preamble,
                                       )
 from dune.perftool.pdelab.driver import (FEM_name_mangling,
-                                         get_formdata,
                                          get_trial_element,
                                          )
 from dune.perftool.pdelab.driver.gridfunctionspace import (name_gfs,
diff --git a/python/dune/perftool/pdelab/driver/error.py b/python/dune/perftool/pdelab/driver/error.py
index b5b616627c111f8eeaef6b29e31c307fc7a40c14..b210d61fec9962f34857861309b0e458b4974b1e 100644
--- a/python/dune/perftool/pdelab/driver/error.py
+++ b/python/dune/perftool/pdelab/driver/error.py
@@ -5,7 +5,7 @@ from dune.perftool.generation import (cached,
                                       preamble,
                                       )
 from dune.perftool.options import get_option
-from dune.perftool.pdelab.driver import (get_formdata,
+from dune.perftool.pdelab.driver import (get_form_ident,
                                          get_trial_element,
                                          preprocess_leaf_data,
                                          )
@@ -65,7 +65,7 @@ def name_discrete_grid_function(gfs, vector_name):
 @preamble
 def typedef_difference_squared_adapter(name, treepath):
     sol = name_exact_solution_gridfunction(treepath)
-    vector = name_vector(get_formdata())
+    vector = name_vector(get_form_ident())
     gfs = name_trial_subgfs(treepath)
     dgf = name_discrete_grid_function(gfs, vector)
 
@@ -82,7 +82,7 @@ def type_difference_squared_adapter(treepath):
 def define_difference_squared_adapter(name, treepath):
     t = type_difference_squared_adapter(treepath)
     sol = name_exact_solution_gridfunction(treepath)
-    vector = name_vector(get_formdata())
+    vector = name_vector(get_form_ident())
     gfs = name_trial_subgfs(treepath)
     dgf = name_discrete_grid_function(gfs, vector)
 
diff --git a/python/dune/perftool/pdelab/driver/gridfunctionspace.py b/python/dune/perftool/pdelab/driver/gridfunctionspace.py
index c14a73dd649efa4517b1a069a0722535a71b30f6..410036b9e8c95b3dae16c8b8eb8aa4a64e44988a 100644
--- a/python/dune/perftool/pdelab/driver/gridfunctionspace.py
+++ b/python/dune/perftool/pdelab/driver/gridfunctionspace.py
@@ -1,7 +1,10 @@
 from dune.perftool.generation import (include_file,
                                       preamble,
                                       )
-from dune.perftool.options import get_option, set_option
+from dune.perftool.options import (get_form_option,
+                                   get_option,
+                                   set_option,
+                                   )
 from dune.perftool.pdelab.driver import (FEM_name_mangling,
                                          get_cell,
                                          get_dimension,
@@ -47,7 +50,7 @@ def type_range():
 def typedef_grid(name):
     dim = get_dimension()
     if isQuadrilateral(get_trial_element().cell()):
-        # For Yasp Grids the jacobi of the transformation is diagonal and constant on each cell
+        # For Yasp Grids the jacobian of the transformation is diagonal and constant on each cell
         set_option('diagonal_transformation_matrix', True)
         set_option('constant_transformation_matrix', True)
 
@@ -121,7 +124,7 @@ def typedef_fem(element, name):
     r = type_range()
     dim = get_dimension()
 
-    if get_option("blockstructured"):
+    if get_form_option("blockstructured"):
         include_file("dune/perftool/blockstructured/blockstructuredqkfem.hh", filetag="driver")
         degree = element.degree() * get_option("number_of_blocks")
         return "using {} = Dune::PDELab::BlockstructuredQkLocalFiniteElementMap<{}, {}, {}, {}>;" \
@@ -340,7 +343,7 @@ def typedef_composite_gfs(element, name, subgfs, root):
 @preamble
 def typedef_vectorbackend(name, element, root):
     include_file("dune/pdelab/backend/istl.hh", filetag="driver")
-    if get_option("fastdg") and root:
+    if get_form_option("fastdg") and root:
         blocking = "Dune::PDELab::ISTL::Blocking::fixed"
         if isinstance(element, MixedElement):
             blocksize = ""
diff --git a/python/dune/perftool/pdelab/driver/gridoperator.py b/python/dune/perftool/pdelab/driver/gridoperator.py
index d0244c107240fb5175329cac47ed145ffbdf02fb..64a6177774beecf3b3438f04f247af307d75542c 100644
--- a/python/dune/perftool/pdelab/driver/gridoperator.py
+++ b/python/dune/perftool/pdelab/driver/gridoperator.py
@@ -2,8 +2,7 @@ from dune.perftool.generation import (get_global_context_value,
                                       include_file,
                                       preamble,
                                       )
-from dune.perftool.pdelab.driver import (form_name_suffix,
-                                         get_cell,
+from dune.perftool.pdelab.driver import (get_cell,
                                          get_dimension,
                                          get_test_element,
                                          get_trial_element,
@@ -23,20 +22,22 @@ from dune.perftool.pdelab.driver.gridfunctionspace import (name_test_gfs,
                                                            )
 from dune.perftool.pdelab.localoperator import localoperator_basename
 from dune.perftool.pdelab.parameter import parameterclass_basename
-from dune.perftool.options import get_option
+from dune.perftool.options import (get_form_option,
+                                   get_option,
+                                   )
 
 
 @preamble
-def typedef_gridoperator(name, formdata):
+def typedef_gridoperator(name, form_ident):
     ugfs = type_trial_gfs()
     vgfs = type_test_gfs()
-    lop = type_localoperator(formdata)
+    lop = type_localoperator(form_ident)
     cc = type_constraintscontainer()
     mb = type_matrixbackend()
     df = type_domainfield()
     r = type_range()
-    if get_option("fastdg"):
-        if not get_option("sumfact"):
+    if get_form_option("fastdg"):
+        if not get_form_option("sumfact"):
             raise PerftoolCodegenError("FastDGGridOperator is only implemented for sumfactorization.")
         include_file("dune/pdelab/gridoperator/fastdg.hh", filetag="driver")
         return "using {} = Dune::PDELab::FastDGGridOperator<{}, {}, {}, {}, {}, {}, {}, {}, {}>;".format(name, ugfs, vgfs, lop, mb, df, r, r, cc, cc)
@@ -45,64 +46,63 @@ def typedef_gridoperator(name, formdata):
         return "using {} = Dune::PDELab::GridOperator<{}, {}, {}, {}, {}, {}, {}, {}, {}>;".format(name, ugfs, vgfs, lop, mb, df, r, r, cc, cc)
 
 
-def type_gridoperator(formdata):
-    name = form_name_suffix("GO", formdata).upper()
-    typedef_gridoperator(name, formdata)
+def type_gridoperator(form_ident):
+    name = "GO_{}".format(form_ident)
+    typedef_gridoperator(name, form_ident)
     return name
 
 
 @preamble
-def define_gridoperator(name, formdata):
-    gotype = type_gridoperator(formdata)
+def define_gridoperator(name, form_ident):
+    gotype = type_gridoperator(form_ident)
     ugfs = name_trial_gfs()
     vgfs = name_test_gfs()
     if ugfs != vgfs:
         raise NotImplementedError("Non-Galerkin methods currently not supported!")
     cc = name_assembled_constraints()
-    lop = name_localoperator(formdata)
+    lop = name_localoperator(form_ident)
     mb = name_matrixbackend()
     return ["{} {}({}, {}, {}, {}, {}, {});".format(gotype, name, ugfs, cc, vgfs, cc, lop, mb),
             "std::cout << \"gfs with \" << {}.size() << \" dofs generated  \"<< std::endl;".format(ugfs),
             "std::cout << \"cc with \" << {}.size() << \" dofs generated  \"<< std::endl;".format(cc)]
 
 
-def name_gridoperator(formdata):
-    name = form_name_suffix("go", formdata)
-    define_gridoperator(name, formdata)
+def name_gridoperator(form_ident):
+    name = "go_{}".format(form_ident)
+    define_gridoperator(name, form_ident)
     return name
 
 
 @preamble
-def typedef_localoperator(name, formdata):
+def typedef_localoperator(name, form_ident):
     ugfs = type_trial_gfs()
     vgfs = type_test_gfs()
-    data = get_global_context_value("data")
-    filename = get_option("operator_file")
+    filename = get_form_option("filename")
     include_file(filename, filetag="driver")
-    lopname = localoperator_basename(formdata, data)
+    lopname = localoperator_basename(form_ident)
     range_type = type_range()
     return "using {} = {}<{}, {}, {}>;".format(name, lopname, ugfs, vgfs, range_type)
 
 
-def type_localoperator(formdata):
-    name = form_name_suffix("LOP", formdata).upper()
-    typedef_localoperator(name, formdata)
+def type_localoperator(form_ident):
+    name = "LOP_{}".format(form_ident.upper())
+    typedef_localoperator(name, form_ident)
     return name
 
 
 @preamble
-def define_localoperator(name, formdata):
+def define_localoperator(name, form_ident):
     trial_gfs = name_trial_gfs()
     test_gfs = name_test_gfs()
-    loptype = type_localoperator(formdata)
+    loptype = type_localoperator(form_ident)
     ini = name_initree()
-    params = name_parameters(formdata)
+    params = name_parameters(form_ident)
     return "{} {}({}, {}, {}, {});".format(loptype, name, trial_gfs, test_gfs, ini, params)
 
 
-def name_localoperator(formdata):
-    name = form_name_suffix("lop", formdata)
-    define_localoperator(name, formdata)
+def name_localoperator(form_ident):
+    name = "lop_{}".format(form_ident)
+    define_localoperator(name, form_ident)
     return name
 
 
@@ -158,20 +158,19 @@ def name_matrixbackend():
     return name
 
 
-def type_parameters(formdata):
-    data = get_global_context_value("data")
-    name = parameterclass_basename(formdata, data)
+def type_parameters(form_ident):
+    name = parameterclass_basename(form_ident)
     return name
 
 
 @preamble
-def define_parameters(name, formdata):
-    partype = type_parameters(formdata)
+def define_parameters(name, form_ident):
+    partype = type_parameters(form_ident)
     return "{} {};".format(partype, name)
 
 
-def name_parameters(formdata, define=True):
-    name = form_name_suffix("params", formdata)
+def name_parameters(form_ident, define=True):
+    name = "params_{}".format(form_ident)
     if define:
-        define_parameters(name, formdata)
+        define_parameters(name, form_ident)
     return name
diff --git a/python/dune/perftool/pdelab/driver/instationary.py b/python/dune/perftool/pdelab/driver/instationary.py
index 3ed2635b91aeb5e3c1af777eea4cbbd2aa0f76aa..11390040db6a1056bc52a00fbf6a0c8f033b7167 100644
--- a/python/dune/perftool/pdelab/driver/instationary.py
+++ b/python/dune/perftool/pdelab/driver/instationary.py
@@ -1,9 +1,7 @@
 from dune.perftool.generation import (include_file,
                                       preamble,
                                       )
-from dune.perftool.pdelab.driver import (get_formdata,
-                                         get_mass_formdata,
-                                         get_trial_element,
+from dune.perftool.pdelab.driver import (get_trial_element,
                                          is_linear,
                                          name_initree,
                                          preprocess_leaf_data,
@@ -33,12 +31,14 @@ from dune.perftool.pdelab.driver.solve import (print_matrix,
 from dune.perftool.pdelab.driver.vtk import (name_vtk_sequence_writer,
                                              visualize_initial_condition,
                                              )
-from dune.perftool.options import get_option
+from dune.perftool.options import (get_form_option,
+                                   get_option,
+                                   )
 
 
 def solve_instationary():
     # Create time loop
-    if get_option('matrix_free'):
+    if get_form_option('matrix_free'):
         raise NotImplementedError("Instationary matrix free not implemented!")
     else:
         time_loop()
diff --git a/python/dune/perftool/pdelab/driver/interpolate.py b/python/dune/perftool/pdelab/driver/interpolate.py
index af2edb556219bbd38fc8e968caa6b95cdda4cae9..73d44a560b16882b024016c1414e38d35cbe1fbc 100644
--- a/python/dune/perftool/pdelab/driver/interpolate.py
+++ b/python/dune/perftool/pdelab/driver/interpolate.py
@@ -5,7 +5,6 @@ from dune.perftool.generation import (cached,
                                       preamble,
                                       )
 from dune.perftool.pdelab.driver import (FEM_name_mangling,
-                                         get_formdata,
                                          get_trial_element,
                                          is_stationary,
                                          preprocess_leaf_data,
diff --git a/python/dune/perftool/pdelab/driver/solve.py b/python/dune/perftool/pdelab/driver/solve.py
index 9b4cad9154f0de3a968c19a500905f849f23628b..44b4cd7e3f87440231399397a7089f6a8a5d875f 100644
--- a/python/dune/perftool/pdelab/driver/solve.py
+++ b/python/dune/perftool/pdelab/driver/solve.py
@@ -1,9 +1,10 @@
 from dune.perftool.generation import (include_file,
                                       preamble,
                                       )
-from dune.perftool.options import get_option
-from dune.perftool.pdelab.driver import (form_name_suffix,
-                                         get_formdata,
+from dune.perftool.options import (get_form_option,
+                                   get_option,
+                                   )
+from dune.perftool.pdelab.driver import (get_form_ident,
                                          is_linear,
                                          name_initree,
                                          )
@@ -21,17 +22,16 @@ from dune.perftool.pdelab.driver.interpolate import interpolate_dirichlet_data
 
 @preamble
 def dune_solve():
+    form_ident = get_form_ident()
     # Test if form is linear in ansatzfunction
     linear = is_linear()
 
     # Test wether we want to do matrix free operator evaluation
-    matrix_free = get_option('matrix_free')
-
+    matrix_free = get_form_option('matrix_free')
     # Get right solve command
     if linear and matrix_free:
-        formdata = get_formdata()
-        go = name_gridoperator(formdata)
-        x = name_vector(formdata)
+        go = name_gridoperator(form_ident)
+        x = name_vector(form_ident)
         include_file("dune/perftool/matrixfree.hh", filetag="driver")
         solve = "solveMatrixFree({},{});".format(go, x)
     elif linear and not matrix_free:
@@ -39,14 +39,13 @@ def dune_solve():
         solve = "{}.apply();".format(slp)
     elif not linear and matrix_free:
         # TODO copy of linear case and obviously broken, used to generate something ;)
-        formdata = get_formdata()
-        go = name_gridoperator(formdata)
-        x = name_vector(formdata)
+        go = name_gridoperator(form_ident)
+        x = name_vector(form_ident)
         include_file("dune/perftool/matrixfree.hh", filetag="driver")
         solve = "solveNonlinearMatrixFree({},{});".format(go, x)
     elif not linear and not matrix_free:
-        go_type = type_gridoperator(get_formdata())
-        go = name_gridoperator(get_formdata())
+        go_type = type_gridoperator(form_ident)
+        go = name_gridoperator(form_ident)
         snp = name_stationarynonlinearproblemsolver(go_type, go)
         solve = "{}.apply();".format(snp)
 
@@ -59,49 +58,43 @@ def dune_solve():
         from dune.perftool.generation import post_include
         post_include("HP_DECLARE_TIMER(solve);", filetag="driver")
 
-        # Print times after solving
-        from dune.perftool.generation import get_global_context_value
-        formdatas = get_global_context_value("formdatas")
-        print_times = []
-        for formdata in formdatas:
-            from dune.perftool.pdelab.driver.gridoperator import name_localoperator
-            lop_name = name_localoperator(formdata)
-            timestream = name_timing_stream()
-            print_times.append("{}.dump_timers({}, {}, true);".format(lop_name, timestream, name_timing_identifier()))
-
         solve = ["HP_TIMER_START(solve);",
                  "{}".format(solve),
                  "HP_TIMER_STOP(solve);",
                  "DUMP_TIMER({}, solve, {}, true);".format(get_option("instrumentation_level"), timestream),
                  ]
+
         if get_option('instrumentation_level') >= 3:
-            solve.extend(print_times)
+            from dune.perftool.pdelab.driver.gridoperator import name_localoperator
+            lop_name = name_localoperator(form_ident)
+            timestream = name_timing_stream()
+            solve.append("{}.dump_timers({}, {}, true);".format(lop_name, timestream, name_timing_identifier()))
 
     return solve
 
 
-def name_vector(formdata):
-    name = form_name_suffix("x", formdata)
-    define_vector(name, formdata)
+def name_vector(form_ident):
+    name = "x_{}".format(form_ident)
+    define_vector(name, form_ident)
     interpolate_dirichlet_data(name)
     return name
 
 
 @preamble
-def typedef_vector(name, formdata):
-    go_type = type_gridoperator(formdata)
+def typedef_vector(name, form_ident):
+    go_type = type_gridoperator(form_ident)
     return "using {} = {}::Traits::Domain;".format(name, go_type)
 
 
-def type_vector(formdata):
-    name = form_name_suffix("V", formdata).upper()
-    typedef_vector(name, formdata)
+def type_vector(form_ident):
+    name = "V_{}".format(form_ident.upper())
+    typedef_vector(name, form_ident)
     return name
 
 
 @preamble
-def define_vector(name, formdata):
-    vtype = type_vector(formdata)
+def define_vector(name, form_ident):
+    vtype = type_vector(form_ident)
     gfs = name_trial_gfs()
     return ["{} {}({});".format(vtype, name, gfs), "{} = 0.0;".format(name)]
 
@@ -155,9 +148,9 @@ def name_reduction():
 @preamble
 def typedef_stationarylinearproblemsolver(name):
     include_file("dune/pdelab/stationary/linearproblem.hh", filetag="driver")
-    gotype = type_gridoperator(get_formdata())
+    gotype = type_gridoperator(get_form_ident())
     lstype = type_linearsolver()
-    xtype = type_vector(get_formdata())
+    xtype = type_vector(get_form_ident())
     return "using {} = Dune::PDELab::StationaryLinearProblemSolver<{}, {}, {}>;".format(name, gotype, lstype, xtype)
 
 
@@ -169,10 +162,9 @@ def type_stationarylinearproblemsolver():
 @preamble
 def define_stationarylinearproblemsolver(name):
     slptype = type_stationarylinearproblemsolver()
-    formdata = get_formdata()
-    go = name_gridoperator(formdata)
+    go = name_gridoperator(get_form_ident())
     ls = name_linearsolver()
-    x = name_vector(formdata)
+    x = name_vector(get_form_ident())
     red = name_reduction()
     return "{} {}({}, {}, {}, {});".format(slptype, name, go, ls, x, red)
 
@@ -186,7 +178,7 @@ def name_stationarylinearproblemsolver():
 def typedef_stationarynonlinearproblemsolver(name, go_type):
     include_file("dune/pdelab/newton/newton.hh", filetag="driver")
     ls_type = type_linearsolver()
-    x_type = type_vector(get_formdata())
+    x_type = type_vector(form_ident)
     return "using {} = Dune::PDELab::Newton<{}, {}, {}>;".format(name, go_type, ls_type, x_type)
 
 
@@ -199,7 +191,7 @@ def type_stationarynonlinearproblemssolver(go_type):
 @preamble
 def define_stationarynonlinearproblemsolver(name, go_type, go):
     snptype = type_stationarynonlinearproblemssolver(go_type)
-    x = name_vector(get_formdata())
+    x = name_vector(get_form_ident())
     ls = name_linearsolver()
     return "{} {}({}, {}, {});".format(snptype, name, go, x, ls)
 
@@ -213,10 +205,9 @@ def name_stationarynonlinearproblemsolver(go_type, go):
 @preamble
 def print_residual():
     ini = name_initree()
-    formdata = get_formdata()
-    n_go = name_gridoperator(formdata)
-    v = name_vector(formdata)
-    t_v = type_vector(formdata)
+    n_go = name_gridoperator(get_form_ident())
+    v = name_vector(get_form_ident())
+    t_v = type_vector(get_form_ident())
     include_file("random", system=True, filetag="driver")
 
     return ["if ({}.get<bool>(\"printresidual\", false)) {{".format(ini),
@@ -236,12 +227,11 @@ def print_residual():
 
 @preamble
 def print_matrix():
-    formdata = get_formdata()
     ini = name_initree()
-    t_go = type_gridoperator(formdata)
-    n_go = name_gridoperator(formdata)
-    v = name_vector(formdata)
-    t_v = type_vector(formdata)
+    t_go = type_gridoperator(get_form_ident())
+    n_go = name_gridoperator(get_form_ident())
+    v = name_vector(get_form_ident())
+    t_v = type_vector(get_form_ident())
 
     return ["if ({}.get<bool>(\"printmatrix\", false)) {{".format(ini),
             "  // Setup random input",
diff --git a/python/dune/perftool/pdelab/driver/visitor.py b/python/dune/perftool/pdelab/driver/visitor.py
index 8ef37f0d8975fee436fc91c7092a62f5a27dc84f..83cece58280c6f8e22404cebed0c74831687477d 100644
--- a/python/dune/perftool/pdelab/driver/visitor.py
+++ b/python/dune/perftool/pdelab/driver/visitor.py
@@ -32,8 +32,7 @@ class DriverUFL2PymbolicVisitor(UFL2LoopyVisitor):
 def ufl_to_code(expr, boundary=True):
     # So far, we only considered this code branch on boundaries!
     assert boundary
-    from dune.perftool.pdelab.driver import get_formdata
-    with global_context(integral_type="exterior_facet", formdata=get_formdata()):
+    with global_context(integral_type="exterior_facet"):
         visitor = DriverUFL2PymbolicVisitor()
         from pymbolic.mapper.c_code import CCodeMapper
         ccm = CCodeMapper()
diff --git a/python/dune/perftool/pdelab/driver/vtk.py b/python/dune/perftool/pdelab/driver/vtk.py
index 6b8b52adc22cf16ec493592dfc79fa5462be8f0b..b131040e12b5cb376acce4c1266504351b82f03f 100644
--- a/python/dune/perftool/pdelab/driver/vtk.py
+++ b/python/dune/perftool/pdelab/driver/vtk.py
@@ -1,8 +1,10 @@
 from dune.perftool.generation import (include_file,
                                       preamble,
                                       )
-from dune.perftool.options import get_option
-from dune.perftool.pdelab.driver import (get_formdata,
+from dune.perftool.options import (get_form_option,
+                                   get_option,
+                                   )
+from dune.perftool.pdelab.driver import (get_form_ident,
                                          get_trial_element,
                                          name_initree,
                                          preprocess_leaf_data,
@@ -45,7 +47,7 @@ def define_subsamplinglevel(name):
     degree = get_trial_element().degree()
     if isinstance(degree, tuple):
         degree = max(degree)
-    if get_option("blockstructured"):
+    if get_form_option("blockstructured"):
         degree *= get_option("number_of_blocks")
     return "Dune::RefinementIntervals {}({}.get<int>(\"vtk.subsamplinglevel\", {}));".format(name, ini, max(degree, 1))
 
@@ -75,7 +77,7 @@ def vtkoutput():
     gfs = name_trial_gfs()
     vtkfile = name_vtkfile()
     predicate = name_predicate()
-    vec = name_vector(get_formdata())
+    vec = name_vector(get_form_ident())
 
     return ["Dune::PDELab::addSolutionToVTKWriter({}, {}, {}, Dune::PDELab::vtk::defaultNameScheme(), {});".format(vtkwriter, gfs, vec, predicate),
             "{}.write({}, Dune::VTK::ascii);".format(vtkwriter, vtkfile)]
diff --git a/python/dune/perftool/pdelab/geometry.py b/python/dune/perftool/pdelab/geometry.py
index 182d8d6ef127967d5ae32ec6e6672a664cb8d8c3..6c400cc4db61472abef9e0d39e953eb9e1e8d1e0 100644
--- a/python/dune/perftool/pdelab/geometry.py
+++ b/python/dune/perftool/pdelab/geometry.py
@@ -13,7 +13,7 @@ from dune.perftool.generation import (backend,
                                       temporary_variable,
                                       valuearg,
                                       )
-from dune.perftool.options import (get_option,
+from dune.perftool.options import (get_form_option,
                                    option_switch,
                                    )
 from dune.perftool.loopy.target import dtype_floatingpoint, type_floatingpoint
@@ -60,8 +60,7 @@ def _component_iname(context, count):
     if context:
         context = '_' + context
     name = 'idim{}{}'.format(context, str(count))
-    formdata = get_global_context_value('formdata')
-    dim = formdata.geometric_dimension
+    dim = world_dimension()
     domain(name, dim)
     return name
 
@@ -232,8 +231,8 @@ def to_cell_coordinates(local, restriction):
 
 
 def world_dimension():
-    formdata = get_global_context_value('formdata')
-    return formdata.geometric_dimension
+    from dune.perftool.pdelab.driver import get_dimension
+    return get_dimension()
 
 
 def intersection_dimension():
@@ -266,7 +265,7 @@ def declare_normal(name, shape, shape_impl):
 
 def pymbolic_unit_outer_normal():
     name = "outer_normal"
-    if not get_option("diagonal_transformation_matrix"):
+    if not get_form_option("diagonal_transformation_matrix"):
         temporary_variable(name, shape=(world_dimension(),), decl_method=declare_normal)
         evaluate_unit_outer_normal(name)
     else:
@@ -470,7 +469,7 @@ def define_cell_volume(name, restriction):
 
 
 def pymbolic_cell_volume(restriction):
-    if get_option("constant_transformation_matrix"):
+    if get_form_option("constant_transformation_matrix"):
         return pymbolic_jacobian_determinant()
     else:
         name = restricted_name("volume", restriction)
@@ -486,7 +485,7 @@ def define_facet_area(name):
 
 
 def pymbolic_facet_area():
-    if get_option("constant_transformation_matrix"):
+    if get_form_option("constant_transformation_matrix"):
         return pymbolic_facet_jacobian_determinant()
     else:
         name = "area"
diff --git a/python/dune/perftool/pdelab/localoperator.py b/python/dune/perftool/pdelab/localoperator.py
index 76d6b0c162e2b5677905ac6c79cfefc8f9f3e01b..05c2adc38ae0904e3265d24c2cdba7373af8e9f6 100644
--- a/python/dune/perftool/pdelab/localoperator.py
+++ b/python/dune/perftool/pdelab/localoperator.py
@@ -3,7 +3,8 @@ from os.path import splitext
 
 import logging
 
-from dune.perftool.options import (get_option,
+from dune.perftool.options import (get_form_option,
+                                   get_option,
                                    option_switch)
 from dune.perftool.generation import (backend,
                                       base_class,
@@ -42,31 +43,6 @@ import loopy as lp
 import cgen
 
 
-def name_form(formdata, data):
-    # Check wether the formdata has a name in UFL
-    try:
-        name = data.object_names[id(formdata.original_form)]
-        return name
-    except:
-        for index, form in enumerate(data.forms):
-            if formdata.preprocessed_form.equals(form):
-                name = str(index)
-                return name
-    # If the form has no name and can not be found in data.forms something went wrong
-    assert False
-
-
-def name_localoperator_file(formdata, data):
-    from dune.perftool.options import get_option
-    if len(data.forms) == 1:
-        filename = get_option("operator_file")
-    else:
-        suffix = '_' + name_form(formdata, data)
-        basename, extension = splitext(get_option("operator_file"))
-        filename = basename + suffix + extension
-    return filename
-
-
 @template_parameter(classtag="operator")
 def lop_template_ansatz_gfs():
     name = "GFSU"
@@ -176,9 +152,8 @@ def name_initree_member():
 
 
 @class_basename(classtag="operator")
-def localoperator_basename(formdata, data):
-    form_name = name_form(formdata, data)
-    return "LocalOperator" + form_name.capitalize()
+def localoperator_basename(form_ident):
+    return get_form_option("classname", form_ident)
 
 
 def class_type_from_cache(classtag):
@@ -244,7 +219,7 @@ def determine_accumulation_space(info, number):
     from loopy.types import NumpyType
     valuearg(lfs, dtype=NumpyType("str"))
 
-    if get_option("blockstructured"):
+    if get_form_option("blockstructured"):
         from dune.perftool.blockstructured.tools import micro_index_to_macro_index
         from dune.perftool.blockstructured.spaces import lfs_inames
         lfsi = micro_index_to_macro_index(subel, lfs_inames(subel, info.restriction, count=number)[0])
@@ -424,10 +399,10 @@ def generate_accumulation_instruction(expr, visitor):
 
 def get_visitor(measure, subdomain_id):
     # Get a transformer instance for this kernel
-    if get_option('sumfact'):
+    if get_form_option('sumfact'):
         from dune.perftool.sumfact import SumFactInterface
         interface = SumFactInterface()
-    elif get_option('blockstructured'):
+    elif get_form_option('blockstructured'):
         from dune.perftool.blockstructured import BlockStructuredInterface
         interface = BlockStructuredInterface()
     else:
@@ -552,8 +527,8 @@ def extract_kernel_from_cache(tag, wrap_in_cgen=True):
     kernel = heuristic_duplication(kernel)
 
     # Maybe apply vectorization strategies
-    if get_option("vectorization_quadloop"):
-        if get_option("sumfact"):
+    if get_form_option("vectorization_quadloop"):
+        if get_form_option("sumfact"):
             from dune.perftool.loopy.transformations.vectorize_quad import vectorize_quadrature_loop
             kernel = vectorize_quadrature_loop(kernel)
         else:
@@ -702,11 +677,14 @@ def cgen_class_from_cache(tag, members=[]):
     return Class(basename, base_classes=base_classes, members=[constructor] + members + pm + decls, tparam_decls=tparams)
 
 
-def generate_localoperator_kernels(formdata, data):
+def generate_localoperator_kernels(operator):
     logger = logging.getLogger(__name__)
 
-    # Extract the relevant attributes of the form data
-    form = formdata.preprocessed_form
+    data = get_global_context_value("data")
+    form = data.object_by_name[get_form_option("form")]
+
+    from dune.perftool.ufl.preprocess import preprocess_form
+    form = preprocess_form(form).preprocessed_form
 
     # Reset the generation cache
     from dune.perftool.generation import delete_cache_items
@@ -726,12 +704,12 @@ def generate_localoperator_kernels(formdata, data):
 
     # Trigger this one once early on to assure that template
     # parameters are set in the right order
-    localoperator_basename(formdata, data)
+    localoperator_basename(operator)
     lop_template_ansatz_gfs()
     lop_template_test_gfs()
     lop_template_range_field()
     from dune.perftool.pdelab.parameter import parameterclass_basename
-    parameterclass_basename(formdata, data)
+    parameterclass_basename(operator)
 
     # Make sure there is always the same constructor arguments (even if parameter class is empty)
     from dune.perftool.pdelab.localoperator import name_initree_member
@@ -769,7 +747,7 @@ def generate_localoperator_kernels(formdata, data):
                     kernel = [k for k in get_backend(interface="generate_kernels_per_integral")(form.integrals_by_type(measure))]
 
                 # Maybe add numerical differentiation
-                if get_option("numerical_jacobian"):
+                if get_form_option("numerical_jacobian"):
                     # Include headers for numerical methods
                     include_file("dune/pdelab/localoperator/defaultimp.hh", filetag="operatorfile")
 
@@ -788,7 +766,7 @@ def generate_localoperator_kernels(formdata, data):
                                      )
 
                     # In the case of matrix free operator evaluation we need jacobian apply methods
-                    if get_option("matrix_free"):
+                    if get_form_option("matrix_free"):
                         from dune.perftool.pdelab.driver import is_linear
                         if is_linear(formdata.original_form):
                             # Numeical jacobian apply base class
@@ -812,7 +790,7 @@ def generate_localoperator_kernels(formdata, data):
                 operator_kernels[(measure, 'residual')] = kernel
 
     # Generate the necessary jacobian methods
-    if not get_option("numerical_jacobian"):
+    if not get_form_option("numerical_jacobian"):
         logger.info("generate_localoperator_kernels: create jacobian methods")
         from ufl import derivative
         jacform = derivative(formdata.original_form, formdata.original_form.coefficients()[0])
@@ -821,7 +799,7 @@ def generate_localoperator_kernels(formdata, data):
         jacform = preprocess_form(jacform).preprocessed_form
 
         with global_context(form_type="jacobian"):
-            if get_option("generate_jacobians"):
+            if get_form_option("generate_jacobians"):
                 for measure in set(i.integral_type() for i in jacform.integrals()):
                     logger.info("generate_localoperator_kernels: measure {}".format(measure))
                     with global_context(integral_type=measure):
@@ -840,7 +818,7 @@ def generate_localoperator_kernels(formdata, data):
                         operator_kernels[(it, 'jacobian')] = [LoopyKernelMethod(assembly_routine_signature(), kernel=None)]
 
         # Jacobian apply methods for matrix-free computations
-        if get_option("matrix_free"):
+        if get_form_option("matrix_free"):
             # The apply vector has reserved index 1 so we directly use Coefficient class from ufl
             from ufl import Coefficient
             apply_coefficient = Coefficient(form.coefficients()[0].ufl_element(), 1)
@@ -871,7 +849,7 @@ def generate_localoperator_kernels(formdata, data):
     return operator_kernels
 
 
-def generate_localoperator_file(formdata, kernels, filename):
+def generate_localoperator_file(kernels, filename):
     operator_methods = []
     for k in kernels.values():
         operator_methods.extend(k)
@@ -888,13 +866,3 @@ def generate_localoperator_file(formdata, kernels, filename):
     # TODO take the name of this thing from the UFL file
     lop = cgen_class_from_cache("operator", members=operator_methods)
     generate_file(filename, "operatorfile", [param, lop])
-
-
-def generate_localoperator_basefile(formdatas, data):
-    filename = get_option("operator_file")
-    for formdata in formdatas:
-        lop_filename = name_localoperator_file(formdata, data)
-        include_file(lop_filename, filetag="operatorbasefile")
-
-    from dune.perftool.file import generate_file
-    generate_file(filename, "operatorbasefile", [])
diff --git a/python/dune/perftool/pdelab/parameter.py b/python/dune/perftool/pdelab/parameter.py
index 36e5cd228fec281a27e17c99fa4006ee95d9e27e..e9e8035d3568b87c991614c52b1da105aecc507a 100644
--- a/python/dune/perftool/pdelab/parameter.py
+++ b/python/dune/perftool/pdelab/parameter.py
@@ -21,13 +21,14 @@ from dune.perftool.pdelab.localoperator import (class_type_from_cache,
                                                 localoperator_basename,
                                                 )
 from dune.perftool.loopy.target import type_floatingpoint
+from dune.perftool.options import get_form_option
 
 from loopy.match import Writes
 
 
 @class_basename(classtag="parameterclass")
-def parameterclass_basename(formdata, data):
-    lopbase = localoperator_basename(formdata, data)
+def parameterclass_basename(form_ident):
+    lopbase = get_form_option("classname", form_ident)
     return "{}Params".format(lopbase)
 
 
diff --git a/python/dune/perftool/pdelab/quadrature.py b/python/dune/perftool/pdelab/quadrature.py
index 02e4a428348d81764ec9672c25370e9be6d14969..bf798fc8392b55a7492a910c1dff586cb43a74d7 100644
--- a/python/dune/perftool/pdelab/quadrature.py
+++ b/python/dune/perftool/pdelab/quadrature.py
@@ -14,7 +14,7 @@ from dune.perftool.generation import (backend,
                                       valuearg,
                                       )
 from dune.perftool.pdelab.localoperator import lop_template_range_field
-from dune.perftool.options import get_option
+from dune.perftool.options import get_form_option
 from dune.perftool.ufl.modified_terminals import Restriction
 
 from pymbolic.primitives import Variable, Subscript
@@ -184,7 +184,8 @@ def _estimate_quadrature_order():
     """Estimate quadrature order using polynomial degree estimation from UFL"""
     # According to UFL documentation estimate_total_polynomial_degree
     # should only be called on preprocessed forms.
-    form = get_global_context_value("formdata").preprocessed_form
+    from dune.perftool.pdelab.driver import get_preprocessed_form
+    form = get_preprocessed_form()
 
     # Estimate polynomial degree of integrals of current type (eg 'Cell')
     integral_type = get_global_context_value("integral_type")
@@ -223,8 +224,8 @@ def quadrature_order():
     - If you use sum factorization and TensorProductElement it is
       possible to use a different quadrature_order per direction.
     """
-    if get_option("quadrature_order"):
-        quadrature_order = tuple(map(int, get_option("quadrature_order").split(',')))
+    if get_form_option("quadrature_order"):
+        quadrature_order = tuple(map(int, get_form_option("quadrature_order").split(',')))
     else:
         quadrature_order = _estimate_quadrature_order()
 
@@ -235,7 +236,7 @@ def quadrature_order():
         if len(quadrature_order) == 1:
             quadrature_order = quadrature_order[0]
     if isinstance(quadrature_order, tuple):
-        if not get_option('sumfact'):
+        if not get_form_option('sumfact'):
             raise NotImplementedError("Different quadrature order per direction is only implemented for kernels using sum factorization.")
         from dune.perftool.pdelab.geometry import world_dimension
         assert(len(quadrature_order) == world_dimension())
diff --git a/python/dune/perftool/sumfact/accumulation.py b/python/dune/perftool/sumfact/accumulation.py
index 3c935161e5d837d2397736a24a996fff4f312384..9566efb71e358dcaf346bd3a24b94858e3a56f26 100644
--- a/python/dune/perftool/sumfact/accumulation.py
+++ b/python/dune/perftool/sumfact/accumulation.py
@@ -16,7 +16,9 @@ from dune.perftool.generation import (backend,
                                       temporary_variable,
                                       transform,
                                       )
-from dune.perftool.options import get_option
+from dune.perftool.options import (get_form_option,
+                                   get_option,
+                                   )
 from dune.perftool.loopy.flatten import flatten_index
 from dune.perftool.loopy.buffer import get_buffer_temporary
 from dune.perftool.sumfact.quadrature import nest_quadrature_loops
@@ -378,7 +380,7 @@ def generate_accumulation_instruction(expr, visitor):
                            (maybe_wrap_subscript(result, prim.Variable(iname)),),
                            )
 
-    if not get_option("fastdg"):
+    if not get_form_option("fastdg"):
         rank = 2 if jacobian_inames else 1
         expr = prim.Call(PDELabAccumulationFunction(accumvar, rank),
                          (test_lfs.get_args() +
diff --git a/python/dune/perftool/sumfact/basis.py b/python/dune/perftool/sumfact/basis.py
index 19278ef7ff5fe6cfbefdd9534f30bbe4ee5ea345..a941b8e0d51e956010f7f11532157993766becf3 100644
--- a/python/dune/perftool/sumfact/basis.py
+++ b/python/dune/perftool/sumfact/basis.py
@@ -33,7 +33,7 @@ from dune.perftool.pdelab.geometry import (local_dimension,
                                            )
 from dune.perftool.loopy.buffer import initialize_buffer, get_buffer_temporary
 from dune.perftool.sumfact.symbolic import SumfactKernel, SumfactKernelInputBase
-from dune.perftool.options import get_option
+from dune.perftool.options import get_form_option
 from dune.perftool.pdelab.driver import FEM_name_mangling
 from dune.perftool.pdelab.restriction import restricted_name
 from dune.perftool.pdelab.spaces import name_lfs, name_lfs_bound, name_leaf_lfs
@@ -93,7 +93,7 @@ class LFSSumfactKernelInput(SumfactKernelInputBase, ImmutableRecord):
 
     @property
     def direct_input(self):
-        if get_option("fastdg"):
+        if get_form_option("fastdg"):
             return self.coeff_func(self.restriction)
         else:
             return None
diff --git a/python/dune/perftool/sumfact/geometry.py b/python/dune/perftool/sumfact/geometry.py
index e17b7aac660718cb424cae73fe2c7c60452fb5f2..0fd86eec08e519960d344cf13cf894f7243f287b 100644
--- a/python/dune/perftool/sumfact/geometry.py
+++ b/python/dune/perftool/sumfact/geometry.py
@@ -20,7 +20,7 @@ from dune.perftool.pdelab.geometry import (local_dimension,
 from dune.perftool.sumfact.switch import get_facedir
 from dune.perftool.sumfact.symbolic import SumfactKernelInputBase
 from dune.perftool.sumfact.vectorization import attach_vectorization_info
-from dune.perftool.options import get_option, option_switch
+from dune.perftool.options import get_form_option, option_switch
 from dune.perftool.ufl.modified_terminals import Restriction
 
 from pytools import ImmutableRecord
@@ -181,7 +181,7 @@ def pymbolic_spatial_coordinate_axiparallel(do_predicates, visitor):
 def pymbolic_unit_outer_normal(visitor_indices):
     index, = visitor_indices
     assert isinstance(index, int)
-    if get_option("diagonal_transformation_matrix"):
+    if get_form_option("diagonal_transformation_matrix"):
         from dune.perftool.sumfact.switch import get_facedir, get_facemod
         if index == get_facedir(Restriction.NEGATIVE):
             if get_facemod(Restriction.NEGATIVE):
@@ -198,7 +198,7 @@ def pymbolic_unit_outer_normal(visitor_indices):
 def pymbolic_unit_inner_normal(visitor_indices):
     index, = visitor_indices
     assert isinstance(index, int)
-    if get_option("diagonal_transformation_matrix"):
+    if get_form_option("diagonal_transformation_matrix"):
         from dune.perftool.sumfact.switch import get_facedir, get_facemod
         if index == get_facedir(Restriction.NEGATIVE):
             if get_facemod(Restriction.NEGATIVE):
@@ -213,7 +213,7 @@ def pymbolic_unit_inner_normal(visitor_indices):
 
 
 def pymbolic_facet_jacobian_determinant():
-    if get_option("constant_transformation_matrix"):
+    if get_form_option("constant_transformation_matrix"):
         return pymbolic_constant_facet_jacobian_determinant()
     else:
         from dune.perftool.pdelab.geometry import pymbolic_facet_jacobian_determinant as _norm
@@ -256,7 +256,7 @@ def define_constant_facet_jacobian_determinant_eval(name):
 
 
 def pymbolic_facet_area():
-    if get_option("constant_transformation_matrix"):
+    if get_form_option("constant_transformation_matrix"):
         return pymbolic_facet_jacobian_determinant()
     else:
         from dune.perftool.pdelab.geometry import pymbolic_facet_area as _norm
diff --git a/python/dune/perftool/sumfact/quadrature.py b/python/dune/perftool/sumfact/quadrature.py
index 8209a41f2890078aba8ffb25b708937fda2c7d01..2b709293a7f2c556c43a3812654c6e178e7d971f 100644
--- a/python/dune/perftool/sumfact/quadrature.py
+++ b/python/dune/perftool/sumfact/quadrature.py
@@ -18,7 +18,7 @@ from dune.perftool.pdelab.argument import name_accumulation_variable
 from dune.perftool.pdelab.geometry import (local_dimension,
                                            world_dimension,
                                            )
-from dune.perftool.options import get_option
+from dune.perftool.options import get_form_option
 from dune.perftool.sumfact.switch import get_facedir
 from dune.perftool.loopy.target import dtype_floatingpoint
 
@@ -142,7 +142,7 @@ def recursive_quadrature_weight(visitor, direction=0):
 
 def quadrature_weight(visitor):
     # Return non-precomputed version
-    if not get_option("precompute_quadrature_info"):
+    if not get_form_option("precompute_quadrature_info"):
         return recursive_quadrature_weight(visitor)
 
     # Quadrature points per (local) direction
@@ -195,7 +195,7 @@ def define_quadrature_position(name, index):
 @backend(interface="quad_pos", name="sumfact")
 def pymbolic_quadrature_position(index, visitor):
     # Return the non-precomputed version
-    if not get_option("precompute_quadrature_info"):
+    if not get_form_option("precompute_quadrature_info"):
         name = 'pos'
         temporary_variable(name, shape=(local_dimension(),), shape_impl=("fv",))
         define_quadrature_position(name, index)
diff --git a/python/dune/perftool/sumfact/realization.py b/python/dune/perftool/sumfact/realization.py
index 703e3e062b72a5615d57d205033a861bb5388ba3..b0e85d43bdaa27409af81596ecf3274b6dbabca7 100644
--- a/python/dune/perftool/sumfact/realization.py
+++ b/python/dune/perftool/sumfact/realization.py
@@ -20,7 +20,9 @@ from dune.perftool.loopy.buffer import (get_buffer_temporary,
 from dune.perftool.pdelab.argument import pymbolic_coefficient
 from dune.perftool.pdelab.basis import shape_as_pymbolic
 from dune.perftool.pdelab.geometry import world_dimension
-from dune.perftool.options import get_option
+from dune.perftool.options import (get_form_option,
+                                   get_option,
+                                   )
 from dune.perftool.pdelab.signatures import assembler_routine_name
 from dune.perftool.sumfact.permutation import (sumfact_permutation_strategy,
                                                permute_backward,
@@ -215,7 +217,7 @@ def _realize_sum_factorization_kernel(sf):
 
         # In case of direct output we directly accumulate the result
         # of the Sumfactorization into some global data structure.
-        if l == len(matrix_sequence) - 1 and get_option('fastdg') and sf.stage == 3:
+        if l == len(matrix_sequence) - 1 and get_form_option('fastdg') and sf.stage == 3:
             ft = get_global_context_value("form_type")
             if sf.test_element_index is None:
                 direct_output = "{}_access".format(sf.accumvar)
diff --git a/python/dune/perftool/sumfact/switch.py b/python/dune/perftool/sumfact/switch.py
index a420850d6bdecb4e479688f7be8c3e1039572351..1fca79729c31dda371a502141eaeabbd7b0a34a2 100644
--- a/python/dune/perftool/sumfact/switch.py
+++ b/python/dune/perftool/sumfact/switch.py
@@ -10,7 +10,7 @@ from dune.perftool.pdelab.signatures import (assembly_routine_args,
                                              assembly_routine_signature,
                                              kernel_name,
                                              )
-from dune.perftool.options import get_option
+from dune.perftool.options import get_form_option
 from dune.perftool.cgen.clazz import ClassMember
 
 
@@ -53,7 +53,7 @@ def get_kernel_name(facedir_s=None, facemod_s=None, facedir_n=None, facemod_n=No
 
 def decide_if_kernel_is_necessary(facedir_s, facemod_s, facedir_n, facemod_n):
     # If we are not using YaspGrid, all variants need to be realized
-    if not get_option("diagonal_transformation_matrix"):
+    if not get_form_option("diagonal_transformation_matrix"):
         return True
 
     # The PDELab machineries visit-once policy combined with Yasp avoids any visits
diff --git a/python/dune/perftool/sumfact/tabulation.py b/python/dune/perftool/sumfact/tabulation.py
index 4009a7bf7a2bb5356f53129f28c80420ff5a8c1a..0fd15bfefdb2b3ac13012ed8ebb188a1bec7cf7f 100644
--- a/python/dune/perftool/sumfact/tabulation.py
+++ b/python/dune/perftool/sumfact/tabulation.py
@@ -1,7 +1,5 @@
 from dune.perftool.ufl.modified_terminals import Restriction
 
-from dune.perftool.options import get_option
-
 from dune.perftool.pdelab.argument import name_coefficientcontainer
 from dune.perftool.pdelab.geometry import world_dimension, local_dimension
 from dune.perftool.generation import (class_member,
diff --git a/python/dune/perftool/sumfact/vectorization.py b/python/dune/perftool/sumfact/vectorization.py
index 4dab904154d80072585453383e482725bfa5da97..38ec9b162e2e604fe485ac5b374feea2944ba696 100644
--- a/python/dune/perftool/sumfact/vectorization.py
+++ b/python/dune/perftool/sumfact/vectorization.py
@@ -19,7 +19,7 @@ from dune.perftool.sumfact.tabulation import (BasisTabulationMatrixArray,
                                               set_quadrature_points,
                                               )
 from dune.perftool.error import PerftoolVectorizationError
-from dune.perftool.options import get_option
+from dune.perftool.options import get_form_option
 from dune.perftool.tools import add_to_frozendict, round_to_multiple, list_diff
 
 from pytools import product
@@ -79,10 +79,10 @@ def fromlist_costmodel(sf):
 def explicit_costfunction(sf):
     # Read the explicitly set values for horizontal and vertical vectorization
     width = get_vcl_type_size(dtype_floatingpoint())
-    horizontal = get_option("vectorization_horizontal")
+    horizontal = get_form_option("vectorization_horizontal")
     if horizontal is None:
         horizontal = width
-    vertical = get_option("vectorization_vertical")
+    vertical = get_form_option("vectorization_vertical")
     if vertical is None:
         vertical = 1
     horizontal = int(horizontal)
@@ -97,7 +97,7 @@ def explicit_costfunction(sf):
 
 def strategy_cost(strategy):
     func = get_backend(interface="vectorization_strategy",
-                       selector=lambda: get_option("vectorization_strategy"))
+                       selector=lambda: get_form_option("vectorization_strategy"))
     keys = set(sf.cache_key for sf in strategy.values())
 
     # Sum over all the sum factorization kernels in the realization
@@ -155,7 +155,7 @@ def decide_vectorization_strategy():
     active_sumfacts = [i for i in all_sumfacts if i.stage == 3 or i in basis_sumfacts]
 
     # If no vectorization is needed, abort now
-    if get_option("vectorization_strategy") == "none":
+    if get_form_option("vectorization_strategy") == "none":
         for sf in all_sumfacts:
             _cache_vectorization_info(sf, sf.copy(buffer=get_counted_variable("buffer")))
         return
@@ -170,7 +170,7 @@ def decide_vectorization_strategy():
     # Optimize over all the possible quadrature point tuples
     #
     quad_points = [quadrature_points_per_direction()]
-    if get_option("vectorization_allow_quadrature_changes"):
+    if get_form_option("vectorization_allow_quadrature_changes"):
         sf = next(iter(active_sumfacts))
         depth = 1
         while depth <= width:
@@ -181,14 +181,14 @@ def decide_vectorization_strategy():
             depth = depth * 2
         quad_points = list(set(quad_points))
 
-    if get_option("vectorization_strategy") == "fromlist":
+    if get_form_option("vectorization_strategy") == "fromlist":
         # This is a bit special and does not follow the minimization procedure at all
 
         def _choose_strategy_from_list(stage1_sumfacts):
             strategy = 0
             for qp in quad_points:
                 for strat in fixed_quad_vectorization_opportunity_generator(frozenset(stage1_sumfacts), width, qp):
-                    if strategy == int(get_option("vectorization_list_index")):
+                    if strategy == int(get_form_option("vectorization_list_index")):
                         set_quadrature_points(qp)
                         # Output the strategy and its cost into a separate file
                         if get_global_context_value("form_type") == "jacobian_apply":
@@ -197,12 +197,12 @@ def decide_vectorization_strategy():
                         return qp, strat
                     strategy = strategy + 1
 
-            raise PerftoolVectorizationError("Specified vectorization list index '{}' was too high!".format(get_option("vectorization_list_index")))
+            raise PerftoolVectorizationError("Specified vectorization list index '{}' was too high!".format(get_form_option("vectorization_list_index")))
 
         s1_sumfacts = frozenset(sf for sf in active_sumfacts if sf.stage == 1)
 
         total = sum(len([s for s in fixed_quad_vectorization_opportunity_generator(frozenset(s1_sumfacts), width, qp)]) for qp in quad_points)
-        print("'fromlist' vectorization is attempting to pick #{} of {} strategies...".format(int(get_option("vectorization_list_index")),
+        print("'fromlist' vectorization is attempting to pick #{} of {} strategies...".format(int(get_form_option("vectorization_list_index")),
                                                                                               total))
         qp, sfdict = _choose_strategy_from_list(s1_sumfacts)
 
diff --git a/python/dune/perftool/ufl/preprocess.py b/python/dune/perftool/ufl/preprocess.py
index 4564fee060d339c93a1f362ccdc7523d258e6856..24d436c9e9dca249c3929f32e704ea2ed72eb804 100644
--- a/python/dune/perftool/ufl/preprocess.py
+++ b/python/dune/perftool/ufl/preprocess.py
@@ -2,7 +2,10 @@
 
 import ufl.classes as uc
 
+from pytools import memoize
 
+
+@memoize
 def preprocess_form(form):
     from ufl.algorithms import compute_form_data
     formdata = compute_form_data(form,
diff --git a/python/dune/perftool/ufl/transformations/__init__.py b/python/dune/perftool/ufl/transformations/__init__.py
index dde8a965177b657a0ed9554302fa3e95c8bf7825..de66173a2a04bf3d5f5b7873d8444d230e665ca0 100644
--- a/python/dune/perftool/ufl/transformations/__init__.py
+++ b/python/dune/perftool/ufl/transformations/__init__.py
@@ -19,10 +19,10 @@ class UFLTransformationWrapper(object):
             return
 
         # Write out a dot file
-        from dune.perftool.options import get_option
-        if get_option("print_transformations"):
+        from dune.perftool.options import get_form_option
+        if get_form_option("print_transformations"):
             import os
-            dir = get_option("print_transformations_dir")
+            dir = get_form_option("print_transformations_dir")
 
             for i, exprtowrite in enumerate(expr):
                 filename = "trafo_{}_{}_{}{}.dot".format(self.name, str(self.counter).zfill(4), "in" if before else "out", "_{}".format(i) if len(expr) > 1 else "")
diff --git a/python/dune/perftool/ufl/visitor.py b/python/dune/perftool/ufl/visitor.py
index e9dab7f919ba5c27ad2a8e93676fcebc0f89a451..bdce4c16afc1ec8d5e5af05ff0792869794459fd 100644
--- a/python/dune/perftool/ufl/visitor.py
+++ b/python/dune/perftool/ufl/visitor.py
@@ -13,7 +13,7 @@ from dune.perftool.ufl.modified_terminals import (ModifiedTerminalTracker,
                                                   Restriction,
                                                   )
 from dune.perftool.tools import maybe_wrap_subscript
-from dune.perftool.options import get_option
+from dune.perftool.options import get_form_option
 from dune.perftool.pdelab.parameter import name_paramclass, name_time
 from loopy import Reduction
 
@@ -61,7 +61,7 @@ class UFL2LoopyVisitor(ModifiedTerminalTracker):
             self.current_info = info
             expr = self._call(o, False)
             if expr != 0:
-                if get_option("simplify"):
+                if get_form_option("simplify"):
                     from dune.perftool.sympy import simplify_pymbolic_expression
                     expr = simplify_pymbolic_expression(expr)
                 self.interface.generate_accumulation_instruction(expr, self)
@@ -471,7 +471,7 @@ class UFL2LoopyVisitor(ModifiedTerminalTracker):
         self.indices = None
 
         # Implement diagonal jacobians for unrolled matrices!
-        if get_option("diagonal_transformation_matrix"):
+        if get_form_option("diagonal_transformation_matrix"):
             if isinstance(i, int) and isinstance(j, int) and i != j:
                 return 0
 
diff --git a/test/laplace/laplace.mini b/test/laplace/laplace.mini
index 1db0ffd82475b75b830f085586bb2d991af2c514..4926306f0716c84a6bbfe1b5a02c0ad4acc6fe51 100644
--- a/test/laplace/laplace.mini
+++ b/test/laplace/laplace.mini
@@ -7,5 +7,5 @@ elements = 4 4
 elementType = simplical
 printmatrix = true
 
-[formcompiler]
+[formcompiler.operator]
 numerical_jacobian = 0, 1 | expand num
diff --git a/test/laplace/laplace.ufl b/test/laplace/laplace.ufl
index 29b6a4bd3d6da839ca622098a222079d716a4f9e..162a1ace0316fb0069e091928f6d6aaeffd29963 100644
--- a/test/laplace/laplace.ufl
+++ b/test/laplace/laplace.ufl
@@ -2,4 +2,4 @@ V = FiniteElement("CG", "triangle", 1)
 u = TrialFunction(V)
 v = TestFunction(V)
 
-forms = [inner(grad(u), grad(v))*dx]
+r = inner(grad(u), grad(v))*dx
diff --git a/test/laplace/laplace_dg.mini b/test/laplace/laplace_dg.mini
index 04a3c3dd0b0c3a5ee4f7d01c073b00f03d498b72..7485e186eb827cd25b60e5fff65202c7d024147d 100644
--- a/test/laplace/laplace_dg.mini
+++ b/test/laplace/laplace_dg.mini
@@ -7,5 +7,5 @@ elements = 2 2
 elementType = simplical
 printmatrix = true
 
-[formcompiler]
+[formcompiler.operator]
 numerical_jacobian = 1, 0 | expand num
diff --git a/test/nonlinear/nonlinear.ufl b/test/nonlinear/nonlinear.ufl
index b22d743777aec41147d727988e020785b7b710bc..425c6cee834c2e940ef2a6b698ac5eb22eaf4bd9 100644
--- a/test/nonlinear/nonlinear.ufl
+++ b/test/nonlinear/nonlinear.ufl
@@ -10,7 +10,6 @@ v = TestFunction(V)
 
 r = (inner(grad(u), grad(v)) + u*u*v - f*v)*dx
 
-forms = [r]
 exact_solution = g
 interpolate_expression = g
 is_dirichlet = 1
\ No newline at end of file
diff --git a/test/nonlinear/nonlinear_dg.ufl b/test/nonlinear/nonlinear_dg.ufl
index 5fba927e4256b2c32d1b4a5766c78046dffbddc1..47b850759e24461a394c230b112a7e5929ae1ea4 100644
--- a/test/nonlinear/nonlinear_dg.ufl
+++ b/test/nonlinear/nonlinear_dg.ufl
@@ -32,5 +32,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/nonlinear/nonlinear_dg_matrix_free.mini b/test/nonlinear/nonlinear_dg_matrix_free.mini
index 8a1de68470a4375d05253337686105bcb35e6fea..7b53a3a4be6d383f2950d242bf38ded4497939fa 100644
--- a/test/nonlinear/nonlinear_dg_matrix_free.mini
+++ b/test/nonlinear/nonlinear_dg_matrix_free.mini
@@ -11,9 +11,11 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
+compare_l2errorsquared = 5e-3
+
+[formcompiler.operator]
 numerical_jacobian = 1, 0 | expand num
 matrix_free = 1
-compare_l2errorsquared = 5e-3
 
 # Disable numerical methods (not working in PDELab?)
 {__exec_suffix} == numdiff | exclude
diff --git a/test/nonlinear/nonlinear_matrix_free.mini b/test/nonlinear/nonlinear_matrix_free.mini
index 644297ab361203ec049c27529904d8b09dd2261a..fdef068b581db4f6d3f786cc5b77e01af7141804 100644
--- a/test/nonlinear/nonlinear_matrix_free.mini
+++ b/test/nonlinear/nonlinear_matrix_free.mini
@@ -11,6 +11,8 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
+compare_l2errorsquared = 6e-4
+
+[formcompiler.operator]
 numerical_jacobian = 1, 0 | expand num
 matrix_free = 1
-compare_l2errorsquared = 6e-4
diff --git a/test/poisson/dimension-grid-variations/poisson_1d_cg_interval.ufl b/test/poisson/dimension-grid-variations/poisson_1d_cg_interval.ufl
index 4125994aadff64d73b54ce2d2e86b9be4462eebf..d584ea4b577deeec3a893e2ca636d0371c36fb5f 100644
--- a/test/poisson/dimension-grid-variations/poisson_1d_cg_interval.ufl
+++ b/test/poisson/dimension-grid-variations/poisson_1d_cg_interval.ufl
@@ -9,7 +9,7 @@ V = FiniteElement("CG", cell, 1)
 u = TrialFunction(V)
 v = TestFunction(V)
 
-forms = [(inner(grad(u), grad(v)) - f*v)*dx]
+r = (inner(grad(u), grad(v)) - f*v)*dx
 exact_solution = g
 is_dirichlet = 1
 interpolate_expression = g
diff --git a/test/poisson/dimension-grid-variations/poisson_1d_dg_interval.ufl b/test/poisson/dimension-grid-variations/poisson_1d_dg_interval.ufl
index 535d102440fa6c5f5405f61ab0ef82285dfeaa7a..15a60efc606078c91fcd91b77c3b21180d73808e 100644
--- a/test/poisson/dimension-grid-variations/poisson_1d_dg_interval.ufl
+++ b/test/poisson/dimension-grid-variations/poisson_1d_dg_interval.ufl
@@ -29,5 +29,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/poisson/dimension-grid-variations/poisson_2d_cg_quadrilateral.ufl b/test/poisson/dimension-grid-variations/poisson_2d_cg_quadrilateral.ufl
index 5b870f088e98633790d8e6f21c67e7cc432f14de..5d1921828127bbe2024d36c8674bbe4d9c868190 100644
--- a/test/poisson/dimension-grid-variations/poisson_2d_cg_quadrilateral.ufl
+++ b/test/poisson/dimension-grid-variations/poisson_2d_cg_quadrilateral.ufl
@@ -9,7 +9,7 @@ V = FiniteElement("CG", cell, 1)
 u = TrialFunction(V)
 v = TestFunction(V)
 
-forms = [(inner(grad(u), grad(v)) - f*v)*dx]
+r = (inner(grad(u), grad(v)) - f*v)*dx
 exact_solution = g
 is_dirichlet = 1
 interpolate_expression = g
diff --git a/test/poisson/dimension-grid-variations/poisson_2d_cg_triangle.ufl b/test/poisson/dimension-grid-variations/poisson_2d_cg_triangle.ufl
index c07533d54897d68192803935ac2b212a46376402..a720e454218bd38347e62246a84e61fb38d5a0ac 100644
--- a/test/poisson/dimension-grid-variations/poisson_2d_cg_triangle.ufl
+++ b/test/poisson/dimension-grid-variations/poisson_2d_cg_triangle.ufl
@@ -8,7 +8,7 @@ V = FiniteElement("CG", "triangle", 1)
 u = TrialFunction(V)
 v = TestFunction(V)
 
-forms = [(inner(grad(u), grad(v)) - f*v)*dx]
+r = (inner(grad(u), grad(v)) - f*v)*dx
 exact_solution = g
 is_dirichlet = 1
 interpolate_expression = g
diff --git a/test/poisson/dimension-grid-variations/poisson_2d_dg_quadrilateral.ufl b/test/poisson/dimension-grid-variations/poisson_2d_dg_quadrilateral.ufl
index 84019c1e87a0a368429d49a9c787344dc6bcdaf2..5b4cf2d71b09181e307ce9dbe835a0d030e05935 100644
--- a/test/poisson/dimension-grid-variations/poisson_2d_dg_quadrilateral.ufl
+++ b/test/poisson/dimension-grid-variations/poisson_2d_dg_quadrilateral.ufl
@@ -29,5 +29,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/poisson/dimension-grid-variations/poisson_2d_dg_triangle.ufl b/test/poisson/dimension-grid-variations/poisson_2d_dg_triangle.ufl
index d518b11815bc40279f08f2a8802f41c02d2f5703..2121ff1b2a635ec4128d521d143f0456c046911a 100644
--- a/test/poisson/dimension-grid-variations/poisson_2d_dg_triangle.ufl
+++ b/test/poisson/dimension-grid-variations/poisson_2d_dg_triangle.ufl
@@ -28,5 +28,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/poisson/dimension-grid-variations/poisson_3d_cg_hexahedron.ufl b/test/poisson/dimension-grid-variations/poisson_3d_cg_hexahedron.ufl
index 7bfd73703228482b235fab077e5a2f874132ae04..9a9b16fc17df5b1d9df5d96b7e7e9e49853a09c5 100644
--- a/test/poisson/dimension-grid-variations/poisson_3d_cg_hexahedron.ufl
+++ b/test/poisson/dimension-grid-variations/poisson_3d_cg_hexahedron.ufl
@@ -8,7 +8,7 @@ V = FiniteElement("CG", cell, 1)
 u = TrialFunction(V)
 v = TestFunction(V)
 
-forms = [(inner(grad(u), grad(v)) - f*v)*dx]
+r = (inner(grad(u), grad(v)) - f*v)*dx
 exact_solution = g
 is_dirichlet = 1
 interpolate_expression = g
diff --git a/test/poisson/dimension-grid-variations/poisson_3d_cg_tetrahedron.ufl b/test/poisson/dimension-grid-variations/poisson_3d_cg_tetrahedron.ufl
index b6e70a98065893ed033907bfea4ef418d39d7363..b9df144c0534c800675cdeaf95343cb504b82bfe 100644
--- a/test/poisson/dimension-grid-variations/poisson_3d_cg_tetrahedron.ufl
+++ b/test/poisson/dimension-grid-variations/poisson_3d_cg_tetrahedron.ufl
@@ -8,7 +8,7 @@ V = FiniteElement("CG", "tetrahedron", 1)
 u = TrialFunction(V)
 v = TestFunction(V)
 
-forms = [(inner(grad(u), grad(v)) - f*v)*dx]
+r = (inner(grad(u), grad(v)) - f*v)*dx
 exact_solution = g
 is_dirichlet = 1
 interpolate_expression = g
diff --git a/test/poisson/dimension-grid-variations/poisson_3d_dg_hexahedron.ufl b/test/poisson/dimension-grid-variations/poisson_3d_dg_hexahedron.ufl
index 8abb5fb0d42197f3b87f38389f222efafeccb30b..3925b02e2f18c7e477aa6109dd4b56fd4ce9b740 100644
--- a/test/poisson/dimension-grid-variations/poisson_3d_dg_hexahedron.ufl
+++ b/test/poisson/dimension-grid-variations/poisson_3d_dg_hexahedron.ufl
@@ -29,5 +29,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/poisson/dimension-grid-variations/poisson_3d_dg_tetrahedron.ufl b/test/poisson/dimension-grid-variations/poisson_3d_dg_tetrahedron.ufl
index 384d582e45cd0d820e0d57ef40abb2d07ada88e2..5cc027cf1fc4ddf0eca4109982c15cb9e340c19c 100644
--- a/test/poisson/dimension-grid-variations/poisson_3d_dg_tetrahedron.ufl
+++ b/test/poisson/dimension-grid-variations/poisson_3d_dg_tetrahedron.ufl
@@ -28,5 +28,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/poisson/opcount_poisson_dg.ufl b/test/poisson/opcount_poisson_dg.ufl
index 1748e05796f0a0c82073c5b7fb24db9631891c71..5db12e14b8a85c7c417534269c18d52e7a3e3596 100644
--- a/test/poisson/opcount_poisson_dg.ufl
+++ b/test/poisson/opcount_poisson_dg.ufl
@@ -35,5 +35,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma_ext*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/poisson/poisson.mini b/test/poisson/poisson.mini
index 3597ac0536fddfa5dcb83f5e6f4fdafe09419360..a851012e8218b44aef8699f795f031acf5695350 100644
--- a/test/poisson/poisson.mini
+++ b/test/poisson/poisson.mini
@@ -12,5 +12,7 @@ reference = poisson_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
 compare_l2errorsquared = 1e-7
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
diff --git a/test/poisson/poisson.ufl b/test/poisson/poisson.ufl
index 025f1e8ca033d97f0ae0399dbf936d72eee9187d..5c6cf421ebc94b7419a8aa1d8d6df29f28cf0bb2 100644
--- a/test/poisson/poisson.ufl
+++ b/test/poisson/poisson.ufl
@@ -11,7 +11,7 @@ u = TrialFunction(V)
 v = TestFunction(V)
 
 
-forms = [(inner(grad(u), grad(v)) - f*v)*dx]
+r = (inner(grad(u), grad(v)) - f*v)*dx
 exact_solution = g
 interpolate_expression = g
-is_dirichlet = 1
\ No newline at end of file
+is_dirichlet = 1
diff --git a/test/poisson/poisson_dg.mini b/test/poisson/poisson_dg.mini
index fd859d6045eb496ab40e5a67a69d14e289bd9d19..5973973e0ac882c6c70d93cdc39990e1fd9001a2 100644
--- a/test/poisson/poisson_dg.mini
+++ b/test/poisson/poisson_dg.mini
@@ -12,5 +12,7 @@ reference = poisson_dg_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
 compare_l2errorsquared = 9e-8
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
\ No newline at end of file
diff --git a/test/poisson/poisson_dg.ufl b/test/poisson/poisson_dg.ufl
index 9ac02bdf061feb9788488ba188e9d3b864ed59e0..fa263e98c6732ba53ca58ab4fd7daa8422d86f36 100644
--- a/test/poisson/poisson_dg.ufl
+++ b/test/poisson/poisson_dg.ufl
@@ -35,5 +35,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma_ext*g*v*ds
 
-forms = [r]
 exact_solution = g
\ No newline at end of file
diff --git a/test/poisson/poisson_dg_matrix_free.mini b/test/poisson/poisson_dg_matrix_free.mini
index 4fa26853c64f1075cc0de65ba1937890066bf54b..ca1f11460caa0bdeedb63e8d74a0b43b3402b860 100644
--- a/test/poisson/poisson_dg_matrix_free.mini
+++ b/test/poisson/poisson_dg_matrix_free.mini
@@ -12,6 +12,8 @@ reference = poisson_dg_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
-matrix_free = 1
 compare_l2errorsquared = 1e-6
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
+matrix_free = 1
\ No newline at end of file
diff --git a/test/poisson/poisson_dg_neumann.mini b/test/poisson/poisson_dg_neumann.mini
index 43157de9980380b50d1a5f8a29df2997309ef7cd..f02d6c821cc5f0799768c6219d2d1465670579e5 100644
--- a/test/poisson/poisson_dg_neumann.mini
+++ b/test/poisson/poisson_dg_neumann.mini
@@ -12,5 +12,7 @@ reference = poisson_dg_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
 compare_l2errorsquared = 9e-8
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
diff --git a/test/poisson/poisson_dg_neumann.ufl b/test/poisson/poisson_dg_neumann.ufl
index 3cedfc0877164e86e6edef343de11784a99d1d0e..0e3a899544f9c2d5941e6d954d2a85e3f65f51ac 100644
--- a/test/poisson/poisson_dg_neumann.ufl
+++ b/test/poisson/poisson_dg_neumann.ufl
@@ -41,5 +41,4 @@ r = inner(grad(u), grad(v))*dx \
   - gamma_ext*g*v*ds(1) \
   - j*v*ds(0)
 
-forms = [r]
 exact_solution = g
diff --git a/test/poisson/poisson_dg_quadrilateral.mini b/test/poisson/poisson_dg_quadrilateral.mini
index 05da536acc9dd1864dd59f5d817f397e131ed541..2a3cef7057e47d3a6c917bacffbf9a21c458362d 100644
--- a/test/poisson/poisson_dg_quadrilateral.mini
+++ b/test/poisson/poisson_dg_quadrilateral.mini
@@ -10,5 +10,7 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
 compare_l2errorsquared = 7e-7
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
diff --git a/test/poisson/poisson_dg_quadrilateral.ufl b/test/poisson/poisson_dg_quadrilateral.ufl
index 30b1eb8bcd225600c9f78801f3498077388444f3..c8e1107df230b5af44cc2262c5e5c20f262ada3f 100644
--- a/test/poisson/poisson_dg_quadrilateral.ufl
+++ b/test/poisson/poisson_dg_quadrilateral.ufl
@@ -29,5 +29,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/poisson/poisson_dg_tensor.mini b/test/poisson/poisson_dg_tensor.mini
index 52df8e1da0f9cadae3d553e00092c053dafac4c1..9db75ec7f7c519e2f315cb234f9099ef94c5d465 100644
--- a/test/poisson/poisson_dg_tensor.mini
+++ b/test/poisson/poisson_dg_tensor.mini
@@ -12,5 +12,7 @@ reference = poisson_dg_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
 compare_l2errorsquared = 4e-6
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
\ No newline at end of file
diff --git a/test/poisson/poisson_dg_tensor.ufl b/test/poisson/poisson_dg_tensor.ufl
index 9409ece67ae08270896caf5adf87ced67ebced01..437e2f708bdcc4a6b455920c0864c20ccd406a8c 100644
--- a/test/poisson/poisson_dg_tensor.ufl
+++ b/test/poisson/poisson_dg_tensor.ufl
@@ -38,5 +38,4 @@ r = (inner(A*grad(u), grad(v)) + c*u*v)*dx \
   - theta*g*inner(A*grad(v), n)*ds \
   - gamma_ext*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/poisson/poisson_matrix_free.mini b/test/poisson/poisson_matrix_free.mini
index 5709ef999b6447810847c80948fb010994fd07fc..8d0e9d1d7b318f83c1498ec358a9ac35cae227b7 100644
--- a/test/poisson/poisson_matrix_free.mini
+++ b/test/poisson/poisson_matrix_free.mini
@@ -11,5 +11,7 @@ reference = poisson_ref
 extension = vtu
 
 [formcompiler]
-matrix_free = 1
 compare_l2errorsquared = 1e-7
+
+[formcompiler.operator]
+matrix_free = 1
\ No newline at end of file
diff --git a/test/poisson/poisson_neumann.mini b/test/poisson/poisson_neumann.mini
index 0c4aa9c7ff2f7470f9f39e146c933dbb4f452ee0..bf836b6fe1cc1a358418e22b6b50524553f9c50e 100644
--- a/test/poisson/poisson_neumann.mini
+++ b/test/poisson/poisson_neumann.mini
@@ -12,5 +12,7 @@ reference = poisson_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
 compare_l2errorsquared = 8e-8
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
diff --git a/test/poisson/poisson_neumann.ufl b/test/poisson/poisson_neumann.ufl
index aa69d6fe42e68476d80da264c625ebead752436d..16eea674b5fd40003af4817ee2c43d461b45e15e 100644
--- a/test/poisson/poisson_neumann.ufl
+++ b/test/poisson/poisson_neumann.ufl
@@ -16,7 +16,7 @@ v = TestFunction(V)
 # Define the boundary measure that knows where we are...
 ds = ds(subdomain_data=bctype)
 
-forms = [(inner(grad(u), grad(v)) - f*v)*dx - j*v*ds(0)]
+r = (inner(grad(u), grad(v)) - f*v)*dx - j*v*ds(0)
 exact_solution = g
 is_dirichlet = bctype
 interpolate_expression = g
\ No newline at end of file
diff --git a/test/poisson/poisson_tensor.mini b/test/poisson/poisson_tensor.mini
index 8711de545698818505fc63e8701206bce8acdf54..11fe49af49ba71a41d2328feb8d92f923485c793 100644
--- a/test/poisson/poisson_tensor.mini
+++ b/test/poisson/poisson_tensor.mini
@@ -12,5 +12,7 @@ reference = poisson_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
-compare_l2errorsquared = 1e-7
\ No newline at end of file
+compare_l2errorsquared = 1e-7
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
\ No newline at end of file
diff --git a/test/poisson/poisson_tensor.ufl b/test/poisson/poisson_tensor.ufl
index 08898c93f560a329072736e1557b2a19f7645646..b527d05258667dae629f608a1a630e5f11f947b8 100644
--- a/test/poisson/poisson_tensor.ufl
+++ b/test/poisson/poisson_tensor.ufl
@@ -12,7 +12,7 @@ V = FiniteElement("CG", cell, 1)
 u = TrialFunction(V)
 v = TestFunction(V)
 
-forms = [(inner(A*grad(u), grad(v)) + c*u*v -f*v)*dx]
+r= (inner(A*grad(u), grad(v)) + c*u*v -f*v)*dx
 exact_solution = g
 is_dirichlet = 1
 interpolate_expression = g
\ No newline at end of file
diff --git a/test/stokes/stokes.mini b/test/stokes/stokes.mini
index c236f3a87ddec3867557ef6056f3096c04f4be1e..61f68cf98fe3ac9e3dd6fdcadaeb00de79d795b6 100644
--- a/test/stokes/stokes.mini
+++ b/test/stokes/stokes.mini
@@ -13,5 +13,7 @@ reference = hagenpoiseuille_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 0, 1 | expand num
 compare_l2errorsquared = 1e-11
+
+[formcompiler.operator]
+numerical_jacobian = 0, 1 | expand num
diff --git a/test/stokes/stokes.ufl b/test/stokes/stokes.ufl
index 6debda650ade492ace562e3f671cefeb8356afc6..99f21bbc943fd51fa1e514b7438e776e8a91e1f0 100644
--- a/test/stokes/stokes.ufl
+++ b/test/stokes/stokes.ufl
@@ -13,7 +13,6 @@ u, p = TrialFunctions(TH)
 
 r = (inner(grad(v), grad(u)) - div(v)*p - q*div(u))*dx
 
-forms = [r]
 is_dirichlet = v_bctype, v_bctype, 0
 interpolate_expression = g_v, None
 exact_solution = g_v, 8.*(1.-x[0])
\ No newline at end of file
diff --git a/test/stokes/stokes_3d_dg_quadrilateral.mini b/test/stokes/stokes_3d_dg_quadrilateral.mini
index d7c82422aab9a3f4522d9ba28c41514da2345f8f..a66a10cd617232e3bb4ac8e3fbbba52570dea53a 100644
--- a/test/stokes/stokes_3d_dg_quadrilateral.mini
+++ b/test/stokes/stokes_3d_dg_quadrilateral.mini
@@ -10,5 +10,7 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
+compare_l2errorsquared = 6e-8
+
+[formcompiler.operator]
 numerical_jacobian = 0, 1 | expand num
-compare_l2errorsquared = 6e-8
\ No newline at end of file
diff --git a/test/stokes/stokes_3d_dg_quadrilateral.ufl b/test/stokes/stokes_3d_dg_quadrilateral.ufl
index 84d1003e16d7f4f36dc0630434d98eb3d633cd3a..3fefa9b763714ca27b7fa7d741bfe94059dae836 100644
--- a/test/stokes/stokes_3d_dg_quadrilateral.ufl
+++ b/test/stokes/stokes_3d_dg_quadrilateral.ufl
@@ -32,5 +32,4 @@ r = inner(grad(u), grad(v))*dx \
   + p*inner(v, n)*ds \
   + q*inner(u-g_v, n)*ds
 
-forms = [r]
 exact_solution = g_v, 8.*(1.-x[0])
diff --git a/test/stokes/stokes_3d_quadrilateral.mini b/test/stokes/stokes_3d_quadrilateral.mini
index 89c4796da75f3212ca59f44cace5a53229b1259a..fb12e781dfa2138ffa0c1a758c7bfe872cbeaf6c 100644
--- a/test/stokes/stokes_3d_quadrilateral.mini
+++ b/test/stokes/stokes_3d_quadrilateral.mini
@@ -11,5 +11,7 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
+compare_l2errorsquared = 1e-10
+
+[formcompiler.operator]
 numerical_jacobian = 1, 0 | expand num
-compare_l2errorsquared = 1e-10
\ No newline at end of file
diff --git a/test/stokes/stokes_3d_quadrilateral.ufl b/test/stokes/stokes_3d_quadrilateral.ufl
index 24d799a8ebea2dc902019efa6d93e289181282bf..f39cdcd42803cc94fd03ca98059707bafc1a844a 100644
--- a/test/stokes/stokes_3d_quadrilateral.ufl
+++ b/test/stokes/stokes_3d_quadrilateral.ufl
@@ -13,7 +13,6 @@ u, p = TrialFunctions(TH)
 
 r = (inner(grad(v), grad(u)) - div(v)*p - q*div(u))*dx
 
-forms = [r]
 exact_solution = g_v, 8.*(1.-x[0])
 is_dirichlet = v_bctype, v_bctype, v_bctype, 0
 interpolate_expression = g_v, None
diff --git a/test/stokes/stokes_dg.mini b/test/stokes/stokes_dg.mini
index 253a347941a128a1f2b38d736a055c6f827a3088..e09fdb1ce0afa8b381f3e2b76bcc2c360b07d11f 100644
--- a/test/stokes/stokes_dg.mini
+++ b/test/stokes/stokes_dg.mini
@@ -15,5 +15,7 @@ zeroThreshold.data_0 = 1e-6
 zeroThreshold.data_1 = 1e-6
 
 [formcompiler]
-numerical_jacobian = 0, 1 | expand num
 compare_l2errorsquared = 1e-9
+
+[formcompiler.operator]
+numerical_jacobian = 0, 1 | expand num
diff --git a/test/stokes/stokes_dg.ufl b/test/stokes/stokes_dg.ufl
index e7176b264efaf2bbb20146d5c19cd8d4b72a2c45..22d239c8166311e74b9dbb489f7b4499a5e9aaa3 100644
--- a/test/stokes/stokes_dg.ufl
+++ b/test/stokes/stokes_dg.ufl
@@ -32,5 +32,4 @@ r = inner(grad(u), grad(v))*dx \
   + p*inner(v, n)*ds \
   + q*inner(u-g_v, n)*ds
 
-forms = [r]
 exact_solution = g_v, 8*(1.-x[0])
\ No newline at end of file
diff --git a/test/stokes/stokes_dg_quadrilateral.mini b/test/stokes/stokes_dg_quadrilateral.mini
index 7f25099677036a8bfe1f715984ee500d9dc4d015..c89a5d3db4596801a454aac2479ff1f4d02fb35c 100644
--- a/test/stokes/stokes_dg_quadrilateral.mini
+++ b/test/stokes/stokes_dg_quadrilateral.mini
@@ -10,5 +10,7 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 0, 1 | expand num
 compare_l2errorsquared = 1e-8
+
+[formcompiler.operator]
+numerical_jacobian = 0, 1 | expand num
diff --git a/test/stokes/stokes_dg_quadrilateral.ufl b/test/stokes/stokes_dg_quadrilateral.ufl
index 8f4415a7a3c30635e9add35faf58336995a820e4..3342fd9a67cc1961817cc9537e50623b012f6ac0 100644
--- a/test/stokes/stokes_dg_quadrilateral.ufl
+++ b/test/stokes/stokes_dg_quadrilateral.ufl
@@ -32,5 +32,4 @@ r = inner(grad(u), grad(v))*dx \
   + p*inner(v, n)*ds \
   + q*inner(u-g_v, n)*ds
 
-forms = [r]
 exact_solution = g_v, 8*(1.-x[0])
diff --git a/test/stokes/stokes_quadrilateral.mini b/test/stokes/stokes_quadrilateral.mini
index e9440771716292bda9664a26b2c8911e00a65a34..bad65475c66d6994474d79d011dde5ce32e90f33 100644
--- a/test/stokes/stokes_quadrilateral.mini
+++ b/test/stokes/stokes_quadrilateral.mini
@@ -11,5 +11,7 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
 compare_l2errorsquared = 1e-10
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
diff --git a/test/stokes/stokes_quadrilateral.ufl b/test/stokes/stokes_quadrilateral.ufl
index 4fa3f5de3749b1e5b9a6f687603378e9049b999f..4411e791138cf85824f24fc6549d52cf6ef1af0d 100644
--- a/test/stokes/stokes_quadrilateral.ufl
+++ b/test/stokes/stokes_quadrilateral.ufl
@@ -13,7 +13,6 @@ u, p = TrialFunctions(TH)
 
 r = (inner(grad(v), grad(u)) - div(v)*p - q*div(u))*dx
 
-forms = [r]
 is_dirichlet = v_bctype, v_bctype, 0
 interpolate_expression = g_v, None
 exact_solution = g_v, 8.*(1.-x[0])
diff --git a/test/stokes/stokes_stress.mini b/test/stokes/stokes_stress.mini
index af72867f443e85aa8b824f66a0c155b9b9f4e9a6..1ac9b06b62e64db652d757673b1f1082978c3490 100644
--- a/test/stokes/stokes_stress.mini
+++ b/test/stokes/stokes_stress.mini
@@ -15,6 +15,7 @@ reference = hagenpoiseuille_ref
 extension = vtu
 
 [formcompiler]
-# numerical_jacobian = 0, 1 | expand num
-numerical_jacobian = 1
 compare_l2errorsquared = 1e-11
+
+[formcompiler.operator]
+numerical_jacobian = 1
diff --git a/test/stokes/stokes_stress.ufl b/test/stokes/stokes_stress.ufl
index 8e8aaf14a831322ef8539a80647f426835b897a3..787e5a232383ad9e0767a6fb41e08cedc988943f 100644
--- a/test/stokes/stokes_stress.ufl
+++ b/test/stokes/stokes_stress.ufl
@@ -14,7 +14,6 @@ u, p, S  = TrialFunctions(TH)
 
 r = (inner(grad(v), S) + inner(grad(u) - S, T) - div(v)*p - q*div(u))*dx
 
-forms = [r]
 is_dirichlet = v_bctype, v_bctype, 0, 0, 0, 0, 0
 interpolate_expression = 4*x[1]*(1.-x[1]), 0.0, None, None, None, None, None
 exact_solution = 4*x[1]*(1.-x[1]), 0.0, 8*(1.-x[0]), 0.0, 0.0, -1.*8*x[1] + 4., 0.0
\ No newline at end of file
diff --git a/test/stokes/stokes_stress_sym.mini b/test/stokes/stokes_stress_sym.mini
index 1aa3d6f087cae99310dd2c7b3a50721f3fd447f6..ec5f2a0d9a45768b87e9157dafff1bed2c9c2e5a 100644
--- a/test/stokes/stokes_stress_sym.mini
+++ b/test/stokes/stokes_stress_sym.mini
@@ -13,5 +13,7 @@ reference = hagenpoiseuille_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1
 compare_l2errorsquared = 1e-6
+
+[formcompiler.operator]
+numerical_jacobian = 1
diff --git a/test/stokes/stokes_stress_sym.ufl b/test/stokes/stokes_stress_sym.ufl
index 012f70dbd5de6dccdf5f88021aaec320b460c970..8e2d55dd4cfca7335932b751d5ae2c2bb71aad59 100644
--- a/test/stokes/stokes_stress_sym.ufl
+++ b/test/stokes/stokes_stress_sym.ufl
@@ -20,7 +20,6 @@ r = (inner(grad(v), S) + inner(2*sym(grad(u)) - S, T) - div(v)*p - q*div(u))*dx
 # \
 #  + inner(S.T*n, v)*ds
 
-forms = [r]
 is_dirichlet = v_bctype, v_bctype, 0, 0, 0, 0
 interpolate_expression = 4*x[1]*(1.-x[1]), 0.0, None, None, None, None
 exact_solution = 4*x[1]*(1.-x[1]), 0.0, 8*(1.-x[0]), 0.0, 0.0, -1.*8*x[1] + 4.
\ No newline at end of file
diff --git a/test/stokes/stokes_sym.mini b/test/stokes/stokes_sym.mini
index 26cc91467a701e8643b040453943d83a2ece3e96..5755676adb0997f4ee61c5f98a282dfa4284fb13 100644
--- a/test/stokes/stokes_sym.mini
+++ b/test/stokes/stokes_sym.mini
@@ -13,5 +13,7 @@ reference = hagenpoiseuille_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 0, 1 | expand num
 compare_l2errorsquared = 1e-10
+
+[formcompiler.operator]
+numerical_jacobian = 0, 1 | expand num
diff --git a/test/stokes/stokes_sym.ufl b/test/stokes/stokes_sym.ufl
index 554101ab123556bfa3affb47fcc6726c88b3d9d8..c7fe07ceafa660906bb33e21ad9603758d23ed3a 100644
--- a/test/stokes/stokes_sym.ufl
+++ b/test/stokes/stokes_sym.ufl
@@ -16,7 +16,6 @@ n = FacetNormal(triangle)('+')
 
 r = (inner(2*sym(grad(u)), grad(v)) - div(v)*p - q*div(u))*dx - inner(grad(u).T*n,v)*ds
 
-forms = [r]
 is_dirichlet = v_bctype, v_bctype, 0
 interpolate_expression = g_v, None
 exact_solution = g_v, 8.*(1.-x[0])
\ No newline at end of file
diff --git a/test/sumfact/mass/mass.mini b/test/sumfact/mass/mass.mini
index 44870439291a6111c18353f2c80e9be2116d3b28..8cb8cb8c83a52ecd42fb8cfa19cee5471b26624e 100644
--- a/test/sumfact/mass/mass.mini
+++ b/test/sumfact/mass/mass.mini
@@ -12,7 +12,7 @@ printmatrix = 1
 name = {__name}
 extension = vtu
 
-[formcompiler]
+[formcompiler.operator]
 numerical_jacobian = 1, 0 | expand num
 vectorization_quadloop = 1, 0 | expand vec
 sumfact = 1
diff --git a/test/sumfact/mass/mass.ufl b/test/sumfact/mass/mass.ufl
index 6434e09f62b28ca9ea003c936997233d31449431..c11e65676418c9472cb77800dd9c8d477ebf7024 100644
--- a/test/sumfact/mass/mass.ufl
+++ b/test/sumfact/mass/mass.ufl
@@ -6,5 +6,3 @@ u = TrialFunction(V)
 v = TestFunction(V)
 
 r = u * v * dx
-
-forms = [r]
diff --git a/test/sumfact/mass/mass_3d.mini b/test/sumfact/mass/mass_3d.mini
index aba93533768a7b5463052c6585f6e45648380b70..1bde4e8cb0c981de2cdc89446f82df0a8e0c3d27 100644
--- a/test/sumfact/mass/mass_3d.mini
+++ b/test/sumfact/mass/mass_3d.mini
@@ -13,7 +13,7 @@ printmatrix = true
 name = {__name}
 extension = vtu
 
-[formcompiler]
+[formcompiler.operator]
 numerical_jacobian = 1, 0 | expand num
 vectorization_quadloop = 1, 0 | expand vec
 sumfact = 1
diff --git a/test/sumfact/mass/mass_3d.ufl b/test/sumfact/mass/mass_3d.ufl
index 5f55103e52f0b84c550e38f68c0acc8d77465793..1336a91db188e0626353294281f109e0d30ad8ae 100644
--- a/test/sumfact/mass/mass_3d.ufl
+++ b/test/sumfact/mass/mass_3d.ufl
@@ -6,5 +6,3 @@ u = TrialFunction(V)
 v = TestFunction(V)
 
 r = u * v * dx
-
-forms = [r]
diff --git a/test/sumfact/mass/sliced.mini b/test/sumfact/mass/sliced.mini
index 90dab43e70b8ddc38830c37afb2dd83b4116f5e7..542548d72ae2bcc2286340f4aeef71f2cbcfb648 100644
--- a/test/sumfact/mass/sliced.mini
+++ b/test/sumfact/mass/sliced.mini
@@ -9,7 +9,7 @@ printmatrix = true
 name = {__name}
 extension = vtu
 
-[formcompiler]
+[formcompiler.operator]
 numerical_jacobian = 1
 vectorization_strategy = explicit
 vectorization_horizontal = 1
diff --git a/test/sumfact/poisson/diagonal.mini b/test/sumfact/poisson/diagonal.mini
index d3744184c52abd1320aa796bc16249c478afe9a7..58abb8c603838a31668a8cc8332bb59f743cbeb8 100644
--- a/test/sumfact/poisson/diagonal.mini
+++ b/test/sumfact/poisson/diagonal.mini
@@ -8,8 +8,10 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
-sumfact = 1
 compare_l2errorsquared = 1e-5
+
+[formcompiler.operator]
+sumfact = 1
 vectorization_quadloop = 1
 vectorization_strategy = explicit
 vectorization_horizontal = 2
diff --git a/test/sumfact/poisson/opcount_poisson_2d_order2.mini b/test/sumfact/poisson/opcount_poisson_2d_order2.mini
index 2350f1137da3df811df6bfa56bb6b373c0b422ab..cc5de56f45a320ca43397be3078622aadd6f9c16 100644
--- a/test/sumfact/poisson/opcount_poisson_2d_order2.mini
+++ b/test/sumfact/poisson/opcount_poisson_2d_order2.mini
@@ -12,11 +12,12 @@ reference = poisson_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 0
 compare_l2errorsquared = 1e-8
-sumfact = 1
 opcounter = 1
 instrumentation_level = 4
 
+[formcompiler.operator]
+sumfact = 1
+
 [formcompiler.ufl_variants]
 degree = 2
diff --git a/test/sumfact/poisson/opcount_sumfact_poisson_dg_2d_vec.mini b/test/sumfact/poisson/opcount_sumfact_poisson_dg_2d_vec.mini
index 063987b96cfb04af6223b426258e412f839582e7..02039a9fbbc9ff213c9d09073f07ba0004ffbde5 100644
--- a/test/sumfact/poisson/opcount_sumfact_poisson_dg_2d_vec.mini
+++ b/test/sumfact/poisson/opcount_sumfact_poisson_dg_2d_vec.mini
@@ -10,11 +10,11 @@ reference = poisson_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 0
-sumfact = 1
 opcounter = 1
 instrumentation_level = 4
 
+[formcompiler.operator]
+sumfact = 1
 
 [formcompiler.ufl_variants]
 degree = 1
diff --git a/test/sumfact/poisson/poisson_2d.mini b/test/sumfact/poisson/poisson_2d.mini
index 9fab490cf7b12a362767520975173927971e2382..a079932bb92ad4812091e29be68acde9cacf1c3e 100644
--- a/test/sumfact/poisson/poisson_2d.mini
+++ b/test/sumfact/poisson/poisson_2d.mini
@@ -14,8 +14,10 @@ reference = poisson_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
 compare_l2errorsquared = 4e-5, 4e-9 | expand deg
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
 sumfact = 1
 vectorization_strategy = explicit, none | expand grad
 quadrature_order = 2, 4
diff --git a/test/sumfact/poisson/poisson_2d.ufl b/test/sumfact/poisson/poisson_2d.ufl
index 97cc99c407283e8e65bda1fdbb98921c96d6a50b..f0cecc18e42eb5ee9b7bc4b488bb8e926f34ea8b 100644
--- a/test/sumfact/poisson/poisson_2d.ufl
+++ b/test/sumfact/poisson/poisson_2d.ufl
@@ -12,7 +12,7 @@ V = TensorProductElement(V_0, V_1, cell=cell)
 u = TrialFunction(V)
 v = TestFunction(V)
 
-forms = [(inner(grad(u), grad(v)) - f*v)*dx]
+r = (inner(grad(u), grad(v)) - f*v)*dx
 exact_solution = g
 is_dirichlet = 1
 interpolate_expression = g
diff --git a/test/sumfact/poisson/poisson_3d.mini b/test/sumfact/poisson/poisson_3d.mini
index 2ddbd626bac79456c3b138e0386a8ac94ee9aa15..31f556f56fd550c252c369d0ea42b7a24e13fd34 100644
--- a/test/sumfact/poisson/poisson_3d.mini
+++ b/test/sumfact/poisson/poisson_3d.mini
@@ -15,8 +15,10 @@ reference = poisson_ref
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 1, 0 | expand num
 compare_l2errorsquared = 1e-4, 1e-8 | expand deg
+
+[formcompiler.operator]
+numerical_jacobian = 1, 0 | expand num
 sumfact = 1
 vectorization_quadloop = 1, 0 | expand quad
 vectorization_strategy = explicit, none | expand grad
diff --git a/test/sumfact/poisson/poisson_3d.ufl b/test/sumfact/poisson/poisson_3d.ufl
index 75f0331a54427d11bcf86dc9613afde38847f881..313cec8ec572013d5604d3a9dde332c13f359b3e 100644
--- a/test/sumfact/poisson/poisson_3d.ufl
+++ b/test/sumfact/poisson/poisson_3d.ufl
@@ -9,7 +9,7 @@ V = FiniteElement("CG", cell, degree)
 u = TrialFunction(V)
 v = TestFunction(V)
 
-forms = [(inner(grad(u), grad(v)) - f*v)*dx]
+r = (inner(grad(u), grad(v)) - f*v)*dx
 exact_solution = g
 is_dirichlet = 1
 interpolate_expression = g
diff --git a/test/sumfact/poisson/poisson_dg_2d.mini b/test/sumfact/poisson/poisson_dg_2d.mini
index 99adc0e31563c6a85ecad18e349ed74adcaf21a0..3b4fd9a55c33e03417705e3fff0729031913dd30 100644
--- a/test/sumfact/poisson/poisson_dg_2d.mini
+++ b/test/sumfact/poisson/poisson_dg_2d.mini
@@ -14,9 +14,11 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
+compare_l2errorsquared = 5e-5, 5e-7 | expand deg
+
+[formcompiler.operator]
 numerical_jacobian = 1, 0 | expand num
 sumfact = 1
-compare_l2errorsquared = 5e-5, 5e-7 | expand deg
 vectorization_quadloop = 1, 0 | expand quad
 vectorization_strategy = explicit, none | expand grad
 
diff --git a/test/sumfact/poisson/poisson_dg_2d.ufl b/test/sumfact/poisson/poisson_dg_2d.ufl
index fefc67d64c9a0b1a6dac9c34423933d149e8c88e..1da5733c19e48a96e931948c2ed5e46214c068d7 100644
--- a/test/sumfact/poisson/poisson_dg_2d.ufl
+++ b/test/sumfact/poisson/poisson_dg_2d.ufl
@@ -29,5 +29,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/sumfact/poisson/poisson_dg_3d.mini b/test/sumfact/poisson/poisson_dg_3d.mini
index b23fda0eba605025076e6a92dfb295ea06525d00..58410c459b76841d6169edf63fc9cfacc61599b1 100644
--- a/test/sumfact/poisson/poisson_dg_3d.mini
+++ b/test/sumfact/poisson/poisson_dg_3d.mini
@@ -14,9 +14,11 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
+compare_l2errorsquared = 1e-4, 5e-6 | expand deg
+
+[formcompiler.operator]
 numerical_jacobian = 1, 0 | expand num
 sumfact = 1
-compare_l2errorsquared = 1e-4, 5e-6 | expand deg
 vectorization_quadloop = 1, 0 | expand quad
 vectorization_strategy = explicit, none | expand grad
 
diff --git a/test/sumfact/poisson/poisson_dg_3d.ufl b/test/sumfact/poisson/poisson_dg_3d.ufl
index 80d78c363b27b6e30476a1189aed521b33b65fa8..20705f2f58a0ba1cdb66dfdcf01f1f0835612ebf 100644
--- a/test/sumfact/poisson/poisson_dg_3d.ufl
+++ b/test/sumfact/poisson/poisson_dg_3d.ufl
@@ -29,5 +29,4 @@ r = inner(grad(u), grad(v))*dx \
   - theta*g*inner(grad(v), n)*ds \
   - gamma*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/sumfact/poisson/poisson_dg_tensor.mini b/test/sumfact/poisson/poisson_dg_tensor.mini
index 4a45e4a1fa469c6fc2753ad230b486bc55bd1b55..aa485d22e7351e5a2cc95bacbaee9984076ad6d0 100644
--- a/test/sumfact/poisson/poisson_dg_tensor.mini
+++ b/test/sumfact/poisson/poisson_dg_tensor.mini
@@ -12,8 +12,10 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
-sumfact = 1
 compare_l2errorsquared = 3e-4
+
+[formcompiler.operator]
+sumfact = 1
 vectorization_quadloop = 1, 0 | expand quad
 vectorization_strategy = explicit, none | expand grad
 
diff --git a/test/sumfact/poisson/poisson_dg_tensor.ufl b/test/sumfact/poisson/poisson_dg_tensor.ufl
index 0d4b7a79ee8e2b7183b8664bb639990960205e0a..8388ed75b40f2108ba56076e974404267d2c7310 100644
--- a/test/sumfact/poisson/poisson_dg_tensor.ufl
+++ b/test/sumfact/poisson/poisson_dg_tensor.ufl
@@ -32,5 +32,4 @@ r = (inner(A*grad(u), grad(v)) + c*u*v)*dx \
   - theta*g*inner(A*grad(v), n)*ds \
   - gamma*g*v*ds
 
-forms = [r]
 exact_solution = g
diff --git a/test/sumfact/poisson/poisson_fastdg_2d.mini b/test/sumfact/poisson/poisson_fastdg_2d.mini
index 541de8712b627f6327076d87586b6d378fa54e78..abd5b214f983243728cfc98f958278b65b1be5a1 100644
--- a/test/sumfact/poisson/poisson_fastdg_2d.mini
+++ b/test/sumfact/poisson/poisson_fastdg_2d.mini
@@ -12,9 +12,11 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
+compare_l2errorsquared = 1e-4
+
+[formcompiler.operator]
 numerical_jacobian = 0
 sumfact = 1
-compare_l2errorsquared = 1e-4
 vectorization_quadloop = 1, 0 | expand quadvec
 vectorization_strategy = explicit, none | expand gradvec
 fastdg = 1
diff --git a/test/sumfact/poisson/poisson_fastdg_3d.mini b/test/sumfact/poisson/poisson_fastdg_3d.mini
index b5974a4fe62609ca038d315656bf3c9f333f3e32..e58d07852e57ddeb2d4ef4f0d5be1b02384030ee 100644
--- a/test/sumfact/poisson/poisson_fastdg_3d.mini
+++ b/test/sumfact/poisson/poisson_fastdg_3d.mini
@@ -12,9 +12,11 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
+compare_l2errorsquared = 1e-4
+
+[formcompiler.operator]
 numerical_jacobian = 0
 sumfact = 1
-compare_l2errorsquared = 1e-4
 vectorization_quadloop = 1, 0 | expand quadvec
 vectorization_strategy = explicit, none | expand gradvec
 fastdg = 1
diff --git a/test/sumfact/poisson/sliced.mini b/test/sumfact/poisson/sliced.mini
index 858b8c6b6b8804f3cede25236ff29ce66bae010b..26652161e2751ea54ed82e326e833f12ae57526d 100644
--- a/test/sumfact/poisson/sliced.mini
+++ b/test/sumfact/poisson/sliced.mini
@@ -8,8 +8,10 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
-sumfact = 1
 compare_l2errorsquared = 1e-5
+
+[formcompiler.operator]
+sumfact = 1
 vectorization_quadloop = 1
 vectorization_strategy = explicit
 vectorization_horizontal = 1
diff --git a/test/sumfact/stokes/stokes.mini b/test/sumfact/stokes/stokes.mini
index 10dca54704f30369881095216a6f5f94c25585e0..203e89e28eaea013e5a5868c9f576f489d98d09f 100644
--- a/test/sumfact/stokes/stokes.mini
+++ b/test/sumfact/stokes/stokes.mini
@@ -12,7 +12,9 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
+compare_l2errorsquared = 1e-12
+
+[formcompiler.operator]
 numerical_jacobian = 1, 0 | expand num
 vectorization_quadloop = 1, 0 | expand quad
-compare_l2errorsquared = 1e-12
 sumfact = 1
diff --git a/test/sumfact/stokes/stokes.ufl b/test/sumfact/stokes/stokes.ufl
index 1286a48d91b8807534345385f8f0b796c72d8012..9c5cb27a59d7662402a0bc6d65f9818dc4761577 100644
--- a/test/sumfact/stokes/stokes.ufl
+++ b/test/sumfact/stokes/stokes.ufl
@@ -14,7 +14,6 @@ u, p = TrialFunctions(TH)
 
 r = (inner(grad(v), grad(u)) - div(v)*p - q*div(u))*dx
 
-forms = [r]
 exact_solution = g_v, 8.*(1.-x[0])
 interpolate_expression = g_v, None
 is_dirichlet = v_bctype, v_bctype, None
\ No newline at end of file
diff --git a/test/sumfact/stokes/stokes_3d_dg.mini b/test/sumfact/stokes/stokes_3d_dg.mini
index b7ec60614d3159d8792b42ab1c416c4926150120..5f4a2bac06330ee40a500e1d04c5656b376bcbfa 100644
--- a/test/sumfact/stokes/stokes_3d_dg.mini
+++ b/test/sumfact/stokes/stokes_3d_dg.mini
@@ -11,8 +11,10 @@ printmatrix = false
 name = {__name}
 extension = vtu
 
+[formcompiler]
+compare_l2errorsquared = 1e-10
+
 [formcompiler]
 numerical_jacobian = 0
 sumfact = 1
-fastdg = 1, 0 | expand fastdg
-compare_l2errorsquared = 1e-10
\ No newline at end of file
+fastdg = 1, 0 | expand fastdg
\ No newline at end of file
diff --git a/test/sumfact/stokes/stokes_3d_dg.ufl b/test/sumfact/stokes/stokes_3d_dg.ufl
index 84d1003e16d7f4f36dc0630434d98eb3d633cd3a..3fefa9b763714ca27b7fa7d741bfe94059dae836 100644
--- a/test/sumfact/stokes/stokes_3d_dg.ufl
+++ b/test/sumfact/stokes/stokes_3d_dg.ufl
@@ -32,5 +32,4 @@ r = inner(grad(u), grad(v))*dx \
   + p*inner(v, n)*ds \
   + q*inner(u-g_v, n)*ds
 
-forms = [r]
 exact_solution = g_v, 8.*(1.-x[0])
diff --git a/test/sumfact/stokes/stokes_dg.mini b/test/sumfact/stokes/stokes_dg.mini
index e3374e4a18e844f6f1356ce45b38e5b5212f015f..9b5202a1a39221aee536f20acd2819ec5d93f9e5 100644
--- a/test/sumfact/stokes/stokes_dg.mini
+++ b/test/sumfact/stokes/stokes_dg.mini
@@ -13,9 +13,11 @@ name = {__name}
 extension = vtu
 
 [formcompiler]
-numerical_jacobian = 0, 1 | expand num
 compare_l2errorsquared = 1e-8
+
+[formcompiler.operator]
+numerical_jacobian = 0, 1 | expand num
 sumfact = 1
 fastdg = 1, 0 | expand fastdg
 
-{formcompiler.fastdg} == 1 and {formcompiler.numerical_jacobian} == 1 | exclude
\ No newline at end of file
+{formcompiler.operator.fastdg} == 1 and {formcompiler.numerical_jacobian} == 1 | exclude
\ No newline at end of file
diff --git a/test/sumfact/stokes/stokes_dg.ufl b/test/sumfact/stokes/stokes_dg.ufl
index 39c243a00c7b16c857f551813b6bf0f4a99ad065..ddefd870a413fe7099150e619d336c6d089da5d6 100644
--- a/test/sumfact/stokes/stokes_dg.ufl
+++ b/test/sumfact/stokes/stokes_dg.ufl
@@ -32,5 +32,4 @@ r = inner(grad(u), grad(v))*dx \
   + p*inner(v, n)*ds \
   + q*inner(u-g_v, n)*ds
 
-forms = [r]
 exact_solution = g_v, 8*(1.-x[0])
\ No newline at end of file