Internal boundaries for thermal conductivity

[andLaing]

I'm currently developing a module to calculate thermal transmittance. The objects that we're modelling are multimaterial with internal gas spaces where there needs to be temperature dependent recalculation of conductivity due to radiation and convection effects. Thanks to the community I've got a simplified implementation using the Newton solver and GridFunctionCoefficient which seems to do what expected. However, if I want to generalize the implementation I think I need to take into account the boundary between solid and gaseous elements rather than putting all the information about emmissivity in the coefficient since it won't be generally true that an element represents the whole gas section.

Reading a bit through the documentation, I think I need to add an interior face integrator to the nonlinearform but I'm not sure and haven't found many examples on how to implement that functionality. Also, there doesn't seem to be a way to limit its use to certain boundaries as with boundary integrators (at least in mfem 4.7), How would I limit the action? Or am I barking up the wrong tree?

I've attached a snippet of what I'm doing at the moment:

mfem::ConstantCoefficient timeVar(0.);
mfem::Vector intExtTransCoeff({ m_thermalParameters.m_intHeatCoeff / 1e6, m_thermalParameters.m_extHeatCoeff / 1e6 });
mfem::Vector transCoeffTemp({ internalT, externalT });
transCoeffTemp *= intExtTransCoeff;
timeDep->AddDomainIntegrator(new mfem::DomainLFIntegrator(timeVar));

mfem::Array<int> rob_bdr(mMesh.bdr_attributes.Max());
rob_bdr = 1;
mfem::PWConstCoefficient PWTransCoeffTemp(transCoeffTemp);
timeDep->AddBoundaryIntegrator(new mfem::BoundaryLFIntegrator(PWTransCoeffTemp), rob_bdr);
timeDep->Assemble();

mfem::Vector epsilon({ 0., m_thermalParameters.m_emissivity });
ConductionCoefficient conductCoeff(tempGrid.get(), conductV, epsilon, m_thermalParameters.m_boltzmann / 1e6);
conduct->AddDomainIntegrator(new ThermalDiffusionIntegrator(conductCoeff));
mfem::PWConstCoefficient PWTransCoeff(intExtTransCoeff);
conduct->AddBoundaryIntegrator(new mfem::MassIntegrator(PWTransCoeff), rob_bdr);
mfem::Array<int> ess_tdof_list(0);
//

ThermalOperator thOper(conduct.get(), fspace->GetVSize());
mfem::Vector B, X;
tempGrid->GetTrueDofs(X);
B.MakeRef(*timeDep, 0, timeDep->Size());

mfem::NewtonSolver newtSolver;
auto J_solver        = new mfem::MINRESSolver;
mfem::Solver* J_prec = new mfem::DSmoother(1);
J_solver->SetPreconditioner(*J_prec);
J_solver->SetRelTol(1e-12);
J_solver->SetAbsTol(0.);
J_solver->SetMaxIter(100);

newtSolver.iterative_mode = false;
newtSolver.SetSolver(*J_solver);
newtSolver.SetOperator(thOper);
newtSolver.SetPrintLevel(1);
newtSolver.Mult(B, X);
//
tempGrid->SetFromTrueDofs(X);

The ThermalDiffusionIntegrator inherits from DiffusionIntegrator and has an override of the AssembleElementGrad method. Am I right in assuming what I need to add is conduct->AddInteriorFaceIntegrator line? Are there examples of the minimum requirement for implementation?