We have used SU2 and SU2-NEMO to predict the heat flux for various shapes in a high-speed flow. Something we have observed is that SU2 requires an extremely small first cell spacing to get an accurate and grid-independent prediction for heat flux. For example, this first cell space is often 1E-9 to 1E-10 m with SU2; whereas with other codes it is a micrometer. The latter is in line with the general rule of thumb of setting the first cell height so that the first cell Reynold's number is between 0.1 and 1.0.

Discussing with @jtneedels and @WallyMaier, we think it may have to do with how the temperature gradient is calculated in CNSSolver.cpp:

/--- Compute the normal gradient in temperature using Twall ---/

su2double dTdn = -(There - Twall)/dist_ij;

/--- Apply a weak boundary condition for the energy equation.
Compute the residual due to the prescribed heat flux. ---/

su2double Res_Conv = 0.0; su2double Res_Visc = thermal_conductivity * dTdn * Area;

which may result in a first-order approximation for the temperature gradient causing this requirement for an extremely small cell size to recover the true gradient. I would be interested to see if calculating this gradient similar to CFVMFlowSolverBase.ini. For example,
for (iDim = 0; iDim < nDim; iDim++) { Grad_Temp[iDim] = nodes->GetGradient_Primitive(iPoint, prim_idx.Temperature(), iDim); }

su2double dTdn = -GeometryToolbox::DotProduct(nDim, Grad_Temp, UnitNormal);

I tried checking this but quickly found out, I am not familiar enough with the code structure to implement.