Webnucleo.org

Mail Lists | Developers
small logo


Warning: This is no longer the latest available version of this module. Please see the releases page for the most recent version. The Webnucleo group strongly recommends the use of the latest version of any of its online modules.

Documentation from the file



libnucnet/0.2/src/Libnucnet__Nuc.h


Table of Contents

Structures
Name: Libnucnet__Species

Description: Libnucnet__Species is a structure that stores nuclear data for a particular nuclear species. Routines act on the structure to retrieve data, update data, or use data to compute functions of the nuclear data, such as the nuclear partition function. The contents of the structure are not made public by the API.

-top-


Name: Libnucnet__Nuc

Description: Libnucnet__Nuc is a structure that stores the collection of nuclear species. Routines act on the structure to retrieve species or add or remove them. The contents of the structure are not made public by the API.

-top-



User-Supplied Routines
Name: Libnucnet__Species__iterateFunction()

Description: User-supplied routine to be applied during an iteration over the species in the species collection.

Syntax:
         int
         Libnucnet__Species__iterateFunction(
           Libnucnet__Species *self,
           void *p_user_data
         );
             
Input:

self: (required) A pointer to a Libnucnet__Species.

p_user_data: (required) A pointer to a user-defined structure containing extra data for the iterate function.

Output:

User's routine must return 1 to continue or 0 to stop. User's routine must not modify self, that is, the pointer to the input Libnucnet__Species structure.


-top-


Name: Libnucnet__Species__nseCorrectionFactorFunction()

Description: Optional user-supplied routine to compute corrections to the NSE factor for a species.

Syntax:
         double
         Libnucnet__Species__nseCorrectionFactorFunction(
           Libnucnet__Species *self,
           double d_t9,
           double d_rho,
           double d_ye,
           void *p_user_data
         );
             
Input:

self: (required) A pointer to the Libnucnet__Species whose NSE correction factor is required.

d_t9: (required) A double giving the temperature in billions of K at which to compute the NSE correction factor.

d_rho: (required) A double giving the density in g/cc at which to compute the NSE correction factor.

d_ye: (required) A double giving the electron-to-baryon ratio at which to compute the NSE correction factor.

p_user_data: (required) A pointer to a user-defined structure containing extra data for the NSE correction factor function.

Output:

User's routine must return a double giving the NSE correction ratio for the input temperature, density, and electron-to-baryon ratio.


-top-



Routines
Name: Libnucnet__Nuc__addSpecies()

Description: Add a species to a species collection.

Syntax:
         void
         Libnucnet__Nuc__addSpecies(
           Libnucnet__Nuc *self,
           Libnucnet__Species *p_species
         );
             
Input:

self: (required) A pointer to the species collection.

p_species: (required) A pointer to the species to be added to the collection.

Output:

For valid input, the routine updates the Libnucnet__Nuc structure to include the new species. If the species already exists, its data are replaced with the new data. If the species cannot be added, Libnucnet__Nuc error handling is invoked.

Example: Add p_species to the species collection p_my_nuclei:

         Libnucnet__Nuc__addSpecies( p_my_nuclei, p_species );
               

-top-


Name: Libnucnet__Nuc__free()

Description: Free the memory allocated for a Libnucnet nuclear species collection.

Syntax:
         void Libnucnet__Nuc__free( Libnucnet__Nuc *self );
             
Input:

self: (required) A pointer to a Libnucnet nuclear species collection.

Output:

Upon successful return, the memory for the collection has been freed.

Example: Free up memory used for Libnucnet__Nuc *p_my_nuclei:

         Libnucnet__Nuc__free( p_my_nuclei );
               

-top-


Name: Libnucnet__Nuc__getNumberOfSpecies()

Description: Retrieves the number of species in the collection of nuclear species.

Syntax:
         size_t
         Libnucnet__Nuc__getNumberOfSpecies( Libnucnet__Nuc *self );
             
Input:

self: (required) A pointer to a Libnucnet__Nuc structure.

Output:

Routine returns an integer representing the number of species in the collection of nuclear species. If the input is invalid, error handling is invoked.

Example: Get number of species in nuclear data collection p_my_nuclei:

         size_t i_species_count;
         i_species_count =
           Libnucnet__Nuc__getNumberOfSpecies( p_my_nuclei );
               

-top-


Name: Libnucnet__Nuc__getSpeciesByName()

Description: Retrieves the specified species.

