# Heat balance multilayer with Neumann condition

I am solving the heat equation (diffusion only) over a 3-layer system (physical properties varying from layer to layer).

``km = 0.128; rhom = 925; Cpm = 1550;  ka = 0.024; rhoa = 1.292; Cpa = 1003;  thickness = 2*0.136;    L1 = N[-thickness/2]; L2 = N[thickness/2]; T01 = 120; T02 = 20; T03 = 120;  Length1 = 0.5;   v= 40/60;  time1 = Length1/v; ``

For describing the properties change in space, I use the following equation:

``slope = 1000; SmoothedStepFunction[fL_, fmax_, fR_, tsL_, tsR_, m_] :=   Function[t, (fL*Exp[tsL*m] + fmax*Exp[m*t])/(Exp[tsL*m] +       Exp[m*t]) - (fR*Exp[tsR*m] + fmax*Exp[m*t])/(Exp[tsR*m] +       Exp[m*t]) + fR];  rhoCp[x_] :=    SmoothedStepFunction[rhoa*Cpa, rhom*Cpm, rhoa*Cpa, L1, L2, slope][x]; k[x_] := SmoothedStepFunction[10^6*ka, 10^6*km, 10^6*ka, L1, L2,      slope][x]; ``

The actual heat balance is solved here:

``heateq = rhoCp[x]*D[u[x, t], t] ==     Inactive[Div][{{k[x]}}.Inactive[Grad][u[x, t], {x}], {x}];  ic[x_] :=    Piecewise[{{T01, x < L1}, {T02, L1 <= x <= L2}, {T03, x > L2}}];  sol1 = First[    NDSolve[{heateq, u[x, 0] == ic[x]},      u, {x, -2*thickness, 2*thickness}, {t, 0, time1},      Method -> {"MethodOfLines",        "SpatialDiscretization" -> {"FiniteElement",          "MeshOptions" -> {"MaxCellMeasure" -> 0.01}}}]]; ``

Despite it works well, I still have a question: next to the Neumann boundary conditions (`Inactive[Div][{{k[x]}}.Inactive[Grad][u[x, t]`) that are applied over the interfaces to connect the solutions in each domain, are there other boundary conditions applied implicitly?