Hello together,

the new timesteppers are great. But wouldn’t it be nice if the `u.dt`

term evaluates to zero in case of classical use of bilinear form and solver.

Consider e.g. following problem:

```
from ngsolve import *
from ngsolve.webgui import Draw
from netgen.occ import *
shape = Rectangle(1,1).Face()
shape.edges.Max(X).name="right"
shape.edges.Min(X).name="left"
shape.edges.Max(Y).name="top"
shape.edges.Min(Y).name="bot"
mesh = Mesh(OCCGeometry(shape, dim=2).GenerateMesh(maxh=0.05))
fes = H1(mesh, order=4, dirichlet="left|right")
u,v = fes.TnT()
eq = ( grad(u) * grad(v) + u.dt * v) * dx
gfu = GridFunction(fes)
def init_gfu(gfu):
t_left = 10
t_right = 20
t_ref = 15
gfu.Interpolate(15)
gfu.Interpolate(mesh.BoundaryCF({"left": t_left, "right": t_right}, default=t_ref),definedon=mesh.Boundaries("left|right"))
```

works as expected in case of timestepper usage:

```
init_gfu(gfu)
ts = timestepping.ImplicitEuler(equation=eq, dt=1e-1)
scene = Draw(gfu)
def callback(t, gfu):
scene.Redraw()
ts.Integrate(gfu, start_time=0, end_time=2, callback=callback)
```

But if we plug `eq`

directly into a bfa and solve the same problem, the `u.dt`

evaluates to `u`

and we get a different result instead of the quasi-static solution:

```
init_gfu(gfu)
a = BilinearForm(fes, symmetric=True)
a += eq.Compile()
solvers.Newton(a=a, u=gfu, inverse="sparsecholesky", printing=False)
scene = Draw (gfu)
import matplotlib.pyplot as plt
import numpy as np
X = np.linspace(0, 1, num=100)
Y = np.ones_like(X) * 0
plt.plot(X, gfu(mesh(X, Y)))
plt.show()
```