Initial commit.

This commit is contained in:
Matt Pharr
2011-06-21 06:23:29 -07:00
commit 18af5226ba
587 changed files with 45117 additions and 0 deletions

26
examples/options/Makefile Normal file
View File

@@ -0,0 +1,26 @@
CXX=g++
CXXFLAGS=-Iobjs/ -g -Wall
ISPC=ispc
ISPCFLAGS=-O2 --target=sse4x2
default: options
.PHONY: dirs clean
dirs:
/bin/mkdir -p objs/
clean:
/bin/rm -rf objs *~ options
options: dirs objs/options.o objs/options_serial.o objs/options_ispc.o
$(CXX) $(CXXFLAGS) -o $@ objs/options.o objs/options_ispc.o objs/options_serial.o -lm
objs/%.o: %.cpp
$(CXX) $< $(CXXFLAGS) -c -o $@
objs/options.o: objs/options_ispc.h options_defs.h
objs/%_ispc.h objs/%_ispc.o: %.ispc options_defs.h
$(ISPC) $(ISPCFLAGS) $< -o objs/$*_ispc.o -h objs/$*_ispc.h

View File

@@ -0,0 +1,151 @@
/*
Copyright (c) 2010-2011, Intel Corporation
All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are
met:
* Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
* 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.
* Neither the name of Intel Corporation nor the names of its
contributors may be used to endorse or promote products derived from
this software without specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS 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 COPYRIGHT OWNER
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 <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <assert.h>
#include <math.h>
#include <algorithm>
#ifndef __APPLE__
#include <malloc.h>
#endif // !__APPLE__
using std::max;
#include "options_defs.h"
#include "../timing.h"
#include "options_ispc.h"
using namespace ispc;
// Allocate memory with 64-byte alignment.
float *AllocFloats(int count) {
int size = count * sizeof(float);
#if defined(_WIN32) || defined(_WIN64)
return (float *)_aligned_malloc(size, 64);
#elif defined (__APPLE__)
// Allocate excess memory to ensure an aligned pointer can be returned
void *mem = malloc(size + (64-1) + sizeof(void*));
char *amem = ((char*)mem) + sizeof(void*);
amem += 64 - (reinterpret_cast<uint64_t>(amem) & (64 - 1));
((void**)amem)[-1] = mem;
return (float *)amem;
#else
return (float *)memalign(64, size);
#endif
}
extern void black_scholes_serial(float Sa[], float Xa[], float Ta[],
float ra[], float va[],
float result[], int count);
extern void binomial_put_serial(float Sa[], float Xa[], float Ta[],
float ra[], float va[],
float result[], int count);
int main() {
// Pointers passed to ispc code must have alignment of the target's
// vector width at minimum.
float *S = AllocFloats(N_OPTIONS);
float *X = AllocFloats(N_OPTIONS);
float *T = AllocFloats(N_OPTIONS);
float *r = AllocFloats(N_OPTIONS);
float *v = AllocFloats(N_OPTIONS);
float *result = AllocFloats(N_OPTIONS);
for (int i = 0; i < N_OPTIONS; ++i) {
S[i] = 100; // stock price
X[i] = 98; // option strike price
T[i] = 2; // time (years)
r[i] = .02; // risk-free interest rate
v[i] = 5; // volatility
}
//
// Binomial options pricing model, ispc implementation
//
reset_and_start_timer();
binomial_put_ispc(S, X, T, r, v, result, N_OPTIONS);
double binomial_ispc = get_elapsed_mcycles();
float sum = 0.f;
for (int i = 0; i < N_OPTIONS; ++i)
sum += result[i];
printf("[binomial ispc]:\t\t[%.3f] million cycles (avg %f)\n",
binomial_ispc, sum / N_OPTIONS);
//
// Binomial options, serial implementation
//
reset_and_start_timer();
binomial_put_serial(S, X, T, r, v, result, N_OPTIONS);
double binomial_serial = get_elapsed_mcycles();
sum = 0.f;
for (int i = 0; i < N_OPTIONS; ++i)
sum += result[i];
printf("[binomial serial]:\t\t[%.3f] million cycles (avg %f)\n",
binomial_serial, sum / N_OPTIONS);
printf("\t\t\t\t(%.2fx speedup from ISPC)\n", binomial_serial / binomial_ispc);
//
// Black-Scholes options pricing model, ispc implementation
//
sum = 0.f;
reset_and_start_timer();
for (int a = 0; a < N_BLACK_SCHOLES_ROUNDS; ++a) {
black_scholes_ispc(S, X, T, r, v, result, N_OPTIONS);
for (int i = 0; i < N_OPTIONS; ++i)
sum += result[i];
}
double bs_ispc = get_elapsed_mcycles();
printf("[black-scholes ispc]:\t\t[%.3f] million cycles (avg %f)\n",
bs_ispc, sum / (N_BLACK_SCHOLES_ROUNDS * N_OPTIONS));
//
// Black-Scholes options pricing model, serial implementation
//
sum = 0.f;
reset_and_start_timer();
for (int a = 0; a < N_BLACK_SCHOLES_ROUNDS; ++a) {
black_scholes_serial(S, X, T, r, v, result, N_OPTIONS);
for (int i = 0; i < N_OPTIONS; ++i)
sum += result[i];
}
double bs_serial = get_elapsed_mcycles();
printf("[black-scholes serial]:\t\t[%.3f] million cycles (avg %f)\n", bs_serial,
sum / (N_BLACK_SCHOLES_ROUNDS * N_OPTIONS));
printf("\t\t\t\t(%.2fx speedup from ISPC)\n", bs_serial / bs_ispc);
return 0;
}

View File

@@ -0,0 +1,103 @@
// -*- mode: c++ -*-
/*
Copyright (c) 2010-2011, Intel Corporation
All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are
met:
* Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
* 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.
* Neither the name of Intel Corporation nor the names of its
contributors may be used to endorse or promote products derived from
this software without specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS 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 COPYRIGHT OWNER
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 "options_defs.h"
// Cumulative normal distribution function
static inline float
CND(float X) {
float L = abs(X);
float k = 1.0 / (1.0 + 0.2316419 * L);
float k2 = k*k;
float k3 = k2*k;
float k4 = k2*k2;
float k5 = k3*k2;
const float invSqrt2Pi = 0.39894228040f;
float w = (0.31938153f * k - 0.356563782f * k2 + 1.781477937f * k3 +
-1.821255978f * k4 + 1.330274429f * k5);
w *= invSqrt2Pi * exp(-L * L * .5f);
if (X > 0.f)
w = 1.0 - w;
return w;
}
export void
black_scholes_ispc(uniform float Sa[], uniform float Xa[], uniform float Ta[],
uniform float ra[], uniform float va[],
uniform float result[], uniform int count) {
for (uniform int i = 0; i < count; i += programCount) {
float S = Sa[i + programIndex], X = Xa[i + programIndex];
float T = Ta[i + programIndex], r = ra[i + programIndex];
float v = va[i + programIndex];
float d1 = (log(S/X) + (r + v * v * .5f) * T) / (v * sqrt(T));
float d2 = d1 - v * sqrt(T);
result[i + programIndex] = S * CND(d1) - X * exp(-r * T) * CND(d2);
}
}
export void
binomial_put_ispc(uniform float Sa[], uniform float Xa[], uniform float Ta[],
uniform float ra[], uniform float va[],
uniform float result[], uniform int count) {
float V[BINOMIAL_NUM];
for (uniform int i = 0; i < count; i += programCount) {
float S = Sa[i + programIndex], X = Xa[i + programIndex];
float T = Ta[i + programIndex], r = ra[i + programIndex];
float v = va[i + programIndex];
float dt = T / BINOMIAL_NUM;
float u = exp(v * sqrt(dt));
float d = 1. / u;
float disc = exp(r * dt);
float Pu = (disc - d) / (u - d);
for (uniform int j = 0; j < BINOMIAL_NUM; ++j) {
float upow = pow(u, (float)(2*j-BINOMIAL_NUM));
V[j] = max(0., X - S * upow);
}
for (uniform int j = BINOMIAL_NUM-1; j >= 0; --j)
for (uniform int k = 0; k < j; ++k)
V[k] = ((1 - Pu) * V[k] + Pu * V[k + 1]) / disc;
result[i + programIndex] = V[0];
}
}

168
examples/options/options.vcxproj Executable file
View File

@@ -0,0 +1,168 @@
<?xml version="1.0" encoding="utf-8"?>
<Project DefaultTargets="Build" ToolsVersion="4.0" xmlns="http://schemas.microsoft.com/developer/msbuild/2003">
<ItemGroup Label="ProjectConfigurations">
<ProjectConfiguration Include="Debug|Win32">
<Configuration>Debug</Configuration>
<Platform>Win32</Platform>
</ProjectConfiguration>
<ProjectConfiguration Include="Debug|x64">
<Configuration>Debug</Configuration>
<Platform>x64</Platform>
</ProjectConfiguration>
<ProjectConfiguration Include="Release|Win32">
<Configuration>Release</Configuration>
<Platform>Win32</Platform>
</ProjectConfiguration>
<ProjectConfiguration Include="Release|x64">
<Configuration>Release</Configuration>
<Platform>x64</Platform>
</ProjectConfiguration>
</ItemGroup>
<PropertyGroup Label="Globals">
<ProjectGuid>{8C7B5D29-1E76-44E6-BBB8-09830E5DEEAE}</ProjectGuid>
<Keyword>Win32Proj</Keyword>
<RootNamespace>options</RootNamespace>
</PropertyGroup>
<Import Project="$(VCTargetsPath)\Microsoft.Cpp.Default.props" />
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'" Label="Configuration">
<ConfigurationType>Application</ConfigurationType>
<UseDebugLibraries>true</UseDebugLibraries>
<CharacterSet>Unicode</CharacterSet>
</PropertyGroup>
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Debug|x64'" Label="Configuration">
<ConfigurationType>Application</ConfigurationType>
<UseDebugLibraries>true</UseDebugLibraries>
<CharacterSet>Unicode</CharacterSet>
</PropertyGroup>
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Release|Win32'" Label="Configuration">
<ConfigurationType>Application</ConfigurationType>
<UseDebugLibraries>false</UseDebugLibraries>
<WholeProgramOptimization>true</WholeProgramOptimization>
<CharacterSet>Unicode</CharacterSet>
</PropertyGroup>
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Release|x64'" Label="Configuration">
<ConfigurationType>Application</ConfigurationType>
<UseDebugLibraries>false</UseDebugLibraries>
<WholeProgramOptimization>true</WholeProgramOptimization>
<CharacterSet>Unicode</CharacterSet>
</PropertyGroup>
<Import Project="$(VCTargetsPath)\Microsoft.Cpp.props" />
<ImportGroup Label="ExtensionSettings">
</ImportGroup>
<ImportGroup Label="PropertySheets" Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">
<Import Project="$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props" Condition="exists('$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props')" Label="LocalAppDataPlatform" />
</ImportGroup>
<ImportGroup Condition="'$(Configuration)|$(Platform)'=='Debug|x64'" Label="PropertySheets">
<Import Project="$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props" Condition="exists('$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props')" Label="LocalAppDataPlatform" />
</ImportGroup>
<ImportGroup Label="PropertySheets" Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">
<Import Project="$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props" Condition="exists('$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props')" Label="LocalAppDataPlatform" />
</ImportGroup>
<ImportGroup Condition="'$(Configuration)|$(Platform)'=='Release|x64'" Label="PropertySheets">
<Import Project="$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props" Condition="exists('$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props')" Label="LocalAppDataPlatform" />
</ImportGroup>
<PropertyGroup Label="UserMacros" />
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">
<LinkIncremental>true</LinkIncremental>
</PropertyGroup>
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Debug|x64'">
<LinkIncremental>true</LinkIncremental>
</PropertyGroup>
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">
<LinkIncremental>false</LinkIncremental>
</PropertyGroup>
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Release|x64'">
<LinkIncremental>false</LinkIncremental>
</PropertyGroup>
<ItemDefinitionGroup Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">
<ClCompile>
<PrecompiledHeader>
</PrecompiledHeader>
<WarningLevel>Level3</WarningLevel>
<Optimization>Disabled</Optimization>
<PreprocessorDefinitions>WIN32;_DEBUG;_CONSOLE;%(PreprocessorDefinitions)</PreprocessorDefinitions>
<DisableSpecificWarnings>4305</DisableSpecificWarnings>
</ClCompile>
<Link>
<SubSystem>Console</SubSystem>
<GenerateDebugInformation>true</GenerateDebugInformation>
</Link>
</ItemDefinitionGroup>
<ItemDefinitionGroup Condition="'$(Configuration)|$(Platform)'=='Debug|x64'">
<ClCompile>
<PrecompiledHeader>
</PrecompiledHeader>
<WarningLevel>Level3</WarningLevel>
<Optimization>Disabled</Optimization>
<PreprocessorDefinitions>WIN32;_DEBUG;_CONSOLE;%(PreprocessorDefinitions)</PreprocessorDefinitions>
<DisableSpecificWarnings>4305</DisableSpecificWarnings>
</ClCompile>
<Link>
<SubSystem>Console</SubSystem>
<GenerateDebugInformation>true</GenerateDebugInformation>
</Link>
</ItemDefinitionGroup>
<ItemDefinitionGroup Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">
<ClCompile>
<WarningLevel>Level3</WarningLevel>
<PrecompiledHeader>
</PrecompiledHeader>
<Optimization>MaxSpeed</Optimization>
<FunctionLevelLinking>true</FunctionLevelLinking>
<IntrinsicFunctions>true</IntrinsicFunctions>
<PreprocessorDefinitions>WIN32;NDEBUG;_CONSOLE;%(PreprocessorDefinitions)</PreprocessorDefinitions>
<DisableSpecificWarnings>4305</DisableSpecificWarnings>
</ClCompile>
<Link>
<SubSystem>Console</SubSystem>
<GenerateDebugInformation>true</GenerateDebugInformation>
<EnableCOMDATFolding>true</EnableCOMDATFolding>
<OptimizeReferences>true</OptimizeReferences>
</Link>
</ItemDefinitionGroup>
<ItemDefinitionGroup Condition="'$(Configuration)|$(Platform)'=='Release|x64'">
<ClCompile>
<WarningLevel>Level3</WarningLevel>
<PrecompiledHeader>
</PrecompiledHeader>
<Optimization>MaxSpeed</Optimization>
<FunctionLevelLinking>true</FunctionLevelLinking>
<IntrinsicFunctions>true</IntrinsicFunctions>
<PreprocessorDefinitions>WIN32;NDEBUG;_CONSOLE;%(PreprocessorDefinitions)</PreprocessorDefinitions>
<DisableSpecificWarnings>4305</DisableSpecificWarnings>
</ClCompile>
<Link>
<SubSystem>Console</SubSystem>
<GenerateDebugInformation>true</GenerateDebugInformation>
<EnableCOMDATFolding>true</EnableCOMDATFolding>
<OptimizeReferences>true</OptimizeReferences>
</Link>
</ItemDefinitionGroup>
<ItemGroup>
<ClCompile Include="options.cpp" />
<ClCompile Include="options_serial.cpp" />
</ItemGroup>
<ItemGroup>
<CustomBuild Include="options.ispc">
<FileType>Document</FileType>
<Command Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">cl /E /TP %(Filename).ispc | ispc -O2 - -o %(Filename).obj -h %(Filename)_ispc.h --arch=x86 --target=sse4x2
</Command>
<Command Condition="'$(Configuration)|$(Platform)'=='Debug|x64'">cl /E /TP %(Filename).ispc | ispc -O2 - -o %(Filename).obj -h %(Filename)_ispc.h --target=sse4x2
</Command>
<Outputs Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">%(Filename).obj;%(Filename)_ispc.h</Outputs>
<Outputs Condition="'$(Configuration)|$(Platform)'=='Debug|x64'">%(Filename).obj;%(Filename)_ispc.h</Outputs>
<Command Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">cl /E /TP %(Filename).ispc | ispc -O2 - -o %(Filename).obj -h %(Filename)_ispc.h --arch=x86 --target=sse4x2
</Command>
<Command Condition="'$(Configuration)|$(Platform)'=='Release|x64'">cl /E /TP %(Filename).ispc | ispc -O2 - -o %(Filename).obj -h %(Filename)_ispc.h --target=sse4x2
</Command>
<Outputs Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">%(Filename).obj;%(Filename)_ispc.h</Outputs>
<Outputs Condition="'$(Configuration)|$(Platform)'=='Release|x64'">%(Filename).obj;%(Filename)_ispc.h</Outputs>
</CustomBuild>
</ItemGroup>
<ItemGroup>
<ClInclude Include="options_defs.h" />
</ItemGroup>
<Import Project="$(VCTargetsPath)\Microsoft.Cpp.targets" />
<ImportGroup Label="ExtensionTargets">
</ImportGroup>
</Project>

View File

@@ -0,0 +1,42 @@
/*
Copyright (c) 2010-2011, Intel Corporation
All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are
met:
* Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
* 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.
* Neither the name of Intel Corporation nor the names of its
contributors may be used to endorse or promote products derived from
this software without specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS 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 COPYRIGHT OWNER
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.
*/
#ifndef OPTIONS_DEFS_H
#define OPTIONS_DEFS_H 1
#define BINOMIAL_NUM 64
#define N_OPTIONS 65536
#define N_BLACK_SCHOLES_ROUNDS 20
#endif // OPTIONS_DEFS_H

