# Contents of /branches/arrexp_2137_win_merge/escript/src/UnaryFuncs.h

Revision 2213 - (show annotations)
Wed Jan 14 00:23:39 2009 UTC (10 years, 6 months ago) by jfenwick
File MIME type: text/plain
File size: 2934 byte(s)
```In preparation for merging to trunk

```
 1 2 /******************************************************* 3 * 4 * Copyright (c) 2003-2008 by University of Queensland 5 * Earth Systems Science Computational Center (ESSCC) 6 * http://www.uq.edu.au/esscc 7 * 8 * Primary Business: Queensland, Australia 9 * Licensed under the Open Software License version 3.0 10 * http://www.opensource.org/licenses/osl-3.0.php 11 * 12 *******************************************************/ 13 14 15 #if !defined escript_UnaryFuncs_20041124_H 16 #define escript_UnaryFuncs_20041124_H 17 #include "system_dep.h" 18 19 namespace escript { 20 21 #ifndef FP_NAN 22 #define FP_NAN IEEE_NaN() 23 #endif 24 25 #ifndef INFINITY 26 #define INFINITY IEEE_Infy() 27 #endif 28 29 //====================================================================== 30 31 inline 32 double log1p (const double x) 33 { 34 volatile double y; 35 y = 1 + x; 36 return log(y) - ((y-1)-x)/y ; 37 } 38 39 //====================================================================== 40 41 inline 42 float IEEE_NaN() 43 { 44 static unsigned char nan[4]={ 0, 0, 0xc0, 0x7f }; 45 return *( float *)nan; 46 } 47 48 //====================================================================== 49 50 inline 51 double IEEE_Infy() 52 { 53 static unsigned char infy[8]={ 0, 0, 0, 0, 0, 0, 0xf0, 0x7f }; 54 return *( double *)infy; 55 } 56 57 58 //====================================================================== 59 60 #if defined (_WIN32) && !defined(__INTEL_COMPILER) 61 inline 62 double 63 acosh_substitute (const double x) 64 { 65 if (x > 1.0 / SQRT_DBL_EPSILON) 66 { 67 return log (x) + M_LN2; 68 } 69 else if (x > 2) 70 { 71 return log (2 * x - 1 / (sqrt (x * x - 1) + x)); 72 } 73 else if (x > 1) 74 { 75 double t = x - 1; 76 return log1p (t + sqrt (2 * t + t * t)); 77 } 78 else if (x == 1) 79 { 80 return 0; 81 } 82 else 83 { 84 return FP_NAN; 85 } 86 } 87 88 89 //====================================================================== 90 91 inline 92 double 93 asinh_substitute (const double x) 94 { 95 double a = fabs (x); 96 double s = (x < 0) ? -1 : 1; 97 98 if (a > 1 / SQRT_DBL_EPSILON) 99 { 100 return s * (log (a) + M_LN2); 101 } 102 else if (a > 2) 103 { 104 return s * log (2 * a + 1 / (a + sqrt (a * a + 1))); 105 } 106 else if (a > SQRT_DBL_EPSILON) 107 { 108 double a2 = a * a; 109 return s * log1p (a + a2 / (1 + sqrt (1 + a2))); 110 } 111 else 112 { 113 return x; 114 } 115 } 116 117 118 //====================================================================== 119 120 inline 121 double 122 atanh_substitute (const double x) 123 { 124 double a = fabs (x); 125 double s = (x < 0) ? -1 : 1; 126 127 if (a > 1) 128 { 129 return FP_NAN; 130 } 131 else if (a == 1) 132 { 133 return (x < 0) ? -INFINITY : INFINITY; 134 } 135 else if (a >= 0.5) 136 { 137 return s * 0.5 * log1p (2 * a / (1 - a)); 138 } 139 else if (a > DBL_EPSILON) 140 { 141 return s * 0.5 * log1p (2 * a + 2 * a * a / (1 - a)); 142 } 143 else 144 { 145 return x; 146 } 147 } 148 #endif // windows substitutes for stupid microsoft compiler. 149 150 151 inline 152 double 153 fsign(double x) 154 { 155 if (x == 0) { 156 return 0; 157 } else { 158 return x/fabs(x); 159 } 160 } 161 162 } 163 164 #endif

## Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

 ViewVC Help Powered by ViewVC 1.1.26