[gnuastro-commits] master debf5daa: Library (cosmology.c): Added sanity checks on arguments

[Mohammad Akhlaghi]

branch: master
commit debf5daabdc83d02f8edc83c23a0eeffbec6945c
Author: Jash Shah <jash28582@gmail.com>
Commit: Mohammad Akhlaghi <mohammad@akhlaghi.org>

    Library (cosmology.c): Added sanity checks on arguments
    
    Until now, there was no internal checking perfomed for the input parameters
    (olambda, omatter, oradiation) in the cosmology library functions that
    needed them.
    
    With this commit, a checking function has been added which ensures that the
    fractional densities are between zero and one (inclusive) and their sum
    does not exceed one. It throws an error with appropriate message if either
    of these constraints are violated.
    
    Also, a typo in the Bootstrapping section of the book has been fixed:
    'README' was written as 'REAME'.
---
 doc/gnuastro.texi |  2 +-
 lib/cosmology.c   | 69 +++++++++++++++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 70 insertions(+), 1 deletion(-)

diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index fc19dd74..e1bca91f 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -7808,7 +7808,7 @@ To get a list of the non-version controlled files instead 
of deleting them, add
 
 Besides the @file{bootstrap} and @file{bootstrap.conf}, the 
@file{bootstrapped/} directory and @file{README-hacking} file are also related 
to the bootstrapping process.
 The former hosts all the imported (bootstrapped) directories.
-Thus, in the version controlled source, it only contains a @file{REAME} file, 
but in the distributed tar-ball it also contains sub-directories filled with 
all bootstrapped files.
+Thus, in the version controlled source, it only contains a @file{README} file, 
but in the distributed tar-ball it also contains sub-directories filled with 
all bootstrapped files.
 @file{README-hacking} contains a summary of the bootstrapping process 
discussed in this section.
 It is a necessary reference when you haven't built this book yet.
 It is thus not distributed in the Gnuastro tarball.
diff --git a/lib/cosmology.c b/lib/cosmology.c
index 3c822469..ad42b794 100644
--- a/lib/cosmology.c
+++ b/lib/cosmology.c
@@ -77,6 +77,67 @@ struct cosmology_integrand_t
 
 
 
+/**************************************************************/
+/************     Constraint Check Function      *************/
+/**************************************************************/
+/* Check if input parameters are witihin the required constraints.
+    i.e All density fractions should be between 0 and 1 AND their
+    sum should not exceed 1. */
+static void
+cosmology_density_check(double o_lambda_0, double o_matter_0,
+                        double o_radiation_0)
+{
+  double sum = o_lambda_0 + o_matter_0 + o_radiation_0;
+
+  /* Check if the density fractions are between 0 and 1. */
+  if(o_lambda_0 > 1 || o_lambda_0 < 0)
+    error(EXIT_FAILURE, 0, "value to option 'olambda' must be between "
+          "zero and one (inclusive), but the given value is '%.8f'. Recall "
+          "that 'olambda' is 'Current cosmological cst. dens. per crit. '"
+          "dens.", o_lambda_0);
+
+  if(o_matter_0 > 1 || o_matter_0 < 0)
+    error(EXIT_FAILURE, 0, "value to option 'omatter' must be between "
+          "zero and one (inclusive), but the given value is '%.8f'. Recall "
+          "that 'omatter' is 'Current matter density per critical density.'",
+          o_matter_0);
+
+  if(o_radiation_0 > 1 || o_radiation_0 < 0)
+    error(EXIT_FAILURE, 0, "value to option 'oradiation' must be between "
+          "zero and one (inclusive), but the given value is '%.8f'. Recall "
+          "that 'oradiation' is 'Current radiation density per critical "
+          "density.", o_radiation_0);
+
+  /* Check if the density fractions add up to 1 (within floating point
+      error). */
+  if( sum > (1+1e-8) || sum < (1-1e-8) )
+    error(EXIT_FAILURE, 0, "sum of fractional densities is not 1, "
+          "but %.8f. The cosmological constant ('olambda'), matter "
+          "('omatter') and radiation ('oradiation') densities are given "
+          "as %.8f, %.8f, %.8f", sum, o_lambda_0, o_matter_0,
+          o_radiation_0);
+
+}
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
 /**************************************************************/
 /************         Integrand functions         *************/
 /**************************************************************/