Syntax:
         Libnucnet__Species *
         Libnucnet__Nuc__getSpeciesByName( 
           Libnucnet__Nuc *self, const char *s_nucname
         );
             
Input:

self: (required) A pointer to a Libnucnet__Nuc structure.

s_nucname: (required) A string giving the name of the particular species to be retrieved.

Output:

Routine returns a pointer to a Libnucnet__Species structure. If the species is not found, routine returns NULL. If the input structure is invalid, error handling is invoked.

Example: Retrieve he4 from nuclear data collection p_my_nuclei:

         p_he4 =
           Libnucnet__Nuc__getSpeciesByName( p_my_nuclei, "he4" );
               

-top-


Name: Libnucnet__Nuc__getSpeciesByZA()

Description: Retrieves a species in the collection of species given the Z, A, and state of the species.

Syntax:
         Libnucnet__Species *
         Libnucnet__Nuc__getSpeciesByZA( 
           Libnucnet__Nuc *self,
           unsigned int i_z,
           unsigned int i_a,
           const char *s_state
         );
             
Input:

self: (required) A pointer to a Libnucnet__Nuc structure.

i_z: (required) The Z of the species desired. i_z is a valid Z of a species in the collection of species.

i_a: (required) The A of the species desired. i_a is a valid A of a species in the collection of species.

s_state: (required) The state of the species desired. s_state is a valid state of a species in the collection of species. If not present, use empty string or NULL.

Output:

Routine returns a pointer to the species retrieved. If no such species is found, the routine returns a NULL. If the input structure is not valid, error handling is invoked.

Examples: Get the pointer to the species al27 in nuclear data structure p_my_nuclei:

         Libnucnet__Species *p_al27;
         p_al27 =
           Libnucnet__Nuc__getSpeciesByZA( p_my_nuclei, 13, 27, NULL ) );
               
Get the pointer to species al26m in a collection of species p_my_nuclei:

         Libnucnet__Species *p_al26m;
         p_al26_m =
           Libnucnet__Nuc__getSpeciesByZA( p_my_nuclei, 13, 26, 'm' ) );
               

-top-


Name: Libnucnet__Nuc__is_valid_input_xml()

Description: Validate an input xml file for Libnucnet__Nuc.

Syntax:
         int
         Libnucnet__Nuc__is_valid_input_xml(
           const char *s_xml_filename
         );
             
Input:

s_xml_filename: (required) A string giving the name of the xml file containing the nuclear data.

Output:

For a valid input nuclear data xml file, the routine returns 1 (true). If the input file is not valid, routine returns 0 (false). For an invalid schema file, or if the schema file cannot be read over the web, routine stops and prints error message.

Example: Validate the input xml file "nuclear_data.xml":

         if( Libnucnet__Nuc__is_valid_input_xml( "nuclear_data.xml" ) ) {
             printf( "Valid xml!\n" );
         }
               

-top-


Name: Libnucnet__Nuc__iterateSpecies()

Description: Iterate through the species and apply the user-supplied iterate function.

Syntax:
         void
         Libnucnet__Nuc__iterateSpecies(
           Libnucnet__Nuc *self,
           Libnucnet__Species__iterateFunction pfFunc,
           void *p_user_data
         );
             
Input:

self: (required) A pointer to a Libnucnet__Nuc structure.

pfFunc: (required) The name of the user-supplied function to apply.

p_user_data: (required) A pointer to the user-supplied data structure.

Output:

The routine iterates through the species collection and applies the user-supplied routine to each species. If any input is invalid, error handling in invoked.

Example: Iterate through the species in p_my_nuclei and apply the function my_iterate_function and the extra data in p_user_data:

         Libnucnet__Nuc__iterateSpecies(
           p_my_nuclei,
           (Libnucnet__Species__iterateFunction) my_iterate_function,
           p_user_data
         );
               

-top-


Name: Libnucnet__Nuc__new()

Description: Create a new Libnucnet__Nuc structure.

Syntax:
         Libnucnet__Nuc *Libnucnet__Nuc__new( );
             
Output:

The routine returns a pointer to a Libnucnet__Nuc structure, that is, a collection of nuclear species. If it is not possible to allocate memory for the new structure, Libnucnet error handling is invoked.

Example: Create a Libnucnet__Nuc structure p_my_nuclei:

         p_my_nuclei = Libnucnet__Nuc__new( );
               

-top-


Name: Libnucnet__Nuc__new_from_xml()

