import changes from `dev' branch of rmottola/Arctic-Fox:

- Bug 1265956 - Assert that no entry is found in HashTable::putNew. (r=terrence) (8770735325)
- Bug 1265483 - Use WeakCache to automate sweeping of ObjectGroupComparment::NewTable; r=jonco (e88fa842ab)
- Bug 1198093 - Part 1: Expose indexedDB to System with [Exposed=System]. r=khuey (2a20a6ecd0)
- Bug 1198093 - Part 2: Set Default Locale Value in ICU Canonicalization Form. r=khuey (1e1cd981ce)
- Bug 1263871 - Fix OOM handling in while resolving function name r=shu (44114a7e8d)
- Bug 1263270 - Sort census reports by smallest node ID counted, rather than number of nodes counted. r=jimb (08d4a431a7)
- Bug 1263218 - Fix possbile race under oomTest involving background threads r=terrence (82c1e3b698)
- Bug 1256488 - Add explicit casts to fix MSVC warning C4365. r=fitzgen (b9bb6b52c5)
- Bug 1235677 - Add assertion to catch unsafe concurrent use of AutoEnterOOMUnsafeRegion r=terrence (2e0876578f)
- Bug 1263902 - check return value from JS_smprintf. r=bbouvier, r=shu (df1d9b5e14)
- Bug 1262208: Generalize the disabled compilation mode message for asm.js; r=luke (9153b2c5ba)
- Bug 1259903: Baldr: unify Select true and false types instead of checking against each other; r=luke (5f89398199)
- Bug 1253344: Remove unused pushPhi/popPhi in WasmIonCompile; r=luke (244967401c)
- Bug 933257 - Part 1: Add a script to import and update fdlibm from FreeBSD. r=jwalden (37c8a85771)
- Bug 933257 - Part 2: Add patches for fdlibm. r=jwalden (bc0dce94a0)
- Bug 933257 - Part 2.1: Import fdlibm from FreeBSD (revision bcea9d50b15e4f0027a5dd526e0e2a612238471e). r=jwalden (223f6d6ce5)
- Bug 933257 - Part 3: Add build scripts for fdlibm. r=jwalden (893f740423)
- Bug 933257 - Part 4: Link fdlibm in SpiderMonkey. r=jwalden (9f1395258a)
- Bug 933257 - Part 5: Use fdlibm in jsmath.cpp. r=jwalden (9d962657ab)
- Bug 933257 - Part 6: Remove unused math polyfill. r=jwalden (cf284ad4e5)
- Bug 1225024 - Allow sloppy tolerance in ecma_6/Math/log10-approx.js. r=jorendorff (7df3bf46dd)
- Bug 933257 - Part 7: Remove or reduce sloppy_tolerance in Math function tests. r=jorendorff (86b978eb14)
- Bug 933257 - Part 8: Add license for k_exp.cpp to about:license. r=gerv (513012fbb9)
- Bug 933257 - Part 9: Use fdlibm in asm.js. r=luke (46bedc10d0)
- Bug 1256490 - Disable C4302 to unblock compilation on VS2015; r=bobowen (7fb6820241)
- Bug 1256499 - Disable C4311 and C4312 to unblock compilation on VS2015; r=bobowen (a9b3b01410)
- Bug 1257036 - Disable C4302 to unblock compilation on VS2015; r=bobowen (307af58682)
- Bug 1124033 - Disable C4311 and C4312 in directories exhibiting warnings; r=ehsan (b6ecd1f8e7)
- Bug 1252931 - Remove INSTALL/PP_TARGETS from js/src/*; r=gps (a1e1d3bb82)
- Bug 1258908: Rename TYPE_MOZILLA_UI to TYPE_MOZILLA_PARENT. r=jld (29aef56b8e)
- Bug 1203835 - Don't ship replace_jemalloc. r=njn (249f927cf5)
- bug 1259753 - fix some C++ unittests to use ScopedXPCOM to init XPCOM. r=ms2ger (a908216277)
This commit is contained in:
2024-04-30 11:53:38 +08:00
parent 4696208ad8
commit b5ee72fe09
126 changed files with 7639 additions and 536 deletions
-3
View File
@@ -100,9 +100,6 @@ endif
ifdef MOZ_SHARED_ICU
DEFINES += -DMOZ_SHARED_ICU
endif
ifdef MOZ_JEMALLOC4
DEFINES += -DMOZ_JEMALLOC4
endif
DEFINES += -DMOZ_ICU_DBG_SUFFIX=$(MOZ_ICU_DBG_SUFFIX)
ifdef MOZ_WIDGET_GTK
-5
View File
@@ -113,11 +113,6 @@
#ifdef MOZ_SHARED_MOZGLUE
@BINPATH@/@DLL_PREFIX@mozglue@DLL_SUFFIX@
#endif
#ifdef MOZ_REPLACE_MALLOC
#ifndef MOZ_JEMALLOC4
@BINPATH@/@DLL_PREFIX@replace_jemalloc@DLL_SUFFIX@
#endif
#endif
#ifdef ANDROID
@RESPATH@/AndroidManifest.xml
@RESPATH@/resources.arsc
-3
View File
@@ -140,9 +140,6 @@ endif
ifdef MOZ_SHARED_ICU
DEFINES += -DMOZ_SHARED_ICU
endif
ifdef MOZ_JEMALLOC4
DEFINES += -DMOZ_JEMALLOC4
endif
DEFINES += -DMOZ_ICU_DBG_SUFFIX=$(MOZ_ICU_DBG_SUFFIX)
ifdef CLANG_CXX
DEFINES += -DCLANG_CXX
-5
View File
@@ -131,11 +131,6 @@
#endif
#endif
#endif
#ifdef MOZ_REPLACE_MALLOC
#ifndef MOZ_JEMALLOC4
@BINPATH@/@DLL_PREFIX@replace_jemalloc@DLL_SUFFIX@
#endif
#endif
#ifdef MOZ_GTK3
@BINPATH@/@DLL_PREFIX@mozgtk@DLL_SUFFIX@
@BINPATH@/gtk2/@DLL_PREFIX@mozgtk@DLL_SUFFIX@
+1
View File
@@ -66,6 +66,7 @@ included_inclnames_to_ignore = set([
'jsautokw.h', # generated in $OBJDIR
'jscustomallocator.h', # provided by embedders; allowed to be missing
'js-config.h', # generated in $OBJDIR
'fdlibm.h', # fdlibm
'pratom.h', # NSPR
'prcvar.h', # NSPR
'prerror.h', # NSPR
+14
View File
@@ -0,0 +1,14 @@
# -*- Mode: python; c-basic-offset: 4; indent-tabs-mode: nil; tab-width: 40 -*-
# vim: set filetype=python:
# This Source Code Form is subject to the terms of the Mozilla Public
# License, v. 2.0. If a copy of the MPL was not distributed with this
# file, You can obtain one at http://mozilla.org/MPL/2.0/.
Library('fdlibm')
with Files('**'):
BUG_COMPONENT = ('Core', 'JavaScript Engine')
DIRS += [
'../../../modules/fdlibm',
]
+3
View File
@@ -35,6 +35,9 @@ if CONFIG['GKMEDIAS_SHARED_LIBRARY']:
if CONFIG['MOZ_SHARED_ICU']:
DEFINES['MOZ_SHARED_ICU'] = True
if CONFIG['MOZ_SYSTEM_ICU']:
DEFINES['MOZ_SYSTEM_ICU'] = True
PYTHON_UNIT_TESTS += [
'tests/test_mozbuild_reading.py',
'tests/unit-expandlibs.py',
+5
View File
@@ -183,3 +183,8 @@ CXXFLAGS += CONFIG['MOZ_CAIRO_CFLAGS']
CXXFLAGS += CONFIG['TK_CFLAGS']
LOCAL_INCLUDES += CONFIG['SKIA_INCLUDES']
if CONFIG['_MSC_VER']:
# This is intended as a temporary workaround to unblock compilation
# on VS2015 in warnings as errors mode.
CXXFLAGS += ['-wd4312']
+39 -3
View File
@@ -424,7 +424,7 @@ IndexedDatabaseManager::Init()
}
if (mLocale.IsEmpty()) {
mLocale.AssignLiteral("en-US");
mLocale.AssignLiteral("en_US");
}
#endif
@@ -584,8 +584,8 @@ IndexedDatabaseManager::CommonPostHandleEvent(EventChainPostVisitor& aVisitor,
// static
bool
IndexedDatabaseManager::DefineIndexedDB(JSContext* aCx,
JS::Handle<JSObject*> aGlobal)
IndexedDatabaseManager::ResolveSandboxBinding(JSContext* aCx,
JS::Handle<JSObject*> aGlobal)
{
MOZ_ASSERT(NS_IsMainThread());
MOZ_ASSERT(js::GetObjectClass(aGlobal)->flags & JSCLASS_DOM_GLOBAL,
@@ -614,6 +614,18 @@ IndexedDatabaseManager::DefineIndexedDB(JSContext* aCx,
return false;
}
return true;
}
// static
bool
IndexedDatabaseManager::DefineIndexedDB(JSContext* aCx,
JS::Handle<JSObject*> aGlobal)
{
MOZ_ASSERT(NS_IsMainThread());
MOZ_ASSERT(js::GetObjectClass(aGlobal)->flags & JSCLASS_DOM_GLOBAL,
"Passed object is not a global object!");
RefPtr<IDBFactory> factory;
if (NS_FAILED(IDBFactory::CreateForMainThreadJS(aCx,
aGlobal,
@@ -722,6 +734,30 @@ IndexedDatabaseManager::ExperimentalFeaturesEnabled()
return gExperimentalFeaturesEnabled;
}
// static
bool
IndexedDatabaseManager::ExperimentalFeaturesEnabled(JSContext* aCx, JSObject* aGlobal)
{
// If, in the child process, properties of the global object are enumerated
// before the chrome registry (and thus the value of |intl.accept_languages|)
// is ready, calling IndexedDatabaseManager::Init will permanently break
// that preference. We can retrieve gExperimentalFeaturesEnabled without
// actually going through IndexedDatabaseManager.
// See Bug 1198093 comment 14 for detailed explanation.
if (IsNonExposedGlobal(aCx, js::GetGlobalForObjectCrossCompartment(aGlobal),
GlobalNames::BackstagePass)) {
MOZ_ASSERT(NS_IsMainThread());
static bool featureRetrieved = false;
if (!featureRetrieved) {
gExperimentalFeaturesEnabled = Preferences::GetBool(kPrefExperimental);
featureRetrieved = true;
}
return gExperimentalFeaturesEnabled;
}
return ExperimentalFeaturesEnabled();
}
// static
bool
IndexedDatabaseManager::IsFileHandleEnabled()
+4 -4
View File
@@ -116,10 +116,7 @@ public:
ExperimentalFeaturesEnabled();
static bool
ExperimentalFeaturesEnabled(JSContext* /* aCx */, JSObject* /* aGlobal */)
{
return ExperimentalFeaturesEnabled();
}
ExperimentalFeaturesEnabled(JSContext* aCx, JSObject* aGlobal);
static bool
IsFileHandleEnabled();
@@ -181,6 +178,9 @@ public:
static nsresult
CommonPostHandleEvent(EventChainPostVisitor& aVisitor, IDBFactory* aFactory);
static bool
ResolveSandboxBinding(JSContext* aCx, JS::Handle<JSObject*> aGlobal);
static bool
DefineIndexedDB(JSContext* aCx, JS::Handle<JSObject*> aGlobal);
+5
View File
@@ -321,4 +321,9 @@ include('/ipc/chromium/chromium-config.mozbuild')
if CONFIG['GNU_CC']:
CXXFLAGS += ['-Wno-error=attributes']
if CONFIG['_MSC_VER']:
# This is intended as a temporary workaround to unblock compilation
# on VS2015 in warnings as errors mode.
CXXFLAGS += ['-wd4312']
FINAL_LIBRARY = 'xul'
+5
View File
@@ -57,6 +57,11 @@ if CONFIG['MOZ_WIDGET_TOOLKIT'] == 'gonk':
if CONFIG['_MSC_VER']:
DEFINES['__PRETTY_FUNCTION__'] = '__FUNCSIG__'
# This is intended as a temporary workaround to enable building with VS2015.
# media\webrtc\trunk\webrtc/base/criticalsection.h(59): warning C4312:
# 'reinterpret_cast': conversion from 'DWORD' to 'HANDLE' of greater size
CXXFLAGS += ['-wd4312']
EXPORTS.mozilla += ['ShmemPool.h',]
EXPORTS.mozilla.media += ['MediaChild.h',
+2 -2
View File
@@ -5,7 +5,7 @@
enum DOMRequestReadyState { "pending", "done" };
[Exposed=(Window,Worker), NoInterfaceObject]
[Exposed=(Window,Worker,System), NoInterfaceObject]
interface DOMRequestShared {
readonly attribute DOMRequestReadyState readyState;
@@ -16,7 +16,7 @@ interface DOMRequestShared {
attribute EventHandler onerror;
};
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface DOMRequest : EventTarget {
// The [TreatNonCallableAsNull] annotation is required since then() should do
// nothing instead of throwing errors when non-callable arguments are passed.
+1 -1
View File
@@ -10,7 +10,7 @@
* liability, trademark and document use rules apply.
*/
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface DOMStringList {
readonly attribute unsigned long length;
getter DOMString? item(unsigned long index);
+2 -2
View File
@@ -14,7 +14,7 @@ enum IDBCursorDirection {
"prevunique"
};
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface IDBCursor {
readonly attribute (IDBObjectStore or IDBIndex) source;
@@ -39,7 +39,7 @@ interface IDBCursor {
IDBRequest delete ();
};
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface IDBCursorWithValue : IDBCursor {
[Throws]
readonly attribute any value;
+1 -1
View File
@@ -10,7 +10,7 @@
* liability, trademark and document use rules apply.
*/
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface IDBDatabase : EventTarget {
readonly attribute DOMString name;
readonly attribute unsigned long long version;
+1 -1
View File
@@ -23,7 +23,7 @@ dictionary IDBOpenDBOptions
* http://dvcs.w3.org/hg/IndexedDB/raw-file/tip/Overview.html#idl-def-IDBFactory
* for more information.
*/
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface IDBFactory {
[Throws]
IDBOpenDBRequest
+1
View File
@@ -8,6 +8,7 @@ dictionary IDBFileMetadataParameters
boolean lastModified = true;
};
[Exposed=(Window,System)]
interface IDBFileHandle : EventTarget
{
readonly attribute IDBMutableFile? mutableFile;
+1
View File
@@ -4,6 +4,7 @@
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, You can obtain one at http://mozilla.org/MPL/2.0/. */
[Exposed=(Window,System)]
interface IDBFileRequest : DOMRequest {
readonly attribute IDBFileHandle? fileHandle;
// this is deprecated due to renaming in the spec
+1 -1
View File
@@ -18,7 +18,7 @@ dictionary IDBIndexParameters {
DOMString? locale = null;
};
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface IDBIndex {
readonly attribute DOMString name;
readonly attribute IDBObjectStore objectStore;
+2 -2
View File
@@ -9,7 +9,7 @@
* liability, trademark and document use rules apply.
*/
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface IDBKeyRange {
[Throws]
readonly attribute any lower;
@@ -33,7 +33,7 @@ interface IDBKeyRange {
static IDBKeyRange bound (any lower, any upper, optional boolean lowerOpen = false, optional boolean upperOpen = false);
};
[Exposed=(Window,Worker),
[Exposed=(Window,Worker,System),
Func="mozilla::dom::IndexedDatabaseManager::ExperimentalFeaturesEnabled"]
interface IDBLocaleAwareKeyRange : IDBKeyRange {
[NewObject, Throws]
+1
View File
@@ -3,6 +3,7 @@
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, You can obtain one at http://mozilla.org/MPL/2.0/. */
[Exposed=(Window,System)]
interface IDBMutableFile : EventTarget {
readonly attribute DOMString name;
readonly attribute DOMString type;
+1 -1
View File
@@ -12,7 +12,7 @@ dictionary IDBObjectStoreParameters {
boolean autoIncrement = false;
};
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface IDBObjectStore {
readonly attribute DOMString name;
+1 -1
View File
@@ -7,7 +7,7 @@
* https://dvcs.w3.org/hg/IndexedDB/raw-file/tip/Overview.html#idl-def-IDBOpenDBRequest
*/
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface IDBOpenDBRequest : IDBRequest {
attribute EventHandler onblocked;
+1 -1
View File
@@ -13,7 +13,7 @@ enum IDBRequestReadyState {
"done"
};
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface IDBRequest : EventTarget {
[Throws]
readonly attribute any result;
+1 -1
View File
@@ -18,7 +18,7 @@ enum IDBTransactionMode {
"versionchange"
};
[Exposed=(Window,Worker)]
[Exposed=(Window,Worker,System)]
interface IDBTransaction : EventTarget {
[Throws]
readonly attribute IDBTransactionMode mode;
+1 -1
View File
@@ -16,7 +16,7 @@ dictionary IDBVersionChangeEventInit : EventInit {
};
[Constructor(DOMString type, optional IDBVersionChangeEventInit eventInitDict),
Exposed=(Window,Worker)]
Exposed=(Window,Worker,System)]
interface IDBVersionChangeEvent : Event {
readonly attribute unsigned long long oldVersion;
readonly attribute unsigned long long? newVersion;
+1 -1
View File
@@ -106,7 +106,7 @@ MessageLoop::MessageLoop(Type type)
get_tls_ptr().Set(this);
switch (type_) {
case TYPE_MOZILLA_UI:
case TYPE_MOZILLA_PARENT:
pump_ = new mozilla::ipc::MessagePump();
return;
case TYPE_MOZILLA_CHILD:
+3 -3
View File
@@ -199,7 +199,7 @@ public:
// This type of ML is used in Mozilla child processes which initialize
// XPCOM and use the gecko event loop.
//
// TYPE_MOZILLA_UI
// TYPE_MOZILLA_PARENT
// This type of ML is used in Mozilla parent processes which initialize
// XPCOM and use the gecko event loop.
//
@@ -216,7 +216,7 @@ public:
TYPE_UI,
TYPE_IO,
TYPE_MOZILLA_CHILD,
TYPE_MOZILLA_UI,
TYPE_MOZILLA_PARENT,
TYPE_MOZILLA_NONMAINTHREAD,
TYPE_MOZILLA_NONMAINUITHREAD
};
@@ -473,7 +473,7 @@ class MessageLoopForUI : public MessageLoop {
return NULL;
Type type = loop->type();
DCHECK(type == MessageLoop::TYPE_UI ||
type == MessageLoop::TYPE_MOZILLA_UI ||
type == MessageLoop::TYPE_MOZILLA_PARENT ||
type == MessageLoop::TYPE_MOZILLA_CHILD);
return static_cast<MessageLoopForUI*>(loop);
}
+6
View File
@@ -136,6 +136,12 @@ SOURCES += [
'URIUtils.cpp',
]
if CONFIG['_MSC_VER']:
# This is intended as a temporary hack to support building with VS2015.
# 'reinterpret_cast': conversion from 'DWORD' to 'HANDLE' of greater size
SOURCES['BackgroundChildImpl.cpp'].flags += ['-wd4312']
SOURCES['BackgroundParentImpl.cpp'].flags += ['-wd4312']
LOCAL_INCLUDES += [
'/caps',
'/dom/broadcastchannel',
+31 -18
View File
@@ -1563,6 +1563,33 @@ class HashTable : private AllocPolicy
// which approach is best.
}
// Note: |l| may be a reference to a piece of |u|, so this function
// must take care not to use |l| after moving |u|.
//
// Prefer to use putNewInfallible; this function does not check
// invariants.
template <typename... Args>
void putNewInfallibleInternal(const Lookup& l, Args&&... args)
{
MOZ_ASSERT(table);
HashNumber keyHash = prepareHash(l);
Entry* entry = &findFreeEntry(keyHash);
MOZ_ASSERT(entry);
if (entry->isRemoved()) {
METER(stats.addOverRemoved++);
removedCount--;
keyHash |= sCollisionBit;
}
entry->setLive(keyHash, mozilla::Forward<Args>(args)...);
entryCount++;
#ifdef JS_DEBUG
mutationCount++;
#endif
}
public:
void clear()
{
@@ -1703,23 +1730,9 @@ class HashTable : private AllocPolicy
template <typename... Args>
void putNewInfallible(const Lookup& l, Args&&... args)
{
MOZ_ASSERT(table);
HashNumber keyHash = prepareHash(l);
Entry* entry = &findFreeEntry(keyHash);
MOZ_ASSERT(entry);
if (entry->isRemoved()) {
METER(stats.addOverRemoved++);
removedCount--;
keyHash |= sCollisionBit;
}
entry->setLive(keyHash, mozilla::Forward<Args>(args)...);
entryCount++;
#ifdef JS_DEBUG
mutationCount++;
#endif
MOZ_ASSERT(!lookup(l).found());
mozilla::ReentrancyGuard g(*this);
putNewInfallibleInternal(l, mozilla::Forward<Args>(args)...);
}
// Note: |l| may be alias arguments in |args|, so this function must take
@@ -1771,7 +1784,7 @@ class HashTable : private AllocPolicy
typename HashTableEntry<T>::NonConstT t(mozilla::Move(*p));
HashPolicy::setKey(t, const_cast<Key&>(k));
remove(*p.entry_);
putNewInfallible(l, mozilla::Move(t));
putNewInfallibleInternal(l, mozilla::Move(t));
}
void rekeyAndMaybeRehash(Ptr p, const Lookup& l, const Key& k)
+3
View File
@@ -31,6 +31,7 @@ class WeakCache : public js::WeakCacheBase<T>,
friend class mozilla::LinkedListElement<WeakCache<T>>;
friend class mozilla::LinkedList<WeakCache<T>>;
WeakCache() = delete;
WeakCache(const WeakCache&) = delete;
using SweepFn = void (*)(T*);
@@ -38,6 +39,8 @@ class WeakCache : public js::WeakCacheBase<T>,
T cache;
public:
using Type = T;
template <typename U>
WeakCache(Zone* zone, U&& initial)
: cache(mozilla::Forward<U>(initial))
+28 -2
View File
@@ -9,6 +9,8 @@
#include "mozilla/Move.h"
#include <algorithm>
#include "jsapi.h"
#include "js/UbiNode.h"
@@ -125,11 +127,31 @@ class JS_FRIEND_API(CountBase) {
~CountBase() { }
public:
explicit CountBase(CountType& type) : type(type), total_(0) { }
explicit CountBase(CountType& type)
: type(type)
, total_(0)
, smallestNodeIdCounted_(SIZE_MAX)
{ }
// Categorize and count |node| as appropriate for this count's type.
bool count(mozilla::MallocSizeOf mallocSizeOf, const Node& node) {
return type.count(*this, mallocSizeOf, node);
total_++;
auto id = node.identifier();
if (id < smallestNodeIdCounted_) {
smallestNodeIdCounted_ = id;
}
#ifdef DEBUG
size_t oldTotal = total_;
#endif
bool ret = type.count(*this, mallocSizeOf, node);
MOZ_ASSERT(total_ == oldTotal,
"CountType::count should not increment total_, CountBase::count handles that");
return ret;
}
// Construct a JavaScript object reporting the counts recorded in this
@@ -147,6 +169,10 @@ class JS_FRIEND_API(CountBase) {
void trace(JSTracer* trc) { type.traceCount(*this, trc); }
size_t total_;
// The smallest JS::ubi::Node::identifier() passed to this instance's
// count() method. This provides a stable way to sort sets.
Node::Id smallestNodeIdCounted_;
};
class RootedCount : JS::CustomAutoRooter {
+14 -18
View File
@@ -8,6 +8,7 @@
#define js_Utility_h
#include "mozilla/Assertions.h"
#include "mozilla/Atomics.h"
#include "mozilla/Attributes.h"
#include "mozilla/Compiler.h"
#include "mozilla/Move.h"
@@ -123,36 +124,26 @@ extern JS_PUBLIC_DATA(uint64_t) maxAllocations;
extern JS_PUBLIC_DATA(uint64_t) counter;
extern JS_PUBLIC_DATA(bool) failAlways;
static inline void
SimulateOOMAfter(uint64_t allocations, uint32_t thread, bool always) {
MOZ_ASSERT(counter + allocations > counter);
MOZ_ASSERT(thread > js::oom::THREAD_TYPE_NONE && thread < js::oom::THREAD_TYPE_MAX);
targetThread = thread;
maxAllocations = counter + allocations;
failAlways = always;
}
extern void
SimulateOOMAfter(uint64_t allocations, uint32_t thread, bool always);
static inline void
ResetSimulatedOOM() {
targetThread = THREAD_TYPE_NONE;
maxAllocations = UINT64_MAX;
failAlways = false;
}
extern void
ResetSimulatedOOM();
static inline bool
inline bool
IsThreadSimulatingOOM()
{
return js::oom::targetThread && js::oom::targetThread == js::oom::GetThreadType();
}
static inline bool
inline bool
IsSimulatedOOMAllocation()
{
return IsThreadSimulatingOOM() &&
(counter == maxAllocations || (counter > maxAllocations && failAlways));
}
static inline bool
inline bool
ShouldFailWithOOM()
{
if (!IsThreadSimulatingOOM())
@@ -166,7 +157,7 @@ ShouldFailWithOOM()
return false;
}
static inline bool
inline bool
HadSimulatedOOM() {
return counter >= maxAllocations;
}
@@ -212,6 +203,7 @@ struct MOZ_RAII AutoEnterOOMUnsafeRegion
oomAfter_(0)
{
if (oomEnabled_) {
MOZ_ALWAYS_TRUE(owner_.compareExchange(nullptr, this));
oomAfter_ = int64_t(oom::maxAllocations) - int64_t(oom::counter);
oom::maxAllocations = UINT64_MAX;
}
@@ -224,10 +216,14 @@ struct MOZ_RAII AutoEnterOOMUnsafeRegion
MOZ_ASSERT(maxAllocations >= 0,
"alloc count + oom limit exceeds range, your oom limit is probably too large");
oom::maxAllocations = uint64_t(maxAllocations);
MOZ_ALWAYS_TRUE(owner_.compareExchange(this, nullptr));
}
}
private:
// Used to catch concurrent use from other threads.
static mozilla::Atomic<AutoEnterOOMUnsafeRegion*> owner_;
bool oomEnabled_;
int64_t oomAfter_;
#endif
+2 -2
View File
@@ -556,7 +556,7 @@ JSVAL_IS_BOOLEAN_IMPL(jsval_layout l)
static inline bool
JSVAL_TO_BOOLEAN_IMPL(jsval_layout l)
{
return l.s.payload.boo;
return bool(l.s.payload.boo);
}
static inline jsval_layout
@@ -564,7 +564,7 @@ BOOLEAN_TO_JSVAL_IMPL(bool b)
{
jsval_layout l;
l.s.tag = JSVAL_TAG_BOOLEAN;
l.s.payload.boo = b;
l.s.payload.boo = uint32_t(b);
return l;
}
+7 -1
View File
@@ -7419,6 +7419,8 @@ CheckBuffer(JSContext* cx, AsmJSModule& module, HandleValue bufferVal,
"valid length is 0x%x",
heapLength,
RoundUpToNextValidAsmJSHeapLength(heapLength)));
if (!msg)
return false;
return LinkFail(cx, msg.get());
}
@@ -7431,6 +7433,8 @@ CheckBuffer(JSContext* cx, AsmJSModule& module, HandleValue bufferVal,
"by const heap accesses).",
heapLength,
module.minHeapLength()));
if (!msg)
return false;
return LinkFail(cx, msg.get());
}
@@ -8131,6 +8135,8 @@ LookupAsmJSModuleInCache(ExclusiveContext* cx, AsmJSParser& parser, bool* loaded
int64_t usecAfter = PRMJ_Now();
int ms = (usecAfter - usecBefore) / PRMJ_USEC_PER_MSEC;
*compilationTimeReport = UniqueChars(JS_smprintf("loaded from cache in %dms", ms));
if (!*compilationTimeReport)
return false;
return true;
}
@@ -8162,7 +8168,7 @@ EstablishPreconditions(ExclusiveContext* cx, AsmJSParser& parser)
switch (parser.options().asmJSOption) {
case AsmJSOption::Disabled:
return Warn(parser, JSMSG_USE_ASM_TYPE_FAIL, "Disabled by javascript.options.asmjs in about:config");
return Warn(parser, JSMSG_USE_ASM_TYPE_FAIL, "Disabled by 'asmjs' runtime option");
case AsmJSOption::DisabledByDebugger:
return Warn(parser, JSMSG_USE_ASM_TYPE_FAIL, "Disabled by debugger");
case AsmJSOption::Enabled:
+3 -6
View File
@@ -459,7 +459,7 @@ DecodeConversionOperator(FunctionDecoder& f, ValType to, ValType argType, ExprTy
}
static bool
DecodeSelectOperator(FunctionDecoder& f, ExprType* type)
DecodeSelect(FunctionDecoder& f, ExprType* type)
{
ExprType trueType;
if (!DecodeExpr(f, &trueType))
@@ -472,9 +472,6 @@ DecodeSelectOperator(FunctionDecoder& f, ExprType* type)
if (!DecodeExpr(f, &falseType))
return false;
if (!CheckType(f, falseType, trueType))
return false;
ExprType condType;
if (!DecodeExpr(f, &condType))
return false;
@@ -482,7 +479,7 @@ DecodeSelectOperator(FunctionDecoder& f, ExprType* type)
if (!CheckType(f, condType, ValType::I32))
return false;
*type = trueType;
*type = Unify(trueType, falseType);
return true;
}
@@ -678,7 +675,7 @@ DecodeExpr(FunctionDecoder& f, ExprType* type)
case Expr::SetLocal:
return DecodeSetLocal(f, type);
case Expr::Select:
return DecodeSelectOperator(f, type);
return DecodeSelect(f, type);
case Expr::Block:
return DecodeBlock(f, /* isLoop */ false, type);
case Expr::Loop:
+9 -22
View File
@@ -1062,27 +1062,6 @@ class FunctionCompiler
return true;
}
bool usesPhiSlot()
{
return curBlock_->stackDepth() == info().firstStackSlot() + 1;
}
void pushPhiInput(MDefinition* def)
{
if (inDeadCode())
return;
MOZ_ASSERT(!usesPhiSlot());
curBlock_->push(def);
}
MDefinition* popPhiOutput()
{
if (inDeadCode())
return nullptr;
MOZ_ASSERT(usesPhiSlot());
return curBlock_->pop();
}
bool startBlock()
{
MOZ_ASSERT_IF(blockDepth_ < blockPatches_.length(), blockPatches_[blockDepth_].empty());
@@ -2287,7 +2266,15 @@ EmitSelect(FunctionCompiler& f, MDefinition** def)
if (!EmitExpr(f, &condExpr))
return false;
*def = f.select(trueExpr, falseExpr, condExpr);
if (trueExpr && falseExpr &&
trueExpr->type() == falseExpr->type() &&
trueExpr->type() != MIRType_None)
{
*def = f.select(trueExpr, falseExpr, condExpr);
} else {
*def = nullptr;
}
return true;
}
+11 -9
View File
@@ -18,6 +18,8 @@
#include "asmjs/WasmTypes.h"
#include "fdlibm.h"
#include "jslibmath.h"
#include "jsmath.h"
@@ -286,23 +288,23 @@ wasm::AddressOf(SymbolicAddress imm, ExclusiveContext* cx)
case SymbolicAddress::TanD:
return FuncCast<double (double)>(tan, Args_Double_Double);
case SymbolicAddress::ASinD:
return FuncCast<double (double)>(asin, Args_Double_Double);
return FuncCast<double (double)>(fdlibm::asin, Args_Double_Double);
case SymbolicAddress::ACosD:
return FuncCast<double (double)>(acos, Args_Double_Double);
return FuncCast<double (double)>(fdlibm::acos, Args_Double_Double);
case SymbolicAddress::ATanD:
return FuncCast<double (double)>(atan, Args_Double_Double);
return FuncCast<double (double)>(fdlibm::atan, Args_Double_Double);
case SymbolicAddress::CeilD:
return FuncCast<double (double)>(ceil, Args_Double_Double);
return FuncCast<double (double)>(fdlibm::ceil, Args_Double_Double);
case SymbolicAddress::CeilF:
return FuncCast<float (float)>(ceilf, Args_Float32_Float32);
return FuncCast<float (float)>(fdlibm::ceilf, Args_Float32_Float32);
case SymbolicAddress::FloorD:
return FuncCast<double (double)>(floor, Args_Double_Double);
return FuncCast<double (double)>(fdlibm::floor, Args_Double_Double);
case SymbolicAddress::FloorF:
return FuncCast<float (float)>(floorf, Args_Float32_Float32);
return FuncCast<float (float)>(fdlibm::floorf, Args_Float32_Float32);
case SymbolicAddress::ExpD:
return FuncCast<double (double)>(exp, Args_Double_Double);
return FuncCast<double (double)>(fdlibm::exp, Args_Double_Double);
case SymbolicAddress::LogD:
return FuncCast<double (double)>(log, Args_Double_Double);
return FuncCast<double (double)>(fdlibm::log, Args_Double_Double);
case SymbolicAddress::PowD:
return FuncCast(ecmaPow, Args_Double_DoubleDouble);
case SymbolicAddress::ATan2D:
-1
View File
@@ -1213,7 +1213,6 @@ ResetOOMFailure(JSContext* cx, unsigned argc, Value* vp)
{
CallArgs args = CallArgsFromVp(argc, vp);
args.rval().setBoolean(js::oom::HadSimulatedOOM());
HelperThreadState().waitForAllThreads();
js::oom::ResetSimulatedOOM();
return true;
}
+29 -15
View File
@@ -65,36 +65,47 @@ class NameResolver
}
/*
* Walk over the given ParseNode, converting it to a stringified name that
* respresents where the function is being assigned to.
* Walk over the given ParseNode, attempting to convert it to a stringified
* name that respresents where the function is being assigned to.
*
* |*foundName| is set to true if a name is found for the expression.
*/
bool nameExpression(ParseNode* n) {
bool nameExpression(ParseNode* n, bool* foundName) {
switch (n->getKind()) {
case PNK_DOT:
return nameExpression(n->expr()) && appendPropertyReference(n->pn_atom);
if (!nameExpression(n->expr(), foundName))
return false;
if (!*foundName)
return true;
return appendPropertyReference(n->pn_atom);
case PNK_NAME:
*foundName = true;
return buf->append(n->pn_atom);
case PNK_THIS:
*foundName = true;
return buf->append("this");
case PNK_ELEM:
return nameExpression(n->pn_left) &&
buf->append('[') &&
nameExpression(n->pn_right) &&
buf->append(']');
if (!nameExpression(n->pn_left, foundName))
return false;
if (!*foundName)
return true;
if (!buf->append('[') || !nameExpression(n->pn_right, foundName))
return false;
if (!*foundName)
return true;
return buf->append(']');
case PNK_NUMBER:
*foundName = true;
return appendNumber(n->pn_dval);
default:
/*
* Technically this isn't an "abort" situation, we're just confused
* on what to call this function, but failures in naming aren't
* treated as fatal.
*/
return false;
/* We're confused as to what to call this function. */
*foundName = false;
return true;
}
}
@@ -212,7 +223,10 @@ class NameResolver
if (assignment) {
if (assignment->isAssignment())
assignment = assignment->pn_left;
if (!nameExpression(assignment))
bool foundName = false;
if (!nameExpression(assignment, &foundName))
return false;
if (!foundName)
return true;
}
+4
View File
@@ -33,3 +33,7 @@ USE_LIBS += [
]
OS_LIBS += CONFIG['MOZ_ZLIB_LIBS']
# This is intended as a temporary workaround to enable VS2015.
if CONFIG['_MSC_VER']:
CXXFLAGS += ['-wd4312']
+8
View File
@@ -0,0 +1,8 @@
if (!('oomTest' in this))
quit();
lfLogBuffer = `this[''] = function() {}`;
loadFile(lfLogBuffer);
loadFile(lfLogBuffer);
function loadFile(lfVarx)
oomTest(function() parseModule(lfVarx))
+7 -1
View File
@@ -441,9 +441,15 @@ wasmEvalText('(module (import $bar "a" "") (func (call_import $bar)) (func $foo
// ----------------------------------------------------------------------------
// select
assertErrorMessage(() => wasmEvalText('(module (func (select (i32.const 0) (f32.const 0) (i32.const 0))))'), TypeError, mismatchError("f32", "i32"));
assertErrorMessage(() => wasmEvalText('(module (func (select (i32.const 0) (i32.const 0) (f32.const 0))))'), TypeError, mismatchError("f32", "i32"));
assertEq(wasmEvalText('(module (func (select (i32.const 0) (f32.const 0) (i32.const 0))) (export "" 0))')(), undefined);
assertEq(wasmEvalText('(module (func (select (block) (i32.const 0) (i32.const 0))) (export "" 0))')(), undefined);
assertEq(wasmEvalText('(module (func (select (return) (i32.const 0) (i32.const 0))) (export "" 0))')(), undefined);
assertEq(wasmEvalText('(module (func (i32.add (i32.const 0) (select (return) (i32.const 0) (i32.const 0)))) (export "" 0))')(), undefined);
assertEq(wasmEvalText('(module (func (select (if (i32.const 1) (i32.const 0) (f32.const 0)) (i32.const 0) (i32.const 0))) (export "" 0))')(), undefined);
assertEq(wasmEvalText('(module (func) (func (select (call 0) (call 0) (i32.const 0))) (export "" 0))')(), undefined);
(function testSideEffects() {
var numT = 0;
+4
View File
@@ -130,3 +130,7 @@ USE_LIBS += [
]
OS_LIBS += CONFIG['MOZ_ZLIB_LIBS']
# This is intended as a temporary workaround to enable VS2015.
if CONFIG['_MSC_VER']:
CXXFLAGS += ['-wd4312']
+47 -310
View File
@@ -22,6 +22,8 @@
# include <unistd.h>
#endif
#include "fdlibm.h"
#ifdef XP_WIN
# include "jswin.h"
#endif
@@ -129,28 +131,18 @@ js::math_abs(JSContext* cx, unsigned argc, Value* vp)
return math_abs_handle(cx, args[0], args.rval());
}
#if defined(SOLARIS) && defined(__GNUC__)
#define ACOS_IF_OUT_OF_RANGE(x) if (x < -1 || 1 < x) return GenericNaN();
#else
#define ACOS_IF_OUT_OF_RANGE(x)
#endif
double
js::math_acos_impl(MathCache* cache, double x)
{
ACOS_IF_OUT_OF_RANGE(x);
return cache->lookup(acos, x, MathCache::Acos);
return cache->lookup(fdlibm::acos, x, MathCache::Acos);
}
double
js::math_acos_uncached(double x)
{
ACOS_IF_OUT_OF_RANGE(x);
return acos(x);
return fdlibm::acos(x);
}
#undef ACOS_IF_OUT_OF_RANGE
bool
js::math_acos(JSContext* cx, unsigned argc, Value* vp)
{
@@ -174,28 +166,18 @@ js::math_acos(JSContext* cx, unsigned argc, Value* vp)
return true;
}
#if defined(SOLARIS) && defined(__GNUC__)
#define ASIN_IF_OUT_OF_RANGE(x) if (x < -1 || 1 < x) return GenericNaN();
#else
#define ASIN_IF_OUT_OF_RANGE(x)
#endif
double
js::math_asin_impl(MathCache* cache, double x)
{
ASIN_IF_OUT_OF_RANGE(x);
return cache->lookup(asin, x, MathCache::Asin);
return cache->lookup(fdlibm::asin, x, MathCache::Asin);
}
double
js::math_asin_uncached(double x)
{
ASIN_IF_OUT_OF_RANGE(x);
return asin(x);
return fdlibm::asin(x);
}
#undef ASIN_IF_OUT_OF_RANGE
bool
js::math_asin(JSContext* cx, unsigned argc, Value* vp)
{
@@ -222,13 +204,13 @@ js::math_asin(JSContext* cx, unsigned argc, Value* vp)
double
js::math_atan_impl(MathCache* cache, double x)
{
return cache->lookup(atan, x, MathCache::Atan);
return cache->lookup(fdlibm::atan, x, MathCache::Atan);
}
double
js::math_atan_uncached(double x)
{
return atan(x);
return fdlibm::atan(x);
}
bool
@@ -257,31 +239,7 @@ js::math_atan(JSContext* cx, unsigned argc, Value* vp)
double
js::ecmaAtan2(double y, double x)
{
#if defined(_MSC_VER)
/*
* MSVC's atan2 does not yield the result demanded by ECMA when both x
* and y are infinite.
* - The result is a multiple of pi/4.
* - The sign of y determines the sign of the result.
* - The sign of x determines the multiplicator, 1 or 3.
*/
if (IsInfinite(y) && IsInfinite(x)) {
double z = js_copysign(M_PI / 4, y);
if (x < 0)
z *= 3;
return z;
}
#endif
#if defined(SOLARIS) && defined(__GNUC__)
if (y == 0) {
if (IsNegativeZero(x))
return js_copysign(M_PI, y);
if (x == 0)
return y;
}
#endif
return atan2(y, x);
return fdlibm::atan2(y, x);
}
bool
@@ -311,11 +269,7 @@ js::math_atan2(JSContext* cx, unsigned argc, Value* vp)
double
js::math_ceil_impl(double x)
{
#ifdef __APPLE__
if (x < 0 && x > -1.0)
return js_copysign(0, -1);
#endif
return ceil(x);
return fdlibm::ceil(x);
}
bool
@@ -401,34 +355,18 @@ js::math_cos(JSContext* cx, unsigned argc, Value* vp)
return true;
}
#ifdef _WIN32
#define EXP_IF_OUT_OF_RANGE(x) \
if (!IsNaN(x)) { \
if (x == PositiveInfinity<double>()) \
return PositiveInfinity<double>(); \
if (x == NegativeInfinity<double>()) \
return 0.0; \
}
#else
#define EXP_IF_OUT_OF_RANGE(x)
#endif
double
js::math_exp_impl(MathCache* cache, double x)
{
EXP_IF_OUT_OF_RANGE(x);
return cache->lookup(exp, x, MathCache::Exp);
return cache->lookup(fdlibm::exp, x, MathCache::Exp);
}
double
js::math_exp_uncached(double x)
{
EXP_IF_OUT_OF_RANGE(x);
return exp(x);
return fdlibm::exp(x);
}
#undef EXP_IF_OUT_OF_RANGE
bool
js::math_exp(JSContext* cx, unsigned argc, Value* vp)
{
@@ -455,7 +393,7 @@ js::math_exp(JSContext* cx, unsigned argc, Value* vp)
double
js::math_floor_impl(double x)
{
return floor(x);
return fdlibm::floor(x);
}
bool
@@ -542,28 +480,18 @@ js::math_fround(JSContext* cx, unsigned argc, Value* vp)
return RoundFloat32(cx, args[0], args.rval());
}
#if defined(SOLARIS) && defined(__GNUC__)
#define LOG_IF_OUT_OF_RANGE(x) if (x < 0) return GenericNaN();
#else
#define LOG_IF_OUT_OF_RANGE(x)
#endif
double
js::math_log_impl(MathCache* cache, double x)
{
LOG_IF_OUT_OF_RANGE(x);
return cache->lookup(math_log_uncached, x, MathCache::Log);
}
double
js::math_log_uncached(double x)
{
LOG_IF_OUT_OF_RANGE(x);
return log(x);
return fdlibm::log(x);
}
#undef LOG_IF_OUT_OF_RANGE
bool
js::math_log_handle(JSContext* cx, HandleValue val, MutableHandleValue res)
{
@@ -843,7 +771,7 @@ js::math_round_impl(double x)
return x;
double add = (x >= 0) ? GetBiggestNumberLessThan(0.5) : 0.5;
return js_copysign(floor(x + add), x);
return js_copysign(fdlibm::floor(x + add), x);
}
float
@@ -858,7 +786,7 @@ js::math_roundf_impl(float x)
return x;
float add = (x >= 0) ? GetBiggestNumberLessThan(0.5f) : 0.5f;
return js_copysign(floorf(x + add), x);
return js_copysign(fdlibm::floorf(x + add), x);
}
bool /* ES5 15.8.2.15. */
@@ -876,12 +804,6 @@ js::math_round(JSContext* cx, unsigned argc, Value* vp)
double
js::math_sin_impl(MathCache* cache, double x)
{
return cache->lookup(math_sin_uncached, x, MathCache::Sin);
}
double
js::math_sin_uncached(double x)
{
#ifdef _WIN64
// Workaround MSVC bug where sin(-0) is +0 instead of -0 on x64 on
@@ -889,6 +811,12 @@ js::math_sin_uncached(double x)
if (IsNegativeZero(x))
return -0.0;
#endif
return cache->lookup(math_sin_uncached, x, MathCache::Sin);
}
double
js::math_sin_uncached(double x)
{
return sin(x);
}
@@ -1046,13 +974,13 @@ static bool math_function(JSContext* cx, unsigned argc, Value* vp)
double
js::math_log10_impl(MathCache* cache, double x)
{
return cache->lookup(log10, x, MathCache::Log10);
return cache->lookup(fdlibm::log10, x, MathCache::Log10);
}
double
js::math_log10_uncached(double x)
{
return log10(x);
return fdlibm::log10(x);
}
bool
@@ -1061,23 +989,16 @@ js::math_log10(JSContext* cx, unsigned argc, Value* vp)
return math_function<math_log10_impl>(cx, argc, vp);
}
#if !HAVE_LOG2
double log2(double x)
{
return log(x) / M_LN2;
}
#endif
double
js::math_log2_impl(MathCache* cache, double x)
{
return cache->lookup(log2, x, MathCache::Log2);
return cache->lookup(fdlibm::log2, x, MathCache::Log2);
}
double
js::math_log2_uncached(double x)
{
return log2(x);
return fdlibm::log2(x);
}
bool
@@ -1086,83 +1007,34 @@ js::math_log2(JSContext* cx, unsigned argc, Value* vp)
return math_function<math_log2_impl>(cx, argc, vp);
}
#if !HAVE_LOG1P
double log1p(double x)
{
if (fabs(x) < 1e-4) {
/*
* Use Taylor approx. log(1 + x) = x - x^2 / 2 + x^3 / 3 - x^4 / 4 with error x^5 / 5
* Since |x| < 10^-4, |x|^5 < 10^-20, relative error less than 10^-16
*/
double z = -(x * x * x * x) / 4 + (x * x * x) / 3 - (x * x) / 2 + x;
return z;
} else {
/* For other large enough values of x use direct computation */
return log(1.0 + x);
}
}
#endif
#ifdef __APPLE__
// Ensure that log1p(-0) is -0.
#define LOG1P_IF_OUT_OF_RANGE(x) if (x == 0) return x;
#else
#define LOG1P_IF_OUT_OF_RANGE(x)
#endif
double
js::math_log1p_impl(MathCache* cache, double x)
{
LOG1P_IF_OUT_OF_RANGE(x);
return cache->lookup(log1p, x, MathCache::Log1p);
return cache->lookup(fdlibm::log1p, x, MathCache::Log1p);
}
double
js::math_log1p_uncached(double x)
{
LOG1P_IF_OUT_OF_RANGE(x);
return log1p(x);
return fdlibm::log1p(x);
}
#undef LOG1P_IF_OUT_OF_RANGE
bool
js::math_log1p(JSContext* cx, unsigned argc, Value* vp)
{
return math_function<math_log1p_impl>(cx, argc, vp);
}
#if !HAVE_EXPM1
double expm1(double x)
{
/* Special handling for -0 */
if (x == 0.0)
return x;
if (fabs(x) < 1e-5) {
/*
* Use Taylor approx. exp(x) - 1 = x + x^2 / 2 + x^3 / 6 with error x^4 / 24
* Since |x| < 10^-5, |x|^4 < 10^-20, relative error less than 10^-15
*/
double z = (x * x * x) / 6 + (x * x) / 2 + x;
return z;
} else {
/* For other large enough values of x use direct computation */
return exp(x) - 1.0;
}
}
#endif
double
js::math_expm1_impl(MathCache* cache, double x)
{
return cache->lookup(expm1, x, MathCache::Expm1);
return cache->lookup(fdlibm::expm1, x, MathCache::Expm1);
}
double
js::math_expm1_uncached(double x)
{
return expm1(x);
return fdlibm::expm1(x);
}
bool
@@ -1171,27 +1043,16 @@ js::math_expm1(JSContext* cx, unsigned argc, Value* vp)
return math_function<math_expm1_impl>(cx, argc, vp);
}
#if !HAVE_SQRT1PM1
/* This algorithm computes sqrt(1+x)-1 for small x */
double sqrt1pm1(double x)
{
if (fabs(x) > 0.75)
return sqrt(1 + x) - 1;
return expm1(log1p(x) / 2);
}
#endif
double
js::math_cosh_impl(MathCache* cache, double x)
{
return cache->lookup(cosh, x, MathCache::Cosh);
return cache->lookup(fdlibm::cosh, x, MathCache::Cosh);
}
double
js::math_cosh_uncached(double x)
{
return cosh(x);
return fdlibm::cosh(x);
}
bool
@@ -1203,13 +1064,13 @@ js::math_cosh(JSContext* cx, unsigned argc, Value* vp)
double
js::math_sinh_impl(MathCache* cache, double x)
{
return cache->lookup(sinh, x, MathCache::Sinh);
return cache->lookup(fdlibm::sinh, x, MathCache::Sinh);
}
double
js::math_sinh_uncached(double x)
{
return sinh(x);
return fdlibm::sinh(x);
}
bool
@@ -1221,13 +1082,13 @@ js::math_sinh(JSContext* cx, unsigned argc, Value* vp)
double
js::math_tanh_impl(MathCache* cache, double x)
{
return cache->lookup(tanh, x, MathCache::Tanh);
return cache->lookup(fdlibm::tanh, x, MathCache::Tanh);
}
double
js::math_tanh_uncached(double x)
{
return tanh(x);
return fdlibm::tanh(x);
}
bool
@@ -1236,47 +1097,16 @@ js::math_tanh(JSContext* cx, unsigned argc, Value* vp)
return math_function<math_tanh_impl>(cx, argc, vp);
}
#if !HAVE_ACOSH
double acosh(double x)
{
const double SQUARE_ROOT_EPSILON = sqrt(std::numeric_limits<double>::epsilon());
if ((x - 1) >= SQUARE_ROOT_EPSILON) {
if (x > 1 / SQUARE_ROOT_EPSILON) {
/*
* http://functions.wolfram.com/ElementaryFunctions/ArcCosh/06/01/06/01/0001/
* approximation by laurent series in 1/x at 0+ order from -1 to 0
*/
return log(x) + M_LN2;
} else if (x < 1.5) {
// This is just a rearrangement of the standard form below
// devised to minimize loss of precision when x ~ 1:
double y = x - 1;
return log1p(y + sqrt(y * y + 2 * y));
} else {
// http://functions.wolfram.com/ElementaryFunctions/ArcCosh/02/
return log(x + sqrt(x * x - 1));
}
} else {
// see http://functions.wolfram.com/ElementaryFunctions/ArcCosh/06/01/04/01/0001/
double y = x - 1;
// approximation by taylor series in y at 0 up to order 2.
// If x is less than 1, sqrt(2 * y) is NaN and the result is NaN.
return sqrt(2 * y) * (1 - y / 12 + 3 * y * y / 160);
}
}
#endif
double
js::math_acosh_impl(MathCache* cache, double x)
{
return cache->lookup(acosh, x, MathCache::Acosh);
return cache->lookup(fdlibm::acosh, x, MathCache::Acosh);
}
double
js::math_acosh_uncached(double x)
{
return acosh(x);
return fdlibm::acosh(x);
}
bool
@@ -1285,59 +1115,16 @@ js::math_acosh(JSContext* cx, unsigned argc, Value* vp)
return math_function<math_acosh_impl>(cx, argc, vp);
}
#if !HAVE_ASINH
// Bug 899712 - gcc incorrectly rewrites -asinh(-x) to asinh(x) when overriding
// asinh.
static double my_asinh(double x)
{
const double SQUARE_ROOT_EPSILON = sqrt(std::numeric_limits<double>::epsilon());
const double FOURTH_ROOT_EPSILON = sqrt(SQUARE_ROOT_EPSILON);
if (x >= FOURTH_ROOT_EPSILON) {
if (x > 1 / SQUARE_ROOT_EPSILON)
// http://functions.wolfram.com/ElementaryFunctions/ArcSinh/06/01/06/01/0001/
// approximation by laurent series in 1/x at 0+ order from -1 to 1
return M_LN2 + log(x) + 1 / (4 * x * x);
else if (x < 0.5)
return log1p(x + sqrt1pm1(x * x));
else
return log(x + sqrt(x * x + 1));
} else if (x <= -FOURTH_ROOT_EPSILON) {
return -my_asinh(-x);
} else {
// http://functions.wolfram.com/ElementaryFunctions/ArcSinh/06/01/03/01/0001/
// approximation by taylor series in x at 0 up to order 2
double result = x;
if (fabs(x) >= SQUARE_ROOT_EPSILON) {
double x3 = x * x * x;
// approximation by taylor series in x at 0 up to order 4
result -= x3 / 6;
}
return result;
}
}
#endif
double
js::math_asinh_impl(MathCache* cache, double x)
{
#ifdef HAVE_ASINH
return cache->lookup(asinh, x, MathCache::Asinh);
#else
return cache->lookup(my_asinh, x, MathCache::Asinh);
#endif
return cache->lookup(fdlibm::asinh, x, MathCache::Asinh);
}
double
js::math_asinh_uncached(double x)
{
#ifdef HAVE_ASINH
return asinh(x);
#else
return my_asinh(x);
#endif
return fdlibm::asinh(x);
}
bool
@@ -1346,44 +1133,16 @@ js::math_asinh(JSContext* cx, unsigned argc, Value* vp)
return math_function<math_asinh_impl>(cx, argc, vp);
}
#if !HAVE_ATANH
double atanh(double x)
{
const double EPSILON = std::numeric_limits<double>::epsilon();
const double SQUARE_ROOT_EPSILON = sqrt(EPSILON);
const double FOURTH_ROOT_EPSILON = sqrt(SQUARE_ROOT_EPSILON);
if (fabs(x) >= FOURTH_ROOT_EPSILON) {
// http://functions.wolfram.com/ElementaryFunctions/ArcTanh/02/
if (fabs(x) < 0.5)
return (log1p(x) - log1p(-x)) / 2;
return log((1 + x) / (1 - x)) / 2;
} else {
// http://functions.wolfram.com/ElementaryFunctions/ArcTanh/06/01/03/01/
// approximation by taylor series in x at 0 up to order 2
double result = x;
if (fabs(x) >= SQUARE_ROOT_EPSILON) {
double x3 = x * x * x;
result += x3 / 3;
}
return result;
}
}
#endif
double
js::math_atanh_impl(MathCache* cache, double x)
{
return cache->lookup(atanh, x, MathCache::Atanh);
return cache->lookup(fdlibm::atanh, x, MathCache::Atanh);
}
double
js::math_atanh_uncached(double x)
{
return atanh(x);
return fdlibm::atanh(x);
}
bool
@@ -1396,16 +1155,7 @@ js::math_atanh(JSContext* cx, unsigned argc, Value* vp)
double
js::ecmaHypot(double x, double y)
{
#ifdef XP_WIN
/*
* Workaround MS hypot bug, where hypot(Infinity, NaN or Math.MIN_VALUE)
* is NaN, not Infinity.
*/
if (mozilla::IsInfinite(x) || mozilla::IsInfinite(y)) {
return mozilla::PositiveInfinity<double>();
}
#endif
return hypot(x, y);
return fdlibm::hypot(x, y);
}
static inline
@@ -1505,13 +1255,13 @@ js::math_hypot_handle(JSContext* cx, HandleValueArray args, MutableHandleValue r
double
js::math_trunc_impl(MathCache* cache, double x)
{
return cache->lookup(trunc, x, MathCache::Trunc);
return cache->lookup(fdlibm::trunc, x, MathCache::Trunc);
}
double
js::math_trunc_uncached(double x)
{
return trunc(x);
return fdlibm::trunc(x);
}
bool
@@ -1546,29 +1296,16 @@ js::math_sign(JSContext* cx, unsigned argc, Value* vp)
return math_function<math_sign_impl>(cx, argc, vp);
}
#if !HAVE_CBRT
double cbrt(double x)
{
if (x > 0) {
return pow(x, 1.0 / 3.0);
} else if (x == 0) {
return x;
} else {
return -pow(-x, 1.0 / 3.0);
}
}
#endif
double
js::math_cbrt_impl(MathCache* cache, double x)
{
return cache->lookup(cbrt, x, MathCache::Cbrt);
return cache->lookup(fdlibm::cbrt, x, MathCache::Cbrt);
}
double
js::math_cbrt_uncached(double x)
{
return cbrt(x);
return fdlibm::cbrt(x);
}
bool
+24
View File
@@ -17,6 +17,8 @@
#include "jstypes.h"
#include "vm/HelperThreads.h"
#ifdef WIN32
# include "jswin.h"
#endif
@@ -31,6 +33,9 @@ using mozilla::PodArrayZero;
#if defined(DEBUG) || defined(JS_OOM_BREAKPOINT)
/* For OOM testing functionality in Utility.h. */
namespace js {
mozilla::Atomic<AutoEnterOOMUnsafeRegion*> AutoEnterOOMUnsafeRegion::owner_;
namespace oom {
JS_PUBLIC_DATA(uint32_t) targetThread = 0;
@@ -54,6 +59,25 @@ GetThreadType(void) {
return threadType.get();
}
void
SimulateOOMAfter(uint64_t allocations, uint32_t thread, bool always) {
MOZ_ASSERT(counter + allocations > counter);
MOZ_ASSERT(thread > js::oom::THREAD_TYPE_NONE && thread < js::oom::THREAD_TYPE_MAX);
targetThread = thread;
maxAllocations = counter + allocations;
failAlways = always;
}
void
ResetSimulatedOOM() {
if (targetThread != THREAD_TYPE_NONE && targetThread != THREAD_TYPE_MAIN)
HelperThreadState().waitForAllThreads();
targetThread = THREAD_TYPE_NONE;
maxAllocations = UINT64_MAX;
failAlways = false;
}
} // namespace oom
} // namespace js
#endif // defined(DEBUG) || defined(JS_OOM_BREAKPOINT)
+2 -1
View File
@@ -622,6 +622,7 @@ if CONFIG['ENABLE_INTL_API']:
]
USE_LIBS += [
'fdlibm',
'nspr',
'zlib',
]
@@ -671,12 +672,12 @@ if CONFIG['_MSC_VER']:
CXXFLAGS += ['-wd4661']
CXXFLAGS += ['-we4067', '-we4258', '-we4275']
CXXFLAGS += ['-wd4146'] # FIXME: unary minus operator applied to unsigned type (bug 1229189)
CXXFLAGS += ['-wd4273'] # FIXME: inconsistent dll linkage (bug 1229666)
# This is intended as a temporary hack to support building with VS2015.
# 'noexcept' used with no exception handling mode specified;
# termination on exception is not guaranteed. Specify /EHsc
CXXFLAGS += ['-wd4577']
CXXFLAGS += ['-wd4312']
if CONFIG['OS_ARCH'] not in ('WINNT', 'HP-UX'):
OS_LIBS += [
+2
View File
@@ -6872,6 +6872,8 @@ SetRuntimeOptions(JSRuntime* rt, const OptionParser& op)
jsCacheDir = JS_smprintf("%s/%u", jsCacheDir, (unsigned)getpid());
else
jsCacheDir = JS_strdup(rt, jsCacheDir);
if (!jsCacheDir)
return false;
jsCacheAsmJSPath = JS_smprintf("%s/asmjs.cache", jsCacheDir);
}
+4
View File
@@ -44,3 +44,7 @@ shellmoduleloader.inputs = [
'../js.msg',
'ModuleLoader.js',
]
# This is intended as a temporary workaround to enable VS2015.
if CONFIG['_MSC_VER']:
CXXFLAGS += ['-wd4312']
+1 -1
View File
@@ -262,7 +262,7 @@ var cosh_data = [
[1875817529344, 28.953212876533797]
];
var sloppy_tolerance = 1000; // FIXME
var sloppy_tolerance = 8; // FIXME
for (var [x, y] of cosh_data)
assertNear(Math.acosh(x), y, sloppy_tolerance);
+8 -10
View File
@@ -282,24 +282,22 @@ var sinh_data = [
[1581915832320, 28.78280496108106]
];
var sloppy_tolerance = 1000; // FIXME
for (var [x, y] of sinh_data)
assertNear(Math.asinh(x), y, sloppy_tolerance);
assertNear(Math.asinh(x), y);
assertNear(Math.asinh(1e300), 691.4686750787737, sloppy_tolerance);
assertNear(Math.asinh(1e-300), 1e-300, sloppy_tolerance);
assertNear(Math.asinh(1e-5), 0.000009999999999833334, sloppy_tolerance);
assertNear(Math.asinh(0.3), 0.29567304756342244, sloppy_tolerance);
assertNear(Math.asinh(1), 0.881373587019543, sloppy_tolerance);
assertNear(Math.asinh(1e300), 691.4686750787737);
assertNear(Math.asinh(1e-300), 1e-300);
assertNear(Math.asinh(1e-5), 0.000009999999999833334);
assertNear(Math.asinh(0.3), 0.29567304756342244);
assertNear(Math.asinh(1), 0.881373587019543);
for (var i = 0; i <= 80; i++) {
var x = (i - 40) / 4;
assertNear(Math.asinh(Math.sinh(x)), x, sloppy_tolerance);
assertNear(Math.asinh(Math.sinh(x)), x);
}
for (var i = -20; i < 20; i++)
assertNear(Math.asinh(Math.sinh(i)), i, sloppy_tolerance);
assertNear(Math.asinh(Math.sinh(i)), i);
reportCompare(0, 0, "ok");
+1 -1
View File
@@ -266,7 +266,7 @@ var tanh_data = [
[1e-10, 1e-10],
];
var sloppy_tolerance = 10; // FIXME
var sloppy_tolerance = 2; // FIXME
for (var [x, y] of tanh_data)
assertNear(Math.atanh(x), y, sloppy_tolerance);
+3 -5
View File
@@ -1,10 +1,8 @@
assertEq(Math.cbrt(1), 1);
assertEq(Math.cbrt(-1), -1);
var sloppy_tolerance = 200; // FIXME
assertNear(Math.cbrt(1e-300), 1e-100, sloppy_tolerance);
assertNear(Math.cbrt(-1e-300), -1e-100, sloppy_tolerance);
assertNear(Math.cbrt(1e-300), 1e-100);
assertNear(Math.cbrt(-1e-300), -1e-100);
var cbrt_data = [
[ Math.E, 1.3956124250860895 ],
@@ -14,6 +12,6 @@ var cbrt_data = [
];
for (var [x, y] of cbrt_data)
assertNear(Math.cbrt(x), y, sloppy_tolerance);
assertNear(Math.cbrt(x), y);
reportCompare(0, 0, "ok");
+4 -6
View File
@@ -1,9 +1,7 @@
var sloppy_tolerance = 100;
assertEq(Math.cosh(1000), Infinity);
assertEq(Math.cosh(Number.MAX_VALUE), Infinity);
assertNear(Math.cosh(1e-30), 1, sloppy_tolerance);
assertNear(Math.cosh(1e-10), 1, sloppy_tolerance);
assertNear(Math.cosh(1e-30), 1);
assertNear(Math.cosh(1e-10), 1);
var cosh_data = [
[0.0016914556651292944, 1.0000014305114746],
@@ -270,9 +268,9 @@ var cosh_data = [
];
for (var [x, y] of cosh_data)
assertNear(Math.cosh(x), y, sloppy_tolerance);
assertNear(Math.cosh(x), y);
for (var i = -20; i < 20; i++)
assertNear(Math.cosh(i), (Math.exp(i) + Math.exp(-i)) / 2, sloppy_tolerance);
assertNear(Math.cosh(i), (Math.exp(i) + Math.exp(-i)) / 2);
reportCompare(0, 0, "ok");
+4 -5
View File
@@ -1,11 +1,10 @@
var sloppy_tolerance = 100;
for (var i = -20; i < 20; i++)
assertNear(Math.sinh(i), (Math.exp(i) - Math.exp(-i)) / 2, sloppy_tolerance);
assertNear(Math.sinh(i), (Math.exp(i) - Math.exp(-i)) / 2);
assertEq(Math.sinh(1000), Infinity);
assertEq(Math.sinh(Number.MAX_VALUE), Infinity);
assertNear(Math.sinh(1e-30), 1e-30, sloppy_tolerance);
assertNear(Math.sinh(1e-10), 1e-10, sloppy_tolerance);
assertNear(Math.sinh(1e-30), 1e-30);
assertNear(Math.sinh(1e-10), 1e-10);
var sinh_data = [
[-6.902103625349695, -497.1816406250001],
@@ -292,7 +291,7 @@ var sinh_data = [
];
for (var [x, y] of sinh_data)
assertNear(Math.sinh(x), y, sloppy_tolerance);
assertNear(Math.sinh(x), y);
reportCompare(0, 0, "ok");
+1 -1
View File
@@ -1,4 +1,4 @@
var sloppy_tolerance = 4;
var sloppy_tolerance = 2;
for (var i = -20; i < 20; i++) {
assertNear(Math.tanh(i),
+2
View File
@@ -6508,6 +6508,8 @@ class DebuggerSourceGetURLMatcher
// end to prevent them from being blacklisted by devtools by having
// the same value as a source mapped URL.
char* buf = JS_smprintf("%s > wasm", wasmModule->module().filename());
if (!buf)
return Nothing();
JSString* str = NewStringCopyZ<CanGC>(cx_, buf);
JS_smprintf_free(buf);
return Some(str);
+21 -15
View File
@@ -405,6 +405,16 @@ struct ObjectGroupCompartment::NewEntry
}
};
class ObjectGroupCompartment::NewTable : public JS::WeakCache<js::GCHashSet<NewEntry, NewEntry,
SystemAllocPolicy>>
{
using Table = js::GCHashSet<NewEntry, NewEntry, SystemAllocPolicy>;
using Base = JS::WeakCache<Table>;
public:
explicit NewTable(Zone* zone) : Base(zone, Table()) {}
};
/* static */ ObjectGroup*
ObjectGroup::defaultNewGroup(ExclusiveContext* cx, const Class* clasp,
TaggedProto proto, JSObject* associated)
@@ -422,7 +432,7 @@ ObjectGroup::defaultNewGroup(ExclusiveContext* cx, const Class* clasp,
ObjectGroupCompartment::NewTable*& table = cx->compartment()->objectGroups.defaultNewTable;
if (!table) {
table = cx->new_<ObjectGroupCompartment::NewTable>();
table = cx->new_<ObjectGroupCompartment::NewTable>(cx->zone());
if (!table || !table->init()) {
js_delete(table);
table = nullptr;
@@ -538,7 +548,7 @@ ObjectGroup::lazySingletonGroup(ExclusiveContext* cx, const Class* clasp, Tagged
ObjectGroupCompartment::NewTable*& table = cx->compartment()->objectGroups.lazyTable;
if (!table) {
table = cx->new_<ObjectGroupCompartment::NewTable>();
table = cx->new_<ObjectGroupCompartment::NewTable>(cx->zone());
if (!table || !table->init()) {
ReportOutOfMemory(cx);
js_delete(table);
@@ -580,8 +590,8 @@ ObjectGroup::setDefaultNewGroupUnknown(JSContext* cx, const Class* clasp, Handle
ObjectGroupCompartment::NewTable* table = cx->compartment()->objectGroups.defaultNewTable;
if (table) {
Rooted<TaggedProto> taggedProto(cx, TaggedProto(obj));
ObjectGroupCompartment::NewTable::Ptr p =
table->lookup(ObjectGroupCompartment::NewEntry::Lookup(clasp, taggedProto, nullptr));
auto lookup = ObjectGroupCompartment::NewEntry::Lookup(clasp, taggedProto, nullptr);
auto p = table->lookup(lookup);
if (p)
MarkObjectGroupUnknownProperties(cx, p->group);
}
@@ -594,8 +604,8 @@ ObjectGroup::hasDefaultNewGroup(JSObject* proto, const Class* clasp, ObjectGroup
ObjectGroupCompartment::NewTable* table = proto->compartment()->objectGroups.defaultNewTable;
if (table) {
ObjectGroupCompartment::NewTable::Ptr p =
table->lookup(ObjectGroupCompartment::NewEntry::Lookup(clasp, TaggedProto(proto), nullptr));
auto lookup = ObjectGroupCompartment::NewEntry::Lookup(clasp, TaggedProto(proto), nullptr);
auto p = table->lookup(lookup);
return p && p->group == group;
}
return false;
@@ -1602,10 +1612,10 @@ void
ObjectGroupCompartment::removeDefaultNewGroup(const Class* clasp, TaggedProto proto,
JSObject* associated)
{
NewTable::Ptr p = defaultNewTable->lookup(NewEntry::Lookup(clasp, proto, associated));
auto p = defaultNewTable->lookup(NewEntry::Lookup(clasp, proto, associated));
MOZ_RELEASE_ASSERT(p);
defaultNewTable->remove(p);
defaultNewTable->get().remove(p);
}
void
@@ -1614,9 +1624,9 @@ ObjectGroupCompartment::replaceDefaultNewGroup(const Class* clasp, TaggedProto p
{
NewEntry::Lookup lookup(clasp, proto, associated);
NewTable::Ptr p = defaultNewTable->lookup(lookup);
auto p = defaultNewTable->lookup(lookup);
MOZ_RELEASE_ASSERT(p);
defaultNewTable->remove(p);
defaultNewTable->get().remove(p);
{
AutoEnterOOMUnsafeRegion oomUnsafe;
if (!defaultNewTable->putNew(lookup, NewEntry(group, associated)))
@@ -1722,10 +1732,6 @@ ObjectGroupCompartment::sweep(FreeOp* fop)
plainObjectTable->sweep();
if (allocationSiteTable)
allocationSiteTable->sweep();
if (defaultNewTable)
defaultNewTable->sweep();
if (lazyTable)
lazyTable->sweep();
}
void
@@ -1782,7 +1788,7 @@ ObjectGroupCompartment::checkNewTableAfterMovingGC(NewTable* table)
clasp = nullptr;
NewEntry::Lookup lookup(clasp, proto, entry.associated);
NewTable::Ptr ptr = table->lookup(lookup);
auto ptr = table->lookup(lookup);
MOZ_RELEASE_ASSERT(ptr.found() && &*ptr == &e.front());
}
}
+1 -1
View File
@@ -540,7 +540,7 @@ class ObjectGroupCompartment
friend class ObjectGroup;
struct NewEntry;
using NewTable = js::GCHashSet<NewEntry, NewEntry, SystemAllocPolicy>;
class NewTable;
// Set of default 'new' or lazy groups in the compartment.
NewTable* defaultNewTable;
+6 -12
View File
@@ -87,7 +87,6 @@ bool
SimpleCount::count(CountBase& countBase, mozilla::MallocSizeOf mallocSizeOf, const Node& node)
{
Count& count = static_cast<Count&>(countBase);
count.total_++;
if (reportBytes)
count.totalBytes_ += node.size(mallocSizeOf);
return true;
@@ -264,7 +263,6 @@ bool
ByCoarseType::count(CountBase& countBase, mozilla::MallocSizeOf mallocSizeOf, const Node& node)
{
Count& count = static_cast<Count&>(countBase);
count.total_++;
switch (node.coarseType()) {
case JS::ubi::CoarseType::Object:
@@ -315,17 +313,17 @@ ByCoarseType::report(JSContext* cx, CountBase& countBase, MutableHandleValue rep
}
// Comparison function for sorting hash table entries by total count. The
// Comparison function for sorting hash table entries by the smallest node ID
// they counted. Node IDs are stable and unique, which ensures ordering of
// results never depends on hash table placement or sort algorithm vagaries. The
// arguments are doubly indirect: they're pointers to elements in an array of
// pointers to table entries.
template<typename Entry>
static int compareEntries(const void* lhsVoid, const void* rhsVoid) {
size_t lhs = (*static_cast<const Entry* const*>(lhsVoid))->value()->total_;
size_t rhs = (*static_cast<const Entry* const*>(rhsVoid))->value()->total_;
auto lhs = (*static_cast<const Entry* const*>(lhsVoid))->value()->smallestNodeIdCounted_;
auto rhs = (*static_cast<const Entry* const*>(rhsVoid))->value()->smallestNodeIdCounted_;
// qsort sorts in "ascending" order, so we should describe entries with
// smaller counts as being "greater than" entries with larger counts. We
// don't want to just subtract the counts, as they're unsigned.
// We don't want to just subtract the values, as they're unsigned.
if (lhs < rhs)
return 1;
if (lhs > rhs)
@@ -459,7 +457,6 @@ bool
ByObjectClass::count(CountBase& countBase, mozilla::MallocSizeOf mallocSizeOf, const Node& node)
{
Count& count = static_cast<Count&>(countBase);
count.total_++;
const char* className = node.jsObjectClassName();
if (!className)
@@ -553,7 +550,6 @@ bool
ByUbinodeType::count(CountBase& countBase, mozilla::MallocSizeOf mallocSizeOf, const Node& node)
{
Count& count = static_cast<Count&>(countBase);
count.total_++;
const char16_t* key = node.typeName();
MOZ_ASSERT(key);
@@ -708,7 +704,6 @@ bool
ByAllocationStack::count(CountBase& countBase, mozilla::MallocSizeOf mallocSizeOf, const Node& node)
{
Count& count = static_cast<Count&>(countBase);
count.total_++;
// If we do have an allocation stack for this node, include it in the
// count for that stack.
@@ -877,7 +872,6 @@ bool
ByFilename::count(CountBase& countBase, mozilla::MallocSizeOf mallocSizeOf, const Node& node)
{
Count& count = static_cast<Count&>(countBase);
count.total_++;
const char* filename = node.scriptFilename();
if (!filename)
+29 -5
View File
@@ -947,13 +947,14 @@ xpc::GlobalProperties::Parse(JSContext* cx, JS::HandleObject obj)
bool
xpc::GlobalProperties::Define(JSContext* cx, JS::HandleObject obj)
{
// Properties will be exposed to System automatically but not to Sandboxes
// if |[Exposed=System]| is specified.
// This function holds common properties not exposed automatically but able
// to be requested either in |Cu.importGlobalProperties| or
// |wantGlobalProperties| of a sandbox.
if (CSS && !dom::CSSBinding::GetConstructorObject(cx, obj))
return false;
if (indexedDB &&
!IndexedDatabaseManager::DefineIndexedDB(cx, obj))
return false;
if (XMLHttpRequest &&
!dom::XMLHttpRequestBinding::GetConstructorObject(cx, obj))
return false;
@@ -1010,6 +1011,29 @@ xpc::GlobalProperties::Define(JSContext* cx, JS::HandleObject obj)
return true;
}
bool
xpc::GlobalProperties::DefineInXPCComponents(JSContext* cx, JS::HandleObject obj)
{
if (indexedDB &&
!IndexedDatabaseManager::DefineIndexedDB(cx, obj))
return false;
return Define(cx, obj);
}
bool
xpc::GlobalProperties::DefineInSandbox(JSContext* cx, JS::HandleObject obj)
{
MOZ_ASSERT(IsSandbox(obj));
if (indexedDB &&
!(IndexedDatabaseManager::ResolveSandboxBinding(cx, obj) &&
IndexedDatabaseManager::DefineIndexedDB(cx, obj)))
return false;
return Define(cx, obj);
}
nsresult
xpc::CreateSandboxObject(JSContext* cx, MutableHandleValue vp, nsISupports* prinOrSop,
SandboxOptions& options)
@@ -1176,7 +1200,7 @@ xpc::CreateSandboxObject(JSContext* cx, MutableHandleValue vp, nsISupports* prin
!JS_DefineFunction(cx, sandbox, "isProxy", SandboxIsProxy, 1, 0)))
return NS_ERROR_XPC_UNEXPECTED;
if (!options.globalProperties.Define(cx, sandbox))
if (!options.globalProperties.DefineInSandbox(cx, sandbox))
return NS_ERROR_XPC_UNEXPECTED;
#ifndef SPIDERMONKEY_PROMISE
+1 -1
View File
@@ -2548,7 +2548,7 @@ nsXPCComponents_Utils::ImportGlobalProperties(HandleValue aPropertyList,
}
if (!options.Parse(cx, propertyList) ||
!options.Define(cx, global))
!options.DefineInXPCComponents(cx, global))
{
return NS_ERROR_FAILURE;
}
+4 -1
View File
@@ -3371,7 +3371,8 @@ struct GlobalProperties {
}
bool Parse(JSContext* cx, JS::HandleObject obj);
bool Define(JSContext* cx, JS::HandleObject obj);
bool DefineInXPCComponents(JSContext* cx, JS::HandleObject obj);
bool DefineInSandbox(JSContext* cx, JS::HandleObject obj);
bool CSS : 1;
bool indexedDB : 1;
bool XMLHttpRequest : 1;
@@ -3388,6 +3389,8 @@ struct GlobalProperties {
bool fetch : 1;
bool caches : 1;
bool fileReader: 1;
private:
bool Define(JSContext* cx, JS::HandleObject obj);
};
// Infallible.
+2
View File
@@ -93,6 +93,7 @@ if CONFIG['MOZ_WEBRTC_SIGNALING']:
'signaling/src/sdp/sipcc/sdp_utils.c',
]
GYP_DIRS['signaling'].sandbox_vars['ALLOW_COMPILER_WARNINGS'] = True
GYP_DIRS['signaling'].non_unified_sources += signaling_non_unified_sources
if CONFIG['MOZ_WIDGET_TOOLKIT'] != 'gonk':
@@ -112,6 +113,7 @@ if CONFIG['MOZ_WIDGET_TOOLKIT'] != 'gonk':
moz_webrtc_mediacodec=0,
build_for_standalone=0
)
GYP_DIRS['signalingtest'].sandbox_vars['ALLOW_COMPILER_WARNINGS'] = True
GYP_DIRS['signalingtest'].non_unified_sources += signaling_non_unified_sources
if CONFIG['MOZ_X11']:
+6
View File
@@ -24,3 +24,9 @@ USE_LIBS += [
'mtransport_s',
]
if CONFIG['_MSC_VER']:
# This is intended as a temporary workaround to enable warning free building
# with VS2015.
# reinterpret_cast': conversion from 'DWORD' to 'HANDLE' of greater size
CXXFLAGS += ['-wd4312']
+17
View File
@@ -0,0 +1,17 @@
This is the fdlibm library imported from
https://github.com/freebsd/freebsd
Upstream code can be viewed at
https://github.com/freebsd/freebsd/tree/master/lib/msun/src
Each file is downloaded separately, as cloning whole repository takes so much
resources.
The in-tree copy is updated by running
sh update.sh
from within the modules/fdlibm directory.
Current version: [commit bcea9d50b15e4f0027a5dd526e0e2a612238471e].
patches 01-14 fixes files to be usable within mozilla-central tree.
See https://bugzilla.mozilla.org/show_bug.cgi?id=933257
+105
View File
@@ -0,0 +1,105 @@
#!/bin/sh
set -e
BASE_URL=https://raw.githubusercontent.com/freebsd/freebsd/master/lib/msun/src
function download_source() {
REMOTE_FILENAME=$1
LOCAL_FILENAME=$2
curl -o "src/${LOCAL_FILENAME}" "${BASE_URL}/${REMOTE_FILENAME}"
}
mkdir -p src
# headers
download_source math.h fdlibm.h
download_source math_private.h math_private.h
# Math.acos
download_source e_acos.c e_acos.cpp
# Math.acosh
download_source e_acosh.c e_acosh.cpp
# Math.asin
download_source e_asin.c e_asin.cpp
# Math.asinh
download_source s_asinh.c s_asinh.cpp
# Math.atan
download_source s_atan.c s_atan.cpp
# Math.atanh
download_source e_atanh.c e_atanh.cpp
# Math.atan2
download_source e_atan2.c e_atan2.cpp
# Math.cbrt
download_source s_cbrt.c s_cbrt.cpp
# Math.ceil
download_source s_ceil.c s_ceil.cpp
download_source s_ceilf.c s_ceilf.cpp
# Math.cos (not used due to poor performance)
# Math.cosh
download_source e_cosh.c e_cosh.cpp
# Math.exp
download_source e_exp.c e_exp.cpp
# Math.expm1
download_source s_expm1.c s_expm1.cpp
# Math.floor and Math.round
download_source s_floor.c s_floor.cpp
# Math.fround
download_source s_floorf.c s_floorf.cpp
# Math.hypot
download_source e_hypot.c e_hypot.cpp
# Math.log
download_source e_log.c e_log.cpp
# Math.log1p
download_source s_log1p.c s_log1p.cpp
# Math.log10
download_source e_log10.c e_log10.cpp
download_source k_log.h k_log.h
# Math.log2
download_source e_log2.c e_log2.cpp
# Math.pow (not used due to poor performance)
# Math.sin (not used due to poor performance)
# Math.sinh
download_source e_sinh.c e_sinh.cpp
# Math.sqrt (not used due to poor performance)
# Math.tan (not used due to poor performance)
# Math.tanh
download_source s_tanh.c s_tanh.cpp
# Math.trunc
download_source s_trunc.c s_trunc.cpp
# dependencies
download_source k_exp.c k_exp.cpp
download_source s_copysign.c s_copysign.cpp
download_source s_fabs.c s_fabs.cpp
download_source s_scalbn.c s_scalbn.cpp
# These are not not used in Math.* functions, but used internally.
download_source e_pow.c e_pow.cpp
download_source e_sqrt.c e_sqrt.cpp
+7
View File
@@ -0,0 +1,7 @@
# -*- Mode: python; c-basic-offset: 4; indent-tabs-mode: nil; tab-width: 40 -*-
# vim: set filetype=python:
# This Source Code Form is subject to the terms of the Mozilla Public
# License, v. 2.0. If a copy of the MPL was not distributed with this
# file, You can obtain one at http://mozilla.org/MPL/2.0/.
DIRS += ['src']
@@ -0,0 +1,505 @@
diff --git a/modules/fdlibm/src/fdlibm.h b/modules/fdlibm/src/fdlibm.h
--- a/modules/fdlibm/src/fdlibm.h
+++ b/modules/fdlibm/src/fdlibm.h
@@ -12,496 +12,45 @@
/*
* from: @(#)fdlibm.h 5.1 93/09/24
* $FreeBSD$
*/
#ifndef _MATH_H_
#define _MATH_H_
-#include <sys/cdefs.h>
-#include <sys/_types.h>
-#include <machine/_limits.h>
-
-/*
- * ANSI/POSIX
- */
-extern const union __infinity_un {
- unsigned char __uc[8];
- double __ud;
-} __infinity;
-
-extern const union __nan_un {
- unsigned char __uc[sizeof(float)];
- float __uf;
-} __nan;
-
-#if __GNUC_PREREQ__(3, 3) || (defined(__INTEL_COMPILER) && __INTEL_COMPILER >= 800)
-#define __MATH_BUILTIN_CONSTANTS
-#endif
-
-#if __GNUC_PREREQ__(3, 0) && !defined(__INTEL_COMPILER)
-#define __MATH_BUILTIN_RELOPS
-#endif
-
-#ifdef __MATH_BUILTIN_CONSTANTS
-#define HUGE_VAL __builtin_huge_val()
-#else
-#define HUGE_VAL (__infinity.__ud)
-#endif
-
-#if __ISO_C_VISIBLE >= 1999
-#define FP_ILOGB0 (-__INT_MAX)
-#define FP_ILOGBNAN __INT_MAX
-
-#ifdef __MATH_BUILTIN_CONSTANTS
-#define HUGE_VALF __builtin_huge_valf()
-#define HUGE_VALL __builtin_huge_vall()
-#define INFINITY __builtin_inff()
-#define NAN __builtin_nanf("")
-#else
-#define HUGE_VALF (float)HUGE_VAL
-#define HUGE_VALL (long double)HUGE_VAL
-#define INFINITY HUGE_VALF
-#define NAN (__nan.__uf)
-#endif /* __MATH_BUILTIN_CONSTANTS */
-
-#define MATH_ERRNO 1
-#define MATH_ERREXCEPT 2
-#define math_errhandling MATH_ERREXCEPT
-
-#define FP_FAST_FMAF 1
-
-/* Symbolic constants to classify floating point numbers. */
-#define FP_INFINITE 0x01
-#define FP_NAN 0x02
-#define FP_NORMAL 0x04
-#define FP_SUBNORMAL 0x08
-#define FP_ZERO 0x10
-
-#if (__STDC_VERSION__ >= 201112L && defined(__clang__)) || \
- __has_extension(c_generic_selections)
-#define __fp_type_select(x, f, d, ld) _Generic((x), \
- float: f(x), \
- double: d(x), \
- long double: ld(x), \
- volatile float: f(x), \
- volatile double: d(x), \
- volatile long double: ld(x), \
- volatile const float: f(x), \
- volatile const double: d(x), \
- volatile const long double: ld(x), \
- const float: f(x), \
- const double: d(x), \
- const long double: ld(x))
-#elif __GNUC_PREREQ__(3, 1) && !defined(__cplusplus)
-#define __fp_type_select(x, f, d, ld) __builtin_choose_expr( \
- __builtin_types_compatible_p(__typeof(x), long double), ld(x), \
- __builtin_choose_expr( \
- __builtin_types_compatible_p(__typeof(x), double), d(x), \
- __builtin_choose_expr( \
- __builtin_types_compatible_p(__typeof(x), float), f(x), (void)0)))
-#else
-#define __fp_type_select(x, f, d, ld) \
- ((sizeof(x) == sizeof(float)) ? f(x) \
- : (sizeof(x) == sizeof(double)) ? d(x) \
- : ld(x))
-#endif
-
-#define fpclassify(x) \
- __fp_type_select(x, __fpclassifyf, __fpclassifyd, __fpclassifyl)
-#define isfinite(x) __fp_type_select(x, __isfinitef, __isfinite, __isfinitel)
-#define isinf(x) __fp_type_select(x, __isinff, __isinf, __isinfl)
-#define isnan(x) \
- __fp_type_select(x, __inline_isnanf, __inline_isnan, __inline_isnanl)
-#define isnormal(x) __fp_type_select(x, __isnormalf, __isnormal, __isnormall)
-
-#ifdef __MATH_BUILTIN_RELOPS
-#define isgreater(x, y) __builtin_isgreater((x), (y))
-#define isgreaterequal(x, y) __builtin_isgreaterequal((x), (y))
-#define isless(x, y) __builtin_isless((x), (y))
-#define islessequal(x, y) __builtin_islessequal((x), (y))
-#define islessgreater(x, y) __builtin_islessgreater((x), (y))
-#define isunordered(x, y) __builtin_isunordered((x), (y))
-#else
-#define isgreater(x, y) (!isunordered((x), (y)) && (x) > (y))
-#define isgreaterequal(x, y) (!isunordered((x), (y)) && (x) >= (y))
-#define isless(x, y) (!isunordered((x), (y)) && (x) < (y))
-#define islessequal(x, y) (!isunordered((x), (y)) && (x) <= (y))
-#define islessgreater(x, y) (!isunordered((x), (y)) && \
- ((x) > (y) || (y) > (x)))
-#define isunordered(x, y) (isnan(x) || isnan(y))
-#endif /* __MATH_BUILTIN_RELOPS */
-
-#define signbit(x) __fp_type_select(x, __signbitf, __signbit, __signbitl)
-
-typedef __double_t double_t;
-typedef __float_t float_t;
-#endif /* __ISO_C_VISIBLE >= 1999 */
-
-/*
- * XOPEN/SVID
- */
-#if __BSD_VISIBLE || __XSI_VISIBLE
-#define M_E 2.7182818284590452354 /* e */
-#define M_LOG2E 1.4426950408889634074 /* log 2e */
-#define M_LOG10E 0.43429448190325182765 /* log 10e */
-#define M_LN2 0.69314718055994530942 /* log e2 */
-#define M_LN10 2.30258509299404568402 /* log e10 */
-#define M_PI 3.14159265358979323846 /* pi */
-#define M_PI_2 1.57079632679489661923 /* pi/2 */
-#define M_PI_4 0.78539816339744830962 /* pi/4 */
-#define M_1_PI 0.31830988618379067154 /* 1/pi */
-#define M_2_PI 0.63661977236758134308 /* 2/pi */
-#define M_2_SQRTPI 1.12837916709551257390 /* 2/sqrt(pi) */
-#define M_SQRT2 1.41421356237309504880 /* sqrt(2) */
-#define M_SQRT1_2 0.70710678118654752440 /* 1/sqrt(2) */
-
-#define MAXFLOAT ((float)3.40282346638528860e+38)
-extern int signgam;
-#endif /* __BSD_VISIBLE || __XSI_VISIBLE */
-
-#if __BSD_VISIBLE
-#if 0
-/* Old value from 4.4BSD-Lite math.h; this is probably better. */
-#define HUGE HUGE_VAL
-#else
-#define HUGE MAXFLOAT
-#endif
-#endif /* __BSD_VISIBLE */
-
-/*
- * Most of these functions depend on the rounding mode and have the side
- * effect of raising floating-point exceptions, so they are not declared
- * as __pure2. In C99, FENV_ACCESS affects the purity of these functions.
- */
-__BEGIN_DECLS
-/*
- * ANSI/POSIX
- */
-int __fpclassifyd(double) __pure2;
-int __fpclassifyf(float) __pure2;
-int __fpclassifyl(long double) __pure2;
-int __isfinitef(float) __pure2;
-int __isfinite(double) __pure2;
-int __isfinitel(long double) __pure2;
-int __isinff(float) __pure2;
-int __isinf(double) __pure2;
-int __isinfl(long double) __pure2;
-int __isnormalf(float) __pure2;
-int __isnormal(double) __pure2;
-int __isnormall(long double) __pure2;
-int __signbit(double) __pure2;
-int __signbitf(float) __pure2;
-int __signbitl(long double) __pure2;
-
-static __inline int
-__inline_isnan(__const double __x)
-{
-
- return (__x != __x);
-}
-
-static __inline int
-__inline_isnanf(__const float __x)
-{
-
- return (__x != __x);
-}
-
-static __inline int
-__inline_isnanl(__const long double __x)
-{
-
- return (__x != __x);
-}
-
-/*
- * Version 2 of the Single UNIX Specification (UNIX98) defined isnan() and
- * isinf() as functions taking double. C99, and the subsequent POSIX revisions
- * (SUSv3, POSIX.1-2001, define it as a macro that accepts any real floating
- * point type. If we are targeting SUSv2 and C99 or C11 (or C++11) then we
- * expose the newer definition, assuming that the language spec takes
- * precedence over the operating system interface spec.
- */
-#if __XSI_VISIBLE > 0 && __XSI_VISIBLE < 600 && __ISO_C_VISIBLE < 1999
-#undef isinf
-#undef isnan
-int isinf(double);
-int isnan(double);
-#endif
-
double acos(double);
double asin(double);
double atan(double);
double atan2(double, double);
-double cos(double);
-double sin(double);
-double tan(double);
double cosh(double);
double sinh(double);
double tanh(double);
double exp(double);
-double frexp(double, int *); /* fundamentally !__pure2 */
-double ldexp(double, int);
double log(double);
double log10(double);
-double modf(double, double *); /* fundamentally !__pure2 */
double pow(double, double);
double sqrt(double);
double ceil(double);
-double fabs(double) __pure2;
+float ceilf(float);
+double fabs(double);
double floor(double);
-double fmod(double, double);
-/*
- * These functions are not in C90.
- */
-#if __BSD_VISIBLE || __ISO_C_VISIBLE >= 1999 || __XSI_VISIBLE
double acosh(double);
double asinh(double);
double atanh(double);
double cbrt(double);
-double erf(double);
-double erfc(double);
-double exp2(double);
double expm1(double);
-double fma(double, double, double);
double hypot(double, double);
-int ilogb(double) __pure2;
-double lgamma(double);
-long long llrint(double);
-long long llround(double);
double log1p(double);
double log2(double);
-double logb(double);
-long lrint(double);
-long lround(double);
-double nan(const char *) __pure2;
-double nextafter(double, double);
-double remainder(double, double);
-double remquo(double, double, int *);
-double rint(double);
-#endif /* __BSD_VISIBLE || __ISO_C_VISIBLE >= 1999 || __XSI_VISIBLE */
-#if __BSD_VISIBLE || __XSI_VISIBLE
-double j0(double);
-double j1(double);
-double jn(int, double);
-double y0(double);
-double y1(double);
-double yn(int, double);
+double copysign(double, double);
+double scalbn(double, int);
+double trunc(double);
-#if __XSI_VISIBLE <= 500 || __BSD_VISIBLE
-double gamma(double);
-#endif
-
-#if __XSI_VISIBLE <= 600 || __BSD_VISIBLE
-double scalb(double, double);
-#endif
-#endif /* __BSD_VISIBLE || __XSI_VISIBLE */
-
-#if __BSD_VISIBLE || __ISO_C_VISIBLE >= 1999
-double copysign(double, double) __pure2;
-double fdim(double, double);
-double fmax(double, double) __pure2;
-double fmin(double, double) __pure2;
-double nearbyint(double);
-double round(double);
-double scalbln(double, long);
-double scalbn(double, int);
-double tgamma(double);
-double trunc(double);
-#endif
-
-/*
- * BSD math library entry points
- */
-#if __BSD_VISIBLE
-double drem(double, double);
-int finite(double) __pure2;
-int isnanf(float) __pure2;
-
-/*
- * Reentrant version of gamma & lgamma; passes signgam back by reference
- * as the second argument; user must allocate space for signgam.
- */
-double gamma_r(double, int *);
-double lgamma_r(double, int *);
-
-/*
- * IEEE Test Vector
- */
-double significand(double);
-#endif /* __BSD_VISIBLE */
-
-/* float versions of ANSI/POSIX functions */
-#if __ISO_C_VISIBLE >= 1999
-float acosf(float);
-float asinf(float);
-float atanf(float);
-float atan2f(float, float);
-float cosf(float);
-float sinf(float);
-float tanf(float);
-
-float coshf(float);
-float sinhf(float);
-float tanhf(float);
-
-float exp2f(float);
-float expf(float);
-float expm1f(float);
-float frexpf(float, int *); /* fundamentally !__pure2 */
-int ilogbf(float) __pure2;
-float ldexpf(float, int);
-float log10f(float);
-float log1pf(float);
-float log2f(float);
-float logf(float);
-float modff(float, float *); /* fundamentally !__pure2 */
-
-float powf(float, float);
-float sqrtf(float);
-
-float ceilf(float);
-float fabsf(float) __pure2;
float floorf(float);
-float fmodf(float, float);
-float roundf(float);
-
-float erff(float);
-float erfcf(float);
-float hypotf(float, float);
-float lgammaf(float);
-float tgammaf(float);
-
-float acoshf(float);
-float asinhf(float);
-float atanhf(float);
-float cbrtf(float);
-float logbf(float);
-float copysignf(float, float) __pure2;
-long long llrintf(float);
-long long llroundf(float);
-long lrintf(float);
-long lroundf(float);
-float nanf(const char *) __pure2;
-float nearbyintf(float);
-float nextafterf(float, float);
-float remainderf(float, float);
-float remquof(float, float, int *);
-float rintf(float);
-float scalblnf(float, long);
-float scalbnf(float, int);
-float truncf(float);
-
-float fdimf(float, float);
-float fmaf(float, float, float);
-float fmaxf(float, float) __pure2;
-float fminf(float, float) __pure2;
-#endif
-
-/*
- * float versions of BSD math library entry points
- */
-#if __BSD_VISIBLE
-float dremf(float, float);
-int finitef(float) __pure2;
-float gammaf(float);
-float j0f(float);
-float j1f(float);
-float jnf(int, float);
-float scalbf(float, float);
-float y0f(float);
-float y1f(float);
-float ynf(int, float);
-
-/*
- * Float versions of reentrant version of gamma & lgamma; passes
- * signgam back by reference as the second argument; user must
- * allocate space for signgam.
- */
-float gammaf_r(float, int *);
-float lgammaf_r(float, int *);
-
-/*
- * float version of IEEE Test Vector
- */
-float significandf(float);
-#endif /* __BSD_VISIBLE */
-
-/*
- * long double versions of ISO/POSIX math functions
- */
-#if __ISO_C_VISIBLE >= 1999
-long double acoshl(long double);
-long double acosl(long double);
-long double asinhl(long double);
-long double asinl(long double);
-long double atan2l(long double, long double);
-long double atanhl(long double);
-long double atanl(long double);
-long double cbrtl(long double);
-long double ceill(long double);
-long double copysignl(long double, long double) __pure2;
-long double coshl(long double);
-long double cosl(long double);
-long double erfcl(long double);
-long double erfl(long double);
-long double exp2l(long double);
-long double expl(long double);
-long double expm1l(long double);
-long double fabsl(long double) __pure2;
-long double fdiml(long double, long double);
-long double floorl(long double);
-long double fmal(long double, long double, long double);
-long double fmaxl(long double, long double) __pure2;
-long double fminl(long double, long double) __pure2;
-long double fmodl(long double, long double);
-long double frexpl(long double value, int *); /* fundamentally !__pure2 */
-long double hypotl(long double, long double);
-int ilogbl(long double) __pure2;
-long double ldexpl(long double, int);
-long double lgammal(long double);
-long long llrintl(long double);
-long long llroundl(long double);
-long double log10l(long double);
-long double log1pl(long double);
-long double log2l(long double);
-long double logbl(long double);
-long double logl(long double);
-long lrintl(long double);
-long lroundl(long double);
-long double modfl(long double, long double *); /* fundamentally !__pure2 */
-long double nanl(const char *) __pure2;
-long double nearbyintl(long double);
-long double nextafterl(long double, long double);
-double nexttoward(double, long double);
-float nexttowardf(float, long double);
-long double nexttowardl(long double, long double);
-long double powl(long double, long double);
-long double remainderl(long double, long double);
-long double remquol(long double, long double, int *);
-long double rintl(long double);
-long double roundl(long double);
-long double scalblnl(long double, long);
-long double scalbnl(long double, int);
-long double sinhl(long double);
-long double sinl(long double);
-long double sqrtl(long double);
-long double tanhl(long double);
-long double tanl(long double);
-long double tgammal(long double);
-long double truncl(long double);
-#endif /* __ISO_C_VISIBLE >= 1999 */
-
-#if __BSD_VISIBLE
-long double lgammal_r(long double, int *);
-#endif
-
-__END_DECLS
#endif /* !_MATH_H_ */
@@ -0,0 +1,35 @@
diff --git a/modules/fdlibm/src/fdlibm.h b/modules/fdlibm/src/fdlibm.h
--- a/modules/fdlibm/src/fdlibm.h
+++ b/modules/fdlibm/src/fdlibm.h
@@ -9,18 +9,18 @@
* ====================================================
*/
/*
* from: @(#)fdlibm.h 5.1 93/09/24
* $FreeBSD$
*/
-#ifndef _MATH_H_
-#define _MATH_H_
+#ifndef mozilla_imported_fdlibm_h
+#define mozilla_imported_fdlibm_h
double acos(double);
double asin(double);
double atan(double);
double atan2(double, double);
double cosh(double);
double sinh(double);
@@ -45,9 +45,9 @@ double log1p(double);
double log2(double);
double copysign(double, double);
double scalbn(double, int);
double trunc(double);
float floorf(float);
-#endif /* !_MATH_H_ */
+#endif /* mozilla_imported_fdlibm_h */
@@ -0,0 +1,34 @@
diff --git a/modules/fdlibm/src/fdlibm.h b/modules/fdlibm/src/fdlibm.h
--- a/modules/fdlibm/src/fdlibm.h
+++ b/modules/fdlibm/src/fdlibm.h
@@ -12,16 +12,18 @@
/*
* from: @(#)fdlibm.h 5.1 93/09/24
* $FreeBSD$
*/
#ifndef mozilla_imported_fdlibm_h
#define mozilla_imported_fdlibm_h
+namespace fdlibm {
+
double acos(double);
double asin(double);
double atan(double);
double atan2(double, double);
double cosh(double);
double sinh(double);
double tanh(double);
@@ -45,9 +47,11 @@ double log1p(double);
double log2(double);
double copysign(double, double);
double scalbn(double, int);
double trunc(double);
float floorf(float);
+} /* namespace fdlibm */
+
#endif /* mozilla_imported_fdlibm_h */
@@ -0,0 +1,631 @@
diff --git a/modules/fdlibm/src/e_acos.cpp b/modules/fdlibm/src/e_acos.cpp
--- a/modules/fdlibm/src/e_acos.cpp
+++ b/modules/fdlibm/src/e_acos.cpp
@@ -35,17 +35,16 @@
* if x is NaN, return x itself;
* if |x|>1, return NaN with invalid signal.
*
* Function needed: sqrt
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double
one= 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
pi = 3.14159265358979311600e+00, /* 0x400921FB, 0x54442D18 */
pio2_hi = 1.57079632679489655800e+00; /* 0x3FF921FB, 0x54442D18 */
static volatile double
pio2_lo = 6.12323399573676603587e-17; /* 0x3C91A626, 0x33145C07 */
diff --git a/modules/fdlibm/src/e_acosh.cpp b/modules/fdlibm/src/e_acosh.cpp
--- a/modules/fdlibm/src/e_acosh.cpp
+++ b/modules/fdlibm/src/e_acosh.cpp
@@ -26,17 +26,16 @@
*
* Special cases:
* acosh(x) is NaN with signal if x<1.
* acosh(NaN) is NaN without signal.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double
one = 1.0,
ln2 = 6.93147180559945286227e-01; /* 0x3FE62E42, 0xFEFA39EF */
double
__ieee754_acosh(double x)
diff --git a/modules/fdlibm/src/e_asin.cpp b/modules/fdlibm/src/e_asin.cpp
--- a/modules/fdlibm/src/e_asin.cpp
+++ b/modules/fdlibm/src/e_asin.cpp
@@ -41,17 +41,16 @@
* Special cases:
* if x is NaN, return x itself;
* if |x|>1, return NaN with invalid signal.
*
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double
one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
huge = 1.000e+300,
pio2_hi = 1.57079632679489655800e+00, /* 0x3FF921FB, 0x54442D18 */
pio2_lo = 6.12323399573676603587e-17, /* 0x3C91A626, 0x33145C07 */
pio4_hi = 7.85398163397448278999e-01, /* 0x3FE921FB, 0x54442D18 */
diff --git a/modules/fdlibm/src/e_atan2.cpp b/modules/fdlibm/src/e_atan2.cpp
--- a/modules/fdlibm/src/e_atan2.cpp
+++ b/modules/fdlibm/src/e_atan2.cpp
@@ -39,17 +39,16 @@
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static volatile double
tiny = 1.0e-300;
static const double
zero = 0.0,
pi_o_4 = 7.8539816339744827900E-01, /* 0x3FE921FB, 0x54442D18 */
pi_o_2 = 1.5707963267948965580E+00, /* 0x3FF921FB, 0x54442D18 */
diff --git a/modules/fdlibm/src/e_atanh.cpp b/modules/fdlibm/src/e_atanh.cpp
--- a/modules/fdlibm/src/e_atanh.cpp
+++ b/modules/fdlibm/src/e_atanh.cpp
@@ -30,17 +30,16 @@
* atanh(x) is NaN if |x| > 1 with signal;
* atanh(NaN) is that NaN with no signal;
* atanh(+-1) is +-INF with signal.
*
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double one = 1.0, huge = 1e300;
static const double zero = 0.0;
double
__ieee754_atanh(double x)
{
diff --git a/modules/fdlibm/src/e_cosh.cpp b/modules/fdlibm/src/e_cosh.cpp
--- a/modules/fdlibm/src/e_cosh.cpp
+++ b/modules/fdlibm/src/e_cosh.cpp
@@ -32,17 +32,16 @@
*
* Special cases:
* cosh(x) is |x| if x is +INF, -INF, or NaN.
* only cosh(0)=1 is exact for finite x.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double one = 1.0, half=0.5, huge = 1.0e300;
double
__ieee754_cosh(double x)
{
double t,w;
diff --git a/modules/fdlibm/src/e_exp.cpp b/modules/fdlibm/src/e_exp.cpp
--- a/modules/fdlibm/src/e_exp.cpp
+++ b/modules/fdlibm/src/e_exp.cpp
@@ -73,17 +73,16 @@
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double
one = 1.0,
halF[2] = {0.5,-0.5,},
o_threshold= 7.09782712893383973096e+02, /* 0x40862E42, 0xFEFA39EF */
u_threshold= -7.45133219101941108420e+02, /* 0xc0874910, 0xD52D3051 */
ln2HI[2] ={ 6.93147180369123816490e-01, /* 0x3fe62e42, 0xfee00000 */
diff --git a/modules/fdlibm/src/e_hypot.cpp b/modules/fdlibm/src/e_hypot.cpp
--- a/modules/fdlibm/src/e_hypot.cpp
+++ b/modules/fdlibm/src/e_hypot.cpp
@@ -43,17 +43,16 @@
*
* Accuracy:
* hypot(x,y) returns sqrt(x^2+y^2) with error less
* than 1 ulps (units in the last place)
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
double
__ieee754_hypot(double x, double y)
{
double a,b,t1,t2,y1,y2,w;
int32_t j,k,ha,hb;
diff --git a/modules/fdlibm/src/e_log.cpp b/modules/fdlibm/src/e_log.cpp
--- a/modules/fdlibm/src/e_log.cpp
+++ b/modules/fdlibm/src/e_log.cpp
@@ -62,17 +62,16 @@
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double
ln2_hi = 6.93147180369123816490e-01, /* 3fe62e42 fee00000 */
ln2_lo = 1.90821492927058770002e-10, /* 3dea39ef 35793c76 */
two54 = 1.80143985094819840000e+16, /* 43500000 00000000 */
Lg1 = 6.666666666666735130e-01, /* 3FE55555 55555593 */
Lg2 = 3.999999999940941908e-01, /* 3FD99999 9997FA04 */
diff --git a/modules/fdlibm/src/e_log10.cpp b/modules/fdlibm/src/e_log10.cpp
--- a/modules/fdlibm/src/e_log10.cpp
+++ b/modules/fdlibm/src/e_log10.cpp
@@ -19,17 +19,16 @@
* comments.
*
* log10(x) = (f - 0.5*f*f + k_log1p(f)) / ln10 + k * log10(2)
* in not-quite-routine extra precision.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
#include "k_log.h"
static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
ivln10hi = 4.34294481878168880939e-01, /* 0x3fdbcb7b, 0x15200000 */
ivln10lo = 2.50829467116452752298e-11, /* 0x3dbb9438, 0xca9aadd5 */
log10_2hi = 3.01029995663611771306e-01, /* 0x3FD34413, 0x509F6000 */
diff --git a/modules/fdlibm/src/e_log2.cpp b/modules/fdlibm/src/e_log2.cpp
--- a/modules/fdlibm/src/e_log2.cpp
+++ b/modules/fdlibm/src/e_log2.cpp
@@ -21,17 +21,16 @@
* This reduces x to {k, 1+f} exactly as in e_log.c, then calls the kernel,
* then does the combining and scaling steps
* log2(x) = (f - 0.5*f*f + k_log1p(f)) / ln2 + k
* in not-quite-routine extra precision.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
#include "k_log.h"
static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
ivln2hi = 1.44269504072144627571e+00, /* 0x3ff71547, 0x65200000 */
ivln2lo = 1.67517131648865118353e-10; /* 0x3de705fc, 0x2eefa200 */
diff --git a/modules/fdlibm/src/e_pow.cpp b/modules/fdlibm/src/e_pow.cpp
--- a/modules/fdlibm/src/e_pow.cpp
+++ b/modules/fdlibm/src/e_pow.cpp
@@ -52,17 +52,16 @@
*
* Constants :
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
-#include "math.h"
#include "math_private.h"
static const double
bp[] = {1.0, 1.5,},
dp_h[] = { 0.0, 5.84962487220764160156e-01,}, /* 0x3FE2B803, 0x40000000 */
dp_l[] = { 0.0, 1.35003920212974897128e-08,}, /* 0x3E4CFDEB, 0x43CFD006 */
zero = 0.0,
one = 1.0,
diff --git a/modules/fdlibm/src/e_sinh.cpp b/modules/fdlibm/src/e_sinh.cpp
--- a/modules/fdlibm/src/e_sinh.cpp
+++ b/modules/fdlibm/src/e_sinh.cpp
@@ -29,17 +29,16 @@
*
* Special cases:
* sinh(x) is |x| if x is +INF, -INF, or NaN.
* only sinh(0)=0 is exact for finite x.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double one = 1.0, shuge = 1.0e307;
double
__ieee754_sinh(double x)
{
double t,h;
diff --git a/modules/fdlibm/src/e_sqrt.cpp b/modules/fdlibm/src/e_sqrt.cpp
--- a/modules/fdlibm/src/e_sqrt.cpp
+++ b/modules/fdlibm/src/e_sqrt.cpp
@@ -81,17 +81,16 @@
* sqrt(NaN) = NaN ... with invalid signal for signaling NaN
*
* Other methods : see the appended file at the end of the program below.
*---------------
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double one = 1.0, tiny=1.0e-300;
double
__ieee754_sqrt(double x)
{
double z;
diff --git a/modules/fdlibm/src/k_exp.cpp b/modules/fdlibm/src/k_exp.cpp
--- a/modules/fdlibm/src/k_exp.cpp
+++ b/modules/fdlibm/src/k_exp.cpp
@@ -24,17 +24,16 @@
* SUCH DAMAGE.
*/
#include <sys/cdefs.h>
__FBSDID("$FreeBSD$");
#include <complex.h>
-#include "math.h"
#include "math_private.h"
static const uint32_t k = 1799; /* constant for reduction */
static const double kln2 = 1246.97177782734161156; /* k * ln2 */
/*
* Compute exp(x), scaled to avoid spurious overflow. An exponent is
* returned separately in 'expt'.
diff --git a/modules/fdlibm/src/math_private.h b/modules/fdlibm/src/math_private.h
--- a/modules/fdlibm/src/math_private.h
+++ b/modules/fdlibm/src/math_private.h
@@ -15,16 +15,18 @@
*/
#ifndef _MATH_PRIVATE_H_
#define _MATH_PRIVATE_H_
#include <sys/types.h>
#include <machine/endian.h>
+#include "fdlibm.h"
+
/*
* The original fdlibm code used statements like:
* n0 = ((*(int*)&one)>>29)^1; * index of high word *
* ix0 = *(n0+(int*)&x); * high word of x *
* ix1 = *((1-n0)+(int*)&x); * low word of x *
* to dig two 32 bit words out of the 64 bit IEEE floating point
* value. That is non-ANSI, and, moreover, the gcc instruction
* scheduler gets it wrong. We instead use the following macros.
diff --git a/modules/fdlibm/src/s_asinh.cpp b/modules/fdlibm/src/s_asinh.cpp
--- a/modules/fdlibm/src/s_asinh.cpp
+++ b/modules/fdlibm/src/s_asinh.cpp
@@ -21,17 +21,16 @@
* asinh(x) := x if 1+x*x=1,
* := sign(x)*(log(x)+ln2)) for large |x|, else
* := sign(x)*log(2|x|+1/(|x|+sqrt(x*x+1))) if|x|>2, else
* := sign(x)*log1p(|x| + x^2/(1 + sqrt(1+x^2)))
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double
one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
ln2 = 6.93147180559945286227e-01, /* 0x3FE62E42, 0xFEFA39EF */
huge= 1.00000000000000000000e+300;
double
diff --git a/modules/fdlibm/src/s_atan.cpp b/modules/fdlibm/src/s_atan.cpp
--- a/modules/fdlibm/src/s_atan.cpp
+++ b/modules/fdlibm/src/s_atan.cpp
@@ -30,17 +30,16 @@
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double atanhi[] = {
4.63647609000806093515e-01, /* atan(0.5)hi 0x3FDDAC67, 0x0561BB4F */
7.85398163397448278999e-01, /* atan(1.0)hi 0x3FE921FB, 0x54442D18 */
9.82793723247329054082e-01, /* atan(1.5)hi 0x3FEF730B, 0xD281F69B */
1.57079632679489655800e+00, /* atan(inf)hi 0x3FF921FB, 0x54442D18 */
};
diff --git a/modules/fdlibm/src/s_cbrt.cpp b/modules/fdlibm/src/s_cbrt.cpp
--- a/modules/fdlibm/src/s_cbrt.cpp
+++ b/modules/fdlibm/src/s_cbrt.cpp
@@ -10,17 +10,16 @@
* ====================================================
*
* Optimized by Bruce D. Evans.
*/
#include <sys/cdefs.h>
__FBSDID("$FreeBSD$");
-#include "math.h"
#include "math_private.h"
/* cbrt(x)
* Return cube root of x
*/
static const u_int32_t
B1 = 715094163, /* B1 = (1023-1023/3-0.03306235651)*2**20 */
B2 = 696219795; /* B2 = (1023-1023/3-54/3-0.03306235651)*2**20 */
diff --git a/modules/fdlibm/src/s_ceil.cpp b/modules/fdlibm/src/s_ceil.cpp
--- a/modules/fdlibm/src/s_ceil.cpp
+++ b/modules/fdlibm/src/s_ceil.cpp
@@ -19,17 +19,16 @@
* Method:
* Bit twiddling.
* Exception:
* Inexact flag raised if x not equal to ceil(x).
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double huge = 1.0e300;
double
ceil(double x)
{
int32_t i0,i1,j0;
diff --git a/modules/fdlibm/src/s_ceilf.cpp b/modules/fdlibm/src/s_ceilf.cpp
--- a/modules/fdlibm/src/s_ceilf.cpp
+++ b/modules/fdlibm/src/s_ceilf.cpp
@@ -11,17 +11,16 @@
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
#include <sys/cdefs.h>
__FBSDID("$FreeBSD$");
-#include "math.h"
#include "math_private.h"
static const float huge = 1.0e30;
float
ceilf(float x)
{
int32_t i0,j0;
diff --git a/modules/fdlibm/src/s_copysign.cpp b/modules/fdlibm/src/s_copysign.cpp
--- a/modules/fdlibm/src/s_copysign.cpp
+++ b/modules/fdlibm/src/s_copysign.cpp
@@ -14,17 +14,16 @@
__FBSDID("$FreeBSD$");
/*
* copysign(double x, double y)
* copysign(x,y) returns a value with the magnitude of x and
* with the sign bit of y.
*/
-#include "math.h"
#include "math_private.h"
double
copysign(double x, double y)
{
u_int32_t hx,hy;
GET_HIGH_WORD(hx,x);
GET_HIGH_WORD(hy,y);
diff --git a/modules/fdlibm/src/s_expm1.cpp b/modules/fdlibm/src/s_expm1.cpp
--- a/modules/fdlibm/src/s_expm1.cpp
+++ b/modules/fdlibm/src/s_expm1.cpp
@@ -105,17 +105,16 @@
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double
one = 1.0,
tiny = 1.0e-300,
o_threshold = 7.09782712893383973096e+02,/* 0x40862E42, 0xFEFA39EF */
ln2_hi = 6.93147180369123816490e-01,/* 0x3fe62e42, 0xfee00000 */
ln2_lo = 1.90821492927058770002e-10,/* 0x3dea39ef, 0x35793c76 */
diff --git a/modules/fdlibm/src/s_fabs.cpp b/modules/fdlibm/src/s_fabs.cpp
--- a/modules/fdlibm/src/s_fabs.cpp
+++ b/modules/fdlibm/src/s_fabs.cpp
@@ -13,17 +13,16 @@
#ifndef lint
static char rcsid[] = "$FreeBSD$";
#endif
/*
* fabs(x) returns the absolute value of x.
*/
-#include "math.h"
#include "math_private.h"
double
fabs(double x)
{
u_int32_t high;
GET_HIGH_WORD(high,x);
SET_HIGH_WORD(x,high&0x7fffffff);
diff --git a/modules/fdlibm/src/s_floor.cpp b/modules/fdlibm/src/s_floor.cpp
--- a/modules/fdlibm/src/s_floor.cpp
+++ b/modules/fdlibm/src/s_floor.cpp
@@ -19,17 +19,16 @@
* Method:
* Bit twiddling.
* Exception:
* Inexact flag raised if x not equal to floor(x).
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double huge = 1.0e300;
double
floor(double x)
{
int32_t i0,i1,j0;
diff --git a/modules/fdlibm/src/s_floorf.cpp b/modules/fdlibm/src/s_floorf.cpp
--- a/modules/fdlibm/src/s_floorf.cpp
+++ b/modules/fdlibm/src/s_floorf.cpp
@@ -20,17 +20,16 @@
* floorf(x)
* Return x rounded toward -inf to integral value
* Method:
* Bit twiddling.
* Exception:
* Inexact flag raised if x not equal to floorf(x).
*/
-#include "math.h"
#include "math_private.h"
static const float huge = 1.0e30;
float
floorf(float x)
{
int32_t i0,j0;
diff --git a/modules/fdlibm/src/s_log1p.cpp b/modules/fdlibm/src/s_log1p.cpp
--- a/modules/fdlibm/src/s_log1p.cpp
+++ b/modules/fdlibm/src/s_log1p.cpp
@@ -75,17 +75,16 @@
* if(u==1.0) return x ; else
* return log(u)*(x/(u-1.0));
*
* See HP-15C Advanced Functions Handbook, p.193.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double
ln2_hi = 6.93147180369123816490e-01, /* 3fe62e42 fee00000 */
ln2_lo = 1.90821492927058770002e-10, /* 3dea39ef 35793c76 */
two54 = 1.80143985094819840000e+16, /* 43500000 00000000 */
Lp1 = 6.666666666666735130e-01, /* 3FE55555 55555593 */
Lp2 = 3.999999999940941908e-01, /* 3FD99999 9997FA04 */
diff --git a/modules/fdlibm/src/s_scalbn.cpp b/modules/fdlibm/src/s_scalbn.cpp
--- a/modules/fdlibm/src/s_scalbn.cpp
+++ b/modules/fdlibm/src/s_scalbn.cpp
@@ -19,17 +19,16 @@ static char rcsid[] = "$FreeBSD$";
* scalbn(x,n) returns x* 2**n computed by exponent
* manipulation rather than by actually performing an
* exponentiation or a multiplication.
*/
#include <sys/cdefs.h>
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
huge = 1.0e+300,
tiny = 1.0e-300;
diff --git a/modules/fdlibm/src/s_tanh.cpp b/modules/fdlibm/src/s_tanh.cpp
--- a/modules/fdlibm/src/s_tanh.cpp
+++ b/modules/fdlibm/src/s_tanh.cpp
@@ -34,17 +34,16 @@
*
* Special cases:
* tanh(NaN) is NaN;
* only tanh(0)=0 is exact for finite argument.
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const volatile double tiny = 1.0e-300;
static const double one = 1.0, two = 2.0, huge = 1.0e300;
double
tanh(double x)
{
diff --git a/modules/fdlibm/src/s_trunc.cpp b/modules/fdlibm/src/s_trunc.cpp
--- a/modules/fdlibm/src/s_trunc.cpp
+++ b/modules/fdlibm/src/s_trunc.cpp
@@ -19,17 +19,16 @@
* Method:
* Bit twiddling.
* Exception:
* Inexact flag raised if x not equal to trunc(x).
*/
#include <float.h>
-#include "math.h"
#include "math_private.h"
static const double huge = 1.0e300;
double
trunc(double x)
{
int32_t i0,i1,j0;
@@ -0,0 +1,21 @@
diff --git a/modules/fdlibm/src/math_private.h b/modules/fdlibm/src/math_private.h
--- a/modules/fdlibm/src/math_private.h
+++ b/modules/fdlibm/src/math_private.h
@@ -12,16 +12,17 @@
/*
* from: @(#)fdlibm.h 5.1 93/09/24
* $FreeBSD$
*/
#ifndef _MATH_PRIVATE_H_
#define _MATH_PRIVATE_H_
+#include <stdint.h>
#include <sys/types.h>
#include <machine/endian.h>
#include "fdlibm.h"
/*
* The original fdlibm code used statements like:
* n0 = ((*(int*)&one)>>29)^1; * index of high word *
@@ -0,0 +1,74 @@
diff --git a/modules/fdlibm/src/math_private.h b/modules/fdlibm/src/math_private.h
--- a/modules/fdlibm/src/math_private.h
+++ b/modules/fdlibm/src/math_private.h
@@ -14,20 +14,21 @@
* $FreeBSD$
*/
#ifndef _MATH_PRIVATE_H_
#define _MATH_PRIVATE_H_
#include <stdint.h>
#include <sys/types.h>
-#include <machine/endian.h>
#include "fdlibm.h"
+#include "mozilla/Endian.h"
+
/*
* The original fdlibm code used statements like:
* n0 = ((*(int*)&one)>>29)^1; * index of high word *
* ix0 = *(n0+(int*)&x); * high word of x *
* ix1 = *((1-n0)+(int*)&x); * low word of x *
* to dig two 32 bit words out of the 64 bit IEEE floating point
* value. That is non-ANSI, and, moreover, the gcc instruction
* scheduler gets it wrong. We instead use the following macros.
@@ -36,27 +37,17 @@
* endianness at run time.
*/
/*
* A union which permits us to convert between a double and two 32 bit
* ints.
*/
-#ifdef __arm__
-#if defined(__VFP_FP__) || defined(__ARM_EABI__)
-#define IEEE_WORD_ORDER BYTE_ORDER
-#else
-#define IEEE_WORD_ORDER BIG_ENDIAN
-#endif
-#else /* __arm__ */
-#define IEEE_WORD_ORDER BYTE_ORDER
-#endif
-
-#if IEEE_WORD_ORDER == BIG_ENDIAN
+#if MOZ_BIG_ENDIAN
typedef union
{
double value;
struct
{
u_int32_t msw;
u_int32_t lsw;
@@ -64,17 +55,17 @@ typedef union
struct
{
u_int64_t w;
} xparts;
} ieee_double_shape_type;
#endif
-#if IEEE_WORD_ORDER == LITTLE_ENDIAN
+#if MOZ_LITTLE_ENDIAN
typedef union
{
double value;
struct
{
u_int32_t lsw;
u_int32_t msw;
@@ -0,0 +1,50 @@
diff --git a/modules/fdlibm/src/math_private.h b/modules/fdlibm/src/math_private.h
--- a/modules/fdlibm/src/math_private.h
+++ b/modules/fdlibm/src/math_private.h
@@ -724,16 +724,46 @@ irintl(long double x)
#define __ieee754_j1f j1f
#define __ieee754_y0f y0f
#define __ieee754_y1f y1f
#define __ieee754_jnf jnf
#define __ieee754_ynf ynf
#define __ieee754_remainderf remainderf
#define __ieee754_scalbf scalbf
+#define acos fdlibm::acos
+#define asin fdlibm::asin
+#define atan fdlibm::atan
+#define atan2 fdlibm::atan2
+#define cosh fdlibm::cosh
+#define sinh fdlibm::sinh
+#define tanh fdlibm::tanh
+#define exp fdlibm::exp
+#define log fdlibm::log
+#define log10 fdlibm::log10
+#define pow fdlibm::pow
+#define sqrt fdlibm::sqrt
+#define ceil fdlibm::ceil
+#define ceilf fdlibm::ceilf
+#define fabs fdlibm::fabs
+#define floor fdlibm::floor
+#define acosh fdlibm::acosh
+#define asinh fdlibm::asinh
+#define atanh fdlibm::atanh
+#define cbrt fdlibm::cbrt
+#define expm1 fdlibm::expm1
+#define hypot fdlibm::hypot
+#define log1p fdlibm::log1p
+#define log2 fdlibm::log2
+#define scalb fdlibm::scalb
+#define copysign fdlibm::copysign
+#define scalbn fdlibm::scalbn
+#define trunc fdlibm::trunc
+#define floorf fdlibm::floorf
+
/* fdlibm kernel function */
int __kernel_rem_pio2(double*,double*,int,int,int);
/* double precision kernel functions */
#ifndef INLINE_REM_PIO2
int __ieee754_rem_pio2(double,double*);
#endif
double __kernel_sin(double,double,int);
@@ -0,0 +1,377 @@
diff --git a/modules/fdlibm/src/e_acos.cpp b/modules/fdlibm/src/e_acos.cpp
--- a/modules/fdlibm/src/e_acos.cpp
+++ b/modules/fdlibm/src/e_acos.cpp
@@ -99,12 +99,8 @@ double
c = (z-df*df)/(s+df);
p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
r = p/q;
w = r*s+c;
return 2.0*(df+w);
}
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(acos, acosl);
-#endif
diff --git a/modules/fdlibm/src/e_acosh.cpp b/modules/fdlibm/src/e_acosh.cpp
--- a/modules/fdlibm/src/e_acosh.cpp
+++ b/modules/fdlibm/src/e_acosh.cpp
@@ -56,12 +56,8 @@ double
} else if (hx > 0x40000000) { /* 2**28 > x > 2 */
t=x*x;
return __ieee754_log(2.0*x-one/(x+sqrt(t-one)));
} else { /* 1<x<2 */
t = x-one;
return log1p(t+sqrt(2.0*t+t*t));
}
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(acosh, acoshl);
-#endif
diff --git a/modules/fdlibm/src/e_asin.cpp b/modules/fdlibm/src/e_asin.cpp
--- a/modules/fdlibm/src/e_asin.cpp
+++ b/modules/fdlibm/src/e_asin.cpp
@@ -105,12 +105,8 @@ double
c = (t-w*w)/(s+w);
r = p/q;
p = 2.0*s*r-(pio2_lo-2.0*c);
q = pio4_hi-2.0*w;
t = pio4_hi-(p-q);
}
if(hx>0) return t; else return -t;
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(asin, asinl);
-#endif
diff --git a/modules/fdlibm/src/e_atan2.cpp b/modules/fdlibm/src/e_atan2.cpp
--- a/modules/fdlibm/src/e_atan2.cpp
+++ b/modules/fdlibm/src/e_atan2.cpp
@@ -117,12 +117,8 @@ double
switch (m) {
case 0: return z ; /* atan(+,+) */
case 1: return -z ; /* atan(-,+) */
case 2: return pi-(z-pi_lo);/* atan(+,-) */
default: /* case 3 */
return (z-pi_lo)-pi;/* atan(-,-) */
}
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(atan2, atan2l);
-#endif
diff --git a/modules/fdlibm/src/e_atanh.cpp b/modules/fdlibm/src/e_atanh.cpp
--- a/modules/fdlibm/src/e_atanh.cpp
+++ b/modules/fdlibm/src/e_atanh.cpp
@@ -56,12 +56,8 @@ double
SET_HIGH_WORD(x,ix);
if(ix<0x3fe00000) { /* x < 0.5 */
t = x+x;
t = 0.5*log1p(t+t*x/(one-x));
} else
t = 0.5*log1p((x+x)/(one-x));
if(hx>=0) return t; else return -t;
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(atanh, atanhl);
-#endif
diff --git a/modules/fdlibm/src/e_cosh.cpp b/modules/fdlibm/src/e_cosh.cpp
--- a/modules/fdlibm/src/e_cosh.cpp
+++ b/modules/fdlibm/src/e_cosh.cpp
@@ -73,12 +73,8 @@ double
/* |x| in [log(maxdouble), overflowthresold] */
if (ix<=0x408633CE)
return __ldexp_exp(fabs(x), -1);
/* |x| > overflowthresold, cosh(x) overflow */
return huge*huge;
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(cosh, coshl);
-#endif
diff --git a/modules/fdlibm/src/e_exp.cpp b/modules/fdlibm/src/e_exp.cpp
--- a/modules/fdlibm/src/e_exp.cpp
+++ b/modules/fdlibm/src/e_exp.cpp
@@ -152,12 +152,8 @@ double
else y = one-((lo-(x*c)/(2.0-c))-hi);
if(k >= -1021) {
if (k==1024) return y*2.0*0x1p1023;
return y*twopk;
} else {
return y*twopk*twom1000;
}
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(exp, expl);
-#endif
diff --git a/modules/fdlibm/src/e_hypot.cpp b/modules/fdlibm/src/e_hypot.cpp
--- a/modules/fdlibm/src/e_hypot.cpp
+++ b/modules/fdlibm/src/e_hypot.cpp
@@ -119,12 +119,8 @@ double
if(k!=0) {
u_int32_t high;
t1 = 1.0;
GET_HIGH_WORD(high,t1);
SET_HIGH_WORD(t1,high+(k<<20));
return t1*w;
} else return w;
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(hypot, hypotl);
-#endif
diff --git a/modules/fdlibm/src/e_log.cpp b/modules/fdlibm/src/e_log.cpp
--- a/modules/fdlibm/src/e_log.cpp
+++ b/modules/fdlibm/src/e_log.cpp
@@ -135,12 +135,8 @@ double
hfsq=0.5*f*f;
if(k==0) return f-(hfsq-s*(hfsq+R)); else
return dk*ln2_hi-((hfsq-(s*(hfsq+R)+dk*ln2_lo))-f);
} else {
if(k==0) return f-s*(f-R); else
return dk*ln2_hi-((s*(f-R)-dk*ln2_lo)-f);
}
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(log, logl);
-#endif
diff --git a/modules/fdlibm/src/e_log10.cpp b/modules/fdlibm/src/e_log10.cpp
--- a/modules/fdlibm/src/e_log10.cpp
+++ b/modules/fdlibm/src/e_log10.cpp
@@ -82,12 +82,8 @@ double
* with some parallelism and it reduces the error for many args.
*/
w = y2 + val_hi;
val_lo += (y2 - w) + val_hi;
val_hi = w;
return val_lo + val_hi;
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(log10, log10l);
-#endif
diff --git a/modules/fdlibm/src/e_log2.cpp b/modules/fdlibm/src/e_log2.cpp
--- a/modules/fdlibm/src/e_log2.cpp
+++ b/modules/fdlibm/src/e_log2.cpp
@@ -105,12 +105,8 @@ double
/* spadd(val_hi, val_lo, y), except for not using double_t: */
w = y + val_hi;
val_lo += (y - w) + val_hi;
val_hi = w;
return val_lo + val_hi;
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(log2, log2l);
-#endif
diff --git a/modules/fdlibm/src/e_sinh.cpp b/modules/fdlibm/src/e_sinh.cpp
--- a/modules/fdlibm/src/e_sinh.cpp
+++ b/modules/fdlibm/src/e_sinh.cpp
@@ -67,12 +67,8 @@ double
/* |x| in [log(maxdouble), overflowthresold] */
if (ix<=0x408633CE)
return h*2.0*__ldexp_exp(fabs(x), -1);
/* |x| > overflowthresold, sinh(x) overflow */
return x*shuge;
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(sinh, sinhl);
-#endif
diff --git a/modules/fdlibm/src/e_sqrt.cpp b/modules/fdlibm/src/e_sqrt.cpp
--- a/modules/fdlibm/src/e_sqrt.cpp
+++ b/modules/fdlibm/src/e_sqrt.cpp
@@ -182,20 +182,16 @@ double
ix0 = (q>>1)+0x3fe00000;
ix1 = q1>>1;
if ((q&1)==1) ix1 |= sign;
ix0 += (m <<20);
INSERT_WORDS(z,ix0,ix1);
return z;
}
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(sqrt, sqrtl);
-#endif
-
/*
Other methods (use floating-point arithmetic)
-------------
(This is a copy of a drafted paper by Prof W. Kahan
and K.C. Ng, written in May, 1986)
Two algorithms are given here to implement sqrt(x)
(IEEE double precision arithmetic) in software.
diff --git a/modules/fdlibm/src/s_asinh.cpp b/modules/fdlibm/src/s_asinh.cpp
--- a/modules/fdlibm/src/s_asinh.cpp
+++ b/modules/fdlibm/src/s_asinh.cpp
@@ -50,12 +50,8 @@ asinh(double x)
t = fabs(x);
w = __ieee754_log(2.0*t+one/(__ieee754_sqrt(x*x+one)+t));
} else { /* 2.0 > |x| > 2**-28 */
t = x*x;
w =log1p(fabs(x)+t/(one+__ieee754_sqrt(one+t)));
}
if(hx>0) return w; else return -w;
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(asinh, asinhl);
-#endif
diff --git a/modules/fdlibm/src/s_atan.cpp b/modules/fdlibm/src/s_atan.cpp
--- a/modules/fdlibm/src/s_atan.cpp
+++ b/modules/fdlibm/src/s_atan.cpp
@@ -112,12 +112,8 @@ atan(double x)
s1 = z*(aT[0]+w*(aT[2]+w*(aT[4]+w*(aT[6]+w*(aT[8]+w*aT[10])))));
s2 = w*(aT[1]+w*(aT[3]+w*(aT[5]+w*(aT[7]+w*aT[9]))));
if (id<0) return x - x*(s1+s2);
else {
z = atanhi[id] - ((x*(s1+s2) - atanlo[id]) - x);
return (hx<0)? -z:z;
}
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(atan, atanl);
-#endif
diff --git a/modules/fdlibm/src/s_cbrt.cpp b/modules/fdlibm/src/s_cbrt.cpp
--- a/modules/fdlibm/src/s_cbrt.cpp
+++ b/modules/fdlibm/src/s_cbrt.cpp
@@ -105,12 +105,8 @@ cbrt(double x)
s=t*t; /* t*t is exact */
r=x/s; /* error <= 0.5 ulps; |r| < |t| */
w=t+t; /* t+t is exact */
r=(r-t)/(w+r); /* r-t is exact; w+r ~= 3*t */
t=t+t*r; /* error <= 0.5 + 0.5/3 + epsilon */
return(t);
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(cbrt, cbrtl);
-#endif
diff --git a/modules/fdlibm/src/s_ceil.cpp b/modules/fdlibm/src/s_ceil.cpp
--- a/modules/fdlibm/src/s_ceil.cpp
+++ b/modules/fdlibm/src/s_ceil.cpp
@@ -65,12 +65,8 @@ ceil(double x)
}
}
i1 &= (~i);
}
}
INSERT_WORDS(x,i0,i1);
return x;
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(ceil, ceill);
-#endif
diff --git a/modules/fdlibm/src/s_expm1.cpp b/modules/fdlibm/src/s_expm1.cpp
--- a/modules/fdlibm/src/s_expm1.cpp
+++ b/modules/fdlibm/src/s_expm1.cpp
@@ -210,12 +210,8 @@ expm1(double x)
SET_HIGH_WORD(t,((0x3ff-k)<<20)); /* 2^-k */
y = x-(e+t);
y += one;
y = y*twopk;
}
}
return y;
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(expm1, expm1l);
-#endif
diff --git a/modules/fdlibm/src/s_floor.cpp b/modules/fdlibm/src/s_floor.cpp
--- a/modules/fdlibm/src/s_floor.cpp
+++ b/modules/fdlibm/src/s_floor.cpp
@@ -66,12 +66,8 @@ floor(double x)
}
}
i1 &= (~i);
}
}
INSERT_WORDS(x,i0,i1);
return x;
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(floor, floorl);
-#endif
diff --git a/modules/fdlibm/src/s_log1p.cpp b/modules/fdlibm/src/s_log1p.cpp
--- a/modules/fdlibm/src/s_log1p.cpp
+++ b/modules/fdlibm/src/s_log1p.cpp
@@ -168,12 +168,8 @@ log1p(double x)
return k*ln2_hi-((R-(k*ln2_lo+c))-f);
}
s = f/(2.0+f);
z = s*s;
R = z*(Lp1+z*(Lp2+z*(Lp3+z*(Lp4+z*(Lp5+z*(Lp6+z*Lp7))))));
if(k==0) return f-(hfsq-s*(hfsq+R)); else
return k*ln2_hi-((hfsq-(s*(hfsq+R)+(k*ln2_lo+c)))-f);
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(log1p, log1pl);
-#endif
diff --git a/modules/fdlibm/src/s_scalbn.cpp b/modules/fdlibm/src/s_scalbn.cpp
--- a/modules/fdlibm/src/s_scalbn.cpp
+++ b/modules/fdlibm/src/s_scalbn.cpp
@@ -53,13 +53,8 @@ scalbn (double x, int n)
if (k <= -54)
if (n > 50000) /* in case integer overflow in n+k */
return huge*copysign(huge,x); /*overflow*/
else return tiny*copysign(tiny,x); /*underflow*/
k += 54; /* subnormal result */
SET_HIGH_WORD(x,(hx&0x800fffff)|(k<<20));
return x*twom54;
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(scalbn, ldexpl);
-__weak_reference(scalbn, scalbnl);
-#endif
diff --git a/modules/fdlibm/src/s_tanh.cpp b/modules/fdlibm/src/s_tanh.cpp
--- a/modules/fdlibm/src/s_tanh.cpp
+++ b/modules/fdlibm/src/s_tanh.cpp
@@ -72,12 +72,8 @@ tanh(double x)
z= -t/(t+two);
}
/* |x| >= 22, return +-1 */
} else {
z = one - tiny; /* raise inexact flag */
}
return (jx>=0)? z: -z;
}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(tanh, tanhl);
-#endif
diff --git a/modules/fdlibm/src/s_trunc.cpp b/modules/fdlibm/src/s_trunc.cpp
--- a/modules/fdlibm/src/s_trunc.cpp
+++ b/modules/fdlibm/src/s_trunc.cpp
@@ -55,12 +55,8 @@ trunc(double x)
i = ((u_int32_t)(0xffffffff))>>(j0-20);
if((i1&i)==0) return x; /* x is integral */
if(huge+x>0.0) /* raise inexact flag */
i1 &= (~i);
}
INSERT_WORDS(x,i0,i1);
return x;
}
-
-#if LDBL_MANT_DIG == 53
-__weak_reference(trunc, truncl);
-#endif
@@ -0,0 +1,753 @@
diff --git a/modules/fdlibm/src/e_acos.cpp b/modules/fdlibm/src/e_acos.cpp
--- a/modules/fdlibm/src/e_acos.cpp
+++ b/modules/fdlibm/src/e_acos.cpp
@@ -6,18 +6,18 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_acos(x)
* Method :
* acos(x) = pi/2 - asin(x)
* acos(-x) = pi/2 + asin(x)
* For |x|<=0.5
* acos(x) = pi/2 - (x + x*x^2*R(x^2)) (see asin.c)
* For x>0.5
diff --git a/modules/fdlibm/src/e_acosh.cpp b/modules/fdlibm/src/e_acosh.cpp
--- a/modules/fdlibm/src/e_acosh.cpp
+++ b/modules/fdlibm/src/e_acosh.cpp
@@ -7,18 +7,18 @@
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_acosh(x)
* Method :
* Based on
* acosh(x) = log [ x + sqrt(x*x-1) ]
* we have
* acosh(x) := log(x)+ln2, if x is large; else
* acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
diff --git a/modules/fdlibm/src/e_asin.cpp b/modules/fdlibm/src/e_asin.cpp
--- a/modules/fdlibm/src/e_asin.cpp
+++ b/modules/fdlibm/src/e_asin.cpp
@@ -6,18 +6,18 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_asin(x)
* Method :
* Since asin(x) = x + x^3/6 + x^5*3/40 + x^7*15/336 + ...
* we approximate asin(x) on [0,0.5] by
* asin(x) = x + x*x^2*R(x^2)
* where
* R(x^2) is a rational approximation of (asin(x)-x)/x^3
diff --git a/modules/fdlibm/src/e_atan2.cpp b/modules/fdlibm/src/e_atan2.cpp
--- a/modules/fdlibm/src/e_atan2.cpp
+++ b/modules/fdlibm/src/e_atan2.cpp
@@ -7,18 +7,18 @@
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_atan2(y,x)
* Method :
* 1. Reduce y to positive by atan2(y,x)=-atan2(-y,x).
* 2. Reduce x to positive by (if x and y are unexceptional):
* ARG (x+iy) = arctan(y/x) ... if x > 0,
* ARG (x+iy) = pi - arctan[y/(-x)] ... if x < 0,
*
diff --git a/modules/fdlibm/src/e_atanh.cpp b/modules/fdlibm/src/e_atanh.cpp
--- a/modules/fdlibm/src/e_atanh.cpp
+++ b/modules/fdlibm/src/e_atanh.cpp
@@ -7,18 +7,18 @@
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_atanh(x)
* Method :
* 1.Reduced x to positive by atanh(-x) = -atanh(x)
* 2.For x>=0.5
* 1 2x x
* atanh(x) = --- * log(1 + -------) = 0.5 * log1p(2 * --------)
* 2 1 - x 1 - x
diff --git a/modules/fdlibm/src/e_cosh.cpp b/modules/fdlibm/src/e_cosh.cpp
--- a/modules/fdlibm/src/e_cosh.cpp
+++ b/modules/fdlibm/src/e_cosh.cpp
@@ -6,18 +6,18 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_cosh(x)
* Method :
* mathematically cosh(x) if defined to be (exp(x)+exp(-x))/2
* 1. Replace x by |x| (cosh(x) = cosh(-x)).
* 2.
* [ exp(x) - 1 ]^2
* 0 <= x <= ln2/2 : cosh(x) := 1 + -------------------
diff --git a/modules/fdlibm/src/e_exp.cpp b/modules/fdlibm/src/e_exp.cpp
--- a/modules/fdlibm/src/e_exp.cpp
+++ b/modules/fdlibm/src/e_exp.cpp
@@ -5,18 +5,18 @@
* Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
*
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_exp(x)
* Returns the exponential of x.
*
* Method
* 1. Argument reduction:
* Reduce x to an r so that |r| <= 0.5*ln2 ~ 0.34658.
* Given x, find r and integer k such that
diff --git a/modules/fdlibm/src/e_hypot.cpp b/modules/fdlibm/src/e_hypot.cpp
--- a/modules/fdlibm/src/e_hypot.cpp
+++ b/modules/fdlibm/src/e_hypot.cpp
@@ -6,18 +6,18 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_hypot(x,y)
*
* Method :
* If (assume round-to-nearest) z=x*x+y*y
* has error less than sqrt(2)/2 ulp, than
* sqrt(z) has error less than 1 ulp (exercise).
*
diff --git a/modules/fdlibm/src/e_log.cpp b/modules/fdlibm/src/e_log.cpp
--- a/modules/fdlibm/src/e_log.cpp
+++ b/modules/fdlibm/src/e_log.cpp
@@ -6,18 +6,18 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_log(x)
* Return the logrithm of x
*
* Method :
* 1. Argument Reduction: find k and f such that
* x = 2^k * (1+f),
* where sqrt(2)/2 < 1+f < sqrt(2) .
diff --git a/modules/fdlibm/src/e_log10.cpp b/modules/fdlibm/src/e_log10.cpp
--- a/modules/fdlibm/src/e_log10.cpp
+++ b/modules/fdlibm/src/e_log10.cpp
@@ -6,32 +6,32 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/*
* Return the base 10 logarithm of x. See e_log.c and k_log.h for most
* comments.
*
* log10(x) = (f - 0.5*f*f + k_log1p(f)) / ln10 + k * log10(2)
* in not-quite-routine extra precision.
*/
#include <float.h>
#include "math_private.h"
#include "k_log.h"
-
+
static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
ivln10hi = 4.34294481878168880939e-01, /* 0x3fdbcb7b, 0x15200000 */
ivln10lo = 2.50829467116452752298e-11, /* 0x3dbb9438, 0xca9aadd5 */
log10_2hi = 3.01029995663611771306e-01, /* 0x3FD34413, 0x509F6000 */
log10_2lo = 3.69423907715893078616e-13; /* 0x3D59FEF3, 0x11F12B36 */
static const double zero = 0.0;
diff --git a/modules/fdlibm/src/e_log2.cpp b/modules/fdlibm/src/e_log2.cpp
--- a/modules/fdlibm/src/e_log2.cpp
+++ b/modules/fdlibm/src/e_log2.cpp
@@ -6,18 +6,18 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/*
* Return the base 2 logarithm of x. See e_log.c and k_log.h for most
* comments.
*
* This reduces x to {k, 1+f} exactly as in e_log.c, then calls the kernel,
* then does the combining and scaling steps
* log2(x) = (f - 0.5*f*f + k_log1p(f)) / ln2 + k
diff --git a/modules/fdlibm/src/e_pow.cpp b/modules/fdlibm/src/e_pow.cpp
--- a/modules/fdlibm/src/e_pow.cpp
+++ b/modules/fdlibm/src/e_pow.cpp
@@ -4,18 +4,18 @@
* Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
*
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_pow(x,y) return x**y
*
* n
* Method: Let x = 2 * (1+f)
* 1. Compute and return log2(x) in two pieces:
* log2(x) = w1 + w2,
* where w1 has 53-24 = 29 bit trailing zeros.
diff --git a/modules/fdlibm/src/e_sinh.cpp b/modules/fdlibm/src/e_sinh.cpp
--- a/modules/fdlibm/src/e_sinh.cpp
+++ b/modules/fdlibm/src/e_sinh.cpp
@@ -6,18 +6,18 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_sinh(x)
* Method :
* mathematically sinh(x) if defined to be (exp(x)-exp(-x))/2
* 1. Replace x by |x| (sinh(-x) = -sinh(x)).
* 2.
* E + E/(E+1)
* 0 <= x <= 22 : sinh(x) := --------------, E=expm1(x)
diff --git a/modules/fdlibm/src/e_sqrt.cpp b/modules/fdlibm/src/e_sqrt.cpp
--- a/modules/fdlibm/src/e_sqrt.cpp
+++ b/modules/fdlibm/src/e_sqrt.cpp
@@ -6,18 +6,18 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* __ieee754_sqrt(x)
* Return correctly rounded sqrt.
* ------------------------------------------
* | Use the hardware sqrt if you have one |
* ------------------------------------------
* Method:
* Bit by bit method using integer arithmetic. (Slow, but portable)
diff --git a/modules/fdlibm/src/k_exp.cpp b/modules/fdlibm/src/k_exp.cpp
--- a/modules/fdlibm/src/k_exp.cpp
+++ b/modules/fdlibm/src/k_exp.cpp
@@ -19,22 +19,22 @@
* DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
* OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
* HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
* LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
* OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
* SUCH DAMAGE.
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
#include <complex.h>
-#include "math_private.h"
+ #include "math_private.h"
static const uint32_t k = 1799; /* constant for reduction */
static const double kln2 = 1246.97177782734161156; /* k * ln2 */
/*
* Compute exp(x), scaled to avoid spurious overflow. An exponent is
* returned separately in 'expt'.
*
diff --git a/modules/fdlibm/src/k_log.h b/modules/fdlibm/src/k_log.h
--- a/modules/fdlibm/src/k_log.h
+++ b/modules/fdlibm/src/k_log.h
@@ -6,18 +6,18 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/*
* k_log1p(f):
* Return log(1+f) - f for 1+f in ~[sqrt(2)/2, sqrt(2)].
*
* The following describes the overall strategy for computing
* logarithms in base e. The argument reduction and adding the final
* term of the polynomial are done by the caller for increased accuracy
diff --git a/modules/fdlibm/src/s_asinh.cpp b/modules/fdlibm/src/s_asinh.cpp
--- a/modules/fdlibm/src/s_asinh.cpp
+++ b/modules/fdlibm/src/s_asinh.cpp
@@ -5,18 +5,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* asinh(x)
* Method :
* Based on
* asinh(x) = sign(x) * log [ |x| + sqrt(x*x+1) ]
* we have
* asinh(x) := x if 1+x*x=1,
* := sign(x)*(log(x)+ln2)) for large |x|, else
diff --git a/modules/fdlibm/src/s_atan.cpp b/modules/fdlibm/src/s_atan.cpp
--- a/modules/fdlibm/src/s_atan.cpp
+++ b/modules/fdlibm/src/s_atan.cpp
@@ -5,18 +5,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* atan(x)
* Method
* 1. Reduce x to positive by atan(x) = -atan(-x).
* 2. According to the integer k=4t+0.25 chopped, t=x, the argument
* is further reduced to one of the following intervals and the
* arctangent of t is evaluated by the corresponding formula:
*
diff --git a/modules/fdlibm/src/s_cbrt.cpp b/modules/fdlibm/src/s_cbrt.cpp
--- a/modules/fdlibm/src/s_cbrt.cpp
+++ b/modules/fdlibm/src/s_cbrt.cpp
@@ -7,18 +7,18 @@
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
* Optimized by Bruce D. Evans.
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
#include "math_private.h"
/* cbrt(x)
* Return cube root of x
*/
static const u_int32_t
B1 = 715094163, /* B1 = (1023-1023/3-0.03306235651)*2**20 */
diff --git a/modules/fdlibm/src/s_ceil.cpp b/modules/fdlibm/src/s_ceil.cpp
--- a/modules/fdlibm/src/s_ceil.cpp
+++ b/modules/fdlibm/src/s_ceil.cpp
@@ -5,18 +5,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/*
* ceil(x)
* Return x rounded toward -inf to integral value
* Method:
* Bit twiddling.
* Exception:
* Inexact flag raised if x not equal to ceil(x).
diff --git a/modules/fdlibm/src/s_ceilf.cpp b/modules/fdlibm/src/s_ceilf.cpp
--- a/modules/fdlibm/src/s_ceilf.cpp
+++ b/modules/fdlibm/src/s_ceilf.cpp
@@ -8,18 +8,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
#include "math_private.h"
static const float huge = 1.0e30;
float
ceilf(float x)
{
diff --git a/modules/fdlibm/src/s_copysign.cpp b/modules/fdlibm/src/s_copysign.cpp
--- a/modules/fdlibm/src/s_copysign.cpp
+++ b/modules/fdlibm/src/s_copysign.cpp
@@ -5,18 +5,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/*
* copysign(double x, double y)
* copysign(x,y) returns a value with the magnitude of x and
* with the sign bit of y.
*/
#include "math_private.h"
diff --git a/modules/fdlibm/src/s_expm1.cpp b/modules/fdlibm/src/s_expm1.cpp
--- a/modules/fdlibm/src/s_expm1.cpp
+++ b/modules/fdlibm/src/s_expm1.cpp
@@ -5,18 +5,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* expm1(x)
* Returns exp(x)-1, the exponential of x minus 1.
*
* Method
* 1. Argument reduction:
* Given x, find r and integer k such that
*
diff --git a/modules/fdlibm/src/s_fabs.cpp b/modules/fdlibm/src/s_fabs.cpp
--- a/modules/fdlibm/src/s_fabs.cpp
+++ b/modules/fdlibm/src/s_fabs.cpp
@@ -1,22 +1,22 @@
-/* @(#)s_fabs.c 5.1 93/09/24 */
+ /* @(#)s_fabs.c 5.1 93/09/24 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
#ifndef lint
-static char rcsid[] = "$FreeBSD$";
+ //static char rcsid[] = "$FreeBSD$";
#endif
/*
* fabs(x) returns the absolute value of x.
*/
#include "math_private.h"
diff --git a/modules/fdlibm/src/s_floor.cpp b/modules/fdlibm/src/s_floor.cpp
--- a/modules/fdlibm/src/s_floor.cpp
+++ b/modules/fdlibm/src/s_floor.cpp
@@ -5,18 +5,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/*
* floor(x)
* Return x rounded toward -inf to integral value
* Method:
* Bit twiddling.
* Exception:
* Inexact flag raised if x not equal to floor(x).
diff --git a/modules/fdlibm/src/s_floorf.cpp b/modules/fdlibm/src/s_floorf.cpp
--- a/modules/fdlibm/src/s_floorf.cpp
+++ b/modules/fdlibm/src/s_floorf.cpp
@@ -8,18 +8,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/*
* floorf(x)
* Return x rounded toward -inf to integral value
* Method:
* Bit twiddling.
* Exception:
* Inexact flag raised if x not equal to floorf(x).
diff --git a/modules/fdlibm/src/s_log1p.cpp b/modules/fdlibm/src/s_log1p.cpp
--- a/modules/fdlibm/src/s_log1p.cpp
+++ b/modules/fdlibm/src/s_log1p.cpp
@@ -5,18 +5,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* double log1p(double x)
*
* Method :
* 1. Argument Reduction: find k and f such that
* 1+x = 2^k * (1+f),
* where sqrt(2)/2 < 1+f < sqrt(2) .
*
diff --git a/modules/fdlibm/src/s_scalbn.cpp b/modules/fdlibm/src/s_scalbn.cpp
--- a/modules/fdlibm/src/s_scalbn.cpp
+++ b/modules/fdlibm/src/s_scalbn.cpp
@@ -6,27 +6,27 @@
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
#ifndef lint
-static char rcsid[] = "$FreeBSD$";
+//static char rcsid[] = "$FreeBSD$";
#endif
/*
* scalbn (double x, int n)
* scalbn(x,n) returns x* 2**n computed by exponent
* manipulation rather than by actually performing an
* exponentiation or a multiplication.
*/
-#include <sys/cdefs.h>
+//#include <sys/cdefs.h>
#include <float.h>
#include "math_private.h"
static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
huge = 1.0e+300,
diff --git a/modules/fdlibm/src/s_tanh.cpp b/modules/fdlibm/src/s_tanh.cpp
--- a/modules/fdlibm/src/s_tanh.cpp
+++ b/modules/fdlibm/src/s_tanh.cpp
@@ -5,18 +5,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/* Tanh(x)
* Return the Hyperbolic Tangent of x
*
* Method :
* x -x
* e - e
* 0. tanh(x) is defined to be -----------
diff --git a/modules/fdlibm/src/s_trunc.cpp b/modules/fdlibm/src/s_trunc.cpp
--- a/modules/fdlibm/src/s_trunc.cpp
+++ b/modules/fdlibm/src/s_trunc.cpp
@@ -5,18 +5,18 @@
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
+//#include <sys/cdefs.h>
+//__FBSDID("$FreeBSD$");
/*
* trunc(x)
* Return x rounded toward 0 to integral value
* Method:
* Bit twiddling.
* Exception:
* Inexact flag raised if x not equal to trunc(x).
@@ -0,0 +1,55 @@
diff --git a/modules/fdlibm/src/k_exp.cpp b/modules/fdlibm/src/k_exp.cpp
--- a/modules/fdlibm/src/k_exp.cpp
+++ b/modules/fdlibm/src/k_exp.cpp
@@ -22,18 +22,16 @@
* LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
* OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
* SUCH DAMAGE.
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
-#include <complex.h>
-
#include "math_private.h"
static const uint32_t k = 1799; /* constant for reduction */
static const double kln2 = 1246.97177782734161156; /* k * ln2 */
/*
* Compute exp(x), scaled to avoid spurious overflow. An exponent is
* returned separately in 'expt'.
@@ -76,32 +74,8 @@ double
double exp_x, scale;
int ex_expt;
exp_x = __frexp_exp(x, &ex_expt);
expt += ex_expt;
INSERT_WORDS(scale, (0x3ff + expt) << 20, 0);
return (exp_x * scale);
}
-
-double complex
-__ldexp_cexp(double complex z, int expt)
-{
- double x, y, exp_x, scale1, scale2;
- int ex_expt, half_expt;
-
- x = creal(z);
- y = cimag(z);
- exp_x = __frexp_exp(x, &ex_expt);
- expt += ex_expt;
-
- /*
- * Arrange so that scale1 * scale2 == 2**expt. We use this to
- * compensate for scalbn being horrendously slow.
- */
- half_expt = expt / 2;
- INSERT_WORDS(scale1, (0x3ff + half_expt) << 20, 0);
- half_expt = expt - half_expt;
- INSERT_WORDS(scale2, (0x3ff + half_expt) << 20, 0);
-
- return (CMPLX(cos(y) * exp_x * scale1 * scale2,
- sin(y) * exp_x * scale1 * scale2));
-}
@@ -0,0 +1,21 @@
diff --git a/modules/fdlibm/src/math_private.h b/modules/fdlibm/src/math_private.h
--- a/modules/fdlibm/src/math_private.h
+++ b/modules/fdlibm/src/math_private.h
@@ -12,16 +12,17 @@
/*
* from: @(#)fdlibm.h 5.1 93/09/24
* $FreeBSD$
*/
#ifndef _MATH_PRIVATE_H_
#define _MATH_PRIVATE_H_
+#include <cfloat>
#include <stdint.h>
#include <sys/types.h>
#include "fdlibm.h"
#include "mozilla/Endian.h"
/*
@@ -0,0 +1,25 @@
diff --git a/modules/fdlibm/src/math_private.h b/modules/fdlibm/src/math_private.h
--- a/modules/fdlibm/src/math_private.h
+++ b/modules/fdlibm/src/math_private.h
@@ -33,16 +33,21 @@
* to dig two 32 bit words out of the 64 bit IEEE floating point
* value. That is non-ANSI, and, moreover, the gcc instruction
* scheduler gets it wrong. We instead use the following macros.
* Unlike the original code, we determine the endianness at compile
* time, not at run time; I don't see much benefit to selecting
* endianness at run time.
*/
+#ifdef WIN32
+#define u_int32_t uint32_t
+#define u_int64_t uint64_t
+#endif
+
/*
* A union which permits us to convert between a double and two 32 bit
* ints.
*/
#if MOZ_BIG_ENDIAN
typedef union
@@ -0,0 +1,31 @@
diff --git a/modules/fdlibm/src/math_private.h b/modules/fdlibm/src/math_private.h
--- a/modules/fdlibm/src/math_private.h
+++ b/modules/fdlibm/src/math_private.h
@@ -285,16 +285,27 @@ do { \
if (sizeof(type) >= sizeof(long double)) \
(lval) = (rval); \
else { \
__lval = (rval); \
(lval) = __lval; \
} \
} while (0)
#endif
+#else
+#define STRICT_ASSIGN(type, lval, rval) do { \
+ volatile type __lval; \
+ \
+ if (sizeof(type) >= sizeof(long double)) \
+ (lval) = (rval); \
+ else { \
+ __lval = (rval); \
+ (lval) = __lval; \
+ } \
+} while (0)
#endif /* FLT_EVAL_METHOD */
/* Support switching the mode to FP_PE if necessary. */
#if defined(__i386__) && !defined(NO_FPSETPREC)
#define ENTERI() \
long double __retval; \
fp_prec_t __oprec; \
\
@@ -0,0 +1,47 @@
diff --git a/modules/fdlibm/src/e_exp.cpp b/modules/fdlibm/src/e_exp.cpp
--- a/modules/fdlibm/src/e_exp.cpp
+++ b/modules/fdlibm/src/e_exp.cpp
@@ -146,14 +146,17 @@ double
if(k >= -1021)
INSERT_WORDS(twopk,0x3ff00000+(k<<20), 0);
else
INSERT_WORDS(twopk,0x3ff00000+((k+1000)<<20), 0);
c = x - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))));
if(k==0) return one-((x*c)/(c-2.0)-x);
else y = one-((lo-(x*c)/(2.0-c))-hi);
if(k >= -1021) {
- if (k==1024) return y*2.0*0x1p1023;
+ if (k==1024) {
+ double const_0x1p1023 = pow(2, 1023);
+ return y*2.0*const_0x1p1023;
+ }
return y*twopk;
} else {
return y*twopk*twom1000;
}
}
diff --git a/modules/fdlibm/src/s_expm1.cpp b/modules/fdlibm/src/s_expm1.cpp
--- a/modules/fdlibm/src/s_expm1.cpp
+++ b/modules/fdlibm/src/s_expm1.cpp
@@ -192,17 +192,20 @@ expm1(double x)
e -= hxs;
if(k== -1) return 0.5*(x-e)-0.5;
if(k==1) {
if(x < -0.25) return -2.0*(e-(x+0.5));
else return one+2.0*(x-e);
}
if (k <= -2 || k>56) { /* suffice to return exp(x)-1 */
y = one-(e-x);
- if (k == 1024) y = y*2.0*0x1p1023;
+ if (k == 1024) {
+ double const_0x1p1023 = pow(2, 1023);
+ y = y*2.0*const_0x1p1023;
+ }
else y = y*twopk;
return y-one;
}
t = one;
if(k<20) {
SET_HIGH_WORD(t,0x3ff00000 - (0x200000>>k)); /* t=1-2^-k */
y = t-(e-x);
y = y*twopk;
+106
View File
@@ -0,0 +1,106 @@
/* @(#)e_acos.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_acos(x)
* Method :
* acos(x) = pi/2 - asin(x)
* acos(-x) = pi/2 + asin(x)
* For |x|<=0.5
* acos(x) = pi/2 - (x + x*x^2*R(x^2)) (see asin.c)
* For x>0.5
* acos(x) = pi/2 - (pi/2 - 2asin(sqrt((1-x)/2)))
* = 2asin(sqrt((1-x)/2))
* = 2s + 2s*z*R(z) ...z=(1-x)/2, s=sqrt(z)
* = 2f + (2c + 2s*z*R(z))
* where f=hi part of s, and c = (z-f*f)/(s+f) is the correction term
* for f so that f+c ~ sqrt(z).
* For x<-0.5
* acos(x) = pi - 2asin(sqrt((1-|x|)/2))
* = pi - 0.5*(s+s*z*R(z)), where z=(1-|x|)/2,s=sqrt(z)
*
* Special cases:
* if x is NaN, return x itself;
* if |x|>1, return NaN with invalid signal.
*
* Function needed: sqrt
*/
#include <float.h>
#include "math_private.h"
static const double
one= 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
pi = 3.14159265358979311600e+00, /* 0x400921FB, 0x54442D18 */
pio2_hi = 1.57079632679489655800e+00; /* 0x3FF921FB, 0x54442D18 */
static volatile double
pio2_lo = 6.12323399573676603587e-17; /* 0x3C91A626, 0x33145C07 */
static const double
pS0 = 1.66666666666666657415e-01, /* 0x3FC55555, 0x55555555 */
pS1 = -3.25565818622400915405e-01, /* 0xBFD4D612, 0x03EB6F7D */
pS2 = 2.01212532134862925881e-01, /* 0x3FC9C155, 0x0E884455 */
pS3 = -4.00555345006794114027e-02, /* 0xBFA48228, 0xB5688F3B */
pS4 = 7.91534994289814532176e-04, /* 0x3F49EFE0, 0x7501B288 */
pS5 = 3.47933107596021167570e-05, /* 0x3F023DE1, 0x0DFDF709 */
qS1 = -2.40339491173441421878e+00, /* 0xC0033A27, 0x1C8A2D4B */
qS2 = 2.02094576023350569471e+00, /* 0x40002AE5, 0x9C598AC8 */
qS3 = -6.88283971605453293030e-01, /* 0xBFE6066C, 0x1B8D0159 */
qS4 = 7.70381505559019352791e-02; /* 0x3FB3B8C5, 0xB12E9282 */
double
__ieee754_acos(double x)
{
double z,p,q,r,w,s,c,df;
int32_t hx,ix;
GET_HIGH_WORD(hx,x);
ix = hx&0x7fffffff;
if(ix>=0x3ff00000) { /* |x| >= 1 */
u_int32_t lx;
GET_LOW_WORD(lx,x);
if(((ix-0x3ff00000)|lx)==0) { /* |x|==1 */
if(hx>0) return 0.0; /* acos(1) = 0 */
else return pi+2.0*pio2_lo; /* acos(-1)= pi */
}
return (x-x)/(x-x); /* acos(|x|>1) is NaN */
}
if(ix<0x3fe00000) { /* |x| < 0.5 */
if(ix<=0x3c600000) return pio2_hi+pio2_lo;/*if|x|<2**-57*/
z = x*x;
p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
r = p/q;
return pio2_hi - (x - (pio2_lo-x*r));
} else if (hx<0) { /* x < -0.5 */
z = (one+x)*0.5;
p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
s = sqrt(z);
r = p/q;
w = r*s-pio2_lo;
return pi - 2.0*(s+w);
} else { /* x > 0.5 */
z = (one-x)*0.5;
s = sqrt(z);
df = s;
SET_LOW_WORD(df,0);
c = (z-df*df)/(s+df);
p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
r = p/q;
w = r*s+c;
return 2.0*(df+w);
}
}
+63
View File
@@ -0,0 +1,63 @@
/* @(#)e_acosh.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_acosh(x)
* Method :
* Based on
* acosh(x) = log [ x + sqrt(x*x-1) ]
* we have
* acosh(x) := log(x)+ln2, if x is large; else
* acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
* acosh(x) := log1p(t+sqrt(2.0*t+t*t)); where t=x-1.
*
* Special cases:
* acosh(x) is NaN with signal if x<1.
* acosh(NaN) is NaN without signal.
*/
#include <float.h>
#include "math_private.h"
static const double
one = 1.0,
ln2 = 6.93147180559945286227e-01; /* 0x3FE62E42, 0xFEFA39EF */
double
__ieee754_acosh(double x)
{
double t;
int32_t hx;
u_int32_t lx;
EXTRACT_WORDS(hx,lx,x);
if(hx<0x3ff00000) { /* x < 1 */
return (x-x)/(x-x);
} else if(hx >=0x41b00000) { /* x > 2**28 */
if(hx >=0x7ff00000) { /* x is inf of NaN */
return x+x;
} else
return __ieee754_log(x)+ln2; /* acosh(huge)=log(2x) */
} else if(((hx-0x3ff00000)|lx)==0) {
return 0.0; /* acosh(1) = 0 */
} else if (hx > 0x40000000) { /* 2**28 > x > 2 */
t=x*x;
return __ieee754_log(2.0*x-one/(x+sqrt(t-one)));
} else { /* 1<x<2 */
t = x-one;
return log1p(t+sqrt(2.0*t+t*t));
}
}
+112
View File
@@ -0,0 +1,112 @@
/* @(#)e_asin.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_asin(x)
* Method :
* Since asin(x) = x + x^3/6 + x^5*3/40 + x^7*15/336 + ...
* we approximate asin(x) on [0,0.5] by
* asin(x) = x + x*x^2*R(x^2)
* where
* R(x^2) is a rational approximation of (asin(x)-x)/x^3
* and its remez error is bounded by
* |(asin(x)-x)/x^3 - R(x^2)| < 2^(-58.75)
*
* For x in [0.5,1]
* asin(x) = pi/2-2*asin(sqrt((1-x)/2))
* Let y = (1-x), z = y/2, s := sqrt(z), and pio2_hi+pio2_lo=pi/2;
* then for x>0.98
* asin(x) = pi/2 - 2*(s+s*z*R(z))
* = pio2_hi - (2*(s+s*z*R(z)) - pio2_lo)
* For x<=0.98, let pio4_hi = pio2_hi/2, then
* f = hi part of s;
* c = sqrt(z) - f = (z-f*f)/(s+f) ...f+c=sqrt(z)
* and
* asin(x) = pi/2 - 2*(s+s*z*R(z))
* = pio4_hi+(pio4-2s)-(2s*z*R(z)-pio2_lo)
* = pio4_hi+(pio4-2f)-(2s*z*R(z)-(pio2_lo+2c))
*
* Special cases:
* if x is NaN, return x itself;
* if |x|>1, return NaN with invalid signal.
*
*/
#include <float.h>
#include "math_private.h"
static const double
one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
huge = 1.000e+300,
pio2_hi = 1.57079632679489655800e+00, /* 0x3FF921FB, 0x54442D18 */
pio2_lo = 6.12323399573676603587e-17, /* 0x3C91A626, 0x33145C07 */
pio4_hi = 7.85398163397448278999e-01, /* 0x3FE921FB, 0x54442D18 */
/* coefficient for R(x^2) */
pS0 = 1.66666666666666657415e-01, /* 0x3FC55555, 0x55555555 */
pS1 = -3.25565818622400915405e-01, /* 0xBFD4D612, 0x03EB6F7D */
pS2 = 2.01212532134862925881e-01, /* 0x3FC9C155, 0x0E884455 */
pS3 = -4.00555345006794114027e-02, /* 0xBFA48228, 0xB5688F3B */
pS4 = 7.91534994289814532176e-04, /* 0x3F49EFE0, 0x7501B288 */
pS5 = 3.47933107596021167570e-05, /* 0x3F023DE1, 0x0DFDF709 */
qS1 = -2.40339491173441421878e+00, /* 0xC0033A27, 0x1C8A2D4B */
qS2 = 2.02094576023350569471e+00, /* 0x40002AE5, 0x9C598AC8 */
qS3 = -6.88283971605453293030e-01, /* 0xBFE6066C, 0x1B8D0159 */
qS4 = 7.70381505559019352791e-02; /* 0x3FB3B8C5, 0xB12E9282 */
double
__ieee754_asin(double x)
{
double t=0.0,w,p,q,c,r,s;
int32_t hx,ix;
GET_HIGH_WORD(hx,x);
ix = hx&0x7fffffff;
if(ix>= 0x3ff00000) { /* |x|>= 1 */
u_int32_t lx;
GET_LOW_WORD(lx,x);
if(((ix-0x3ff00000)|lx)==0)
/* asin(1)=+-pi/2 with inexact */
return x*pio2_hi+x*pio2_lo;
return (x-x)/(x-x); /* asin(|x|>1) is NaN */
} else if (ix<0x3fe00000) { /* |x|<0.5 */
if(ix<0x3e500000) { /* if |x| < 2**-26 */
if(huge+x>one) return x;/* return x with inexact if x!=0*/
}
t = x*x;
p = t*(pS0+t*(pS1+t*(pS2+t*(pS3+t*(pS4+t*pS5)))));
q = one+t*(qS1+t*(qS2+t*(qS3+t*qS4)));
w = p/q;
return x+x*w;
}
/* 1> |x|>= 0.5 */
w = one-fabs(x);
t = w*0.5;
p = t*(pS0+t*(pS1+t*(pS2+t*(pS3+t*(pS4+t*pS5)))));
q = one+t*(qS1+t*(qS2+t*(qS3+t*qS4)));
s = sqrt(t);
if(ix>=0x3FEF3333) { /* if |x| > 0.975 */
w = p/q;
t = pio2_hi-(2.0*(s+s*w)-pio2_lo);
} else {
w = s;
SET_LOW_WORD(w,0);
c = (t-w*w)/(s+w);
r = p/q;
p = 2.0*s*r-(pio2_lo-2.0*c);
q = pio4_hi-2.0*w;
t = pio4_hi-(p-q);
}
if(hx>0) return t; else return -t;
}
+124
View File
@@ -0,0 +1,124 @@
/* @(#)e_atan2.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_atan2(y,x)
* Method :
* 1. Reduce y to positive by atan2(y,x)=-atan2(-y,x).
* 2. Reduce x to positive by (if x and y are unexceptional):
* ARG (x+iy) = arctan(y/x) ... if x > 0,
* ARG (x+iy) = pi - arctan[y/(-x)] ... if x < 0,
*
* Special cases:
*
* ATAN2((anything), NaN ) is NaN;
* ATAN2(NAN , (anything) ) is NaN;
* ATAN2(+-0, +(anything but NaN)) is +-0 ;
* ATAN2(+-0, -(anything but NaN)) is +-pi ;
* ATAN2(+-(anything but 0 and NaN), 0) is +-pi/2;
* ATAN2(+-(anything but INF and NaN), +INF) is +-0 ;
* ATAN2(+-(anything but INF and NaN), -INF) is +-pi;
* ATAN2(+-INF,+INF ) is +-pi/4 ;
* ATAN2(+-INF,-INF ) is +-3pi/4;
* ATAN2(+-INF, (anything but,0,NaN, and INF)) is +-pi/2;
*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
#include <float.h>
#include "math_private.h"
static volatile double
tiny = 1.0e-300;
static const double
zero = 0.0,
pi_o_4 = 7.8539816339744827900E-01, /* 0x3FE921FB, 0x54442D18 */
pi_o_2 = 1.5707963267948965580E+00, /* 0x3FF921FB, 0x54442D18 */
pi = 3.1415926535897931160E+00; /* 0x400921FB, 0x54442D18 */
static volatile double
pi_lo = 1.2246467991473531772E-16; /* 0x3CA1A626, 0x33145C07 */
double
__ieee754_atan2(double y, double x)
{
double z;
int32_t k,m,hx,hy,ix,iy;
u_int32_t lx,ly;
EXTRACT_WORDS(hx,lx,x);
ix = hx&0x7fffffff;
EXTRACT_WORDS(hy,ly,y);
iy = hy&0x7fffffff;
if(((ix|((lx|-lx)>>31))>0x7ff00000)||
((iy|((ly|-ly)>>31))>0x7ff00000)) /* x or y is NaN */
return x+y;
if((hx-0x3ff00000|lx)==0) return atan(y); /* x=1.0 */
m = ((hy>>31)&1)|((hx>>30)&2); /* 2*sign(x)+sign(y) */
/* when y = 0 */
if((iy|ly)==0) {
switch(m) {
case 0:
case 1: return y; /* atan(+-0,+anything)=+-0 */
case 2: return pi+tiny;/* atan(+0,-anything) = pi */
case 3: return -pi-tiny;/* atan(-0,-anything) =-pi */
}
}
/* when x = 0 */
if((ix|lx)==0) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny;
/* when x is INF */
if(ix==0x7ff00000) {
if(iy==0x7ff00000) {
switch(m) {
case 0: return pi_o_4+tiny;/* atan(+INF,+INF) */
case 1: return -pi_o_4-tiny;/* atan(-INF,+INF) */
case 2: return 3.0*pi_o_4+tiny;/*atan(+INF,-INF)*/
case 3: return -3.0*pi_o_4-tiny;/*atan(-INF,-INF)*/
}
} else {
switch(m) {
case 0: return zero ; /* atan(+...,+INF) */
case 1: return -zero ; /* atan(-...,+INF) */
case 2: return pi+tiny ; /* atan(+...,-INF) */
case 3: return -pi-tiny ; /* atan(-...,-INF) */
}
}
}
/* when y is INF */
if(iy==0x7ff00000) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny;
/* compute y/x */
k = (iy-ix)>>20;
if(k > 60) { /* |y/x| > 2**60 */
z=pi_o_2+0.5*pi_lo;
m&=1;
}
else if(hx<0&&k<-60) z=0.0; /* 0 > |y|/x > -2**-60 */
else z=atan(fabs(y/x)); /* safe to do y/x */
switch (m) {
case 0: return z ; /* atan(+,+) */
case 1: return -z ; /* atan(-,+) */
case 2: return pi-(z-pi_lo);/* atan(+,-) */
default: /* case 3 */
return (z-pi_lo)-pi;/* atan(-,-) */
}
}
+63
View File
@@ -0,0 +1,63 @@
/* @(#)e_atanh.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_atanh(x)
* Method :
* 1.Reduced x to positive by atanh(-x) = -atanh(x)
* 2.For x>=0.5
* 1 2x x
* atanh(x) = --- * log(1 + -------) = 0.5 * log1p(2 * --------)
* 2 1 - x 1 - x
*
* For x<0.5
* atanh(x) = 0.5*log1p(2x+2x*x/(1-x))
*
* Special cases:
* atanh(x) is NaN if |x| > 1 with signal;
* atanh(NaN) is that NaN with no signal;
* atanh(+-1) is +-INF with signal.
*
*/
#include <float.h>
#include "math_private.h"
static const double one = 1.0, huge = 1e300;
static const double zero = 0.0;
double
__ieee754_atanh(double x)
{
double t;
int32_t hx,ix;
u_int32_t lx;
EXTRACT_WORDS(hx,lx,x);
ix = hx&0x7fffffff;
if ((ix|((lx|(-lx))>>31))>0x3ff00000) /* |x|>1 */
return (x-x)/(x-x);
if(ix==0x3ff00000)
return x/zero;
if(ix<0x3e300000&&(huge+x)>zero) return x; /* x<2**-28 */
SET_HIGH_WORD(x,ix);
if(ix<0x3fe00000) { /* x < 0.5 */
t = x+x;
t = 0.5*log1p(t+t*x/(one-x));
} else
t = 0.5*log1p((x+x)/(one-x));
if(hx>=0) return t; else return -t;
}
+80
View File
@@ -0,0 +1,80 @@
/* @(#)e_cosh.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_cosh(x)
* Method :
* mathematically cosh(x) if defined to be (exp(x)+exp(-x))/2
* 1. Replace x by |x| (cosh(x) = cosh(-x)).
* 2.
* [ exp(x) - 1 ]^2
* 0 <= x <= ln2/2 : cosh(x) := 1 + -------------------
* 2*exp(x)
*
* exp(x) + 1/exp(x)
* ln2/2 <= x <= 22 : cosh(x) := -------------------
* 2
* 22 <= x <= lnovft : cosh(x) := exp(x)/2
* lnovft <= x <= ln2ovft: cosh(x) := exp(x/2)/2 * exp(x/2)
* ln2ovft < x : cosh(x) := huge*huge (overflow)
*
* Special cases:
* cosh(x) is |x| if x is +INF, -INF, or NaN.
* only cosh(0)=1 is exact for finite x.
*/
#include <float.h>
#include "math_private.h"
static const double one = 1.0, half=0.5, huge = 1.0e300;
double
__ieee754_cosh(double x)
{
double t,w;
int32_t ix;
/* High word of |x|. */
GET_HIGH_WORD(ix,x);
ix &= 0x7fffffff;
/* x is INF or NaN */
if(ix>=0x7ff00000) return x*x;
/* |x| in [0,0.5*ln2], return 1+expm1(|x|)^2/(2*exp(|x|)) */
if(ix<0x3fd62e43) {
t = expm1(fabs(x));
w = one+t;
if (ix<0x3c800000) return w; /* cosh(tiny) = 1 */
return one+(t*t)/(w+w);
}
/* |x| in [0.5*ln2,22], return (exp(|x|)+1/exp(|x|)/2; */
if (ix < 0x40360000) {
t = __ieee754_exp(fabs(x));
return half*t+half/t;
}
/* |x| in [22, log(maxdouble)] return half*exp(|x|) */
if (ix < 0x40862E42) return half*__ieee754_exp(fabs(x));
/* |x| in [log(maxdouble), overflowthresold] */
if (ix<=0x408633CE)
return __ldexp_exp(fabs(x), -1);
/* |x| > overflowthresold, cosh(x) overflow */
return huge*huge;
}
+162
View File
@@ -0,0 +1,162 @@
/* @(#)e_exp.c 1.6 04/04/22 */
/*
* ====================================================
* Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
*
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_exp(x)
* Returns the exponential of x.
*
* Method
* 1. Argument reduction:
* Reduce x to an r so that |r| <= 0.5*ln2 ~ 0.34658.
* Given x, find r and integer k such that
*
* x = k*ln2 + r, |r| <= 0.5*ln2.
*
* Here r will be represented as r = hi-lo for better
* accuracy.
*
* 2. Approximation of exp(r) by a special rational function on
* the interval [0,0.34658]:
* Write
* R(r**2) = r*(exp(r)+1)/(exp(r)-1) = 2 + r*r/6 - r**4/360 + ...
* We use a special Remes algorithm on [0,0.34658] to generate
* a polynomial of degree 5 to approximate R. The maximum error
* of this polynomial approximation is bounded by 2**-59. In
* other words,
* R(z) ~ 2.0 + P1*z + P2*z**2 + P3*z**3 + P4*z**4 + P5*z**5
* (where z=r*r, and the values of P1 to P5 are listed below)
* and
* | 5 | -59
* | 2.0+P1*z+...+P5*z - R(z) | <= 2
* | |
* The computation of exp(r) thus becomes
* 2*r
* exp(r) = 1 + -------
* R - r
* r*R1(r)
* = 1 + r + ----------- (for better accuracy)
* 2 - R1(r)
* where
* 2 4 10
* R1(r) = r - (P1*r + P2*r + ... + P5*r ).
*
* 3. Scale back to obtain exp(x):
* From step 1, we have
* exp(x) = 2^k * exp(r)
*
* Special cases:
* exp(INF) is INF, exp(NaN) is NaN;
* exp(-INF) is 0, and
* for finite argument, only exp(0)=1 is exact.
*
* Accuracy:
* according to an error analysis, the error is always less than
* 1 ulp (unit in the last place).
*
* Misc. info.
* For IEEE double
* if x > 7.09782712893383973096e+02 then exp(x) overflow
* if x < -7.45133219101941108420e+02 then exp(x) underflow
*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
#include <float.h>
#include "math_private.h"
static const double
one = 1.0,
halF[2] = {0.5,-0.5,},
o_threshold= 7.09782712893383973096e+02, /* 0x40862E42, 0xFEFA39EF */
u_threshold= -7.45133219101941108420e+02, /* 0xc0874910, 0xD52D3051 */
ln2HI[2] ={ 6.93147180369123816490e-01, /* 0x3fe62e42, 0xfee00000 */
-6.93147180369123816490e-01,},/* 0xbfe62e42, 0xfee00000 */
ln2LO[2] ={ 1.90821492927058770002e-10, /* 0x3dea39ef, 0x35793c76 */
-1.90821492927058770002e-10,},/* 0xbdea39ef, 0x35793c76 */
invln2 = 1.44269504088896338700e+00, /* 0x3ff71547, 0x652b82fe */
P1 = 1.66666666666666019037e-01, /* 0x3FC55555, 0x5555553E */
P2 = -2.77777777770155933842e-03, /* 0xBF66C16C, 0x16BEBD93 */
P3 = 6.61375632143793436117e-05, /* 0x3F11566A, 0xAF25DE2C */
P4 = -1.65339022054652515390e-06, /* 0xBEBBBD41, 0xC5D26BF1 */
P5 = 4.13813679705723846039e-08; /* 0x3E663769, 0x72BEA4D0 */
static volatile double
huge = 1.0e+300,
twom1000= 9.33263618503218878990e-302; /* 2**-1000=0x01700000,0*/
double
__ieee754_exp(double x) /* default IEEE double exp */
{
double y,hi=0.0,lo=0.0,c,t,twopk;
int32_t k=0,xsb;
u_int32_t hx;
GET_HIGH_WORD(hx,x);
xsb = (hx>>31)&1; /* sign bit of x */
hx &= 0x7fffffff; /* high word of |x| */
/* filter out non-finite argument */
if(hx >= 0x40862E42) { /* if |x|>=709.78... */
if(hx>=0x7ff00000) {
u_int32_t lx;
GET_LOW_WORD(lx,x);
if(((hx&0xfffff)|lx)!=0)
return x+x; /* NaN */
else return (xsb==0)? x:0.0; /* exp(+-inf)={inf,0} */
}
if(x > o_threshold) return huge*huge; /* overflow */
if(x < u_threshold) return twom1000*twom1000; /* underflow */
}
/* argument reduction */
if(hx > 0x3fd62e42) { /* if |x| > 0.5 ln2 */
if(hx < 0x3FF0A2B2) { /* and |x| < 1.5 ln2 */
hi = x-ln2HI[xsb]; lo=ln2LO[xsb]; k = 1-xsb-xsb;
} else {
k = (int)(invln2*x+halF[xsb]);
t = k;
hi = x - t*ln2HI[0]; /* t*ln2HI is exact here */
lo = t*ln2LO[0];
}
STRICT_ASSIGN(double, x, hi - lo);
}
else if(hx < 0x3e300000) { /* when |x|<2**-28 */
if(huge+x>one) return one+x;/* trigger inexact */
}
else k = 0;
/* x is now in primary range */
t = x*x;
if(k >= -1021)
INSERT_WORDS(twopk,0x3ff00000+(k<<20), 0);
else
INSERT_WORDS(twopk,0x3ff00000+((k+1000)<<20), 0);
c = x - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))));
if(k==0) return one-((x*c)/(c-2.0)-x);
else y = one-((lo-(x*c)/(2.0-c))-hi);
if(k >= -1021) {
if (k==1024) {
double const_0x1p1023 = pow(2, 1023);
return y*2.0*const_0x1p1023;
}
return y*twopk;
} else {
return y*twopk*twom1000;
}
}
+126
View File
@@ -0,0 +1,126 @@
/* @(#)e_hypot.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_hypot(x,y)
*
* Method :
* If (assume round-to-nearest) z=x*x+y*y
* has error less than sqrt(2)/2 ulp, than
* sqrt(z) has error less than 1 ulp (exercise).
*
* So, compute sqrt(x*x+y*y) with some care as
* follows to get the error below 1 ulp:
*
* Assume x>y>0;
* (if possible, set rounding to round-to-nearest)
* 1. if x > 2y use
* x1*x1+(y*y+(x2*(x+x1))) for x*x+y*y
* where x1 = x with lower 32 bits cleared, x2 = x-x1; else
* 2. if x <= 2y use
* t1*y1+((x-y)*(x-y)+(t1*y2+t2*y))
* where t1 = 2x with lower 32 bits cleared, t2 = 2x-t1,
* y1= y with lower 32 bits chopped, y2 = y-y1.
*
* NOTE: scaling may be necessary if some argument is too
* large or too tiny
*
* Special cases:
* hypot(x,y) is INF if x or y is +INF or -INF; else
* hypot(x,y) is NAN if x or y is NAN.
*
* Accuracy:
* hypot(x,y) returns sqrt(x^2+y^2) with error less
* than 1 ulps (units in the last place)
*/
#include <float.h>
#include "math_private.h"
double
__ieee754_hypot(double x, double y)
{
double a,b,t1,t2,y1,y2,w;
int32_t j,k,ha,hb;
GET_HIGH_WORD(ha,x);
ha &= 0x7fffffff;
GET_HIGH_WORD(hb,y);
hb &= 0x7fffffff;
if(hb > ha) {a=y;b=x;j=ha; ha=hb;hb=j;} else {a=x;b=y;}
a = fabs(a);
b = fabs(b);
if((ha-hb)>0x3c00000) {return a+b;} /* x/y > 2**60 */
k=0;
if(ha > 0x5f300000) { /* a>2**500 */
if(ha >= 0x7ff00000) { /* Inf or NaN */
u_int32_t low;
/* Use original arg order iff result is NaN; quieten sNaNs. */
w = fabs(x+0.0)-fabs(y+0.0);
GET_LOW_WORD(low,a);
if(((ha&0xfffff)|low)==0) w = a;
GET_LOW_WORD(low,b);
if(((hb^0x7ff00000)|low)==0) w = b;
return w;
}
/* scale a and b by 2**-600 */
ha -= 0x25800000; hb -= 0x25800000; k += 600;
SET_HIGH_WORD(a,ha);
SET_HIGH_WORD(b,hb);
}
if(hb < 0x20b00000) { /* b < 2**-500 */
if(hb <= 0x000fffff) { /* subnormal b or 0 */
u_int32_t low;
GET_LOW_WORD(low,b);
if((hb|low)==0) return a;
t1=0;
SET_HIGH_WORD(t1,0x7fd00000); /* t1=2^1022 */
b *= t1;
a *= t1;
k -= 1022;
} else { /* scale a and b by 2^600 */
ha += 0x25800000; /* a *= 2^600 */
hb += 0x25800000; /* b *= 2^600 */
k -= 600;
SET_HIGH_WORD(a,ha);
SET_HIGH_WORD(b,hb);
}
}
/* medium size a and b */
w = a-b;
if (w>b) {
t1 = 0;
SET_HIGH_WORD(t1,ha);
t2 = a-t1;
w = sqrt(t1*t1-(b*(-b)-t2*(a+t1)));
} else {
a = a+a;
y1 = 0;
SET_HIGH_WORD(y1,hb);
y2 = b - y1;
t1 = 0;
SET_HIGH_WORD(t1,ha+0x00100000);
t2 = a - t1;
w = sqrt(t1*y1-(w*(-w)-(t1*y2+t2*b)));
}
if(k!=0) {
u_int32_t high;
t1 = 1.0;
GET_HIGH_WORD(high,t1);
SET_HIGH_WORD(t1,high+(k<<20));
return t1*w;
} else return w;
}
+142
View File
@@ -0,0 +1,142 @@
/* @(#)e_log.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_log(x)
* Return the logrithm of x
*
* Method :
* 1. Argument Reduction: find k and f such that
* x = 2^k * (1+f),
* where sqrt(2)/2 < 1+f < sqrt(2) .
*
* 2. Approximation of log(1+f).
* Let s = f/(2+f) ; based on log(1+f) = log(1+s) - log(1-s)
* = 2s + 2/3 s**3 + 2/5 s**5 + .....,
* = 2s + s*R
* We use a special Reme algorithm on [0,0.1716] to generate
* a polynomial of degree 14 to approximate R The maximum error
* of this polynomial approximation is bounded by 2**-58.45. In
* other words,
* 2 4 6 8 10 12 14
* R(z) ~ Lg1*s +Lg2*s +Lg3*s +Lg4*s +Lg5*s +Lg6*s +Lg7*s
* (the values of Lg1 to Lg7 are listed in the program)
* and
* | 2 14 | -58.45
* | Lg1*s +...+Lg7*s - R(z) | <= 2
* | |
* Note that 2s = f - s*f = f - hfsq + s*hfsq, where hfsq = f*f/2.
* In order to guarantee error in log below 1ulp, we compute log
* by
* log(1+f) = f - s*(f - R) (if f is not too large)
* log(1+f) = f - (hfsq - s*(hfsq+R)). (better accuracy)
*
* 3. Finally, log(x) = k*ln2 + log(1+f).
* = k*ln2_hi+(f-(hfsq-(s*(hfsq+R)+k*ln2_lo)))
* Here ln2 is split into two floating point number:
* ln2_hi + ln2_lo,
* where n*ln2_hi is always exact for |n| < 2000.
*
* Special cases:
* log(x) is NaN with signal if x < 0 (including -INF) ;
* log(+INF) is +INF; log(0) is -INF with signal;
* log(NaN) is that NaN with no signal.
*
* Accuracy:
* according to an error analysis, the error is always less than
* 1 ulp (unit in the last place).
*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
#include <float.h>
#include "math_private.h"
static const double
ln2_hi = 6.93147180369123816490e-01, /* 3fe62e42 fee00000 */
ln2_lo = 1.90821492927058770002e-10, /* 3dea39ef 35793c76 */
two54 = 1.80143985094819840000e+16, /* 43500000 00000000 */
Lg1 = 6.666666666666735130e-01, /* 3FE55555 55555593 */
Lg2 = 3.999999999940941908e-01, /* 3FD99999 9997FA04 */
Lg3 = 2.857142874366239149e-01, /* 3FD24924 94229359 */
Lg4 = 2.222219843214978396e-01, /* 3FCC71C5 1D8E78AF */
Lg5 = 1.818357216161805012e-01, /* 3FC74664 96CB03DE */
Lg6 = 1.531383769920937332e-01, /* 3FC39A09 D078C69F */
Lg7 = 1.479819860511658591e-01; /* 3FC2F112 DF3E5244 */
static const double zero = 0.0;
static volatile double vzero = 0.0;
double
__ieee754_log(double x)
{
double hfsq,f,s,z,R,w,t1,t2,dk;
int32_t k,hx,i,j;
u_int32_t lx;
EXTRACT_WORDS(hx,lx,x);
k=0;
if (hx < 0x00100000) { /* x < 2**-1022 */
if (((hx&0x7fffffff)|lx)==0)
return -two54/vzero; /* log(+-0)=-inf */
if (hx<0) return (x-x)/zero; /* log(-#) = NaN */
k -= 54; x *= two54; /* subnormal number, scale up x */
GET_HIGH_WORD(hx,x);
}
if (hx >= 0x7ff00000) return x+x;
k += (hx>>20)-1023;
hx &= 0x000fffff;
i = (hx+0x95f64)&0x100000;
SET_HIGH_WORD(x,hx|(i^0x3ff00000)); /* normalize x or x/2 */
k += (i>>20);
f = x-1.0;
if((0x000fffff&(2+hx))<3) { /* -2**-20 <= f < 2**-20 */
if(f==zero) {
if(k==0) {
return zero;
} else {
dk=(double)k;
return dk*ln2_hi+dk*ln2_lo;
}
}
R = f*f*(0.5-0.33333333333333333*f);
if(k==0) return f-R; else {dk=(double)k;
return dk*ln2_hi-((R-dk*ln2_lo)-f);}
}
s = f/(2.0+f);
dk = (double)k;
z = s*s;
i = hx-0x6147a;
w = z*z;
j = 0x6b851-hx;
t1= w*(Lg2+w*(Lg4+w*Lg6));
t2= z*(Lg1+w*(Lg3+w*(Lg5+w*Lg7)));
i |= j;
R = t2+t1;
if(i>0) {
hfsq=0.5*f*f;
if(k==0) return f-(hfsq-s*(hfsq+R)); else
return dk*ln2_hi-((hfsq-(s*(hfsq+R)+dk*ln2_lo))-f);
} else {
if(k==0) return f-s*(f-R); else
return dk*ln2_hi-((s*(f-R)-dk*ln2_lo)-f);
}
}
+89
View File
@@ -0,0 +1,89 @@
/* @(#)e_log10.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/*
* Return the base 10 logarithm of x. See e_log.c and k_log.h for most
* comments.
*
* log10(x) = (f - 0.5*f*f + k_log1p(f)) / ln10 + k * log10(2)
* in not-quite-routine extra precision.
*/
#include <float.h>
#include "math_private.h"
#include "k_log.h"
static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
ivln10hi = 4.34294481878168880939e-01, /* 0x3fdbcb7b, 0x15200000 */
ivln10lo = 2.50829467116452752298e-11, /* 0x3dbb9438, 0xca9aadd5 */
log10_2hi = 3.01029995663611771306e-01, /* 0x3FD34413, 0x509F6000 */
log10_2lo = 3.69423907715893078616e-13; /* 0x3D59FEF3, 0x11F12B36 */
static const double zero = 0.0;
static volatile double vzero = 0.0;
double
__ieee754_log10(double x)
{
double f,hfsq,hi,lo,r,val_hi,val_lo,w,y,y2;
int32_t i,k,hx;
u_int32_t lx;
EXTRACT_WORDS(hx,lx,x);
k=0;
if (hx < 0x00100000) { /* x < 2**-1022 */
if (((hx&0x7fffffff)|lx)==0)
return -two54/vzero; /* log(+-0)=-inf */
if (hx<0) return (x-x)/zero; /* log(-#) = NaN */
k -= 54; x *= two54; /* subnormal number, scale up x */
GET_HIGH_WORD(hx,x);
}
if (hx >= 0x7ff00000) return x+x;
if (hx == 0x3ff00000 && lx == 0)
return zero; /* log(1) = +0 */
k += (hx>>20)-1023;
hx &= 0x000fffff;
i = (hx+0x95f64)&0x100000;
SET_HIGH_WORD(x,hx|(i^0x3ff00000)); /* normalize x or x/2 */
k += (i>>20);
y = (double)k;
f = x - 1.0;
hfsq = 0.5*f*f;
r = k_log1p(f);
/* See e_log2.c for most details. */
hi = f - hfsq;
SET_LOW_WORD(hi,0);
lo = (f - hi) - hfsq + r;
val_hi = hi*ivln10hi;
y2 = y*log10_2hi;
val_lo = y*log10_2lo + (lo+hi)*ivln10lo + lo*ivln10hi;
/*
* Extra precision in for adding y*log10_2hi is not strictly needed
* since there is no very large cancellation near x = sqrt(2) or
* x = 1/sqrt(2), but we do it anyway since it costs little on CPUs
* with some parallelism and it reduces the error for many args.
*/
w = y2 + val_hi;
val_lo += (y2 - w) + val_hi;
val_hi = w;
return val_lo + val_hi;
}
+112
View File
@@ -0,0 +1,112 @@
/* @(#)e_log10.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/*
* Return the base 2 logarithm of x. See e_log.c and k_log.h for most
* comments.
*
* This reduces x to {k, 1+f} exactly as in e_log.c, then calls the kernel,
* then does the combining and scaling steps
* log2(x) = (f - 0.5*f*f + k_log1p(f)) / ln2 + k
* in not-quite-routine extra precision.
*/
#include <float.h>
#include "math_private.h"
#include "k_log.h"
static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
ivln2hi = 1.44269504072144627571e+00, /* 0x3ff71547, 0x65200000 */
ivln2lo = 1.67517131648865118353e-10; /* 0x3de705fc, 0x2eefa200 */
static const double zero = 0.0;
static volatile double vzero = 0.0;
double
__ieee754_log2(double x)
{
double f,hfsq,hi,lo,r,val_hi,val_lo,w,y;
int32_t i,k,hx;
u_int32_t lx;
EXTRACT_WORDS(hx,lx,x);
k=0;
if (hx < 0x00100000) { /* x < 2**-1022 */
if (((hx&0x7fffffff)|lx)==0)
return -two54/vzero; /* log(+-0)=-inf */
if (hx<0) return (x-x)/zero; /* log(-#) = NaN */
k -= 54; x *= two54; /* subnormal number, scale up x */
GET_HIGH_WORD(hx,x);
}
if (hx >= 0x7ff00000) return x+x;
if (hx == 0x3ff00000 && lx == 0)
return zero; /* log(1) = +0 */
k += (hx>>20)-1023;
hx &= 0x000fffff;
i = (hx+0x95f64)&0x100000;
SET_HIGH_WORD(x,hx|(i^0x3ff00000)); /* normalize x or x/2 */
k += (i>>20);
y = (double)k;
f = x - 1.0;
hfsq = 0.5*f*f;
r = k_log1p(f);
/*
* f-hfsq must (for args near 1) be evaluated in extra precision
* to avoid a large cancellation when x is near sqrt(2) or 1/sqrt(2).
* This is fairly efficient since f-hfsq only depends on f, so can
* be evaluated in parallel with R. Not combining hfsq with R also
* keeps R small (though not as small as a true `lo' term would be),
* so that extra precision is not needed for terms involving R.
*
* Compiler bugs involving extra precision used to break Dekker's
* theorem for spitting f-hfsq as hi+lo, unless double_t was used
* or the multi-precision calculations were avoided when double_t
* has extra precision. These problems are now automatically
* avoided as a side effect of the optimization of combining the
* Dekker splitting step with the clear-low-bits step.
*
* y must (for args near sqrt(2) and 1/sqrt(2)) be added in extra
* precision to avoid a very large cancellation when x is very near
* these values. Unlike the above cancellations, this problem is
* specific to base 2. It is strange that adding +-1 is so much
* harder than adding +-ln2 or +-log10_2.
*
* This uses Dekker's theorem to normalize y+val_hi, so the
* compiler bugs are back in some configurations, sigh. And I
* don't want to used double_t to avoid them, since that gives a
* pessimization and the support for avoiding the pessimization
* is not yet available.
*
* The multi-precision calculations for the multiplications are
* routine.
*/
hi = f - hfsq;
SET_LOW_WORD(hi,0);
lo = (f - hi) - hfsq + r;
val_hi = hi*ivln2hi;
val_lo = (lo+hi)*ivln2lo + lo*ivln2hi;
/* spadd(val_hi, val_lo, y), except for not using double_t: */
w = y + val_hi;
val_lo += (y - w) + val_hi;
val_hi = w;
return val_lo + val_hi;
}
+305
View File
@@ -0,0 +1,305 @@
/* @(#)e_pow.c 1.5 04/04/22 SMI */
/*
* ====================================================
* Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
*
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_pow(x,y) return x**y
*
* n
* Method: Let x = 2 * (1+f)
* 1. Compute and return log2(x) in two pieces:
* log2(x) = w1 + w2,
* where w1 has 53-24 = 29 bit trailing zeros.
* 2. Perform y*log2(x) = n+y' by simulating multi-precision
* arithmetic, where |y'|<=0.5.
* 3. Return x**y = 2**n*exp(y'*log2)
*
* Special cases:
* 1. (anything) ** 0 is 1
* 2. (anything) ** 1 is itself
* 3. (anything) ** NAN is NAN except 1 ** NAN = 1
* 4. NAN ** (anything except 0) is NAN
* 5. +-(|x| > 1) ** +INF is +INF
* 6. +-(|x| > 1) ** -INF is +0
* 7. +-(|x| < 1) ** +INF is +0
* 8. +-(|x| < 1) ** -INF is +INF
* 9. +-1 ** +-INF is 1
* 10. +0 ** (+anything except 0, NAN) is +0
* 11. -0 ** (+anything except 0, NAN, odd integer) is +0
* 12. +0 ** (-anything except 0, NAN) is +INF
* 13. -0 ** (-anything except 0, NAN, odd integer) is +INF
* 14. -0 ** (odd integer) = -( +0 ** (odd integer) )
* 15. +INF ** (+anything except 0,NAN) is +INF
* 16. +INF ** (-anything except 0,NAN) is +0
* 17. -INF ** (anything) = -0 ** (-anything)
* 18. (-anything) ** (integer) is (-1)**(integer)*(+anything**integer)
* 19. (-anything except 0 and inf) ** (non-integer) is NAN
*
* Accuracy:
* pow(x,y) returns x**y nearly rounded. In particular
* pow(integer,integer)
* always returns the correct integer provided it is
* representable.
*
* Constants :
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
#include "math_private.h"
static const double
bp[] = {1.0, 1.5,},
dp_h[] = { 0.0, 5.84962487220764160156e-01,}, /* 0x3FE2B803, 0x40000000 */
dp_l[] = { 0.0, 1.35003920212974897128e-08,}, /* 0x3E4CFDEB, 0x43CFD006 */
zero = 0.0,
one = 1.0,
two = 2.0,
two53 = 9007199254740992.0, /* 0x43400000, 0x00000000 */
huge = 1.0e300,
tiny = 1.0e-300,
/* poly coefs for (3/2)*(log(x)-2s-2/3*s**3 */
L1 = 5.99999999999994648725e-01, /* 0x3FE33333, 0x33333303 */
L2 = 4.28571428578550184252e-01, /* 0x3FDB6DB6, 0xDB6FABFF */
L3 = 3.33333329818377432918e-01, /* 0x3FD55555, 0x518F264D */
L4 = 2.72728123808534006489e-01, /* 0x3FD17460, 0xA91D4101 */
L5 = 2.30660745775561754067e-01, /* 0x3FCD864A, 0x93C9DB65 */
L6 = 2.06975017800338417784e-01, /* 0x3FCA7E28, 0x4A454EEF */
P1 = 1.66666666666666019037e-01, /* 0x3FC55555, 0x5555553E */
P2 = -2.77777777770155933842e-03, /* 0xBF66C16C, 0x16BEBD93 */
P3 = 6.61375632143793436117e-05, /* 0x3F11566A, 0xAF25DE2C */
P4 = -1.65339022054652515390e-06, /* 0xBEBBBD41, 0xC5D26BF1 */
P5 = 4.13813679705723846039e-08, /* 0x3E663769, 0x72BEA4D0 */
lg2 = 6.93147180559945286227e-01, /* 0x3FE62E42, 0xFEFA39EF */
lg2_h = 6.93147182464599609375e-01, /* 0x3FE62E43, 0x00000000 */
lg2_l = -1.90465429995776804525e-09, /* 0xBE205C61, 0x0CA86C39 */
ovt = 8.0085662595372944372e-0017, /* -(1024-log2(ovfl+.5ulp)) */
cp = 9.61796693925975554329e-01, /* 0x3FEEC709, 0xDC3A03FD =2/(3ln2) */
cp_h = 9.61796700954437255859e-01, /* 0x3FEEC709, 0xE0000000 =(float)cp */
cp_l = -7.02846165095275826516e-09, /* 0xBE3E2FE0, 0x145B01F5 =tail of cp_h*/
ivln2 = 1.44269504088896338700e+00, /* 0x3FF71547, 0x652B82FE =1/ln2 */
ivln2_h = 1.44269502162933349609e+00, /* 0x3FF71547, 0x60000000 =24b 1/ln2*/
ivln2_l = 1.92596299112661746887e-08; /* 0x3E54AE0B, 0xF85DDF44 =1/ln2 tail*/
double
__ieee754_pow(double x, double y)
{
double z,ax,z_h,z_l,p_h,p_l;
double y1,t1,t2,r,s,t,u,v,w;
int32_t i,j,k,yisint,n;
int32_t hx,hy,ix,iy;
u_int32_t lx,ly;
EXTRACT_WORDS(hx,lx,x);
EXTRACT_WORDS(hy,ly,y);
ix = hx&0x7fffffff; iy = hy&0x7fffffff;
/* y==zero: x**0 = 1 */
if((iy|ly)==0) return one;
/* x==1: 1**y = 1, even if y is NaN */
if (hx==0x3ff00000 && lx == 0) return one;
/* y!=zero: result is NaN if either arg is NaN */
if(ix > 0x7ff00000 || ((ix==0x7ff00000)&&(lx!=0)) ||
iy > 0x7ff00000 || ((iy==0x7ff00000)&&(ly!=0)))
return (x+0.0)+(y+0.0);
/* determine if y is an odd int when x < 0
* yisint = 0 ... y is not an integer
* yisint = 1 ... y is an odd int
* yisint = 2 ... y is an even int
*/
yisint = 0;
if(hx<0) {
if(iy>=0x43400000) yisint = 2; /* even integer y */
else if(iy>=0x3ff00000) {
k = (iy>>20)-0x3ff; /* exponent */
if(k>20) {
j = ly>>(52-k);
if((j<<(52-k))==ly) yisint = 2-(j&1);
} else if(ly==0) {
j = iy>>(20-k);
if((j<<(20-k))==iy) yisint = 2-(j&1);
}
}
}
/* special value of y */
if(ly==0) {
if (iy==0x7ff00000) { /* y is +-inf */
if(((ix-0x3ff00000)|lx)==0)
return one; /* (-1)**+-inf is 1 */
else if (ix >= 0x3ff00000)/* (|x|>1)**+-inf = inf,0 */
return (hy>=0)? y: zero;
else /* (|x|<1)**-,+inf = inf,0 */
return (hy<0)?-y: zero;
}
if(iy==0x3ff00000) { /* y is +-1 */
if(hy<0) return one/x; else return x;
}
if(hy==0x40000000) return x*x; /* y is 2 */
if(hy==0x3fe00000) { /* y is 0.5 */
if(hx>=0) /* x >= +0 */
return sqrt(x);
}
}
ax = fabs(x);
/* special value of x */
if(lx==0) {
if(ix==0x7ff00000||ix==0||ix==0x3ff00000){
z = ax; /*x is +-0,+-inf,+-1*/
if(hy<0) z = one/z; /* z = (1/|x|) */
if(hx<0) {
if(((ix-0x3ff00000)|yisint)==0) {
z = (z-z)/(z-z); /* (-1)**non-int is NaN */
} else if(yisint==1)
z = -z; /* (x<0)**odd = -(|x|**odd) */
}
return z;
}
}
/* CYGNUS LOCAL + fdlibm-5.3 fix: This used to be
n = (hx>>31)+1;
but ANSI C says a right shift of a signed negative quantity is
implementation defined. */
n = ((u_int32_t)hx>>31)-1;
/* (x<0)**(non-int) is NaN */
if((n|yisint)==0) return (x-x)/(x-x);
s = one; /* s (sign of result -ve**odd) = -1 else = 1 */
if((n|(yisint-1))==0) s = -one;/* (-ve)**(odd int) */
/* |y| is huge */
if(iy>0x41e00000) { /* if |y| > 2**31 */
if(iy>0x43f00000){ /* if |y| > 2**64, must o/uflow */
if(ix<=0x3fefffff) return (hy<0)? huge*huge:tiny*tiny;
if(ix>=0x3ff00000) return (hy>0)? huge*huge:tiny*tiny;
}
/* over/underflow if x is not close to one */
if(ix<0x3fefffff) return (hy<0)? s*huge*huge:s*tiny*tiny;
if(ix>0x3ff00000) return (hy>0)? s*huge*huge:s*tiny*tiny;
/* now |1-x| is tiny <= 2**-20, suffice to compute
log(x) by x-x^2/2+x^3/3-x^4/4 */
t = ax-one; /* t has 20 trailing zeros */
w = (t*t)*(0.5-t*(0.3333333333333333333333-t*0.25));
u = ivln2_h*t; /* ivln2_h has 21 sig. bits */
v = t*ivln2_l-w*ivln2;
t1 = u+v;
SET_LOW_WORD(t1,0);
t2 = v-(t1-u);
} else {
double ss,s2,s_h,s_l,t_h,t_l;
n = 0;
/* take care subnormal number */
if(ix<0x00100000)
{ax *= two53; n -= 53; GET_HIGH_WORD(ix,ax); }
n += ((ix)>>20)-0x3ff;
j = ix&0x000fffff;
/* determine interval */
ix = j|0x3ff00000; /* normalize ix */
if(j<=0x3988E) k=0; /* |x|<sqrt(3/2) */
else if(j<0xBB67A) k=1; /* |x|<sqrt(3) */
else {k=0;n+=1;ix -= 0x00100000;}
SET_HIGH_WORD(ax,ix);
/* compute ss = s_h+s_l = (x-1)/(x+1) or (x-1.5)/(x+1.5) */
u = ax-bp[k]; /* bp[0]=1.0, bp[1]=1.5 */
v = one/(ax+bp[k]);
ss = u*v;
s_h = ss;
SET_LOW_WORD(s_h,0);
/* t_h=ax+bp[k] High */
t_h = zero;
SET_HIGH_WORD(t_h,((ix>>1)|0x20000000)+0x00080000+(k<<18));
t_l = ax - (t_h-bp[k]);
s_l = v*((u-s_h*t_h)-s_h*t_l);
/* compute log(ax) */
s2 = ss*ss;
r = s2*s2*(L1+s2*(L2+s2*(L3+s2*(L4+s2*(L5+s2*L6)))));
r += s_l*(s_h+ss);
s2 = s_h*s_h;
t_h = 3.0+s2+r;
SET_LOW_WORD(t_h,0);
t_l = r-((t_h-3.0)-s2);
/* u+v = ss*(1+...) */
u = s_h*t_h;
v = s_l*t_h+t_l*ss;
/* 2/(3log2)*(ss+...) */
p_h = u+v;
SET_LOW_WORD(p_h,0);
p_l = v-(p_h-u);
z_h = cp_h*p_h; /* cp_h+cp_l = 2/(3*log2) */
z_l = cp_l*p_h+p_l*cp+dp_l[k];
/* log2(ax) = (ss+..)*2/(3*log2) = n + dp_h + z_h + z_l */
t = (double)n;
t1 = (((z_h+z_l)+dp_h[k])+t);
SET_LOW_WORD(t1,0);
t2 = z_l-(((t1-t)-dp_h[k])-z_h);
}
/* split up y into y1+y2 and compute (y1+y2)*(t1+t2) */
y1 = y;
SET_LOW_WORD(y1,0);
p_l = (y-y1)*t1+y*t2;
p_h = y1*t1;
z = p_l+p_h;
EXTRACT_WORDS(j,i,z);
if (j>=0x40900000) { /* z >= 1024 */
if(((j-0x40900000)|i)!=0) /* if z > 1024 */
return s*huge*huge; /* overflow */
else {
if(p_l+ovt>z-p_h) return s*huge*huge; /* overflow */
}
} else if((j&0x7fffffff)>=0x4090cc00 ) { /* z <= -1075 */
if(((j-0xc090cc00)|i)!=0) /* z < -1075 */
return s*tiny*tiny; /* underflow */
else {
if(p_l<=z-p_h) return s*tiny*tiny; /* underflow */
}
}
/*
* compute 2**(p_h+p_l)
*/
i = j&0x7fffffff;
k = (i>>20)-0x3ff;
n = 0;
if(i>0x3fe00000) { /* if |z| > 0.5, set n = [z+0.5] */
n = j+(0x00100000>>(k+1));
k = ((n&0x7fffffff)>>20)-0x3ff; /* new k for n */
t = zero;
SET_HIGH_WORD(t,n&~(0x000fffff>>k));
n = ((n&0x000fffff)|0x00100000)>>(20-k);
if(j<0) n = -n;
p_h -= t;
}
t = p_l+p_h;
SET_LOW_WORD(t,0);
u = t*lg2_h;
v = (p_l-(t-p_h))*lg2+t*lg2_l;
z = u+v;
w = v-(z-u);
t = z*z;
t1 = z - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))));
r = (z*t1)/(t1-two)-(w+z*w);
z = one-(r-z);
GET_HIGH_WORD(j,z);
j += (n<<20);
if((j>>20)<=0) z = scalbn(z,n); /* subnormal output */
else SET_HIGH_WORD(z,j);
return s*z;
}
+74
View File
@@ -0,0 +1,74 @@
/* @(#)e_sinh.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_sinh(x)
* Method :
* mathematically sinh(x) if defined to be (exp(x)-exp(-x))/2
* 1. Replace x by |x| (sinh(-x) = -sinh(x)).
* 2.
* E + E/(E+1)
* 0 <= x <= 22 : sinh(x) := --------------, E=expm1(x)
* 2
*
* 22 <= x <= lnovft : sinh(x) := exp(x)/2
* lnovft <= x <= ln2ovft: sinh(x) := exp(x/2)/2 * exp(x/2)
* ln2ovft < x : sinh(x) := x*shuge (overflow)
*
* Special cases:
* sinh(x) is |x| if x is +INF, -INF, or NaN.
* only sinh(0)=0 is exact for finite x.
*/
#include <float.h>
#include "math_private.h"
static const double one = 1.0, shuge = 1.0e307;
double
__ieee754_sinh(double x)
{
double t,h;
int32_t ix,jx;
/* High word of |x|. */
GET_HIGH_WORD(jx,x);
ix = jx&0x7fffffff;
/* x is INF or NaN */
if(ix>=0x7ff00000) return x+x;
h = 0.5;
if (jx<0) h = -h;
/* |x| in [0,22], return sign(x)*0.5*(E+E/(E+1))) */
if (ix < 0x40360000) { /* |x|<22 */
if (ix<0x3e300000) /* |x|<2**-28 */
if(shuge+x>one) return x;/* sinh(tiny) = tiny with inexact */
t = expm1(fabs(x));
if(ix<0x3ff00000) return h*(2.0*t-t*t/(t+one));
return h*(t+t/(t+one));
}
/* |x| in [22, log(maxdouble)] return 0.5*exp(|x|) */
if (ix < 0x40862E42) return h*__ieee754_exp(fabs(x));
/* |x| in [log(maxdouble), overflowthresold] */
if (ix<=0x408633CE)
return h*2.0*__ldexp_exp(fabs(x), -1);
/* |x| > overflowthresold, sinh(x) overflow */
return x*shuge;
}
+446
View File
@@ -0,0 +1,446 @@
/* @(#)e_sqrt.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/* __ieee754_sqrt(x)
* Return correctly rounded sqrt.
* ------------------------------------------
* | Use the hardware sqrt if you have one |
* ------------------------------------------
* Method:
* Bit by bit method using integer arithmetic. (Slow, but portable)
* 1. Normalization
* Scale x to y in [1,4) with even powers of 2:
* find an integer k such that 1 <= (y=x*2^(2k)) < 4, then
* sqrt(x) = 2^k * sqrt(y)
* 2. Bit by bit computation
* Let q = sqrt(y) truncated to i bit after binary point (q = 1),
* i 0
* i+1 2
* s = 2*q , and y = 2 * ( y - q ). (1)
* i i i i
*
* To compute q from q , one checks whether
* i+1 i
*
* -(i+1) 2
* (q + 2 ) <= y. (2)
* i
* -(i+1)
* If (2) is false, then q = q ; otherwise q = q + 2 .
* i+1 i i+1 i
*
* With some algebric manipulation, it is not difficult to see
* that (2) is equivalent to
* -(i+1)
* s + 2 <= y (3)
* i i
*
* The advantage of (3) is that s and y can be computed by
* i i
* the following recurrence formula:
* if (3) is false
*
* s = s , y = y ; (4)
* i+1 i i+1 i
*
* otherwise,
* -i -(i+1)
* s = s + 2 , y = y - s - 2 (5)
* i+1 i i+1 i i
*
* One may easily use induction to prove (4) and (5).
* Note. Since the left hand side of (3) contain only i+2 bits,
* it does not necessary to do a full (53-bit) comparison
* in (3).
* 3. Final rounding
* After generating the 53 bits result, we compute one more bit.
* Together with the remainder, we can decide whether the
* result is exact, bigger than 1/2ulp, or less than 1/2ulp
* (it will never equal to 1/2ulp).
* The rounding mode can be detected by checking whether
* huge + tiny is equal to huge, and whether huge - tiny is
* equal to huge for some floating point number "huge" and "tiny".
*
* Special cases:
* sqrt(+-0) = +-0 ... exact
* sqrt(inf) = inf
* sqrt(-ve) = NaN ... with invalid signal
* sqrt(NaN) = NaN ... with invalid signal for signaling NaN
*
* Other methods : see the appended file at the end of the program below.
*---------------
*/
#include <float.h>
#include "math_private.h"
static const double one = 1.0, tiny=1.0e-300;
double
__ieee754_sqrt(double x)
{
double z;
int32_t sign = (int)0x80000000;
int32_t ix0,s0,q,m,t,i;
u_int32_t r,t1,s1,ix1,q1;
EXTRACT_WORDS(ix0,ix1,x);
/* take care of Inf and NaN */
if((ix0&0x7ff00000)==0x7ff00000) {
return x*x+x; /* sqrt(NaN)=NaN, sqrt(+inf)=+inf
sqrt(-inf)=sNaN */
}
/* take care of zero */
if(ix0<=0) {
if(((ix0&(~sign))|ix1)==0) return x;/* sqrt(+-0) = +-0 */
else if(ix0<0)
return (x-x)/(x-x); /* sqrt(-ve) = sNaN */
}
/* normalize x */
m = (ix0>>20);
if(m==0) { /* subnormal x */
while(ix0==0) {
m -= 21;
ix0 |= (ix1>>11); ix1 <<= 21;
}
for(i=0;(ix0&0x00100000)==0;i++) ix0<<=1;
m -= i-1;
ix0 |= (ix1>>(32-i));
ix1 <<= i;
}
m -= 1023; /* unbias exponent */
ix0 = (ix0&0x000fffff)|0x00100000;
if(m&1){ /* odd m, double x to make it even */
ix0 += ix0 + ((ix1&sign)>>31);
ix1 += ix1;
}
m >>= 1; /* m = [m/2] */
/* generate sqrt(x) bit by bit */
ix0 += ix0 + ((ix1&sign)>>31);
ix1 += ix1;
q = q1 = s0 = s1 = 0; /* [q,q1] = sqrt(x) */
r = 0x00200000; /* r = moving bit from right to left */
while(r!=0) {
t = s0+r;
if(t<=ix0) {
s0 = t+r;
ix0 -= t;
q += r;
}
ix0 += ix0 + ((ix1&sign)>>31);
ix1 += ix1;
r>>=1;
}
r = sign;
while(r!=0) {
t1 = s1+r;
t = s0;
if((t<ix0)||((t==ix0)&&(t1<=ix1))) {
s1 = t1+r;
if(((t1&sign)==sign)&&(s1&sign)==0) s0 += 1;
ix0 -= t;
if (ix1 < t1) ix0 -= 1;
ix1 -= t1;
q1 += r;
}
ix0 += ix0 + ((ix1&sign)>>31);
ix1 += ix1;
r>>=1;
}
/* use floating add to find out rounding direction */
if((ix0|ix1)!=0) {
z = one-tiny; /* trigger inexact flag */
if (z>=one) {
z = one+tiny;
if (q1==(u_int32_t)0xffffffff) { q1=0; q += 1;}
else if (z>one) {
if (q1==(u_int32_t)0xfffffffe) q+=1;
q1+=2;
} else
q1 += (q1&1);
}
}
ix0 = (q>>1)+0x3fe00000;
ix1 = q1>>1;
if ((q&1)==1) ix1 |= sign;
ix0 += (m <<20);
INSERT_WORDS(z,ix0,ix1);
return z;
}
/*
Other methods (use floating-point arithmetic)
-------------
(This is a copy of a drafted paper by Prof W. Kahan
and K.C. Ng, written in May, 1986)
Two algorithms are given here to implement sqrt(x)
(IEEE double precision arithmetic) in software.
Both supply sqrt(x) correctly rounded. The first algorithm (in
Section A) uses newton iterations and involves four divisions.
The second one uses reciproot iterations to avoid division, but
requires more multiplications. Both algorithms need the ability
to chop results of arithmetic operations instead of round them,
and the INEXACT flag to indicate when an arithmetic operation
is executed exactly with no roundoff error, all part of the
standard (IEEE 754-1985). The ability to perform shift, add,
subtract and logical AND operations upon 32-bit words is needed
too, though not part of the standard.
A. sqrt(x) by Newton Iteration
(1) Initial approximation
Let x0 and x1 be the leading and the trailing 32-bit words of
a floating point number x (in IEEE double format) respectively
1 11 52 ...widths
------------------------------------------------------
x: |s| e | f |
------------------------------------------------------
msb lsb msb lsb ...order
------------------------ ------------------------
x0: |s| e | f1 | x1: | f2 |
------------------------ ------------------------
By performing shifts and subtracts on x0 and x1 (both regarded
as integers), we obtain an 8-bit approximation of sqrt(x) as
follows.
k := (x0>>1) + 0x1ff80000;
y0 := k - T1[31&(k>>15)]. ... y ~ sqrt(x) to 8 bits
Here k is a 32-bit integer and T1[] is an integer array containing
correction terms. Now magically the floating value of y (y's
leading 32-bit word is y0, the value of its trailing word is 0)
approximates sqrt(x) to almost 8-bit.
Value of T1:
static int T1[32]= {
0, 1024, 3062, 5746, 9193, 13348, 18162, 23592,
29598, 36145, 43202, 50740, 58733, 67158, 75992, 85215,
83599, 71378, 60428, 50647, 41945, 34246, 27478, 21581,
16499, 12183, 8588, 5674, 3403, 1742, 661, 130,};
(2) Iterative refinement
Apply Heron's rule three times to y, we have y approximates
sqrt(x) to within 1 ulp (Unit in the Last Place):
y := (y+x/y)/2 ... almost 17 sig. bits
y := (y+x/y)/2 ... almost 35 sig. bits
y := y-(y-x/y)/2 ... within 1 ulp
Remark 1.
Another way to improve y to within 1 ulp is:
y := (y+x/y) ... almost 17 sig. bits to 2*sqrt(x)
y := y - 0x00100006 ... almost 18 sig. bits to sqrt(x)
2
(x-y )*y
y := y + 2* ---------- ...within 1 ulp
2
3y + x
This formula has one division fewer than the one above; however,
it requires more multiplications and additions. Also x must be
scaled in advance to avoid spurious overflow in evaluating the
expression 3y*y+x. Hence it is not recommended uless division
is slow. If division is very slow, then one should use the
reciproot algorithm given in section B.
(3) Final adjustment
By twiddling y's last bit it is possible to force y to be
correctly rounded according to the prevailing rounding mode
as follows. Let r and i be copies of the rounding mode and
inexact flag before entering the square root program. Also we
use the expression y+-ulp for the next representable floating
numbers (up and down) of y. Note that y+-ulp = either fixed
point y+-1, or multiply y by nextafter(1,+-inf) in chopped
mode.
I := FALSE; ... reset INEXACT flag I
R := RZ; ... set rounding mode to round-toward-zero
z := x/y; ... chopped quotient, possibly inexact
If(not I) then { ... if the quotient is exact
if(z=y) {
I := i; ... restore inexact flag
R := r; ... restore rounded mode
return sqrt(x):=y.
} else {
z := z - ulp; ... special rounding
}
}
i := TRUE; ... sqrt(x) is inexact
If (r=RN) then z=z+ulp ... rounded-to-nearest
If (r=RP) then { ... round-toward-+inf
y = y+ulp; z=z+ulp;
}
y := y+z; ... chopped sum
y0:=y0-0x00100000; ... y := y/2 is correctly rounded.
I := i; ... restore inexact flag
R := r; ... restore rounded mode
return sqrt(x):=y.
(4) Special cases
Square root of +inf, +-0, or NaN is itself;
Square root of a negative number is NaN with invalid signal.
B. sqrt(x) by Reciproot Iteration
(1) Initial approximation
Let x0 and x1 be the leading and the trailing 32-bit words of
a floating point number x (in IEEE double format) respectively
(see section A). By performing shifs and subtracts on x0 and y0,
we obtain a 7.8-bit approximation of 1/sqrt(x) as follows.
k := 0x5fe80000 - (x0>>1);
y0:= k - T2[63&(k>>14)]. ... y ~ 1/sqrt(x) to 7.8 bits
Here k is a 32-bit integer and T2[] is an integer array
containing correction terms. Now magically the floating
value of y (y's leading 32-bit word is y0, the value of
its trailing word y1 is set to zero) approximates 1/sqrt(x)
to almost 7.8-bit.
Value of T2:
static int T2[64]= {
0x1500, 0x2ef8, 0x4d67, 0x6b02, 0x87be, 0xa395, 0xbe7a, 0xd866,
0xf14a, 0x1091b,0x11fcd,0x13552,0x14999,0x15c98,0x16e34,0x17e5f,
0x18d03,0x19a01,0x1a545,0x1ae8a,0x1b5c4,0x1bb01,0x1bfde,0x1c28d,
0x1c2de,0x1c0db,0x1ba73,0x1b11c,0x1a4b5,0x1953d,0x18266,0x16be0,
0x1683e,0x179d8,0x18a4d,0x19992,0x1a789,0x1b445,0x1bf61,0x1c989,
0x1d16d,0x1d77b,0x1dddf,0x1e2ad,0x1e5bf,0x1e6e8,0x1e654,0x1e3cd,
0x1df2a,0x1d635,0x1cb16,0x1be2c,0x1ae4e,0x19bde,0x1868e,0x16e2e,
0x1527f,0x1334a,0x11051,0xe951, 0xbe01, 0x8e0d, 0x5924, 0x1edd,};
(2) Iterative refinement
Apply Reciproot iteration three times to y and multiply the
result by x to get an approximation z that matches sqrt(x)
to about 1 ulp. To be exact, we will have
-1ulp < sqrt(x)-z<1.0625ulp.
... set rounding mode to Round-to-nearest
y := y*(1.5-0.5*x*y*y) ... almost 15 sig. bits to 1/sqrt(x)
y := y*((1.5-2^-30)+0.5*x*y*y)... about 29 sig. bits to 1/sqrt(x)
... special arrangement for better accuracy
z := x*y ... 29 bits to sqrt(x), with z*y<1
z := z + 0.5*z*(1-z*y) ... about 1 ulp to sqrt(x)
Remark 2. The constant 1.5-2^-30 is chosen to bias the error so that
(a) the term z*y in the final iteration is always less than 1;
(b) the error in the final result is biased upward so that
-1 ulp < sqrt(x) - z < 1.0625 ulp
instead of |sqrt(x)-z|<1.03125ulp.
(3) Final adjustment
By twiddling y's last bit it is possible to force y to be
correctly rounded according to the prevailing rounding mode
as follows. Let r and i be copies of the rounding mode and
inexact flag before entering the square root program. Also we
use the expression y+-ulp for the next representable floating
numbers (up and down) of y. Note that y+-ulp = either fixed
point y+-1, or multiply y by nextafter(1,+-inf) in chopped
mode.
R := RZ; ... set rounding mode to round-toward-zero
switch(r) {
case RN: ... round-to-nearest
if(x<= z*(z-ulp)...chopped) z = z - ulp; else
if(x<= z*(z+ulp)...chopped) z = z; else z = z+ulp;
break;
case RZ:case RM: ... round-to-zero or round-to--inf
R:=RP; ... reset rounding mod to round-to-+inf
if(x<z*z ... rounded up) z = z - ulp; else
if(x>=(z+ulp)*(z+ulp) ...rounded up) z = z+ulp;
break;
case RP: ... round-to-+inf
if(x>(z+ulp)*(z+ulp)...chopped) z = z+2*ulp; else
if(x>z*z ...chopped) z = z+ulp;
break;
}
Remark 3. The above comparisons can be done in fixed point. For
example, to compare x and w=z*z chopped, it suffices to compare
x1 and w1 (the trailing parts of x and w), regarding them as
two's complement integers.
...Is z an exact square root?
To determine whether z is an exact square root of x, let z1 be the
trailing part of z, and also let x0 and x1 be the leading and
trailing parts of x.
If ((z1&0x03ffffff)!=0) ... not exact if trailing 26 bits of z!=0
I := 1; ... Raise Inexact flag: z is not exact
else {
j := 1 - [(x0>>20)&1] ... j = logb(x) mod 2
k := z1 >> 26; ... get z's 25-th and 26-th
fraction bits
I := i or (k&j) or ((k&(j+j+1))!=(x1&3));
}
R:= r ... restore rounded mode
return sqrt(x):=z.
If multiplication is cheaper then the foregoing red tape, the
Inexact flag can be evaluated by
I := i;
I := (z*z!=x) or I.
Note that z*z can overwrite I; this value must be sensed if it is
True.
Remark 4. If z*z = x exactly, then bit 25 to bit 0 of z1 must be
zero.
--------------------
z1: | f2 |
--------------------
bit 31 bit 0
Further more, bit 27 and 26 of z1, bit 0 and 1 of x1, and the odd
or even of logb(x) have the following relations:
-------------------------------------------------
bit 27,26 of z1 bit 1,0 of x1 logb(x)
-------------------------------------------------
00 00 odd and even
01 01 even
10 10 odd
10 00 even
11 01 even
-------------------------------------------------
(4) Special cases (see (4) of Section A).
*/
+60
View File
@@ -0,0 +1,60 @@
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/*
* from: @(#)fdlibm.h 5.1 93/09/24
* $FreeBSD$
*/
#ifndef mozilla_imported_fdlibm_h
#define mozilla_imported_fdlibm_h
namespace fdlibm {
double acos(double);
double asin(double);
double atan(double);
double atan2(double, double);
double cosh(double);
double sinh(double);
double tanh(double);
double exp(double);
double log(double);
double log10(double);
double pow(double, double);
double sqrt(double);
double ceil(double);
float ceilf(float);
double fabs(double);
double floor(double);
double acosh(double);
double asinh(double);
double atanh(double);
double cbrt(double);
double expm1(double);
double hypot(double, double);
double log1p(double);
double log2(double);
double copysign(double, double);
double scalbn(double, int);
double trunc(double);
float floorf(float);
} /* namespace fdlibm */
#endif /* mozilla_imported_fdlibm_h */
+81
View File
@@ -0,0 +1,81 @@
/*-
* Copyright (c) 2011 David Schultz <das@FreeBSD.ORG>
* All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions
* are met:
* 1. Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* 2. Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in the
* documentation and/or other materials provided with the distribution.
*
* THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
* ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
* ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
* FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
* DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
* OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
* HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
* LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
* OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
* SUCH DAMAGE.
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
#include "math_private.h"
static const uint32_t k = 1799; /* constant for reduction */
static const double kln2 = 1246.97177782734161156; /* k * ln2 */
/*
* Compute exp(x), scaled to avoid spurious overflow. An exponent is
* returned separately in 'expt'.
*
* Input: ln(DBL_MAX) <= x < ln(2 * DBL_MAX / DBL_MIN_DENORM) ~= 1454.91
* Output: 2**1023 <= y < 2**1024
*/
static double
__frexp_exp(double x, int *expt)
{
double exp_x;
uint32_t hx;
/*
* We use exp(x) = exp(x - kln2) * 2**k, carefully chosen to
* minimize |exp(kln2) - 2**k|. We also scale the exponent of
* exp_x to MAX_EXP so that the result can be multiplied by
* a tiny number without losing accuracy due to denormalization.
*/
exp_x = exp(x - kln2);
GET_HIGH_WORD(hx, exp_x);
*expt = (hx >> 20) - (0x3ff + 1023) + k;
SET_HIGH_WORD(exp_x, (hx & 0xfffff) | ((0x3ff + 1023) << 20));
return (exp_x);
}
/*
* __ldexp_exp(x, expt) and __ldexp_cexp(x, expt) compute exp(x) * 2**expt.
* They are intended for large arguments (real part >= ln(DBL_MAX))
* where care is needed to avoid overflow.
*
* The present implementation is narrowly tailored for our hyperbolic and
* exponential functions. We assume expt is small (0 or -1), and the caller
* has filtered out very large x, for which overflow would be inevitable.
*/
double
__ldexp_exp(double x, int expt)
{
double exp_x, scale;
int ex_expt;
exp_x = __frexp_exp(x, &ex_expt);
expt += ex_expt;
INSERT_WORDS(scale, (0x3ff + expt) << 20, 0);
return (exp_x * scale);
}
+100
View File
@@ -0,0 +1,100 @@
/* @(#)e_log.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
//#include <sys/cdefs.h>
//__FBSDID("$FreeBSD$");
/*
* k_log1p(f):
* Return log(1+f) - f for 1+f in ~[sqrt(2)/2, sqrt(2)].
*
* The following describes the overall strategy for computing
* logarithms in base e. The argument reduction and adding the final
* term of the polynomial are done by the caller for increased accuracy
* when different bases are used.
*
* Method :
* 1. Argument Reduction: find k and f such that
* x = 2^k * (1+f),
* where sqrt(2)/2 < 1+f < sqrt(2) .
*
* 2. Approximation of log(1+f).
* Let s = f/(2+f) ; based on log(1+f) = log(1+s) - log(1-s)
* = 2s + 2/3 s**3 + 2/5 s**5 + .....,
* = 2s + s*R
* We use a special Reme algorithm on [0,0.1716] to generate
* a polynomial of degree 14 to approximate R The maximum error
* of this polynomial approximation is bounded by 2**-58.45. In
* other words,
* 2 4 6 8 10 12 14
* R(z) ~ Lg1*s +Lg2*s +Lg3*s +Lg4*s +Lg5*s +Lg6*s +Lg7*s
* (the values of Lg1 to Lg7 are listed in the program)
* and
* | 2 14 | -58.45
* | Lg1*s +...+Lg7*s - R(z) | <= 2
* | |
* Note that 2s = f - s*f = f - hfsq + s*hfsq, where hfsq = f*f/2.
* In order to guarantee error in log below 1ulp, we compute log
* by
* log(1+f) = f - s*(f - R) (if f is not too large)
* log(1+f) = f - (hfsq - s*(hfsq+R)). (better accuracy)
*
* 3. Finally, log(x) = k*ln2 + log(1+f).
* = k*ln2_hi+(f-(hfsq-(s*(hfsq+R)+k*ln2_lo)))
* Here ln2 is split into two floating point number:
* ln2_hi + ln2_lo,
* where n*ln2_hi is always exact for |n| < 2000.
*
* Special cases:
* log(x) is NaN with signal if x < 0 (including -INF) ;
* log(+INF) is +INF; log(0) is -INF with signal;
* log(NaN) is that NaN with no signal.
*
* Accuracy:
* according to an error analysis, the error is always less than
* 1 ulp (unit in the last place).
*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
static const double
Lg1 = 6.666666666666735130e-01, /* 3FE55555 55555593 */
Lg2 = 3.999999999940941908e-01, /* 3FD99999 9997FA04 */
Lg3 = 2.857142874366239149e-01, /* 3FD24924 94229359 */
Lg4 = 2.222219843214978396e-01, /* 3FCC71C5 1D8E78AF */
Lg5 = 1.818357216161805012e-01, /* 3FC74664 96CB03DE */
Lg6 = 1.531383769920937332e-01, /* 3FC39A09 D078C69F */
Lg7 = 1.479819860511658591e-01; /* 3FC2F112 DF3E5244 */
/*
* We always inline k_log1p(), since doing so produces a
* substantial performance improvement (~40% on amd64).
*/
static inline double
k_log1p(double f)
{
double hfsq,s,z,R,w,t1,t2;
s = f/(2.0+f);
z = s*s;
w = z*z;
t1= w*(Lg2+w*(Lg4+w*Lg6));
t2= z*(Lg1+w*(Lg3+w*(Lg5+w*Lg7)));
R = t2+t1;
hfsq=0.5*f*f;
return s*(hfsq+R);
}

Some files were not shown because too many files have changed in this diff Show More