C ALGORITHM 764, COLLECTED ALGORITHMS FROM ACM. C THIS WORK PUBLISHED IN TRANSACTIONS ON MATHEMATICAL SOFTWARE, C VOL. 23, NO. 1, March, 1997, P. 1--15. C #! /bin/sh # This is a shell archive, meaning: # 1. Remove everything above the #! /bin/sh line. # 2. Save the resulting text in a file. # 3. Execute the file with /bin/sh (not csh) to create the files: # Src # Doc # This archive created: Sat Apr 19 22:37:40 1997 export PATH; PATH=/bin:$PATH if test ! -d 'Src' then mkdir 'Src' fi cd 'Src' if test ! -d 'Cubpack++' then mkdir 'Cubpack++' fi cd 'Cubpack++' if test ! -d 'Drivers' then mkdir 'Drivers' fi cd 'Drivers' if test -f 'CONTENTS' then echo shar: will not over-write existing file "'CONTENTS'" else cat << \SHAR_EOF > 'CONTENTS' What you find in the current directory -------------------------------------- CONTENTS : this file Small examples of main programs illustrating the use of different geometries. You can start playing with Cubpack++ by modifying these. ---------------------- ------------------------------ source region output ------ ------ ------ vb1.c : TRIANGLE : vb1.out vb2.c : CIRCLE : vb2.out vb3.c : PLANE_SECTOR (wedge) : vb3.out vb4.c : GENERALIZED_RECTANGLE (parabola) : vb4.out vb5.c : GENERALIZED_RECTANGLE (lemniscate) : vb5.out vb6.c : SEMI_INFINITE_STRIP : vb6.out vb7.c : SEMI_INFINITE_STRIP : vb7.out vb8.c : INFINITE_STRIP : vb8.out vb9.c : RECTANGLE : vb9.out vb10.c : TRIANGLE : vb10.out vb11.c : PLANE_SECTOR (half plane) : vb11.out vb12.c : PLANE_SECTOR (quadrant) : vb12.out vb13.c : House = TRIANGLE + RECTANGLE : vb13.out vb14.c : SEMI_INFINITE_STRIP (see vb6) : vb14.out vb15.c : House (see vb13) : vb15.out vb16.c : House (see vb13) : vb16.out vb17.c : GENERALIZED_RECTANGLE (see vb4) : vb17.out vb18.c : PLANE_SECTOR (Bivariate normal distribution) vb_hankel.c: call of Fortran function from C++, see User Manual. The output files where generated using the xlC compiler on an IBM RS6000 running AIX version 3.2.5 SHAR_EOF fi # end of overwriting check if test -f 'EXAMPLES' then echo shar: will not over-write existing file "'EXAMPLES'" else cat << \SHAR_EOF > 'EXAMPLES' Nr File Integration region -- ----- ---------------------------------- 0 interactive choice 1 vb1.c TRIANGLE 2 vb2.c CIRCLE 3 vb3.c PLANE_SECTOR (wedge) 4 vb4.c GENERALIZED_RECTANGLE (parabola) 5 vb5.c GENERALIZED_RECTANGLE (lemniscate) 6 vb6.c SEMI_INFINITE_STRIP 7 vb7.c SEMI_INFINITE_STRIP 8 vb8.c INFINITE_STRIP 9 vb9.c RECTANGLE 10 vb10.c TRIANGLE 11 vb11.c PLANE_SECTOR (half plane) 12 vb12.c PLANE_SECTOR (quadrant) 13 vb13.c House = TRIANGLE + RECTANGLE 14 vb14.c SEMI_INFINITE_STRIP (variant of vb6) 15 vb15.c House (variant of vb13) 16 vb16.c House (variant of vb13) 17 vb17.c GENERALIZED_RECTANGLE (variant of vb4) 18 vb18.c PLANE_SECTOR (Bivariate normal distribution) SHAR_EOF fi # end of overwriting check if test -f 'vb1.c' then echo shar: will not over-write existing file "'vb1.c'" else cat << \SHAR_EOF > 'vb1.c' #include #include real f(const Point& p) { real x=p.X(); return x*x; } int main() { Point p1(0,0), p2(1,1), p3(2,0); TRIANGLE T(p1,p2,p3); cout << "The integral is " << Integrate(f,T) << endl; return 0; } SHAR_EOF fi # end of overwriting check if test -f 'vb10.c' then echo shar: will not over-write existing file "'vb10.c'" else cat << \SHAR_EOF > 'vb10.c' // Example 20 from ditamo // exact value = 4/15 #include #include real f(const Point& p) { real x=p.X() , y=p.Y() ; return sqrt(1-x-y); } int main () { Point origin(0,0),p1(1,0),p2(0,1); TRIANGLE ditamo20(origin,p1,p2); EvaluationCounter count; count.Start(); cout <<"The integral is " << Integrate(f,ditamo20,0,0.5e-4,10000); cout <<" with absolute error " << ditamo20.AbsoluteError() << endl; count.Stop(); cout << count.Read() << " function evaluations were used." << endl; return 0; } SHAR_EOF fi # end of overwriting check if test -f 'vb11.c' then echo shar: will not over-write existing file "'vb11.c'" else cat << \SHAR_EOF > 'vb11.c' // This example assumes that math.h contains the function erfc. // The erfc function is only used to compute the exact value. // Problem posted in sci.math.num-analysis on 9 Feb 1995 // by N. Sukumar (E-mail: n-sukumar@nwu.edu), // Theoretical & Applied Mechanics, Northwestern U, Evanston IL 60208. // Exact value = 1/2 - (1/2) * erf (h/sqrt(2)) // where h = b / sqrt (1 + a^2) #include #include real f(const Point& p) { real x=p.X() , y=p.Y(); return exp(-(x*x+y*y)/2.0)/(2*M_PI); } int main () { real a=1,b=1,h=b/sqrt(1.0+a*a),t=h/sqrt(2.0); Point A(0,b),B(1,a+b),D=2*A-B; PLANE_SECTOR halfplane(A,B,D); EvaluationCounter count; count.Start(); cout <<"The estimated integral is " << Integrate(f,halfplane,0,0.5e-6)< 'vb12.c' // Example 34 from DITAMO // The function is rescaled such that the solution is 1. #include #include real f(const Point& p) { real x=p.X() , y=p.Y(),z , a=1.0/3.0 , b=2.0/3.0 ; if (x*y < 1e-15) return 0.0; z = x + y + 1/(x*y); if ((z > 100) || (z < -100)) { return 0.0;} else { return exp(-z) / (pow(x,b)*pow(y,a)*0.1806075059054343159);} } int main () { // Using the first constructor: Point origin(0,0),p1(1,0),p2(0,1); PLANE_SECTOR quadrant(origin,p1,p2); // Using the second constructor: // Point origin(0,0); // PLANE_SECTOR quadrant(origin,0.0, 0.0, M_PI/2); EvaluationCounter count; count.Start(); cout <<"The integral is " << Integrate(f,quadrant,0,0.5e-4); cout <<" with absolute error " << quadrant.AbsoluteError() << endl; count.Stop(); cout << count.Read() << " function evaluations were used." << endl; return 0; } SHAR_EOF fi # end of overwriting check if test -f 'vb13.c' then echo shar: will not over-write existing file "'vb13.c'" else cat << \SHAR_EOF > 'vb13.c' // Example suggested by a referee, who called this `reasonable'. // One can easily see, e.g. by substituting y=1, that this // integral involves severe difficulties such as fast oscillations // in the neighbourhood of the line x+y=2, and discontinuous derivatives. #include #include #define sqr(x) ((x)*(x)) static const real sqrt_PI = sqrt(M_PI), p35_PI = pow(M_PI,0.35e0); real f(const Point& p) { real x=p.X() , y=p.Y() ; return ( sqrt( sqr(x-sqrt_PI) + sqr(y-1) + 0.00001) + pow(pow(y-p35_PI,4) + sqr(x-1) + 0.00002,1.0/3.0) ) * sin( (x-y-sqrt_PI+p35_PI)/(sqr(x+y-2) + 0.01) ); } int main () { TRIANGLE Roof( Point(0,1), Point(2,1), Point(1,2)); RECTANGLE Walls( Point(0,0), Point(0,1), Point(2,0)); REGION_COLLECTION House = Walls + Roof; EvaluationCounter count; count.Start(); cout <<"The integral is " << Integrate(f,House,0,0.5e-1,1000000); cout <<" with absolute error " << House.AbsoluteError() << endl; count.Stop(); cout << count.Read() << " function evaluations were used." << endl; return 0; } SHAR_EOF fi # end of overwriting check if test -f 'vb14.c' then echo shar: will not over-write existing file "'vb14.c'" else cat << \SHAR_EOF > 'vb14.c' // Example 28 from ditamo #include #include real f(const Point& p) { real y=p.Y(); return 1/(y*y); } int main () { Point A(2,-1), B(1,-1); SEMI_INFINITE_STRIP langelat(A,B); real Integral_est, Error_est; Boolean Success; EvaluationCounter TikTak; TikTak.Start(); cout.setf(ios::scientific,ios::floatfield); do { Integrate(f,langelat,Integral_est,Error_est,Success,0,0.5e-7,10000); cout <<"The integral is " << Integral_est; cout <<" with estimated absolute error " << Error_est << endl; cout <<"The real error is " << Integral_est - 1 < 'vb15.c' // The `House' example with on each subregion a different integrand. #include #include real f1(const Point& p) { return ( p.X()*p.X() ); } real f2(const Point& p) { return ( p.Y()*p.Y() ); } int main () { TRIANGLE Roof( Point(0,1), Point(2,1), Point(1,2)); RECTANGLE Walls( Point(0,0), Point(0,1), Point(2,0)); REGION_COLLECTION House = Walls + Roof; Roof.LocalIntegrand(f1); Walls.LocalIntegrand(f2); EvaluationCounter count; count.Start(); cout <<"The integral is " << Integrate(House); cout <<" with estimated absolute error " << House.AbsoluteError() << endl; cout <<"The contributions of the subregions are:"< 'vb16.c' // Example suggested by a referee, who called this `reasonable'. // One can easily see, e.g. by substituting y=1, that this // integral involves severe difficulties such as fast oscillations // in the neighbourhood of the line x+y=2, and discontinuous derivatives. #include #include #define sqr(x) ((x)*(x)) static const real sqrt_PI = sqrt(M_PI), p35_PI = pow(M_PI,0.35e0); real f(const Point& p) { real x=p.X() , y=p.Y() ; return ( sqrt( sqr(x-sqrt_PI) + sqr(y-1) + 0.00001) + pow(pow(y-p35_PI,4) + sqr(x-1) + 0.00002,1.0/3.0) ) * sin( (x-y-sqrt_PI+p35_PI)/(sqr(x+y-2) + 0.01) ); } int main () { TRIANGLE T1( Point(0,0), Point(1,0), Point(0,1)), T2( Point(1,0), Point(2,0), Point(1.5,0.5)), T3( Point(2,0), Point(2,1), Point(1.5,0.5)); RECTANGLE R1( Point(1,0), Point(0,1), Point(2,1)); REGION_COLLECTION House1 = T1 + T2 + T3 + R1; RECTANGLE R3( Point(1,0), Point(0,1), Point(2,1)), R2( Point(1.8,0), Point(1.4,0.4), Point(2,0.2)); TRIANGLE T4( Point(1,0), Point(1.8,0), Point(1.4,0.4)), T7( Point(0,0), Point(1,0), Point(0,1)), T5( Point(1.8,0), Point(2,0), Point(2,0.2)), T6( Point(2,0.2), Point(2,1), Point(1.6,0.6)); REGION_COLLECTION House2 = T7 + R3 + R2 + T4 + T5 + T6; EvaluationCounter count; count.Start(); cout <<"The integral is " << Integrate(f,House1,0,0.5e-1,1000000); cout <<" with absolute error " << House1.AbsoluteError() << endl; count.Stop(); cout << count.Read() << " function evaluations were used." << endl; cout << "Contribution of the different parts:" << endl; cout << "R1: integral and error are " << R1.Integral() << " and " << R1.AbsoluteError() < 'vb17.c' // Example 23 from ditamo // The function is rescaled such that the solution is 1. #include #include real f(const Point& p) { real x=p.X() , y=p.Y(); return 0.25/sqrt(x*(2.0-x)); } real Height(const Point& p) { return sqrt(4-2*p.X()); } int main () { real Integral_est, Error_est; Boolean Success; Point A(0,0), B(2,0); GENERALIZED_RECTANGLE parabola(Height,A,B); do { Integrate(f,parabola,Integral_est,Error_est,Success,0,0.5e-1,10000); cout <<"The integral is " << Integral_est < 'vb18.c' // Bivariate normal distribution // Note that special purpose methods exist for this type of integral. #include #include // The integral depends on the following 3 parameters: real rho, a, b; real N; // This is used to save work real f(const Point& p) { real x=p.X() , y=p.Y(); return exp(N*(-x*x+2*rho*x*y-y*y)) ; } int main () { EvaluationCounter count; real c; // Read the 3 parameters from standard input: cout << "Give a,b and rho ( |rho| < 1 ): "; cin >> a >> b >> rho; // Define the region of integration: Point Center(a,b); PLANE_SECTOR quadrant(Center,0.0, M_PI , 3*M_PI/2); c = 1.0/(2.0*M_PI*sqrt(1.0 - rho*rho)); N = 1.0/(2.0 - 2.0*rho*rho); count.Start(); cout <<"The integral is " << Integrate(f,quadrant,0,0.5e-6)*c; cout <<" with absolute error " << quadrant.AbsoluteError()*c << endl; count.Stop(); cout << count.Read() << " function evaluations were used." << endl; return 0; } SHAR_EOF fi # end of overwriting check if test -f 'vb2.c' then echo shar: will not over-write existing file "'vb2.c'" else cat << \SHAR_EOF > 'vb2.c' #include #include real f(const Point& p) { real r=p.Length(); return 1.0/(1.0+r*r*r*r); } int main () { Point origin(0,0); real radius=1; CIRCLE cir(origin,radius); cout <<"The integral is " << Integrate(f, cir, 0, 1.0e-6, 10000); cout <<" with absolute error " << cir.AbsoluteError() << endl; return 0; } SHAR_EOF fi # end of overwriting check if test -f 'vb3.c' then echo shar: will not over-write existing file "'vb3.c'" else cat << \SHAR_EOF > 'vb3.c' // Example 33 from ditamo // Exact value = arctan(2) // Note that the "exact" value in the original paper is wrong. #include #include #include real f(const Point& p) { real x=p.X() , y=p.Y(); return exp(0.5*(-x*x -y*y)); } int main () { Point origin(0,0); real innerradius=0, alfa=0, beta=atan(2.0); PLANE_SECTOR wedge(origin,innerradius,alfa,beta); EvaluationCounter count; cout.setf(ios::scientific,ios::floatfield); cout<<"req. rel. error est integral est error abs error evaluations" <1e-12; req_err/=10) { cout << setprecision(1) << " <" << req_err <<" " << setprecision(10) << Integrate(f,wedge,0,req_err) << " "; cout << setprecision(1) << wedge.AbsoluteError() << " " << fabs(wedge.Integral() - atan(2.0)) << " " << count.Read() << endl; } count.Stop(); return 0; } SHAR_EOF fi # end of overwriting check if test -f 'vb4.c' then echo shar: will not over-write existing file "'vb4.c'" else cat << \SHAR_EOF > 'vb4.c' // Example 23 from ditamo // The function is rescaled such that the solution is 1. #include #include real f(const Point& p) { real x=p.X() , y=p.Y(); return 0.25/sqrt(x*(2.0-x)); } real Height(const Point& p) { return sqrt(4-2*p.X()); } int main () { Point A(0,0), B(2,0); GENERALIZED_RECTANGLE parabola(Height,A,B); Chrono TikTak; TikTak.Start(); cout <<"The integral is " << Integrate(f,parabola,0,0.5e-4); TikTak.Stop(); cout <<" with absolute error " << parabola.AbsoluteError() << endl; cout <<"Elapsed Time: " << TikTak.Read()/1000<<" seconds" << endl; return 0; } SHAR_EOF fi # end of overwriting check if test -f 'vb5.c' then echo shar: will not over-write existing file "'vb5.c'" else cat << \SHAR_EOF > 'vb5.c' // Example 5 from ditamo #include #include real f(const Point& p) { real x=p.X() , y=p.Y(); return p.Y(); } real Height(const Point& p) { real x=p.X(); return sqrt(4*cos(x)*cos(x)/9 - 16*sin(x)*sin(x)/25); } int main () { Point A(0,0), B(atan(5.0/6.0),0); GENERALIZED_RECTANGLE lemniscate(Height,A,B); EvaluationCounter TikTak; TikTak.Start(); // cout.setf(ios::left,ios::adjustfield); cout.setf(ios::scientific,ios::floatfield); cout <<"The integral is " << Integrate(f,lemniscate,0,0.5e-7,10000); cout <<" with estimated absolute error " << lemniscate.AbsoluteError() << endl; cout <<"The real error is " << lemniscate.Integral() - (2-11*atan(5.0/6.0)/15)/15 < 'vb6.c' // Example 28 from ditamo #include #include real f(const Point& p) { real y=p.Y(); return 1/(y*y); } int main () { Point A(2,-1), B(1,-1); SEMI_INFINITE_STRIP langelat(A,B); EvaluationCounter TikTak; TikTak.Start(); // cout.setf(ios::left,ios::adjustfield); cout.setf(ios::scientific,ios::floatfield); cout <<"The integral is " << Integrate(f,langelat,0,0.5e-7,20000); cout <<" with estimated absolute error " << langelat.AbsoluteError() << endl; cout <<"The real error is " << langelat.Integral() - 1 < 'vb7.c' // Example 29 from ditamo #include #include real f(const Point& p) { return fabs(p.Y())*exp(-p.Y()*p.Y()/2); } int main () { Point A(0,0), B(1.107148717794090503,0); SEMI_INFINITE_STRIP langelat(A,B); EvaluationCounter TikTak; TikTak.Start(); // cout.setf(ios::left,ios::adjustfield); cout.setf(ios::scientific,ios::floatfield); cout <<"The integral is " << Integrate(f,langelat,0,0.5e-4,10000); cout <<" with estimated absolute error " << langelat.AbsoluteError() << endl; TikTak.Stop(); cout<<"Number of evaluations = "< 'vb8.c' // Example #include #include #include real f(const Point& p) { return exp(-p.Y()*p.Y()); } int main () { Point A(0,0), B(1,0); INFINITE_STRIP langelat(A,B); EvaluationCounter TikTak; TikTak.Start(); // cout.setf(ios::left,ios::adjustfield); cout.setf(ios::scientific,ios::floatfield); cout <<"The integral is " << Integrate(f,langelat,0,0.5e-3,10000); cout <<" with estimated absolute error " << langelat.AbsoluteError() << endl; cout <<"The real error is " << langelat.Integral() - sqrt(M_PI) < 'vb9.c' // This example comes from Problem 94-7 by D.E. Loper // Siam Review, Vol. 36 (1994) page 277 #include #include #include static real q; real A(const Point& p) { real m=p.X() , f=p.Y(); real t=(1-m*m)*(1-m*m), s = sin(f); s = s*s; return (t*s)/(t*s*s+q*q*m*m); } real B(const Point& p) { real m=p.X() , f=p.Y(); real t=(1-m*m), s = sin(f); s = s*s; return (t*s*(t*s+m*m))/(t*t*s*s+q*q*m*m); } int main () { Point Origin(0,0), mu(1,0),phi(0,M_PI/2); cout<<"q A(q) B(q) A(q)/B(q)"< 'vb_hankel.c' #include #include extern "C" {void zbesh_ (real&, real&, real&, int&, int&, int&, real[], real[], int&, int&);} real AbsHankel1 ( const Point& z) { real x=z.X(), y=z.Y(), cr[10], ci[10], fnu=0; int kode=1, m=1, n=1, nz, ierr; zbesh_(x,y,fnu,kode,m,n,cr,ci,nz,ierr); return sqrt(cr[0]*cr[0]+ci[0]*ci[0]); } real AbsHankel2 ( const Point& z) { real x=(z.X()*0.8+z.Y()*0.6), y=(0.8*z.Y()-z.X()*0.6), cr[10], ci[10], fnu=0; int kode=1, m=1, n=1, nz, ierr; zbesh_(x,y,fnu,kode,m,n,cr,ci,nz,ierr); return sqrt(cr[0]*cr[0]+ci[0]*ci[0]); } int main() { Point O(0,0), A(-1,0),B(-0.8,-0.6); EvaluationCounter count; cout << " Omega = 1"< 'vb1.out' The integral is 1.16667 SHAR_EOF fi # end of overwriting check if test -f 'vb10.out' then echo shar: will not over-write existing file "'vb10.out'" else cat << \SHAR_EOF > 'vb10.out' The integral is 0.266677 with absolute error 1.16825e-05 925 function evaluations were used. SHAR_EOF fi # end of overwriting check if test -f 'vb11.out' then echo shar: will not over-write existing file "'vb11.out'" else cat << \SHAR_EOF > 'vb11.out' The estimated integral is 0.23975 The exact value is 0.23975 17620 function evaluations were used. SHAR_EOF fi # end of overwriting check if test -f 'vb12.out' then echo shar: will not over-write existing file "'vb12.out'" else cat << \SHAR_EOF > 'vb12.out' The integral is 1 with absolute error 4.27146e-05 5091 function evaluations were used. SHAR_EOF fi # end of overwriting check if test -f 'vb13.out' then echo shar: will not over-write existing file "'vb13.out'" else cat << \SHAR_EOF > 'vb13.out' The integral is -0.40791 with absolute error 0.0203924 266548 function evaluations were used. SHAR_EOF fi # end of overwriting check if test -f 'vb14.out' then echo shar: will not over-write existing file "'vb14.out'" else cat << \SHAR_EOF > 'vb14.out' The integral is 1.000000e+00 with estimated absolute error 3.035122e-05 The real error is 1.779739e-07 ------------------------------------------------- The integral is 1.000000e+00 with estimated absolute error 4.073745e-07 The real error is 5.862995e-10 ------------------------------------------------- The integral is 1.000000e+00 with estimated absolute error 6.832923e-08 The real error is 2.741229e-11 ------------------------------------------------- The integral is 1.000000e+00 with estimated absolute error 4.986300e-08 The real error is 2.816547e-11 ------------------------------------------------- Total number of evaluations = 32067 SHAR_EOF fi # end of overwriting check if test -f 'vb15.out' then echo shar: will not over-write existing file "'vb15.out'" else cat << \SHAR_EOF > 'vb15.out' The integral is 1.83333 with estimated absolute error 2.03541e-14 The contributions of the subregions are: Walls: integral = 0.666667, error = 7.40149e-15 Roof: integral = 1.16667, error = 1.29526e-14 The exact value is 2/3 + 7/6 = 11/6 74 function evaluations were used. SHAR_EOF fi # end of overwriting check if test -f 'vb16.out' then echo shar: will not over-write existing file "'vb16.out'" else cat << \SHAR_EOF > 'vb16.out' The integral is -0.407393 with absolute error 0.0203507 105968 function evaluations were used. Contribution of the different parts: R1: integral and error are -0.283031 and 0.00816722 T1: integral and error are -0.221416 and 1.09496e-06 T2: integral and error are 0.068406 and 0.0055473 T3: integral and error are 0.0286485 and 0.00663513 The integral is -0.407853 with absolute error 0.0203784 94572 function evaluations were used. SHAR_EOF fi # end of overwriting check if test -f 'vb17.out' then echo shar: will not over-write existing file "'vb17.out'" else cat << \SHAR_EOF > 'vb17.out' The integral is 0.996521 with estimated absolute error 0.151082 ------------------------------------------------- The integral is 0.997535 with estimated absolute error 0.107081 ------------------------------------------------- The integral is 0.997926 with estimated absolute error 0.0900689 ------------------------------------------------- The integral is 0.998255 with estimated absolute error 0.0758062 ------------------------------------------------- The integral is 0.998393 with estimated absolute error 0.0697915 ------------------------------------------------- The integral is 0.998532 with estimated absolute error 0.0637768 ------------------------------------------------- The integral is 0.99867 with estimated absolute error 0.0577621 ------------------------------------------------- The integral is 0.998765 with estimated absolute error 0.0536343 ------------------------------------------------- The integral is 0.998814 with estimated absolute error 0.0515078 ------------------------------------------------- The integral is 0.998851 with estimated absolute error 0.0499129 ------------------------------------------------- SHAR_EOF fi # end of overwriting check if test -f 'vb18.out' then echo shar: will not over-write existing file "'vb18.out'" else cat << \SHAR_EOF > 'vb18.out' Give a,b and rho ( |rho| < 1 ): 1.0 1.0 0.5 The integral is 0.745204 with absolute error 3.68313e-07 17194 function evaluations were used. SHAR_EOF fi # end of overwriting check if test -f 'vb2.out' then echo shar: will not over-write existing file "'vb2.out'" else cat << \SHAR_EOF > 'vb2.out' The integral is 2.4674 with absolute error 5.77166e-08 SHAR_EOF fi # end of overwriting check if test -f 'vb3.out' then echo shar: will not over-write existing file "'vb3.out'" else cat << \SHAR_EOF > 'vb3.out' req. rel. error est integral est error abs error evaluations <5.0e-02 1.1096504402e+00 3.4e-02 2.5e-03 292 <5.0e-03 1.1071935478e+00 4.1e-03 4.5e-05 1181 <5.0e-04 1.1071470438e+00 5.1e-04 1.7e-06 3105 <5.0e-05 1.1071488562e+00 5.3e-05 1.4e-07 5657 <5.0e-06 1.1071486939e+00 5.3e-06 2.4e-08 11239 <5.0e-07 1.1071487183e+00 5.3e-07 5.5e-10 16377 <5.0e-08 1.1071487178e+00 5.5e-08 1.0e-12 25553 <5.0e-09 1.1071487178e+00 5.5e-09 1.1e-11 38577 <5.0e-10 1.1071487178e+00 5.4e-10 5.2e-13 53081 <5.0e-11 1.1071487178e+00 5.5e-11 2.6e-13 73949 <5.0e-12 1.1071487178e+00 5.5e-12 1.8e-15 117794 SHAR_EOF fi # end of overwriting check if test -f 'vb4.out' then echo shar: will not over-write existing file "'vb4.out'" else cat << \SHAR_EOF > 'vb4.out' The integral is 0.99886 with absolute error 0.0495064 Elapsed Time: 0 seconds SHAR_EOF fi # end of overwriting check if test -f 'vb5.out' then echo shar: will not over-write existing file "'vb5.out'" else cat << \SHAR_EOF > 'vb5.out' The integral is 9.936835e-02 with estimated absolute error 7.230197e-14 The real error is 1.387779e-17 Number of evaluations = 37 SHAR_EOF fi # end of overwriting check if test -f 'vb6.out' then echo shar: will not over-write existing file "'vb6.out'" else cat << \SHAR_EOF > 'vb6.out' The integral is 1.000000e+00 with estimated absolute error 4.173926e-07 The real error is 5.863554e-10 Number of evaluations = 20079 SHAR_EOF fi # end of overwriting check if test -f 'vb7.out' then echo shar: will not over-write existing file "'vb7.out'" else cat << \SHAR_EOF > 'vb7.out' The integral is 1.107149e+00 with estimated absolute error 5.263967e-05 Number of evaluations = 5657 SHAR_EOF fi # end of overwriting check if test -f 'vb8.out' then echo shar: will not over-write existing file "'vb8.out'" else cat << \SHAR_EOF > 'vb8.out' The integral is 1.772455e+00 with estimated absolute error 6.976425e-04 The real error is 1.067348e-06 Number of evaluations = 1183 SHAR_EOF fi # end of overwriting check if test -f 'vb9.out' then echo shar: will not over-write existing file "'vb9.out'" else cat << \SHAR_EOF > 'vb9.out' q A(q) B(q) A(q)/B(q) 1 1.513682 0.7568409 2 2 0.9107223 0.4553611 2 3 0.6587322 0.3293661 2 4 0.5176153 0.2588077 2 5 0.4268288 0.2134144 2 6 0.3633539 0.181677 2 7 0.3164155 0.1582078 2 8 0.2802695 0.1401347 2 9 0.2515646 0.1257823 2 10 0.2282109 0.1141055 2 11 0.2088359 0.104418 2 12 0.1925005 0.09625026 2 13 0.1785401 0.08927006 2 14 0.166471 0.08323551 2 15 0.1559327 0.07796635 2 16 0.1466509 0.07332547 2 17 0.1384134 0.06920668 2 18 0.131053 0.06552649 2 19 0.1244366 0.06221831 2 SHAR_EOF fi # end of overwriting check cd .. if test ! -d 'Documentation' then mkdir 'Documentation' fi cd 'Documentation' if test -f 'CONTENTS' then echo shar: will not over-write existing file "'CONTENTS'" else cat << \SHAR_EOF > 'CONTENTS' What you find in the current directory -------------------------------------- CONTENTS : this file COMPILERS : overview of compilers used to test Cubpack++ MACHINE_DEPEND : list of machine dependable source files paper.dvi : copy of submitted paper in dvi format paper.ps : copy of submitted paper in Postscript format userman.dvi : copy of User Manual in dvi format userman.ps : copy of User Manual in Postscript format SHAR_EOF fi # end of overwriting check if test -f 'COMPILERS' then echo shar: will not over-write existing file "'COMPILERS'" else cat << \SHAR_EOF > 'COMPILERS' Systems for which Cubpack++ has proven it works =============================================== Compiler Version Operating System Machine type -------------------------------------------------------------- g++/gcc 2.5.8 Ultrix 4.4 DECstation 5000 SunOS 5.3 Sparcenter 1000 Linux 1.0 i586 2.6.0 HP-UX 2.6.2 SunOS 5.3 Sparcenter 1000 2.6.3* MSDOS 2.7.2 Ultrix 4.4 DECstation 5000 Solaris2.3 Sparcenter 1000 -------------------------------------------------------------- xlC AIX Version 3.2 IBM RS6000 -------------------------------------------------------------- cxx OSF/1 V1.3 DECalpha -------------------------------------------------------------- CC SC3.0.1 SunOS 5.3 Sparcenter 1000 -------------------------------------------------------------- Borland C++ 4.0 Windows -------------------------------------------------------------- Turbo C++ 3.0 MSDOS 5.0 386SX -------------------------------------------------------------- *) DJGPP 1.12m4 Systems for which a workaround is provided ========================================== Compiler Version Operating System Machine type -------------------------------------------------------------- CC SC3.0 SunOS 5.3 Sparcenter 1000 (15 Dec 1993) Workaround activated with compiler option -DSOLARIS Note that this bug was removed from CC in version SC3.0.1. -------------------------------------------------------------- Systems for which Cubpack++ is too complex ========================================== Compiler Version Operating System Machine type -------------------------------------------------------------- CC SC2.0.1 SunOS 4. Sun4 -------------------------------------------------------------- The answer from vendors to these problems is: upgrade. SHAR_EOF fi # end of overwriting check if test -f 'MACHINE_DEPEND' then echo shar: will not over-write existing file "'MACHINE_DEPEND'" else cat << \SHAR_EOF > 'MACHINE_DEPEND' MACHINE DEPENDENT FILES: ------------------------ chrono.c real.h templist.h -> For details we refer to the User Manual, section 11. Users might need to change these. deccxxtp.c -> This file is only used by the cxx compiler on DECalpha OSF/1 V1.3 Users should not change this. gnu262tp.c -> This file is only used by the gnu compiler version at least 2.6.2 Users should not change this. SHAR_EOF fi # end of overwriting check cd .. if test ! -d 'Tools' then mkdir 'Tools' fi cd 'Tools' if test -f 'CONTENTS' then echo shar: will not over-write existing file "'CONTENTS'" else cat << \SHAR_EOF > 'CONTENTS' What you find in the current directory -------------------------------------- CONTENTS : this file c2C : change from suffix .c to .C c2cpp : change from suffix .c to .cpp The following makefile-s are used by quick_install and quick_run. mkfile.decalpha : makefile for cxx compiler on DECalpha OSF mkfile.gn1 : makefile for gnu compiler at most v2.6.0 mkfile.gn2 : makefile for gnu compiler v2.6.2 mkfile.gn3 : makefile for gnu compiler v2.7.0 mkfile.rs6000 : makefile for xlC compiler on RS6000 AIX mkfile.solaris : makefile for CC compiler on Solaris 2.* mkfile.tc : makefile for Turbo C++ under MS-DOS The following Makefile-s are models which can be modified by users. Mk.decalpha : makefile for cxx compiler on DECalpha OSF Mk.gn1 : makefile for gnu compiler at most v2.6.0 Mk.gn2 : makefile for gnu compiler v2.6.2 Mk.gn3 : makefile for gnu compiler v2.7.0 Mk.solaris : makefile for CC compiler on Solaris 2.* Mk.tc : makefile for Turbo C++ under MS-DOS SHAR_EOF fi # end of overwriting check if test -f 'c2C' then echo shar: will not over-write existing file "'c2C'" else cat << \SHAR_EOF > 'c2C' #!/bin/sh SRC_DIR=Src cd ../$SRC_DIR for i in s_adapt.h set.h stack.h vector.h vstack.h heap.h atomic.h div.h passbuck.h pointer.h regproc.h rule.h samediv.h userint.h do cat $i | sed 's/\.c>/.C>/' > tmp mv tmp $i done # for i in `ls *.c` do a=`echo $i | sed "s/\.c$/.C/"` mv $i $a done SHAR_EOF chmod +x 'c2C' fi # end of overwriting check if test -f 'c2cpp' then echo shar: will not over-write existing file "'c2cpp'" else cat << \SHAR_EOF > 'c2cpp' #!/bin/sh SRC_DIR=Src cd ../$SRC_DIR for i in s_adapt.h set.h stack.h vector.h vstack.h heap.h atomic.h div.h passbuck.h pointer.h regproc.h rule.h samediv.h userint.h do cat $i | sed 's/\.c>/.cpp>/' > tmp mv tmp $i done # for i in `ls *.c` do a=`echo $i | sed "s/\.c$/.cpp/"` mv $i $a done SHAR_EOF chmod +x 'c2cpp' fi # end of overwriting check if test -f 'mkfile.decalpha' then echo shar: will not over-write existing file "'mkfile.decalpha'" else cat << \SHAR_EOF > 'mkfile.decalpha' .SUFFIXES: .o .C .C.o: ; $(CC) $(CFLAGS) -c $*.C CC = cxx CFLAGS = -I. libdir = -L. libs = -lcubpack -lm MYTESTPROG = main.o OBJS = \ deccxxtp.o \ point2D.o \ C2.o \ C2dv2.o \ C2dv4.o \ C2interf.o \ C2prc.o \ C2rule13.o \ C2toS2.o \ C2togr.o \ E2tostrp.o \ E2.o \ E2adapt.o \ E2interf.o \ E2sec.o \ E2secitf.o \ E2secprc.o \ S2.o \ S2adapt.o \ S2interf.o \ S2rule13.o \ T2.o \ T2dv4.o \ T2interf.o \ T2rule13.o \ T2tops.o \ atomreg.o \ boolean.o \ chrono.o \ compreg.o \ eval_ctr.o \ geometry.o \ gr.o \ gritf.o \ gs.o \ gsitf.o \ gsprc.o \ integran.o \ integrat.o \ invert.o \ outS2.o \ outS2itf.o \ polC2.o \ polC2itf.o \ polC2prc.o \ ps.o \ psitf.o \ refcount.o \ regcoll.o \ reginfo.o \ region.o \ semstitf.o \ semistrp.o \ sttosmst.o \ strip.o \ stripitf.o \ translat.o \ trnsfrm.o all:: libcubpack.a libcubpack.a: $(OBJS) ar rus libcubpack.a $(OBJS) tst: $(MYTESTPROG) $(CC) $(MYTESTPROG) $(libdir) $(libs) -o tst clean: /bin/rm -f *.o tst core SHAR_EOF fi # end of overwriting check if test -f 'mkfile.gn1' then echo shar: will not over-write existing file "'mkfile.gn1'" else cat << \SHAR_EOF > 'mkfile.gn1' CC = gcc CFLAGS = -x c++ -I. # COPTIONS = -DEGTRUSAGE COPTIONS = libdir = -L. libs = -lcubpack -lm -lg++ MYTESTPROG = main.o .c.o: ; $(CC) $(CFLAGS) -c $*.c chrono.o: chrono.c $(CC) $(CFLAGS) $(COPTIONS) -c chrono.c OBJS = \ point2D.o \ C2.o \ C2dv2.o \ C2dv4.o \ C2interf.o \ C2prc.o \ C2rule13.o \ C2toS2.o \ C2togr.o \ E2tostrp.o \ E2.o \ E2adapt.o \ E2interf.o \ E2sec.o \ E2secitf.o \ E2secprc.o \ S2.o \ S2adapt.o \ S2interf.o \ S2rule13.o \ T2.o \ T2dv4.o \ T2interf.o \ T2rule13.o \ T2tops.o \ atomreg.o \ boolean.o \ chrono.o \ compreg.o \ eval_ctr.o \ geometry.o \ gr.o \ gritf.o \ gs.o \ gsitf.o \ gsprc.o \ integran.o \ integrat.o \ invert.o \ outS2.o \ outS2itf.o \ polC2.o \ polC2itf.o \ polC2prc.o \ ps.o \ psitf.o \ refcount.o \ regcoll.o \ reginfo.o \ region.o \ semstitf.o \ semistrp.o \ sttosmst.o \ strip.o \ stripitf.o \ translat.o \ trnsfrm.o all:: libcubpack.a libcubpack.a: $(OBJS) ar ru libcubpack.a $(OBJS) ranlib libcubpack.a tst: $(MYTESTPROG) $(CC) $(MYTESTPROG) $(libdir) $(libs) -o tst clean: /bin/rm -f *.o tst core SHAR_EOF fi # end of overwriting check if test -f 'mkfile.gn2' then echo shar: will not over-write existing file "'mkfile.gn2'" else cat << \SHAR_EOF > 'mkfile.gn2' CC = gcc CFLAGS = -x c++ -I. -fno-implicit-templates # COPTIONS = -DEGTRUSAGE COPTIONS = libdir = -L. libs = -lcubpack -lm -lg++ MYTESTPROG = main.o .c.o: ; $(CC) $(CFLAGS) -c $*.c chrono.o: chrono.c $(CC) $(CFLAGS) $(COPTIONS) -c chrono.c OBJS = \ point2D.o \ C2.o \ C2dv2.o \ C2dv4.o \ C2interf.o \ C2prc.o \ C2rule13.o \ C2toS2.o \ C2togr.o \ E2tostrp.o \ E2.o \ E2adapt.o \ E2interf.o \ E2sec.o \ E2secitf.o \ E2secprc.o \ S2.o \ S2adapt.o \ S2interf.o \ S2rule13.o \ T2.o \ T2dv4.o \ T2interf.o \ T2rule13.o \ T2tops.o \ atomreg.o \ boolean.o \ chrono.o \ compreg.o \ eval_ctr.o \ geometry.o \ gr.o \ gritf.o \ gs.o \ gsitf.o \ gsprc.o \ integran.o \ integrat.o \ invert.o \ outS2.o \ outS2itf.o \ polC2.o \ polC2itf.o \ polC2prc.o \ ps.o \ psitf.o \ refcount.o \ regcoll.o \ reginfo.o \ region.o \ semstitf.o \ semistrp.o \ sttosmst.o \ strip.o \ stripitf.o \ translat.o \ trnsfrm.o \ gnu262tp.o all:: libcubpack.a libcubpack.a: $(OBJS) ar rus libcubpack.a $(OBJS) tst: $(MYTESTPROG) $(CC) $(MYTESTPROG) $(libdir) $(libs) -o tst clean: /bin/rm -f *.o tst core SHAR_EOF fi # end of overwriting check if test -f 'mkfile.gn3' then echo shar: will not over-write existing file "'mkfile.gn3'" else cat << \SHAR_EOF > 'mkfile.gn3' CC = g++ CFLAGS = -I. -fno-implicit-templates # COPTIONS = -DEGTRUSAGE COPTIONS = libdir = -L. libs = -lcubpack -lm MYTESTPROG = main.o .c.o: ; $(CC) $(CFLAGS) -c $*.c chrono.o: chrono.c $(CC) $(CFLAGS) $(COPTIONS) -c chrono.c OBJS = \ point2D.o \ C2.o \ C2dv2.o \ C2dv4.o \ C2interf.o \ C2prc.o \ C2rule13.o \ C2toS2.o \ C2togr.o \ E2tostrp.o \ E2.o \ E2adapt.o \ E2interf.o \ E2sec.o \ E2secitf.o \ E2secprc.o \ S2.o \ S2adapt.o \ S2interf.o \ S2rule13.o \ T2.o \ T2dv4.o \ T2interf.o \ T2rule13.o \ T2tops.o \ atomreg.o \ boolean.o \ chrono.o \ compreg.o \ eval_ctr.o \ geometry.o \ gr.o \ gritf.o \ gs.o \ gsitf.o \ gsprc.o \ integran.o \ integrat.o \ invert.o \ outS2.o \ outS2itf.o \ polC2.o \ polC2itf.o \ polC2prc.o \ ps.o \ psitf.o \ refcount.o \ regcoll.o \ reginfo.o \ region.o \ semstitf.o \ semistrp.o \ sttosmst.o \ strip.o \ stripitf.o \ translat.o \ trnsfrm.o \ gnu262tp.o all:: libcubpack.a libcubpack.a: $(OBJS) ar rus libcubpack.a $(OBJS) tst: $(MYTESTPROG) $(CC) $(MYTESTPROG) $(libdir) $(libs) -o tst clean: /bin/rm -f *.o tst core SHAR_EOF fi # end of overwriting check if test -f 'mkfile.rs6000' then echo shar: will not over-write existing file "'mkfile.rs6000'" else cat << \SHAR_EOF > 'mkfile.rs6000' CC = xlC -+ #CFLAGS = -O CFLAGS = -I. libdir = libs = -lm MYTESTPROG = main.o .c.o: ; $(CC) $(CFLAGS) -c $*.c OBJS = \ point2D.o \ C2.o \ C2dv2.o \ C2dv4.o \ C2interf.o \ C2prc.o \ C2rule13.o \ C2toS2.o \ C2togr.o \ E2tostrp.o \ E2.o \ E2adapt.o \ E2interf.o \ E2sec.o \ E2secitf.o \ E2secprc.o \ S2.o \ S2adapt.o \ S2interf.o \ S2rule13.o \ T2.o \ T2dv4.o \ T2interf.o \ T2rule13.o \ T2tops.o \ atomreg.o \ boolean.o \ chrono.o \ compreg.o \ eval_ctr.o \ geometry.o \ gr.o \ gritf.o \ gs.o \ gsitf.o \ gsprc.o \ integran.o \ integrat.o \ invert.o \ outS2.o \ outS2itf.o \ polC2.o \ polC2itf.o \ polC2prc.o \ ps.o \ psitf.o \ refcount.o \ regcoll.o \ reginfo.o \ region.o \ semstitf.o \ semistrp.o \ sttosmst.o \ strip.o \ stripitf.o \ translat.o \ trnsfrm.o all:: $(OBJS) tst: $(MYTESTPROG) $(CC) $(CFLAGS) $(MYTESTPROG) $(OBJS) $(libdir) $(libs) -o tst clean: /bin/rm -f *.o run tst core SHAR_EOF fi # end of overwriting check if test -f 'mkfile.solaris' then echo shar: will not over-write existing file "'mkfile.solaris'" else cat << \SHAR_EOF > 'mkfile.solaris' CC = CC # If you have CC version SC3.0 of 15 Dec 1993, then use the following # to activated the workaround # CFLAGS = -I. -DSOLARIS # If you have CC version SC3.0.1 of 13 Jul 1994, then use CFLAGS = -I. libdir = -L. libs = -lcubpack -lm MYTESTPROG = main.o .c.o: ; $(CC) $(CFLAGS) -c $*.c OBJS = \ point2D.o \ C2.o \ C2dv2.o \ C2dv4.o \ C2interf.o \ C2prc.o \ C2rule13.o \ C2toS2.o \ C2togr.o \ E2tostrp.o \ E2.o \ E2adapt.o \ E2interf.o \ E2sec.o \ E2secitf.o \ E2secprc.o \ S2.o \ S2adapt.o \ S2interf.o \ S2rule13.o \ T2.o \ T2dv4.o \ T2interf.o \ T2rule13.o \ T2tops.o \ atomreg.o \ boolean.o \ chrono.o \ compreg.o \ eval_ctr.o \ geometry.o \ gr.o \ gritf.o \ gs.o \ gsitf.o \ gsprc.o \ integran.o \ integrat.o \ invert.o \ outS2.o \ outS2itf.o \ polC2.o \ polC2itf.o \ polC2prc.o \ ps.o \ psitf.o \ refcount.o \ regcoll.o \ reginfo.o \ region.o \ semstitf.o \ semistrp.o \ sttosmst.o \ strip.o \ stripitf.o \ translat.o \ trnsfrm.o all:: libcubpack.a libcubpack.a: $(OBJS) $(CC) -xar -o libcubpack.a $(OBJS) tst: $(MYTESTPROG) $(CC) $(MYTESTPROG) $(libdir) $(libs) -o tst clean: /bin/rm -f *.o run tst core /bin/rm -r -f Templates.DB SHAR_EOF fi # end of overwriting check if test -f 'mkfile.tc' then echo shar: will not over-write existing file "'mkfile.tc'" else cat << \SHAR_EOF > 'mkfile.tc' # Makefile for Cubpack++ installation using Turbo C++ under MS-DOS # File suffixes C = .c # C++ code H = .h # C++ header O = .obj # Compiled code X = .exe # Executable # File delete utility RM = del # Compiler name and flags CC = tcc # Must be 3.0 or later CFLAGS = -v -P -I. -ml # This combination works. Many others don't. # COPTIONS = -DEGTRUSAGE COPTIONS = MYTESTPROG = main$(O) $(C)$(O): $(CC) $(CFLAGS) -c $*$(C) chrono$(O): chrono$(C) $(CC) $(CFLAGS) $(COPTIONS) -c chrono$(C) main$(O): main$(C) $(CC) $(CFLAGS) -c main$(C) OBJS = \ point2D$(O) \ C2$(O) \ C2dv2$(O) \ C2dv4$(O) \ C2interf$(O) \ C2prc$(O) \ C2rule13$(O) \ C2toS2$(O) \ C2togr$(O) \ E2tostrp$(O) \ E2$(O) \ E2adapt$(O) \ E2interf$(O) \ E2sec$(O) \ E2secitf$(O) \ E2secprc$(O) \ S2$(O) \ S2adapt$(O) \ S2interf$(O) \ S2rule13$(O) \ T2$(O) \ T2dv4$(O) \ T2interf$(O) \ T2rule13$(O) \ T2tops$(O) \ atomreg$(O) \ boolean$(O) \ chrono$(O) \ compreg$(O) \ eval_ctr$(O) \ geometry$(O) \ gr$(O) \ gritf$(O) \ gs$(O) \ gsitf$(O) \ gsprc$(O) \ integran$(O) \ integrat$(O) \ invert$(O) \ outS2$(O) \ outS2itf$(O) \ polC2$(O) \ polC2itf$(O) \ polC2prc$(O) \ ps$(O) \ psitf$(O) \ refcount$(O) \ regcoll$(O) \ reginfo$(O) \ region$(O) \ semstitf$(O) \ semistrp$(O) \ sttosmst$(O) \ strip$(O) \ stripitf$(O) \ translat$(O) \ trnsfrm$(O) all: cubpack.lib cubpack.lib: $(OBJS) tlib.rsp tlib cubpack @tlib.rsp tst$(X): $(MYTESTPROG) $(CC) $(CFLAGS) -etst$(X) cubpack.lib $(MYTESTPROG) tst.out: tst$(X) tst > tst.out type tst.out clean: $(RM) *$(O) $(RM) tst$(X) $(RM) cubpack.lib $(RM) tlib.rsp $(RM) *.out tlib.rsp: # TLIB response file copy &&| /C & /P64 & +point2D$(O) & +C2$(O) & +C2dv2$(O) & +C2dv4$(O) & +C2interf$(O) & +C2prc$(O) & +C2rule13$(O) & +C2toS2$(O) & +C2togr$(O) & +E2tostrp$(O) & +E2$(O) & +E2adapt$(O) & +E2interf$(O) & +E2sec$(O) & +E2secitf$(O) & +E2secprc$(O) & +S2$(O) & +S2adapt$(O) & +S2interf$(O) & +S2rule13$(O) & +T2$(O) & +T2dv4$(O) & +T2interf$(O) & +T2rule13$(O) & +T2tops$(O) & +atomreg$(O) & +boolean$(O) & +chrono$(O) & +compreg$(O) & +eval_ctr$(O) & +geometry$(O) & +gr$(O) & +gritf$(O) & +gs$(O) & +gsitf$(O) & +gsprc$(O) & +integran$(O) & +integrat$(O) & +invert$(O) & +outS2$(O) & +outS2itf$(O) & +polC2$(O) & +polC2itf$(O) & +polC2prc$(O) & +ps$(O) & +psitf$(O) & +refcount$(O) & +regcoll$(O) & +reginfo$(O) & +region$(O) & +semstitf$(O) & +semistrp$(O) & +sttosmst$(O) & +strip$(O) & +stripitf$(O) & +translat$(O) & +trnsfrm$(O) | tlib.rsp SHAR_EOF fi # end of overwriting check if test -f 'Mk.decalpha' then echo shar: will not over-write existing file "'Mk.decalpha'" else cat << \SHAR_EOF > 'Mk.decalpha' CUBPACK = /home/ronald/Research/Cubpack++ MYTESTPROG = tmp.o TARGET = tst ######################################################### LIB_DIR=$(CUBPACK)/Src INCLUDE_DIR=$(CUBPACK)/Src CC = cxx CFLAGS = -I$(INCLUDE_DIR) libdir = -L$(LIB_DIR) libs = -lcubpack -lm .SUFFIXES: .o .C .C.o: ; $(CC) $(CFLAGS) -c $*.C $(TARGET): $(MYTESTPROG) $(CC) $(MYTESTPROG) $(libdir) $(libs) -o $(TARGET) clean: /bin/rm -f *.o $(TARGET) core SHAR_EOF fi # end of overwriting check if test -f 'Mk.gn1' then echo shar: will not over-write existing file "'Mk.gn1'" else cat << \SHAR_EOF > 'Mk.gn1' CUBPACK = /home/ronald/Research/Cubpack++ MYTESTPROG = vb5.o TARGET = tst ######################################################### LIB_DIR=$(CUBPACK)/Src INCLUDE_DIR=$(CUBPACK)/Src CC = gcc CFLAGS = -x c++ -I$(INCLUDE_DIR) libdir = -L$(LIB_DIR) libs = -lcubpack -lm -lg++ .c.o: ; $(CC) $(CFLAGS) -c $*.c $(TARGET): $(MYTESTPROG) $(CC) $(MYTESTPROG) $(libdir) $(libs) -o $(TARGET) clean: /bin/rm -f *.o $(TARGET) core SHAR_EOF fi # end of overwriting check if test -f 'Mk.gn2' then echo shar: will not over-write existing file "'Mk.gn2'" else cat << \SHAR_EOF > 'Mk.gn2' CUBPACK = /home/ronald/Research/Cubpack++ MYTESTPROG = vb5.o TARGET = tst ######################################################### LIB_DIR=$(CUBPACK)/Src INCLUDE_DIR=$(CUBPACK)/Src CC = gcc CFLAGS = -x c++ -I$(INCLUDE_DIR) -fno-implicit-templates libdir = -L$(LIB_DIR) libs = -lcubpack -lm -lg++ .c.o: ; $(CC) $(CFLAGS) -c $*.c $(TARGET): $(MYTESTPROG) $(CC) $(MYTESTPROG) $(libdir) $(libs) -o $(TARGET) clean: /bin/rm -f *.o $(TARGET) core SHAR_EOF fi # end of overwriting check if test -f 'Mk.gn3' then echo shar: will not over-write existing file "'Mk.gn3'" else cat << \SHAR_EOF > 'Mk.gn3' CUBPACK = /home/ronald/Research/Cubpack++ MYTESTPROG = vb5.o TARGET = tst ######################################################### LIB_DIR=$(CUBPACK)/Src INCLUDE_DIR=$(CUBPACK)/Src CC = g++ CFLAGS = -I$(INCLUDE_DIR) -fno-implicit-templates libdir = -L$(LIB_DIR) libs = -lcubpack -lm .c.o: ; $(CC) $(CFLAGS) -c $*.c $(TARGET): $(MYTESTPROG) $(CC) $(MYTESTPROG) $(libdir) $(libs) -o $(TARGET) clean: /bin/rm -f *.o $(TARGET) core SHAR_EOF fi # end of overwriting check if test -f 'Mk.solaris' then echo shar: will not over-write existing file "'Mk.solaris'" else cat << \SHAR_EOF > 'Mk.solaris' CUBPACK = /home/ronald/Cubpack++ MYTESTPROG = vb7.o TARGET = tst ######################################################### LIB_DIR=$(CUBPACK)/Src INCLUDE_DIR=$(CUBPACK)/Src CC = CC CFLAGS = -I$(INCLUDE_DIR) libdir = -L$(LIB_DIR) libs = -lcubpack -lm .c.o: ; $(CC) $(CFLAGS) -c $*.c $(TARGET): $(MYTESTPROG) $(CC) $(MYTESTPROG) $(libdir) $(libs) -o $(TARGET) clean: /bin/rm -f *.o $(TARGET) core SHAR_EOF fi # end of overwriting check if test -f 'Mk.tc' then echo shar: will not over-write existing file "'Mk.tc'" else cat << \SHAR_EOF > 'Mk.tc' # Makefile for running Cubpack++ applications using Turbo C++ under MS-DOS CUBPACK = D:\CUBPACK MYTESTPROG = C2toquad.obj quadril.obj quadritf.obj testquad.obj TARGET = tst.exe # You should not need to change anything after this ############################################################################# LIB_DIR=$(CUBPACK)/Src INCLUDE_DIR=$(CUBPACK)/Src # File suffixes C = .c # C++ code H = .h # C++ header O = .obj # Compiled code X = .exe # Executable # File delete utility RM = del # Compiler name and flags. The given compiler options work with # Turbo C++ 3.0. Many others don't. CC = tcc CFLAGS = -v -P -I$(INCLUDE_DIR) -ml libdir = -L$(LIB_DIR) $(C)$(O): $(CC) $(CFLAGS) -c $*$(C) >> compiler.log $(TARGET): $(MYTESTPROG) $(CC) $(CFLAGS) -e$(TARGET) $(libdir) cubpack.lib $(MYTESTPROG) clean: $(RM) *.$(O) $(TARGET) SHAR_EOF fi # end of overwriting check cd .. if test -f 'README.1st' then echo shar: will not over-write existing file "'README.1st'" else cat << \SHAR_EOF > 'README.1st' ////////////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools (ronald@cs.kuleuven.ac.be) // // Dirk Laurie (dirk@calvyn.puk.ac.za) // // Luc Pluym // // // ////////////////////////////////////////////////////////////// // // // Last modification of source code: April 4, 1996 // // Last modification of text : May 23, 1996 // // // ////////////////////////////////////////////////////////////// This version is available by anonymous ftp from ftp.cs.kuleuven.ac.be in /pub/NumAnal-ApplMath/Cubpack/all.tar.gz or /pub/NumAnal-ApplMath/Cubpack/all.zip and also by using a World Wide Web browser such as xmosaic via "http://www.cs.kuleuven.ac.be/~ronald/". If you tell us who you are, we promise to keep you informed about corrections, updates and extensions. If you need assistance, please always mention the two dates above. The next file to read is named INSTALL. Enjoy Cubpack++ ! SHAR_EOF fi # end of overwriting check if test -f 'README.DOS' then echo shar: will not over-write existing file "'README.DOS'" else cat << \SHAR_EOF > 'README.DOS' [ If you are using UNIX, read the file README instead of this one. ] 23 May 1996 The current version of Cubpack++ is Version 0.1h = 1.0 ====================================================== |The complete Cubpack++ distribution is in the file all.zip . | |Switch your ftp to binary mode before you get this file. | |Once you have this file on your machine, you have to execute | unzip -aL all.zip |You now have a directory named CUBPACK . --------------------------------------------------------------------- Start by reading CUBPACK\DOCUMENT\PAPER.* and CUBPACK\DOCUMENT\USERMAN.* and CUBPACK\README.1ST Enjoy Cubpack++! SHAR_EOF fi # end of overwriting check if test -f 'README' then echo shar: will not over-write existing file "'README'" else cat << \SHAR_EOF > 'README' [ If you are using MSDOS, read the file README.DOS instead of this one. ] 23 May 1996 The current version of Cubpack++ is Version 0.1h = 1.0 ====================================================== First check which decompression program is available on your machine. IF you have gunzip and tar THEN |The complete Cubpack++ distribution is in the file all.tar.gz . |(Yes, all patches are already included in this.) | |Switch your ftp to binary mode before you get this file. | |Once you have this file on your machine, you have to execute | gunzip all.tar.gz | tar xf all.tar |You now have a directory named Cubpack++ . --------------------------------------------------------------------- IF you have unzip THEN |The complete Cubpack++ distribution is in the file all.zip . | |Switch your ftp to binary mode before you get this file. | |Once you have this file on your machine, you have to execute | unzip all.zip |You now have a directory named Cubpack . --------------------------------------------------------------------- Start by reading Cubpack{++}/Documentation/paper.* and Cubpack{++}/Documentation/userman.* and Cubpack{++}/README.1st Enjoy Cubpack++! SHAR_EOF fi # end of overwriting check if test -f 'CONTENTS' then echo shar: will not over-write existing file "'CONTENTS'" else cat << \SHAR_EOF > 'CONTENTS' What you find in the current directory -------------------------------------- README.1st : information you should first look at CONTENTS : this file INSTALL : hints to get Cubpack++ up and running Makefile : primitive Makefile for unix systems INSTALL.BAT : installation for MS-DOS systems TESTCASE.BAT : running test cases on MS-DOS systems quick_install : script invoked by `make install' quick_run : script invoked by `make run' Src : directory with all source code Documentation : directory with documentation See Documentation/CONTENTS for details. Drivers : directory with example programs See Documentation/CONTENTS for details. Tools : directory with unix tools to change suffixes and makefiles for several systems See Documentation/CONTENTS for details. Note that on MSDOS systems, the file names are all in the same case and truncated to conform to the 8.3 naming restriction. We have taken care to ensure that no conflicts arise. SHAR_EOF fi # end of overwriting check if test -f 'INSTALL' then echo shar: will not over-write existing file "'INSTALL'" else cat << \SHAR_EOF > 'INSTALL' How to get Cubpack++ up and running? ------------------------------------ 1) Getting the suffix right --------------------------- At this stage of C++, a uniform installation procedure is impossible (one of the reasons is that templates are treated differently by different systems). With this file we provide some help. At the moment all source code is located in one directory, named `Src'. All C++ source files have suffix .c on the distribution. Start with verifying what suffix your compiler expects. - If this is .c, you have nothing to change. (This is the case with e.g. gcc (gnu), xlC (IBM), CC (Sun). With some other compilers (e.g. Turbo C++) this is not the default, but you can specify it with a compiler option.) - If there is absolutely no way to persuade your compiler to accept .c files as C++ code, you have to change not only the suffixes of the files itself but also inside the header files of template classes where .c files are included. For those working on Unix systems, we provided some shell scripts to do this. These are located in the directory `Tools'. Change your current working directory to `Tools', learn from the file `CONTENTS' which script you need and execute if from where you are. Other users can learn from these scripts (e.g. c2C) which files they have to modify. 2) Compile the package ---------------------- Those working on a Unix system can now try to execute make install in the parent directory. Then you can choose between the several compilers we used. Those working on MSDOS can now try to execute INSTALL in the parent directory, and choose a compiler. Now you can go for coffee or lunch, depending on your platform. Some timings to give you an idea: IBM RS6000 AIX 3.2.5, xlC compiler - 181.0u 34.0s 6:09 58% i586 90 MHz Linux1.2.4, g++ compiler (gcc v2.7.2) - 208.79u 36.69s 4:29.19 91% DECalpha OSF/1 V1.3, cxx compiler - 272.08u 62.45s 13:00 42% Sparccenter 1000, SunOS 5.3, g++ compiler ( gcc v2.7.2) - 296.0u 49.0s 7:28 76% DECstation 5000/240 ULTRIX V4.4 , g++ compiler (gcc v2.7.2) - 366.3u 85.6s 15:31 48% Sparccenter 1000, SunOS 5.3, CC compiler - 286.0u 186.0s 28:49 27% i486 33MHz MS-DOS Turbo C++ 3.0 compiler - 7min38s elapsed time HP, HP-UX, g++ compiler (gcc v2.6.0) - 499.2u 41.3s 10:36 84% DECstation 5000/120 ULTRIX V4.4 , g++ compiler (gcc v2.7.2) - 747.8u 148.6s 19:45 75% --------------- | If you use one of the platforms on which we tested Cubpack++ no | problem should occur. If it does, please let us know. | | If you use another Unix platform, writing a specific makefile (by | modifying one of those we provided in the directory Tools) and | extending the script quick_install, should be straitforward. If | you are not feeling that familiar with makefiles, please ask and | we will help. | | If you are not using a Unix platform, we hope you know enough | about your system to get Cubpack++ at work. Probably you can | still find useful information in the makefiles we provided in the | directory Tools. --------------- This was the time consuming part which never has to be done again if you continue to use the same platform. 3) Trying some examples ----------------------- Those working on a Unix system can now try to execute make run in the parent directory. You are requested to select a main program. You can choose among the examples we provided or enter the name of your own main program. When the compilation finished successfully, an executable named `tst' is available in the parent directory. Those working on MSDOS can now try to execute TESTCASE n in the parent directory, where n is a number in the range 1 to 11. (Executing TESTCASE without a number will give a menu.) The corresponding main program will be copied into the CODE directory, compiled, linked and executed. The output will be written to the file CODE\VBn.OUT and may be compared with EXAMPLES\VBn.OUT. --------------- | If you have problems linking with libcubpack.a you can replace | this with $(OBJS) in the lines that define tst. --------------- 4) Using Cubpack++ ------------------ It is of course far from ideal that the user program has to be copied into the directory where the Cubpack++ source code is located. Because templates are treated differently on different platforms, we cannot (yet) provide a good alternative for all platforms. The ideal solution we have in mind can be realized for those using the gnu gcc compiler, the cxx compiler on DECalpha OSF or the CC compiler on Solaris. In the directory Tools we provided a Makefile which a user can copy into his working directory. He should modify the values of the constants CUBPACK, MYTESTPROG and TARGET in the beginning of these files, depending on his situation and wishes. 5) Cleaning up -------------- To remove all object files, archive, the template data base, execute make clean in the parent directory. 6) In case of problems ---------------------- In case of installation problems, please don't wait to contact one of the authors: ----------------------------------------------------------------------------- | Name : Ronald Cools | | Email : ronald.cools@cs.kuleuven.ac.be Katholieke Universiteit Leuven | | ronald@cs.kuleuven.ac.be Department of Computer Science | | Celestijnenlaan 200 A | | Fax : +(32) 16 32 79 96 B-3001 HEVERLEE | | Phone : +(32) 16 32 75 62 BELGIUM | ----------------------------------------------------------------------------- For MSDOS installation, contact: ----------------------------------------------------------------------------- | Name : Dirk Laurie | | Email : dirk@calvyn.puk.ac.za Potchefstroomse Universiteit vir CHO | | na.dlaurie@na-net.ornl.gov Department of Mathematics | | P.O.Box 1174 | | Fax : +(27) 16 807 3614 1900 VANDERBIJLPARK | | Phone : +(27) 16 807 3600 South Africa | ----------------------------------------------------------------------------- SHAR_EOF fi # end of overwriting check if test -f 'INSTALL.BAT' then echo shar: will not over-write existing file "'INSTALL.BAT'" else cat << \SHAR_EOF > 'INSTALL.BAT' @echo off if "%1"=="" goto 1 rem Quick installation on MS-DOS systems using Turbo C++ copy tools\mkfile.%1 code\makefile echo --------------------------------------------------- echo Compilation will now start. This will take a while. echo --------------------------------------------------- chdir code make clean rem The REM'd lines below record start/finish times, and require an rem empty file called CR in the root directory rem type ..\cr|time>time make all rem type ..\cr|time>>time chdir .. goto 2 :1 echo 旼컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴커 echo � Available makefiles under MSDOS � echo � � echo � tc Turbo C++ version 3.00 and later � echo � gnu Gnu C++ version 2.6.0 and earlier � echo � gn2 Gnu C++ version 2.6.2 and later � echo � � echo � Type 'install OPTION' at the DOS prompt to install � echo � the chosen OPTION; e.g. 'install tc'. � echo 읕컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴켸 :2 SHAR_EOF fi # end of overwriting check if test -f 'TESTCASE.BAT' then echo shar: will not over-write existing file "'TESTCASE.BAT'" else cat << \SHAR_EOF > 'TESTCASE.BAT' @echo off if "%1"=="" goto 1 copy examples\vb%1.c code\main.c chdir code del main.obj del tst.* make tst.out rename tst.out vb%1.out chdir .. goto 2 :1 echo 旼컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴� echo � Available examples for Cubpack++ � echo � � echo � vb1 CIRCLE � echo � vb2 PLANE_SECTOR (quadrant) � echo � vb3 PLANE_SECTOR (wedge) � echo � vb4 GENERALIZED_RECTANGLE (parabola) � echo � vb5 GENERALIZED_RECTANGLE (lemniscate) � echo � vb6 SEMI_INFINITE_STRIP � echo � vb7 SEMI_INFINITE_STRIP � echo � vb8 INFINITE_STRIP � echo � vb9 RECTANGLE � echo � vb10 TRIANGLE � echo � vb10 TRIANGLE � echo � vb11 PLANE_SECTOR (wedge) � echo � � echo � Type 'testcase 1' etc. at the DOS � echo � prompt to activate one of these options. � echo 읕컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴� :2 SHAR_EOF fi # end of overwriting check if test -f 'quick_install' then echo shar: will not over-write existing file "'quick_install'" else cat << \SHAR_EOF > 'quick_install' #!/bin/sh # # Constants that must be adjusted # SRC_DIR=Src TOOL_DIR=Tools # ####################################################################### echo "This script should make it easier for you to install Cubpack++." echo "An appropriate makefile is selected based on the brand of unix" echo "you use." echo "Available are" echo " (1) IBM RS6000 - AIX - XlC compiler" echo " (2) DEC alpha - DEC OSF - cxx compiler" echo " (3) Unix gnu - gcc or g++ compiler" echo " (4) Solaris 2.* - CC compiler version 4.*" echo echo "Make your choice. \c" read x case $x in 1) cp $TOOL_DIR/mkfile.rs6000 $SRC_DIR/makefile;; 2) cp $TOOL_DIR/mkfile.decalpha $SRC_DIR/makefile;; 3) OLDVERSION1=2.6.0 OLDVERSION2=2.7.0 gcc -v 2> tmp$$ v=`tail -1 tmp$$ | sed "s/^.*version //" | sed "s/ .*$//"` /bin/rm tmp$$ echo "You are using gcc version $v" a=`expr "$v" ">" "$OLDVERSION2"` if test $a = "1" then cp $TOOL_DIR/mkfile.gn3 $SRC_DIR/makefile else a=`expr "$v" ">" "$OLDVERSION1"` if test $a = "1" then cp $TOOL_DIR/mkfile.gn2 $SRC_DIR/makefile else cp $TOOL_DIR/mkfile.gn1 $SRC_DIR/makefile fi fi;; 4) cp $TOOL_DIR/mkfile.solaris $SRC_DIR/makefile;; *) echo "This was an illegal choice!" exit 1 ;; esac ####################################################################### if test -f .lastinstall then y=`cat .lastinstall` /bin/rm -f .lastinstall if test $x != $y then echo "Removing the remains of a previous compilation which might interfere ..." cd $SRC_DIR make clean cd .. echo "... Done." fi fi echo $x > .lastinstall ####################################################################### echo "Do you want to edit the makefile first? [n/y] \c" read a case $a in [Yy]*) vi $SRC_DIR/makefile ;; esac ####################################################################### # echo echo "Compilation will now start. This will take a while." echo cd $SRC_DIR make all SHAR_EOF chmod +x 'quick_install' fi # end of overwriting check if test -f 'quick_run' then echo shar: will not over-write existing file "'quick_run'" else cat << \SHAR_EOF > 'quick_run' #!/bin/sh # # Constants that must be adjusted # SRC_DIR=Src # ####################################################################### if test -f .lastinstall then y=`cat .lastinstall` case $y in 1|3|4) SUFFIX=c;; 2) SUFFIX=C;; *) echo "Unknown compiler used for installation." exit 1;; esac else echo "make install must be executed first." exit 1 fi ####################################################################### echo "Choose a main program from the following list." echo cat Drivers/EXAMPLES echo echo "Make your choice. \c" read a case $a in 0) echo "Give the name of your example file." echo "(Full pathname from current position.)" read a cp $a $SRC_DIR/main.$SUFFIX ;; [1-9]|1[0-8]) cp Drivers/vb${a}.c $SRC_DIR/main.$SUFFIX ;; *) echo "That was an illegal choice!" exit 1 ;; esac ####################################################################### # echo echo "Compilation will now start." echo cd $SRC_DIR make tst SHAR_EOF chmod +x 'quick_run' fi # end of overwriting check if test -f 'Makefile' then echo shar: will not over-write existing file "'Makefile'" else cat << \SHAR_EOF > 'Makefile' SRC_DIR = Src install: ./quick_install run: ./quick_run mv $(SRC_DIR)/tst . echo 'executable "tst" is available.' clean: cd $(SRC_DIR) ; /bin/rm -f *.o tst core makefile cd $(SRC_DIR) ; /bin/rm -f main.c libcubpack.a cd $(SRC_DIR) ; /bin/rm -r -f ptrepository tempinc Templates.DB SHAR_EOF fi # end of overwriting check if test ! -d 'Src' then mkdir 'Src' fi cd 'Src' if test -f 'C2.c' then echo shar: will not over-write existing file "'C2.c'" else cat << \SHAR_EOF > 'C2.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ///////////////////////////////////////////////// //File C2.c // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) //////////////////////////////////////////////// #include #include ////////////////////////////////////////////// Parallelogram::Parallelogram (const Point& p1,const Point& p2,const Point& p3) : Geometry(2),Vertices(3),TheVolumeKnown(False) { Vertices[0] = p1; Vertices[1] = p2; Vertices[2] = p3; } ////////////////////////////////////////////// void Parallelogram::ComputeVolume() { TheVolume = fabs( ((Vertices[0][0]*Vertices[1][1]-Vertices[1][0]*Vertices[0][1]) +(Vertices[1][0]*Vertices[2][1]-Vertices[2][0]*Vertices[1][1]) +(Vertices[2][0]*Vertices[0][1]-Vertices[0][0]*Vertices[2][1])) ); TheVolumeKnown = True; } /////////////////////////////////////////////// const Point& Parallelogram::Vertex(int i) const { return Vertices[i]; } ////////////////////////////////////////////////// real Parallelogram::Volume() const { Parallelogram* P = (Parallelogram*) this; if (!TheVolumeKnown) { P->ComputeVolume(); }; return P->TheVolume; } /////////////////////////////////////////////// void Parallelogram::Volume(real v) { TheVolume = v; TheVolumeKnown = True; } //////////////////////////////////////////////// //Processor* //Parallelogram::DefaultProcessor() //const //{ //return new SimpleAdaptive //(new Parallelogram_Rule13, //new Parallelogram_Divide4); //} /////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'C2dv2.c' then echo shar: will not over-write existing file "'C2dv2.c'" else cat << \SHAR_EOF > 'C2dv2.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ////////////////////////////////////////////////// //File C2dv2.cpp // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) // 25 Jan 1996 V0.1f(long lines split) ///////////////////////////////////////////////// #include #include #include #include ////////////////////////////////////////////////// void Parallelogram_Divide2::Apply (const Parallelogram& t, Stack& Offspring, const Vector& DiffOrder) { Point m12 = (t.Vertex(0)+t.Vertex(DiffOrder[0]))/2; Point m34 = t.Vertex(DiffOrder[1])+(t.Vertex(DiffOrder[0])-t.Vertex(0))/2; Parallelogram* t1 = new Parallelogram(t.Vertex(0),m12,t.Vertex(DiffOrder[1])); Parallelogram* t2 = new Parallelogram(m12,t.Vertex(DiffOrder[0]),m34); Offspring.Push(t1); Offspring.Push(t2); } ////////////////////////////////////////////////// Parallelogram_Divide2::Parallelogram_Divide2() :SameShapeDivisor() { } ////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'C2dv4.c' then echo shar: will not over-write existing file "'C2dv4.c'" else cat << \SHAR_EOF > 'C2dv4.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ////////////////////////////////////////////////// //File C2dv4.cpp // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) // 25 Jan 1996 V0.1f(long lines split) ///////////////////////////////////////////////// #include #include #include ////////////////////////////////////////////////// void Parallelogram_Divide4::Apply (const Parallelogram& t, Stack& Offspring, const Vector& DiffOrder) { Point m12 = (t.Vertex(0)+t.Vertex(1))/2; Point m13 = (t.Vertex(0)+t.Vertex(2))/2; Point m23 = (t.Vertex(2)+t.Vertex(1))/2; Point m24 = t.Vertex(1)+(t.Vertex(2)-t.Vertex(0))/2; Point m34 = t.Vertex(2)+(t.Vertex(1)-t.Vertex(0))/2; Parallelogram* t1 = new Parallelogram(t.Vertex(0),m12,m13); Parallelogram* t2 = new Parallelogram(m12,t.Vertex(1),m23); Parallelogram* t3 = new Parallelogram(m13,m23,t.Vertex(2)); Parallelogram* t4 = new Parallelogram(m23,m24,m34); Offspring.Push(t1); Offspring.Push(t2); Offspring.Push(t3); Offspring.Push(t4); } ////////////////////////////////////////////////// Parallelogram_Divide4::Parallelogram_Divide4() :SameShapeDivisor() { } ////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'C2interf.c' then echo shar: will not over-write existing file "'C2interf.c'" else cat << \SHAR_EOF > 'C2interf.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// //////////////////////////////////////////////// //File C2interf.c // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) // 25 Jan 1996 V0.1f(typedef introduced) //////////////////////////////////////////////// #include #include #include #include #include #include #include /////////////////////////////////////////////// typedef Rule RuleParallelogram; typedef SameShapeDivisor SameShapeDivisorParallelogram; PARALLELOGRAM::PARALLELOGRAM(const Point& p1, const Point& p2, const Point& p3) :USERINTERFACE() {Pointer R13(new Parallelogram_Rule13); Pointer D4(new Parallelogram_Divide4); Processor *SAP(new SimpleAdaptive(R13,D4)); StoreAtomic(new Parallelogram(p1,p2,p3),SAP); } /////////////////////////////////////////////// RECTANGLE::RECTANGLE(const Point& p1, const Point& p2, const Point& p3) :USERINTERFACE() { Error( fabs((p2-p1)*(p3-p1)) > 100*REAL_EPSILON, "Sides of RECTANGLE not orthogonal"); StoreAtomic(new Parallelogram(p1,p2,p3), //new SimpleAdaptive( //new Parallelogram_Rule13, //new Parallelogram_Divide4)); new Parallelogram_Processor); } /////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'C2prc.c' then echo shar: will not over-write existing file "'C2prc.c'" else cat << \SHAR_EOF > 'C2prc.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// //File C2prc.c // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) // 25 Jan 1996 V0.1f(typedef introduced and long lines split) //////////////////////////////////////////////////////////// #include #include #include #include #include #include //////////////////////////////////////////////////////// Pointer < Parallelogram_Processor::RuleParallelogram > Parallelogram_Processor::TheRule = new Parallelogram_Rule13; Pointer < Parallelogram_Processor::SameShapeDivisorParallelogram > Parallelogram_Processor::TheDivisor2 = new Parallelogram_Divide2; Pointer < Parallelogram_Processor::SameShapeDivisorParallelogram > Parallelogram_Processor::TheDivisor4 = new Parallelogram_Divide4; //////////////////////////////////////////////////////// Parallelogram_Processor::Parallelogram_Processor() :TimesCalled(0), Diffs(2) { } ///////////////////////////////////////////////////////////// Parallelogram_Processor* Parallelogram_Processor::Descendant() const { Parallelogram_Processor* r = new Parallelogram_Processor; return r; } //////////////////////////////////////////////// void Parallelogram_Processor::Process( Stack& Offspring) { TimesCalled ++; if (TimesCalled == 1) { TheRule->ApplyWithDiffs(LocalIntegrand(),Geometry(),Integral(), AbsoluteError(),Diffs); Offspring.MakeEmpty(); return; }; if(TimesCalled == 2) { real NewVolume = Geometry().Volume()/2; Stack Parts; Vector DiffOrder(Diffs.Size()); const real difffac = 1.81818, difftreshold = 0.00004; if (max(Diffs[0],Diffs[1]) < difftreshold) { TheDivisor4->Apply(Geometry(),Parts,DiffOrder); NewVolume /=2; } else if (Diffs[0]>difffac*Diffs[1]) { DiffOrder[0] = 1 ; DiffOrder[1] = 2; TheDivisor2->Apply (Geometry(),Parts,DiffOrder); } else if (Diffs[1]>difffac*Diffs[0]) { DiffOrder[0] = 2; DiffOrder[1] =1; TheDivisor2->Apply (Geometry(),Parts,DiffOrder); } else { TheDivisor4->Apply(Geometry(),Parts,DiffOrder); NewVolume /=2; }; unsigned int N = Parts.Size(); for (unsigned int i =0;iVolume(NewVolume); Processor* p = Descendant(); Atomic* a = new Atomic(g,p); a->LocalIntegrand(&LocalIntegrand()); Offspring.Push(a); }; return; }; Error(TimesCalled > 2, "Parallelogram_Processor : more than two calls of Process()"); } /////////////////////////////////////////////// Processor* Parallelogram_Processor::NewCopy() const { return new Parallelogram_Processor(*this); } /////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'C2rule13.c' then echo shar: will not over-write existing file "'C2rule13.c'" else cat << \SHAR_EOF > 'C2rule13.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// // //File C2rule13.cpp // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) // 25 Jan 1996 V0.1f(long lines split) ///////////////////////////////////////////////////////// // // A cubature formula of degree 13 with 37 points and associated // null rules. // The cubature formula, from Rabinowitz & Richter, was recomputed // to obtain higher accuracy. // ///////////////////////////////////////////////////////// #include static int K[]={1,2,3,2}; const int NumberOfPoints=K[0]+4*K[1]+4*K[2]+8*K[3]; const int Orbits=K[0]+K[1]+K[2]+K[3]; static real Type1[]={ 0.9909890363004326469792722978603e+00, 0.6283940712305315063814483471116e+00}; static real l1sdl2s = Type1[1]*Type1[1]/Type1[0]/Type1[0]; static real Type2[]={ 0.9194861553393073086142137772149e+00, 0.6973201917871173078084506730937e+00, 0.3805687186904854497424188074662e+00}; static real Type3[2][2]={{0.9708504361720225062147290554088e+00, 0.6390348393207252159077623446225e+00}, {0.8623637916722781475018696425693e+00, 0.3162277660168700033875075593701e+00}}; static real Weight[8][8]={{0.2995235559387052215463143056692e+00, 0.3311006686692356205977471655046e-01, 0.1802214941550624038355347399683e+00, 0.3916727896035153300761243260674e-01, 0.1387748348777288706306435595057e+00, 0.2268881207335707037147066705814e+00, 0.3657395765508995601240002438981e-01, 0.1169047000557533546701746277951e+00}, {7.610781847149629154716409791983e-2, 1.486101247399760261471935168346e-1, -2.077685631717747007172983323970e-1, 6.850758313011924198538315395405e-2, 2.024205813317813585572881715385e-1, 1.108627473745508429879249169864e-1, -1.187411393304862640859204217487e-1, -5.208857468077715683772080394959e-2}, {4.016494861405949747097510013162e-2, -1.093132962444079541048635452881e-1, -2.270251673633777452624380129694e-1, 1.231674163356097016086203579325e-2, -1.420402526499201540699111172200e-1, 1.189080551229557928776504129312e-1, -4.482039658150474743804189300793e-3, 1.730383808319875827592824151609e-1}, {-5.643905795781771973971259866415e-1, 2.878418073676293225652331648545e-2, 1.159354231997583294689565314470e-1, 1.376081498690624477894043101438e-1, -7.909780225340130915490382973570e-2, 1.174335441429478112778176601234e-1, -1.107251942334134124782600707843e-1, 2.094226883312045633400182488252e-2}, {-2.269001713589584730602581694579e-1, 2.976190892690301120078774620049e-2, -7.440193483272787588251423144751e-2, -1.224665989043784131260454301280e-1, -4.857910454732976198562745578156e-2, 2.228157325962656425537280474671e-1, 1.459764751457503859063666414952e-1, -1.211789553452468781539987084682e-1}, {-3.326760468009974589269134283992e-1, 1.796655319904795478676993902115e-1, -4.389976396805911868560791966472e-2, -2.295841771339316497310760908889e-1, 6.182618387692816082856552878852e-2, -1.202703885325137746461829140891e-1, 5.109536580363550180208564374234e-3, 1.126062761533095493689566169969e-1}, {2.290638530086106512999345512401e-1, 2.702070398116919449911037051753e-1, -9.078047988731123605988441792069e-3, 4.618480310858703283999169489655e-2, -2.598231009547631799096616255056e-1, -2.518433931146441037986247681820e-2, -1.257796993152456033984707367389e-2, -2.720818902721190304043617320910e-2}, {2.746908885094872977794932213372e-1, -1.149427039769738298032807785523e-2, 1.596178537820019535731955591283e-1, -2.180626972663360443142752377527e-1, -8.711748038292630173597899697063e-3, 1.902786182960269617633915869710e-1, -1.189840649092108827784089292890e-1, 2.883382565767354162177931122471e-2}}; // // The error estimator constants // static const real crival=0.4 , facmed=8.0; static const real facopt = facmed/(crival*crival); #include #include #include #include #define sqr(x) ((x)*(x)) void Parallelogram_Rule13::ApplyWithDiffs(Integrand& F,Parallelogram& R, real& TheResult,real& TheError, Vector& DiffOrder) { int i,j,p,Type,nr,number; real Null[7]; real Tres = -1.0; real noise,r,r1,r2,r3,Deg5,Deg3,Deg1,Deg7,z1,z2,sumval; real D1,D2; real F0; real FE[2][2][2]; // [l1,l2][vertex1,vertex2][plus,minus] Point x[8]; if (Tres <= 0) { Tres=50*REAL_EPSILON; } TheResult = 0.0; for (i=0;i<7;i++){Null[i]=0;} p=0; for (Type=0;Type <=3;Type++){ for (nr=0;nr= 1) { TheError =facmed*Deg7; } else { if (r>=crival) { TheError = facmed*r*Deg7; } else { TheError = facopt *r*r*r*Deg7; } } TheError = max(noise,TheError); } TheError *= R.Volume()/4; TheResult *=R.Volume()/4; } ////////////////////////////////////////////////// Parallelogram_Rule13::Parallelogram_Rule13() :Rule() { } ////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'C2toS2.c' then echo shar: will not over-write existing file "'C2toS2.c'" else cat << \SHAR_EOF > 'C2toS2.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ////////////////////////////////////////// //File C2toS2.cpp // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) ////////////////////////////////////////// #include #include #include ////////////////////////////////////////// void PolarToRectangular::Transform(real& w, Point& p) { Point P(p.R()*cos(p.Theta()),p.R()*sin(p.Theta())); w *= p.R(); p = P; } /////////////////////////////////////////// PolarToRectangular::PolarToRectangular() : Transformation() { } /////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'C2togr.c' then echo shar: will not over-write existing file "'C2togr.c'" else cat << \SHAR_EOF > 'C2togr.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ///////////////////////////////////////////////// //File gr.cpp // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) //////////////////////////////////////////////// #include #include //////////////////////////////////////////// void C2toGR::Transform(real& w, Point& p) { GeneralizedRectangle& G = *GR_ptr; Point M = G.B()-G.A(); Point C = G.A() + p.X()*M; real ml = M.Length(); real dist = G.Boundary(C), ratio=dist/ml; Point P(-ratio*M.Y(),ratio*M.X()); w *= dist*ml; p = C+p.Y()*P; } /////////////////////////////////////////////////// C2toGR::C2toGR(GeneralizedRectangle* g) : Transformation(),GR_ptr(g) { } ////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'E2.c' then echo shar: will not over-write existing file "'E2.c'" else cat << \SHAR_EOF > 'E2.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ///////////////////////////////////////////////////////// //File E2.cpp // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) // 25 Nov 1994 V0.1d(re-ordering of member initializers) ///////////////////////////////////////////////////////// #include #include #include #include ///////////////////////////////////////////////////////// Plane::Plane () : Geometry(2), xscale(1.0), yscale(1.0), TheCenter(0,0) { } ///////////////////////////////////////////////////////// Plane::Plane (const Point& center) : Geometry(2), xscale(1.0), yscale(1.0), TheCenter(center) { } ///////////////////////////////////////////////////////// Plane::Plane(const Point& center, real x, real y) : Geometry(2), xscale(x), yscale(y), TheCenter(center) { } ///////////////////////////////////////////////////////// real Plane::ScaleX() const { return xscale; } ///////////////////////////////////////////////////////// real Plane::ScaleY() const { return yscale; } ///////////////////////////////////////////////////////// //Processor* //Plane::DefaultProcessor() //const //{ //return new PlaneAdaptive; //} ///////////////////////////////////////////////////////// const Point& Plane::Center() const { return TheCenter; } ///////////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'E2adapt.c' then echo shar: will not over-write existing file "'E2adapt.c'" else cat << \SHAR_EOF > 'E2adapt.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ///////////////////////////////////////////////// //File E2.cpp // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) //////////////////////////////////////////////// #include #include #include #include #include #include #include #include /////////////////////////////////////////////// #define sqr(x) ((x)*(x)) static int K[]={0,3,2,2}; const int NumberOfPoints=K[0]+4*K[1]+4*K[2]+8*K[3]; const int Orbits=K[0]+K[1]+K[2]+K[3]; static real Type1[]= {0.5460844152585403666984417e+00 , 0.1376486679696350642323050e+01 , 0.2358932261980680967828022e+01}; static real Type2[]= {0.8462499884480199068814807e+00 , 0.1985531339917884990890422e+01}; static real Type3[2][2]= {{0.1779246738905319006374908e+01 , 0.9318240227761998500475049e+00}, {0.3030436065368579151664045e+01 , 0.9962534261727632104611162e+00}}; static real OrigWeight[7]= {0.1418036987406628465242303e+00, 0.3807850175899473939094384e-01, 0.1290856883945041937607418e-02, 0.5772571826901623901546648e-01, 0.1620590503146543183518239e-03, 0.5445472363593369147167489e-02, 0.2411028493987025953260209e-04}; static int d[]={0,3,1,5,2,4,6}; static real Weight[7][7]= {{0.6002713316445080080142126, 0.7955990316249109629216728, 1.0584888443837658614605899, 0.7595383160846588054989698, 1.3523576929325600784468914, 0.9663531601396869556259704, 1.9895758469520118113052969}, {-0.1714120818975139, -1.4261615934658296, 20.9491066359021694, 1.7029699666146943, 57.2013795938251143, -4.3332227154742827, -175.0624620496226842}, {-0.2883800695213331, 2.5450105556734605, 2.8699527274413806, 0.304671504643251, 138.1913876850488233, -8.5573830821954439, -18.7930973279112966}, {-0.250456836419066, 0.3273926337592694, -83.9551965137636766, 0.4609375440071548, -365.1974538111883157, 10.4127651350263225, 1537.3605596124529503}, {-0.2595143739988262, 0.1255865390202144, 85.4652862386879934, -0.2489288730337109, -712.5513360902265439, 6.3018929219370941, -645.6559424743654782}, {-0.24820543230912, -0.2143965818640334, 8.5337065366470588, -0.2238393816932916, 867.6556565044116888, 7.6084037588845301, -4712.804729866979117}, {0.3052968227487744, 1.0011209710514362, -6.3357177267870899, 0.7227067892680122, -358.0249035190603089, 2.7469614313645364, -6407.2910382311609777}}; // Radii of rings associated with weights of published formula static real RingRadius[]= {0.9151577658134376, 1.2649388389767648, 1.7333214530567136, 2.2616892861784437, 2.6609731353556634, 2.9246248487342039, 6.0979478618219725}; //This is `infinity'. Never used actually // // The error estimator constants // const real crival=0.5 , facmed=5.0; const real facopt = facmed/(crival*crival); ////////////////////////////////////////////////////// void PlaneAdaptive::Rule(Integrand& F,Plane& R,real& TheResult,real& TheError,real& HalfValueRadius) { int i,j,p,Type,nr,number; real Null[6],OrbitContrib[7],below,above,fraction; real Tres = 50*REAL_EPSILON; real noise,r,r1,r2,r3,Deg5,Deg3,Deg7,Deg1,sumval; real A=1.0/R.ScaleX() , B=1.0/R.ScaleY(); Point x[8]; TheResult = 0.0; for (i=0;i<6;i++){Null[i]=0;} p=0; for (Type=0;Type <=3;Type++){ for (nr=0;nr= 1) { TheError =facmed*Deg7; } else { if (r>=crival) { TheError = facmed*r*Deg7; } else { TheError = facopt *r*r*r*Deg7; } } TheError = max(noise,TheError); } TheResult *= (A*B); TheError *= (A*B); // cout <<"=> Error = "<& Offspring) { TimesCalled++; if (TimesCalled==1) { Rule(LocalIntegrand(),Geometry(),Integral(),AbsoluteError(),HalfValueRadius); } else { const Point& C=Geometry().Center(); Offspring.Push(( AtomicRegion*) CIRCLE(C,HalfValueRadius)); Offspring.Push(( AtomicRegion*) OUT_CIRCLE(C,HalfValueRadius)); Offspring.IteratorReset(); while (!Offspring.IteratorAtEnd()) { Offspring.IteratorNext()->LocalIntegrand(&LocalIntegrand()); }; } } /////////////////////////////////////////////////// PlaneAdaptive::PlaneAdaptive() :Processor(),TimesCalled(0),HalfValueRadius(0) { } //////////////////////////////////////////////////// Processor* PlaneAdaptive::NewCopy() const { return new PlaneAdaptive(*this); } ////////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'E2interf.c' then echo shar: will not over-write existing file "'E2interf.c'" else cat << \SHAR_EOF > 'E2interf.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// /////////// //File E2interf.c // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) //////////////// #include #include /////////////////////// PLANE::PLANE() { StoreAtomic(new Plane,new PlaneAdaptive); } ////////////////////////////////////////////// PLANE::PLANE(const Point& Center) { StoreAtomic(new Plane(Center),new PlaneAdaptive); } ////////////////////////////////////////////// PLANE::PLANE(const Point& Center, real ScaleX, real ScaleY) { StoreAtomic(new Plane(Center,ScaleX, ScaleY),new PlaneAdaptive); } /////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'E2sec.c' then echo shar: will not over-write existing file "'E2sec.c'" else cat << \SHAR_EOF > 'E2sec.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ////////////////////////////////////////// //File E2sec.c // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) // 14 Feb 1995 V0.1e(bug fix in transformation) ////////////////////////////////////////// #include #include #include ////////////////////////////////////////// PlaneSector::PlaneSector( const Point& A, const Point& B, const Point& C) :Geometry(2) { if(C==B) { TheCenter=A; } else { real dp = (C-B)*(A-B); Error(dp==0, "Error: trying to construct a plane sector with centre at infinity"); TheCenter = (A+(A-B)*((C-B)*((C+B)/2-A))/dp); }; TheInnerRadius=(TheCenter-A).Length(); Point BB(B-TheCenter), CC(C-TheCenter); TheSmallAngle=atan2(BB.Y(),BB.X()); TheBigAngle=atan2(CC.Y(),CC.X()); if (TheBigAngle<=TheSmallAngle) TheBigAngle += 2*M_PI; } /////////////////////////////////////////////////// PlaneSector::PlaneSector(const Point& O,real r, real theta1, real theta2) : Geometry(2) { TheCenter = O; TheInnerRadius = r; TheSmallAngle = theta1; TheBigAngle = theta2; if (TheSmallAngle == TheBigAngle) TheBigAngle += 2*M_PI; } ///////////////////////////////////////////////////// real PlaneSector::InnerRadius() const { return TheInnerRadius; } //////////////////////////////////////////////////// real PlaneSector::SmallAngle() const { return TheSmallAngle; } ///////////////////////////////////////////////////// real PlaneSector::BigAngle() const { return TheBigAngle; } //////////////////////////////////////////////////// const Point& PlaneSector::Center() const { return TheCenter; } //////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'E2secitf.c' then echo shar: will not over-write existing file "'E2secitf.c'" else cat << \SHAR_EOF > 'E2secitf.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ////////////////////////////////////////////// //File E2secitf.c // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) ////////////////////////////////////////////// #include #include ////////////////////////////////////////////// PLANE_SECTOR::PLANE_SECTOR (const Point& A, const Point& B,const Point& C) { StoreAtomic(new PlaneSector(A,B,C),new PlaneSector_Processor); } ////////////////////////////////////////////// PLANE_SECTOR::PLANE_SECTOR (const Point& O, real r , real theta1, real theta2) { StoreAtomic(new PlaneSector(O,r,theta1,theta2), new PlaneSector_Processor); } /////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'E2secprc.c' then echo shar: will not over-write existing file "'E2secprc.c'" else cat << \SHAR_EOF > 'E2secprc.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// //////////////////////////////////////////////////// //File E2secprc.c // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) ///////////////////////////////////////////////////////// #include #include #include #include #include #include #include //////////////////////////////////////////////////////// void PlaneSector_Processor::Process( Stack& Offspring) { PlaneSector& G = Geometry(); Point P1(G.InnerRadius(),G.BigAngle()), P2(G.InnerRadius(),G.SmallAngle()); AtomicRegion* A ; A= (AtomicRegion*) SEMI_INFINITE_STRIP(P1,P2); Integrand I1(LocalIntegrand(),new Translation(G.Center())); A->LocalIntegrand(new Integrand(I1, new PolarToRectangular)); Offspring.Push(A); } ///////////////////////////////////////////////// PlaneSector_Processor::PlaneSector_Processor() :Processor() { } ///////////////////////////////////////////////// Processor* PlaneSector_Processor::NewCopy() const { return new PlaneSector_Processor(*this); } //////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'E2tostrp.c' then echo shar: will not over-write existing file "'E2tostrp.c'" else cat << \SHAR_EOF > 'E2tostrp.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ///////////////////////////////////////////////// //File E2tostrp.c // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) //////////////////////////////////////////////// #include #include #include //////////////////////////////////////////// void E2toIS::Transform(real& w, Point& p) { InfiniteStrip& s = *IS_ptr; Point D = s.B()- s.A(); Point C(-D.Y(),D.X()); C = C/C.Length(); w *= D.Length()*.5*(1-tanh(p.X())*tanh(p.X())); p = s.A() + p.Y()*C + (.5*(1+tanh(p.X())))*D; } /////////////////////////////////////////////////// E2toIS::E2toIS(InfiniteStrip* g) : Transformation(),IS_ptr(g) { } ////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'S2.c' then echo shar: will not over-write existing file "'S2.c'" else cat << \SHAR_EOF > 'S2.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ////////////////////////////////////////// //File S2.cpp // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) ////////////////////////////////////////// #include #include #include #include #include #include #include #include #include //////////////////////////////////////////// Circle::Circle(const Point& Center, real Radius) : Geometry(2), TheCenter(Center), TheRadius(fabs(Radius)) {} /////////////////////////////////////////////// Circle::Circle(const Point& Center, const Point& Boundary) : Geometry(2), TheCenter(Center), TheRadius(sqrt((Center.X()-Boundary.X())*(Center.X()-Boundary.X()) +(Center.Y()-Boundary.Y())*(Center.Y()-Boundary.Y()))) {} /////////////////////////////////////////////// real Circle::Volume() const { return M_PI*TheRadius*TheRadius; } ////////////////////////////////////////////// const Point& Circle::Center() const { return TheCenter; } ///////////////////////////////////////////// real Circle::Radius() const { return TheRadius; } ////////////////////////////////////////////// //Processor* //Circle::DefaultProcessor() //const //{ //return new CircleAdaptive(new Circle_Rule13); //} ///////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'S2adapt.c' then echo shar: will not over-write existing file "'S2adapt.c'" else cat << \SHAR_EOF > 'S2adapt.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ////////////////////////////////////////// //File S2adapt.c // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) // 30 Jan 1996 V0.1g(no longer compare signed and unsigned) ////////////////////////////////////////// #include #include #include #include #include //////////////////////////////////////////////// void CircleAdaptive::Process(Stack& Offspring) { Circle& C = Geometry(); static unsigned int mxgen=2; //if (mxgen == 0) //{ //cerr<<"How many times may I apply the circle rule? "; //cin>>mxgen; //}; TimesCalled++; if(TimesCalled ==1) { if (GenerationNumber < mxgen) { TheRule->Apply(LocalIntegrand(),C,Integral(),AbsoluteError()); } else { Point rv (C.Radius(),0); Point bp = C.Center()+rv; AtomicRegion* A ; A= (AtomicRegion*) POLAR_RECTANGLE(C.Center(),bp,bp); A->LocalIntegrand(&(LocalIntegrand())); Offspring.Push(A); }; } else { const real CircleRatio = 0.44721359549995793928; Point m1(C.Radius()*CircleRatio ,0.); Point m2(C.Radius() , 0.); Point m3(0.,C.Radius()); Point m4(0.,C.Radius()*CircleRatio); Point origin=C.Center(); AtomicRegion* t1=(AtomicRegion*) POLAR_RECTANGLE(origin+m1,origin+m2,origin+m3); AtomicRegion* t2=(AtomicRegion*) POLAR_RECTANGLE(origin+m4,origin+m3,origin-m2); AtomicRegion* t3=(AtomicRegion*) POLAR_RECTANGLE(origin-m1,origin-m2,origin-m3); AtomicRegion* t4=(AtomicRegion*) POLAR_RECTANGLE(origin-m4,origin-m3,origin+m2); Offspring.Push(t1); Offspring.Push(t2); Offspring.Push(t3); Offspring.Push(t4); Offspring.IteratorReset(); while (!Offspring.IteratorAtEnd()) { Offspring.IteratorNext()->LocalIntegrand(&LocalIntegrand()); }; if (CircleRatio != 0.0) { Circle* tmp = new Circle(origin,C.Radius()*CircleRatio); Atomic* a =new Atomic(tmp,Descendant()); a->LocalIntegrand(&LocalIntegrand()); Offspring.Push(a); }; }; } ////////////////////////////////////////////////// CircleAdaptive::CircleAdaptive(Rule * R) :Processor(),TheRule(R),TimesCalled(0),GenerationNumber(0) { } ///////////////////////////////////////////////////// CircleAdaptive* CircleAdaptive::Descendant() const { CircleAdaptive * CA = new CircleAdaptive (&(*(TheRule))); CA->GenerationNumber = GenerationNumber+1; return CA; } ////////////////////////////////////////////////////// Processor* CircleAdaptive::NewCopy() const { return new CircleAdaptive(*this); } ////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'S2interf.c' then echo shar: will not over-write existing file "'S2interf.c'" else cat << \SHAR_EOF > 'S2interf.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// ///////////////////////////////////////////////////////// //File S2interf.c // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) ///////////////////////////////////////////////////////// #include #include #include ///////////////////////////////////////////////////////// CIRCLE::CIRCLE(const Point& c,const Point& b) { StoreAtomic(new Circle(c,b),new CircleAdaptive(new Circle_Rule13)); } ///////////////////////////////////////////////////////// CIRCLE::CIRCLE(const Point& c, real Radius) { StoreAtomic(new Circle(c,Radius),new CircleAdaptive(new Circle_Rule13)); } ///////////////////////////////////////////////////////// SHAR_EOF fi # end of overwriting check if test -f 'S2rule13.c' then echo shar: will not over-write existing file "'S2rule13.c'" else cat << \SHAR_EOF > 'S2rule13.c' ///////////////////////////////////////////////////////// // // // Cubpack++ // // // // A Package For Automatic Cubature // // // // Authors : Ronald Cools // // Dirk Laurie // // Luc Pluym // // // ///////////////////////////////////////////////////////// //File S2rule13.cpp // History: // (date) (version) // 19 Aug 1994 V0.1 (first limited distribution) // // A cubature formula of degree 13 with 36 points and associated // null rules (Cools & Haegemans). // #include static real facmed = 1.5; static real crival = 0.3; static int K[]={0,3,2,2}; const int NumberOfPoints=K[0]+4*K[1]+4*K[2]+8*K[3]; const int Orbits=K[0]+K[1]+K[2]+K[3]; static real Type1[]= {0.283402832348840340260055531925e+00, 0.649007230578023857142954047869e+00, 0.920575474036255978410400945265e+00}; static real Type2[]= {0.677355106028069632559870394597e+00, 0.412672216763278882147845015245e+00}; static real Type3[2][2]= {{0.335767600829543545504321849367e+00, 0.938694583835106683089433405007e+00}, {0.752042776803937943743749637665e+00, 0.379717011170079696840056017956e+00}}; static real Weight[7][7]= {{0.497326951852944153115255529183e-01, 0.439732090591084987635001867639e-01, 0.221057969395291114959209696751e-01, 0.138071954801617944635783199052e-01, 0.478191639416850908138653549756e-01, 0.579302409092634509107760741653e-02, 0.304879456061841994847272004644e-01}, // Nullrule exact for the first 6 basispolynomials, {-0.8073879159066197e-02, -0.2498315518214109e-01, 0.5349170539028515e-01, 0.3595547716318644e-01, 0.4320518775944201e-01, -0.1967607360410405e-01, -0.3012159438174910e-01}, // Nullrule exact for the first 5 basispolynomials, {0.2175100120594341e-01, -0.7117820786684134e-01, 0.1801191272610644e-01, -0.1775642346529020e-01, -0.2754333027079005e-02, -0.1170796597217573e-01, 0.3767099118575607e-01}, // Nullrule exact for the first 4 basispolynomials, {0.2413156434309756e-01, 0.2327141313816948e-01, 0.3438009174891747e-01, 0.3084063431284037e-01, -0.6858211678704531e-01, -0.2561913343310581e-01, 0.3598340055116032e-02}, // Nullrule exact for the first 3 basispolynomials, {-0.4853627844062673e-01, 0.3645578330411178e-01, 0.4003943562394422e-01, -0.4494496660283691e-01, 0.1266357696444992e-01, -0.2028566592818438e-01, 0.2244689050366323e-01}, // Nullrule exact for the first 2 basispolynomials, {-0.2979266039117054e-01, 0.3977474716810517e-02, -0.4396076050818155e-01, 0.5463380585542593e-01, 0.1777067704103634e-01, -0.2876745680547196e-01, 0.2745318844851159e-01}, // Nullrule exact for the first 1 basispolynomials, {0.6246499865107039e-01, 0.2624836264288300e-01, -0.1903648724133519e-01, -0.2649038633167001e-01, 0.3481279195257287e-01, -0.3459307949232926e-01, -0.4406560344431272e-02}}; // // The error estimator constants // #include #include #include #include //#include #define sqr(x) ((x)*(x)) void Circle_Rule13::Apply(Integrand& F,Circle& R,real& TheResult,real& TheError) { //if (crival<0 || facmed <0) // { // cerr << " crival,facmed? " << endl; // cin >> crival; // cin >> facmed; // }; const real facopt = facmed/(crival*crival); int i,j,p,Type,nr,number; real Null[6]; real Tres = -1.0; real noise,r,r1,r2,Deg5,Deg3,Deg7,sumval; Point x[8]; if (Tres <= 0) { Tres=50*REAL_EPSILON; } TheResult = 0.0; for (i=0;i<6;i++){Null[i]=0;} p=0; for (Type=0;Type <=3;Type++){ for (nr=0;nr= 1) { TheError =10*Deg7; } else { if (r>=crival) { TheError = facmed*r*Deg7; } else { TheError = facopt *r*r*r*Deg7; } } TheError = max(noise,TheError); } TheError *= R.Volume(); TheResult *= R.Volume(); // cout <<"=> Error = "<() { } /////////////////////