doubleutil.c File Reference

Working with doubles. More...

#include "tpcclibConfig.h"
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <string.h>
#include <ctype.h>
#include <unistd.h>
#include <math.h>
#include "tpcextensions.h"

Go to the source code of this file.

Functions
int	doubleMatch (const double v1, const double v2, const double lim)
int	doubleArrayMatch (const double *a1, const double *a2, const unsigned int n, const double lim)
int	doubleMatchRel (const double v1, const double v2, const double lim)
double	doubleMachEps ()
void	doubleCopy (double *t, double *s, const unsigned int n)
unsigned int	doubleCopyFinite (double *t, double *s, const unsigned int n)
unsigned int	doubleNaNs (double *a, const unsigned int n)
unsigned int	doubleRange (double *a, const unsigned int n, double *amin, double *amax)
double	doubleSum (double *a, const unsigned int n)
double	doubleMean (double *a, const unsigned int n)
double	doubleWMean (double *a, double *w, const unsigned int n)
int	doubleGetWithUnit (const char *s, double *v, int *u)
int	doubleSpanPositives (double *a, const int n)
int	doubleCSpanPositives (double *a, const int n)
unsigned int	doubleNonzeroes (double *a, const unsigned int n)
unsigned int	doubleMaxIndex (double *a, const unsigned int n)
unsigned int	doubleAbsMaxIndex (double *a, const unsigned int n)
unsigned int	doubleMinIndex (double *a, const unsigned int n)
unsigned int	doubleAbsMinIndex (double *a, const unsigned int n)
unsigned int	doubleGEIndex (double *a, const unsigned int n, double lim)
unsigned int	doubleGTIndex (double *a, const unsigned int n, double lim)
double	inverfc (double x)

Detailed Description

Working with doubles.

Author

Vesa Oikonen

Definition in file doubleutil.c.

Function Documentation

doubleAbsMaxIndex()

unsigned int doubleAbsMaxIndex	(	double *	a,
		const unsigned int	n )

Find the maximum absolute value in given double array.

Returns

The index [0..n-1] of maximum absolute value in array; 0 is returned also in case of errors.

See also

doubleMaxIndex, doubleGTIndex, doubleGEIndex, doubleRange

Parameters

a	Pointer to double array.
n	Length of array.

Definition at line 377 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(0);
  unsigned int i, mi=0;
  double mv=nan("");
  for(i=0; i<n; i++) if(isnan(mv) || fabs(a[i])>mv) {mv=fabs(a[i]); mi=i;}
  return(mi);
}

doubleAbsMinIndex()

unsigned int doubleAbsMinIndex	(	double *	a,
		const unsigned int	n )

Find the minimum absolute value in given double array.

Returns

The index [0..n-1] of minimum absolute value in array; 0 is returned also in case of errors.

See also

doubleAbsMaxIndex, doubleMinIndex, doubleGTIndex, doubleGEIndex

Parameters

a	Pointer to double array.
n	Length of array.

Definition at line 416 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(0);
  unsigned int i, mi=0;
  double mv=nan("");
  for(i=0; i<n; i++) if(isnan(mv) || fabs(a[i])<mv) {mv=fabs(a[i]); mi=i;}
  return(mi);
}

Referenced by spectralBFExtract().

doubleArrayMatch()

int doubleArrayMatch	(	const double *	a1,
		const double *	a2,
		const unsigned int	n,
		const double	lim )

Verifies that given two double arrays have the same values inside given limits. Values are considered to match also if both are NaNs.

Returns

1 if match is found, and 0 if not.

Author

Vesa Oikonen

See also

doubleMatch, doubleMatchRel, doubleMachEps, statMeanSD

Parameters

a1	First array.
a2	Second array.
n	Length of arrays.
lim	Limit for absolute difference (if <0 then test will fail every time).

Definition at line 51 of file doubleutil.c.

  {
  if(a1==NULL || a2==NULL || n<1) return 0;
  for(unsigned int i=0; i<n; i++) if(!doubleMatch(a1[i], a2[i], lim)) return 0;
  return 1;
}

