diff --git a/Makefile b/Makefile
index 9553e53ce898295ede7325d88d15c46871f8f997..52aa79cf0dccf2801debede4a1ca18ae1a4ee718 100644
--- a/Makefile
+++ b/Makefile
@@ -40,10 +40,17 @@ DEBUG = -ggdb3
RANLIB = ranlib
#OPTIM = -O2
-INCLUDE = -I.
+INCLUDE = -I include
CXXFLAGS += -std=c++17 -fPIC -D_REENTRANT $(DEBUG) $(OPTIM) $(WARN) $(INCLUDE)
CFLAGS = $(CXXFLAGS)
PROJECTHOME = .
+
+SRC_DIR = src
+OBJ_DIR = obj
+TEST_DIR = test
+TEST_OBJ_DIR = test_obj
+TEST_RUNNER = test_runner
+
##############################################
# support for shared libray versioning
#
@@ -84,9 +91,12 @@ LOCALDESTDIRHEADERS = $(LOCALPREFIX)/include
# header files to be installed for distribution
INSTHEADERS = periodicspline.h spline.h multipolynomial.h interpolatingpolynomial.h interpolator.h
################################################
+SRCS = $(wildcard $(SRC_DIR)/*.cpp)
+OBJS = $(patsubst $(SRC_DIR)/%.cpp, $(OBJ_DIR)/%.o, $(SRCS))
-OBJS = periodicspline.o spline.o multipolynomial.o interpolator.o
+TEST_SRCS = test_src/test_runner.cpp
+TEST_OBJS = test_obj/test_runner.o
########################################################################################
# compiler warnings
@@ -97,48 +107,59 @@ WARN+= -pedantic -Wno-long-long # test with mysql due to long long variables!
#######################################
#options for coverage and performance analysis
-#DEBUG += -fprofile-arcs -ftest-coverage -pg
-#LDFLAGS += -fprofile-arcs -ftest-coverage -pg
+DEBUG += -fprofile-arcs -ftest-coverage -pg
+LDFLAGS += -fprofile-arcs -ftest-coverage -pg
#####################################
+# Compile source files to object files
+$(OBJ_DIR)/%.o: $(SRC_DIR)/%.cpp
+ @mkdir -p $(OBJ_DIR)
+ $(CXX) $(CXXFLAGS) -c $< -o $@
+
-all: lib libstatic test
-
-lib: $(PROJECTHOME)/lib/$(DT_SHLIB)
-
-libstatic: $(PROJECTHOME)/lib/$(LIBRARY)
-
-check: check_spline check_multipoly
+# Compile test source files to object files
+$(TEST_OBJ_DIR)/%.o: $(TEST_DIR)/%.cpp
+ @mkdir -p $(TEST_OBJ_DIR)
+ $(CXX) $(CXXFLAGS) -c $< -o $@
+$(TEST_OBJS): $(TEST_SRCS)
+ $(CXX) $(CXXFLAGS) -c $< -o $@
+
+# Test target: compile and link the test executable
+test: $(TEST_OBJS) $(OBJS)
+ $(info [$(TEST_OBJS)])
+ $(info [$(TEST_SRCS)])
+ $(info [$(OBJS)])
+ $(CXX) $(CXXFLAGS) $(TEST_OBJS) $(OBJS) -o $(TEST_RUNNER)
+ ./$(TEST_RUNNER)
-test: test_spline test_inverse test_freespline test_multipoly test_periodicspline spltest
-test_spline: test_spline.o spline.o interpolator.o
+all: lib libstatic
-test_inverse: test_inverse.o spline.o interpolator.o
+lib: $(PROJECTHOME)/lib/$(DT_SHLIB)
-test_freespline: test_freespline.o spline.o interpolator.o
+libstatic: $(PROJECTHOME)/lib/$(LIBRARY)
-test_periodicspline: test_periodicspline.o periodicspline.o interpolator.o
+check: test
+ gcov test_src/test_runner.cpp
+ gcov src/spline.cpp
+ gcov src/interpolator.cpp
+ gcov src/multipolynomial.cpp
+ gcov src/peridocspline.cpp
+ valgrind --num-callers=20 --leak-check=yes --leak-resolution=high --leak-check=full --show-reachable=yes --log-file=test_runner.grind ./test_runner
-spltest: spltest.o spline.o interpolator.o
-test_multipoly: test_multipoly.o multipolynomial.o interpolator.o
-check_spline: test_spline
+check_exec: test_runner
./test_spline
- gcov test_spline.cpp
- gcov spline.cpp
- gcov interpolator.cpp
+ gcov test_src/test_runner.cpp
+ gcov src/spline.cpp
+ gcov src/interpolator.cpp
+ gcov src/multipolynomial.cpp
+ gcov src/peridocspline.cpp
valgrind --num-callers=20 --leak-check=yes --leak-resolution=high --leak-check=full --show-reachable=yes --log-file=test_spline.grind ./test_spline
-check_multipoly: test_multipoly
- gcov test_multipoly.cpp
- gcov multipoly.cpp
- gcov interpolator.cpp
- valgrind --num-callers=20 --leak-check=yes --leak-resolution=high --leak-check=full --show-reachable=yes --log-file=test_multipoly.grind ./test_multipoly
-
$(PROJECTHOME)/lib/$(LIBRARY):$(OBJS) $(MAKEFILE)
# @echo "Archiving $(LIBRARY) ... $(OBJS)"
mkdir -p -m 755 $(PROJECTHOME)/lib
@@ -160,7 +181,10 @@ doc:
docclean:
rm -rf docs doxy.log
clean:
- rm -f *.o core* *.gcov *.gcno *.gcda *.grind *grind.core.* *.valgrind dump* gmon.out doxygen_log.txt dump* *.csv
+ rm -rf $(OBJ_DIR)
+ rm -rf $(TEST_RUNNER)
+ rm -rf $(TEST_OBJS)
+ rm -f *.o core* *.gcov *.gcno *.gcda *.grind *grind.core.* *.valgrind dump* gmon.out doxygen_log.txt dump* *.csv perf.data
rm -f test_spline test_inverse test_periodicspline test_freespline test_multipoly test_magnets test_dipole test_quadrupole test_sextupole test_kicker test_solenoid spltest
rm -f *.csv
rm -f doxy.log
diff --git a/interpolatingpolynomial.h b/include/interpolatingpolynomial.h
similarity index 100%
rename from interpolatingpolynomial.h
rename to include/interpolatingpolynomial.h
diff --git a/interpolator.h b/include/interpolator.h
similarity index 100%
rename from interpolator.h
rename to include/interpolator.h
diff --git a/multipolynomial.h b/include/multipolynomial.h
similarity index 100%
rename from multipolynomial.h
rename to include/multipolynomial.h
diff --git a/periodicspline.h b/include/periodicspline.h
similarity index 100%
rename from periodicspline.h
rename to include/periodicspline.h
diff --git a/spline.h b/include/spline.h
similarity index 100%
rename from spline.h
rename to include/spline.h
diff --git a/spltest.cpp b/spltest.cpp
deleted file mode 100644
index 3bd43a21d7b03442e7a865455c7cb2303f9c6b03..0000000000000000000000000000000000000000
--- a/spltest.cpp
+++ /dev/null
@@ -1,138 +0,0 @@
-
-
-/* $Author: claudio $
- *
- * $Revision: 1.1.1.1 $
- *
- * $Log: spltest.cpp,v $
- * Revision 1.1.1.1 2012-09-05 08:17:22 claudio
- * Interpolator utility classes
- *
- */
-
-#include <cstdlib>
-#include <iostream>
-#include <iomanip>
-#include "spline.h"
-
-
-using namespace std;
-using namespace Interpolator;
-
-int main(){
-
- try{
- double i_qd[]={0.0,14.8056,74.028,148.056,222.084,296.112,333.126,370.14};
- double g_qd[]={0.0,0.7745,3.9519,7.9231,11.8187,15.7069,17.5649,19.2667};
-// doubleVector iqd;
-// doubleVector fieldqd;
- vector<double> iqd;
- vector<double> fieldqd;
-
- //controle printout of table
- int PRECISION=16;
- int WIDTH=20;
-
- for (unsigned int i=0;i<8;i++){
- iqd.push_back(i_qd[i]);
- fieldqd.push_back(g_qd[i]);
-
- cout << setprecision(PRECISION) << setw(WIDTH) << right;
- cout << iqd[i];
-
- cout << setprecision(PRECISION) << setw(WIDTH) << right;
- cout << fieldqd[i] << endl;
- }
- Spline* dir= new Spline(iqd,fieldqd);
- Spline* inv= new Spline(fieldqd,iqd);
-
-
- double min=0.0,max=i_qd[7];
- unsigned int steps=200;
- double mincur=min;
- double inccur=(max-min)/steps;
- double I,g,inv_I;
- cout << setw(WIDTH) << right;
- cout << "I";
-
- cout << setw(WIDTH) << right;
- cout << "g" ;
-
- cout << setw(WIDTH) << right;
- cout << "inv_I" ;
-
-
-
- cout << endl;
- for(unsigned int i=0;i<=steps;i++){
- I=mincur+ i*inccur;
- g=dir->evaluate(I);
- inv_I=inv->evaluate(g);
-
-
-
- cout << setprecision(PRECISION) << setw(WIDTH) << right;
- cout << I;
-
- cout << setprecision(PRECISION) << setw(WIDTH) << right;
- cout << g;
-
- cout << setprecision(PRECISION) << setw(WIDTH) << right;
- cout << inv_I;
-
-
- if(I != 0.0) {
- cout << setprecision(PRECISION) << setw(WIDTH) << right;
- cout << 10000.0*((I-inv_I)/I);
- }
- else{
- cout << setprecision(PRECISION) << setw(WIDTH) << right;
- cout << 0.0;
- }
-
-// if(dir->bracket(g,I,1.0e-4,xmin,xmax)){
-// cout << setprecision(PRECISION) << setw(WIDTH) << right;
-// cout << xmin;
-//
-// cout << setprecision(PRECISION) << setw(WIDTH) << right;
-// cout << xmax;
-//
-// double I2=dir->solve(g,xmin,xmax,1.0e-8,100);
-//
-// cout << setprecision(PRECISION) << setw(WIDTH) << right;
-// cout << I2;
-//
-// }
-// else{
-// cout << setprecision(PRECISION) << setw(WIDTH) << right;
-// cout << "------";
-//
-// cout << setprecision(PRECISION) << setw(WIDTH) << right;
-// cout << "++++++";
-//
-// cout << setprecision(PRECISION) << setw(WIDTH) << right;
-// cout << "//////";
-// }
- {
-
-
- double I2=dir->inverse_evaluate(g,I);
-
- cout << setprecision(PRECISION) << setw(WIDTH) << right;
- cout << I2;
-
- }
-
-
- cout << endl;
- }
-
- delete inv;
- delete dir;
- }
- catch(std::invalid_argument& ex)
- {
- cout<<ex.what()<<endl;
- }
- return 0;
-}
diff --git a/interpolator.cpp b/src/interpolator.cpp
similarity index 100%
rename from interpolator.cpp
rename to src/interpolator.cpp
diff --git a/multipolynomial.cpp b/src/multipolynomial.cpp
similarity index 100%
rename from multipolynomial.cpp
rename to src/multipolynomial.cpp
diff --git a/periodicspline.cpp b/src/periodicspline.cpp
similarity index 100%
rename from periodicspline.cpp
rename to src/periodicspline.cpp
diff --git a/spline.cpp b/src/spline.cpp
similarity index 100%
rename from spline.cpp
rename to src/spline.cpp
diff --git a/test_freespline.cpp b/test_freespline.cpp
deleted file mode 100644
index 185684042fa571342cf6806b1a2bc5033f8f294b..0000000000000000000000000000000000000000
--- a/test_freespline.cpp
+++ /dev/null
@@ -1,263 +0,0 @@
-/*
- * program to test spline class
- *
- * $Author: claudio $
- *
- * $Revision: 1.1.1.1 $
- *
- * $Log: test_freespline.cpp,v $
- * Revision 1.1.1.1 2012-09-05 08:17:22 claudio
- * Interpolator utility classes
- *
- *
- */
-
-#include <iostream>
-#include <fstream>
-#include "spline.h"
-#include <cmath>
-
-
-using namespace std;
-using namespace Interpolator;
-
-int main(){
- cout<< "test of Spline interpolator class" << endl;
- doubleVector x0,y0;
- const int NP=30;
- for(int i=0; i < NP; i++){
- x0.push_back(i*i*.02);
- y0.push_back(sin(i*i*.02));
- }
- Spline* test_spline = new Spline(x0,y0,NAN,0.0);
- {
- doubleVector a;
- vector<double> b;
- test_spline->get_base_points(a,b);
- cout<<"test get_base_points()"<<endl;
- for(int i=0; i < NP; i++ )
- if((x0[i] != a[i]) ||(y0[i] != b[i])){
- cout << i << " " << x0[i] << " " << a[i] << " | "<< y0[i] << " " << b[i] << endl;
- return -1;
- }
-
- }
- {
- cout << "test get_range()" << endl;
- double minx,maxx;
- test_spline->get_range(minx,maxx);
- if( (minx != x0[0]) || (maxx != x0[NP-1]) ){
- cout << minx << " " << x0[0] << " | " << maxx << " " << x0[NP-1] << endl;
- return -1;
- }
- }
- {
- cout << "test evaluate()" << endl;
- double xn=0.02;
- double yn=test_spline->evaluate(xn);
- cout << xn << " " << yn << endl;
- yn=test_spline->evaluate(x0[5]);
- cout << xn << " " << yn << " " << y0[5] << endl;
- }
- {
- cout << "test evaluate()" << endl;
- double xn=0.02;
- double yn=test_spline->evaluate(xn);
- cout << xn << " " << yn << endl;
- yn=test_spline->evaluate(x0[5]);
- cout << xn << " " << yn << " " << y0[5] << endl;
- }
- { cout << "test evaluate()for plot" << endl;
- doubleVector X;
- doubleVector Y;
- double range=x0[NP-1];
- double nsteps=300;
- for(int i=0;i<=nsteps;i++){
- X.push_back((range*i)/nsteps);
- Y.push_back(test_spline->evaluate(X[i]));
- }
- //print out data
- cout<<"---------------------------"<<endl;
- for(int i=0; i < NP; i++) cout<<x0[i]<<"\t"<<y0[i]<<endl;
- cout<<"---------------------------"<<endl;
- for(int i=0; i < nsteps; i++) cout<<X[i]<<"\t"<<Y[i]<<endl;
- cout<<"+++++++++++++++++++++++++++"<<endl;
- }
- {
- cout << "test evaluate()for doubleVector" << endl;
- doubleVector X;
- doubleVector Y;
- double range=x0[NP-1];
- double nsteps=300;
- for(int i=0;i<=nsteps;i++)
- X.push_back((range*i)/nsteps);
- Y=test_spline->evaluate(X);
- //print out data
- cout<<"###########################"<<endl;
- for(int i=0; i < nsteps; i++) cout<<X[i]<<"\t"<<Y[i]<<endl;
- cout<<"***************************"<<endl;
- }
- {
- cout << "test range exception" << endl;
- try{
- test_spline->evaluate(-10.0);
- cerr << "range_error not thrown!" << endl;
- }
- catch(std::range_error& ex){
- cout << "OK, catched expected exception: "<< ex.what() << endl;
- }
- catch(...){
- cerr << "unexpected exception thrown!" << endl;
- }
- }
- {
- cout << "test constructor with unequal lengths" << endl;
- try{
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(3.0);
-
- y.push_back(4.0);
- y.push_back(5.0);
- Spline spl(x,y,NAN,0.0);
- cerr << "length_error not thrown!" << endl;
- }
- catch(std::length_error& ex){
- cout << "OK, catched expected exception: "<< ex.what() << endl;
- }
- catch(...){
- cerr << "unexpected exception thrown!" << endl;
- }
- }
- {
- cout << "test constructor with out of order data" << endl;
- try{
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(1.5);
-
- y.push_back(4.0);
- y.push_back(5.0);
- y.push_back(6.0);
- Spline spl(x,y,NAN,0.0);
- cerr << "domain_error not thrown!" << endl;
- }
- catch(std::domain_error& ex){
- cout << "OK, catched expected exception: "<< ex.what() << endl;
- }
- catch(...){
- cerr << "unexpected exception thrown!" << endl;
- }
- }
- {
- cout << "Test copy constructor" << endl;
- Spline copy_spline(*test_spline);
- double minx,maxx,copyminx,copymaxx;
- test_spline->get_range(minx,maxx);
- copy_spline.get_range(copyminx,copymaxx);
- if( (minx != copyminx) || (maxx != copymaxx) ){
- cerr << "COPY FAILED!" << endl;
- return -1;
- }
- double range=x0[NP-1];
- double nsteps=100;
- for(int i=0;i<nsteps;i++){
- double val=(range*i)/nsteps;
- double y,ycopy;
- ycopy=copy_spline.evaluate(val);
- y=test_spline->evaluate(val);
- if(y != ycopy){
- cerr << "COPY FAILED " << val << " " << y << " " << ycopy << endl;
- }
- }
- }
- {
- cout << "Test copy assignment operator" << endl;
- Spline oper_spline=*test_spline;
- double minx,maxx,operminx,opermaxx;
- test_spline->get_range(minx,maxx);
- oper_spline.get_range(operminx,opermaxx);
- if( (minx != operminx) || (maxx != opermaxx) ){
- cerr << "ASSIGN COPY FAILED!" << endl;
- return -1;
- }
- double range=x0[NP-1];
- double nsteps=100;
- for(int i=0;i<nsteps;i++){
- double val=(range*i)/nsteps;
- double y,yoper;
- yoper=oper_spline.evaluate(val);
- y=test_spline->evaluate(val);
- if(y != yoper){
- cerr << "COPY FAILED " << val << " " << y << " " << yoper << endl;
- }
- }
- }
- cout << "test delete" << endl;
- delete test_spline;
- {
- cout << "test leaks - new/delete" << endl;
- Spline *dynspl;
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(3.5);
-
- y.push_back(4.0);
- y.push_back(5.0);
- y.push_back(6.0);
- for(int i=0; i < 100 ; i++){
- dynspl = new Spline(x,y,NAN,0.0);
- delete dynspl;
- }
- }
- {
- cout << "test leaks - stack" << endl;
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(3.5);
-
- y.push_back(4.0);
- y.push_back(5.0);
- y.push_back(6.0);
- for(int i=0; i < 100 ; i++){
- Spline stspl(x,y,NAN,0.0);
- }
- }
-
- {
- cout <<"test end conditions" << endl;
- std::ofstream outfile;
- outfile.open("splineval.csv",ios::out); //CHECK success
- if ( (outfile.rdstate() & ifstream::failbit ) != 0 ) return -1;
- doubleVector x,y;
- x.push_back(0.0);y.push_back(0.0);
- x.push_back(1.0);y.push_back(1.0);
- x.push_back(1.5);y.push_back(1.0);
- x.push_back(2.0);y.push_back(1.0);
- x.push_back(3.0);y.push_back(0.0);
- Spline nsp(x,y); //natural spline
- Spline zsp(x,y,0.0,0.0); // spline with 0 derivate start and stop
- Spline nzsp(x,y,NAN,0.0); // spline with 0 derivate at stop
- Spline znsp(x,y,0.0,NAN); // spline with 0 derivate at start
- Spline ddsp(x,y,-0.5,0.5); // spline with derivatives at start and stop
- double X;
- double Yn,Yz,Ynz,Yzn,Ydd;
- double range=3.0;
- int nsteps=300;
- for(int i=0;i<=nsteps;i++){
- X=((range*i)/nsteps);
- Yn=nsp.evaluate(X);
- Yz=zsp.evaluate(X);
- Ynz=nzsp.evaluate(X);
- Yzn=znsp.evaluate(X);
- Ydd=ddsp.evaluate(X);
- outfile<< X << "\t\t" << Yn << "\t\t" << Yz << "\t\t" << Ynz << "\t\t" << Yzn << "\t\t" << Ydd << std::endl;
- }
- outfile.close();
- }
- return 0;
-}
diff --git a/test_inverse.cpp b/test_inverse.cpp
deleted file mode 100644
index 0fa48ce0b5f087194e339875280f590d7587b1b5..0000000000000000000000000000000000000000
--- a/test_inverse.cpp
+++ /dev/null
@@ -1,67 +0,0 @@
-/*
- * program to test spline classes
- *
- * $Author: claudio $
- *
- * $Revision: 1.1 $
- *
- * $Log: test_inverse.cpp,v $
- * Revision 1.1 2017-03-03 15:32:40 claudio
- * correctetd spline::inverse_calculate if returned value is at range limit
- *
- * Revision 1.1.1.1 2012-09-05 08:17:22 claudio
- * Interpolator utility classes
- *
- *
- */
-
-#include <iostream>
-#include "spline.h"
-#include <cmath>
-//#include <stdexcept>
-
-using namespace std;
-using namespace Interpolator;
-int main(){
- cout << "test spline::inverse_evaluate" << endl;
- {
-
- doubleVector x,y;
- x.push_back(10.0); y.push_back(83.6);
- x.push_back(12); y.push_back(68.3);
- x.push_back( 14); y.push_back( 55.7);
- x.push_back( 16); y.push_back(45.5);
- x.push_back( 18); y.push_back(37.3);
- x.push_back( 20); y.push_back(30.7);
- x.push_back( 22); y.push_back( 25.4);
- x.push_back( 24); y.push_back( 21.2);
- x.push_back( 26); y.push_back( 17.9);
- x.push_back( 28); y.push_back( 15.2);
- x.push_back( 30); y.push_back( 13.1);
- x.push_back( 32); y.push_back( 11.4);
- x.push_back( 34); y.push_back( 10.1);
- x.push_back( 36); y.push_back( 9);
- x.push_back( 40); y.push_back( 7.5);
- x.push_back( 50); y.push_back( 5.77);
- Spline *lambdasplie = new Spline(x,y);
- double gmin,gmax;
-
- lambdasplie->get_range(gmin,gmax);
- double max = lambdasplie->evaluate(gmin);
- double min = lambdasplie->evaluate(gmax);
- double initial=min +((max-min)/2);
-
- double l=min;
- double gap, linv;
- cout << "l\tgap\tlinverse" << endl;
- cout.precision(8);
- while (l<max){
- gap = lambdasplie->inverse_evaluate(l,initial);
- linv = lambdasplie->evaluate(gap);
- cout << l <<"\t" << gap << "\t" << linv <<endl;
- l+=0.05;
- }
- delete lambdasplie;
- }
- return 0;
-}
diff --git a/test_multipoly.cpp b/test_multipoly.cpp
deleted file mode 100644
index d67022b86af83f73e0b22dbc1311184d72a77f45..0000000000000000000000000000000000000000
--- a/test_multipoly.cpp
+++ /dev/null
@@ -1,282 +0,0 @@
-/*
- * program to test multipolynomial class
- *
- * $Author: claudio $
- *
- * $Revision: 1.1.1.1 $
- *
- * $Log: test_multipoly.cpp,v $
- * Revision 1.1.1.1 2012-09-05 08:17:22 claudio
- * Interpolator utility classes
- *
- */
-
-#include <iostream>
-#include "multipolynomial.h"
-
-using namespace std;
-using namespace Interpolator;
-
-int main(){
- cout << "test of multiPolynomial interpolator class" << endl;
- try{
- //create interpolating polynomials
-
- InterpolatingPolynomial p0,p1,p2;
- p0.xMin=0.5;
- p0.xMax=1.5;
- p0.coefficient.push_back(1.1);
-
- p1.xMin=1.5;
- p1.xMax=2.5;
- p1.coefficient.push_back(1.1);
- p1.coefficient.push_back(2.2);
-
- p2.xMin=2.5;
- p2.xMax=3.5;
- p2.coefficient.push_back(1.1);
- p2.coefficient.push_back(2.2);
- p2.coefficient.push_back(3.3);
-
- //stuff polynomial into vector
-
- InterpolatingPolynoamialVector pvector;
- pvector.push_back(p0);
- pvector.push_back(p1);
- pvector.push_back(p2);
-
- multiPolynomial *multipoly = new multiPolynomial(pvector);
- double minX,maxX;
- multipoly->get_range(minX,maxX);
- if( (minX != p0.xMin) || ( maxX != p2.xMax) ){
- cout << "multiPoly.get_range(): error " << minX << " " << maxX << " " << p0.xMin << " " << p2.xMax << endl;
- return -1; //stop further testing
- }
- delete multipoly;
- }
- catch(...){
- cerr << " unexpected excpetion!" << endl;
- return -1;
- }
-
- {
- cout << "test a single 2nd degree polynomial " << endl;
- try{
- InterpolatingPolynomial p0;
- p0.xMin=0.0;
- p0.xMax=4.0;
- p0.coefficient.push_back(0.0); //a0
- p0.coefficient.push_back(0.0); //a1
- p0.coefficient.push_back(1.0); //a2 --> Y=X*X
- InterpolatingPolynoamialVector pvector;
- pvector.push_back(p0);
- multiPolynomial *multipoly = new multiPolynomial(pvector);
- double v=0.0,y;
- y=multipoly->evaluate(v); if( y != 0.0 ){ cerr << "value not correct (0.0): " << y << endl; return -1; }
- v=1.0;
- y=multipoly->evaluate(v); if( y != 1.0 ){ cerr << "value not correct (1.0): " << y << endl; return -1; }
- v=4.0;
- y=multipoly->evaluate(v); if( y != 16.0 ){ cerr << "value not correct (16.0): " << y << endl; return -1; }
- doubleVector X;
- for (int i=0; i<=400;i++) X.push_back(i/100.0); //create range for testing
- doubleVector Y=multipoly->evaluate(X);
- cout << "----------------------------------" << endl;
- for (unsigned int i=0; i < X.size(); i++) cout << X[i] << "\t" << Y[i] << endl;
- cout << "----------------------------------" << endl;
- delete multipoly;
- }
- catch(std::range_error& ex){
- std::cerr << ex.what() << std::endl;
- return -1;
- }
- catch(std::length_error& ex){
- std::cerr << ex.what() << std::endl;
- return -1;
- }
- try{
- //create interpolating polynomials
- cout << "test non contiguos" << endl;
- InterpolatingPolynomial p0,p1,p2;
- p0.xMin=0.5;
- p0.xMax=1.5;
- p0.coefficient.push_back(1.1);
-
- p1.xMin=1.6; // <-------
- p1.xMax=2.5;
- p1.coefficient.push_back(1.1);
- p1.coefficient.push_back(2.2);
-
- p2.xMin=2.5;
- p2.xMax=3.5;
- p2.coefficient.push_back(1.1);
- p2.coefficient.push_back(2.2);
- p2.coefficient.push_back(3.3);
-
- //stuff polynomial into vector
-
- InterpolatingPolynoamialVector pvector;
- pvector.push_back(p0);
- pvector.push_back(p1);
- pvector.push_back(p2);
-
- multiPolynomial *multipoly = new multiPolynomial(pvector);
- delete multipoly;
- }
- catch(std::domain_error& ex){
- cout << "OK expected exception: " << ex.what() << endl;
-
- }
- catch(...){
- cerr << " unexpected excpetion!" << endl;
- return -1;
- }
-
- try{
- //create interpolating polynomials
- cout << "test 3 polynomials" << endl;
- InterpolatingPolynomial p0,p1,p2;
- p0.xMin=0.0;
- p0.xMax=1.5;
- p0.coefficient.push_back(1.0); //c[0]
-
- p1.xMin=1.5;
- p1.xMax=2.5;
- p1.coefficient.push_back(0.0); //c[0]
- p1.coefficient.push_back(0.0); //c[1]
- p1.coefficient.push_back(1.0); //c[2]
-
- p2.xMin=2.5;
- p2.xMax=3.5;
- p2.coefficient.push_back(5.0); //c[0]
- p2.coefficient.push_back(-2.0); //c[1]
- //stuff polynomial into vector
-
- InterpolatingPolynoamialVector pvector;
- pvector.push_back(p0);
- pvector.push_back(p1);
- pvector.push_back(p2);
-
- multiPolynomial *multipoly = new multiPolynomial(pvector);
- doubleVector X;
- for (int i=0; i<=350;i++) X.push_back(i/100.0); //create range for testing
- doubleVector Y=multipoly->evaluate(X);
- cout << "-+++++++++++++++++++++++++++++++++" << endl;
- for (unsigned int i=0; i < X.size(); i++) cout << X[i] << "\t" << Y[i] << endl;
- cout << "++++++++++++++++++++++++++++++++++" << endl;
- delete multipoly;
- }
- catch(std::domain_error& ex){
- cout << "unexpected exception: " << ex.what() << endl;
- return -1;
- }
- catch(std::range_error& ex){
- std::cerr << ex.what() << std::endl;
- return -1;
- }
- catch(std::length_error& ex){
- std::cerr << ex.what() << std::endl;
- return -1;
- }
- catch(...){
- cerr << " unexpected excpetion!" << endl;
- return -1;
- }
-
- try{
- //create interpolating polynomials
- cout << "test 3 identical polynomials" << endl;
- InterpolatingPolynomial p0,p1,p2;
- p0.xMin=0.0;
- p0.xMax=1.5;
- p0.coefficient.push_back(1.0); //c[0]
-
-
- InterpolatingPolynoamialVector pvector;
- pvector.push_back(p0);
- pvector.push_back(p0);
- pvector.push_back(p0);
-
- multiPolynomial *multipoly = new multiPolynomial(pvector);
- doubleVector X;
- for (int i=0; i<=150;i++) X.push_back(i/100.0); //create range for testing
- doubleVector Y=multipoly->evaluate(X);
- cout << "**********************************" << endl;
- for (unsigned int i=0; i < X.size(); i++) cout << X[i] << "\t" << Y[i] << endl;
- cout << "**********************************" << endl;
- delete multipoly;
- }
- catch(std::domain_error& ex){
- cout << "OK, exception caught: " << ex.what() << endl;
- }
- catch(std::range_error& ex){
- std::cerr << ex.what() << std::endl;
- return -1;
- }
- catch(std::length_error& ex){
- std::cerr << ex.what() << std::endl;
- return -1;
- }
- catch(...){
- cerr << " unexpected excpetion!" << endl;
- return -1;
- }
- try{
- //create interpolating polynomials
- cout << "test copy constructor" << endl;
- InterpolatingPolynomial p0,p1,p2;
- p0.xMin=0.0;
- p0.xMax=1.5;
- p0.coefficient.push_back(1.0); //c[0]
-
- p1.xMin=1.5;
- p1.xMax=2.5;
- p1.coefficient.push_back(0.0); //c[0]
- p1.coefficient.push_back(0.0); //c[1]
- p1.coefficient.push_back(1.0); //c[2]
-
- p2.xMin=2.5;
- p2.xMax=3.5;
- p2.coefficient.push_back(5.0); //c[0]
- p2.coefficient.push_back(-2.0); //c[1]
- //stuff polynomial into vector
-
- InterpolatingPolynoamialVector pvector;
- pvector.push_back(p0);
- pvector.push_back(p1);
- pvector.push_back(p2);
-
- multiPolynomial *multipoly = new multiPolynomial(pvector);
- multiPolynomial copypoly=*multipoly;
- doubleVector X;
- for (int i=0; i<=350;i++) X.push_back(i/100.0); //create range for testing
- doubleVector Y=multipoly->evaluate(X);
- doubleVector Yc=copypoly.evaluate(X);
-
- for (unsigned int i=0; i < X.size(); i++){
- if (Y[i] != Yc[i])
- cerr << "COPY ERROR: " << X[i] << " " << Y[i] << " " << Yc[i] <<endl;
- }
- delete multipoly;
- }
- catch(...){
- cerr << " test copy constructor: unexpected excpetion!" << endl;
- return -1;
- }
- }
- {
- cout << "test interpolatore for constat value on simmetric range" << endl;
- InterpolatingPolynomial p0;
- InterpolatingPolynoamialVector pvector;
- p0.xMin=-5.0;
- p0.xMax=5.0;
- p0.coefficient.push_back(0.52); //c[0]
- pvector.push_back(p0);
- multiPolynomial *multipoly = new multiPolynomial(pvector);
- double a0=-4.99;
- cout << "for "<< a0 << " "<< multipoly->evaluate(a0);
- a0=4.99;
- cout << "for "<< a0 << " "<< multipoly->evaluate(a0);
- delete multipoly;
- }
- return 0;
-}
diff --git a/test_obj/test_runner.gcda b/test_obj/test_runner.gcda
new file mode 100644
index 0000000000000000000000000000000000000000..4efd016f8520e9fa7d1cf9bf6795fd8490652704
Binary files /dev/null and b/test_obj/test_runner.gcda differ
diff --git a/test_obj/test_runner.gcno b/test_obj/test_runner.gcno
new file mode 100644
index 0000000000000000000000000000000000000000..4b985fd5a481894c197ae972ac0804f07a7a4f90
Binary files /dev/null and b/test_obj/test_runner.gcno differ
diff --git a/test_periodicspline.cpp b/test_periodicspline.cpp
deleted file mode 100644
index 2d8af4e338944262aaeae7d968faba30e40947b1..0000000000000000000000000000000000000000
--- a/test_periodicspline.cpp
+++ /dev/null
@@ -1,242 +0,0 @@
-/*
- * program to test peridic spline class
- *
- * $Author: claudio $
- *
- * $Revision: 1.2 $
- *
- * $Log: test_periodicspline.cpp,v $
- * Revision 1.2 2017-03-06 13:13:20 claudio
- * use unsigned for comparison
- *
- * Revision 1.1.1.1 2012-09-05 08:17:22 claudio
- * Interpolator utility classes
- *
- *
- */
-
-#include <iostream>
-#include "periodicspline.h"
-#include <cmath>
-
-using namespace std;
-using namespace Interpolator;
-
-int main(){
- cout<< "test of periodicSpline interpolator class" << endl;
- doubleVector x0,y0;
- const int NP=15;
- const double k = (2*M_PI)/NP;
- for(int i=0; i <= NP; i++){
- x0.push_back(i*k);
- y0.push_back(sin(i*k));
- }
- periodicSpline* test_spline = new periodicSpline(x0,y0);
- {
- doubleVector a;
- doubleVector b;
- test_spline->get_base_points(a,b);
- cout<<"test get_base_points()"<<endl;
- for(unsigned int i=0; i < a.size()-1; i++ )
- if((x0[i] != a[i]) ||(y0[i] != b[i])){
- cout << i << " " << x0[i] << " " << a[i] << " | "<< y0[i] << " " << b[i] << endl;
- return -1;
- }
-
- }
- {
- cout << "test get_range()" << endl;
- double minx,maxx;
- test_spline->get_range(minx,maxx);
- int I=x0.size()-1;
- if( (minx != x0[0]) || (maxx != x0[I]) ){
- cout << minx << " " << x0[0] << " | " << maxx << " " << x0[I] << endl;
- return -1;
- }
- }
- {
- cout << "test evaluate()" << endl;
- double xn=0.02;
- double yn=test_spline->evaluate(xn);
- cout << xn << " " << yn << endl;
- yn=test_spline->evaluate(x0[5]);
- cout << xn << " " << yn << " " << y0[5] << endl;
- }
- {
- cout << "test evaluate()" << endl;
- double xn=0.02;
- double yn=test_spline->evaluate(xn);
- cout << xn << " " << yn << endl;
- yn=test_spline->evaluate(x0[5]);
- cout << xn << " " << yn << " " << y0[5] << endl;
- }
- { cout << "test evaluate()for plot" << endl;
- doubleVector X;
- doubleVector Y;
- double minx;
- double range;
- test_spline->get_range(minx,range);
- double nsteps=100;
- for(int i=0;i<=nsteps;i++){
- X.push_back((range*i)/nsteps);
- Y.push_back(test_spline->evaluate(X[i]));
- }
- //print out data
- cout<<"---------------------------"<<endl;
- for(unsigned int i=0; i < x0.size(); i++) cout<<x0[i]<<"\t"<<y0[i]<<endl;
- cout<<"---------------------------"<<endl;
- for(int i=0; i <= nsteps; i++) cout<<X[i]<<"\t"<<Y[i]<<endl;
- cout<<"+++++++++++++++++++++++++++"<<endl;
- }
- {
- cout << "test evaluate()for doubleVector" << endl;
- doubleVector X;
- doubleVector Y;
- double minx;
- double range;
- test_spline->get_range(minx,range);
- double nsteps=100;
- for(int i=0;i<=nsteps;i++)
- X.push_back((range*i)/nsteps);
- Y=test_spline->evaluate(X);
- //print out data
- cout<<"###########################"<<endl;
- for(int i=0; i <= nsteps; i++) cout<<X[i]<<"\t"<<Y[i]<<endl;
- cout<<"***************************"<<endl;
- }
- {
- cout << "test range exception" << endl;
- try{
- test_spline->evaluate(-10.0);
- cerr << "range_error not thrown!" << endl;
- }
- catch(std::range_error& ex){
- cout << "OK, catched expected exception: "<< ex.what() << endl;
- }
- catch(...){
- cerr << "unexpected exception thrown!" << endl;
- }
- }
- {
- cout << "test constructor with unequal lengths" << endl;
- try{
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(3.0);
-
- y.push_back(4.0);
- y.push_back(5.0);
- periodicSpline spl(x,y);
- cerr << "length_error not thrown!" << endl;
- }
- catch(std::length_error& ex){
- cout << "OK, catched expected exception: "<< ex.what() << endl;
- }
- catch(...){
- cerr << "unexpected exception thrown!" << endl;
- }
- }
- {
- cout << "test constructor with out of order data" << endl;
- try{
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(1.5);
-
- y.push_back(4.0);
- y.push_back(5.0);
- y.push_back(6.0);
- periodicSpline spl(x,y);
- cerr << "domain_error not thrown!" << endl;
- }
- catch(std::domain_error& ex){
- cout << "OK, catched expected exception: "<< ex.what() << endl;
- }
- catch(...){
- cerr << "unexpected exception thrown!" << endl;
- }
- }
- {
- cout << "Test copy constructor" << endl;
- periodicSpline copy_spline(*test_spline);
- double minx,maxx,copyminx,copymaxx;
- test_spline->get_range(minx,maxx);
- copy_spline.get_range(copyminx,copymaxx);
- if( (minx != copyminx) || (maxx != copymaxx) ){
- cerr << "COPY FAILED!" << endl;
- return -1;
- }
- double range;
- test_spline->get_range(minx,range);
- double nsteps=100;
- for(int i=0;i<nsteps;i++){
- double val=(range*i)/nsteps;
- double y,ycopy;
- ycopy=copy_spline.evaluate(val);
- y=test_spline->evaluate(val);
- if(y != ycopy){
- cerr << "COPY FAILED " << val << " " << y << " " << ycopy << endl;
- }
- }
- }
- {
- cout << "Test copy operator" << endl;
- periodicSpline copy_spline(*test_spline);
- copy_spline=*test_spline;
- double minx,maxx,copyminx,copymaxx;
- test_spline->get_range(minx,maxx);
- copy_spline.get_range(copyminx,copymaxx);
- if( (minx != copyminx) || (maxx != copymaxx) ){
- cerr << "COPY FAILED!" << endl;
- return -1;
- }
- double range;
- test_spline->get_range(minx,range);
- double nsteps=100;
- for(int i=0;i<nsteps;i++){
- double val=(range*i)/nsteps;
- double y,ycopy;
- ycopy=copy_spline.evaluate(val);
- y=test_spline->evaluate(val);
- if(y != ycopy){
- cerr << "COPY FAILED " << val << " " << y << " " << ycopy << endl;
- }
- }
- }
- cout << "test delete" << endl;
- delete test_spline;
- {
- cout << "test leaks - new/delete" << endl;
- periodicSpline *dynspl;
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(3.5);
-
- y.push_back(4.0);
- y.push_back(5.0);
- y.push_back(6.0);
- for(int i=0; i < 100 ; i++){
- dynspl = new periodicSpline(x,y);
- delete dynspl;
- }
- }
- {
- cout << "test leaks - stack" << endl;
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(3.5);
-
- y.push_back(4.0);
- y.push_back(5.0);
- y.push_back(6.0);
- for(int i=0; i < 100 ; i++){
- periodicSpline stspl(x,y);
- }
- }
-
- return 0;
-}
diff --git a/test_spline.cpp b/test_spline.cpp
deleted file mode 100644
index b2f0a97fccf92507e2f263af3fa33fc4bb569fad..0000000000000000000000000000000000000000
--- a/test_spline.cpp
+++ /dev/null
@@ -1,209 +0,0 @@
-/*
- * program to test spline classes
- *
- * $Author: claudio $
- *
- * $Revision: 1.1.1.1 $
- *
- * $Log: test_spline.cpp,v $
- * Revision 1.1.1.1 2012-09-05 08:17:22 claudio
- * Interpolator utility classes
- *
- *
- */
-
-#include <iostream>
-#include "spline.h"
-#include <cmath>
-//#include <stdexcept>
-
-using namespace std;
-using namespace Interpolator;
-int main(){
- cout<< "test of Spline interpolator class" << endl;
- doubleVector x0,y0;
- const int NP=30;
- for(int i=0; i < NP; i++){
- x0.push_back(i*i*.02);
- y0.push_back(sin(i*i*.02));
- }
- Spline* test_spline = new Spline(x0,y0);
- {
- doubleVector a;
- doubleVector b;
- test_spline->get_base_points(a,b);
- cout<<"test get_base_points()"<<endl;
- for(int i=0; i < NP; i++ )
- if((x0[i] != a[i]) ||(y0[i] != b[i])){
- cout << i << " " << x0[i] << " " << a[i] << " | "<< y0[i] << " " << b[i] << endl;
- return -1;
- }
-
- }
- {
- cout << "test get_range()" << endl;
- double minx,maxx;
- test_spline->get_range(minx,maxx);
- if( (minx != x0[0]) || (maxx != x0[NP-1]) ){
- cout << minx << " " << x0[0] << " | " << maxx << " " << x0[NP-1] << endl;
- return -1;
- }
- }
- {
- cout << "test evaluate()" << endl;
- double xn=0.02;
- double yn=test_spline->evaluate(xn);
- cout << xn << " " << yn << endl;
- yn=test_spline->evaluate(x0[5]);
- cout << xn << " " << yn << " " << y0[5] << endl;
- }
- {
- cout << "test evaluate()" << endl;
- double xn=0.02;
- double yn=test_spline->evaluate(xn);
- cout << xn << " " << yn << endl;
- yn=test_spline->evaluate(x0[5]);
- cout << xn << " " << yn << " " << y0[5] << endl;
- }
- { cout << "test evaluate()for plot" << endl;
- doubleVector X;
- doubleVector Y;
- double range=x0[NP-1];
- double nsteps=100;
- for(int i=0;i<nsteps;i++){
- X.push_back((range*i)/nsteps);
- Y.push_back(test_spline->evaluate(X[i]));
- }
- //print out data
- cout<<"---------------------------"<<endl;
- for(int i=0; i < NP; i++) cout<<x0[i]<<"\t"<<y0[i]<<endl;
- cout<<"---------------------------"<<endl;
- for(int i=0; i < nsteps; i++) cout<<X[i]<<"\t"<<Y[i]<<endl;
- cout<<"+++++++++++++++++++++++++++"<<endl;
- }
- {
- cout << "test evaluate()for doubleVector" << endl;
- doubleVector X;
- doubleVector Y;
- double range=x0[NP-1];
- double nsteps=100;
- for(int i=0;i<nsteps;i++)
- X.push_back((range*i)/nsteps);
- Y=test_spline->evaluate(X);
- //print out data
- cout<<"###########################"<<endl;
- for(int i=0; i < nsteps; i++) cout<<X[i]<<"\t"<<Y[i]<<endl;
- cout<<"***************************"<<endl;
- }
- {
- cout << "test range exception" << endl;
- try{
- test_spline->evaluate(-10.0);
- cerr << "range_error not thrown!" << endl;
- }
- catch(std::range_error& ex){
- cout << "OK, catched expected exception: "<< ex.what() << endl;
- }
- catch(...){
- cerr << "unexpected exception thrown!" << endl;
- }
- }
- {
- cout << "test constructor with unequal lengths" << endl;
- try{
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(3.0);
-
- y.push_back(4.0);
- y.push_back(5.0);
- Spline spl(x,y);
- cerr << "length_error not thrown!" << endl;
- }
- catch(std::length_error& ex){
- cout << "OK, catched expected exception: "<< ex.what() << endl;
- }
- catch(...){
- cerr << "unexpected exception thrown!" << endl;
- }
- }
- {
- cout << "test constructor with out of order data" << endl;
- try{
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(1.5);
-
- y.push_back(4.0);
- y.push_back(5.0);
- y.push_back(6.0);
- Spline spl(x,y);
- cerr << "domain_error not thrown!" << endl;
- }
- catch(std::domain_error& ex){
- cout << "OK, catched expected exception: "<< ex.what() << endl;
- }
- catch(...){
- cerr << "unexpected exception thrown!" << endl;
- }
- }
- {
- cout << "Test copy constructor" << endl;
- //Spline copy_spline=*test_spline;
- Spline copy_spline(*test_spline);
- double minx,maxx,copyminx,copymaxx;
- test_spline->get_range(minx,maxx);
- copy_spline.get_range(copyminx,copymaxx);
- if( (minx != copyminx) || (maxx != copymaxx) ){
- cerr << "COPY FAILED!" << endl;
- return -1;
- }
- double range=x0[NP-1];
- double nsteps=100;
- for(int i=0;i<nsteps;i++){
- double val=(range*i)/nsteps;
- double y,ycopy;
- ycopy=copy_spline.evaluate(val);
- y=test_spline->evaluate(val);
- if(y != ycopy){
- cerr << "COPY FAILED " << val << " " << y << " " << ycopy << endl;
- }
- }
- }
- cout << "test delete" << endl;
- delete test_spline;
- {
- cout << "test leaks - new/delete" << endl;
- Spline *dynspl;
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(3.5);
-
- y.push_back(4.0);
- y.push_back(5.0);
- y.push_back(6.0);
- for(int i=0; i < 100 ; i++){
- dynspl = new Spline(x,y);
- delete dynspl;
- }
- }
- {
- cout << "test leaks - stack" << endl;
- doubleVector x,y;
- x.push_back(1.0);
- x.push_back(2.0);
- x.push_back(3.5);
-
- y.push_back(4.0);
- y.push_back(5.0);
- y.push_back(6.0);
- for(int i=0; i < 100 ; i++){
- Spline stspl(x,y);
- }
- }
-
- return 0;
-}
diff --git a/test_src/test_runner.cpp b/test_src/test_runner.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..346ac106e4ee2b4d9378747849224f67bbe6d5cc
--- /dev/null
+++ b/test_src/test_runner.cpp
@@ -0,0 +1,1064 @@
+/*
+ * program to test interpolator library
+ *
+ */
+
+#include <iostream>
+#include <fstream>
+#include <iomanip>
+#include "spline.h"
+#include "multipolynomial.h"
+#include "periodicspline.h"
+#include <cmath>
+
+using namespace std;
+using namespace Interpolator;
+
+int main() {
+ {
+ cout << "test of multiPolynomial interpolator class" << endl;
+ try {
+ //create interpolating polynomials
+
+ InterpolatingPolynomial p0, p1, p2;
+ p0.xMin = 0.5;
+ p0.xMax = 1.5;
+ p0.coefficient.push_back(1.1);
+
+ p1.xMin = 1.5;
+ p1.xMax = 2.5;
+ p1.coefficient.push_back(1.1);
+ p1.coefficient.push_back(2.2);
+
+ p2.xMin = 2.5;
+ p2.xMax = 3.5;
+ p2.coefficient.push_back(1.1);
+ p2.coefficient.push_back(2.2);
+ p2.coefficient.push_back(3.3);
+
+ //stuff polynomial into vector
+
+ InterpolatingPolynoamialVector pvector;
+ pvector.push_back(p0);
+ pvector.push_back(p1);
+ pvector.push_back(p2);
+
+ multiPolynomial *multipoly = new multiPolynomial(pvector);
+ double minX, maxX;
+ multipoly->get_range(minX, maxX);
+ if ((minX != p0.xMin) || (maxX != p2.xMax)) {
+ cout << "multiPoly.get_range(): error " << minX << " " << maxX
+ << " " << p0.xMin << " " << p2.xMax << endl;
+ return -1; //stop further testing
+ }
+ delete multipoly;
+ } catch (...) {
+ cerr << " unexpected excpetion!" << endl;
+ return -1;
+ }
+
+ {
+ cout << "test a single 2nd degree polynomial " << endl;
+ try {
+ InterpolatingPolynomial p0;
+ p0.xMin = 0.0;
+ p0.xMax = 4.0;
+ p0.coefficient.push_back(0.0); //a0
+ p0.coefficient.push_back(0.0); //a1
+ p0.coefficient.push_back(1.0); //a2 --> Y=X*X
+ InterpolatingPolynoamialVector pvector;
+ pvector.push_back(p0);
+ multiPolynomial *multipoly = new multiPolynomial(pvector);
+ double v = 0.0, y;
+ y = multipoly->evaluate(v);
+ if (y != 0.0) {
+ cerr << "value not correct (0.0): " << y << endl;
+ return -1;
+ }
+ v = 1.0;
+ y = multipoly->evaluate(v);
+ if (y != 1.0) {
+ cerr << "value not correct (1.0): " << y << endl;
+ return -1;
+ }
+ v = 4.0;
+ y = multipoly->evaluate(v);
+ if (y != 16.0) {
+ cerr << "value not correct (16.0): " << y << endl;
+ return -1;
+ }
+ doubleVector X;
+ for (int i = 0; i <= 400; i++)
+ X.push_back(i / 100.0); //create range for testing
+ doubleVector Y = multipoly->evaluate(X);
+ cout << "----------------------------------" << endl;
+ for (unsigned int i = 0; i < X.size(); i++)
+ cout << X[i] << "\t" << Y[i] << endl;
+ cout << "----------------------------------" << endl;
+ delete multipoly;
+ } catch (std::range_error &ex) {
+ std::cerr << ex.what() << std::endl;
+ return -1;
+ } catch (std::length_error &ex) {
+ std::cerr << ex.what() << std::endl;
+ return -1;
+ }
+ try {
+ //create interpolating polynomials
+ cout << "test non contiguos" << endl;
+ InterpolatingPolynomial p0, p1, p2;
+ p0.xMin = 0.5;
+ p0.xMax = 1.5;
+ p0.coefficient.push_back(1.1);
+
+ p1.xMin = 1.6; // <-------
+ p1.xMax = 2.5;
+ p1.coefficient.push_back(1.1);
+ p1.coefficient.push_back(2.2);
+
+ p2.xMin = 2.5;
+ p2.xMax = 3.5;
+ p2.coefficient.push_back(1.1);
+ p2.coefficient.push_back(2.2);
+ p2.coefficient.push_back(3.3);
+
+ //stuff polynomial into vector
+
+ InterpolatingPolynoamialVector pvector;
+ pvector.push_back(p0);
+ pvector.push_back(p1);
+ pvector.push_back(p2);
+
+ multiPolynomial *multipoly = new multiPolynomial(pvector);
+ delete multipoly;
+ } catch (std::domain_error &ex) {
+ cout << "OK expected exception: " << ex.what() << endl;
+
+ } catch (...) {
+ cerr << " unexpected excpetion!" << endl;
+ return -1;
+ }
+
+ try {
+ //create interpolating polynomials
+ cout << "test 3 polynomials" << endl;
+ InterpolatingPolynomial p0, p1, p2;
+ p0.xMin = 0.0;
+ p0.xMax = 1.5;
+ p0.coefficient.push_back(1.0); //c[0]
+
+ p1.xMin = 1.5;
+ p1.xMax = 2.5;
+ p1.coefficient.push_back(0.0); //c[0]
+ p1.coefficient.push_back(0.0); //c[1]
+ p1.coefficient.push_back(1.0); //c[2]
+
+ p2.xMin = 2.5;
+ p2.xMax = 3.5;
+ p2.coefficient.push_back(5.0); //c[0]
+ p2.coefficient.push_back(-2.0); //c[1]
+ //stuff polynomial into vector
+
+ InterpolatingPolynoamialVector pvector;
+ pvector.push_back(p0);
+ pvector.push_back(p1);
+ pvector.push_back(p2);
+
+ multiPolynomial *multipoly = new multiPolynomial(pvector);
+ doubleVector X;
+ for (int i = 0; i <= 350; i++)
+ X.push_back(i / 100.0); //create range for testing
+ doubleVector Y = multipoly->evaluate(X);
+ cout << "-+++++++++++++++++++++++++++++++++" << endl;
+ for (unsigned int i = 0; i < X.size(); i++)
+ cout << X[i] << "\t" << Y[i] << endl;
+ cout << "++++++++++++++++++++++++++++++++++" << endl;
+ delete multipoly;
+ } catch (std::domain_error &ex) {
+ cout << "unexpected exception: " << ex.what() << endl;
+ return -1;
+ } catch (std::range_error &ex) {
+ std::cerr << ex.what() << std::endl;
+ return -1;
+ } catch (std::length_error &ex) {
+ std::cerr << ex.what() << std::endl;
+ return -1;
+ } catch (...) {
+ cerr << " unexpected excpetion!" << endl;
+ return -1;
+ }
+
+ try {
+ //create interpolating polynomials
+ cout << "test 3 identical polynomials" << endl;
+ InterpolatingPolynomial p0, p1, p2;
+ p0.xMin = 0.0;
+ p0.xMax = 1.5;
+ p0.coefficient.push_back(1.0); //c[0]
+
+ InterpolatingPolynoamialVector pvector;
+ pvector.push_back(p0);
+ pvector.push_back(p0);
+ pvector.push_back(p0);
+
+ multiPolynomial *multipoly = new multiPolynomial(pvector);
+ doubleVector X;
+ for (int i = 0; i <= 150; i++)
+ X.push_back(i / 100.0); //create range for testing
+ doubleVector Y = multipoly->evaluate(X);
+ cout << "**********************************" << endl;
+ for (unsigned int i = 0; i < X.size(); i++)
+ cout << X[i] << "\t" << Y[i] << endl;
+ cout << "**********************************" << endl;
+ delete multipoly;
+ } catch (std::domain_error &ex) {
+ cout << "OK, exception caught: " << ex.what() << endl;
+ } catch (std::range_error &ex) {
+ std::cerr << ex.what() << std::endl;
+ return -1;
+ } catch (std::length_error &ex) {
+ std::cerr << ex.what() << std::endl;
+ return -1;
+ } catch (...) {
+ cerr << " unexpected excpetion!" << endl;
+ return -1;
+ }
+ try {
+ //create interpolating polynomials
+ cout << "test copy constructor" << endl;
+ InterpolatingPolynomial p0, p1, p2;
+ p0.xMin = 0.0;
+ p0.xMax = 1.5;
+ p0.coefficient.push_back(1.0); //c[0]
+
+ p1.xMin = 1.5;
+ p1.xMax = 2.5;
+ p1.coefficient.push_back(0.0); //c[0]
+ p1.coefficient.push_back(0.0); //c[1]
+ p1.coefficient.push_back(1.0); //c[2]
+
+ p2.xMin = 2.5;
+ p2.xMax = 3.5;
+ p2.coefficient.push_back(5.0); //c[0]
+ p2.coefficient.push_back(-2.0); //c[1]
+ //stuff polynomial into vector
+
+ InterpolatingPolynoamialVector pvector;
+ pvector.push_back(p0);
+ pvector.push_back(p1);
+ pvector.push_back(p2);
+
+ multiPolynomial *multipoly = new multiPolynomial(pvector);
+ multiPolynomial copypoly = *multipoly;
+ doubleVector X;
+ for (int i = 0; i <= 350; i++)
+ X.push_back(i / 100.0); //create range for testing
+ doubleVector Y = multipoly->evaluate(X);
+ doubleVector Yc = copypoly.evaluate(X);
+
+ for (unsigned int i = 0; i < X.size(); i++) {
+ if (Y[i] != Yc[i])
+ cerr << "COPY ERROR: " << X[i] << " " << Y[i] << " "
+ << Yc[i] << endl;
+ }
+ delete multipoly;
+ } catch (...) {
+ cerr << " test copy constructor: unexpected excpetion!" << endl;
+ return -1;
+ }
+ }
+ {
+ cout << "test interpolatore for constat value on simmetric range"
+ << endl;
+ InterpolatingPolynomial p0;
+ InterpolatingPolynoamialVector pvector;
+ p0.xMin = -5.0;
+ p0.xMax = 5.0;
+ p0.coefficient.push_back(0.52); //c[0]
+ pvector.push_back(p0);
+ multiPolynomial *multipoly = new multiPolynomial(pvector);
+ double a0 = -4.99;
+ cout << "for " << a0 << " " << multipoly->evaluate(a0);
+ a0 = 4.99;
+ cout << "for " << a0 << " " << multipoly->evaluate(a0);
+ delete multipoly;
+ }
+ }
+ //################################################################
+ {
+ cout << "test of Spline interpolator class" << endl;
+ doubleVector x0, y0;
+ const int NP = 30;
+ for (int i = 0; i < NP; i++) {
+ x0.push_back(i * i * .02);
+ y0.push_back(sin(i * i * .02));
+ }
+ Spline *test_spline = new Spline(x0, y0, NAN, 0.0);
+ {
+ doubleVector a;
+ vector<double> b;
+ test_spline->get_base_points(a, b);
+ cout << "test get_base_points()" << endl;
+ for (int i = 0; i < NP; i++)
+ if ((x0[i] != a[i]) || (y0[i] != b[i])) {
+ cout << i << " " << x0[i] << " " << a[i] << " | " << y0[i]
+ << " " << b[i] << endl;
+ return -1;
+ }
+
+ }
+ {
+ cout << "test get_range()" << endl;
+ double minx, maxx;
+ test_spline->get_range(minx, maxx);
+ if ((minx != x0[0]) || (maxx != x0[NP - 1])) {
+ cout << minx << " " << x0[0] << " | " << maxx << " "
+ << x0[NP - 1] << endl;
+ return -1;
+ }
+ }
+ {
+ cout << "test evaluate()" << endl;
+ double xn = 0.02;
+ double yn = test_spline->evaluate(xn);
+ cout << xn << " " << yn << endl;
+ yn = test_spline->evaluate(x0[5]);
+ cout << xn << " " << yn << " " << y0[5] << endl;
+ }
+ {
+ cout << "test evaluate()" << endl;
+ double xn = 0.02;
+ double yn = test_spline->evaluate(xn);
+ cout << xn << " " << yn << endl;
+ yn = test_spline->evaluate(x0[5]);
+ cout << xn << " " << yn << " " << y0[5] << endl;
+ }
+ {
+ cout << "test evaluate()for plot" << endl;
+ doubleVector X;
+ doubleVector Y;
+ double range = x0[NP - 1];
+ double nsteps = 300;
+ for (int i = 0; i <= nsteps; i++) {
+ X.push_back((range * i) / nsteps);
+ Y.push_back(test_spline->evaluate(X[i]));
+ }
+ //print out data
+ cout << "---------------------------" << endl;
+ for (int i = 0; i < NP; i++)
+ cout << x0[i] << "\t" << y0[i] << endl;
+ cout << "---------------------------" << endl;
+ for (int i = 0; i < nsteps; i++)
+ cout << X[i] << "\t" << Y[i] << endl;
+ cout << "+++++++++++++++++++++++++++" << endl;
+ }
+ {
+ cout << "test evaluate()for doubleVector" << endl;
+ doubleVector X;
+ doubleVector Y;
+ double range = x0[NP - 1];
+ double nsteps = 300;
+ for (int i = 0; i <= nsteps; i++)
+ X.push_back((range * i) / nsteps);
+ Y = test_spline->evaluate(X);
+ //print out data
+ cout << "###########################" << endl;
+ for (int i = 0; i < nsteps; i++)
+ cout << X[i] << "\t" << Y[i] << endl;
+ cout << "***************************" << endl;
+ }
+ {
+ cout << "test range exception" << endl;
+ try {
+ test_spline->evaluate(-10.0);
+ cerr << "range_error not thrown!" << endl;
+ } catch (std::range_error &ex) {
+ cout << "OK, catched expected exception: " << ex.what() << endl;
+ } catch (...) {
+ cerr << "unexpected exception thrown!" << endl;
+ }
+ }
+ {
+ cout << "test constructor with unequal lengths" << endl;
+ try {
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(3.0);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ Spline spl(x, y, NAN, 0.0);
+ cerr << "length_error not thrown!" << endl;
+ } catch (std::length_error &ex) {
+ cout << "OK, catched expected exception: " << ex.what() << endl;
+ } catch (...) {
+ cerr << "unexpected exception thrown!" << endl;
+ }
+ }
+ {
+ cout << "test constructor with out of order data" << endl;
+ try {
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(1.5);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ y.push_back(6.0);
+ Spline spl(x, y, NAN, 0.0);
+ cerr << "domain_error not thrown!" << endl;
+ } catch (std::domain_error &ex) {
+ cout << "OK, catched expected exception: " << ex.what() << endl;
+ } catch (...) {
+ cerr << "unexpected exception thrown!" << endl;
+ }
+ }
+ {
+ cout << "Test copy constructor" << endl;
+ Spline copy_spline(*test_spline);
+ double minx, maxx, copyminx, copymaxx;
+ test_spline->get_range(minx, maxx);
+ copy_spline.get_range(copyminx, copymaxx);
+ if ((minx != copyminx) || (maxx != copymaxx)) {
+ cerr << "COPY FAILED!" << endl;
+ return -1;
+ }
+ double range = x0[NP - 1];
+ double nsteps = 100;
+ for (int i = 0; i < nsteps; i++) {
+ double val = (range * i) / nsteps;
+ double y, ycopy;
+ ycopy = copy_spline.evaluate(val);
+ y = test_spline->evaluate(val);
+ if (y != ycopy) {
+ cerr << "COPY FAILED " << val << " " << y << " " << ycopy
+ << endl;
+ }
+ }
+ }
+ {
+ cout << "Test copy assignment operator" << endl;
+ Spline oper_spline = *test_spline;
+ double minx, maxx, operminx, opermaxx;
+ test_spline->get_range(minx, maxx);
+ oper_spline.get_range(operminx, opermaxx);
+ if ((minx != operminx) || (maxx != opermaxx)) {
+ cerr << "ASSIGN COPY FAILED!" << endl;
+ return -1;
+ }
+ double range = x0[NP - 1];
+ double nsteps = 100;
+ for (int i = 0; i < nsteps; i++) {
+ double val = (range * i) / nsteps;
+ double y, yoper;
+ yoper = oper_spline.evaluate(val);
+ y = test_spline->evaluate(val);
+ if (y != yoper) {
+ cerr << "COPY FAILED " << val << " " << y << " " << yoper
+ << endl;
+ }
+ }
+ }
+ cout << "test delete" << endl;
+ delete test_spline;
+ {
+ cout << "test leaks - new/delete" << endl;
+ Spline *dynspl;
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(3.5);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ y.push_back(6.0);
+ for (int i = 0; i < 100; i++) {
+ dynspl = new Spline(x, y, NAN, 0.0);
+ delete dynspl;
+ }
+ }
+ {
+ cout << "test leaks - stack" << endl;
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(3.5);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ y.push_back(6.0);
+ for (int i = 0; i < 100; i++) {
+ Spline stspl(x, y, NAN, 0.0);
+ }
+ }
+
+ {
+ cout << "test end conditions" << endl;
+ std::ofstream outfile;
+ outfile.open("splineval.csv", ios::out); //CHECK success
+ if ((outfile.rdstate() & ifstream::failbit) != 0)
+ return -1;
+ doubleVector x, y;
+ x.push_back(0.0);
+ y.push_back(0.0);
+ x.push_back(1.0);
+ y.push_back(1.0);
+ x.push_back(1.5);
+ y.push_back(1.0);
+ x.push_back(2.0);
+ y.push_back(1.0);
+ x.push_back(3.0);
+ y.push_back(0.0);
+ Spline nsp(x, y); //natural spline
+ Spline zsp(x, y, 0.0, 0.0); // spline with 0 derivate start and stop
+ Spline nzsp(x, y, NAN, 0.0); // spline with 0 derivate at stop
+ Spline znsp(x, y, 0.0, NAN); // spline with 0 derivate at start
+ Spline ddsp(x, y, -0.5, 0.5); // spline with derivatives at start and stop
+ double X;
+ double Yn, Yz, Ynz, Yzn, Ydd;
+ double range = 3.0;
+ int nsteps = 300;
+ for (int i = 0; i <= nsteps; i++) {
+ X = ((range * i) / nsteps);
+ Yn = nsp.evaluate(X);
+ Yz = zsp.evaluate(X);
+ Ynz = nzsp.evaluate(X);
+ Yzn = znsp.evaluate(X);
+ Ydd = ddsp.evaluate(X);
+ outfile << X << "\t\t" << Yn << "\t\t" << Yz << "\t\t" << Ynz
+ << "\t\t" << Yzn << "\t\t" << Ydd << std::endl;
+ }
+ outfile.close();
+ }
+ }
+ //######################################################################
+ {
+ cout << "test of Spline interpolator class" << endl;
+ doubleVector x0, y0;
+ const int NP = 30;
+ for (int i = 0; i < NP; i++) {
+ x0.push_back(i * i * .02);
+ y0.push_back(sin(i * i * .02));
+ }
+ Spline *test_spline = new Spline(x0, y0, NAN, 0.0);
+ {
+ doubleVector a;
+ vector<double> b;
+ test_spline->get_base_points(a, b);
+ cout << "test get_base_points()" << endl;
+ for (int i = 0; i < NP; i++)
+ if ((x0[i] != a[i]) || (y0[i] != b[i])) {
+ cout << i << " " << x0[i] << " " << a[i] << " | " << y0[i]
+ << " " << b[i] << endl;
+ return -1;
+ }
+
+ }
+ {
+ cout << "test get_range()" << endl;
+ double minx, maxx;
+ test_spline->get_range(minx, maxx);
+ if ((minx != x0[0]) || (maxx != x0[NP - 1])) {
+ cout << minx << " " << x0[0] << " | " << maxx << " "
+ << x0[NP - 1] << endl;
+ return -1;
+ }
+ }
+ {
+ cout << "test evaluate()" << endl;
+ double xn = 0.02;
+ double yn = test_spline->evaluate(xn);
+ cout << xn << " " << yn << endl;
+ yn = test_spline->evaluate(x0[5]);
+ cout << xn << " " << yn << " " << y0[5] << endl;
+ }
+ {
+ cout << "test evaluate()" << endl;
+ double xn = 0.02;
+ double yn = test_spline->evaluate(xn);
+ cout << xn << " " << yn << endl;
+ yn = test_spline->evaluate(x0[5]);
+ cout << xn << " " << yn << " " << y0[5] << endl;
+ }
+ {
+ cout << "test evaluate()for plot" << endl;
+ doubleVector X;
+ doubleVector Y;
+ double range = x0[NP - 1];
+ double nsteps = 300;
+ for (int i = 0; i <= nsteps; i++) {
+ X.push_back((range * i) / nsteps);
+ Y.push_back(test_spline->evaluate(X[i]));
+ }
+ //print out data
+ cout << "---------------------------" << endl;
+ for (int i = 0; i < NP; i++)
+ cout << x0[i] << "\t" << y0[i] << endl;
+ cout << "---------------------------" << endl;
+ for (int i = 0; i < nsteps; i++)
+ cout << X[i] << "\t" << Y[i] << endl;
+ cout << "+++++++++++++++++++++++++++" << endl;
+ }
+ {
+ cout << "test evaluate()for doubleVector" << endl;
+ doubleVector X;
+ doubleVector Y;
+ double range = x0[NP - 1];
+ double nsteps = 300;
+ for (int i = 0; i <= nsteps; i++)
+ X.push_back((range * i) / nsteps);
+ Y = test_spline->evaluate(X);
+ //print out data
+ cout << "###########################" << endl;
+ for (int i = 0; i < nsteps; i++)
+ cout << X[i] << "\t" << Y[i] << endl;
+ cout << "***************************" << endl;
+ }
+ {
+ cout << "test range exception" << endl;
+ try {
+ test_spline->evaluate(-10.0);
+ cerr << "range_error not thrown!" << endl;
+ } catch (std::range_error &ex) {
+ cout << "OK, catched expected exception: " << ex.what() << endl;
+ } catch (...) {
+ cerr << "unexpected exception thrown!" << endl;
+ }
+ }
+ {
+ cout << "test constructor with unequal lengths" << endl;
+ try {
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(3.0);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ Spline spl(x, y, NAN, 0.0);
+ cerr << "length_error not thrown!" << endl;
+ } catch (std::length_error &ex) {
+ cout << "OK, catched expected exception: " << ex.what() << endl;
+ } catch (...) {
+ cerr << "unexpected exception thrown!" << endl;
+ }
+ }
+ {
+ cout << "test constructor with out of order data" << endl;
+ try {
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(1.5);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ y.push_back(6.0);
+ Spline spl(x, y, NAN, 0.0);
+ cerr << "domain_error not thrown!" << endl;
+ } catch (std::domain_error &ex) {
+ cout << "OK, catched expected exception: " << ex.what() << endl;
+ } catch (...) {
+ cerr << "unexpected exception thrown!" << endl;
+ }
+ }
+ {
+ cout << "Test copy constructor" << endl;
+ Spline copy_spline(*test_spline);
+ double minx, maxx, copyminx, copymaxx;
+ test_spline->get_range(minx, maxx);
+ copy_spline.get_range(copyminx, copymaxx);
+ if ((minx != copyminx) || (maxx != copymaxx)) {
+ cerr << "COPY FAILED!" << endl;
+ return -1;
+ }
+ double range = x0[NP - 1];
+ double nsteps = 100;
+ for (int i = 0; i < nsteps; i++) {
+ double val = (range * i) / nsteps;
+ double y, ycopy;
+ ycopy = copy_spline.evaluate(val);
+ y = test_spline->evaluate(val);
+ if (y != ycopy) {
+ cerr << "COPY FAILED " << val << " " << y << " " << ycopy
+ << endl;
+ }
+ }
+ }
+ {
+ cout << "Test copy assignment operator" << endl;
+ Spline oper_spline = *test_spline;
+ double minx, maxx, operminx, opermaxx;
+ test_spline->get_range(minx, maxx);
+ oper_spline.get_range(operminx, opermaxx);
+ if ((minx != operminx) || (maxx != opermaxx)) {
+ cerr << "ASSIGN COPY FAILED!" << endl;
+ return -1;
+ }
+ double range = x0[NP - 1];
+ double nsteps = 100;
+ for (int i = 0; i < nsteps; i++) {
+ double val = (range * i) / nsteps;
+ double y, yoper;
+ yoper = oper_spline.evaluate(val);
+ y = test_spline->evaluate(val);
+ if (y != yoper) {
+ cerr << "COPY FAILED " << val << " " << y << " " << yoper
+ << endl;
+ }
+ }
+ }
+ cout << "test delete" << endl;
+ delete test_spline;
+ {
+ cout << "test leaks - new/delete" << endl;
+ Spline *dynspl;
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(3.5);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ y.push_back(6.0);
+ for (int i = 0; i < 100; i++) {
+ dynspl = new Spline(x, y, NAN, 0.0);
+ delete dynspl;
+ }
+ }
+ {
+ cout << "test leaks - stack" << endl;
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(3.5);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ y.push_back(6.0);
+ for (int i = 0; i < 100; i++) {
+ Spline stspl(x, y, NAN, 0.0);
+ }
+ }
+
+ {
+ cout << "test end conditions" << endl;
+ std::ofstream outfile;
+ outfile.open("splineval.csv", ios::out); //CHECK success
+ if ((outfile.rdstate() & ifstream::failbit) != 0)
+ return -1;
+ doubleVector x, y;
+ x.push_back(0.0);
+ y.push_back(0.0);
+ x.push_back(1.0);
+ y.push_back(1.0);
+ x.push_back(1.5);
+ y.push_back(1.0);
+ x.push_back(2.0);
+ y.push_back(1.0);
+ x.push_back(3.0);
+ y.push_back(0.0);
+ Spline nsp(x, y); //natural spline
+ Spline zsp(x, y, 0.0, 0.0); // spline with 0 derivate start and stop
+ Spline nzsp(x, y, NAN, 0.0); // spline with 0 derivate at stop
+ Spline znsp(x, y, 0.0, NAN); // spline with 0 derivate at start
+ Spline ddsp(x, y, -0.5, 0.5); // spline with derivatives at start and stop
+ double X;
+ double Yn, Yz, Ynz, Yzn, Ydd;
+ double range = 3.0;
+ int nsteps = 300;
+ for (int i = 0; i <= nsteps; i++) {
+ X = ((range * i) / nsteps);
+ Yn = nsp.evaluate(X);
+ Yz = zsp.evaluate(X);
+ Ynz = nzsp.evaluate(X);
+ Yzn = znsp.evaluate(X);
+ Ydd = ddsp.evaluate(X);
+ outfile << X << "\t\t" << Yn << "\t\t" << Yz << "\t\t" << Ynz
+ << "\t\t" << Yzn << "\t\t" << Ydd << std::endl;
+ }
+ outfile.close();
+ }
+ }
+ //#######################################################################
+ {
+ cout << "test spline::inverse_evaluate" << endl;
+ {
+
+ doubleVector x, y;
+ x.push_back(10.0);
+ y.push_back(83.6);
+ x.push_back(12);
+ y.push_back(68.3);
+ x.push_back(14);
+ y.push_back(55.7);
+ x.push_back(16);
+ y.push_back(45.5);
+ x.push_back(18);
+ y.push_back(37.3);
+ x.push_back(20);
+ y.push_back(30.7);
+ x.push_back(22);
+ y.push_back(25.4);
+ x.push_back(24);
+ y.push_back(21.2);
+ x.push_back(26);
+ y.push_back(17.9);
+ x.push_back(28);
+ y.push_back(15.2);
+ x.push_back(30);
+ y.push_back(13.1);
+ x.push_back(32);
+ y.push_back(11.4);
+ x.push_back(34);
+ y.push_back(10.1);
+ x.push_back(36);
+ y.push_back(9);
+ x.push_back(40);
+ y.push_back(7.5);
+ x.push_back(50);
+ y.push_back(5.77);
+ Spline *lambdasplie = new Spline(x, y);
+ double gmin, gmax;
+
+ lambdasplie->get_range(gmin, gmax);
+ double max = lambdasplie->evaluate(gmin);
+ double min = lambdasplie->evaluate(gmax);
+ double initial = min + ((max - min) / 2);
+
+ double l = min;
+ double gap, linv;
+ cout << "l\tgap\tlinverse" << endl;
+ cout.precision(8);
+ while (l < max) {
+ gap = lambdasplie->inverse_evaluate(l, initial);
+ linv = lambdasplie->evaluate(gap);
+ cout << l << "\t" << gap << "\t" << linv << endl;
+ l += 0.05;
+ }
+ delete lambdasplie;
+ }
+ }
+ //###############################################################
+ {
+ cout << "test of periodicSpline interpolator class" << endl;
+ doubleVector x0, y0;
+ const int NP = 15;
+ const double k = (2 * M_PI) / NP;
+ for (int i = 0; i <= NP; i++) {
+ x0.push_back(i * k);
+ y0.push_back(sin(i * k));
+ }
+ periodicSpline *test_spline = new periodicSpline(x0, y0);
+ {
+ doubleVector a;
+ doubleVector b;
+ test_spline->get_base_points(a, b);
+ cout << "test get_base_points()" << endl;
+ for (unsigned int i = 0; i < a.size() - 1; i++)
+ if ((x0[i] != a[i]) || (y0[i] != b[i])) {
+ cout << i << " " << x0[i] << " " << a[i] << " | " << y0[i]
+ << " " << b[i] << endl;
+ return -1;
+ }
+
+ }
+ {
+ cout << "test get_range()" << endl;
+ double minx, maxx;
+ test_spline->get_range(minx, maxx);
+ int I = x0.size() - 1;
+ if ((minx != x0[0]) || (maxx != x0[I])) {
+ cout << minx << " " << x0[0] << " | " << maxx << " " << x0[I]
+ << endl;
+ return -1;
+ }
+ }
+ {
+ cout << "test evaluate()" << endl;
+ double xn = 0.02;
+ double yn = test_spline->evaluate(xn);
+ cout << xn << " " << yn << endl;
+ yn = test_spline->evaluate(x0[5]);
+ cout << xn << " " << yn << " " << y0[5] << endl;
+ }
+ {
+ cout << "test evaluate()" << endl;
+ double xn = 0.02;
+ double yn = test_spline->evaluate(xn);
+ cout << xn << " " << yn << endl;
+ yn = test_spline->evaluate(x0[5]);
+ cout << xn << " " << yn << " " << y0[5] << endl;
+ }
+ {
+ cout << "test evaluate()for plot" << endl;
+ doubleVector X;
+ doubleVector Y;
+ double minx;
+ double range;
+ test_spline->get_range(minx, range);
+ double nsteps = 100;
+ for (int i = 0; i <= nsteps; i++) {
+ X.push_back((range * i) / nsteps);
+ Y.push_back(test_spline->evaluate(X[i]));
+ }
+ //print out data
+ cout << "---------------------------" << endl;
+ for (unsigned int i = 0; i < x0.size(); i++)
+ cout << x0[i] << "\t" << y0[i] << endl;
+ cout << "---------------------------" << endl;
+ for (int i = 0; i <= nsteps; i++)
+ cout << X[i] << "\t" << Y[i] << endl;
+ cout << "+++++++++++++++++++++++++++" << endl;
+ }
+ {
+ cout << "test evaluate()for doubleVector" << endl;
+ doubleVector X;
+ doubleVector Y;
+ double minx;
+ double range;
+ test_spline->get_range(minx, range);
+ double nsteps = 100;
+ for (int i = 0; i <= nsteps; i++)
+ X.push_back((range * i) / nsteps);
+ Y = test_spline->evaluate(X);
+ //print out data
+ cout << "###########################" << endl;
+ for (int i = 0; i <= nsteps; i++)
+ cout << X[i] << "\t" << Y[i] << endl;
+ cout << "***************************" << endl;
+ }
+ {
+ cout << "test range exception" << endl;
+ try {
+ test_spline->evaluate(-10.0);
+ cerr << "range_error not thrown!" << endl;
+ } catch (std::range_error &ex) {
+ cout << "OK, catched expected exception: " << ex.what() << endl;
+ } catch (...) {
+ cerr << "unexpected exception thrown!" << endl;
+ }
+ }
+ {
+ cout << "test constructor with unequal lengths" << endl;
+ try {
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(3.0);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ periodicSpline spl(x, y);
+ cerr << "length_error not thrown!" << endl;
+ } catch (std::length_error &ex) {
+ cout << "OK, catched expected exception: " << ex.what() << endl;
+ } catch (...) {
+ cerr << "unexpected exception thrown!" << endl;
+ }
+ }
+ {
+ cout << "test constructor with out of order data" << endl;
+ try {
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(1.5);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ y.push_back(6.0);
+ periodicSpline spl(x, y);
+ cerr << "domain_error not thrown!" << endl;
+ } catch (std::domain_error &ex) {
+ cout << "OK, catched expected exception: " << ex.what() << endl;
+ } catch (...) {
+ cerr << "unexpected exception thrown!" << endl;
+ }
+ }
+ {
+ cout << "Test copy constructor" << endl;
+ periodicSpline copy_spline(*test_spline);
+ double minx, maxx, copyminx, copymaxx;
+ test_spline->get_range(minx, maxx);
+ copy_spline.get_range(copyminx, copymaxx);
+ if ((minx != copyminx) || (maxx != copymaxx)) {
+ cerr << "COPY FAILED!" << endl;
+ return -1;
+ }
+ double range;
+ test_spline->get_range(minx, range);
+ double nsteps = 100;
+ for (int i = 0; i < nsteps; i++) {
+ double val = (range * i) / nsteps;
+ double y, ycopy;
+ ycopy = copy_spline.evaluate(val);
+ y = test_spline->evaluate(val);
+ if (y != ycopy) {
+ cerr << "COPY FAILED " << val << " " << y << " " << ycopy
+ << endl;
+ }
+ }
+ }
+ {
+ cout << "Test copy operator" << endl;
+ periodicSpline copy_spline(*test_spline);
+ copy_spline = *test_spline;
+ double minx, maxx, copyminx, copymaxx;
+ test_spline->get_range(minx, maxx);
+ copy_spline.get_range(copyminx, copymaxx);
+ if ((minx != copyminx) || (maxx != copymaxx)) {
+ cerr << "COPY FAILED!" << endl;
+ return -1;
+ }
+ double range;
+ test_spline->get_range(minx, range);
+ double nsteps = 100;
+ for (int i = 0; i < nsteps; i++) {
+ double val = (range * i) / nsteps;
+ double y, ycopy;
+ ycopy = copy_spline.evaluate(val);
+ y = test_spline->evaluate(val);
+ if (y != ycopy) {
+ cerr << "COPY FAILED " << val << " " << y << " " << ycopy
+ << endl;
+ }
+ }
+ }
+ cout << "test delete" << endl;
+ delete test_spline;
+ {
+ cout << "test leaks - new/delete" << endl;
+ periodicSpline *dynspl;
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(3.5);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ y.push_back(6.0);
+ for (int i = 0; i < 100; i++) {
+ dynspl = new periodicSpline(x, y);
+ delete dynspl;
+ }
+ }
+ {
+ cout << "test leaks - stack" << endl;
+ doubleVector x, y;
+ x.push_back(1.0);
+ x.push_back(2.0);
+ x.push_back(3.5);
+
+ y.push_back(4.0);
+ y.push_back(5.0);
+ y.push_back(6.0);
+ for (int i = 0; i < 100; i++) {
+ periodicSpline stspl(x, y);
+ }
+ }
+ }
+return 0;
+}