Description: Creates a Libnucnet__Nuc structure from an input xml nuclear data file.

Syntax:
         Libnucnet__Nuc *
         Libnucnet__Nuc__new_from_xml(
           const char *s_xml_filename,
           const char *s_xpath
         );
             
Input:

s_xml_filename: (required) A string giving the name of the xml file containing the nuclear data. This may be the name of a local file or a URL.

s_xpath: (required) A string giving an xpath expression to apply.

Output:

For a valid input nuclear data xml file, the routine returns a pointer to a Libnucnet__Nuc structure. If the routine cannot allocate enough memory or if the xpath expression is invalid, Libnucnet__Nuc error handling is invoked.

Examples: Store the nuclear data in "nuclear_data.xml" to p_my_nuclei:

         p_my_nuclei =
           Libnucnet__Nuc__new_from_xml( "nuclear_data.xml", NULL );
               
Store the nuclear data for neutrons and protons and nuclei with z >= 6 in "nuclear_data.xml" to p_my_nuclei:

         p_my_nuclei =
           Libnucnet__Nuc__new_from_xml(
             "nuclear_data.xml", "[ a = 1 or z >= 6]"
           );
               

-top-


Name: Libnucnet__Nuc__removeSpecies()

Description: Remove a species from a Libnucnet__Nuc structure.

Syntax:
         void
         Libnucnet__Nuc__removeSpecies(
           Libnucnet__Nuc *self, Libnucnet__Species *p_species
         );
             
Input:

self: (required) A pointer to a Libnucnet nuclear species collection.

p_species: (required) A pointer to a Libnucnet species structure to be removed from self.

Output:

Upon successful return, the species has been removed from the Libnucnet__Nuc structure and the species indices have been updated. If the input is invalid, Libnucnet__Nuc error handling is invoked.

Example: Remove he4 from Libnucnet__Nuc *p_my_nuclei:

         Libnucnet__Species *p_he4 =
           Libnucnet__Nuc__getSpeciesByName( p_my_nuclei, "he4" );
         Libnucnet__Nuc__removeSpecies( p_my_nuclei, p_he4 );
               

-top-


Name: Libnucnet__Nuc__sortSpecies()

Description: Sort the species in a Libnucnet__Nuc structure in increasing Z and A.

Syntax:
         void Libnucnet__Nuc__sortSpecies( Libnucnet__Nuc *self );
             
Input:

self: (required) A pointer to a Libnucnet__Nuc structure.

Output:

On successful return, the species in the structure have been sorted so that successive species are higher in Z and/or A. If the input structure is not valid or memory cannot be allocated for the sorting, error handling is invoked.

Example: Sort the species in p_my_nuclei:

         Libnucnet__Nuc__sortSpecies( p_my_nuclei );
               

-top-


Name: Libnucnet__Nuc__updateFromXml()

Description: Updates the data in a Libnucnet__Nuc structure from an input xml nuclear data file.

Syntax:
         void
         Libnucnet__Nuc__updateFromXml(
           Libnucnet__Nuc *self,
           const char *s_xml_filename,
           const char *s_xpath
         );
             
Input:

self: (required) A pointer to a pre-existing Libnucnet__Nuc structure.

s_xml_filename: (required) A string giving the name of the xml file containing the nuclear data. This may be the name of a local file or a URL.

s_xpath: (required) A string giving an xpath expression to apply to the data in s_xml_filename.

Output:

For a valid Libnucnet__Nuc structure and a valid input nuclear data xml file, the routine updates the structure with the data in the input data file. If the routine cannot allocate enough memory or if the xpath expression is invalid, Libnucnet__Nuc error handling is invoked.

Example: Update the data in p_my_nuclei with data in "new_nuclear_data.xml":

         Libnucnet__Nuc__updateFromXml(
           p_my_nuclei, "new_nuclear_data.xml", NULL
         );
               

-top-


Name: Libnucnet__Nuc__writeToXmlFile()

Description: Output a Libnucnet__Nuc structure to an xml file.

Syntax:
         void
         Libnucnet__Nuc__writeToXmlFile(
           Libnucnet__Nuc *self,
           const char *s_xml_filename
         );
             
Input:

self: (required) A pointer to a Libnucnet nuclear species collection.

s_xml_filename: (required) The name of the output xml file.

Output:

Upon successful return, the contents of the collection of species have been written to an xml file. If the input is invalid, Libnucnet__Nuc error handling is invoked.

