[Daniel_Khazaei]

Originally Posted by gridley2

Thanks very much for the explanation Daniel! This makes sense. So to obtain the true divergence, I have to multiply "mygrad" by face areas first? To be honest, I'm surprised that ::div doesn't automatically do this. I wonder why not. I cannot thank you enough for helping me figure this out... I was certainly scratching my head over this for a few hours this past week!

well, what I think about and I'm not sure is:
surfaceScalarFields are variables defined on the faces of each cell and it should be used for variables that actually have physical meaning when defined on the surface, e.g. mass flux(rho*U*A)

fvc::div(anything) will give a result as long as you follow its rules and the provided argument doesn't violates any mathematical or C++ methods defined for it!

Whether the given argument and its divergence actually have any physical meaning or not is a different story and of-course non of fvc::div() business.