doubleCopy()

void doubleCopy	(	double *	t,
		double *	s,
		const unsigned int	n )

Copy double values from the 2nd array to the first.

See also

floatCopy, doubleCopyFinite

Parameters

t	Target array.
s	Source array.
n	Length of arrays.

Definition at line 117 of file doubleutil.c.

  {
  unsigned int i;
  if(t==NULL || s==NULL || n<1) return;
  for(i=0; i<n; i++) t[i]=s[i];
}

Referenced by nloptIATGO(), nloptMPSO(), and nloptSimplexARRS().

doubleCopyFinite()

unsigned int doubleCopyFinite	(	double *	t,
		double *	s,
		const unsigned int	n )

Copy finite double values from the 2nd array to the first.

See also

doubleCopy

Returns

The number of finite values that were copied.

Parameters

t	Target array.
s	Source array.
n	Length of arrays.

Definition at line 136 of file doubleutil.c.

  {
  unsigned int i, tn=0;
  if(t==NULL || s==NULL || n<1) return(0);
  for(i=0; i<n; i++) if(isfinite(s[i])) t[tn++]=s[i];
  return(tn);
}

doubleCSpanPositives()

int doubleCSpanPositives	(	double *	a,
		const int	n )

Returns the length of array consisting of non-positive (<=0 and NaN) values.

Returns

Returns the index of first positive value in the given array.

See also

floatCSpanPositives, doubleSpanPositives, atofList

Parameters

a	Pointer to the array.
n	Length of the array.

Definition at line 320 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(0);
  int i=0;
  for(i=0; i<n; i++) if(a[i]>0.0) break;
  return(i);
}

Referenced by spectralBFExtract(), and spectralBFNr().

doubleGEIndex()

unsigned int doubleGEIndex	(	double *	a,
		const unsigned int	n,
		double	lim )

Get the index of first array value greater or equal than given limit.

Returns

The index, or n, if not found.

See also

doubleGTIndex, doubleMaxIndex, doubleMinIndex, fitLine, statMeanSD

Parameters

a	Pointer to an array of values.
n	The number of samples (length of a[]).
lim	The limit value.

Definition at line 435 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(n);
  for(unsigned int i=0; i<n; i++) if(a[i]>=lim) return(i);
  return(n);
}

doubleGetWithUnit()

int doubleGetWithUnit	(	const char *	s,
		double *	v,
		int *	u )

Reads numerical value and its unit from argument string.

Both decimal point and comma are accepted. Optional result double value is set to NaN if string was not valid value. Optional unit is set to UNIT_UNKNOWN, if not valid or not found.

Returns

0 if numerical value successfully read, and 1 in case of an error.

See also

floatGetWithUnit, atofVerified, atofCheck

Parameters

s	String from where double and unit is read; string must not contain any extra space characters.
v	Pointer to the double float; enter NULL, if not needed.
u	Pointer to int for unit code; enter NULL, if not needed.

Definition at line 272 of file doubleutil.c.

  {
  if(v!=NULL) *v=nan(""); 
  if(u!=NULL) *u=UNIT_UNKNOWN; 
  if(s==NULL) return 1;
  int n=strTokenNr(s, " \t"); if(n==0 || n>2) return 1;
  if(n==1) return(atofCheck(s, v)); // no unit
  /* separate number and unit */
  char buf[64];
  n=strTokenNCpy(s, " \t", 1, buf, 64);
  if(atofCheck(buf, v)) return 1;
  if(u==NULL) return 0;
  n=strTokenNCpy(s, " \t", 2, buf, 64);
  *u=unitIdentify(buf);
  return 0;
}

Referenced by iftCheckKeyValue(), and iftCheckKeyValues().

doubleGTIndex()

unsigned int doubleGTIndex	(	double *	a,
		const unsigned int	n,
		double	lim )

Get the index of first array value greater than given limit.

Returns

The index, or n, if not found.

See also

