diff --git a/.gitignore b/.gitignore index 1bf4bd66e386359aef0b89f2338533a28bf35c0c..609777ae3ae95ca9fbb815beab4a9f2bbab562aa 100644 --- a/.gitignore +++ b/.gitignore @@ -35,6 +35,7 @@ .vscode/ .history/ build/ +srcGen/ main *.txt __pycache__/ diff --git a/Makefile b/Makefile index 0a535438f5d6a48d4527faa530325b7946a91d54..291ac1ead82e5a63e404e1ce81cbc65c65d90851 100644 --- a/Makefile +++ b/Makefile @@ -4,11 +4,11 @@ CC = gcc CPP = g++ INCLUDE = -Iinclude -I/usr/include/python3.8 # LIBS= -L/usr/lib/python3.8/config-3.8-x86_64-linux-gnu -lpython3.8 # may need -L to assign the Python lobrary path -LIBS = -lpython3.8 -lfmt +LIBS = -lm -lpython3.8 -lfmt ECHO = printf # $(info $(CFLAGS) ) -override CFLAGS += -g -Wall -Wextra -Wpedantic -Wno-unused-function -fdiagnostics-color=always +override CFLAGS += -g -Wall -Wextra -Wpedantic -Wno-unused-function -fdiagnostics-color=always --std=c++2a # $(info $(CFLAGS) ) EXPRAUTO_ALL_SRCS_CPP = $(wildcard src/*.cpp) diff --git a/benchMarkInterval.txt b/benchMarkInterval.txt index 16b257f7a79a6e7e2478854444946e01c875c3ab..0715ee06c51ecf68ca8cc7ca826f1e55e5423968 100644 --- a/benchMarkInterval.txt +++ b/benchMarkInterval.txt @@ -4,7 +4,7 @@ 0 999 -8 8 0.01 100 --1 1 +-0.5 0.5 0.01 3 0.01 100 0.01 100 diff --git a/benchMarkThreshold.txt b/benchMarkThreshold.txt new file mode 100644 index 0000000000000000000000000000000000000000..cd0b7e420a5c923ace6a3b89eae6d2e62100a5a4 --- /dev/null +++ b/benchMarkThreshold.txt @@ -0,0 +1,44 @@ +0.8 +1 +1 +1 +10 +2 +2 +2.1 +1 +2 +0.5 +10 +5 +1 +1.5 +200 +1 +1 +2 +1 +2 +2 +2 +2 +1 +3 +0.6 +1.5 +1.2 +1.2 +2 +350,4 +1 +2.6 +2.5 +2.8 +3 +3 +3.2 +2.6 +2.6 +2.6 +2 +2 diff --git a/detectErrorOneFPEDParallel.sh b/detectErrorOneFPEDParallel.sh index d203a6fc044092aff39b8da7732e69e87f5cd07b..5d69f167d7392863f81c2ec572b67f043f8a44d4 100755 --- a/detectErrorOneFPEDParallel.sh +++ b/detectErrorOneFPEDParallel.sh @@ -25,13 +25,23 @@ fi testFileName=test1paramFPEDParallel numProcs=32 -echo "Detecting error: ${uniqueLabel} ${x0Start} ${x0End} ${x0Size} ${prefix} ${middle} ${suffix}" +# echo "Detecting error: ${uniqueLabel} ${x0Start} ${x0End} ${x0Size} ${prefix} ${middle} ${suffix}" directory="./srcTest"/${uniqueLabel} # echo "${CC} ${testFileName}.c ${prefix}_${suffix}.c ${prefix}_mpfr.c computeULP.c -IincludeTEST -DEXPRESSION=${prefix}_ -DSUFFIX=${suffix} -lmpfr -lm -O3 -o ${testFileName}.exe" -${CC} ./srcTest/${testFileName}.c ${directory}/${prefix}_${suffix}.c ${directory}/${prefix}_mpfr.c ./srcTest/computeULP.c -IincludeTEST -DEXPRESSION=${prefix}_ -DSUFFIX=${suffix} -lmpfr -lm -o ${testFileName}.exe +${CC} ./srcTest/${testFileName}.c ${directory}/${prefix}_${suffix}.c ${directory}/${prefix}_mpfr.c ./srcTest/computeULP.c -IincludeTEST -IincludeDD -DEXPRESSION=${prefix}_ -DSUFFIX=${suffix} -Llibs -lTGen -lmpfr -lm -lqd -o ${testFileName}.exe # echo "mpirun -n ${numProcs} ./${testFileName}.exe ${x0Start} ${x0End} ${x0Size} ${prefix}__${middle}_${suffix}" mpirun -n ${numProcs} ./${testFileName}.exe ${x0Start} ${x0End} ${x0Size} ${prefix}__${middle}_${suffix} ${uniqueLabel} # mv outputs/${prefix}__${middle}_${suffix}_error.txt ./outputs/${uniqueLabel}/${prefix}__${middle}_${suffix}_error.txt rm ${testFileName}.exe -echo "end detecting ${uniqueLabel}" + +# combine files +cd ./outputs/${uniqueLabel} +findWord="${prefix}__${middle}_${suffix}_sample_*.txt" +# echo "For suffix = ${suffix}, Find and combine by shell command cat: ${findWord}" +find . -name "${findWord}" | sort -h | xargs cat > sample_${uniqueLabel}_${suffix}.txt +# echo "sample file: `pwd`/sample_${uniqueLabel}_${suffix}.txt" +rm ${findWord} +cd - > /dev/null + +echo "end detecting error ${uniqueLabel}" echo diff --git a/detectErrorThreeFPEDParallel.sh b/detectErrorThreeFPEDParallel.sh index 313c9ba2a264fcaede737c3b568fd2c258a4dfa3..d90d24615bb6bfb981b19cc4652b70ac1507e2a3 100755 --- a/detectErrorThreeFPEDParallel.sh +++ b/detectErrorThreeFPEDParallel.sh @@ -36,7 +36,7 @@ numProcs=32 echo "Detecting error: ${uniqueLabel} ${x0Start} ${x0End} ${x1Start} ${x1End} ${x2Start} ${x2End} ${x0Size} ${x1Size} ${x2Size} ${prefix} ${middle} ${suffix}" directory="./srcTest"/${uniqueLabel} # echo "${CC} ${testFileName}.c ${prefix}_${suffix}.c ${prefix}_mpfr.c computeULP.c -IincludeTEST -DEXPRESSION=${prefix}_ -DSUFFIX=${suffix} -lmpfr -lm -O3 -o ${testFileName}.exe" -${CC} ./srcTest/${testFileName}.c ${directory}/${prefix}_${suffix}.c ${directory}/${prefix}_mpfr.c ./srcTest/computeULP.c -IincludeTEST -DEXPRESSION=${prefix}_ -DSUFFIX=${suffix} -lmpfr -lm -o ${testFileName}.exe +${CC} ./srcTest/${testFileName}.c ${directory}/${prefix}_${suffix}.c ${directory}/${prefix}_mpfr.c ./srcTest/computeULP.c -IincludeTEST -IincludeDD -DEXPRESSION=${prefix}_ -DSUFFIX=${suffix} -lmpfr -lm -lqd -o ${testFileName}.exe # echo "mpirun -n ${numProcs} ./${testFileName}.exe ${x0Start} ${x0End} ${x1Start} ${x1End} ${x2Start} ${x2End} ${x0Size} ${x1Size} ${x2Size} ${prefix}__${middle}_${suffix}" mpirun -n ${numProcs} ./${testFileName}.exe ${x0Start} ${x0End} ${x1Start} ${x1End} ${x2Start} ${x2End} ${x0Size} ${x1Size} ${x2Size} ${prefix}__${middle}_${suffix} ${uniqueLabel} rm ${testFileName}.exe diff --git a/detectErrorTwoFPEDParallel.sh b/detectErrorTwoFPEDParallel.sh index 0828a8a48b1f0a2776696e4e8961568b4a69bf51..dcb50691049e8411c0aeacca75ff33bd32fe5d06 100755 --- a/detectErrorTwoFPEDParallel.sh +++ b/detectErrorTwoFPEDParallel.sh @@ -34,7 +34,7 @@ numProcs=32 echo "Detecting error: ${uniqueLabel} ${x0Start} ${x0End} ${x1Start} ${x1End} ${x0Size} ${x1Size} ${prefix} ${middle} ${suffix}" directory="./srcTest"/${uniqueLabel} # echo "${CC} ${testFileName}.c ${prefix}_${suffix}.c ${prefix}_mpfr.c computeULP.c -IincludeTEST -DEXPRESSION=${prefix}_ -DSUFFIX=${suffix} -lmpfr -lm -O3 -o ${testFileName}.exe" -${CC} ./srcTest/${testFileName}.c ${directory}/${prefix}_${suffix}.c ${directory}/${prefix}_mpfr.c ./srcTest/computeULP.c -IincludeTEST -DEXPRESSION=${prefix}_ -DSUFFIX=${suffix} -lmpfr -lm -o ${testFileName}.exe +${CC} ./srcTest/${testFileName}.c ${directory}/${prefix}_${suffix}.c ${directory}/${prefix}_mpfr.c ./srcTest/computeULP.c -IincludeTEST -IincludeDD -DEXPRESSION=${prefix}_ -DSUFFIX=${suffix} -lmpfr -lm -lqd -o ${testFileName}.exe # echo "mpirun -n ${numProcs} ./${testFileName}.exe ${x0Start} ${x0End} ${x1Start} ${x1End} ${x2Start} ${x2End} ${x0Size} ${x1Size} ${x2Size} ${prefix}__${middle}_${suffix}" mpirun -n ${numProcs} ./${testFileName}.exe ${x0Start} ${x0End} ${x1Start} ${x1End} ${x0Size} ${x1Size} ${prefix}__${middle}_${suffix} ${uniqueLabel} diff --git a/include/basic.hpp b/include/basic.hpp index fbe60ea0812f2409ce4ba416fa7fdf101833b19c..e665a42ba82480882a77dfa8a89d39590c01972d 100644 --- a/include/basic.hpp +++ b/include/basic.hpp @@ -36,6 +36,8 @@ extern const size_t promtTimes; /// ExprAST - Base class for all expression nodes. class ExprAST { + int order = 0; + string opType = "double"; public: virtual ~ExprAST() = default; @@ -44,6 +46,14 @@ public: virtual string type() { return "Base"; } virtual std::unique_ptr Clone() { return makePtr(); } + + virtual int getOrder() { return order; } + + virtual void setOrder(int input) { order = input; } + + virtual string getOpType() { return opType; } + + virtual void setOpType(string input) { opType = input; } }; /// NumberExprAST - Expression class for numeric literals like "1.0". diff --git a/include/devideUpEdgeData.hpp b/include/devideUpEdgeData.hpp new file mode 100644 index 0000000000000000000000000000000000000000..2f3b48673e6d81a3ccc140cf4ea2ffdcda381396 --- /dev/null +++ b/include/devideUpEdgeData.hpp @@ -0,0 +1,5 @@ +#ifndef _DEVIDEUPEDGEDATA +#define _DEVIDEUPEDGEDATA +// void devideUpEdgeData(string upEdgeFileName, string intervalFileName, double threshold = 2, double thresholdCombine = 101); +vector devideUpEdgeData(string upEdgeFileName, double threshold = 2, double thresholdCombine = 0); +#endif \ No newline at end of file diff --git a/include/exprAuto.hpp b/include/exprAuto.hpp index 387098581d8739eab1443de4bfe996486167fbb3..2a9208a8b16e45cfe8409fccd38d7c6d81b870a4 100644 --- a/include/exprAuto.hpp +++ b/include/exprAuto.hpp @@ -31,7 +31,7 @@ vector dealWithBinOpKernel(const ast_ptr &expr); vector dealWithBinOp(vector &exprs, const char &op); -// TODO: delete the same expr in exprs using the function 'isEqual' +// TODO: (not sure) delete the same expr in exprs using the function 'isEqual' void deleteTheSame(vector &exprs); vector dealWithCalls(const ast_ptr &expr); @@ -44,6 +44,8 @@ void sortExpr(ast_ptr &expr); vector tryRewrite(ast_ptr expr, bool addSelf = true); +string pickTheBest(string uniqueLabel, vector testSet, vector intervals, vector scales); + void geneSampleData(); vector createAll(vector &numerators, vector &denominators); diff --git a/include/geneCode.hpp b/include/geneCode.hpp index 0a44bf0f259ddc5cfe7403d66efbd8fc10a72711..39e399e7a507aa3b9c69e156a72bed4f35dab9f9 100644 --- a/include/geneCode.hpp +++ b/include/geneCode.hpp @@ -19,11 +19,15 @@ vector getVariablesFromExpr(const ast_ptr &expr); bool getVariablesFromExpr(const ast_ptr &expr, vector &vars); -string geneOriginCodeKernel(string exprStr, vector vars, string uniqueLabel,string tail); +string geneExprCodeKernel(string exprStr, vector vars, string uniqueLabel,string tail); -string geneOriginCode(string exprStr, string uniqueLabel, string tail); +string geneExprCode(string exprStr, string uniqueLabel, string tail); -void geneHerbieCode(string exprstr, vector cs, string exprname,double v[], double u[]); +string geneTGenCode(string exprStr, vector vars, string uniqueLabel, string tail); + +void geneHerbieCode(string exprstr, vector cs, string exprname, double v[], double u[]); + +string geneHerbieCode(string uniqueLabel); void geneDaisyCode(string exprStr); @@ -35,4 +39,8 @@ string geneFinalCodeKernel(string exprStr, string uniqueLabel, vector string geneFinalCode(string exprStr, string uniqueLabel, vector exprInfoVector); +void getDepth(ast_ptr &expr, int &depth); + +int codegenWrapper(ast_ptr &expr, vector &vars, const string uniqueLabel, string tail); + #endif \ No newline at end of file diff --git a/include/tools.hpp b/include/tools.hpp index 1159fdcc06c92ade20242a4c3bf19c842c6a9605..6175c2f2b2f39f9f0d7fa67421217af78bee0124 100644 --- a/include/tools.hpp +++ b/include/tools.hpp @@ -35,16 +35,18 @@ exprInfo testError(string uniqueLabel, string suffix, double x0Start, double x0E exprInfo testError(string uniqueLabel, string suffix, double x0Start, double x0End, double x1Start, double x1End, double x2Start, double x2End, int x0Size, int x1Size, int x2Size); -void geneBoundaryData(string uniqueLabel, string suffix); +double testPerformance(string uniqueLabel, string suffix, const vector &intervals); -void geneIntervalData(string uniqeuLabel, vector &intervals, vector &threholds); +string geneBoundaryData(string uniqueLabel, string suffix); -vector getIntervalData(); +vector geneIntervalData(vector upEdgeFileNames, string uniqueLabel, vector &thresholds); vector> getIntervalData(string filename); -vector rewrite(string exprSr, string uniqueLabel); +vector> getIntervalData(vector upEdgeFileNames, vector &thresholds); -vector> zuhe(vector> vec); +vector rewrite(string exprStr, string uniqueLabel, vector> &intervalData); + +vector> permuteMultiVec(vector> vec); #endif \ No newline at end of file diff --git a/includeDD/dd.h b/includeDD/dd.h new file mode 100644 index 0000000000000000000000000000000000000000..91a1c94be9e1941f4abc0c882de6ba01bdfbed1c --- /dev/null +++ b/includeDD/dd.h @@ -0,0 +1,130 @@ +#ifndef _MULTITYPE +#define _MULTITYPE + +#define xTmp(x) x##Tmp + +// =========================================== ADD ======================================= +#define add_d_d_d(a, b, c) double c = (a) + (b); +#define add_ld_d_d(a, b, c) long double xTmp(c) = (a) + (b); double c = xTmp(c); +#define add_dd_d_d(a, b, c) double xTmp(c)[2]; c_dd_add_dd_d(a, b, xTmp(c)); double c = xTmp(c)[0]; +#define add_d_dd_d(a, b, c) double xTmp(c)[2]; c_dd_add_d_dd(a, b, xTmp(c)); double c = xTmp(c)[0]; +#define add_ld_dd_d(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double xTmp(c)[2]; c_dd_add(xTmp(a), b, xTmp(c)); c = xTmp(c)[0]; // check c = ((a) + (b[0])) + b[1]; +#define add_dd_dd_d(a, b, c) double xTmp(c)[2]; c_dd_add(a, b, xTmp(c)); double c; c = xTmp(c)[0]; +#define add_d_ld_d(a, b, c) double c = (a) + (b); +#define add_ld_ld_d(a, b, c) double c = (a) + (b); +#define add_dd_ld_d(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double xTmp(c)[2]; c_dd_add(a, xTmp(b), xTmp(c)); c = xTmp(c)[0]; // check double c = ((a[0]) + (b[0])) + a[1]; + +#define add_d_d_dd(a, b, c) double c[2]; c[0] = two_add_mine(a, b, &(c[1])); // check c[0] = (a) + (b) +#define add_ld_d_dd(a, b, c) long double xTmp(c) = (a) + (b); double c[2]; c[0] = xTmp(c); c[1] = xTmp(c) - c[0]; +#define add_dd_d_dd(a, b, c) double c[2]; c_dd_add_dd_d(a, b, c); +#define add_d_dd_dd(a, b, c) double c[2]; c_dd_add_d_dd(a, b, c); +#define add_ld_dd_dd(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double c[2]; c_dd_add(xTmp(a), b, c); +#define add_dd_dd_dd(a, b, c) double c[2]; c_dd_add(a, b, c); +#define add_d_ld_dd(a, b, c) long double xTmp(c) = (a) + (b); double c[2]; c[0] = xTmp(c); c[1] = c - xTmp(c)[0]; +#define add_ld_ld_dd(a, b, c) long double xTmp(c) = (a) + (b); double c[2]; c[0] = xTmp(c); c[1] = c - xTmp(c)[0]; +#define add_dd_ld_dd(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double c[2]; c_dd_add(a, xTmp(b), c); + +#define add_d_d_ld(a, b, c) long double c = (a) + (b); +#define add_ld_d_ld(a, b, c) long double c = (a) + (b); +#define add_dd_d_ld(a, b, c) double xTmp(c)[2]; c_dd_add_dd_d(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define add_d_dd_ld(a, b, c) double xTmp(c)[2]; c_dd_add_d_dd(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define add_dd_dd_ld(a, b, c) double xTmp(c)[2]; c_dd_add(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define add_ld_dd_ld(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double xTmp(c)[2]; c_dd_add(xTmp(a), b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define add_d_ld_ld(a, b, c) long double c; c = (a) + (b); +#define add_ld_ld_ld(a, b, c) long double c; c = (a) + (b); +#define add_dd_ld_ld(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double xTmp(c)[2]; c_dd_add(a, xTmp(b), xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; + +// =========================================== SUB ======================================= +#define sub_d_d_d(a, b, c) double c = (a) - (b); +#define sub_ld_d_d(a, b, c) long double xTmp(c) = (a) - (b); double c = xTmp(c); +#define sub_dd_d_d(a, b, c) double xTmp(c)[2]; c_dd_sub_dd_d(a, b, xTmp(c)); double c = xTmp(c)[0]; +#define sub_d_dd_d(a, b, c) double xTmp(c)[2]; c_dd_sub_d_dd(a, b, xTmp(c)); double c = xTmp(c)[0]; +#define sub_ld_dd_d(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double xTmp(c)[2]; c_dd_sub(xTmp(a), b, xTmp(c)); c = xTmp(c)[0]; // check c = ((a) - (b[0])) - b[1]; +#define sub_dd_dd_d(a, b, c) double xTmp(c)[2]; c_dd_sub(a, b, xTmp(c)); double c; c = xTmp(c)[0]; +#define sub_d_ld_d(a, b, c) double c = (a) - (b); +#define sub_ld_ld_d(a, b, c) double c = (a) - (b); +#define sub_dd_ld_d(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double xTmp(c)[2]; c_dd_sub(a, xTmp(b), xTmp(c)); c = xTmp(c)[0]; // check double c = ((a[0]) - (b[0])) - a[1]; + +#define sub_d_d_dd(a, b, c) double c[2]; c[0] = two_sub_mine(a, b, &(c[1])); // check c[0] = (a) - (b) +#define sub_ld_d_dd(a, b, c) long double xTmp(c) = (a) - (b); double c[2]; c[0] = xTmp(c); c[1] = xTmp(c) - c[0]; +#define sub_dd_d_dd(a, b, c) double c[2]; c_dd_sub_dd_d(a, b, c); +#define sub_d_dd_dd(a, b, c) double c[2]; c_dd_sub_d_dd(a, b, c); +#define sub_ld_dd_dd(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double c[2]; c_dd_sub(xTmp(a), b, c); +#define sub_dd_dd_dd(a, b, c) double c[2]; c_dd_sub(a, b, c); +#define sub_d_ld_dd(a, b, c) long double xTmp(c) = (a) - (b); double c[2]; c[0] = xTmp(c); c[1] = c - xTmp(c)[0]; +#define sub_ld_ld_dd(a, b, c) long double xTmp(c) = (a) - (b); double c[2]; c[0] = xTmp(c); c[1] = c - xTmp(c)[0]; +#define sub_dd_ld_dd(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double c[2]; c_dd_sub(a, xTmp(b), c); + +#define sub_d_d_ld(a, b, c) long double c = (a) - (b); +#define sub_ld_d_ld(a, b, c) long double c = (a) - (b); +#define sub_dd_d_ld(a, b, c) double xTmp(c)[2]; c_dd_sub_dd_d(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define sub_d_dd_ld(a, b, c) double xTmp(c)[2]; c_dd_sub_d_dd(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define sub_dd_dd_ld(a, b, c) double xTmp(c)[2]; c_dd_sub(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define sub_ld_dd_ld(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double xTmp(c)[2]; c_dd_sub(xTmp(a), b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define sub_d_ld_ld(a, b, c) long double c; c = (a) - (b); +#define sub_ld_ld_ld(a, b, c) long double c; c = (a) - (b); +#define sub_dd_ld_ld(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double xTmp(c)[2]; c_dd_sub(a, xTmp(b), xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; + +// =========================================== MUL ======================================= +#define mul_d_d_d(a, b, c) double c = (a) * (b); +#define mul_ld_d_d(a, b, c) long double xTmp(c) = (a) * (b); double c = xTmp(c); +#define mul_dd_d_d(a, b, c) double xTmp(c)[2]; c_dd_mul_dd_d(a, b, xTmp(c)); double c = xTmp(c)[0]; +#define mul_d_dd_d(a, b, c) double xTmp(c)[2]; c_dd_mul_d_dd(a, b, xTmp(c)); double c = xTmp(c)[0]; +#define mul_ld_dd_d(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double xTmp(c)[2]; c_dd_mul(xTmp(a), b, xTmp(c)); c = xTmp(c)[0]; // TODO: check c = ((a) + (b[0])) + b[1]; +#define mul_dd_dd_d(a, b, c) double xTmp(c)[2]; c_dd_mul(a, b, xTmp(c)); double c; c = xTmp(c)[0]; +#define mul_d_ld_d(a, b, c) double c = (a) * (b); +#define mul_ld_ld_d(a, b, c) double c = (a) * (b); +#define mul_dd_ld_d(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double xTmp(c)[2]; c_dd_mul(a, xTmp(b), xTmp(c)); c = xTmp(c)[0]; // TODO: check double c = ((a[0]) + (b)) + a[1]; + +#define mul_d_d_dd(a, b, c) double c[2]; c[0] = two_mul_mine(a, b, &(c[1])); // check c[0] = (a) * (b) +#define mul_ld_d_dd(a, b, c) long double xTmp(c) = (a) * (b); double c[2]; c[0] = xTmp(c); c[1] = xTmp(c) - c[0]; +#define mul_dd_d_dd(a, b, c) double c[2]; c_dd_mul_dd_d(a, b, c); +#define mul_d_dd_dd(a, b, c) double c[2]; c_dd_mul_d_dd(a, b, c); +#define mul_ld_dd_dd(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double c[2]; c_dd_mul(xTmp(a), b, c); +#define mul_dd_dd_dd(a, b, c) double c[2]; c_dd_mul(a, b, c); +#define mul_d_ld_dd(a, b, c) long double xTmp(c) = (a) * (b); double c[2]; c[0] = xTmp(c); c[1] = c - xTmp(c)[0]; +#define mul_ld_ld_dd(a, b, c) long double xTmp(c) = (a) * (b); double c[2]; c[0] = xTmp(c); c[1] = c - xTmp(c)[0]; +#define mul_dd_ld_dd(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double c[2]; c_dd_mul(a, xTmp(b), c); + +#define mul_d_d_ld(a, b, c) long double c = (a) * (b); +#define mul_ld_d_ld(a, b, c) long double c = (a) * (b); +#define mul_dd_d_ld(a, b, c) double xTmp(c)[2]; c_dd_mul_dd_d(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define mul_d_dd_ld(a, b, c) double xTmp(c)[2]; c_dd_mul_d_dd(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define mul_dd_dd_ld(a, b, c) double xTmp(c)[2]; c_dd_mul(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define mul_ld_dd_ld(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double xTmp(c)[2]; c_dd_mul(xTmp(a), b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define mul_d_ld_ld(a, b, c) long double c; c = (a) * (b); +#define mul_ld_ld_ld(a, b, c) long double c; c = (a) * (b); +#define mul_dd_ld_ld(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double xTmp(c)[2]; c_dd_mul(a, xTmp(b), xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; + +// =========================================== DIV ======================================= +#define div_d_d_d(a, b, c) double c = (a) / (b); +#define div_ld_d_d(a, b, c) long double xTmp(c) = (a) / (b); double c = xTmp(c); +#define div_dd_d_d(a, b, c) double xTmp(c)[2]; c_dd_div_dd_d(a, b, xTmp(c)); double c = xTmp(c)[0]; +#define div_d_dd_d(a, b, c) double xTmp(c)[2]; c_dd_div_d_dd(a, b, xTmp(c)); double c = xTmp(c)[0]; +#define div_ld_dd_d(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double xTmp(c)[2]; c_dd_div(xTmp(a), b, xTmp(c)); c = xTmp(c)[0]; // TODO: check c = ((a) + (b[0])) + b[1]; +#define div_dd_dd_d(a, b, c) double xTmp(c)[2]; c_dd_div(a, b, xTmp(c)); double c; c = xTmp(c)[0]; +#define div_d_ld_d(a, b, c) double c = (a) / (b); +#define div_ld_ld_d(a, b, c) double c = (a) / (b); +#define div_dd_ld_d(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double xTmp(c)[2]; c_dd_div(a, xTmp(b), xTmp(c)); c = xTmp(c)[0]; // TODO: check double c = ((a[0]) + (b)) + a[1]; + +#define div_d_d_dd(a, b, c) double c[2]; c[0] = div_mine(a, b, &(c[1])); // check c[0] = (a) / (b) +#define div_ld_d_dd(a, b, c) long double xTmp(c) = (a) / (b); double c[2]; c[0] = xTmp(c); c[1] = xTmp(c) - c[0]; +#define div_dd_d_dd(a, b, c) double c[2]; c_dd_div_dd_d(a, b, c); +#define div_d_dd_dd(a, b, c) double c[2]; c_dd_div_d_dd(a, b, c); +#define div_ld_dd_dd(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double c[2]; c_dd_div(xTmp(a), b, c); +#define div_dd_dd_dd(a, b, c) double c[2]; c_dd_div(a, b, c); +#define div_d_ld_dd(a, b, c) long double xTmp(c) = (a) / (b); double c[2]; c[0] = xTmp(c); c[1] = c - xTmp(c)[0]; +#define div_ld_ld_dd(a, b, c) long double xTmp(c) = (a) / (b); double c[2]; c[0] = xTmp(c); c[1] = c - xTmp(c)[0]; +#define div_dd_ld_dd(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double c[2]; c_dd_div(a, xTmp(b), c); + +#define div_d_d_ld(a, b, c) long double c = (a) / (b); +#define div_ld_d_ld(a, b, c) long double c = (a) / (b); +#define div_dd_d_ld(a, b, c) double xTmp(c)[2]; c_dd_div_dd_d(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define div_d_dd_ld(a, b, c) double xTmp(c)[2]; c_dd_div_d_dd(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define div_dd_dd_ld(a, b, c) double xTmp(c)[2]; c_dd_div(a, b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define div_ld_dd_ld(a, b, c) double xTmp(a)[2]; xTmp(a)[0] = a; xTmp(a)[1] = a - xTmp(a)[0]; double xTmp(c)[2]; c_dd_div(xTmp(a), b, xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; +#define div_d_ld_ld(a, b, c) long double c; c = (a) / (b); +#define div_ld_ld_ld(a, b, c) long double c; c = (a) / (b); +#define div_dd_ld_ld(a, b, c) double xTmp(b)[2]; xTmp(b)[0] = b; xTmp(b)[1] = b - xTmp(b)[0]; double xTmp(c)[2]; c_dd_div(a, xTmp(b), xTmp(c)); long double c; c = xTmp(c)[0]; c += xTmp(c)[1]; + +#endif // _MULTITYPE diff --git a/includeDD/inline.h b/includeDD/inline.h new file mode 100644 index 0000000000000000000000000000000000000000..294dfc8c43a8a56d5d8613b66e65ad39d9f0398c --- /dev/null +++ b/includeDD/inline.h @@ -0,0 +1,69 @@ +#include +#include + +/* Computes fl(a+b) and err(a+b). */ +static inline double two_add_mine(double a, double b, double *err) +{ + double s = a + b; + double bb = s - a; + *err = (a - (s - bb)) + (b - bb); + return s; +} + +/* Computes fl(a*b) and err(a*b). */ +static inline double two_mul_mine(double a, double b, double *err) +{ + double p = a * b; + *err = fma(a, b, -p); + return p; +} + +static inline double two_sub_mine(double a, double b, double *err) +{ + double s = a - b; + double bb = s - a; + *err = (a - (s - bb)) - (b + bb); + return s; +} + +static inline double quick_two_sum_mine(double a, double b, double *err) +{ + double s = a + b; + *err = b - (s - a); + return s; +} + +static inline double div_mine(double a, double b, double *err) +{ + double q1, q2; + double p1, p2; + double s; + + q1 = a / b; + + /* Compute a - q1 * b */ + p1 = two_mul_mine(q1, b, &p2); + s = two_sub_mine(a, p1, err); + *err -= p2; + + /* get next approximation */ + q2 = (s + *err) / b; + + s = quick_two_sum_mine(q1, q2, err); + return s; +} + +// pow = exp(b * log(a)) +static inline void c_dd_pow_mine(double *a, double *b, double *result) { + if(a[0] == 0 && a[1] == 0) + { + result[0] = 0; + result[1] = 0; + return; + } + double tmp1[2]; + double tmp2[2]; + c_dd_log(a, tmp1); + c_dd_mul(tmp1, b, tmp2); + c_dd_exp(tmp2, result); +} \ No newline at end of file diff --git a/includeTEST/binary.h b/includeTEST/binary.h new file mode 100644 index 0000000000000000000000000000000000000000..77871a7cda55922fe857544715de50a0e3cb0429 --- /dev/null +++ b/includeTEST/binary.h @@ -0,0 +1,36 @@ +#ifndef MYHEAD +#define MYHEAD + +#define DATA_COUNT 1000000 + +#define MAXEXP 0x7ffu +#define MINEXP 0x001u +#define DEXPBIAS 0x3ff +#define DMANTISSA 0x000fffffffffffffll +#define DMANTWIDTH 52 +#define DEXPWIDTH 11 + +typedef unsigned long int _UL; +typedef long int _L; +typedef union +{ + double d; + unsigned long int l; +} DL; + +#define D(h,l) l,h +typedef union +{ + struct + { + unsigned int lo; + unsigned int hi; + } word; + + double d; +} du; +#define QNANHI 0xfff80000 +#define QNANLO 0x00000000 +#define QNANF 0xffc00000 + +#endif diff --git a/includeTEST/tgen.h b/includeTEST/tgen.h new file mode 100644 index 0000000000000000000000000000000000000000..61cf97e43da3fb4818187fb5b2c349d51926aa65 --- /dev/null +++ b/includeTEST/tgen.h @@ -0,0 +1,18 @@ +double sin_gen(double); +double cos_gen(double); +double tan_gen(double); +double asin_gen(double); +double acos_gen(double); +double atan_gen(double); +double sinh_gen(double); +double cosh_gen(double); +double tanh_gen(double); +double asinh_gen(double); +double acosh_gen(double); +double atanh_gen(double); +double exp_gen(double); +double exp2_gen(double); +double exp10_gen(double); +double log_gen(double); +double log2_gen(double); +double log10_gen(double); \ No newline at end of file diff --git a/src/devideUpEdgeData.cpp b/src/devideUpEdgeData.cpp new file mode 100644 index 0000000000000000000000000000000000000000..2e5827a9d0023aa9966588358ffa3668776ee872 --- /dev/null +++ b/src/devideUpEdgeData.cpp @@ -0,0 +1,325 @@ +#include +#include +#include +#include +#include +#include +#include +#include +clock_t start; +clock_t end; + +using namespace std; +#define DEBUG +#undef DEBUG + +int read_scanf(const string &filename, const int &cols, vector &_vector, const double &ulp) { + FILE *fp = fopen(filename.c_str(), "r"); + bool flag = true; + double tmp; + if (!fp) { + cout << "File open error!\n"; + return 0; + } + while (flag) { + //将txt文本文件中的一行数据存入rowArray中，并将rowArray存入vector中 + double *rowArray = new double[cols]; // new一个double类型的动态数组 + if (EOF == fscanf(fp, "%le", &rowArray[0])) { //读取数据，存在_vector[cols]中 + break; + } + if (EOF == fscanf(fp, "%le", &rowArray[1])) { + break; + } + if (EOF == fscanf(fp, "%le", &tmp)) { + break; + } + // if (tmp < ulp) { + // tmp = ulp; + // } + rowArray[2] = tmp; + _vector.push_back(rowArray); + } + fclose(fp); + return 1; +} + +//将采样点对应的行号/序号、坐标值、误差、距离存入vector中 +void storeData(double rows, double x, double ulpErr, double distance, vector &points) { + double *point = new double[4]; + point[0] = rows; + point[1] = x; + point[2] = ulpErr; + point[3] = distance; + points.push_back(point); +} + +//打印区间数据 +void printRegime(const int &num, vector &vector, string fintervalName, string fdistanceName, string fwidthName) { + ofstream finterval, fdistance, fwidth; + finterval.open(fintervalName); + fdistance.open(fdistanceName); + fwidth.open(fwidthName); + for (int i = 0; i < num; i += 2) { + // if (vector[i][1] == vector[i + 1][1]) { + // continue; // 注释掉后，单个突出也会记录 + // } + // finterval << "lineNo: [" << vector[i][0] << ", " << vector[i + 1][0] << "]; input: " << "[" << vector[i][1] << ", " << vector[i + 1][1] << "]" << endl; + finterval << "[" << vector[i][1] << ", " << vector[i + 1][1] << "]\n"; // print intervals data to file fintervalName + // fdistance << vector[i][3] << endl; // print intervals' distance data to file "distanceTmp.txt" + // fwidth << vector[i+1][3] << endl; // print intervals' width data to file "widthTmp.txt" + } + finterval.close(); + fdistance.close(); + fwidth.close(); +} + +void ifRegimeulp(const int &num, vector &vector, std::vector &store_vector, const double &threshold) { + double distance; + double k, b, tmp; + + for (int i = 0; i < num - 1; ++i) { + // find the interval start point and set the distance value. The distance value is the distance to the previous interval + if (vector[i][2] > threshold) { + if (i == 0) { + // if 1st interval, the distance value is 0 + distance = 0; + storeData(0, vector[i][1], vector[i][2], distance, store_vector); + } else if (vector[i - 1][2] <= threshold) { + if(store_vector.size() == 0) { + distance = 0; + } else { + distance = i - store_vector.back()[0]; + } + k = (vector[i][1] - vector[i - 1][1]) / (vector[i][2] - vector[i - 1][2]); + b = vector[i][1] - k * vector[i][2]; + tmp = k * threshold + b; + // storeData(i, vector[i][1], vector[i][2], distance, store_vector); // 保存的是误差大于阈值的最左边点 + storeData(i - 0.5, tmp, threshold, distance, store_vector); // 保存的左端点连线与阈值横线的交点，这也是为何i-0.5 + } + } + // find the interval end point and set the distance value. The distance value is the width of the current interval + if (vector[i][2] > threshold) { + if (vector[i + 1][2] <= threshold) { + distance = i - store_vector.back()[0]; + k = (vector[i + 1][1] - vector[i][1]) / (vector[i + 1][2] - vector[i][2]); + b = vector[i][1] - k * vector[i][2]; + tmp = k * threshold + b; + // storeData(i, vector[i][1], vector[i][2], distance, store_vector); // 保存的是误差大于阈值的最右边点 + storeData(i + 0.5, tmp, threshold, distance, store_vector); // 保存的是最右边点与阈值横线的交点，这也是为何i+0.5 + } + else if (i == num - 2) { + // the last point's ULP error > threshold, so store the last point, that is vector[num - 1] + distance = i - store_vector.back()[0]; + storeData(num - 1, vector[num - 1][1], vector[num - 1][2], distance, store_vector); + } + } + } +} + +// 函数的轮廓线数据，即每个轮廓点对应的输入、误差，是基于扩展后的投影面误差采样数据，利用matlab的boundary函数得到的，保存在matlab_vector向量中。 +// 基于轮廓线数据，调用ifRegimeulp函数，初步得到精度优化区间的起止点信息，保存在store_vector向量中。 +// store_vector 保存的是精度优化区间的起止点信息，即其在matlab_vector向量中对应的序号、输入、误差、距离（matlab序号）。 +// store_vector 每相邻的2个元素，分别对应区间的起点和止点信息，构成一个精度优化区间。每个元素中又有3个信息，即其在matlab_vector向量中对应的序号、输入、误差。 +// output_vector 保存的是精度优化区间的起止点信息，即其在matlab_vector向量中对应的序号、输入、误差、距离（真实输入）。 +int findMinMax(vector &matlab_vector, vector &store_vector, vector &output_vector) { + int vectorNum = store_vector.size(); + double start, end, startInput, endInput, currentInput, startError, endError, distance; + double maxInput, minInput, lastmaxInput, minIdx, maxIdx; + + for (int i = 0; i < vectorNum; i += 2) { + start = store_vector[i][0]; + end = store_vector[i+1][0]; + startInput = store_vector[i][1]; + endInput = store_vector[i+1][1]; + startError = store_vector[i][2]; + endError = store_vector[i+1][2]; + minIdx = start; + maxIdx = end; + minInput = startInput; + maxInput = endInput; + + // Considering that the boundary line may protrude on both sides, the start and end points stored in 'store_vector' that intersect with the threshold are not necessarily the endpoints (leftmost and rightmost points) of the interval. + // Therefore, it is necessary to traverse the coordinate information of all points within the interval and update the endpoints of the interval to the leftmost and rightmost points found. + // 考虑到轮廓线可能向两边突出去，所以store_vector中存储的与阈值相交的起止点并不一定是区间的端点（最左最右点）。因此，需要遍历区间内所有点的坐标信息，将区间端点更新为找到的最小值和最大值。 + for (int j = ceil(start); j <= floor(end); j++) { + currentInput = matlab_vector[j][1]; + if(currentInput < minInput) { + minIdx = j; + minInput = currentInput; + } else if(currentInput > maxInput) { + maxIdx = j; + maxInput = currentInput; + } + } + // store the start point to output vector + if(minInput != startInput) { + #ifdef DEBUG + printf("\tFor No.%d: minInput != startInput. startInput: %g, minInput: %g.\n", i / 2, startInput, minInput); + #endif + startInput = matlab_vector[minIdx][1]; + startError = matlab_vector[minIdx][2]; + } + // 更新distance，即当前区间左端点到其右侧区间右端点的距离 + if(i != 0) { + distance = minInput - lastmaxInput; + } else { + distance = 0; + } + storeData(minIdx, startInput, startError, distance, output_vector); + // store the end point to output vector + if(maxInput != endInput) { + #ifdef DEBUG + printf("\tFor No.%d: maxInput != endInput. endInput: %g, maxInput: %g.\n", i / 2, endInput, maxInput); + #endif + endInput = matlab_vector[maxIdx][1]; + endError = matlab_vector[maxIdx][2]; + } + distance = endInput - startInput; + storeData(maxIdx, endInput, endError, distance, output_vector); + + // update the lastmaxInput + lastmaxInput = maxInput; + + // cout << i / 2 << ": " << minIdx << " " << minInput << " " << startInput << " " << maxIdx << " " << maxInput << " " << endInput << endl; + // cout << "No." << i/2+1 << ": the index is [" << minInput << ", " << maxInput << "], the corresponding input is [" << startInput << ", " << endInput << "]" << endl; + } + + // the output vector is from edge, a Matlab variable which stores the index, the input and the error of the sample data + // minIdx is the index of the element in the output vector. + // minInput is the index of the element in the sample data. + // startInput is the corresponding input value in the sample data. + // cout << "i: minIdx minInput startInput maxIdx maxInput endInput" << endl; + // cout << i << ": " << minIdx << " " << minInput << " " << startInput << " " << maxIdx << " " << maxInput << " " << endInput << endl; + // cout << i << ": the output vector index is " << minIdx << " " << maxIdx << endl; + // cout << i << ": the sampleData index is " << minInput << " " << maxInput << endl; + // cout << i << ": the corresponding input to the new sampleData index is " << startInput << " " << endInput << endl; + // cout << i << ": the width is " << maxInput - minInput + 1 << "\n" << endl; + + return 1; +} + +// add_zzy interval_merge +// TODO-done: Just compare the input value but not the idx of input. Because the idx is for the boundary line which is useless. +void interval_merge(vector in_vector, vector &out_vector, const double &thresholdCombine) { + if(in_vector.size() == 0) { + cout << "NOTE: interval_merge: none intervals" << endl; + return; + } + if (in_vector.size() == 2) { + double *temp = in_vector[0]; + out_vector.push_back(temp); + double *temp1 = in_vector[1]; + out_vector.push_back(temp1); + return; + } + double *temp = in_vector.front(); + out_vector.push_back(temp); + // cout << out_vector.size() << endl; + for (int i = 2; i < in_vector.size(); i += 2) { + // cout << i << " " << out_vector.size() << endl; + // if (in_vector.at(i)[3] < thresholdCombine && in_vector.at(i - 1)[3] > thresholdCombine && in_vector.at(i + 1)[3] > thresholdCombine) { + if (in_vector.at(i)[3] < thresholdCombine) { + if (i == in_vector.size() - 2) { + temp = in_vector.at(i + 1); + out_vector.push_back(temp); + } else { + continue; + } + } else { + temp = in_vector.at(i-1); + out_vector.push_back(temp); + temp = in_vector.at(i); + out_vector.push_back(temp); + // temp = in_vector.at(i+1); + // out_vector.push_back(temp); + if (i == in_vector.size() - 2) { + temp = in_vector.at(i + 1); + out_vector.push_back(temp); + } + } + } + // cout << out_vector.size() << endl; + // for (int i = 0; i < out_vector.size(); i += 2) { + // cout << "["; + // cout << out_vector.at(i)[1] << "," << out_vector.at(i + 1)[1] << "]" < extractInfo(vector &origin) { + vector destination; + for (size_t i = 0; i < origin.size(); i++) { + destination.push_back(origin[i][1]); + } + return destination; +} + +vector devideUpEdgeData(string upEdgeFileName, double threshold, double thresholdCombine) { + vector intervalData; + // thresholdCombine = 0; // TODO: set thresholdCombine + // upEdgeFileName 文件中有3列，分别是序号、输入、误差。该文件是由matlab生成，故文件中的序号从1开始，而非0。 + // matlab_vector可视为二维数组，保存的是 upEdgeFileName 中的边界点数据，包括序号、输入、误差;输出数组元素： + int columns = 3; + vector matlab_vector; + if (!read_scanf(upEdgeFileName, columns, matlab_vector, threshold)) { + fprintf(stderr, "ERROR: devideUpEdgeData: read file error\n"); + exit(EXIT_FAILURE); + } + + // 初步获取区间信息 + vector store_vector; + ifRegimeulp(matlab_vector.size(), matlab_vector, store_vector, threshold); + // printRegime(store_vector.size(), store_vector, "intervalTmp.txt", "distanceTmp.txt", "widthTmp.txt"); + + // 进一步优化区间信息 + vector output_vector; + findMinMax(matlab_vector, store_vector, output_vector); + // printRegime(output_vector.size(), output_vector, "intervalTmp_beforeCombine.txt", "distanceTmp_beforeCombine.txt", "widthTmp_beforeCombine.txt"); + + // Merge intervals that are close in distance. + vector merge_vector; + interval_merge(output_vector, merge_vector, thresholdCombine); + // printRegime(merge_vector.size(), merge_vector, intervalFileName, "distanceTmp_afterCombine.txt", "widthTmp_afterCombine.txt"); + + // generate the interval data vector + intervalData = extractInfo(merge_vector); + return intervalData; +} + +// void devideKernel(const string &filePath, const string &intervalFileName, const double &threshold, const double &thresholdCombine) { +// // filePath文件中有3列，分别是序号、输入、误差。该文件是由matlab生成，故文件中的序号从1开始，而非0。 +// // matlab_vector可视为二维数组，保存的是filePath中的边界点数据，包括序号、输入、误差;输出数组元素： +// vector matlab_vector; +// int columns = 3; +// if (!read_scanf(filePath, columns, matlab_vector, threshold)) { +// return; +// } +// // 初步获取区间信息 +// vector store_vector; +// ifRegimeulp(matlab_vector.size(), matlab_vector, store_vector, threshold); +// // printRegime(store_vector.size(), store_vector, "intervalTmp.txt", "distanceTmp.txt", "widthTmp.txt"); +// // 进一步优化区间信息 +// vector output_vector; +// findMinMax(matlab_vector, store_vector, output_vector); +// // printRegime(output_vector.size(), output_vector, "intervalTmp_beforeCombine.txt", "distanceTmp_beforeCombine.txt", "widthTmp_beforeCombine.txt"); +// // Merge intervals that are close in distance. +// vector merge_vector; +// interval_merge(output_vector, merge_vector, thresholdCombine); +// printRegime(merge_vector.size(), merge_vector, intervalFileName, "distanceTmp_afterCombine.txt", "widthTmp_afterCombine.txt"); +// } +// void devideUpEdgeData(string upEdgeFileName, string intervalFileName, double threshold, double thresholdCombine) { +// // cout << upEdgeFileName << " " << threshold << "\n"; +// // start = clock(); +// devideKernel(upEdgeFileName, intervalFileName, threshold, thresholdCombine); +// // end = clock(); +// // double endtime = (double)(end - start)/CLOCKS_PER_SEC; +// // cout << "Total time:" << endtime << endl; +// // cout << "Total time:" << endtime * 1000 << "ms" << endl; +// } \ No newline at end of file diff --git a/src/expandAST.cpp b/src/expandAST.cpp index 84064638747ed2db4810b93b24e0034719a54c7e..8a828911516cabc96513bbb66697abc4a2025a25 100644 --- a/src/expandAST.cpp +++ b/src/expandAST.cpp @@ -135,7 +135,7 @@ ast_ptr expandExpr(const ast_ptr &expr) return exprFinal; } - else if((opL != '+') && (opR == '+'))// TODO: opL or opR is / + else if((opL != '+') && (opR == '+')) { ast_ptr &lhsR = binOpR->getLHS(); ast_ptr &rhsR = binOpR->getRHS(); @@ -156,7 +156,7 @@ ast_ptr expandExpr(const ast_ptr &expr) // { // } - // else if((opL != '/') && (opR == '/'))// TODO: opL or opR is / + // else if((opL != '/') && (opR == '/')) // { // } diff --git a/src/exprAuto.cpp b/src/exprAuto.cpp index 9762b596bc78f43f78961762f6aef32f34650ab3..3b874d90120cd1c05229028d470d8636f3d2bde5 100644 --- a/src/exprAuto.cpp +++ b/src/exprAuto.cpp @@ -1007,6 +1007,32 @@ size_t computeRandomNum(size_t sum) return randomNum; } +// for main +string pickTheBest(string uniqueLabel, vector testSet, vector intervals, vector scales) +{ + string bestExpr = "origin"; + double maxError = INFINITY; + double aveError = INFINITY; + for (size_t i = 0; i < testSet.size(); i++) + { + string suffixTmp = testSet.at(i); + cout << "*-*-*-pickTheBest: for item No." << i << ": type = " << suffixTmp << endl; + + // generate function code and test error + auto tempError = testError(uniqueLabel, suffixTmp, intervals, scales); + // cout << "pickTheBest: for item No." << i << ": maxError: " << tempError.maxError << "\n"; + // cout << "pickTheBest: for item No." << i << ": aveError: " << tempError.aveError << "\n"; + + if ((tempError.maxError < maxError) || ((tempError.maxError == maxError) && (tempError.aveError < aveError))) + { + bestExpr = suffixTmp; + maxError = tempError.maxError; + aveError = tempError.aveError; + } + } + return bestExpr; +} + size_t pickTheBest(vector &items, ast_ptr &originExpr) { auto uniqueLabel = getUniqueLabel(); @@ -1016,8 +1042,8 @@ size_t pickTheBest(vector &items, ast_ptr &originExpr) getVariablesFromExpr(originExpr, vars); auto funcNameMpfr = geneMpfrCode(originExpr, uniqueLabel, vars); - string filename = "./intervalData.txt"; // TODO: get the filename from uniqueLabel - auto intervalData = getIntervalData(filename); + string filename = "./intervalData.txt"; // TODO: (useful for multi parameters) get the filename from uniqueLabel + auto intervalData = getIntervalData(filename); // TODO: (useful for multi parameters) update to getInterval(filename, dimension) size_t scale; auto intervalTmp = intervalData.at(0); auto dimension = intervalTmp.size() / 2; @@ -1044,7 +1070,7 @@ size_t pickTheBest(vector &items, ast_ptr &originExpr) // generate function code and test error string suffixTmp = suffix + std::to_string(i); - geneOriginCodeKernel(item, vars, uniqueLabel, suffixTmp); + geneExprCodeKernel(item, vars, uniqueLabel, suffixTmp); // auto timeTmp1 = std::chrono::high_resolution_clock::now(); auto tempError = testError(uniqueLabel, suffixTmp, intervalTmp, scales); // auto timeTmp2 = std::chrono::high_resolution_clock::now(); @@ -1083,6 +1109,7 @@ size_t pickTheBest(vector &items, ast_ptr &originExpr) return maxIdx; } +// be useful for multi parameters test void geneSampleData() { // static size_t callCount = 0; @@ -1092,8 +1119,8 @@ void geneSampleData() // prompt += "geneSampleData: "; cout << "geneSampleData start--------" < parameters; @@ -1366,9 +1393,9 @@ vector exprAutoWrapper(ast_ptr &expr, const std::vector &interv getVariablesFromExpr(expr, vars); auto exprStr = PrintExpression(expr); - auto funcNameOrigin = geneOriginCodeKernel(exprStr, vars, uniqueLabel, "origin"); + auto funcNameOrigin = geneExprCodeKernel(exprStr, vars, uniqueLabel, "origin"); auto expr1Str = PrintExpression(expr1); - auto funcNameSympy = geneOriginCodeKernel(expr1Str, vars, uniqueLabel, "sympy"); + auto funcNameSympy = geneExprCodeKernel(expr1Str, vars, uniqueLabel, "sympy"); auto funcNameMpfr = geneMpfrCode(exprStr, uniqueLabel, vars); // int scale = 256; diff --git a/src/geneCode.cpp b/src/geneCode.cpp index 115af756a6d738542529c562d9b5bbba5190b6a8..c792f1f92e4b75103c6884a5115c50e8413da27c 100644 --- a/src/geneCode.cpp +++ b/src/geneCode.cpp @@ -3,6 +3,8 @@ #include "parserASTLY.hpp" #include +#include +#include using std::cout; using std::endl; @@ -115,7 +117,7 @@ bool getVariablesFromExpr(const ast_ptr &expr, vector &vars) return true; } -string geneOriginCodeKernel(string exprStr, vector vars, string uniqueLabel, string tail) +string geneExprCodeKernel(string exprStr, vector vars, string uniqueLabel, string tail) { // print input info // cout << "expression : " << exprStr << endl; @@ -152,13 +154,56 @@ string geneOriginCodeKernel(string exprStr, vector vars, string uniqueLa return funcName; } -string geneOriginCode(string exprStr, string uniqueLabel, string tail) +string geneExprCode(string exprStr, string uniqueLabel, string tail) { auto originExpr = ParseExpressionFromString(exprStr); vector vars; getVariablesFromExpr(originExpr, vars); - auto funcName = geneOriginCodeKernel(exprStr, vars, uniqueLabel, tail); + auto funcName = geneExprCodeKernel(exprStr, vars, uniqueLabel, tail); + + return funcName; +} + +string geneTGenCode(string exprStr, vector vars, string uniqueLabel, string tail) +{ + // replace func to func_gen + string exprNewStr(exprStr); + vector supportFuncList = {"sin", "cos", "tan", "exp", "log", "sinh", "cosh", "tanh"}; + for(auto replaceFrom : supportFuncList) + { + string replaceTo(replaceFrom); + replaceTo.append("_gen("); + replaceFrom += "("; + int pos = exprNewStr.find(replaceFrom); + while(pos != -1) + { + exprNewStr.replace(pos, replaceFrom.length(), replaceTo); + pos = exprNewStr.find(replaceFrom); + } + } + // generate TGen code + std::ofstream fout; + string directory = "srcTest/" + uniqueLabel + "/"; + string funcName = "expr_" + uniqueLabel + "_" + tail; + string fileName = directory + funcName + ".c"; + fout.open(fileName); + fout << "#include \n"; + fout << "#include \"tgen.h\"\n"; // this one is important + fout << "double " << funcName << "("; + for (size_t i = 0; i < vars.size(); ++i) + { + if (i != vars.size() - 1) + fout << "double" << " " << vars.at(i) << ", "; + else + fout << "double" << " " << vars.at(i); + } + fout << ")\n"; + fout << "{\n"; + fout << "\t" << "double result = " << exprNewStr << ";\n"; + fout << "\t" << "return result;\n"; + fout << "}\n"; + fout.close(); return funcName; } @@ -189,6 +234,62 @@ void geneHerbieCode(string exprstr, vector cs, string exprname, double v fout.close(); } +// TODO: not real implementation +string geneHerbieCode(string uniqueLabel) +{ + map benchmarkHerbie = { + {"Bsplines3", "-pow(x, 3.0) /6.0"}, + {"exp1x", "expm1(x) / x"}, + {"exp1x_log", "expm1(x) / x"}, + {"intro_example", "(x / (1.0 + pow(x, 3.0))) * fma(x, x, (1.0 - x))"}, + {"logexp", "log1p(exp(x))"}, + {"NMSEexample31", "(x + (1.0 - x)) / fma(pow(x, 0.25), pow(x, 0.25), sqrt((x + 1.0)))"}, + {"NMSEexample310", "log1p(-x) / log1p(x)"}, + {"NMSEexample34", "tan((x / 2.0))"}, + {"NMSEexample35", "atan2((x + (1.0- x)), (1.0+fma(sqrt(x), sqrt(x), (x * x))))"}, + {"NMSEexample36", ""}, + {"NMSEexample37", "expm1(x)"}, + {"NMSEexample38", "fma(x, log((1.0 + (1.0 / x))), log1p(x)) + -1.0"}, + {"NMSEexample39", ""}, + {"NMSEproblem331", "-1.0 / fma(x, x, x)"}, + {"NMSEproblem333", "((1.0 / (x + -1.0)) + (1.0 / (1.0 + x))) + (-2.0 / x)"}, + {"NMSEproblem334", ""}, + {"NMSEproblem336", "log1p((1.0 / x))"}, + {"NMSEproblem337", ""}, + {"NMSEproblem341", "(sin(x) / x) * (tan((x / 2.0)) / x)"}, + {"NMSEproblem343", ""}, + {"NMSEproblem344", "sqrt((1.0 + exp(x)))"}, + {"NMSEproblem345", ""}, + {"NMSEsection311", "(1.0 + expm1(x)) / expm1(x)"}, + {"predatorPrey", ""}, + {"sine", ""}, + {"sineorder3", ""}, + {"sqroot", ""}, + {"sqrt_add", ""}, + {"test05_nonlin1_r4", ""}, + {"test05_nonlin1_test2", ""}, + {"verhulst", ""}, + }; + + auto pos = benchmarkHerbie.find(uniqueLabel); + if (pos != benchmarkHerbie.end()) + { + string herbieExpr = pos->second; + if(herbieExpr != "") + { + geneExprCode(herbieExpr, uniqueLabel, "herbie"); + return herbieExpr; + } + else + { + fprintf(stderr, "ERROR: geneHerbieCode: we can not handle %s\n", uniqueLabel.c_str()); + exit(EXIT_FAILURE); + } + } + fprintf(stderr, "ERROR: geneHerbieCode: we can not handle %s\n", uniqueLabel.c_str()); + exit(EXIT_FAILURE); +} + void geneDaisyCode(string exprStr) { cout << exprStr << endl; @@ -208,14 +309,19 @@ string geneMpfrCode(const ast_ptr &exprAst, const string uniqueLabel, vector vars; // getVariablesFromExpr(exprAst, vars); size_t mpfr_variables = 0; @@ -274,11 +380,12 @@ string geneFinalCodeKernel(string exprStr, string uniqueLabel, std::vector" << std::endl; + fout << "#include \n"; fout << "double " << funcName << "("; for (size_t i = 0; i < vars.size(); ++i) { @@ -299,7 +406,7 @@ string geneFinalCodeKernel(string exprStr, string uniqueLabel, std::vector exprIn auto funcName = geneFinalCodeKernel(exprStr, uniqueLabel, exprInfoVector, vars); return funcName; -} \ No newline at end of file +} + +void getDepth(ast_ptr &expr, int &depth) +{ + static int depthNow = 0; + depthNow++; + auto type = expr->type(); + if(type == "Number" || type == "Variable") + { + // nothing + } + else if(type == "Call") + { + CallExprAST *callExpr = dynamic_cast(expr.get()); + auto &args = callExpr->getArgs(); + for(auto& arg : args) + { + getDepth(arg, depth); + } + } + else if(type == "Binary") + { + BinaryExprAST *binPtr = dynamic_cast(expr.get()); + ast_ptr &lhs = binPtr->getLHS(); + ast_ptr &rhs = binPtr->getRHS(); + getDepth(lhs, depth); + getDepth(rhs, depth); + } + else + { + fprintf(stderr, "ERROR: unknown type %s", type.c_str()); + exit(EXIT_FAILURE); + } + if(depthNow > depth) + { + depth = depthNow; + } + depthNow--; +} + +void setTypeKernel(ast_ptr &expr, const vector &opTypes) +{ + static size_t depthNow = 0; + auto type = expr->type(); + auto opType = opTypes.at(depthNow); + depthNow++; + if(type == "Number" || type == "Variable") + { + expr->setOpType(opType); + } + else if(type == "Call") + { + CallExprAST *callExpr = dynamic_cast(expr.get()); + auto &args = callExpr->getArgs(); + for(auto& arg : args) + { + setTypeKernel(arg, opTypes); + } + expr->setOpType(opType); + } + else if(type == "Binary") + { + BinaryExprAST *binPtr = dynamic_cast(expr.get()); + ast_ptr &lhs = binPtr->getLHS(); + ast_ptr &rhs = binPtr->getRHS(); + setTypeKernel(lhs, opTypes); + setTypeKernel(rhs, opTypes); + expr->setOpType(opType); + } + else + { + fprintf(stderr, "ERROR: unknown type %s", type.c_str()); + exit(EXIT_FAILURE); + } + depthNow--; +} + +// set opTypes for different tree levels by middle +void setType(ast_ptr &expr, int depth, int middle) +{ + vector opTypes; + for(int i = 0; i < middle; i++) + { + opTypes.push_back("DD"); + } + for(int i = middle; i < depth; i++) + { + opTypes.push_back("double"); + } + + // start setting the operators' type by vector opTypes + setTypeKernel(expr, opTypes); +} + +// set opTypes for different tree nodes by opTypes +void setType(ast_ptr &expr, vector opTypes) +{ + auto type = expr->type(); + if(type == "Number" || type == "Variable") + { + return; + } + + auto order = expr->getOrder(); + auto opType = opTypes.at(order); + // cout << "setType: NO." << order << ", type = " << type << ", opType = " << opType << "\n"; + expr->setOpType(opType); + if(type == "Call") + { + CallExprAST *callExpr = dynamic_cast(expr.get()); + auto &args = callExpr->getArgs(); + for(auto& arg : args) + { + setType(arg, opTypes); + } + } + else if(type == "Binary") + { + BinaryExprAST *binPtr = dynamic_cast(expr.get()); + ast_ptr &lhs = binPtr->getLHS(); + ast_ptr &rhs = binPtr->getRHS(); + setType(lhs, opTypes); + setType(rhs, opTypes); + } + else + { + fprintf(stderr, "ERROR: unknown type %s", type.c_str()); + exit(EXIT_FAILURE); + } +} + +// set opTypes for different tree nodes by opTypes +void setType(ast_ptr &expr, map opTypes) +{ + auto type = expr->type(); + if(type == "Number" || type == "Variable") + { + return; + } + + auto order = expr->getOrder(); + auto opType = opTypes.at(order); + // cout << "setType: NO." << order << ", type = " << type << ", opType = " << opType << "\n"; + expr->setOpType(opType); + if(type == "Call") + { + CallExprAST *callExpr = dynamic_cast(expr.get()); + auto &args = callExpr->getArgs(); + for(auto& arg : args) + { + setType(arg, opTypes); + } + } + else if(type == "Binary") + { + BinaryExprAST *binPtr = dynamic_cast(expr.get()); + ast_ptr &lhs = binPtr->getLHS(); + ast_ptr &rhs = binPtr->getRHS(); + setType(lhs, opTypes); + setType(rhs, opTypes); + } + else + { + fprintf(stderr, "ERROR: unknown type %s", type.c_str()); + exit(EXIT_FAILURE); + } +} + +void setOrdersKernel(ast_ptr &expr, int &orderNow, vector &opOrder) +{ + auto type = expr->type(); + if(type == "Number") + { + expr->setOrder(orderNow); + orderNow++; + } + else if(type == "Variable") + { + expr->setOrder(orderNow); + orderNow++; + } + else if(type == "Call") + { + CallExprAST *callPtr = dynamic_cast(expr.get()); + auto &args = callPtr->getArgs(); + vector paramOrders; + for(auto& arg : args) + { + setOrdersKernel(arg, orderNow, opOrder); + } + expr->setOrder(orderNow); + opOrder.push_back(orderNow); + orderNow++; + } + else if(type == "Binary") + { + BinaryExprAST *binPtr = dynamic_cast(expr.get()); + ast_ptr &lhs = binPtr->getLHS(); + ast_ptr &rhs = binPtr->getRHS(); + setOrdersKernel(lhs, orderNow, opOrder); + setOrdersKernel(rhs, orderNow, opOrder); + expr->setOrder(orderNow); + opOrder.push_back(orderNow); + orderNow++; + } + else + { + fprintf(stderr, "ERROR: unknowntype %s\n", type.c_str()); + exit(EXIT_FAILURE); + } +} + +// set orders for nodes in expr and return the operation nodes' order (aka binary & call) +vector setOrder(ast_ptr &expr) +{ + int order = 0; + vector opOrder; + setOrdersKernel(expr, order, opOrder); + return opOrder; +} + +static std::map opMap = { + {'+', "add"}, + {'-', "sub"}, + {'*', "mul"}, + {'/', "div"}}; + +static std::map callMap = { + {"sin", "c_dd_sin"}, + {"cos", "c_dd_cos"}, + {"tan", "c_dd_tan"}, + {"exp", "c_dd_exp"}, + // {"exp2", "c_dd_exp2"}, + // {"exp10", "c_dd_exp10"}, + {"log", "c_dd_log"}, + // {"log2", "c_dd_log2"}, + {"log10", "c_dd_log10"}, + {"asin", "c_dd_asin"}, + {"acos", "c_dd_acos"}, + {"atan", "c_dd_atan"}, + {"sinh", "c_dd_sinh"}, + {"cosh", "c_dd_cosh"}, + {"tanh", "c_dd_tanh"}, + {"asinh", "c_dd_asinh"}, + {"acosh", "c_dd_acosh"}, + {"atanh", "c_dd_atanh"}, + {"sqrt", "c_dd_sqrt"}, + {"pow", "c_dd_pow_mine"} +}; + +int codegenKernel(ofstream &ofs, const ast_ptr &expr) +{ + auto type = expr->type(); + auto opType = expr->getOpType(); + auto order = expr->getOrder(); + if(type == "Number") + { + NumberExprAST *numPtr = dynamic_cast(expr.get()); + auto num = numPtr->getNumber(); + if(opType == "double") + { + ofs << "\t" << "double tmp" << order << " = " << num << ";\n"; + } + else if(opType == "ld") + { + ofs << "\t" << "long double tmp" << order << " = " << num << ";\n"; + } + else if(opType == "DD") + { + ofs << "\t" << "double tmp" << order << "[2];\n"; + ofs << "\t" << "tmp" << order << "[0] = " << num << ";\n"; + ofs << "\t" << "tmp" << order << "[1] = 0.0;\n"; + } + else + { + fprintf(stderr, "opType %s of number %lf is error\n", opType.c_str(), num); + exit(EXIT_FAILURE); + } + return order; + } + else if(type == "Variable") + { + VariableExprAST *varPtr = dynamic_cast(expr.get()); + auto var = varPtr->getVariable(); + if(opType == "double") + { + ofs << "\t" << "double tmp" << order << " = " << var << ";\n"; + } + else if(opType == "ld") + { + ofs << "\t" << "long double tmp" << order << " = " << var << ";\n"; + } + else if(opType == "DD") + { + ofs << "\t" << "double tmp" << order << "[2];\n"; + ofs << "\t" << "tmp" << order << "[0] = " << var << ";\n"; + ofs << "\t" << "tmp" << order << "[1] = 0.0;\n"; + } + else + { + fprintf(stderr, "opType %s of variable %s is error\n", opType.c_str(), var.c_str()); + exit(EXIT_FAILURE); + } + return order; + } + else if(type == "Call") // the parameters' type should be the same to the return value's type + { + CallExprAST *callExpr = dynamic_cast(expr.get()); + auto callee = (callExpr->getCallee()); + auto &args = callExpr->getArgs(); + vector paramOpTypes; + vector paramOrders; + for(const auto& arg : args) + { + auto opType = arg->getOpType(); + paramOpTypes.push_back(opType); + auto paramOrder = codegenKernel(ofs, arg); + paramOrders.push_back(paramOrder); + } + + if(opType == "double") + { + for(size_t i = 0; i < paramOpTypes.size(); i++) + { + auto ¶mOpType = paramOpTypes.at(i); + if(paramOpType == "double") + { + ofs << "\t" << "double p" << paramOrders.at(i) << " = tmp" << paramOrders.at(i) << ";\n"; + } + else if(paramOpType == "ld") // ld to double is straight. + { + ofs << "\t" << "double p" << paramOrders.at(i) << " = tmp" << paramOrders.at(i) << ";\n"; + } + else if(paramOpType == "DD") + { + ofs << "\t" << "double p" << paramOrders.at(i) << " = tmp" << paramOrders.at(i) << "[0];\n"; + } + else + { + fprintf(stderr, "ERROR: Unknown paramOpType %s\n", paramOpType.c_str()); + exit(EXIT_FAILURE); + } + } + ofs << "\t" << "double tmp" << order << " = " << callee << "("; + for(size_t i = 0; i < paramOrders.size() - 1; i++) + { + ofs << "p" << paramOrders.at(i) << ", "; + } + ofs << "p" << paramOrders.back(); + ofs << ");\n"; + } + else if(opType == "ld") + { + for(size_t i = 0; i < paramOpTypes.size(); i++) + { + auto ¶mOpType = paramOpTypes.at(i); + if(paramOpType == "double") + { + ofs << "\t" << "long double p" << paramOrders.at(i) << " = tmp" << paramOrders.at(i) << ";\n"; + } + else if(paramOpType == "ld") + { + ofs << "\t" << "long double p" << paramOrders.at(i) << " = tmp" << paramOrders.at(i) << ";\n"; + } + else if(paramOpType == "DD") + { + ofs << "\t" << "long double p" << paramOrders.at(i) << " = tmp" << paramOrders.at(i) << "[0];\n"; + ofs << "\t" << "p" << paramOrders.at(i) << " += tmp" << paramOrders.at(i) << "[1];\n"; + } + else + { + fprintf(stderr, "ERROR: Unknown paramOpType %s\n", paramOpType.c_str()); + exit(EXIT_FAILURE); + } + } + ofs << "\t" << "long double tmp" << order << " = " << callee << "("; + for(size_t i = 0; i < paramOrders.size() - 1; i++) + { + ofs << "p" << paramOrders.at(i) << ", "; + } + ofs << "p" << paramOrders.back(); + ofs << ");\n"; + } + else if(opType == "DD") + { + for(size_t i = 0; i < paramOpTypes.size(); i++) + { + auto ¶mOpType = paramOpTypes.at(i); + if(paramOpType == "double") + { + ofs << "\t" << "double p" << paramOrders.at(i) << "[2];\n"; + ofs << "\t" << "p" << paramOrders.at(i) << "[0] = tmp" << paramOrders.at(i) << ";\n"; + ofs << "\t" << "p" << paramOrders.at(i) << "[1] = 0.0;\n"; + } + else if(paramOpType == "ld") + { + ofs << "\t" << "double p" << paramOrders.at(i) << "[2];\n"; + ofs << "\t" << "p" << paramOrders.at(i) << "[0] = tmp" << paramOrders.at(i) << ";\n"; + ofs << "\t" << "p" << paramOrders.at(i) << "[1] = tmp" << paramOrders.at(i) << " - p" << paramOrders.at(i) << "[0];\n"; + } + else if(paramOpType == "DD") + { + ofs << "\t" << "double p" << paramOrders.at(i) << "[2];\n"; + ofs << "\t" << "p" << paramOrders.at(i) << "[0] = tmp" << paramOrders.at(i) << "[0];\n"; + ofs << "\t" << "p" << paramOrders.at(i) << "[1] = tmp" << paramOrders.at(i) << "[1];\n"; + } + else + { + fprintf(stderr, "ERROR: Unknown paramOpType %s\n", paramOpType.c_str()); + exit(EXIT_FAILURE); + } + } + ofs << "\t" << "double tmp" << order << "[2];\n"; + auto it = callMap.find(callee); + auto callStr = it->second; + ofs << "\t" << callStr << "("; + for(size_t i = 0; i < paramOrders.size(); i++) + { + ofs << "p" << paramOrders.at(i) << ", "; + } + ofs << "tmp" << order; + ofs << ");\n"; + } + else + { + fprintf(stderr, "ERROR: Unknown opType %s\n", opType.c_str()); + exit(EXIT_FAILURE); + } + + return order; + } + else if(type == "Binary") + { + BinaryExprAST *binPtr = dynamic_cast(expr.get()); + char op = binPtr->getOp(); + auto it = opMap.find(op); + auto opStr = it->second; + ast_ptr &lhs = binPtr->getLHS(); + ast_ptr &rhs = binPtr->getRHS(); + auto orderL = codegenKernel(ofs, lhs); + auto orderR = codegenKernel(ofs, rhs); + auto opTypeL = lhs->getOpType(); + auto opTypeR = rhs->getOpType(); + std::map opType2funcType = { + {"double", "d"}, + {"ld", "ld"}, + {"DD", "dd"}, + }; + auto funcType = opType2funcType.at(opType); + auto funcTypeL = opType2funcType.at(opTypeL); + auto funcTypeR = opType2funcType.at(opTypeR); + ofs << "\t" << opStr << "_" << funcTypeL << "_" << funcTypeR << "_" << funcType << "(tmp" << orderL << ", tmp" << orderR << ", tmp" << order << ");\n"; + + return order; + } + else + { + fprintf(stderr, "ERROR: unsupported type: type = %s\n", type.c_str()); + exit(EXIT_FAILURE); + } +} + +void codegen(ast_ptr &expr, vector &vars, const string funcName, ofstream &ofs) +{ + auto opType = expr->getOpType(); + auto order = expr->getOrder(); + ofs << "#include \n"; + ofs << "#include \n"; + ofs << "#include \"dd.h\"\n"; + ofs << "#include \"inline.h\"\n"; + ofs << "double " << funcName << "("; + for (size_t i = 0; i < vars.size() - 1; i++) + { + ofs << "double" << " " << vars.at(i) << ", "; + } + ofs << "double" << " " << vars.back(); + ofs << ") {\n"; + codegenKernel(ofs, expr); + if(opType == "double" || opType == "ld") + { + ofs << "\t" << "double result = tmp" << order << ";\n"; + } + else if(opType == "DD") + { + ofs << "\t" << "double result = tmp" << order << "[0];\n"; + } + else + { + fprintf(stderr, "ERROR: unsupported opType: opType = %s\n", opType.c_str()); + exit(EXIT_FAILURE); + } + ofs << "\t" << "return result;\n"; + ofs << "}" << endl; +} + +// Note: Each operator node has two states, double or double-double. if there are n node in that expression, there are 2^n states. +// Generate double-double implementations in all states and return the number of generated codes +int codegenWrapper(ast_ptr &expr, vector &vars, const string uniqueLabel, string tail) +{ + // AST init + auto opOrder = setOrder(expr); + auto lenOp = opOrder.size(); + + // set opTypes: method NO.1 + // int depth = 0; + // getDepth(expr, depth); + // cout << "codegen: depth: " << depth << endl; + // int middle = (depth - 3) > 1 ? (depth - 3) : 1; // !!! + // int middle = depth - 1; + // setType(expr, depth, middle); + // set opTypes: method NO.2 + // auto order = expr->getOrder(); + // vector opTypes; + // for(int j = 0; j <= order; j++) + // { + // auto tmp = currentNum % 2; + // if(tmp == 0) + // { + // opTypes.push_back("DD"); + // } + // else + // { + // opTypes.push_back("double"); + // } + // currentNum = currentNum >> 1; + // } + + vector dataTypes = {"DD", "ld", "double"}; + int lenDataTypes = dataTypes.size(); + int maxNum = pow(lenDataTypes, lenOp); // for dataTypes = {"ld", "double"}, maxNum = 1 << lenOp; + for(int num = 0; num < maxNum; num++) + { + // set opTypes: method NO.3 + auto currentNum = num; + std::map opTypes; + for(size_t i = 0; i < lenOp; i++) + { + auto id = opOrder[i]; + auto tmp = currentNum % lenDataTypes; + opTypes[id] = dataTypes[tmp]; + currentNum = currentNum / lenDataTypes; + } + setType(expr, opTypes); + + // init to generate code + string directory = "srcTest/" + uniqueLabel + "/"; + string funcName = "expr_" + uniqueLabel + "_" + tail + "_" + to_string(num); + string fileName = directory + funcName + ".c"; + cout << "fileName: " << fileName << "\topTypes: "; + for(map::iterator it = opTypes.begin(); it != opTypes.end(); it++) + { + cout << it->first << " " << it->second << ", "; + } + cout << "\n"; + ofstream file_clean(fileName, ios_base::out); + ofstream ofs(fileName, ios::app); + + // call codegen to generate code + codegen(expr, vars, funcName, ofs); + } + + return maxNum; +} diff --git a/src/getUpEdge.m b/src/getUpEdge.m new file mode 100644 index 0000000000000000000000000000000000000000..c45ef64ecf17465b5c1dde69b6583635ddde8087 --- /dev/null +++ b/src/getUpEdge.m @@ -0,0 +1,201 @@ +% disp(sampleFileName); +sampleData = importdata(sampleFileName); +lenTmp = length(sampleData); +maxTmp = max(sampleData(:,2)); +maxIndex = find(sampleData(:,2)>=maxTmp); +meanTmp = mean(sampleData(:,2)); +upperTmp = (maxTmp + meanTmp) / 10; +domainTmp = 0.01; +leftTmp = int64(max(1, maxIndex - lenTmp * domainTmp * 0.5)); +rightTmp = int64(min(lenTmp, maxIndex + lenTmp * domainTmp * 0.5)); + +for j = 1:1 % 右端点需要更正为maxIndex的长度 + leftTmp(j) = int64(max(1, maxIndex(j) - lenTmp * domainTmp * 0.5)); + rightTmp(j) = int64(min(lenTmp, maxIndex(j) + lenTmp * domainTmp * 0.5)); + localTmp = find(sampleData(leftTmp(j):rightTmp(j),2) > upperTmp); + wholeTmp = find(sampleData(:,2) > upperTmp); + rateTmp = length(localTmp) / length(wholeTmp); + edgeIndex = boundary(sampleData(:,1), sampleData(:,2), 1.0); + edge = [edgeIndex, sampleData(edgeIndex,1), sampleData(edgeIndex,2)]; + + % simple way + % upEdgeIndex = edge(:,3)>=2; + % upEdge = edge(upEdgeIndex,:); + + % upEdge: tradional way + edge = flipud(edge); % 水平轴上下翻转数组 + if isempty(edge) + fprintf("\tedge is empty!\n"); + break; + end + upEdgeStart = 1; + upEdgeEnd = length(edge); + edgeTmp = edge; + thresholdTrash = 5; + minEdgeInput = min(edgeTmp(:,2)); % 记录的是误差对应的输入 + maxEdgeInput = max(edgeTmp(:,2)); % 记录的是误差对应的输入 + + threshold = 2; + edgeTmp(edge(:,3) thresholdTrash) + upEdgeStart = max(1,k-1-thresholdTrash); + break + end + else + trash = 0; + end + end + % ready for upEdge's tail + trash = 0; + for k = length(edgeTmp):-1:1 + if(edgeTmp(k,2) == maxEdgeInput) + upEdgeEnd = min(k, length(edgeTmp)); + break + end + if(edgeTmp(k,end) ~= threshold) + trash = trash + 1; + if(trash > thresholdTrash) + upEdgeEnd = min(k+1+thresholdTrash, length(edgeTmp)); + break + end + else + trash = 0; + end + end + % after getting head & tail, we get upEdge + upEdge = edge(upEdgeStart:upEdgeEnd,:); % upEdge是保留了上半段轮廓线中小于阈值的部分 + + if rateTmp > 0.5 % 做判断是因为rateTmp大于0.5即说明其图形unstable，出现了畸变 + % fprintf("%s is unstable, rate = %g, maxError = %g\n", sampleFileName, rateTmp, maxTmp); + % draw pic and save pic + % figure(1); + % scatter(sampleData(:,1), sampleData(:,2), '.'); + % hold on; plot(edge(:,2), edge(:,3), '-r', 'LineWidth', 1); hold off; + % set(gca, 'Fontname', 'Times newman', 'Fontsize', 55); + % nameTmp = strcat("unstable_scatterEdge_", sampleFileName); + % % title(nameTmp, 'Interpreter', 'none'); + % % grid on; + % img =gcf; print(img, '-djpeg', '-r100', nameTmp); + + % update boundary error to [2, 10]。这一段其实可以放到前边去 + radio = (maxTmp - 2) / 8; + indexTmp = sampleData(:,2) >= 2; + sampleData(indexTmp, 2) = (sampleData(indexTmp,2) - 2) / radio + 2; + maxTmp = max(sampleData(indexTmp,2)); + edgeIndex = boundary(sampleData(:,1), sampleData(:,2), 1.0); + edge = [edgeIndex, sampleData(edgeIndex,1), sampleData(edgeIndex,2)]; + + % upEdge: tradional way + edge = flipud(edge); + if isempty(edge) + fprintf("\tedge is empty!\n"); + break; + end + upEdgeStart = 1; + upEdgeEnd = length(edge); + edgeTmp = edge; + thresholdTrash = 5; + minEdgeInput = min(edgeTmp(:,2)); + maxEdgeInput = max(edgeTmp(:,2)); + + threshold = 2; + edgeTmp(edge(:,3) thresholdTrash) + upEdgeStart = max(1,k-1-thresholdTrash); + break + end + else + trash = 0; + end + end + % tail + trash = 0; + for k = length(edgeTmp):-1:1 + if(edgeTmp(k,2) == maxEdgeInput) + upEdgeEnd = min(k, length(edgeTmp)); + break + end + if(edgeTmp(k,end) ~= threshold) + trash = trash + 1; + if(trash > thresholdTrash) + upEdgeEnd = min(k+1+thresholdTrash, length(edgeTmp)); + break + end + else + trash = 0; + end + end + % after head & tail + upEdge = edge(upEdgeStart:upEdgeEnd,:); % upEdge是保留了上半段轮廓线中小于阈值的部分 + + % draw pic and save pic + % figure(1); + % scatter(sampleData(:,1), sampleData(:,2), '.'); + % hold on; plot(upEdge(:,2), upEdge(:,3), '-r', 'LineWidth', 1); hold off; + % set(gca, 'Fontname', 'Times newman', 'Fontsize', 55); + % nameTmp = strcat("unstable_1_scatterEdge_", sampleFileName); + % % title(nameTmp, 'Interpreter', 'none'); + % % grid on; + % img =gcf; print(img, '-djpeg', '-r100', nameTmp); + else % do nothing, just for debug + % fprintf("%s is stable, rate = %g, maxError = %g\n", sampleFileName, rateTmp, maxTmp); + % draw pic and save pic + % figure(1); + % scatter(sampleData(:,1), sampleData(:,2), '.'); + % hold on; plot(upEdge(:,2), upEdge(:,3), '-r', 'LineWidth', 1); hold off; + % set(gca, 'Fontname', 'Times newman', 'Fontsize', 55); + % nameTmp = strcat("stable_scatterEdge_", sampleFileName); + % % title(nameTmp, 'Interpreter', 'none'); + % % grid on; + % img =gcf; print(img, '-djpeg', '-r100', nameTmp); + end + + % save upEdge data to file whose format like upEdge_Bspline3.txt + if(~isempty(edge)) + % fprintf("\tupEdgeFileName: %s with %d lines\n", upEdgeFileName, length(upEdge)); + save(upEdgeFileName, 'upEdge', '-ascii', '-double'); % the 2nd parameter mean the data to save to file + else + fprintf("\tedge is empty\n"); + end + + % if(length(maxIndex) > 1) + % fprintf("\tleftTmp(%d) = %d, rightTmp(%d) = %d\n", j, leftTmp(j), j, rightTmp(j)); + % end +end +if maxTmp > 2 + indexTmp = sampleData(:,2) >= 2; + meanTmp = mean(sampleData(indexTmp, 2)); +else + meanTmp = mean(sampleData(:, 2)); +end + +% compute the threshold for deviding intervals +% thresholdTmp = (maxTmp + meanTmp) / 2; +% fprintf("%s: max = %g, average = %g, threshold = %d\n", sampleFileName, maxTmp, meanTmp, thresholdTmp); + +% 去除误差小于0.5的点，然后绘图。感觉区别不大。下面这段代码就是示例代码，仅仅是绘图，并没有保存数据或图片 +% indexTmp = sampleData(:,2)>0.5; +% sampleDataTmp1 = (sampleDataTmp(indexTmp,:)); +% edgeIndex = boundary(sampleDataTmp1(:,1), sampleDataTmp1(:,2), 1.0); +% edge = [edgeIndex, sampleDataTmp1(edgeIndex,1), sampleDataTmp1(edgeIndex,2)]; +% figure(2); +% scatter(sampleDataTmp1(:,1), sampleDataTmp1(:,2), '.'); +% hold on; plot(edge(:,2), edge(:,3), '-r', 'LineWidth', 1); hold off; diff --git a/src/main.cpp b/src/main.cpp index 1008febc1d9f66c355ee3028699ca1e5f716d49d..a8a622ad4c2ec9bdf7d5257fc124624daf06be83 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -22,6 +22,41 @@ using std::vector; #ifndef DOUBLE_PRECISION #define DOUBLE_PRECISION 17 #endif + +map> benchmarkThresholds = { + {"Bsplines3", {0.8}}, + {"exp1x", {1}}, + {"exp1x_log", {1}}, + {"intro_example", {1}}, + {"logexp", {10}}, + {"NMSEexample31", {2}}, + {"NMSEexample310", {2}}, + {"NMSEexample34", {2.1}}, + {"NMSEexample35", {1}}, + {"NMSEexample36", {2}}, + {"NMSEexample37", {0.5}}, + {"NMSEexample38", {10}}, + {"NMSEexample39", {5}}, + {"NMSEproblem331", {1}}, + {"NMSEproblem333", {1.5}}, + {"NMSEproblem334", {200}}, + {"NMSEproblem336", {1}}, + {"NMSEproblem337", {1}}, + {"NMSEproblem341", {2}}, + {"NMSEproblem343", {1}}, + {"NMSEproblem344", {2}}, + {"NMSEproblem345", {2}}, + {"NMSEsection311", {2}}, + {"predatorPrey", {2}}, + {"sine", {1}}, + {"sineorder3", {3}}, + {"sqroot", {0.6}}, + {"sqrt_add", {1.5}}, + {"test05_nonlin1_r4", {1.2}}, + {"test05_nonlin1_test2", {1.2}}, + {"verhulst", {2}}, +}; + //===----------------------------------------------------------------------===// // Main driver code. //===----------------------------------------------------------------------===// @@ -46,7 +81,7 @@ int main() cout << "open tmpResult.txt failed"; exit(EXIT_FAILURE); } - int getlineCount = 0; + // int getlineCount = 0; fprintf(stderr, GREEN "ready> " RESET); string inputStr = ""; @@ -54,9 +89,9 @@ int main() while (getline(cin, inputStr)) // read line from keyboard input { // only rewrite - getlineCount++; - if (getlineCount == 35 || getlineCount == 36 || getlineCount < 0) - continue; + // getlineCount++; + // if (getlineCount == 35 || getlineCount == 36 || getlineCount < 0) + // continue; if (inputStr == "exit;" || inputStr == "quit;" || inputStr == "exit" || inputStr == "quit") { @@ -79,11 +114,16 @@ int main() auto pos = benchMarkData.find(inputStr); bool isBenchMark = false; string name; + vector thresholds; if (pos != benchMarkData.end()) { name = inputStr; inputStr = pos->second.begin()->first; - cout << inputStr << endl; + cout << "expression := " << inputStr << endl; + + thresholds = benchmarkThresholds.at(name); + // fmt::print("thresholds := {}\n", thresholds); + isBenchMark = true; } @@ -93,8 +133,8 @@ int main() auto originExpr = ParseExpressionFromString(inputStr); vector vars; getVariablesFromExpr(originExpr, vars); + int dimension = vars.size(); // cout << vars.size() << endl; - fprintf(stderr, GREEN "interval> " RESET); // after input completed, start timing @@ -114,20 +154,21 @@ int main() } // For temporary use only . It will be replaced by geneBoundaryData and geneIntervalData - ofstream ofs; - ofs.open("intervalData.txt", ios::out); - for (size_t i = 0; i < intervals.size(); i++) - { - if (i == intervals.size() - 1) - { - ofs << intervals.at(i); - } - else - { - ofs << intervals.at(i) << " "; - } - } - ofs.close(); + // ofstream ofs; + // ofs.open("intervalData.txt", ios::out); + // TODO-DONE: (but no need do it): judge if intervals.size() == dimension + // for (int i = 0; i < dimension; i++) + // { + // if (i == dimension - 1) + // { + // ofs << intervals.at(i); + // } + // else + // { + // ofs << intervals.at(i) << " "; + // } + // } + // ofs.close(); // Read scale from keyboard // fprintf(stderr, GREEN "scale> " RESET); @@ -137,15 +178,15 @@ int main() // Default scale setting according to the variables' size int sampleScale; - if (vars.size() == 1) + if (dimension == 1) { sampleScale = 500000; } - else if (vars.size() == 2) + else if (dimension == 2) { sampleScale = 1024; } - else if (vars.size() == 3) + else if (dimension == 3) { sampleScale = 256; } @@ -154,7 +195,7 @@ int main() sampleScale = 10; } vector scales; - for (size_t i = 0; i < vars.size(); i++) + for (int i = 0; i < dimension; i++) { scales.push_back(sampleScale); } @@ -171,43 +212,68 @@ int main() string mkdirCommand = "mkdir -p srcTest/" + uniqueLabel + " outputs/" + uniqueLabel; system(mkdirCommand.c_str()); - auto funcNameOrigin = geneOriginCodeKernel(inputStr, vars, uniqueLabel, "origin"); - // auto funcNameOrigin = geneOriginCode(inputStr, uniqueLabel, "origin"); + auto funcNameOrigin = geneExprCodeKernel(inputStr, vars, uniqueLabel, "origin"); + // auto funcNameOrigin = geneExprCode(inputStr, uniqueLabel, "origin"); // auto funcNameHerbie = geneHerbieCode(inputStr, uniqueLabel, "herbie"); + auto funcNameHerbie = geneHerbieCode(uniqueLabel); // auto funcNameDaisy = geneDaisyCode(inputStr, uniqueLabel, "daisy"); auto funcNameMpfr = geneMpfrCode(inputStr, uniqueLabel, vars); - // TODO: pick the best from origin, herbie, daisy + // TODO: improve pickTheBest to support more suffix // pickTheBest(uniqueLabel, 0, 1, 100); - - cout << BLUE << "main: start testError for origin: " << inputStr << RESET << endl; + vector suffixSet = {"origin", "herbie"}; + auto suffix = pickTheBest(uniqueLabel, suffixSet, intervals, scales); + if (suffix == "origin") + { + inputStr = funcNameOrigin; + } + else if (suffix == "herbie") + { + inputStr = funcNameHerbie; + } + else + { + fprintf(stderr, "ERROR: main: we do not support the suffix %s now\n", suffix.c_str()); + exit(EXIT_FAILURE); + } + // cout << BLUE << "main: start testError for origin: " << inputStr << RESET << endl; // auto timeTmp1 = std::chrono::high_resolution_clock::now(); // geneSampleData(); - auto infoTmp = testError(uniqueLabel, "origin", intervals, scales); + // auto infoTmp = testError(uniqueLabel, "origin", intervals, scales); // auto timeTmp2 = std::chrono::high_resolution_clock::now(); - cout << BLUE << "main: ending testError for origin: " << inputStr << RESET << endl; + // cout << BLUE << "main: ending testError for origin: " << inputStr << RESET << endl; // std::chrono::duration testError_seconds = timeTmp2 - timeTmp1; // cout << BLUE << "testError time: " << testError_seconds.count() << "s" << RESET << endl; // fprintf(stderr, GREEN "ready> " RESET); // continue; - - geneBoundaryData(uniqueLabel, "origin"); // matlab - - vector intervalsTemp; - vector threholdsTemp; - geneIntervalData(uniqueLabel, intervalsTemp, threholdsTemp); + + // auto infoTmp = testError(uniqueLabel, "origin", intervals, scales); // TODO-done: 完善origin误差测试 // WHY do it? Because no pickTheBest before. So, have to use testError to get the sample data. + auto upEdgeFileName = geneBoundaryData(uniqueLabel, suffix); // sample data == matlab ==> upEdge data // TODO: support multiple dimension + cout << "upEdgeFileName: " << upEdgeFileName << "\n"; + vector upEdgeFileNames; + upEdgeFileNames.push_back(upEdgeFileName); + auto intervalData = getIntervalData(upEdgeFileNames, thresholds); + fmt::print("after regime, we have {} intervals: {}\n", intervalData.size(), intervalData); + // fprintf(stderr, GREEN "ready> " RESET); + // continue; cout << "=-=-=-=-=-=-=-=-=-=-=-=-= rewrite start =-=-=-=-=-=-=-=-=-=-=-=-=" << endl; // auto timeTmp3 = std::chrono::high_resolution_clock::now(); - auto exprInfoVector = rewrite(inputStr, uniqueLabel); + auto exprInfoVector = rewrite(inputStr, uniqueLabel, intervalData); // auto timeTmp4 = std::chrono::high_resolution_clock::now(); // std::chrono::duration rewrite_seconds = timeTmp4 - timeTmp3; // cout << BLUE << "rewrite time: " << rewrite_seconds.count() << " s" << RESET << endl; // fprintf(stderr, GREEN "ready> " RESET); // continue; cout << "=-=-=-=-=-=-=-=-=-=-=-=-= rewrite end =-=-=-=-=-=-=-=-=-=-=-=-=" << endl; + auto funcNameFinal = geneFinalCodeKernel(inputStr, uniqueLabel, exprInfoVector, vars); + cout << "=-=-=-=-=-=-=-=-=-=-=-=-= test final code's error and performance start =-=-=-=-=-=-=-=-=-=-=-=-=\n"; + auto infoTmp = testError(uniqueLabel, "final", intervals, scales); + infoTmp.performance = testPerformance(uniqueLabel, "final", intervals); + cout << "performance: " << infoTmp.performance << "\n\n"; + cout << "=-=-=-=-=-=-=-=-=-=-=-=-= test final code's error and performance end =-=-=-=-=-=-=-=-=-=-=-=-=\n"; // write the results to file // auto results = exprAutoWrapper(inputStr, intervals, scales); @@ -233,7 +299,7 @@ int main() // continue; // write the results to file - fout << "-------------------------------------NO." << getlineCount << ": " << inputStr << endl; + // fout << "-------------------------------------NO." << getlineCount << ": " << inputStr << endl; for (size_t i = 0; i < results.size(); i++) { auto &result = results.at(i); diff --git a/src/tools.cpp b/src/tools.cpp index ac0dc82e6b81e1accce98d91adeaf3e3a153628d..6b2b892d307fdab3342d1812b7a60a16f0849884 100644 --- a/src/tools.cpp +++ b/src/tools.cpp @@ -4,9 +4,11 @@ #include "basic.hpp" #include "geneCode.hpp" #include "color.hpp" +#include "devideUpEdgeData.hpp" #include #include #include +#include #include #include #include @@ -17,6 +19,7 @@ #include #include #include +#include using std::cout; using std::endl; @@ -387,7 +390,9 @@ exprInfo testError(string uniqueLabel, string suffix, const vector &inte } } string fileNameKernel = prefix + "__" + middle + "_" + suffix; - string testName = "./outputs/" + uniqueLabel + "/" + fileNameKernel + "_error.txt"; + namespace fs = std::filesystem; + string currentPath = fs::current_path(); + string testName = currentPath + "/outputs/" + uniqueLabel + "/" + fileNameKernel + "_error.txt"; string number[3] = {"One", "Two", "Three"}; string scriptName = "./detectError" + number[size - 1] + "FPEDParallel.sh"; string commandStr = scriptName + " " + uniqueLabel; @@ -397,7 +402,9 @@ exprInfo testError(string uniqueLabel, string suffix, const vector &inte } commandStr = commandStr + " " + prefix + " " + middle + " " + suffix; - cout << "fileNameKernel: " << fileNameKernel << "\ncommand: " << commandStr << "\ntestName: " << testName << endl; + // cout << "fileNameKernel: " << fileNameKernel << "\n"; + cout << "command: " << commandStr << "\n"; + // cout << "testName: " << testName << "\n"; char command[200] = {0}; strcat(command, commandStr.c_str()); system(command); @@ -405,7 +412,7 @@ exprInfo testError(string uniqueLabel, string suffix, const vector &inte std::ifstream ifs(testName, std::ios::in); if(ifs.eof()) { - std::cout << "is null" << std::endl; + cout << "is null" << "\n"; ifs.close(); tempError.intervals = intervals; @@ -437,53 +444,164 @@ exprInfo testError(string uniqueLabel, string suffix, const vector &inte } else { - fprintf(stderr, "WRONG: rewrite: the intervalTmp's dimension is %ld, which we don't support now.\n", size); + fprintf(stderr, "WRONG: testError: the intervalTmp's dimension is %ld, which we don't support now.\n", size); exit(EXIT_FAILURE); } return tempError; } -// TODO: implement -// call matlab to generate the boundaryData to file -void geneBoundaryData(string uniqueLabel, string suffix) +double testPerformance(string uniqueLabel, string suffix, const vector &intervals) { - suffix = "origin"; - string filename = "./outputs/" + uniqueLabel + "/expr_" + uniqueLabel + "_" + suffix + "_boudary.txt"; - std::ofstream ofs; - ofs.open(filename, std::ios::out); - ofs << "123" << std::endl; - ofs.close(); - std::cout << "generate boundaryData file:" << filename << std::endl; -} + size_t size = intervals.size() / 2; + string prefix = "expr_" + uniqueLabel; + string fileNameKernel = prefix + "_" + suffix; + if (size < 4) + { + vector params; + for(const auto &interval : intervals) + { + auto paraTmp = fmt::format("{}", interval); + params.push_back(paraTmp); + } -// TODO: improve -// generate the interval info to file -void geneIntervalData(string uniqueLabel, vector &intervals, vector &threholds) -{ - // filename = "expr_" + uniqueLabel + "_interval.txt"; - for (size_t j = 1; j < threholds.size(); j++) + string number[3] = {"One", "Two", "Three"}; + string scriptName = "./testPerformance" + number[size - 1] + ".sh"; + string commandStr = "taskset -c 0 " +scriptName + " " + uniqueLabel + " " + suffix; + for(const auto & param : params) + { + commandStr = commandStr + " " + param; + } + commandStr = "cd srcTest; " + commandStr + "; cd - > /dev/null;"; + // cout << "fileNameKernel: " << fileNameKernel << "\n"; + cout << "command: " << commandStr << "\n"; + char command[200] = {0}; + strcat(command, commandStr.c_str()); + system(command); + } + else if (size == 4) { - intervals.push_back(to_string(j * 2)); - intervals.push_back(to_string(j * 2 + 1)); - threholds.push_back(j + 0.0002); + fprintf(stderr, "WRONG: testPerformance: the intervalTmp's dimension is %ld, which we don't support now.\n", size); + exit(EXIT_FAILURE); } - for (size_t i = 0; i < threholds.size(); i++) + else { - string fileName; - string prefix = "./outputs/" + uniqueLabel + "/expr_" + uniqueLabel; - string middle = "_" + intervals.at(i * 2) + "_" + intervals.at(i * 2 + 1) + "_" + to_string(threholds.at(i)); - string suffix = "_interval.txt"; - fileName = prefix + middle + suffix; - std::ofstream ofs; - ofs.open(fileName, std::ios::out); - ofs << "123" << std::endl; - ofs.close(); - std::cout << "generate boudaryDate file:" << fileName << std::endl; + fprintf(stderr, "WRONG: testPerformance: the intervalTmp's dimension is %ld, which we don't support now.\n", size); + exit(EXIT_FAILURE); } + + namespace fs = std::filesystem; + string currentPath = fs::current_path(); + string testName = currentPath + "/outputs/" + uniqueLabel + "/" + fileNameKernel + "_performance.txt"; + std::ifstream ifs(testName, std::ios::in); + double tempCycles; + if(ifs.eof()) + { + std::cout << "is null" << std::endl; + ifs.close(); + tempCycles = 0; + } + else + { + std::string lineStr; + std::getline(ifs, lineStr); + std::istringstream ss(lineStr); + ss >> tempCycles; + } + + return tempCycles; } -//生成区间数据 +// TODO: support multiple dimension +// call matlab to generate the upEdgeData to file +string geneBoundaryData(string uniqueLabel, string suffix) +{ + char *currentPath; + if((currentPath = getcwd(NULL, 0)) == NULL) + { + perror("getcwd error"); + } + string currentPathStr(currentPath); + string outputFilesPath = currentPathStr + "/outputs/" + uniqueLabel + "/"; + string sampleFileName = outputFilesPath + "sample_" + uniqueLabel + "_" + suffix + ".txt"; + string upEdgeFileName = outputFilesPath + "upEdge_" + uniqueLabel + "_" + suffix + ".txt"; + string matlabExecutable = "/home/xyy/helloBro/softwares/matlab/R2020a/bin/matlab"; // set to your own matlab executable + string matlabCommands = " -batch \"sampleFileName=\'" + sampleFileName + "\'; upEdgeFileName = \'" + upEdgeFileName + "\'; run " + currentPathStr + "/src/getUpEdge\""; // DO NOT need to add the suffix ".m" for matlab file! + string commandStr = matlabExecutable + matlabCommands; + cout << "Matlab command: " << commandStr << "\n"; + // cout << "length of commandStr: " << commandStr.size() << "\n"; + char command[512] = {0}; + strcat(command, commandStr.c_str()); + system(command); + free(currentPath); + // std::cout << "generate upEdge file:" << upEdgeFileName << std::endl; + + return upEdgeFileName; +} + +// string geneIntervalData1D(string upEdgeFileName, string uniqueLabel, double threshold) // generate the interval info to file at the target dimension +// { +// // init +// double thresholdCombine = 101; +// string prefix = "./outputs/" + uniqueLabel + "/expr_" + uniqueLabel; +// string suffix = "_interval.txt"; +// // call devideUpEdgeData for each dimension +// string middle = "_" + to_string(threshold); +// string intervalFileName = prefix + middle + suffix; +// cout << "generate interval file:" << intervalFileName << endl; +// devideUpEdgeData(upEdgeFileName, intervalFileName, threshold, thresholdCombine); +// return intervalFileName; +// } + +// vector getIntervalData1D(string filename) // read filename to get vector +// { +// vector intervalData; +// std::ifstream ifs; +// ifs.open(filename, std::ios::in); +// if (!ifs.is_open()) +// { +// std::cout << "ERROR: getIntervalData: fail to open " << filename << std::endl; +// exit(EXIT_FAILURE); +// } +// string currentLine; +// while (getline(ifs, currentLine, '\n')) +// { +// // cout << currentLine << endl; +// std::istringstream iss(currentLine); +// string tmp; +// while(getline(iss, tmp, ' ')) +// { +// intervalData.push_back(atof(tmp.c_str())); +// } +// } +// ifs.close(); +// return intervalData; +// } + +// vector geneIntervalData(vector upEdgeFileNames, string uniqueLabel, vector &thresholds) // according the upEdgeFile data, generate the interval info into file +// { +// // init +// double thresholdCombine = 101; +// vector intervalFilenames; +// string prefix = "./outputs/" + uniqueLabel + "/expr_" + uniqueLabel; +// string suffix = "_interval.txt"; +// // call devideUpEdgeData for each dimension +// for (size_t i = 0; i < thresholds.size(); i++) +// { +// auto &upEdgeFileName = upEdgeFileNames.at(i); +// auto &threshold = thresholds.at(i); +// string middle = "_" + to_string(threshold); +// string intervalFileName = prefix + middle + suffix; +// cout << "generate interval file:" << intervalFileName << endl; +// intervalFilenames.push_back(intervalFileName); +// devideUpEdgeData(upEdgeFileName, intervalFileName, threshold, thresholdCombine); +// } +// return intervalFilenames; +// } + +// may useless, ready to comment it +// read filename to get vector>. +// A line is a vector. lines are vector<~> (~ is vector) vector> getIntervalData(string filename) { vector> intervalData; @@ -495,52 +613,48 @@ vector> getIntervalData(string filename) exit(EXIT_FAILURE); } - string buf; - while (getline(ifs, buf)) + string currentLine; + while (getline(ifs, currentLine, '\n')) { + // cout << currentLine << endl; vector lineData; - // cout << buf << endl; - int flag = 1; - vector index; - for (size_t i = 0; i < buf.length(); i++) + std::istringstream iss(currentLine); + string tmp; + while(getline(iss, tmp, ' ')) { - if (buf[i] == ' ') - { - flag += 1; - index.push_back(i); - } - } - for (size_t i = 0; i < index.size(); i++) - { - if (i == 0) - { - string data = buf.substr(0, index.at(i)); - lineData.push_back(atof(data.c_str())); - } - else - { - int n = index.at(i) - index.at(i - 1) - 1; - string data = buf.substr(index.at(i - 1) + 1, n); - lineData.push_back(atof(data.c_str())); - } + lineData.push_back(atof(tmp.c_str())); } - int pos = index.back() + 1; - int n = buf.length() - 1 - index.at(index.size() - 1); - string last = buf.substr(pos, n); - lineData.push_back(atof(last.c_str())); intervalData.push_back(lineData); } ifs.close(); return intervalData; } -// call exprAuto to rewrite the input expression -vector rewrite(string exprStr, string uniqueLabel) +// according to the corresponding intervalData file in each dimension, generate the interval data, and then generate all permutation of all dimensions. +vector> getIntervalData(vector upEdgeFileNames, vector &thresholds) { - // string filename = "expr_" + uniqueLabel + "_interval.txt"; - string filename = "./intervalData.txt"; // TODO: get the filename from uniqueLabel + vector> intervalDataMultiDim; + int dimension = thresholds.size(); + // auto thresholdCombine = 0; + for (int i = 0; i < dimension; i++) // iterate all the dimensions + { + // call devideUpEdgeData for each dimension: matlab upEdge ==> interval vector + auto &upEdgeFileName = upEdgeFileNames.at(i); + auto &threshold = thresholds.at(i); + auto intervalData1D = devideUpEdgeData(upEdgeFileName, threshold); // TODO: add the parameter thresholdCombine + // push_back interval vector + intervalDataMultiDim.push_back(intervalData1D); + } + // call permuteMultiVec to get all permutation + auto results = permuteMultiVec(intervalDataMultiDim); + return results; +} - auto intervalData = getIntervalData(filename); +// call exprAuto to rewrite the input expression +vector rewrite(string exprStr, string uniqueLabel, vector> &intervalData) +{ + // string filename = "expr_" + uniqueLabel + "_interval.txt"; // useless + // string filename = "./intervalData.txt"; // TODO-done: get the filename from uniqueLabel. Now we get intervalData from input parameters directly. auto tempExpr = ParseExpressionFromString(exprStr); @@ -586,7 +700,7 @@ vector rewrite(string exprStr, string uniqueLabel) // generate function code and test error string suffixTmp = suffix + std::to_string(j); - geneOriginCode(newExpr, uniqueLabel, suffixTmp); + geneExprCode(newExpr, uniqueLabel, suffixTmp); // auto timeTmp1 = std::chrono::high_resolution_clock::now(); auto tempError = testError(uniqueLabel, suffixTmp, intervalTmp, scales); // auto timeTmp2 = std::chrono::high_resolution_clock::now(); @@ -633,7 +747,8 @@ vector rewrite(string exprStr, string uniqueLabel) return exprInfoVector; } -vector> zuhe(vector> vec) +// permute multiple vector +vector> permuteMultiVec (vector> vec) { vector> suzu; size_t row = vec.size(); diff --git a/srcTest/Makefile b/srcTest/Makefile new file mode 100644 index 0000000000000000000000000000000000000000..03466e1355db238e66eee4c671d2dba7e40ea39a --- /dev/null +++ b/srcTest/Makefile @@ -0,0 +1,33 @@ +export PROJECT_NAME = exprAuto + +CC = gcc +CPP = g++ +INCLUDE = -I../include -I../includeTEST +LIBS = +ECHO = printf + +# $(info $(CFLAGS) ) +override CFLAGS += -O3 -g -Wall -Wextra -Wpedantic -Wno-unused-function -fdiagnostics-color=always +# $(info $(CFLAGS) ) + +# default: dataClean.exe + +dataClean.exe: dataClean.o + @$(ECHO) "\033[1;32mBuilding $@ \n\033[0m" + $(CC) -o $@ $^ $(LIBS) + +dataClean.o: dataClean.c + @$(ECHO) "\033[1;32mBuilding $@ \n\033[0m" + $(CC) -o $@ -c $< $(CFLAGS) $(INCLUDE) + +gccPerformanceTest.o: gccPerformanceTest.c + @$(ECHO) "\033[1;32mBuilding $@ \n\033[0m" + $(CC) -o $@ -c $< $(CFLAGS) $(INCLUDE) + +binary.o: binary.c + @$(ECHO) "\033[1;32mBuilding $@ \n\033[0m" + $(CC) -o $@ -c $< $(CFLAGS) $(INCLUDE) + +.PHONY: clean +clean: + rm -f *.o *.exe \ No newline at end of file diff --git a/srcTest/binary.c b/srcTest/binary.c new file mode 100644 index 0000000000000000000000000000000000000000..d0f68d151ee390e8b499dbf38f8efd4a7aa1d656 --- /dev/null +++ b/srcTest/binary.c @@ -0,0 +1,210 @@ +#include +#include + +#include "binary.h" + +/* copy: copy 'from' into 'to'; assume to is big enough */ +void copy(char to[], char from[]) +{ + int i; + i = 0; + while((to[i] = from[i]) != '\0') + ++i; +} + +_UL ftoid(double x) { return *((_UL *)(&x)); } + +double itofd(_UL x) { return *((double *)(&x)); } + +void getbinary(_UL x, int *result) +{ + int i; + _UL temp; + + temp = x; + for(i = 0; i <= 63; i++) + { + result[i] = temp & 0x0000000000000001; + temp = temp >> 1; + } +} + +_UL getUL(int *result) +{ + int i; + _UL temp; + + temp = 0; + for(i = 63; i >= 0; i--) + { + // result[i] = temp & 0x0000000000000001; + // temp = temp >> 1; + temp = temp << 1; + temp = temp + result[i]; + } + + return temp; +} + +void printbinary(int *result) +{ + int i; + + for(i = 63; i >= 0; i--) + { + printf("%d", result[i]); + } + printf("\n"); +} + +void fprintbinary(FILE *file, int *result) +{ + int i; + + for(i = 63; i >= 0; i--) + { + fprintf(file, "%d", result[i]); + } + fprintf(file, "\n"); +} + +void binaryshow(double x) +{ + _UL temp; + int sign, e, i, result[64]; + + temp = ftoid(x); + getbinary(temp, result); + sign = result[63]; + e = 0; + for(i = 10; i >= 0; i--) + { + e = (e << 1) + result[52 + i]; + } + + if(e == 0) + { + printf(" 0."); + } + else + { + e = e - (1 << 10) + 1; + if(sign == 0) + { + printf("+1."); + } + else + { + printf("-1."); + } + } + for(i = 51; i >= 0; i--) + { + printf("%d", result[i]); + } + printf("b%d\n", e); +} + +_UL computeULPDiff(double x1, double x2) +{ + _UL temp1, temp2; + _UL t; + + temp1 = ftoid(x1); + temp2 = ftoid(x2); + + temp1 = temp1 & 0x000fffffffffffff; + temp2 = temp2 & 0x000fffffffffffff; + + if(temp1 > temp2) + { + t = temp1 - temp2; + } + else + { + t = temp2 - temp1; + } + + return t; +} + +int computeAccurateBit(double x1, double x2) +{ + int y1[64], y2[64], i, count; + _UL temp1, temp2; + + temp1 = ftoid(x1); + temp2 = ftoid(x2); + getbinary(temp1, y1); + getbinary(temp2, y2); + + count = 0; + for(i = 51; i >= 0; i--) + { + if(y1[i] != y2[i]) + { + break; + } + count++; + } + + return count; +} + +int jiou(int x) +{ + if(x % 2 == 0) + return 1; + else + return 0; +} + +double computeULP(double y) +{ + double x = 0, res = 0, powermin = 0, powermax = 0, powermiddle = 0; + int expmin = 0, expmax = 0, expmiddle = 0, jioupanduan = 0; + x = fabs(y); + + // printf("res=%.8f\n", pow(2, logb(x)+1-53)); + + if(x < pow(2, -1021)) + res = pow(2, -1074); + else if(x > (1 - pow(2, -53)) * pow(2, 1024)) + res = pow(2, 971); + else + { + powermin = pow(2, -1021); + expmin = -1021; + powermax = pow(2, 1024); + expmax = 1024; + } + + while(expmax - expmin > 1) + { + jioupanduan = jiou(expmin + expmax); + + if(jioupanduan == 1) + expmiddle = (expmax + expmin) / 2; + else + expmiddle = (expmax + expmin + 1) / 2; + + powermiddle = pow(2, expmiddle); + + if(x >= powermiddle) + { + powermin = powermiddle; + expmin = expmiddle; + } + else + { + powermax = powermiddle; + expmax = expmiddle; + } + + if(x == powermin) + res = pow(2, expmin - 53); + else + res = pow(2, expmin - 52); + } + return res; +} diff --git a/srcTest/dataClean.c b/srcTest/dataClean.c new file mode 100644 index 0000000000000000000000000000000000000000..7c3a3c453706e41e3baf2e07f32114c7230b7fb2 --- /dev/null +++ b/srcTest/dataClean.c @@ -0,0 +1,55 @@ +#include +#include + +#define MAX_LEN 255 +#define PERCENT 0.0005 +#ifndef RUNTIME +#define RUNTIME 10000 +#endif + +int main(int argc, char *argv[]) +{ + int row = 0, i = 0; + long cycle = 0; + FILE *data; + long time_base[RUNTIME]; + long time_base_sum = 0; + int start, end; + + if(argc == 2) + { + data = fopen(argv[1], "r"); + } + else + { + printf("ERROR: Invalid arguments. argc should be 2, not %d\n", argc); + exit(EXIT_FAILURE); + } + + if(!data) + { + printf("Error in opening file.\r\n"); + exit(1); + } + + // printf("----Starting data clean: %d----\n", RUNTIME); + row = 0; + while(fscanf(data, "%ld\n", &time_base[row]) != EOF) + { + row++; + } + start = row * PERCENT; + end = row - start; + // printf("the data size is %d, we choose: %d to %d\n", row, start, end); + for(i = start; i < end; i++) + { + time_base_sum += time_base[i]; + } + // printf("----End data clean----\n"); + + cycle = time_base_sum / (end - start); + // printf("the performance result (cycle): %ld\n", cycle); + printf("%ld\n", cycle); + + fclose(data); +} diff --git a/srcTest/gccPerformanceTest.c b/srcTest/gccPerformanceTest.c new file mode 100644 index 0000000000000000000000000000000000000000..fdf7f9ba8b5ab792682f5692b2f1aa793ed97ec2 --- /dev/null +++ b/srcTest/gccPerformanceTest.c @@ -0,0 +1,105 @@ +#include +#include +#include +// #include +#include +#include + +#include "common.h" + +#ifndef FUNCNAME +#define FUNCNAME sin +#endif +#ifndef RUNTIME +#define RUNTIME 10000 +#endif + +double FUNCNAME(double); +/* copy: copy 'from' into 'to'; assume to is big enough */ +void copy(char to[], char from[]); + +// for performance test +extern inline __attribute__((always_inline)) unsigned long rdtscp() +{ + unsigned long a, d, c; + + __asm__ volatile("rdtscp" : "=a"(a), "=d"(d), "=c"(c)); + + return (a | (d << 32)); +} + +int main(int argc, char *argv[]) +{ + int i; + unsigned long times[RUNTIME]; + unsigned long sum; + uint64_t start, end; + double inputdata[RUNTIME]; + double result[RUNTIME], sumResult = 0, test = 1; + double data_start = 0, data_end = 1; + FILE *stream_time; + char outputFile[64] = "FUNCNAME_time.txt"; // the name of the outputFile. + + srand((unsigned)time(NULL)); + if(argc == 2) + { + copy(outputFile, argv[1]); + // printf("Performance test function: %s\n", outputFile); + // printf("outputFile is %s.\ntest from %f, end at %f by default\n", outputFile, data_start, data_end); + } + if(argc == 4) + { + copy(outputFile, argv[1]); + data_start = atof(argv[2]); + data_end = atof(argv[3]); + // printf("Performance test function: %s\n", outputFile); + // printf("outputFile is %s.\ntest from %f, end at %f\n", outputFile, data_start, data_end); + } + stream_time = fopen(outputFile, "w"); + if(stream_time == (FILE *)0) + { + printf("open file error!\n"); + exit(1); + } + + // printf("----Starting generating %d test data----\n", RUNTIME); + for(i = 0; i < RUNTIME; i++) + { + inputdata[i] = (double)rand() / ((double)RAND_MAX + 1) * ((data_end) - (data_start)) + (data_start); + } + // printf("----End generating %d test data----\n", RUNTIME); + + // printf("----Starting performance test----\n"); + for(i = 0; i < 100; i++) + { + FUNCNAME(inputdata[i]); + } + for(i = 0; i < RUNTIME; i++) + { + start = rdtscp(); + result[i] = FUNCNAME(inputdata[i]); + end = rdtscp(); + // printf("the %d result is %f\n", i, result[i]); + result[i] = result[i] + test; + times[i] = end - start; + // printf("time is %llu clock cycles\n", times[i]); + } + + // printf("----End performance test----\n"); + + sum = 0; + for(i = 0; i < RUNTIME; i++) + { + fprintf(stream_time, "%lu\n", times[i]); + sum = sum + times[i]; + sumResult = sumResult + result[i]; + } + sum = sum / RUNTIME; + printf("sumResult is %f\n", sumResult); + // printf("the performance result to be clean: %lu\n", sum); + if(fclose(stream_time) == EOF) + { + printf("close file error!\n"); + exit(1); + } +} diff --git a/srcTest/test1paramFPEDParallel.c b/srcTest/test1paramFPEDParallel.c index 01627c536fc04dbcd0690692e9a520e59e73c239..1b8c18aeda1211f07c34039512074cd884212b58 100644 --- a/srcTest/test1paramFPEDParallel.c +++ b/srcTest/test1paramFPEDParallel.c @@ -21,6 +21,7 @@ struct errorInfo { #define TESTNUMX0 500000 // #define FP // #define DEBUG +#define ERRORFILE double EXPRESSIONMPFR(double, mpfr_t); double EXPRESSIONMINE(double); @@ -39,7 +40,7 @@ int computeResult1param(double x0, mpfr_t mpfrResult) return status; } -struct errorInfo test1FPEDparamParallel(DL x0Start, DL x0End, unsigned long int testNumX0, const char* fileNameKernel, int myid, int i0StartLocal, int i0EndLocal) { +struct errorInfo test1FPEDparamParallel(DL x0Start, DL x0End, unsigned long int testNumX0, const char* uniqueLabel, const char* fileNameKernel, int myid, int i0StartLocal, int i0EndLocal) { DL maxInputX0; int i0; // int flag; @@ -49,24 +50,26 @@ struct errorInfo test1FPEDparamParallel(DL x0Start, DL x0End, unsigned long int mpfr_t mpfrOrcle, mpfrResult; mpfr_inits2(PRECISION, mpfrOrcle, mpfrResult, (mpfr_ptr) 0); - // file - // char *directory = "outputs"; - // char *suffix = "sample"; - // char *fileNameSample; - // FILE *f; - // fileNameSample = (char *) calloc(strlen(fileNameKernel) + strlen(suffix) + 128, sizeof(char)); - // sprintf(fileNameSample, "./%s/%s_%s_%d.txt", directory, fileNameKernel, suffix, myid); - // // printf("%s\n", fileNameSample); - // if ((f = fopen(fileNameSample, "w")) == NULL) { - // printf("Error opening file %s.\n", fileNameSample); - // exit(0); - // } + // write error data to file + #ifdef ERRORFILE + char *directory = "./outputs"; + char *suffix = "sample"; + char *fileNameSample; + FILE *f; + fileNameSample = (char *) calloc(strlen(directory) + strlen(uniqueLabel) + strlen(fileNameKernel) + strlen(suffix) + 128, sizeof(char)); + sprintf(fileNameSample, "%s/%s/%s_%s_%02d.txt", directory, uniqueLabel, fileNameKernel, suffix, myid); + // printf("%s\n", fileNameSample); + if ((f = fopen(fileNameSample, "w")) == NULL) { + printf("Error opening file %s.\n", fileNameSample); + exit(0); + } // printf("test expression: %s\n", STR2(EXPRESSION)); // printf("x0Start : %lg 0x%016lx\nx0End : %lg 0x%016lx\n", x0Start.d, x0Start.l, x0End.d, x0End.l); // fprintf(f, "x0Start : %lg 0x%016lx\nx0End : %lg 0x%016lx\n", x0Start.d, x0Start.l, x0End.d, x0End.l); // fprintf(f, "\nxInput\t\txInput (Hex)\t\tresult\t\tresult (Hex)\t\torcle\t\torcle (Hex)\t\tulp error\n"); // printf("testNum : %lu 0x%lx\n", testNumX0, testNumX0); + #endif // loop boundary lenX0 = x0End.d - x0Start.d; @@ -96,7 +99,9 @@ struct errorInfo test1FPEDparamParallel(DL x0Start, DL x0End, unsigned long int // printf("happen to NaN or inf\n"); // exit(1); // } - // fprintf(f, "%le\t%e\n", x0, reUlp); + #ifdef ERRORFILE + fprintf(f, "%le\t%le\n", x0, reUlp); + #endif if (isnormal(reUlp) != 0) { sumError += reUlp; @@ -107,7 +112,6 @@ struct errorInfo test1FPEDparamParallel(DL x0Start, DL x0End, unsigned long int } } } - // fprintf(f, "\n"); // aveReUlp = aveReUlp / (testNumX0); // if(flag == 1) { // printf("all error are 0!!\n"); @@ -116,10 +120,13 @@ struct errorInfo test1FPEDparamParallel(DL x0Start, DL x0End, unsigned long int // printf("%le\t%le\n", aveReUlp, maxReUlp); // printf("\naveReUlp = %lg\nmaxInputX0 = 0x%016lx %lg, maxReUlp = %lg\n", aveReUlp, maxInputX0.l, maxInputX0.d, maxReUlp); // } + #ifdef ERRORFILE + // fprintf(f, "\n"); // clear - // fclose(f); - // free(fileNameSample); + fclose(f); + free(fileNameSample); + #endif mpfr_clears(mpfrOrcle, mpfrResult, (mpfr_ptr) 0); struct errorInfo err; err.sumError = sumError; @@ -178,10 +185,10 @@ int main(int argc, char **argv) { printf("Usage: \tthe fixed inputs [%g %g %lu] will be used\n", x0Start.d, x0End.d, testNumX0); } - if(myid == 0) { - printf("\n---------------------------------------------------start test1paramFPEDParallel\n"); - printf("Parameters: x0Start: %lg, x0End: %lg, testNumX0 = %lu, fileNameKernel: %s\n", x0Start.d, x0End.d, testNumX0, fileNameKernel); - } + // if(myid == 0) { + // printf("\n---------------------------------------------------start test1paramFPEDParallel\n"); + // printf("Parameters: x0Start: %lg, x0End: %lg, testNumX0 = %lu, fileNameKernel: %s\n", x0Start.d, x0End.d, testNumX0, fileNameKernel); + // } // local parameters init int lenX0Local = testNumX0 / numProcs; @@ -195,7 +202,7 @@ int main(int argc, char **argv) { } // call the error test function - struct errorInfo err = test1FPEDparamParallel(x0Start, x0End, testNumX0, fileNameKernel, myid, i0StartLocal, i0EndLocal); + struct errorInfo err = test1FPEDparamParallel(x0Start, x0End, testNumX0, uniqueLabel, fileNameKernel, myid, i0StartLocal, i0EndLocal); // gather errors and find the max struct errorInfo *errs; @@ -234,7 +241,7 @@ int main(int argc, char **argv) { } printf("average ulp\tmax ulp\n"); printf("%lg\t%lg\n", aveError, maxError); - printf("\naveReUlp = %lg\nmaxInputX0 = 0x%016lx %lg, maxReUlp = %lg\n", aveError, maxInputX0.l, maxInputX0.d, maxError); + // printf("\naveReUlp = %lg\nmaxInputX0 = 0x%016lx %lg, maxReUlp = %lg\n", aveError, maxInputX0.l, maxInputX0.d, maxError); fprintf(fErr, "average ulp\tmax ulp\n"); fprintf(fErr, "%lg\t%lg\n", aveError, maxError); fprintf(fErr, "\naveReUlp = %lg\nmaxInputX0 = 0x%016lx %lg, maxReUlp = %lg\n", aveError, maxInputX0.l, maxInputX0.d, maxError); diff --git a/srcTest/testPerformanceOne.sh b/srcTest/testPerformanceOne.sh new file mode 100755 index 0000000000000000000000000000000000000000..fa7ee700a1c2bf51d734c39ddb19ffc62352b129 --- /dev/null +++ b/srcTest/testPerformanceOne.sh @@ -0,0 +1,45 @@ +#example: cd path/to/srcTest; taskset -c 0 ./testPerformanceOne.sh Bsplines3 origin 0 100 + +if [ $# == 1 ]; then + uniqueLabel=${1} + suffix=origin + start=0.01 + end=100 + runtime=10000 +elif [ $# == 4 ]; then + uniqueLabel=${1} + suffix=${2} + start=${3} + end=${4} + runtime=10000 +elif [ $# == 5 ]; then + uniqueLabel=${1} + suffix=${2} + start=${3} + end=${4} + runtime=${5} +else + echo "please input 1, 4 or 5 parameters" + exit +fi + +resultFileName=result_${uniqueLabel}_${start}_${end}.txt +funcName=expr_${uniqueLabel}_${suffix} + +# gcc gccPerformanceTest.c binary.c ${uniqueLabel}/${funcName}.c -DFUNCNAME=${funcName} -DRUNTIME=${runtime} -I../includeTEST -lm -o gccPerformanceTest_${uniqueLabel}.exe -O3 + +# make gccPerformanceTest.o -s CFLAGS="-DFUNCNAME=${funcName} -DRUNTIME=${runtime}" +gcc gccPerformanceTest.c -DFUNCNAME=${funcName} -DRUNTIME=${runtime} -I../includeTEST -c -O3 +make binary.o -s # gcc binary.c -I../includeTEST -c -O3 +gcc ${uniqueLabel}/${funcName}.c -I../includeDD -c -O3 +gcc gccPerformanceTest.o binary.o ${funcName}.o -lm -lqd -o gccPerformanceTest_${uniqueLabel}.exe -O3 +./gccPerformanceTest_${uniqueLabel}.exe ${resultFileName} ${start} ${end} > /dev/null +# rm -rf gccPerformanceTest_${uniqueLabel}.exe + +### compute the cycles of the function ${uniqueLabel} +make dataClean.exe -s CFLAGS="-DRUNTIME=${runtime}" # gcc dataClean.c -DRUNTIME=${runtime} -o dataClean.exe +sort -n ${resultFileName} > ${resultFileName}.sort +./dataClean.exe ${resultFileName}.sort > ../outputs/${uniqueLabel}/${funcName}_performance.txt +rm -rf ${resultFileName} ${resultFileName}.sort +# rm dataClean.exe + diff --git a/srcTest/testPerformanceOneManual.sh b/srcTest/testPerformanceOneManual.sh new file mode 100755 index 0000000000000000000000000000000000000000..ff23c00b8a0de661684cf0cf05418a482ef75f8d --- /dev/null +++ b/srcTest/testPerformanceOneManual.sh @@ -0,0 +1,44 @@ +#example: cd path/to/srcTest; taskset -c 0 ./testPerformanceOneManual.sh Bsplines3 origin 0 100 + +if [ $# == 1 ]; then + uniqueLabel=${1} + suffix=origin + start=0.01 + end=100 + runtime=10000 +elif [ $# == 4 ]; then + uniqueLabel=${1} + suffix=${2} + start=${3} + end=${4} + runtime=10000 +elif [ $# == 5 ]; then + uniqueLabel=${1} + suffix=${2} + start=${3} + end=${4} + runtime=${5} +else + echo "please input 1, 4 or 5 parameters" + exit +fi + +resultFileName=result_${uniqueLabel}_${start}_${end}.txt +funcName=expr_${uniqueLabel}_${suffix} + +# gcc gccPerformanceTest.c binary.c ${uniqueLabel}/${funcName}.c -DFUNCNAME=${funcName} -DRUNTIME=${runtime} -I../includeTEST -lm -o gccPerformanceTest_${uniqueLabel}.exe -O3 + +# make gccPerformanceTest.o -s CFLAGS="-DFUNCNAME=${funcName} -DRUNTIME=${runtime}" +gcc gccPerformanceTest.c -DFUNCNAME=${funcName} -DRUNTIME=${runtime} -I../includeTEST -c -O3 +make binary.o -s # gcc binary.c -I../includeTEST -c -O3 +gcc ${uniqueLabel}/${funcName}.c -I../includeDD -c -O3 +gcc gccPerformanceTest.o binary.o ${funcName}.o -lmpfr -lm -lqd -o gccPerformanceTest_${uniqueLabel}.exe -O3 +./gccPerformanceTest_${uniqueLabel}.exe ${resultFileName} ${start} ${end} +# rm -rf gccPerformanceTest_${uniqueLabel}.exe + +### compute the cycles of the function ${uniqueLabel} +make dataClean.exe -s CFLAGS="-DRUNTIME=${runtime}" # gcc dataClean.c -DRUNTIME=${runtime} -o dataClean.exe +sort -n ${resultFileName} > ${resultFileName}.sort +./dataClean.exe ${resultFileName}.sort +rm -rf ${resultFileName} ${resultFileName}.sort +# rm dataClean.exe