View File

@@ -0,0 +1,114 @@
/*
Copyright (c) 2010-2011, Intel Corporation
All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are
met:
* Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
* 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.
* Neither the name of Intel Corporation nor the names of its
contributors may be used to endorse or promote products derived from
this software without specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS 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 COPYRIGHT OWNER
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.
*/
#ifdef _MSC_VER
#define _CRT_SECURE_NO_WARNINGS
#define NOMINMAX
#pragma warning (disable: 4244)
#pragma warning (disable: 4305)
#endif
#include "options_defs.h"
#include <math.h>
#include <algorithm>
// Cumulative normal distribution function
static inline float
CND(float X) {
float L = fabsf(X);
float k = 1.0 / (1.0 + 0.2316419 * L);
float k2 = k*k;
float k3 = k2*k;
float k4 = k2*k2;
float k5 = k3*k2;
const float invSqrt2Pi = 0.39894228040f;
float w = (0.31938153f * k - 0.356563782f * k2 + 1.781477937f * k3 +
-1.821255978f * k4 + 1.330274429f * k5);
w *= invSqrt2Pi * expf(-L * L * .5f);
if (X > 0.f)
w = 1.0 - w;
return w;
}
void
black_scholes_serial(float Sa[], float Xa[], float Ta[],
float ra[], float va[],
float result[], int count) {
for (int i = 0; i < count; ++i) {
float S = Sa[i], X = Xa[i];
float T = Ta[i], r = ra[i];
float v = va[i];
float d1 = (logf(S/X) + (r + v * v * .5f) * T) / (v * sqrtf(T));
float d2 = d1 - v * sqrtf(T);
result[i] = S * CND(d1) - X * expf(-r * T) * CND(d2);
}
}
void
binomial_put_serial(float Sa[], float Xa[], float Ta[],
float ra[], float va[],
float result[], int count) {
float V[BINOMIAL_NUM];
for (int i = 0; i < count; ++i) {
float S = Sa[i], X = Xa[i];
float T = Ta[i], r = ra[i];
float v = va[i];
float dt = T / BINOMIAL_NUM;
float u = expf(v * sqrtf(dt));
float d = 1. / u;
float disc = expf(r * dt);
float Pu = (disc - d) / (u - d);
for (int j = 0; j < BINOMIAL_NUM; ++j) {
float upow = powf(u, (float)(2*j-BINOMIAL_NUM));
V[j] = std::max(0.f, X - S * upow);
}
for (int j = BINOMIAL_NUM-1; j >= 0; --j)
for (int k = 0; k < j; ++k)
V[k] = ((1 - Pu) * V[k] + Pu * V[k + 1]) / disc;
result[i] = V[0];
}
}