Couple subdomains with interiour boundary condition

Question

Dear Fenics developers,
I have asked before, how to couple different PDEs in two separate subdomains ( https://fenicsproject.org/qa/11060/fenics-handle-individual-coupled-regions-conjugate-transfer ), unfortunately with no reply. The aim is e.g. to solve turbulent, non-isothermal flow (subdomain 1) in thermal contact with a heat conducting solid (subdomain 2) with the interior boundary condition that the temperature and the heat flux on both sides of the boundary are the same.

Now I have found a similar question, that could help to use a variable on subdomain 1 as boundary conditions on a subdomain 2 on the common interior boundary:
https://fenicsproject.org/qa/2532/how-to-set-up-interior-neumann-boundaries and
https://fenicsproject.org/qa/2744/how-to-make-interior-neumann-bc-dependent-on-fe-function

Unfortunately, the second link tells that u, which I'd need to couple values from subdomain 1 to 2, must be a Function, not a TrialFunction. However, u must be solved on subdomain 1 and thus has to be a TrialFunction.

Is it possible to form a single equation system from one PDE on subdomain 1, another PDE on subdomain 2 and a linear relation for the solution and its derivatives on the interior boundary?

Comments

I'm not 100% sure about the implications of this, but as you haven't got any answers so far I'll comment on a part of your question.
$u$ can be a Function instead of a TrialFunction, if handled in the following manner (like done for nonlinear forms):

u = Function(V)
v = TestFunction(V)
a = inner(grad(u), grad(v))*dx
J = derivative(a, u)
A = assemble(J)

Should give the same as if $u$ was a TrialFunction and the derivative was omitted.