Example: Dump the contents of Libnucnet__Nuc *p_my_nuclei to the xml file my.xml:

         Libnucnet__Nuc__writeToXmlFile( p_my_nuclei, "my.xml" );
               

-top-


Name: Libnucnet__Species__computePartitionFunction()

Description: Compute the partition function of a species in the collection of nuclear species given the temperature.

Syntax:
         double
         Libnucnet__Species__computePartitionFunction(
           Libnucnet__Species *self, double d_t9
         );
             
Input:

self: (required) A pointer to a Libnucnet species structure.

d_t9: (required) The temperature (in billions of K) at which the partition function will be calculated.

Output:

Routine returns a double representing the partition function of the species retrieved. If input is invalid, error handling is invoked.

Example: Compute the partition function of he4 in the collection of nuclear species p_my_nuclei at a t9 of 0.2:

         if(
             (
               p_species =
                 Libnucnet__Nuc__getSpeciesByName( p_my_nuclei, "he4" )
             )
         ) {
             d_partf =
               Libnucnet__Species__computePartitionFunction(
                  p_species, 0.2
               );
         }
               

-top-


Name: Libnucnet__Species__getA()

Description: Retrieves the A of a species in the collection of nuclear species.

Syntax:
         unsigned int
         Libnucnet__Species__getA( Libnucnet__Species *self );
             
Input:

self: (required) A pointer to a Libnucnet__Species structure.

Output:

Routine returns the A of the input species. If the species is invalid, Libnucnet__Nuc error handling is invoked.

Example: Print the A of au197 in p_my_nuclei (should of course be 197!):

         p_species =
           Libnucnet__Nuc__getSpeciesByName( p_my_nuclei, "au197" );
         if( p_species )
           printf(
             "A = %d\n",
             Libnucnet__Species__getA( p_species )
           );
               

-top-


Name: Libnucnet__Species__getIndex()

Description: Routine to return the index of a species in a Libnucnet__Nuc structure.

Syntax:
         size_t
         Libnucnet__Species__getIndex( Libnucnet__Species *self );
             
Input:

self: (required) A pointer to a Libnucnet species structure.

Output:

Returns the index of the species in the nuclear species collection. If the species is invalid, Libnucnet__Nuc error handling is invoked.

Example: Print the index of he4 in the Libnucnet structure p_my_nuclei:

         if( 
             p_he4 = 
               Libnucnet__Nuc__getSpeciesByName(
                 p_my_nuclei, "he4"
               )
         ) {
           printf(
             "Index of he4 = %d\n",
             Libnucnet__Species__getIndex( p_he4 )
           );
         }
               

-top-


Name: Libnucnet__Species__getMassExcess()

Description: Retrieves the mass excess of a species in the nuclear species collection.

Syntax:
         double
         Libnucnet__Species__getMassExcess( Libnucnet__Species *self );
             
Input:

self: (required) A pointer to a Libnucnet__Species structure.

Output:

Routine returns a double giving the mass excess of the requested species. If the input species is invalid, error handling is invoked.

Example: Print the mass excess of sn120 in the structure Libnucnet__Nuc *p_my_nuclei:

         p_species =
           Libnucnet__Nuc__getSpeciesByName( p_my_nuclei, "sn120" );
         if( p_species ) {
           printf(
             "Mass excess (MeV) = %lf\n,
             Libnucnet__Species__getMassExcess( p_species )
           );
         }
               

-top-


Name: Libnucnet__Species__getName()

Description: Retrieves the name of a species in the collection of nuclear species.

Syntax:
         const char *
         Libnucnet__Species__getName( Libnucnet__Species *self );
             
Input:

self: (required) A pointer to a Libnucnet species structure.

Output:

Routine returns a string representing the name of the requested species. If the species is not valid, Libnucnet__Nuc error handling is invoked.

Example: Print the name of the Libnucnet__Species *p_species:

         printf(
           "Name of species = %s\n",
           Libnucnet__Species__getName( p_species )
         );
               

-top-


Name: Libnucnet__Species__getSource()

Description: Retrieves the message about the source of data for a species in the collection of nuclear species.

Syntax:
         const char *
         Libnucnet__Species__getSource( Libnucnet__Species *self );
             
Input:

self: (required) A pointer to a Libnucnet species structure.

Output:

Routine returns a string containing information about the data for the requested species. If there are no data, an empty string is returned. If the species is not valid, Libnucnet__Nuc error handling is invoked.