doubleGEIndex, doubleMaxIndex, doubleMinIndex, fitLine, statMeanSD, doubleMatch

Parameters

a	Pointer to an array of values.
n	The number of samples (length of a[]).
lim	The limit value.

Definition at line 452 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(n);
  for(unsigned int i=0; i<n; i++) if(a[i]>lim) return(i);
  return(n);
}

doubleMachEps()

double doubleMachEps

(

)

Estimates the machine epsilon, the upper bound on the relative error due to rounding in floating point arithmetic, within one order of magnitude of the true machine epsilon.

Standard C library should also have DBL_EPSILON in float.h.

Returns

Estimate of machine epsilon.

Author

Vesa Oikonen

See also

floatMachEps, doubleMatchRel, doubleMatch, factorial

Definition at line 105 of file doubleutil.c.

{
  double macheps=1.0;
  do {macheps/=2.0;} while((1.0+macheps/2.0)!=1.0);
  return(macheps);
}

Referenced by householder_transform(), nlopt1D(), nloptFixedNr(), nloptLimitFixedNr(), and nloptPowellBrent().

doubleMatch()

int doubleMatch	(	const double	v1,
		const double	v2,
		const double	lim )

Verifies that given two doubles have the same value inside given limits. Values are considered to match also if both are NaNs.

Returns

1 if match is found, and 0 if not.

Author

Vesa Oikonen

See also

floatMatch, doubleMatchRel, doubleMachEps, doubleArrayMatch, atofCheck, statMeanSD

Parameters

v1	First value.
v2	Second value.
lim	Limit for absolute difference (if <0 then test will fail every time).

Definition at line 27 of file doubleutil.c.

  {
  if(isnan(v1) && isnan(v2)) return 1;
  if(isnan(v1) || isnan(v2)) return 0;
  if(v1==v2) return 1;
  if(isnan(lim) || lim<0.0) return 0;
  if(fabs(v1-v2)<=lim) return 1;
  return 0;
}

Referenced by doubleArrayMatch(), parCompareParameters(), parCompareWSS(), parIsFitRange(), tacCheckX1X2X(), tacCompareConc(), tacCompareTimes(), tacCompareWeights(), and tacReadQView().

doubleMatchRel()

int doubleMatchRel	(	const double	v1,
		const double	v2,
		const double	lim )

Verifies that given two doubles have the same value inside given relative limit |2*(v1-v2)/(v1+v2)|.

Values are considered to match also if both are NaNs, but not if mean is zero (unless both are exactly zero).

Returns

1 if match is found, and 0 if not.

Author

Vesa Oikonen

See also

floatMatchRel, doubleMatch, doubleMachEps

Parameters

v1	First value.
v2	Second value.
lim	Limit for relative difference (if <0 then test will fail every time).

Definition at line 77 of file doubleutil.c.

  {
  if(isnan(v1) && isnan(v2)) return 1;
  if(isnan(v1) || isnan(v2)) return 0;
  if(v1==v2) return 1;
  if(isnan(lim)) return 0;
  double mean;
  mean=0.5*(v1+v2); if(!isnormal(mean)) return 0;
  if(fabs((v1-v2)/mean)<=lim) return 1;
  return 0;
}

Referenced by parCompareParameters(), parCompareWSS(), tacCompareConc(), tacCompareTimes(), and tacCompareWeights().

doubleMaxIndex()

unsigned int doubleMaxIndex	(	double *	a,
		const unsigned int	n )

Find the maximum value in given double array.

Returns

The index [0..n-1] of maximum value in array; 0 is returned also in case of errors.

See also

floatMaxIndex, doubleSpanPositives, doubleCSpanPositives, doubleGTIndex, doubleRange

Parameters

a	Pointer to double array.
n	Length of array.

Definition at line 357 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(0);
  unsigned int i, mi=0;
  double mv=nan("");
  for(i=0; i<n; i++) if(isnan(mv) || a[i]>mv) {mv=a[i]; mi=i;}
  return(mi);
}

doubleMean()

double doubleMean	(	double *	a,
		const unsigned int	n )

Calculate the mean of values in given double array.

Returns

The mean of array values, or NaN in case of an error.

See also

floatMean, doubleWMean, doubleSum, doubleMaxIndex, statMeanSD, statMedian

Parameters

a	Pointer to double array.
n	Length of array.

Definition at line 224 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(nan(""));
  double s=0.0;
  unsigned int i, ci=0;
  for(i=0; i<n; i++) if(!isnan(a[i])) {ci++; s+=a[i];}
  if(ci<1) return(nan(""));
  return(s/(double)ci);
}

doubleMinIndex()

unsigned int doubleMinIndex	(	double *	a,
		const unsigned int	n )

Find the minimum value in given double array.

Returns

The index [0..n-1] of minimum value in array; 0 is returned also in case of errors.

See also

doubleMaxIndex, doubleAbsMinIndex, doubleRange

Parameters

a	Pointer to double array.
n	Length of array.

Definition at line 396 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(0);
  unsigned int i, mi=0;
  double mv=nan("");
  for(i=0; i<n; i++) if(isnan(mv) || a[i]<mv) {mv=a[i]; mi=i;}
  return(mi);
}

doubleNaNs()

unsigned int doubleNaNs	(	double *	a,
		const unsigned int	n )

Check the number of NaNs in a double array.

See also

doubleNonzeroes, doubleCopyFinite, doubleRange, doubleMean

Returns

The number of NaNs.

Parameters

a	Pointer to double array.
n	Length of array.

Definition at line 156 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(0);
  unsigned int i, nn=0;
  for(i=0; i<n; i++) if(isnan(a[i])) nn++;
  return(nn);
}

doubleNonzeroes()

unsigned int doubleNonzeroes	(	double *	a,
		const unsigned int	n )

Find the number of non-zero values in given double array. Negative values are counted, but NaNs are not.

Returns

The number of non-zero values in array; 0 is returned also in case of errors.

See also

floatNonzeroes, doubleSpanPositives, doubleNaNs

Parameters

a	Pointer to double array.
n	Length of array.

Definition at line 339 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(0);
  unsigned int i, m=0;
  for(i=0; i<n; i++) if(fabs(a[i])>0.0) m++;
  return(m);
}

doubleRange()

unsigned int doubleRange	(	double *	a,
		const unsigned int	n,
		double *	amin,
		double *	amax )

Get the range of finite values in a double array.

See also

doubleMaxIndex, doubleMinIndex, doubleNonzeroes, doubleNaNs, doubleMean, statMedian

Returns

The number of finite values in the array, or 0 in case of an error.

Parameters

a	Pointer to double array.
n	Length of array.
amin	Pointer to the double float where minimum value is stored; enter NULL, if not needed.
amax	Pointer to the double float where maximum value is stored; enter NULL, if not needed.

Definition at line 174 of file doubleutil.c.

  {
  if(amin!=NULL) *amin=nan("");
  if(amax!=NULL) *amax=nan("");
  if(a==NULL || n<1) return(0);
 
  double mi=nan(""), ma=nan("");
  unsigned int i, nn=0;
  for(i=0; i<n; i++) if(isfinite(a[i])) {
    nn++;
    if(!isfinite(mi) || mi>a[i]) mi=a[i];
    if(!isfinite(ma) || ma<a[i]) ma=a[i];
  }
  if(amin!=NULL) *amin=mi;
  if(amax!=NULL) *amax=ma;
  return(nn);
}

doubleSpanPositives()

int doubleSpanPositives	(	double *	a,
		const int	n )

Returns the length of array consisting of only positive (>0 and not NaN) values.

Returns

Returns the index of first non-positive value in the given array.

See also

floatSpanPositives, doubleCSpanPositives, atofList, doubleNonzeroes

Parameters

a	Pointer to the array.
n	Length of the array.

Definition at line 302 of file doubleutil.c.

  {
  if(a==NULL || n<1) return(0);
  int i=0;
  for(i=0; i<n; i++) if(!(a[i]>0.0)) break;
  return(i);
}

