COSMOS v7.655  COSMOSv7655
(AirShowerMC)
Functions/Subroutines
mneig.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine mneig (A, NDIMA, N, MITS, WORK, PRECIS, IFAULT)
 

Function/Subroutine Documentation

◆ mneig()

subroutine mneig ( dimension(ndima,*)  A,
  NDIMA,
  N,
  MITS,
dimension(*)  WORK,
  PRECIS,
  IFAULT 
)

Definition at line 10 of file mneig.f.

References a, b, c, e, f, h, i, j, m, n, o, parameter(), pt, r, and z.

Referenced by mnpsdf().

10 *
11 * $Id: d506dp.inc,v 1.1.1.1 1996/03/07 14:31:32 mclareni Exp $
12 *
13 * $Log: d506dp.inc,v $
14 * Revision 1.1.1.1 1996/03/07 14:31:32 mclareni
15 * Minuit
16 *
17 *
18 *
19 *
20 * d506dp.inc
21 *
22 C ************ DOUBLE PRECISION VERSION *************
23  IMPLICIT DOUBLE PRECISION (a-h,o-z)
24 C
25  parameter(zero=0.0, one=1.0, two=2.0)
26  parameter(tol=1.0e-35)
27  dimension a(ndima,*),work(*)
28 C PRECIS is the machine precision EPSMAC
29  ifault = 1
30 C
31  i = n
32  DO 70 i1 = 2,n
33  l = i-2
34  f = a(i,i-1)
35  gl = zero
36 C
37  IF(l .LT. 1) GO TO 25
38 C
39  DO 20 k = 1,l
40  20 gl = gl+a(i,k)**2
41  25 h = gl + f**2
42 C
43  IF(gl .GT. tol) GO TO 30
44 C
45  work(i) = zero
46  work(n+i) = f
47  GO TO 65
48  30 l = l+1
49 C
50  gl = sqrt(h)
51 C
52  IF(f .GE. zero) gl = -gl
53 C
54  work(n+i) = gl
55  h = h-f*gl
56  a(i,i-1) = f-gl
57  f = zero
58  DO 50 j = 1,l
59  a(j,i) = a(i,j)/h
60  gl = zero
61  DO 40 k = 1,j
62  40 gl = gl+a(j,k)*a(i,k)
63 C
64  IF(j .GE. l) GO TO 47
65 C
66  j1 = j+1
67  DO 45 k = j1,l
68  45 gl = gl+a(k,j)*a(i,k)
69  47 work(n+j) = gl/h
70  f = f+gl*a(j,i)
71  50 CONTINUE
72  hh = f/(h+h)
73  DO 60 j = 1,l
74  f = a(i,j)
75  gl = work(n+j)-hh*f
76  work(n+j) = gl
77  DO 60 k = 1,j
78  a(j,k) = a(j,k)-f*work(n+k)-gl*a(i,k)
79  60 CONTINUE
80  work(i) = h
81  65 i = i-1
82  70 CONTINUE
83  work(1) = zero
84  work(n+1) = zero
85  DO 110 i = 1,n
86  l = i-1
87 C
88  IF(work(i) .EQ. zero .OR. l .EQ. 0) GO TO 100
89 C
90  DO 90 j = 1,l
91  gl = zero
92  DO 80 k = 1,l
93  80 gl = gl+a(i,k)*a(k,j)
94  DO 90 k = 1,l
95  a(k,j) = a(k,j)-gl*a(k,i)
96  90 CONTINUE
97  100 work(i) = a(i,i)
98  a(i,i) = one
99 C
100  IF(l .EQ. 0) GO TO 110
101 C
102  DO 105 j = 1,l
103  a(i,j) = zero
104  a(j,i) = zero
105  105 CONTINUE
106  110 CONTINUE
107 C
108 C
109  n1 = n-1
110  DO 130 i = 2,n
111  i0 = n+i-1
112  130 work(i0) = work(i0+1)
113  work(n+n) = zero
114  b = zero
115  f = zero
116  DO 210 l = 1,n
117  j = 0
118  h = precis*(abs(work(l))+abs(work(n+l)))
119 C
120  IF(b .LT. h) b = h
121 C
122  DO 140 m1 = l,n
123  m = m1
124 C
125  IF(abs(work(n+m)) .LE. b) GO TO 150
126 C
127  140 CONTINUE
128 C
129  150 IF(m .EQ. l) GO TO 205
130 C
131  160 IF(j .EQ. mits) RETURN
132 C
133  j = j+1
134  pt = (work(l+1)-work(l))/(two*work(n+l))
135  r = sqrt(pt*pt+one)
136  pr = pt+r
137 C
138  IF(pt .LT. zero) pr=pt-r
139 C
140  h = work(l)-work(n+l)/pr
141  DO 170 i=l,n
142  170 work(i) = work(i)-h
143  f = f+h
144  pt = work(m)
145  c = one
146  s = zero
147  m1 = m-1
148  i = m
149  DO 200 i1 = l,m1
150  j = i
151  i = i-1
152  gl = c*work(n+i)
153  h = c*pt
154 C
155  IF(abs(pt) .GE. abs(work(n+i))) GO TO 180
156 C
157  c = pt/work(n+i)
158  r = sqrt(c*c+one)
159  work(n+j) = s*work(n+i)*r
160  s = one/r
161  c = c/r
162  GO TO 190
163  180 c = work(n+i)/pt
164  r = sqrt(c*c+one)
165  work(n+j) = s*pt*r
166  s = c/r
167  c = one/r
168  190 pt = c*work(i)-s*gl
169  work(j) = h+s*(c*gl+s*work(i))
170  DO 200 k = 1,n
171  h = a(k,j)
172  a(k,j) = s*a(k,i)+c*h
173  a(k,i) = c*a(k,i)-s*h
174  200 CONTINUE
175  work(n+l) = s*pt
176  work(l) = c*pt
177 C
178  IF(abs(work(n+l)) .GT. b) GO TO 160
179 C
180  205 work(l) = work(l)+f
181  210 CONTINUE
182  DO 240 i=1,n1
183  k = i
184  pt = work(i)
185  i1 = i+1
186  DO 220 j = i1,n
187 C
188  IF(work(j) .GE. pt) GO TO 220
189 C
190  k = j
191  pt = work(j)
192  220 CONTINUE
193 C
194  IF(k .EQ. i) GO TO 240
195 C
196  work(k) = work(i)
197  work(i) = pt
198  DO 230 j=1,n
199  pt = a(j,i)
200  a(j,i) = a(j,k)
201  a(j,k) = pt
202  230 CONTINUE
203  240 CONTINUE
204  ifault = 0
205 C
206  RETURN
parameter
integer npitbl real *nx parameter(n=101, npitbl=46, nx=n-1) real *8 uconst
z
nodes z
Definition: cstdatmos0.fNew2.h:93
e
dE dx *! Nuc Int sampling table e
Definition: cblkMuInt.h:130
r
*Zfirst p fm *Zfirst p Zfirst p Zfirst p *Zfirst p *Zfirst pos xyz r
Definition: ZavoidUnionMap.h:1
i
nodes i
Definition: cstdatmos0.fNew2.h:93
pt
! structure defining a particle at production ! Basic idea of what is to be contained in ! the particle structue is that dynamical ones should be included those derivable from the particle code ! is not included ******************************************************type fmom momentum sequence union map real e endmap map real t endmap map ! pt before pz is set real pt
Definition: Zptcl.h:21
nrl_atmos::h
real(4), dimension(:), allocatable, save h
Definition: cNRLAtmos.f:28
j
********************block data cblkHeavy ********************integer j data *HeavyG2symbol *data *HeavyG2code kiron data j
Definition: cblkHeavy.h:36
o
struct ob o[NpMax]
Definition: Zprivate.h:34
m
block data include Zlatfit h c fitting region data data data data data d0 data data d0 data data m
Definition: ZlatfitBD.h:35
nrl_atmos::a
real(4), save a
Definition: cNRLAtmos.f:20
nrl_atmos::b
real(4), save b
Definition: cNRLAtmos.f:21
n
integer n
Definition: Zcinippxc.h:1
f
dE dx *! Nuc Int sampling table f
Definition: cblkMuInt.h:130
c
dE dx *! Nuc Int sampling table c
Definition: cblkMuInt.h:130
Here is the call graph for this function:
Here is the caller graph for this function: