[aspect-devel] questions regarding accuracy of no-heat flux boundary condition and help needed for the right entropy viscosity parameters
Wolfgang Bangerth
bangerth at tamu.edu
Mon Jul 11 12:34:45 PDT 2016
>>> What is the size of the heat flux through the left/right and how does
>>> it compare to the top/bottom heat flux? If this is non-zero but
>>> several orders of magnitude smaller, then that is as good as it can
>>> get. If not, then I don't know.
>> Hi Timo, thanks for you response. Yes, the left and right flux (O(1))
>> are comparable small to the top and bottom (O(10^4)). However, I still
>> think if we assigned the zero non-flux boundary, then it as least should
>> be in some lower order accuracy near zero, say (10^{-1})... But perhaps
>> you are right, considering the error from the matrix iteration solver.
As well as the discretization. In the finite element context, the heat
flux is not zero in a pointwise sense, but only in the sense of some
weighted integral. A heat flux of ~1 on the sides is relatively small
compared to the heat flux of ~10^4 at the top/bottom boundaries. If you
want to know whether this is a systematic error or just a discretization
artifact, you can run on a sequence of finer and finer meshes. The
top/bottom heat flux should converge to something that does not change
with the mesh size, whereas the heat flux at the left/right should go to
zero with mesh refinement.
>> Also I am more concerned about the right formula to compute the heat
>> flux in the statistics postprocessor. Numerically, we are solving
>> transport equation
>>
>> \rho dT/dt - \grad\cdot (k \grad T) + \grad\cdot (\nu \grad T) + v\cdot
>> \grad T = f
>>
>> where \grad\cdot (\nu \grad T) comes from the artificial entropy
>> viscosity. Therefore, in the weak formulation, we assume
>>
>> (-k+\nu) dT/dn = 0
>>
> I just realized that the correctly formula should be
>
> \rho dT/dt - \grad\cdot (k \grad T) - \grad\cdot (\nu \grad T) + v\cdot
> \grad T = f
>
> therefore we imposed
>
> -(k+\nu) dT/dn = 0
>
> and we only assume k and \nu positive numbers, so it could be the
> magnitude of (k +\nu) affecting the final flux calculation. I'll check
> the value dT/dn instead.
Conceptually, the boundary conditions we impose is indeed
-(k+\nu) dT/dn = 0
In practice, you of course only want
-k dT/dn = 0
but that is not what is imposed. If you want, the boundary condition you
get is really
-k dT/dn = nu dT/dn
and so the total heat flux across the boundary is not going to be zero
but proportional to the artificial viscosity. However, because nu goes
to zero as the mesh is refined, you will converge to the correct
solution. This is all you can probably hope.
>> which is the no-heat flux condition. However, in the statistics
>> postprocessor, we only compute (-k) dT/dn, so perhaps this is the
>> reason why we couldn't have a zero heat flux?
Yes, this seems correct.
Best
W.
--
------------------------------------------------------------------------
Wolfgang Bangerth email: bangerth at math.tamu.edu
www: http://www.math.tamu.edu/~bangerth/
More information about the Aspect-devel
mailing list