_apply_stiffness_transpose_operator: uses einsum in a way that causes high cost

[inducer]

Since our einsum is not smart (it doesn't identify possible common subexpressions), this:

grudge/grudge/op.py

Lines 296 to 306 in 1c44c4b

	actx.einsum("dij,ej,ej,dej->ei",
	reference_stiffness_transpose_matrix(
	actx,
	out_element_group=out_grp,
	in_element_group=in_grp
	),
	ae_i,
	vec_i,
	inv_jac_t_i,
	arg_names=("ref_stiffT_mat", "jac", "vec", "inv_jac_t"),
	tagged=(FirstAxisIsElementsTag(),))

amounts to computing the rst derivatives once per dimension. (Together with #141, we're computing them 9 times in 3D, when once would suffice.) The problem is that the sum is computed once per d, when it is actually independent of d.

This is something I should have caught during review of #74.

cc @thomasgibson @lukeolson

[inducer]

Again, not totally sure what I was smoking here. For general curvi, the inverse metric derivatives depend on the point, so the best one could conceivably do is precompute the sum over d, but that would lead to a matrix-per-element size temporary, which is entirely unaffordable.

More can be done for affine, but that's the domain of #161.