Referenced by spectralBFExtract(), and spectralBFNr().

doubleSum()

double doubleSum	(	double *	a,
		const unsigned int	n )

Calculate the sum of values in given double array.

Returns

The sum of array values.

See also

floatSum, doubleMean, doubleMaxIndex

Parameters

a	Pointer to double array.
n	Length of array.

Definition at line 206 of file doubleutil.c.

  {
  double s=0.0;
  if(a==NULL || n<1) return(s);
  for(unsigned int i=0; i<n; i++) if(!isnan(a[i])) s+=a[i];
  return(s);
}

Referenced by spectralBFExtract().

doubleWMean()

double doubleWMean	(	double *	a,
		double *	w,
		const unsigned int	n )

Calculate the weighted mean of values in given double array.

Returns

The mean of array values, or NaN in case of an error.

See also

doubleMean, floatMean, doubleSum, doubleMaxIndex, statMeanSD

Parameters

a	Pointer to double array to calculated weighted mean from.
w	Pointer to double array containing the weights. Negative or NaN weight is interpreted as zero weight.
n	Length of arrays.

Definition at line 244 of file doubleutil.c.

  {
  if(a==NULL || w==NULL || n<1) return(nan(""));
  double s=0.0, ws=0.0;
  unsigned int i, ci=0;
  for(i=0; i<n; i++) if(isnormal(a[i]) && isnormal(w[i]) && w[i]>0.0) {
    ci++; s+=w[i]*a[i]; ws+=w[i];}
  if(ci<1 || !(ws>0.0)) return(nan(""));
  return(s/ws);
}

Referenced by spectralBFExtract().

inverfc()

double inverfc

(

double

x

)

Inverse complementary error function erfc^{-1}(x).

Based on the code by Acklam PJ, 2010.

Returns

Returns an approximate value of erfc^{-1}(x).

See also

doubleMachEps, igam, igamc

Parameters

x	Parameter for the inverse complementary error function.

Definition at line 473 of file doubleutil.c.

  {
  static const double a[] = {
    -3.969683028665376E+01, 2.209460984245205E+02, -2.759285104469687E+02,
    1.383577518672690E+02, -3.066479806614716E+01,  2.506628277459239
  };
  static const double b[] = {
    -5.447609879822406E+01, 1.615858368580409E+02, -1.556989798598866E+02, 
    6.680131188771972E+01, -1.328068155288572E+01
  };
  static const double c[] = {
    -7.784894002430293E-03, -3.223964580411365E-01, -2.400758277161838, 
    -2.549732539343734, 4.374664141464968, 2.938163982698783
  };
  static const double d[] = {
    7.784695709041462E-03, 3.224671290700398E-01, 2.445134137142996, 3.754408661907416
  };
  double y, e, u;
  x/=2.0;
  if(x<0.02425) {
    double q=sqrt(-2.0 * log(x));
    y = (((((c[0]*q + c[1])*q + c[2])*q + c[3])*q + c[4])*q + c[5])
         / ((((d[0]*q + d[1])*q + d[2])*q + d[3])*q + 1.0);
  } else if(x<=0.97575) {
    double q=x-0.5;
    double r=q*q;
    y = (((((a[0]*r + a[1])*r + a[2])*r + a[3])*r + a[4])*r + a[5])*q
         / (((((b[0]*r + b[1])*r + b[2])*r + b[3])*r + b[4])*r + 1.0);
  } else {
    double q=sqrt(-2.0 * log(1.0-x));
    y = -(((((c[0]*q + c[1])*q + c[2])*q + c[3])*q + c[4])*q + c[5])
          / ((((d[0]*q + d[1])*q + d[2])*q + d[3])*q + 1.0);
  }
  e=0.5*erfc(-y/M_SQRT2)-x;
  u=(e)*M_SQRT2*M_PI*exp(0.5*y*y);
  y-=u/(1.0 + 0.5*y*u);
  return(fabs(y/M_SQRT2));
}