/*Translated by FOR_C, v3.4.2 (-), on 07/09/115 at 08:32:07 */ /*FOR_C Options SET: ftn=u io=c no=p op=aimnv s=dbov str=l x=f - prototypes */ #include #include "fcrt.h" #include "sintm.h" #include /* PARAMETER translations */ #define KINT 29 #define KLOG 11 #define KREAL 169 /* end of PARAMETER translations */ /* COMMON translations */ struct t_sintnc { float ainit, binit, fncval, s, tp, fer, fer1, relobt, tps, xj, xjp; long int fea, fea1, kdim, inc, inc2, istop[2][2], jprint, iprint, kk, kmaxf, ndim, nfindx, nfmax, nfmaxm, reltol, reverm, revers, wherem; LOGICAL32 needh; } sintnc; struct t_sintc { double acum, pacum, result[2]; float aacum, local[4], abscis, ta, delta, delmin, diff, discx[2], end[2], errina, errinb, fat[2], fsave, funct[24], f2, paacum, pf1, pf2, phisum, phtsum, px, space[6], step[2], start[2], sum, t, tasave, tb, tend, worry[2], x1, x2, x, f1, count, xt[17], ft[17], phi[34], absdif, edue2a, edue2b, ep, epnoiz, epsmax, epso, epsr, epss, errat[2], errc, errf, errt[2], esold, extra, pepsmn, releps, rep, rndc, tlen, xjump, erri, err, epsmin, eps, re, reprod; long int discf, dischk, endpts, inew, iold, ip, ixkdim, j, j1, j1old, j2, j2old, kmax, kmin, l, lendt, nfjump, nsubsv, nxkdim, taloc, where2, i, k, kaimt, nsub, part, search, where, nfeval; LOGICAL32 did1, fail, fats[2], fsaved, havdif, iend, init, roundf, xcdobt[2], pad[7]; } sintc; struct t_sintec { float emeps, eepsm8, edelm2, edelm3, esqeps, ersqep, ersqe6, eminf, esmall, enzer, edelm1, eninf; } sintec; /* end of COMMON translations */ void /*FUNCTION*/ sintm( long ndimi, float *answer, float work[], long nwork, long iopt[]) { long int itemp; long int *const imove = (long*)&sintc.discf; LOGICAL32 *const lmove = (LOGICAL32*)&sintc.did1; float *const rmove = (float*)&sintc.aacum; /* OFFSET Vectors w/subscript range: 1 to dimension */ long *const Imove = &imove[0] - 1; long *const Iopt = &iopt[0] - 1; LOGICAL32 *const Lmove = &lmove[0] - 1; float *const Rmove = &rmove[0] - 1; float *const Work = &work[0] - 1; /* end of OFFSET VECTORS */ /* Copyright (c) 1996 California Institute of Technology, Pasadena, CA. * ALL RIGHTS RESERVED. * Based on Government Sponsored Research NAS7-03001. *>> 2009-11-03 SINTM Krogh Initialized lots of variables. *>> 1996-03-31 SINTM Krogh Removed unused variable in common. *>> 1995-11-20 SINTM Krogh Converted from SFTRAN to Fortran 77. *>> 1994-11-14 SINTM Krogh Declared all vars. *>> 1994-10-19 SINTM Krogh Changes to use M77CON *>> 1994-07-07 SINTM Snyder set up for CHGTYP. *>> 1993-05-18 SINTM Krogh -- Changed "END" to "END PROGRAM" *>> 1994-05-02 SINTM Snyder corrected some comments *>> 1991-09-20 SINTM Krogh converted '(1)' dimensioning to '(*)'. *>> 1987-11-19 SINTM Snyder Initial code. * *--S replaces "?": ?INTA,?intc,?intec,?INTF,?INTM,?INTMA,?INTNC,?INTOP * * ****************************************************************** * * THIS SUBROUTINE ATTEMPTS TO CALCULATE THE INTEGRAL OF A FUNCTION * PROVIDED BY THE USER. THE FUNCTION IS PROVIDED BY THE USER VIA A * SUBROUTINE REFERENCED BY CALL SINTF(F,X,IFLAG) OR BY REVERSE * COMMUNICATION. ALL ABSCISSAE ARE WITHIN THE REGION OF * INTEGRATION. * * THE RESULT IS OBTAINED USING QUADRATURE FORMULAE DUE TO * T. N. L. PATTERSON, MATHEMATICS OF COMPUTATION, VOLUME 22, * PAGES 847-856, 1968. * * ***** WARNING *********************************************** * * THE RELIABILITY AND EFFICIENCY OF THIS PROGRAM ARE STRONGLY * INFLUENCED BY DISCONTINUITIES IN THE FUNCTION OR IT'S * DERIVATIVES, INTEGRABLE SINGULARITIES IN THE REGION OF * INTEGRATION, AND NON-INTEGRABLE SINGULARITIES NEAR THIS REGION. * (INCLUDING COMPLEX VALUES). THE EFFICIENCY AND RELIABILITY * OF INTEGRATING SUCH FUNCTIONS MAY BE GREATLY IMPROVED BY MANUALLY * SUBDIVIDING THE REGION OF INTEGRATION AT THE DISCONTINUITY, * SINGULARITY OR CLOSEST POINT TO A POLE, AND SUMMING THE ANSWERS. * A CHANGE OF VARIABLE TO ELIMINATE OR REDUCE THE STRENGTH OF THE * SINGULARITY WILL SIGNIFICANTLY IMPROVE PERFORMANCE. * * ***** FORMAL ARGUMENTS ************************************* * * NDIMI IS THE NUMBER OF DIMENSIONS OF INTEGRATION. * ANSWER IS THE ESTIMATE OF THE INTEGRAL. WHEN REVERSE COMMUNICATION * IS SPECIFIED IT IS USED TO PASS FUNCTION VALUES FROM * THE USER PROGRAM TO SINTA OR SINTMA. */ /* WORK CONTAINS THE LIMITS, WORKING SPACE AND MAY BE REFERENCED BY * THE OPTION VECTOR (SEE IOPT BELOW). WHEN THE INTEGRATION * IS COMPLETE, WORK(1) CONTAINS THE ESTIMATED ABSOLUTE ERROR. * WHEN REVERSE COMMUNICATION IS SPECIFIED, WORK(1) IS USED TO * PASS ABSCISSAE FROM SINTA OR SINTMA TO THE USER PROGRAM. * WORK(1), ..., WORK(NDIMI) ARE COMPONENTS OF THE VECTOR OF * ABSCISSAE. WORK(NDIMI+1), ..., WORK(2*NDIMI) ARE LOWER * LIMITS. WORK(2*NDIMI+1), ..., WORK(3*NDIMI) ARE UPPER LIMITS. * THE ABSCISSAE AND LIMITS ARE STORED INNERMOST FIRST. SINTF * IS CALLED TO REQUEST THE LIMITS OF EVERY DIMENSION, BUT * CONSTANT LIMITS (THE LIMITS OF THE OUTER DIMENSION ARE ALWAYS * CONSTANT) MAY INSTEAD BE STORED BEFORE INTEGRATION BEGINS. * WORK(3*NDIMI+1), ..., WORK(3*NDIMI+KWORK*(NDIMI-1)) ARE * WORKING STORAGE, WHERE KWORK DEPENDS ON THE MACHINE AND THE * PRECISION OF THE PROGRAM VERSION. KWORK MUST SPECIFY ENOUGH * STORAGE FOR 4 VARIABLES WHICH ARE ALWAYS DOUBLE PRECISION, * 139 VARIABLES WHICH ARE DOUBLE PRECISION IN THE DOUBLE * PRECISION VERSION, 30 VARIABLES WHICH ARE DOUBLE PRECISION * IN THE DOUBLE PRECISION PROGRAM IF THE LARGEST REPRESENTABLE * NUMBERS IN SINGLE AND DOUBLE PRECISION ARE DIFFERENT BY MORE * THAN A FACTOR OF 4, 29 INTEGER VARIABLES AND 11 LOGICAL * VARIABLES. REMEMBER TO CONSIDER THE DATA TYPE OF WORK WHEN * SPECIFYING THE VALUE OF KWORK. EXAMPLES ARE *-- Begin mask code changes * IF PRECISION IS REAL, AND INTEGER AND LOGICAL VARIABLES * OCCUPY THE SAME SPACE AS REAL VARIABLES, KWORK=217. * IF PRECISION IS REAL, AND INTEGER AND LOGICAL VARIABLES * OCCUPY HALF AS MUCH SPACE AS REAL VARIABLES, KWORK=? * IF PRECISION IS PARTLY DOUBLE, AND INTEGER AND LOGICAL * VARIABLES OCCUPY THE SAME SPACE AS REAL VARIABLES, * KWORK=? * IF PRECISION IS PARTLY DOUBLE, AND INTEGER AND LOGICAL * VARIABLES OCCUPY HALF AS MUCH SPACE AS REAL VARIABLES, * KWORK=? * IF PRECISION IS DOUBLE, AND INTEGER AND LOGICAL VARIABLES * OCCUPY THE SAME SPACE AS REAL VARIABLES, KWORK=193. * IF PRECISION IS DOUBLE, AND INTEGER AND LOGICAL VARIABLES * OCCUPY HALF AS MUCH SPACE AS REAL VARIABLES, KWORK=? * FOR PORTABILITY, ASSUME THE WORST CASE. THAT IS, IF THE * PRECISION IS SINGLE, KWORK=217. IF THE PRECISION IS DOUBLE, * KWORK=193. * WORK(3*NDIMI+KWORK*(NDIMI-1)+1), ... ARE AVAILABLE FOR * REFERENCE THROUGH THE OPTION VECTOR. *-- End mask code changes */ /* NWORK IS THE NUMBER OF ELEMENTS ALLOCATED FOR WORK. * NWORK MUST BE GREATER THAN 3*NDIMI+KWORK*(NDIMI-1). * IOPT IS A VECTOR OF INTEGERS USED TO RETURN THE STATUS AND TO * SPECIFY OPTIONS. */ /* IOPT(1) IS USED TO RETURN A STATUS INDICATOR TO THE USER. * VALUES OF THIS INDICATOR ARE * -NDIM - NORMAL TERMINATION, EITHER THE ABSOLUTE OR THE * RELATIVE ERROR CRITERIA ARE SATISFIED. * -NDIM-1 - NORMAL TERMINATION, NEITHER THE ABSOLUTE NOR THE * RELATIVE ERROR CRITERIA ARE SATISFIED, BUT THE * ERROR RELATIVE TO THE OBTAINABLE PRECISION * PRECISION CRITERION IS SATISFIED. * -NDIM-2 - NORMAL TERMINATION, BUT NONE OF THE ERROR CRITERIA * ARE SATISFIED. * -NDIM-3 - NOT ENOUGH SPACE IN WORK, NWORK IS TOO SMALL. * -NDIM-4 - BAD IOPT VALUE. * -NDIM-5 - TOO MANY FUNCTION VALUES NEEDED. * -NDIM-KDIM-5 - APPARENT NON-INTEGRABLE SINGULARITY IN * DIMENSION KDIM. ANSWER CONTAINS THE APPROXIMATE * ABSCISSA OF THE SINGULARITY IN THE KDIM-TH * DIMENSION, WORK(KDIM+1) THROUGH WORK(NDIM) * CONTAIN THE ABSCISSAE OF EXTERIOR INTEGRALS. * -2*NDIM-KDIM-5 - INCORRECT INNER INTEGRAL DIMENSIONALITY * SPECIFIED DURING EVALUATION OF INTEGRAL AT * DIMENSION KDIM, THAT IS, IABS(IFLAG(1)) = KDIM. * ENTRIES IN IOPT STARTING WITH IOPT(2) ARE DESCRIBED BELOW. * IOPT(I) ENTRY IN IOPT(I) MEANS * 0 NO MORE OPTIONS, BEGIN INTEGRATION. * 1 IOPT(I+1) CONTAINS THE UNIT NUMBER FOR OUTPUT. * 2 IOPT(I+1) IS AN NDIMI DIGIT INTEGER, WHERE * EACH DIGIT IS THE DIAGNOSTIC PRINT LEVEL FOR * ONE DIMENSION. THE LOW ORDER DIGIT IS THE PRINT * LEVEL FOR THE INNER DIMENSION. * 0 - NO PRINTING * 1 - MINIMUM PRINTING - ERROR MESSAGES (DEFAULT) * 2 - PANEL BOUNDARIES AND ANSWERS * 3 - ERROR ESTIMATES FOR EACH QUADRATURE FORMULA * 4 - DETAILED OUTPUT (DIFFERENCE LINES, ETC). * 3 WORK(IOPT(I+1)) CONTAINS THE ABSOLUTE TOLERANCE, * WORK(IOPT(I+1)+1) CONTAINS TOLERANCE RELATIVE TO * THE LOCALLY OBTAINABLE PRECISION, AND * WORK(IOPT(I+1)+2) CONTAINS TOLERANCE RELATIVE TO * THE VALUE OF THE INTEGRAL. THE TOLERANCE RELATIVE * TO THE LOCALLY OBTAINABLE PRECISION IS SPECIFIED AS * THE FRACTION OF LOCALLY OBTAINABLE DIGITS THAT ARE * PERMITTED TO BE WRONG, AND IS INTERNALLY BOUNDED * BETWEEN 0.0 AND 1.0. IF ANY OF THE ERROR CONTROL * CRITERIA ARE SATISFIED, THE ANSWER IS ACCEPTED. IF * THIS OPTION IS NOT USED, THE ABSOLUTE AND RELATIVE * TOLERANCES ARE SET TO ZERO, AND THE TOLERANCE * RELATIVE TO THE LOCALLY OBTAINABLE PRECISION IS SET * TO 0.25. * 4 WORK(IOPT(I+1)) CONTAINS ABSOLUTE ERROR COMMITTED * COMPUTING F. WORK(IOPT(I+1)) MAY BE CHANGED DURING * THE INTEGRATION. * 5 WORK(IOPT(I+1)) CONTAINS RELATIVE ERROR EXPECTED TO * BE COMMITTED COMPUTING F. CHANGES TO THIS VALUE * DURING INTEGRATION WILL NOT BE DETECTED. * 6 USE REVERSE COMMUNICATION - * CALL SINTM (NDIMI,ANSWER,WORK,NWORK,IOPT) * 1 CALL SINTMA (ANSWER,WORK,IOPT(1)) * IF (IOPT(1)) 3,2,4 * 2 COMPUTE THE INNERMOST INTEGRAND - * ANSWER=F(WORK(K)), K=1,NDIMI * THEN EITHER * GO TO 1 * OR (FOR BETTER EFFICIENCY) * CALL SINTA (ANSWER,WORK,IOPT(1)) * IF (IOPT(1)) 1,2,5 * VALUES OF IOPT(1) PRODUCED BY SINTA AND SINTMA * HAVE DIFFERENT MEANINGS. * 3 IF (IOPT(1)+NDIMI.LE.0) GO TO 6 * THE INTEGRAL OVER THE (-IOPT(1))TH DIMENSION IS * CONTAINED IN ANSWER. IF A TRANSFORMATION OF * THE INTEGRAL IS NOT NECESSARY BEFORE IT IS USED * AS AN INTEGRAND FOR THE (1-IOPT(1))TH DIMENSION * GO TO 1 * ELSE COMPUTE ANSWER=G(ANSWER,WORK(K)), * K=1-IOPT(1),NDIMI AND G IS THE TRANSFORMATION. * COMPUTE ALSO WORK(1)=PARTIAL DERIVATIVE OF G * WITH RESPECT TO ANSWER THEN * GO TO 1 * 4 COMPUTE WORK(NDIMI+IOPT(1))=A(WORK(K)), * K=IOPT(1)+1,NDIMI, WORK(2*NDIMI+IOPT(1))=B(WORK(K)) * K=IOPT(1)+1,NDIM, WHERE A IS THE LOWER LIMIT OF * INTEGRATION FOR THE (IOPT(1))TH DIMENSION, AND B IS * THE UPPER LIMIT. IF A TRANSFORMATION AS DISCUSSED * ABOVE IS TO BE APPLIED WHEN THE INTEGRATION OF THIS * DIMENSION IS COMPLETE, COMPUTE WORK(1)=AN ESTIMATE * OF THE UPPER BOUND OF THE PARTIAL DERIVATIVE OF G * WITH RESPECT TO THE INTEGRAL. IF G IS LINEAR IN * THE INTEGRAL, THIS DERIVATIVE CAN BE COMPUTED * WITHOUT KNOWING THE INTEGRAL, RATHER THAN * ESTIMATED. IF THERE IS A SINGULARITY IN THE * INTEGRAND FOR THIS DIMENSION (PERHAPS INTRODUCED * BY THE BOUNDARY OR THE TRANSFORMATION G), SET * IOPT(1) TO THE LOCATION IN WORK OF THE ABSCISSA IN * THIS DIMENSION OF THE SINGULARITY, AND WITH THE * SIGN SET AS DESCRIBED FOR OPTION 11 BELOW. THE * MAGNITUDE OF IOPT(1) MUST BE GREATER THAN NDIM FOR * THIS REQUEST TO BE DETECTED. IN PARTICULAR, IF * THERE IS A SINGULARITY AT ONE OF THE LIMITS, SET * IOPT(1) TO IOPT(1)+NDIM OR -(IOPT(1)+2*NDIM). THEN * GO TO 1. */ /* 5 IOPT(1)=-(IOPT(1)+NDIM) * 6 CONTINUE * AT THIS POINT, IERR CONTAINS THE STATUS DESCRIBED * BELOW, FOR OPTION 8. * 7 SET MINIMUM INDEX OF QUADRATURE FORMULA TO USE * BEFORE SUBDIVISION TO IOPT(I+1). * 8 NO EFFECT. IOPT(I+1) MAY BE USED TO PASS * INFORMATION TO SINTF. INCREMENT I BY 2. * 9 IOPT(I+1) IS THE MAXIMUM NUMBER OF FUNCTION VALUES * TO USE TO INTEGRATE THE ENTIRE PROBLEM. IF * IOPT(I+1) .LE. 0, THE NUMBER OF FUNCTION VALUES * IS NOT CONTROLLED. * 10 IOPT(I+1) IS USED TO RETURN THE NUMBER OF FUNCTION * VALUES USED TO INTEGRATE THE ENTIRE PROBLEM. * 11 WORK(IABS(IOPT(I+1))) IS THE LOCATION OF A * SINGULARITY OR DISCONTINUITY IN THE OUTERMOST * INTEGRAL. IF THE LOCATION IS INSIDE THE INTERVAL, * THE INTERVAL IS IMMEDIATELY SUBDIVIDED. IF * IOPT(I+1) .GT. 0, A T**2 SUBSTITUTION BASED AT * WORK(...) WILL BE USED. IF IOPT(I+1) .LT. 0, * A T**4 SUBSTITUTION WILL BE USED. IF IOPT(I+1) * .EQ. 0, NO SUBSTITUTION WILL BE USED. TO NOTIFY * THE PROGRAM OF SINGULARITIES OR DISCONTINUITIES IN * INTERIOR INTEGRANDS, SET IFLAG NEGATIVE WHEN * ASKED TO COMPUTE LIMITS. SEE THE DESCRIPTION * OF REVERSE COMMUNICATION ABOVE, OR THE DESCRIPTION * OF SINTF BELOW. * 12 WORK(IOPT(I+1)) CONTAINS THE ABSOLUTE ERROR IN THE * LOWER LIMIT OF THE CURRENT DIMENSION. THE ERROR IN * THE UPPER LIMIT IS IN WORK(IOPT(I+1)+1). * 13 IOPT(I+1) IS THE LOCATION IN IOPT IN WHICH THE USER * IS TO STORE NON-STANDARD CHANGES TO THE DIMENSION OF * INNER INTEGRALS. THE DEFAULT VALUE IS 1, ALLOWING * SUCH NONSTANDARD CHANGES TO BE STORED IN THE CELL OF * IFLAG (=IOPT) USUALLY USED FOR COMMUNICATION WITH * SINTF, IN THE CASE WHEN IT IS NOT SIMULTANEOUSLY * NECESSARY TO INDICATE A SINGULARITY. SEE THE * DESCRIPTION OF THE INTERFACE TO SINTF BELOW. * * ALL OPTIONS ARE SET TO THEIR NOMINAL VALUES BEFORE THE OPTION * VECTOR IS PROCESSED. * * ***** SINTF *********************************************** * * SINTF IS REFERENCED VIA CALL SINTF (ANSWER,WORK,IFLAG). * WORK AND ANSWER ARE AS DESCRIBED ABOVE. VALUES OF IFLAG ARE * IFLAG .LT. 0 * THE INTEGRAL OVER THE (-IFLAG)TH DIMENSION IS * CONTAINED IN ANSWER. IF A TRANSFORMATION OF * THE INTEGRAL IS NOT NECESSARY BEFORE IT IS USED * AS AN INTEGRAND FOR THE (1-IFLAG)TH DIMENSION * RETURN * ELSE IF THE INTEGRAND DEPENDS ON ONLY ONE INNER * INTEGRAL COMPUTE ANSWER=G(ANSWER,WORK(K)), * ALSO MULTIPLY WORK(1) BY THE PARTIAL DERIVATIVE * OF G WITH RESPECT TO ANSWER THEN * RETURN * ELSE IF THE INTEGRAND DEPENDS ON MORE INTEGRALS * SAVE ANSWER AND WORK(1), AND SET IFLAG(IXKDIM) TO * THE DIMENSION OF INTEGRAL NEXT TO BE COMPUTED * (IXKDIM IS THE VALUE SPECIFIED BY OPTION 12, OR 1 * IF OPTION 12 IS NOT SELECTED), THEN * RETURN * ELSE (ALL INNER INTEGRALS UPON WHICH THE INNER * INTEGRAND DEPENDS HAVE BEEN EVALUATED) CALCULATE * ANSWER=G(ANSWER,WORK(K),SAVED INTEGRALS) * K=1-IFLAG, NDIMI (G IS THE TRANSFORMATION), * MULTIPLY WORK(1) BY THE PARTIAL DERIVATIVE OF G * WITH RESPECT TO THE LAST INTEGRAL, AND ADD THE * PRODUCT OF EACH OF THE SAVED VALUES OF WORK(1) AND * PARTIAL DERIVATIVE OF G WITH RESPECT TO THE * CORRESPONDING SAVED INTEGRAL ONTO WORK(1). THAT * IS, WORK(1) IS THE TOTAL ERROR IN G, CALCULATED AS * THE INNER PRODUCT OF THE ERRORS IN THE INDIVIDUAL * ARGUMENTS, AND THE PARTIAL DERIVATIVES OF G WITH * RESPECT TO THOSE ARGUMENTS. * IFLAG .EQ. 0 * COMPUTE THE INNERMOST INTEGRAND - * ANSWER=F(WORK(K)), K=1,NDIMI * AND * RETURN * IFLAG .GT. 0 * COMPUTE WORK(NDIMI+IFLAG)=A(WORK(K)), * K=IFLAG+1,NDIMI, WORK(2*NDIMI+IFLAG)=B(WORK(K)), * K=IFLAG+1,NDIM, WHERE A IS THE LOWER LIMIT OF * INTEGRATION FOR THE (IFLAG)TH DIMENSION, AND B IS * THE UPPER LIMIT. IF A TRANSFORMATION AS DISCUSSED * ABOVE IS TO BE APPLIED WHEN THE INTEGRATION OF THIS * DIMENSION IS COMPLETE, COMPUTE WORK(1)=AN ESTIMATE * OF THE UPPER BOUND OF THE PARTIAL DERIVATIVE * OF G WITH RESPECT TO THE INTEGRAL. IF G IS * LINEAR IN THE INTEGRAL, THIS DERIVATIVE MAY BE * COMPUTED WITHOUT KNOWING THE INTEGRAL, RATHER * THAN ESTIMATED. IF THERE IS A SINGULARITY IN * THE INTEGRAND FOR THIS DIMENSION (PERHAPS * INTRODUCED BY THE BOUNDARY OR THE TRANSFORMATION), * SET IFLAG TO THE LOCATION IN WORK OF THE ABSCISSA * IN THIS DIMENSION OF THE SINGULARITY, AND SET THE * SIGN AS FOR OPTION 11 ABOVE. NOTE THAT THE * MAGNITUDE OF IFLAG MUST BE GREATER THAN NDIM FOR */ /* THIS REQUEST TO BE DETECTED. IN PARTICULAR, IF * THERE IS A SINGULARITY AT ONE OF THE LIMITS, SET * THE MAGNITUDE OF IFLAG TO IFLAG+NDIM OR * IFLAG+2*NDIM, WITH THE SIGN AS DESCRIBED FOR * OPTION 11. THEN * RETURN. * IF THE FIRST INTEGRAL TO BE COMPUTED DOES NOT HAVE * DIMENSION IFLAG(1)-1, SET IFLAG(IXKDIM) TO THE * DIMENSION OF INTEGRAL TO BE NEXT COMPUTED. SEE * OPTION 12 AND THE DISCUSSION OF IFLAG .LT. 0 ABOVE * FOR INFORMATION ON IXKDIM. * * ***** EXTERNAL REFERENCES ********************************* * * SINTMA TO DO THE INTEGRATION. * SINTOP TO SET OPTIONS. * * ***** LOCAL VARIABLES ************************************* * */ /* ***** COMMON STORAGE ************************************** * * COMMON /SINTNC/ CONTAINS VARIABLES NOT SEPARATELY SAVED FOR * EACH DIMENSION OF A MULTIPLE QUADRATURE. COMMON /SINTC/ * CONTAINS VARIABLES THAT MUST BE SAVED FOR EACH DIMENSION OF THE * QUADRATURE. THE VARIABLES IN EACH COMMON BLOCK ARE STORED IN THE * ORDER - ALWAYS DOUBLE, DOUBLE IF DOUBLE PRECISION PROGRAM, DOUBLE * IF DOUBLE PRECISION PROGRAM AND EXPONENT RANGE OF DOUBLE AND * SINGLE VERY DIFFERENT, SINGLE, INTEGER, LOGICAL. A PAD OF LOGICAL * VARIABLES IS INCLUDED AT THE END OF /SINTC/. THE DIMENSION OF * THE PAD MAY NEED TO BE VARIED SO THAT NO VARIABLES BEYOND THE END * OF THE COMMON BLOCK ARE ALTERED. * * DECLARATIONS OF COMMON /SINTNC/ VARIABLES. * */ /* DECLARATIONS OF COMMON /SINTC/ VARIABLES. * *--D Next line special: S => D, X => Q, D => D, P => D */ /* 139 $.TYPE.$ VARIABLES */ /* Note XT, FT, and PHI above are last, because they must be in adjacent * locations in SINTC. * 30 $DSTYP$ VARIABLES */ /* 29 INTEGER VARIABLES */ /* 11 TO 18 LOGICALS (7 ARE PADDING). */ /* THE COMMON BLOCKS. * */ /* 1 2 3 4 5 6 7 8 * 9 10 11 12 13 1 2 3 * 4 (2,2) 8 9 10 11 12 13 14 * 15 16 17 18 19 20 */ /* 1 2 (4) 6 7 8 9 10 11 (2) * 13 (2) 15 16 17 (2) 19 20 (24) 44 * 45 46 47 48 49 50 51 (6) * 57 (2) 59 (2) 61 62 63 64 65 * 66 (2) 68 69 70 71 72 * 73 (17) 90 (17) 107 (34) */ /* 141 142 143 144 145 146 * 147 148 149 150 (2) 152 153 * 154 (2) 156 157 158 159 160 * 161 162 163 * 164 165 166 167 168 169 */ /* 170 171 172 * 1 2 3 4 5 6 7 8 */ /* THE VARIABLES HERE DEFINE THE MACHINE ENVIRONMENT. ALL ARE SET * IN SINTOP. THE MEANING ATTACHED TO THESE VARIABLES CAN BE * FOUND BY LOOKING AT THE DEFINITIONS IN SINTOP. */ /* For initialization of common blocks */ /* ***** PROCEDURES ****************************************** * */ /* Initialize common blocks to avoid references to undefined variables. */ for (itemp = 1; itemp <= KLOG; itemp++) { Lmove[itemp] = FALSE; } for (itemp = 1; itemp <= KREAL; itemp++) { Rmove[itemp] = 0.e0; } for (itemp = 1; itemp <= KINT; itemp++) { Imove[itemp] = 0; } sintnc.wherem = 0; sintc.nfeval = 0; sintnc.ndim = ndimi; sintnc.kdim = nwork; sintc.x = 0.e0; /* KDIM IS TEMPORARILY USED IN SINTMA TO CHECK THE DIMENSION OF * WORK. */ sintop( iopt, work ); if (Iopt[1] == 0) { /* CALL SINTMA TO DO THE INTEGRATION. * */ if (sintnc.reverm == 0) sintma( answer, work, iopt ); } else { Iopt[1] = -sintnc.ndim - 4; } return; } /* end of function */