]>
Commit | Line | Data |
---|---|---|
fe4da5cc | 1 | * |
2 | * $Id$ | |
3 | * | |
4 | * $Log$ | |
5 | * Revision 1.1.1.1 1996/04/01 15:02:48 mclareni | |
6 | * Mathlib gen | |
7 | * | |
8 | * | |
9 | #include "gen/pilot.h" | |
10 | SUBROUTINE VVISET(RKA,BE2,MODE,XL,XU) | |
11 | ||
12 | COMMON /G116C1/ H(7),T0,T1,T,OMEGA,A(155),B(155),X0 | |
13 | EXTERNAL G116F1,G116F2 | |
14 | CHARACTER NAME*(*) | |
15 | CHARACTER*80 ERRTXT | |
16 | PARAMETER (NAME = 'VVISET') | |
17 | DIMENSION XP(8),XQ(6) | |
18 | ||
19 | C H1 = 5 ln 10 - ln(2/PI**2), h2 = 3 ln 10, h0 = Ln 0.0005... | |
20 | ||
21 | PARAMETER (PI = 3.14159 265, EU = 0.57721 566) | |
22 | PARAMETER (PI2 = 2*PI, RPI = 1/PI) | |
23 | PARAMETER (H1 = 5*2.30258 509-1.59631 259, H2 = 3*2.30258 509) | |
24 | PARAMETER (H0 =-7.6, PIH = PI/2, EPS = 1E-5) | |
25 | ||
26 | DATA XP /9.29, 2.47, 0.89, 0.36, 0.15, 0.07, 0.03, 0.02/ | |
27 | DATA XQ /0.012, 0.03, 0.08, 0.26, 0.87, 3.83/ | |
28 | ||
29 | IF(RKA .LT. 0.01 .OR. RKA .GT. 10) THEN | |
30 | WRITE(ERRTXT,101) RKA | |
31 | CALL MTLPRT(NAME,'G116.1',ERRTXT) | |
32 | ELSEIF(BE2 .LT. 0 .OR. BE2 .GT. 1) THEN | |
33 | WRITE(ERRTXT,102) BE2 | |
34 | CALL MTLPRT(NAME,'G116.2',ERRTXT) | |
35 | ELSE | |
36 | H(5)=1-BE2*(1-EU)-H0/RKA | |
37 | H(6)=BE2 | |
38 | H(7)=1-BE2 | |
39 | H4=H0/RKA-(1+BE2*EU) | |
40 | H5=LOG(RKA) | |
41 | H6=1/RKA | |
42 | T0=(H4-H(5)*H5-(H(5)+BE2)*(LOG(H(5))+REXPIN(H(5)))+EXP(-H(5)))/ | |
43 | 1 H(5) | |
44 | DO 1 LP = 1,8 | |
45 | IF(RKA .GE. XP(LP)) GO TO 11 | |
46 | 1 CONTINUE | |
47 | LP=9 | |
48 | 11 DO 2 LQ = 1,6 | |
49 | IF(RKA .LE. XQ(LQ)) GO TO 22 | |
50 | 2 CONTINUE | |
51 | LQ=7 | |
52 | 22 CALL RZERO(-LP-0.5,LQ-7.5,U,XX,EPS,1000,G116F2) | |
53 | Q=1/U | |
54 | T1=H4*Q-H5-(1+BE2*Q)*(LOG(ABS(U))+REXPIN(U))+EXP(-U)*Q | |
55 | ||
56 | T=T1-T0 | |
57 | OMEGA=PI2/T | |
58 | H(1)=RKA*(2+BE2*EU)+H1 | |
59 | IF(RKA .GE. 0.07) H(1)=H(1)+H2 | |
60 | H(2)=BE2*RKA | |
61 | H(3)=H6*OMEGA | |
62 | H(4)=PIH*OMEGA | |
63 | CALL RZERO(5.,155.,X0,XX,EPS,1000,G116F1) | |
64 | N=X0+1 | |
65 | ||
66 | D=RPI*EXP(RKA*(1+BE2*(EU-H5))) | |
67 | A(N)=0 | |
68 | IF(MODE .EQ. 0) A(N)=RPI*OMEGA | |
69 | Q=-1 | |
70 | Q2=2 | |
71 | DO 3 K = 1,N-1 | |
72 | L=N-K | |
73 | X=OMEGA*K | |
74 | X1=H6*X | |
75 | C1=LOG(X)-COSINT(X1) | |
76 | C2=SININT(X1) | |
77 | C3=SIN(X1) | |
78 | C4=COS(X1) | |
79 | XF1=RKA*(BE2*C1-C4)-X*C2 | |
80 | XF2=X*C1+RKA*(C3+BE2*C2)+T0*X | |
81 | IF(MODE .EQ. 0) THEN | |
82 | D1=Q*D*OMEGA*EXP(XF1) | |
83 | A(L)=D1*COS(XF2) | |
84 | B(L)=-D1*SIN(XF2) | |
85 | ELSE | |
86 | D1=Q*D*EXP(XF1)/K | |
87 | A(L)=D1*SIN(XF2) | |
88 | B(L)=D1*COS(XF2) | |
89 | A(N)=A(N)+Q2*A(L) | |
90 | ENDIF | |
91 | Q=-Q | |
92 | 3 Q2=-Q2 | |
93 | ENDIF | |
94 | XL=T0 | |
95 | XU=T1 | |
96 | RETURN | |
97 | 101 FORMAT('KAPPA = ',E10.3,' - OUT OF RANGE') | |
98 | 102 FORMAT('BETA**2 = ',E10.3,' - OUT OF RANGE') | |
99 | END |