C OPTABJ SOURCE CB215821 23/10/18 21:15:09 11760 SUBROUTINE OPTABj(NBTHR ,ITHR,IOPERA,NTABEN, & XVAL0,XVAL1,XVAL2, & NN0 ,NN1 ,NN2 ,IARG ,I1I ,X1I ,IRETOU) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C Cette subroutine effectue des operations elementaires ainsi que C les fonctions sur des tableaux FORTRAN de REAL*8 C Elle est prevue pour etre executee en parallele CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C NBTHR : Nombre de Thread disponibles C ITHR : Numero du Thread courant C IRETOU : Entier contenant le code d'erreur C IOPERA : Type d''operation a realiser (Voir ci-dessous) C NTABEN : Nombre de tableaux constituant l''entree (Exemple : 2 pour ATAN2) C XVAL0 : Tableau de valeur d''entree C XVAL1 : Tableau de valeur d''entree (Deuxieme argument pour ATAN2) C XVAL2 : Tableau de valeur de sortie C NN0 : Taille du tableau XVAL0 C NN1 : Taille du tableau XVAL1 C NN2 : Taille du tableau XVAL2 C C IARG = 0 ==> ARGUMENT I1I ET X1I INUTILISES C IARG = 1 ==> ARGUMENT I1I UTILISE C IARG = 11 ==> ARGUMENT I1I UTILISE MAIS COMMUTE AVEC LE TABLEAU (PUISSANCE, SOUSTRACTION, DIVISION : POSITIONNEL) C IARG = 2 ==> ARGUMENT X1I UTILISE C IARG = 21 ==> ARGUMENT X1I UTILISE MAIS COMMUTE AVEC LE TABLEAU (PUISSANCE, SOUSTRACTION, DIVISION : POSITIONNEL) C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C Elle realise les operations suivantes : C Operations elementaires entre un TABLEAU et un ENTIER ou FLOTTANT C IOPERA= 1 PUISSANCE C = 2 PRODUIT C = 3 ADDITION C = 4 SOUSTRACTION C = 5 DIVISION C C Fonctions sur un TABLEAU C = 6 COSINUS C = 7 SINUS C = 8 TANGENTE C = 9 ARCOSINUS C = 10 ARCSINUS C = 11 ARCTANGENTE C = 12 EXPONENTIELLE C = 13 LOGARITHME C = 14 VALEUR ABSOLUE C = 15 COSINUS HYPERBOLIQUE C = 16 SINUS HYPERBOLIQUE C = 17 TANGENTE HYPERBOLIQUE C = 18 ERF FONCTION D''ERRREUR DE GAUSS C = 19 ERFC FONCTION D''ERRREUR COMPLEMENTAIRE DE GAUSS (1-ERF(X)) C = 20 ARGCH (FONCTION RECIPROQUE DE COSH) C = 21 ARGSH (FONCTION RECIPROQUE DE SINH) C = 22 ARGTH (FONCTION RECIPROQUE DE TANH) C = 23 SIGN (renvoie -1 ou +1) C = 24 BESSEL J0 C = 25 BESSEL J1 C = 26 BESSEL Y0 C = 27 BESSEL Y1 C = 28 FRESNEL CX C = 29 FRESNEL SX C = 30 GAMMA (Fonction Gamma d'Euler) C = 31 BESSEL JN (Fonction BESSEL de type J d'ordre N) C = 32 BESSEL YN (Fonction BESSEL de type Y d'ordre N) C HISTORIQUE : C - CB215821 31/08/2016 --> Creation C - CB215821 05/06/2018 --> Ajout de la fonction SIGN a un argument C - CB215821 17/10/2023 --> Ajout des fonctions BESSEL, FRESNEL et GAMMA C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC IMPLICIT INTEGER(I-N) IMPLICIT REAL*8 (A-H,O-Z) -INC CCREEL -INC PPARAM -INC CCOPTIO SEGMENT ITOTO(0) INTEGER NTABEN REAL*8 XNOR,XINV,XFLOT,XFLOT1,XFLOT2,XIINV,XPREC,XTRA,X2,XF1,UN PARAMETER (XNOR = XPI / 180.D0) PARAMETER (XINV = 180.D0 / XPI) PARAMETER (UN = 1.D0) REAL*8 XVAL0(NN0),XVAL1(NN1),XVAL2(NN2) INTEGER IRETOU I2 = I1I X2 = X1I IARG2 = IARG IRETOU = 0 C On assure le travail contigu en memoire NNC = MAX(NN0,NN1,NN2) IF(NBTHR .EQ. 1)THEN IDEB = 1 IFIN = NNC ELSE IRES = MOD(NNC,NBTHR) IF(IRES .EQ. 0)THEN ILON = NNC / NBTHR IDEB = (ithr -1)* ILON + 1 ELSE IF (ithr .LE. IRES) THEN ILON = (NNC / NBTHR) + 1 IDEB = (ithr -1)* ILON + 1 ELSE ILON = NNC / NBTHR IDEB = (IRES * (ILON+1)) + (ithr-IRES-1)* ILON + 1 ENDIF ENDIF IFIN = IDEB + ILON - 1 ENDIF C PRINT *, 'OPTABJ:',ITHR,IOPERA,NN0,NTABEN IF (NTABEN .EQ. 2) GOTO 5000 C======================================================================C C OPERATIONS ENTRE UN TABLEAU ET UN FLOTTANT / ENTIER C======================================================================C C IF (IARG2 .EQ. 1 .OR. IARG2 .EQ. 11) THEN C X2 = REAL(I1I) C ELSEIF(IARG2 .NE. 0 .AND. IARG2 .NE. 2 .AND. IARG2 .NE. 21) THEN C Surveillance de la validite des PARAMETRES d'entree C IRETOU = 21 C RETURN C ENDIF GOTO ( 1, 2, 3, 4, 5, 6, 7, 8, 9,10,11,12,13,14,15,16,17, & 18,19,20,21,22,23,24,25,26,27,28,29,30,31,32 ),IOPERA C Erreur si l''operation demandee n''est pas dans la liste IRETOU = 21 RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C PUISSANCE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 1 CONTINUE IF (IARG2 .EQ. 1) THEN C PRINT *,'TABLEAU ** ENTIER',ITHR IF (I2 .EQ. 0) THEN C PRINT *,' Cas TABLEAU ** 0' DO 101 IA = IDEB,IFIN XVAL2(IA)= REAL(1.D0) 101 CONTINUE RETURN ELSEIF(I2 .EQ. 1)THEN C PRINT *,' Cas TABLEAU ** 1' DO 102 IA = IDEB,IFIN XVAL2(IA)= XVAL0(IA) 102 CONTINUE RETURN ELSE DO 103 IA = IDEB,IFIN XTRA=XVAL0(IA) IF(ABS(XTRA).LE.XPETIT .AND. I2.LT.0)THEN IRETOU = 213 RETURN ELSE XVAL2(IA)= XTRA ** I2 ENDIF 103 CONTINUE RETURN ENDIF RETURN ELSEIF(IARG2 .EQ. 2) THEN C PRINT *,'TABLEAU ** FLOTTANT',ITHR I2 = NINT(X2) XFLOT = ABS(X2 - REAL(I2)) XPREC = XZPREC*ABS(X2) C Verification si puissance ENTIERE possible IF ( XFLOT .LE. XPREC) THEN IARG2=1 GOTO 1 ENDIF C Verification si SQRT possible XF1 = X2 - REAL(0.5D0) I2 = NINT(XF1) I3 = (I2 * 2) + 1 XFLOT = ABS(XF1 - REAL(I2)) IF (XFLOT .LE. XPREC) THEN IF (I2 .EQ. 0) THEN C PRINT *,' Cas SQRT simple' DO 104 IA = IDEB,IFIN IF (XVAL0(IA) .LT. REAL(0.D0)) THEN IRETOU = 213 RETURN ELSE XVAL2(IA)= SQRT(XVAL0(IA)) ENDIF 104 CONTINUE RETURN ELSE C PRINT *,' Nouveau cas SQRT ** I3',I3 DO 105 IA = IDEB,IFIN IF (XVAL0(IA) .LT. REAL(0.D0)) THEN IRETOU = 213 RETURN ELSE XVAL2(IA)= (SQRT(XVAL0(IA))) ** I3 ENDIF 105 CONTINUE RETURN ENDIF RETURN ELSE C Verification si racine Nieme possible IF (X2 .GT. XPETIT) THEN XIINV=UN/X2 IINV = NINT(XIINV) XFLOT= ABS(XIINV - REAL(IINV)) XPREC= XZPREC*ABS(XIINV) IF (XFLOT .LE. XPREC .AND. MOD(IINV,2).NE. 0) THEN C PRINT *,' Racine Nieme' DO 106 IA = IDEB,IFIN XFLOT = XVAL0(IA) XVAL2(IA)=SIGN(UN,XFLOT)*(ABS(XFLOT)**X2) 106 CONTINUE RETURN ENDIF ENDIF C PRINT *,' Cas general' DO 107 IA = IDEB,IFIN IF ((ABS(XVAL0(IA)) .LE. XPETIT) .AND. & (X2 .LT. REAL(0.D0))) THEN IRETOU = 213 RETURN ELSEIF (XVAL0(IA) .LT. REAL(0.D0)) THEN IRETOU = 213 RETURN ELSE XVAL2(IA)= XVAL0(IA) ** X2 ENDIF 107 CONTINUE RETURN ENDIF RETURN ELSEIF(IARG2 .EQ. 11 .OR. IARG2 .EQ. 21) THEN C PRINT *,'ENTIER ** TABLEAU ou FLOTTANT ** TABLEAU' IF(IARG2 .EQ. 11) X2 = REAL(I1I) DO 108 IA = IDEB,IFIN I2 = NINT(XVAL0(IA)) XFLOT1= ABS(XVAL0(IA) - REAL(I2 )) XFLOT2= ABS(XVAL0(IA) - REAL(0.5D0)) XPREC = XZPREC*ABS(XVAL0(IA)) IF (((ABS(X2) .LE. XPETIT) .AND. & (XVAL0(IA).LT.REAL(0.D0) )) .OR. & (X2 .LT. REAL(0.D0))) THEN IRETOU = 213 RETURN ELSEIF ( XFLOT1 .LE. XPREC ) THEN C PRINT *,' Puissance Entiere Possible' XVAL2(IA)= X2 ** I2 ELSEIF ( XFLOT2 .LE. XPREC) THEN C PRINT *,' SQRT Possible' XVAL2(IA)= SQRT(X2) ELSE C PRINT *,' Cas general' XVAL2(IA)= X2 ** XVAL0(IA) ENDIF 108 CONTINUE RETURN ENDIF RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C PRODUIT CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 2 CONTINUE IF(IARG2 .EQ. 1 .OR. IARG2 .EQ. 11) THEN IF (I1I .EQ. 0) RETURN X2 = REAL(I1I) ENDIF IF(ABS(1.D0-X2) .LE. XZPREC)THEN IF(SIGN(1.D0,X2) .LT. 0.D0)THEN DO 201 IA = IDEB,IFIN XVAL2(IA)=-XVAL0(IA) 201 CONTINUE ELSE DO 202 IA = IDEB,IFIN XVAL2(IA)= XVAL0(IA) 202 CONTINUE ENDIF ELSE DO 203 IA = IDEB,IFIN XVAL2(IA)= XVAL0(IA) * X2 203 CONTINUE ENDIF RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ADDITION CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 3 CONTINUE IF ((IARG2 .EQ. 1 .OR. IARG2 .EQ. 11) .AND. (I1I .EQ. 0)) THEN C PRINT *,'ADDITION Cas 1' DO 301 IA=IDEB,IFIN XVAL2(IA) = XVAL0(IA) 301 CONTINUE RETURN ENDIF IF(IARG2 .EQ. 1 .OR. IARG2 .EQ. 11) X2 = REAL(I1I) C PRINT *,'ADDITION Cas 2' DO 303 IA=IDEB,IFIN XVAL2(IA) = XVAL0(IA) + X2 303 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C SOUSTRACTION CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 4 CONTINUE IF (IARG2 .EQ. 1 .AND. I1I .EQ. 0) THEN C PRINT *,'SOUSTRACTION Cas 1' DO 401 IA=IDEB,IFIN XVAL2(IA) = XVAL0(IA) 401 CONTINUE RETURN ENDIF IF (IARG2 .EQ. 11 .AND. I1I .EQ. 0) THEN C PRINT *,'SOUSTRACTION Cas 2' DO 403 IA=IDEB,IFIN XVAL2(IA) = -XVAL0(IA) 403 CONTINUE RETURN ENDIF IF (IARG2 .EQ. 1) THEN C Cas TABLEAU - I1I C PRINT *,'SOUSTRACTION Cas 3' X2 = REAL(I1I) DO 405 IA = IDEB,IFIN XVAL2(IA)= XVAL0(IA) - X2 405 CONTINUE RETURN ELSEIF(IARG2 .EQ. 2) THEN C Cas TABLEAU - X2 C PRINT *,'SOUSTRACTION Cas 4' DO 406 IA = IDEB,IFIN XVAL2(IA)= XVAL0(IA) - X2 406 CONTINUE RETURN ELSEIF(IARG2 .EQ. 11) THEN C Cas I1I - TABLEAU C PRINT *,'SOUSTRACTION Cas 5' X2 = REAL(I1I) DO 407 IA = IDEB,IFIN XVAL2(IA)= X2 - XVAL0(IA) 407 CONTINUE RETURN ELSEIF(IARG2 .EQ. 21) THEN C Cas X2 - TABLEAU C PRINT *,'SOUSTRACTION Cas 6' DO 408 IA = IDEB,IFIN XVAL2(IA)= X2 - XVAL0(IA) 408 CONTINUE RETURN ENDIF RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C DIVISION CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 5 CONTINUE IF(IARG2 .EQ. 11 .AND. I1I .EQ. 0) THEN C PRINT *,'DIVISION Cas 1' RETURN C ELSEIF(IARG2 .EQ. 21 .AND. (ABS(X2) .LE. XPETIT)) THEN CC PRINT *,'DIVISION Cas 2' XPETIT divise par qqc n'est pas forcement negligeable !!! C RETURN ENDIF IF (IARG2 .EQ. 1) THEN IF (I1I .NE. 0) THEN C PRINT *,'DIVISION Cas 3' X2 = REAL(I1I) DO 501 IA = IDEB,IFIN C Cas TABLEAU / ENTIER XVAL2(IA)= XVAL0(IA) / X2 501 CONTINUE RETURN ELSE IRETOU = 835 RETURN ENDIF ELSEIF (IARG2 .EQ. 2) THEN IF (ABS(X2) .GT. XPETIT) THEN C PRINT *,'DIVISION Cas 4' X3 = 1.D0 / X2 DO 502 IA = IDEB,IFIN C Cas TABLEAU / FLOTTANT XVAL2(IA)= XVAL0(IA) * X3 502 CONTINUE RETURN ELSE IRETOU = 835 RETURN ENDIF ELSEIF(IARG2 .EQ. 11 .OR. IARG2 .EQ. 21) THEN C PRINT *,'DIVISION Cas 5' IF (IARG2 .EQ. 11) X2 = REAL(I1I) DO 503 IA = IDEB,IFIN C Cas FLOTTANT / TABLEAU ou ENTIER / TABLEAU (terme a terme) IF (ABS(XVAL0(IA)) .GT. XPETIT) THEN XVAL2(IA)= X2 / XVAL0(IA) ELSE IRETOU = 835 RETURN ENDIF 503 CONTINUE RETURN ENDIF RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C COSINUS CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 6 CONTINUE DO 601 IA = IDEB,IFIN XVAL2(IA)= COS(XNOR * XVAL0(IA)) 601 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C SINUS CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 7 CONTINUE DO 701 IA = IDEB,IFIN XVAL2(IA)= SIN(XNOR * XVAL0(IA)) 701 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C TANGENTE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 8 CONTINUE DO 801 IA = IDEB,IFIN XVAL2(IA)= TAN(XNOR * XVAL0(IA)) 801 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ARCCOS CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 9 CONTINUE DO 901 IA = IDEB,IFIN X2 = XVAL0(IA) IF (ABS(X2) .LE. UN) THEN XVAL2(IA)= XINV * ACOS(X2) ELSE IRETOU = 21 RETURN ENDIF 901 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ARCSIN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 10 CONTINUE DO 1001 IA = IDEB,IFIN X2 = XVAL0(IA) IF (ABS(X2) .LE. UN) THEN XVAL2(IA)= XINV * ASIN(X2) ELSE IRETOU = 21 RETURN ENDIF 1001 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ARCTANGENTE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 11 CONTINUE DO 1101 IA = IDEB,IFIN XVAL2(IA)= XINV * ATAN(XVAL0(IA)) 1101 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C EXPONENTIELLE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 12 CONTINUE DO 1201 IA = IDEB,IFIN XVAL2(IA)= EXP(XVAL0(IA)) 1201 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C LOGARITHME CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 13 CONTINUE DO 1301 IA = IDEB,IFIN X2 = XVAL0(IA) IF (X2 .GE. REAL(0.D0) .AND. X2 .LE. XPETIT) THEN XVAL2(IA)=-1.D0*XGRAND ELSEIF(X2 .GT. XPETIT) THEN XVAL2(IA)= LOG(X2) ELSE IRETOU = 21 RETURN ENDIF 1301 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C VALEUR ABSOLUE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 14 CONTINUE DO 1401 IA = IDEB,IFIN XVAL2(IA)= ABS(XVAL0(IA)) 1401 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C COSINUS HYPERBOLIQUE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 15 CONTINUE DO 1501 IA = IDEB,IFIN XVAL2(IA)= COSH(XVAL0(IA)) 1501 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C SINUS HYPERBOLIQUE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 16 CONTINUE DO 1601 IA = IDEB,IFIN XVAL2(IA)= SINH(XVAL0(IA)) 1601 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C TANGENTE HYPERBOLIQUE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 17 CONTINUE DO 1701 IA = IDEB,IFIN XVAL2(IA)= TANH(XVAL0(IA)) 1701 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ERF (Fonction Erreur) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 18 CONTINUE DO 1801 IA = IDEB,IFIN XVAL2(IA)= ERF(XVAL0(IA)) 1801 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ERFC (Fonction Erreur Complementaire 1-ERF(x)) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 19 CONTINUE DO 1901 IA = IDEB,IFIN XVAL2(IA)= ERFC(XVAL0(IA)) 1901 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ARCOSH CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 20 CONTINUE DO 2001 IA = IDEB,IFIN X2 = XVAL0(IA) IF (X2 .GE. UN) THEN XVAL2(IA)= LOG(X2 + SQRT((X2**2) - UN)) ELSE IRETOU = 21 RETURN ENDIF 2001 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ARSINH CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 21 CONTINUE DO 2101 IA = IDEB,IFIN X2 = XVAL0(IA) XVAL2(IA)= LOG(X2 + SQRT((X2**2) + UN)) 2101 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ARTANH CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 22 CONTINUE DO 2201 IA = IDEB,IFIN X2 = XVAL0(IA) IF (ABS(X2) .LT. UN) THEN XVAL2(IA)=REAL(0.5D0)*LOG((UN+X2) / (UN - X2)) ELSE IRETOU = 21 RETURN ENDIF 2201 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C SIGN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 23 CONTINUE DO 2301 IA = IDEB,IFIN XVAL2(IA)= SIGN(UN,XVAL0(IA)) 2301 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C BESSEL J0 (Fortran 2008) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 24 CONTINUE DO 2401 IA = IDEB,IFIN XVAL2(IA)= BESSEL_J0(XVAL0(IA)) 2401 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C BESSEL J1 (Fortran 2008) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 25 CONTINUE DO 2501 IA = IDEB,IFIN XVAL2(IA)= BESSEL_J1(XVAL0(IA)) 2501 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C BESSEL Y0 (Fortran 2008) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 26 CONTINUE DO 2601 IA = IDEB,IFIN X0 = XVAL0(IA) IF (X0 .GT. XZERO) THEN XVAL2(IA)= BESSEL_Y0(XVAL0(IA)) ELSE IRETOU = 21 RETURN ENDIF 2601 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C BESSEL Y1 (Fortran 2008) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 27 CONTINUE DO 2701 IA = IDEB,IFIN X0 = XVAL0(IA) IF (X0 .GT. XZERO) THEN XVAL2(IA)= BESSEL_Y1(XVAL0(IA)) ELSE IRETOU = 21 RETURN ENDIF 2701 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C FRESNEL CX CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 28 CONTINUE DO 2801 IA = IDEB,IFIN C XVAL2(IA)= SIGN(UN,XVAL0(IA)) 2801 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C FRESNEL SX CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 29 CONTINUE DO 2901 IA = IDEB,IFIN C XVAL2(IA)= SIGN(UN,XVAL0(IA)) 2901 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C GAMMA d'Euler (Fortran 2008) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 30 CONTINUE DO 3001 IA = IDEB,IFIN X0 = XVAL0(IA) C X0 ne peut pas etre egal a ZERO ni un entier negatif I0 = NINT(X0) XI0= REAL(I0) IF(A_EGALE_B(XZERO,X0) .OR. & ((I0 .LT. 0) .AND. A_EGALE_B(X0,XI0)))THEN IRETOU = 21 RETURN ELSE XVAL2(IA) = GAMMA(X0) ENDIF 3001 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C BESSEL JN (Fortran 2008) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 31 CONTINUE DO 3101 IA = IDEB,IFIN XVAL2(IA)= BESSEL_JN(I1I,XVAL0(IA)) 3101 CONTINUE RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C BESSEL YN (Fortran 2008) CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 32 CONTINUE DO 3201 IA = IDEB,IFIN X0 = XVAL0(IA) IF (X0 .GT. XZERO) THEN XVAL2(IA)= BESSEL_YN(I1I,X0) ELSE IRETOU = 21 RETURN ENDIF 3201 CONTINUE RETURN C======================================================================C C OPERATIONS TERMES A TERMES DE 2 TABLEAUX C======================================================================C 5000 CONTINUE IF (NN0 .GT. NN1)THEN NSAUT1 = NN0/NN1 NSAUT2 = 0 ELSEIF(NN1 .GT. NN0)THEN NSAUT1 = 0 NSAUT2 = NN1/NN0 ELSE NSAUT1 = 1 NSAUT2 = 1 ENDIF GOTO ( 5001,5002,5003,5004,5005,9999,9999,9999,9999,9999,5011 ), & IOPERA C Erreur si l''operation demandee n''est pas dans la liste 9999 CONTINUE IRETOU = 21 RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C PUISSANCE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 5001 CONTINUE IF (NN0 .EQ. NN1 ) THEN C OPERATION TERME A TERME DO 5101 IA = IDEB,IFIN XFLO = XVAL0(IA) X2 = XVAL1(IA) IF ( ((ABS(XFLO) .LE. XPETIT) .AND. (X2 .LT. REAL(0.D0))) .OR. & (XFLO .LT. REAL(0.D0))) THEN IRETOU = 213 RETURN ELSE I2 = NINT(X2) XFLOT1 = ABS (X2 - REAL(I2)) XFLOT = ABS (X2 - REAL(0.5D0)) IF ( XFLOT1 .LE. XZPREC*ABS(X2)) THEN C PUISSANCE ENTIERE XVAL2(IA)= XFLO ** I2 ELSEIF (XFLOT .LE. XZPREC*ABS(X2) ) THEN C RACINE CARREE SQRT XVAL2(IA)= SQRT(XFLO) ELSE C CAS GENERAL XVAL2(IA)= XFLO ** X2 ENDIF ENDIF 5101 CONTINUE ELSEIF(NN0 .GT. NN1) THEN C OPERATION EN AVANCANT 1 FOIS SUR NSAUT1 POUR LE TABLEAU 2 DO 5201 IA = IDEB,IFIN IB = (IA-1) / NSAUT1 + 1 XFLO = XVAL0(IB) X2 = XVAL1(IA) IF ( ((ABS(XFLO) .LE. XPETIT) .AND. (X2 .LT. REAL(0.D0))) .OR. & (XFLO .LT. REAL(0.D0))) THEN IRETOU = 213 RETURN ELSE I2 = NINT(X2) XFLOT1 = ABS(X2 - REAL(I2)) XFLOT = ABS(X2 - REAL(0.5D0)) IF ( XFLOT1 .LE. XZPREC*ABS(X2)) THEN C PUISSANCE ENTIERE XVAL2(IA)= XFLO ** I2 ELSEIF (XFLOT .LE. XZPREC*ABS(X2) ) THEN C RACINE CARREE SQRT XVAL2(IA)= SQRT(XFLO) ELSE C CAS GENERAL XVAL2(IA)= XFLO ** X2 ENDIF ENDIF 5201 CONTINUE ELSE C OPERATION EN AVANCANT 1 FOIS SUR NSAUT2 POUR LE TABLEAU 1 DO 5301 IA = IDEB,IFIN IB = (IA-1) / NSAUT2 + 1 XFLO = XVAL0(IA) X2 = XVAL1(IB) IF ( ((ABS(XFLO) .LE. XPETIT) .AND. (X2 .LT. REAL(0.D0))) .OR. & (XFLO .LT. REAL(0.D0))) THEN IRETOU = 213 RETURN ELSE I2 = NINT(X2) XFLOT1 = ABS(X2 - REAL(I2)) XFLOT = ABS(X2 - REAL(0.5D0)) IF ( XFLOT1 .LE. XZPREC*ABS(X2)) THEN C PUISSANCE ENTIERE XVAL2(IA)= XFLO ** I2 ELSEIF (XFLOT .LE. XZPREC*ABS(X2) ) THEN C RACINE CARREE SQRT XVAL2(IA)= SQRT(XFLO) ELSE C CAS GENERAL XVAL2(IA)= XFLO ** X2 ENDIF ENDIF 5301 CONTINUE ENDIF RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C PRODUIT CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 5002 CONTINUE IF (NSAUT1 .EQ. 1 ) THEN C OPERATION TERME A TERME DO 5102 IA = IDEB,IFIN XVAL2(IA)= XVAL0(IA) * XVAL1(IA) 5102 CONTINUE ELSEIF(NSAUT1 .GT. 0 ) THEN C OPERATION EN AVANCANT 1 FOIS SUR NSAUT1 POUR LE TABLEAU 2 DO 5202 IA = IDEB,IFIN IB = (IA-1) / NSAUT1 + 1 XVAL2(IA)= XVAL0(IA) * XVAL1(IB) 5202 CONTINUE ELSE C OPERATION EN AVANCANT 1 FOIS SUR NSAUT2 POUR LE TABLEAU 1 DO 5302 IA = IDEB,IFIN IB = (IA-1) / NSAUT2 + 1 XVAL2(IA)= XVAL0(IB) * XVAL1(IA) 5302 CONTINUE ENDIF RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ADDITION CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 5003 CONTINUE IF (NSAUT1 .EQ. 1 ) THEN C OPERATION TERME A TERME DO 5103 IA=IDEB,IFIN XVAL2(IA) = XVAL0(IA) + XVAL1(IA) 5103 CONTINUE ELSEIF(NSAUT1 .GT. 0 ) THEN C OPERATION EN AVANCANT 1 FOIS SUR NSAUT1 POUR LE TABLEAU 2 DO 5203 IA=IDEB,IFIN IB = (IA-1) / NSAUT1 + 1 XVAL2(IA) = XVAL0(IA) + XVAL1(IB) 5203 CONTINUE ELSE C OPERATION EN AVANCANT 1 FOIS SUR NSAUT2 POUR LE TABLEAU 1 DO 5303 IA=IDEB,IFIN IB = (IA-1) / NSAUT2 + 1 XVAL2(IA) = XVAL0(IB) + XVAL1(IA) 5303 CONTINUE ENDIF RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C SOUSTRACTION CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 5004 CONTINUE IF (NSAUT1 .EQ. 1 ) THEN C OPERATION TERME A TERME DO 5104 IA = IDEB,IFIN XVAL2(IA)= XVAL0(IA) - XVAL1(IA) 5104 CONTINUE ELSEIF(NSAUT1 .GT. 0 ) THEN C OPERATION EN AVANCANT 1 FOIS SUR NSAUT1 POUR LE TABLEAU 2 DO 5204 IA = IDEB,IFIN IB = (IA-1) / NSAUT1 + 1 XVAL2(IA)= XVAL0(IA) - XVAL1(IB) 5204 CONTINUE ELSE C OPERATION EN AVANCANT 1 FOIS SUR NSAUT2 POUR LE TABLEAU 1 DO 5304 IA = IDEB,IFIN IB = (IA-1) / NSAUT2 + 1 XVAL2(IA)= XVAL0(IB) - XVAL1(IA) 5304 CONTINUE ENDIF RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C DIVISION CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 5005 CONTINUE IF (NSAUT1 .EQ. 1 ) THEN C OPERATION TERME A TERME DO 5105 IA = IDEB,IFIN X2 = XVAL1(IA) IF (ABS(X2) .GT. XPETIT) THEN XVAL2(IA)= XVAL0(IA) / X2 ELSE IRETOU = 835 RETURN ENDIF 5105 CONTINUE ELSEIF(NSAUT1 .GT. 0 ) THEN C OPERATION EN AVANCANT 1 FOIS SUR NSAUT1 POUR LE TABLEAU 2 DO 5205 IA = IDEB,IFIN IB = (IA-1) / NSAUT1 + 1 X2 = XVAL1(IB) IF (ABS(X2) .GT. XPETIT) THEN XVAL2(IA)= XVAL0(IA) / X2 ELSE IRETOU = 835 RETURN ENDIF 5205 CONTINUE ELSE C OPERATION EN AVANCANT 1 FOIS SUR NSAUT2 POUR LE TABLEAU 1 DO 5305 IA = IDEB,IFIN IB = (IA-1) / NSAUT2 + 1 X2 = XVAL1(IA) IF (ABS(X2) .GT. XPETIT) THEN XVAL2(IA)= XVAL0(IB) / X2 ELSE IRETOU = 835 RETURN ENDIF 5305 CONTINUE ENDIF RETURN CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ARCTANGENTE A 2 ARGUMENTS CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 5011 CONTINUE IF (NSAUT1 .EQ. 1 ) THEN C OPERATION TERME A TERME DO 5111 IA = IDEB,IFIN XVAL2(IA)= XINV * ATAN2(XVAL0(IA),XVAL1(IA)) 5111 CONTINUE ELSEIF(NSAUT1 .GT. 0 ) THEN C OPERATION EN AVANCANT 1 FOIS SUR NSAUT1 POUR LE TABLEAU 2 DO 5211 IA = IDEB,IFIN IB = (IA-1) / NSAUT1 + 1 XVAL2(IA)= XINV * ATAN2(XVAL0(IA),XVAL1(IB)) 5211 CONTINUE ELSE C OPERATION EN AVANCANT 1 FOIS SUR NSAUT2 POUR LE TABLEAU 1 DO 5311 IA = IDEB,IFIN IB = (IA-1) / NSAUT2 + 1 XVAL2(IA)= XINV * ATAN2(XVAL0(IB),XVAL1(IA)) 5311 CONTINUE ENDIF RETURN END