[mesa-users] Total Energy Integral
Max Katz
maxpkatz at gmail.com
Wed Mar 11 00:23:28 EDT 2015
Slightly related: for a self-gravitating system, the total gravitational
energy is (1/2) int(phi dm). Noted bugs aside, the missing factor of 1/2
might explain why "total_energy_integral" is positive everywhere in the
figure Michael sent (whereas we expect negative energy for a bound system).
Max Katz
Department of Physics and Astronomy
Vice President, Graduate Student Organization
Stony Brook University
http://astro.sunysb.edu/mkatz/
On Tue, Mar 10, 2015 at 10:34 PM, Bill Paxton <paxton at kitp.ucsb.edu> wrote:
>
> On Mar 10, 2015, at 4:43 PM, Michael Medford wrote:
>
> Bill,
>
> Thanks as always for your help, and for guiding me to understand the
> answer instead of simply having it. I have worked my way through the code
> and these definitions.
>
> For total_integral_energy, both dm(k) and energy(k) are averages for the
> cell and therefore should be thought of as occurring at the cell's center.
> The multiplication of these two factors should yield the internal
> energy occurring "within" this zone. Add them all up to get the total
> internal energy that occurs between this zone and the surface (where k=1).
> Total_internal_energy should be near-zero at the surface and a larger
> positive value as the integral penetrates into the star.
>
> For total_gravitational_energy, m_grav(k) and r(k) are defined at the
> cell's outer edge. I have checked the profile to ensure that there is no
> difference between baryonic and gravitational mass in my models. Indeed
> these two values in my circumstances are identical. I see why dm_bar(k) is
> the more appropriate term to multiply into m_grav(k) and r(k), since
> dm_bar(k) is an average of adjacent zone-centered masses from dm(k) and
> ends up lying on each zone's outer edge. Multiplying these factors together
> gives you the gravitational energy occurring "at" the outer boundary of the
> zone. Add them all up to get the total gravitational energy that occurs
> between this zone's outer edge and the surface (where k=1).
> Total_gravitational_energy should be near-zero at the surface and a larger
> negative value as the integral penetrates into the star.
>
>
> great.
>
>
> It seems to me that my volume_term (4*pi*r^2*rho*dr) will always be equal
> to dm(k) instead of dm_bar(k), since it only uses values that are defined
> as at the center of each cell. (I checked and indeed the volume_term very
> nearly equals dm(k). There is a slight discrepancy at the center and
> surface due to different boundary conditions, but that should not
> dramatically shift the overall shape of the integral). To make sure each
> value occurred at the same place in using my volume_term, I shifted each of
> my radii and mass values to occur at the center of each cell (I created my
> own rmid and averaged the masses between zones).
>
>
> there's a small point concerning the right way to find r at the center of
> the cell. if we consider the center to be by mass, and treat rho as
> average density, the what would you use for center radius?
>
> Essentially instead of shifting a dm(k) to a dm_bar(k) (averaging center
> values to get outer edge values), I shifted r(k) and m(k) to the centers
> (averaging outer edge values to get center values). By placing both the
> energy_integral and gravitational_integral at the center of each zone, I
> can then add them together in a simple fashion to get my own
> total_energy_integral.
>
> It looks to me that the code calculates total_energy_integral without
> taking into account the differences between an
> energy_integral(k) occurring at the center and a
> gravitational_integral(k) occurring at the outer boundary. From
> profile_getval.f:
>
> case (p_total_energy_integral)
> val = 0
> do kk=1,k
> val = val + s% dm(kk)*( &
> exp_cr(s% lnE(kk)) - &
> s% cgrav(kk)*s% m_grav(kk)/s% r(kk) + &
> 0.5d0*s% velocity(kk)*s% velocity(kk))
> end do
>
> What am I missing here?
>
>
> you are only missing some changes I've made since your version (and, yes,
> that's a bug in the code you show)!. here's what I have now:
>
> case (p_total_energy_integral)
> val = get_total_energy_integral(s,k)
>
> the routine get_total_energy_integral is defined in star_utils. it does
> a sum from the surface to k of the following:
>
> total_internal_energy = dot_product(s% dm(1:k), s% energy(1:k))
> total_gravitational_energy = -dot_product(s% dm_bar(1:k), &
> s% cgrav(1:k)*s% m_grav(1:k)/s% r(1:k))
> if (s% v_flag) then
> total_linear_kinetic_energy = &
> 0.5d0*dot_product(s% dm_bar(1:k), s% v(1:k)*s% v(1:k))
> else
> total_linear_kinetic_energy = &
> 0.5d0*dot_product(s% dm_bar(1:k), &
> s% dlnR_dt(1:k)*s% dlnR_dt(1:k)*s% r(1:k)*s% r(1:k))
> end if
> if (s% rotation_flag) then
> total_rotational_kinetic_energy = 0.5d0*dot_product( &
> s% i_rot(1:k), s% omega(1:k)*s% omega(1:k))
> else
> total_rotational_kinetic_energy = 0d0
> end if
>
>
> Hope that helps.
>
> Cheers,
> Bill
>
>
>
>
>
>
>
>
> It looks like dm(k) is being multiplied into both the internal energy and
> the gravitational term (and the kinetic term) without regard to where the
> value occurs. I am also unable to find the eval_section_total_energy_integrals
> that Richard describes as occurring in star/private/star_utils.f90. I am
> running Mesa Version 7393.
>
> Thanks as always for the help!
> ~Michael
>
> On Mon, Mar 9, 2015 at 6:25 PM, Bill Paxton <paxton at kitp.ucsb.edu> wrote:
>
>> Hi,
>>
>> I'd suggest a careful reading of the following two lines of code -- think
>> about what the various terms mean. check where things are defined (cell
>> center vs cell boundary). what's the difference between dm and dm_bar, and
>> why do we use dm for total energy and dm_bar for gravitational energy?
>> what's the definition of s% energy(k)? what disagreements could arise
>> between dm or dm_bar vs 4*pi*r^2*rho*dr given that we are considering
>> finite volumes rather than differential equations? hopefully enlightenment
>> will follow. ;D
>>
>> we sum the following for k = 1, s% nz
>> total_internal_energy = dot_product(s% dm(1:k), s% energy(1:k))
>> total_gravitational_energy = -dot_product(s% dm_bar(1:k), &
>> s% cgrav(1:k)*s% m_grav(1:k)/s% r(1:k))
>>
>> Cheers,
>> Bill
>>
>>
>>
>> On Mar 9, 2015, at 1:22 PM, Michael Medford wrote:
>>
>> Hello everyone,
>>
>> Over the last few weeks I have been investigating the profile column
>> variable labelled total_energy_integral. The description of this variable
>> is "sum from surface inwards of dm*(v^2/2 - G m / r + E) (ergs)". I was
>> confused as to whether this sum includes the necessary volume elements to
>> properly sum the spherical shells within this integral. So I set about
>> calculating my own energy integral using the necessary data from each cell
>> within the profile. Here is how I calculated my integral:
>>
>> -Routines to shift each of the mass, radii and density terms to the
>> center of each cell by taking averages of the outer and inner boundaries
>> -Define "dr" as the radii width of each shell
>> -Defining an "energy element" = (-Gm/r + E) * (4*pi*r^2*density*dr) where
>> all of the terms are in cgs and the velocity of all of the cells are 0
>> -Sum from the surface inwards of this energy element
>>
>> Attached you will see the results of my variable versus the output from
>> the profile's total_energy_integral. As the energy integral claims to have
>> a dm multiplied into each term, I did not multiply the
>> total_energy_integral times anything, since this should be equivalent to
>> (4*pi*r^2*density*dr).
>>
>> Does anyone know why these two values would be different, and by such
>> large orders of magnitude? What am I missing about the
>> total_energy_integral and what it is recording? Or am I manually summing
>> the integral incorrectly?
>>
>> Thanks so much for the help.
>> ~Michael
>>
>> ----
>> Michael Medford
>> Northwestern University
>> <Energy_Integrals_2.png><Energy_Integrals_1.png>
>> ------------------------------------------------------------------------------
>> Dive into the World of Parallel Programming The Go Parallel Website,
>> sponsored
>> by Intel and developed in partnership with Slashdot Media, is your hub
>> for all
>> things parallel software development, from weekly thought leadership
>> blogs to
>> news, videos, case studies, tutorials and more. Take a look and join the
>> conversation now.
>> http://goparallel.sourceforge.net/_______________________________________________
>> mesa-users mailing list
>> mesa-users at lists.sourceforge.net
>> https://lists.sourceforge.net/lists/listinfo/mesa-users
>>
>>
>>
>
>
> --
> Michael Medford
> Northwestern University
>
>
>
>
> ------------------------------------------------------------------------------
> Dive into the World of Parallel Programming The Go Parallel Website,
> sponsored
> by Intel and developed in partnership with Slashdot Media, is your hub for
> all
> things parallel software development, from weekly thought leadership blogs
> to
> news, videos, case studies, tutorials and more. Take a look and join the
> conversation now. http://goparallel.sourceforge.net/
> _______________________________________________
> mesa-users mailing list
> mesa-users at lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/mesa-users
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20150311/af629ab4/attachment.html>
More information about the Mesa-users
mailing list