n_ph of dust torus about agnpy HOT 8 CLOSED

Hi Julian,

the multiplication is by 2 pi, while I think it should be by 4 pi.

I thought about it for some time and I came to the conclusion it should be 2 pi, this is the reasoning I did:

This is the formula Finke 2016 gives us for the energy density of the dust torus.
It is differential in solid angle, so we have to integrate from 0 to 2 pi in dphi (-> a factor 2 pi) and from -1 to 1 in mu (-> a factor 2).
But here we have a Dirac delta for mu, so our integration in mu vanishes, right?
I mean that the integral of a Dirac delta without any function (or with a function constant to 1) is 1.
https://www.wolframalpha.com/input/?i=int+from+-1+to+1+delta%28x-a%29+dx

Same goes for the BLR due to the other Dirac delta in mu

the integration in mu is eliminated by the Dirac delta on mu.

Did I get it right? Sorry if my math is not fresh enough.

It is basically this function multiplied by Gamma^2 (1-mu * beta)^2, do you prefer it in the target.py or in emission_region.py? It needs variables from both blob and DT.

You sure it's not cubed, the u is invariant with the energy cubed (Chapter 5 of Dermer 2009)

I thought about inserting the same function to help you with B2 1420, I thought to something like

target.u_ph(r, frame="blob")
target.u_ph(r, frame="galaxy")

but then in the first case the blob has to be passed. Or at least the quantities relevant for the transformation (Gamma, Beta).
Anyhow I think it is better to keep such function in targets.

Cheers

from agnpy.

about 1:
I still think it should be 4 pi
imagine a simple example with the blob sitting in the middle of the ring with radius R_DT. Such a ring will work effectively the same as a point source at the same distance.
So if you have luminosity L_DT=xi_dt * L_disk you would expect energy density at the center of the ring of L_DT / (4 *pi * R_DT * c), and to reproduce this number you would need 4 pi in the formula.
About how to reconcile with your argument: I think the delta factor is not in mu, but in Omega, and then to get from delta factor in Omega to delta in mu and phi you need to multiply by 2 for mu and by 2 pi for phi
Look at formula 9 in here: https://iopscience.iop.org/article/10.1086/341431/fulltext/
and the formula with mu just under in the text, when moving from Omega to mu there is an extra factor 2 appearing.

About 2:
I've spend some time today computing those things, and I'm pretty sure it should be square. Indeed I started with the same formula with u'(epsilon', Omega') scalling with cube, but then if you integrate the epsilon out (with a delta in epsilon) you get one extra power from the transformation from epsilon' to epsilon, and then when you integrate the Omega' out you lose two powers ending up with a power of two.
In the same paper as I cited above there is a similar integration done (but not for a monodirectional radiation, but isotropic) and you also end up with a power of 2, but different numerical factor.
In other words you can look at the original formula with a power of 3 like this: beaming makes: emission region contract along the motion due to Lorentz contraction (1 power), directions in solid angles also contract (you get 2 powers here), nothing changes with energy because you have energy flux, but per unit energy, so 2+1 powers = 3.
And in this particular case you have dE/dV, so beaming gives you boost in energy (1 power) and the same Lorentz contraction (1 power), so you end up with 1+1=2.

I also modified this function computing the energy density, I did it before seeing your response, so I used a somewhat different scheme, passing blob as an optional parameter, if it is not there, the calculations are done in the frame of the galaxy, otherwise in the frame of the blob. It is a simpler way of making it I think, but if you prefer a tag instead I can change it (but then we need to include in the code a check that blob is given when this frame is selected). I committed it to my fork:
jsitarek@4c3999a
and commited together with it a macro to test the two functions to get the gamma break and gamma max from the EC (and this change from 2 pi to 4 pi).
The same macro does also some comparisons between SSC and EC code, and as far as I can see it works pretty well. With the same energy density of synchrotron photons as EC in the blob I get ~40% larger EC emission then SSC, but I think this is because of the angular distribution of the radiation (the difference goes in the correct direction, and I also checked some extreme de-beamed case that the factor makes sense).

If you like the changes let me know and I will make a PR.

from agnpy.

Dear Julian,