@@ -160,6 +221,7 @@ double
 gal_cosmology_age(double z, double H0, double o_lambda_0, double o_matter_0,
                   double o_radiation_0)
 {
+  cosmology_density_check(o_lambda_0, o_matter_0, o_radiation_0);
   gsl_function F;
   double result, error;
   double o_curv_0 = 1.0 - ( o_lambda_0 + o_matter_0 + o_radiation_0 );
@@ -186,6 +248,7 @@ double
 gal_cosmology_proper_distance(double z, double H0, double o_lambda_0,
                               double o_matter_0, double o_radiation_0)
 {
+  cosmology_density_check(o_lambda_0, o_matter_0, o_radiation_0);
   size_t neval;
   gsl_function F;
   double result, error, c=GSL_CONST_MKSA_SPEED_OF_LIGHT;
@@ -215,6 +278,7 @@ double
 gal_cosmology_comoving_volume(double z, double H0, double o_lambda_0,
                               double o_matter_0, double o_radiation_0)
 {
+  cosmology_density_check(o_lambda_0, o_matter_0, o_radiation_0);
   size_t neval;
   gsl_function F;
   double result, error;
@@ -246,6 +310,7 @@ double
 gal_cosmology_critical_density(double z, double H0, double o_lambda_0,
                                double o_matter_0, double o_radiation_0)
 {
+  cosmology_density_check(o_lambda_0, o_matter_0, o_radiation_0);
   double H;
   double H0s=H0/1000/GSL_CONST_MKSA_PARSEC;     /* H0 in units of seconds. */
   double o_curv_0 = 1.0 - ( o_lambda_0 + o_matter_0 + o_radiation_0 );
@@ -266,6 +331,7 @@ double
 gal_cosmology_angular_distance(double z, double H0, double o_lambda_0,
                                double o_matter_0, double o_radiation_0)
 {
+  cosmology_density_check(o_lambda_0, o_matter_0, o_radiation_0);
   return gal_cosmology_proper_distance(z, H0, o_lambda_0, o_matter_0,
                                        o_radiation_0) / (1+z);
 }
@@ -279,6 +345,7 @@ double
 gal_cosmology_luminosity_distance(double z, double H0, double o_lambda_0,
                                   double o_matter_0, double o_radiation_0)
 {
+  cosmology_density_check(o_lambda_0, o_matter_0, o_radiation_0);
   return gal_cosmology_proper_distance(z, H0, o_lambda_0, o_matter_0,
                                        o_radiation_0) * (1+z);
 }
@@ -292,6 +359,7 @@ double
 gal_cosmology_distance_modulus(double z, double H0, double o_lambda_0,
                                double o_matter_0, double o_radiation_0)
 {
+  cosmology_density_check(o_lambda_0, o_matter_0, o_radiation_0);
   double ld=gal_cosmology_luminosity_distance(z, H0, o_lambda_0, o_matter_0,
                                               o_radiation_0);
   return 5*(log10(ld*1000000)-1);
@@ -306,6 +374,7 @@ double
 gal_cosmology_to_absolute_mag(double z, double H0, double o_lambda_0,
                               double o_matter_0, double o_radiation_0)
 {
+  cosmology_density_check(o_lambda_0, o_matter_0, o_radiation_0);
   double dm=gal_cosmology_distance_modulus(z, H0, o_lambda_0, o_matter_0,
                                            o_radiation_0);
   return dm-2.5*log10(1.0+z);