Example: Print the source information about the Libnucnet__Species *p_species:

         printf(
           "Data source for species = %s\n",
           Libnucnet__Species__getSource( p_species )
         );
               

-top-


Name: Libnucnet__Species__getSpin()

Description: Retrieves the spin of a species in the collection of nuclear species.

Syntax:
         double
         Libnucnet__Species__getSpin( Libnucnet__Species *self );
             
Input:

self: (required) A pointer to a Libnucnet__Species structure.

Output:

Routine returns a double giving the spin of the requested species. If the species is invalid, Libnucnet__Nuc error handling is invoked.

Example: Print the spin of bi209 in the Libnucnet__Nuc structure *p_my_nuclei:

         if(
             ( p_species =
                 Libnucnet__Nuc__getSpeciesByName( p_my_nuclei, "bi209 )
           )
         {
            printf(
              "Spin = "%lf\n",
              Libnucnet__Species__getSpin( p_species )
            );
         }
               

-top-


Name: Libnucnet__Species__getZ()

Description: Retrieves the Z of a species in the collection of nuclear species.

Syntax:
         unsigned int
         Libnucnet__Species__getZ( Libnucnet__Species *self );
             
Input:

self: (required) A pointer to a Libnucnet__Species structure.

Output:

Routine returns the Z of the input species. If the species is invalid, Libnucnet__Nuc error handling is invoked.

Example: Print the Z of au197 in p_my_nuclei:

         p_species =
           Libnucnet__Nuc__getSpeciesByName( p_my_nuclei, "au197" );
         if( p_species )
           printf(
             "Z = %d\n",
             Libnucnet__Species__getZ( p_species )
           );
               

-top-


Name: Libnucnet__Species__new()

Description: Create a new species from the input data.

Syntax:
         Libnucnet__Species *
         Libnucnet__Species__new(
             unsigned int i_z,
             unsigned int i_a,
             const char *s_source,
             int i_state_flag,
             const char *s_state,
             double d_mass_excess,
             double d_spin,
             gsl_vector *p_t9,
             gsl_vector *p_log10_partf
           );
             
Input:

i_z: (required) An unsigned int giving the atomic number of the species.

i_a: (required) An unsigned int giving the atomic mass of the species.

s_source: (required) A string giving a message containing information about the source of the data or NULL.

i_state_flag: (required) An int indicating whether the nuclear species has only one state included in the collection (0) or a ground state and at least one excited state included in the collection (1).

s_state: (required) A string giving the state id of the species. This string should be empty (or NULL) for no state suffix.

d_spin: (required) A double giving the spin of the species.

d_mass_excess: (required) A double giving the mass excess (in MeV) of the species.

p_t9: (required) A gsl_vector giving the temperatures in 10^9 K at which the partition function of the species is evaluated.

p_log10_partf: (required) A gsl_vector giving the log10 of the partition function factor of the species. The partition function is 10 raised to this factor times (2*spin + 1) of the state.

Output:

For valid input, the routine returns a pointer to a new Libnucnet__Species structure. If the species cannot be created, error handling is invoked.

Examples: Create li7 (Z = 3, A = 7, spin = 3/2, mass excess = 14.908 MeV), and partition function data stored in gsl_vectors p_t9 and p_log10_partf:

         p_li7 =
           Libnucnet__Species__new(
             3,
             7,
             "example",
             0,
             "",
             14.908,
             1.5,
             p_t9,
             p_log10_partf
           );
               
Create ground state (g) of al26 (Z = 13, A = 26, spin = 5, mass excess = -12.210 MeV), and partition function data stored in p_t9_g and p_log10_partf_g and metastable state (m) of al26 (Z = 13, A = 26, spin = 0, mass excess = -11.982), and partition function data stored in p_t9_m and p_log10_partf_m:

         p_al26g =
           Libnucnet__Species__new(
             13,
             26,
             "example data",
             1,
             "g",
             -12.210,
             5.,
             p_t9_g,
             p_log10_partf_g
           );
        
         p_al26m =
           Libnucnet__Species__new(
             13,
             26,
             "example data",
             1,
             "m",
             -11.982,
             0.,
             p_t9_m,
             p_log10_partf_m
           );
               

-top-




Valid XHTML 1.1        Copyright © 2001-2012, Clemson University. All rights reserved.        Valid CSS!
Page last modified on 2008/08/01 14:45