Models with viscous heating

If you put parameter alpha_vis to something >0 in Parameter.in, the gas will be heated by viscous energy dissipation.

1.E-2    ! alpha_vis

The formula for the local gas heating rate applied [erg/cm3/s] is given by Eqs.(104) and (105), and text, in the first ProDiMo-paper (Woitke et al.2009). It scales as

Gamma_vis = 9/4 alpha press Omega

where alpha~0.01 is the viscosity parameter, press the gas pressure, and Omega the Keplerian angular velocity.

There is a big conceptional problem with this description of viscous heating, which in unsolved in ProDiMo, as follows. This heating scales roughly as ~rho, whereas, in the upper layers, all cooling processes scale as rho2, which means that for low densities the viscous heating always dominates, cannot be stopped by any cooling, so using alpha>0 will result in unbound gas temperatures in all tenuous layers.

There is another parameter called dust_nonRE in ProDiMo, which you can set to .true. in Parameter.in. In that case, the thermal accomodation cooling rate of the gas (energy transfer gas=>dust [erg/cm3/s]) will be stored in heatcoolgas_it.f, and used in the next global iteration to heat the dust, additionally to the radiative heating, in function FUNCTION RADEQUIL(Tdold,Td,ix,iz) in tdust_from_RT.f. Thus there is an additional heating (or cooling!) process for the dust in the disc, the dust is not anymore assumed to be in radiative equilibrium, therefore the name of that switch dust_nonRE.

If you put dust_nonRE=.true. (also possible with alpha=0 of course), ProDiMo needs to do global iterations, because the dust non-radiative heating can only be calculated when Tgas has been calculated. Thus, even for given disc density structure, ProDiMo will start iterating until the change in total (volume-integrated) dust non-RE heating rate [erg/s] becomes small, i.e. <5%.

These two quantities, the total viscous heating rate Lvisc [erg/s] and the total non-radiative dust heating rate LnonRE [erg/s], are calculated, and written to stdout, in routine dust_heating.f. Only a part of Lvisc will be transferred to the dust, because in the more optically thin regions, the gas will simply radiate away that excess energy, so one should expect LnonRELvisc, because there are other processes that lead to Tgas>Tdust, causing non-radiative dust heating.

I find it extremely hard to explain near-IR excess with this viscous heating, because (i) Lvisc is often too low, even for alpha=0.1 and Mdisk=0.1 Msun, (ii) only a part of Lvisc ends up heating the dust, (iii) for epsilon=1.0, must of the viscous energy is created in the outer disk regions, resulting in a far-IR (non near-IR) excess. An alternative seems to be ChemHeatFac>0 in combination with dust_nonRE=.true., which, because of rho2-behavior, does produce a near-IR excess.

Another option is to specify the mass accretion rate

1.E-7     ! Mdot [Msun/yr]

In that case we apply Eq.(2) of D'Alessio et al. (1998), stating the half column (z=0 in oo) heating rate [erg/cm2/s] as

F_vis(r) = 3 G Mstar Mdot / (8pi r3) * [1-sqrt(Rstar/r)] ,

assuming that Mdot is a constant throughout the disk, and assuming that (a certain part of) the gravitational energy released when the disk "shrinks" is converted into heat. To turn this heating rate per column into a heating rate per volume, we need to make an additional assumption about how that total heating rate is distributed within the column as function of height z. We choose

Gamma_vis(r,z) = F_vis(r) * rhop(r,z) / \int rhop(r,z') dz'

with p=2 in order to avoid the strange effects caused in the uppermost layers.

If you use, in addition, MCMAX_LIKE = .true., this viscous heating is applied to both the gas and the dust. That option does produce a near-IR excess.