Indeed, derivative doesn’t work on a tuple.
The system you are solving is block shaped, and derivative should be used per block, w.r.t to only one variable.
When using MixedNonlinearVariationalProblem, both the form and the Jacobian are defined as lists corresponding to the decomposition per block.
To obtain the corresponding list of subforms, you can do
L_blocks = extract_blocks(L)
and derivative can then be used on each subforms to build the list of blocks for the Jacobian :
Js = []
for Li in L_blocks :
for zi in z.split():
d = derivative(Li, zi)
Js.append(d)