Supplying your own screening function

In my last post, I showed how to compute electron screening and include it in a network calculation. In this post, I will how a user could change the screening to his or her own parameterization.

The parameterization I will use here is that of DeWitt, Graboske, and Cooper. I first write files my_screen.cpp and my_screen.h, based on screen.cpp and screen.h in the nucnet-tools-code/user directory. The file my_screen.cpp, like the parent screen.cpp, has five key functions:

  1. set_screening_function(). This routine sets the screening function for a zone and the zone function that retrieves the screening data for the zone. The only change a user would make here is to the name of the zone screening function.
  2. get_screening_data(). This routine returns a user-defined structure giving the data for the screening function. The user would probably not change this function except perhaps to pass in other data. The user-defined structure is defined in my_screen.h, and it is here the user would make changes.
  3. screening_function(). This routine computes the screening factor for a reaction and applies it to the (already computed) forward and reverse reaction rates. In screen.cpp, this function applies the larger of the screening factors for the forward and reverse reactions to both the forward and reverse rates. It also computes and applies the reverse ratio correction factor (see libnucnet documentation), if present. In my_screen.cpp, the routine applies the screening for the forward rate to both rates and the reverse ratio correcton factor only to the reverse rate. The user would change how screening factors are applied to forward and reverse rates here.
  4. reaction_screening_function(). This routine computes the screening factor for the forward reaction direction and the reverse reaction direction. The user could change how screening is applied for more than two reacting species here. The default for a reaction a + b + c is to compute the screening for a + b and then the screening for (a + b) + c; that is, the reaction is imagined to proceed through an intermediate step a + b.
  5. pair_screening_function(). This routine computes the screening factor between two reactants.

To compute the new screening, I place my_screen.cpp and my_screen.h in the my_users directory, as well as the file I next put the files compute_screening_factors_in_zones.cpp, compare_screened_rates_in_zones.cpp, and Makefile in my_examples/screen. Finally, I place the new files run_single_zone.cpp and Makefile in my_examples/network. These files are appropriately modifed versions of their counterparts in the examples and user directories. Note that the codes in my_examples had to be modified to use the routines in the my_user namespace.

For convenience, I have created a tar ball 2016-06-08.tar. Rather than separately download the individual files, I can place this tar file in nucnet-tools-code and type

tar xvf 2016-06-08.tar

However I managed to download the new files, once they are available, I can compile and run them. To do so, I type

cd my_examples/network

make run_single_zone

I then run the network calculation I previously ran by typing

./run_single_zone ../../data_pub/my_net.xml ../../data_pub/zone.xml my_output_screening.xml "[z <= 10]"

Were I to plot the 4He mass fraction versus time for this screened rates calculation versus the previous unscreened calculation, I would get a figure very similar to that in the previous post.

I can now compute the screening factors and screened rates for this screening treatment. I type

cd ../screen

make all_screen

I then compute the PP reaction screening factors by typing

./compute_screening_factors_in_zones ../network/my_output_screening.xml --zone_xpath "[last()]" --reac_xpath "[reactant[position() = 1] = 'h1' and reactant[position() = 2] = 'h1']"

which yields

time(s) = 3.1500e+17  t9 = 0.0150  rho(g/cc) = 1.5000e+02  Ye = 0.5816

        Reaction              Forward Screening  Reverse Screening

=======================================================   ================ ===========

h1 + h1 + electron -> h2 + neutrino_e                     1.05111e+00   1.00000e+00
h1 + h1 -> h2 + positron + neutrino_e                     1.05111e+00   1.00000e+00

This shows this treatment of screening gives slightly smaller screening factors for these reactions that the default treatment (1.05111 vs. 1.05765). I can then compute the screened rates as before by typing, for example,

./compare_screened_rates_in_zones ../network/my_output.xml --zone_xpath "[last()]" --reac_xpath "[reactant[position() = 1] = 'h1' and reactant[position() = 2] = 'h1']"

The screened rates values reflect this different screening treatment.

With these instructions, you should be able to employ your own screening directly or within a network calculation. The simplest approach would be to modify my_screen.cpp and my_screen.h for your purposes. Of course, you should use subversion, git, or similar to keep track of your own versions of these files--recall that the my_user and my_examples directories are not version controlled and are ignored in NucNet Tools.

Posted by Bradley S. Meyer 2016-06-08


Cancel  Add attachments

Get latest updates about Open Source Projects, Conferences and News.

Sign up for the SourceForge newsletter:

No, thanks