About 1.
Sorry but I still disagree, the delta in mu just removes the integration and there is no additional factor 2 provided by the integration in mu.
To convince myself (and yourself) in this document I have reworked from scratch the calculations on u, both in the galaxy and comoving frame, for some of the photon fields
https://github.com/cosimoNigro/radiation_fields_agn/blob/master/tex/main.pdf
I started this document to have all these formulas ready for implementation in agnpy.

The expression for the BLR - cancelling the integration with the delta in mu - is correct!
Indeed, as you can see, the u for the BLR reduces to the case of a point-like source behind the jet (both in the static and comoving frame) for very large distances (Section 3.3.1 and 3.3.2).
You are right that the expression for the torus does not do the same.

What I realised reworking the calculations is that actually the expression for the differential u of the torus

might be wrong, this is why it does not gives us back a simple

luminosity / 4 pi c r^2

when integrated.

It is indeed strange that the emissivity (which is used to compute u) is almost the same for the BLR and the Dust Torus despite they have different geometries,

for the BLR emissivity I understand the denumerator is (4 \pi R^2), as (4 \pi R^2 dR) is the differential (infinitesimal) volume of the shell. But in the ring emissivity the differential volume emitting should be the one of the ring, so something proportional to 2 \pi... right?

from agnpy.

What I mean is that equations 9 and 12 in the documents should be the limit toward which the energy of these radiation fields tends for large distances.
Eliminating the integration in mu with the delta (as I told you I believe it is correct) I can get the energy density of the BLR (in both frames) to reduce to equations 9 and 12.
So I think my implementation in the code of u for the BLR was correct.

For the torus equation it cannot be reduced to equations 9 and 12, using the same integration procedure used for the BLR. Therefore I think that there is a factor 2 missing in the expression of the differential energy density of the torus u(\epsilon, \Omega).
It is not clear to me how the emissivities that Finke uses are calculated. But I think the error must be in Eq. 90 of Finke 2016.

Hope I have made myself more clear.

from agnpy.

I need to go over your document carefully (too many integrals for this hour ;-) ), but I think in many things we actually agree, I also started from Eq 9 and I did the calculations only for the DT, where I was missing this factor of 2, and I think you ended up also with the same factor of 2.
The question is how to get to it from the emissivities, I suspect that there might be some non-standard definition of 'j' there w.r.t. mu or Omega, but let me come back to this after going carefully through the .pdf file.

from agnpy.

Hi,

I went through the calculations, more or less I agree with everything (some comments are embedded in the attached .pdf file)
main2.pdf

So let me try to summarize: there is most probably a mistake in Eq 20 in your pdf file which is Eq 91 in Finke paper, and probably the mistake is already in Eq 90 in that paper. This later on affects the u_ph function in the DT class. The current code (after you merged) has a factor of 2 added by me, which should make the u_ph agree again with our calculation.
In fact in one of the test notebooks: simple_ec.ipynb in Test 1 I was comparing the BLR and DT EC code, and was getting a factor of 2 difference in the energy density of radiation that is also solved now.

What I want to understand is if this possibly also affects the EC calculations itself. I did a test in the past (again in simple_ec file, Test 1) when I run DT and BLR EC with similar parameters. I did also checks comparing SED obtained from BLR EC with DT EC for the same energy density of the radiation field (in the frame of the source). The test that I did it did not look bad, EC DT had 1.3 times the flux of EC BLR. But this test was done with strong beaming, that changes with the angular distriution considerably and thus makes it more difficult to interpret the results. I redone the test also with no beaming (delta=1.02, Gamma=1.01) and I get then that EC DT is 0.65 of EC BLR. If we are missing the factor of 2 in EC emissivity then we would have a factor of 0.5 here, but again there should be some effect also from the angular distribution of the radiation.

In this issue:
#13
When I was comparing the EC radiation with back of the envelope calculations I was also missing about a factor of 2 which is probably the same factor of 2.

Also in the first tests that I was doing, comparing with some old modellings I got a factor of 3 difference which again goes in the same direction.

Julian

from agnpy.

Great, thanks a lot for checking the equations.

Your summary is right and I think all the hints you have found were pointing in the direction of this issue.
And if Eq. 90 is bugged then yes, also the EC computation is.
I obtained the formula for the EC SED starting from Eq. 90.

I'll try to write to Finke, I cannot figure now how to compute the emissivity myself...

Cheers!

from agnpy.

I think this issue is addressed by PR #30, all the radiation energy densities are crosschecked by the case of point source behind the jet.

from agnpy.