I am trying to solve the axisymmetric stationary heat transfer problem of a wire having a positive constant internal heat source term and flux boundary conditions (BC) in the top, low and lateral surfaces using the pde toolbox. I set the BCs in the top and low surfaces simply as convective outgoing heat fluxes, and the BC at the lateral surface as convective plus radiative outgoing heat fluxes. When analyzing my results I encountered an unexpected behavior. When I set the infrared emissivity coefficient of the lateral surface equal to one, the temperature values at the surface of the wire are unchanged or higher than when I make it equal to zero. Emissivity equal to one means more heat loses through the lateral surface of the wire, therefore I would expect lower stationary temperature values than when I set the emissivity equal to zero (no radiative heat fluxes being considered).
I tried to use both the createpde() and the createpde('thermal') functions to define my thermal model, with both giving me unexpected results. When using the generic pde model definition (i.e. createpde()), I was careful to define the boundary conditions in such a way that the solver treats the problem as nonlinear, following the suggestions given in an older thread. However, even if it seems that the solver is considering the problem as nonlinear, it is dismissing the nonlinear part of my radiative heat flux definition, which results in unchanged temperature values with and without the radiative flux definition. With the createpde('thermal') function, there is no much one can play with, I simply changed the emissivity parameter in the thermalBC function and found out that the temperature with the radiative BC was higher than without it. I managed to get more realistic (expected) results by defining the problem as 'transient' and not as 'steady-state'. This requires that I define an evaluation time long enough so that the temperature reaches the stationary state. However, for computation time reasons, I would appreciate if the nonlinear stationary solver could also give a good solution to the problem.
I am using MATLAB R2017a. In the attachment I joined the AXheatTransferPDE.m file, containing the code I used for defining stationary and transient problems using the createpde() and the createpde('thermal') functions. The .m file loads the workspace workspacePDE.mat, also attached. Besides, I added .fig files with the results I got from my PC.
Thanks for any explanations about this issue.
