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

Go to the source code of this file.

Functions/Subroutines

subroutine mnimpr (FCN, FUTIL)
 

Function/Subroutine Documentation

◆ mnimpr()

subroutine mnimpr ( external  FCN,
external  FUTIL 
)

Definition at line 10 of file mnimpr.f.

References a, b, c, d, e, f, g, h, softenpik::half, i, j, m, maxp, mnamin(), mncalf(), mninex(), mnprin(), mnrazz(), mnrn15(), mnrset(), mnsimp(), mnvert(), n, o, p, parameter(), to, true, up(), x, xs, y, and z.

Referenced by mnexcm().

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 CC Attempts to improve on a good local minimum by finding a
25 CC better one. The quadratic part of FCN is removed by MNCALF
26 CC and this transformed function is minimized using the simplex
27 CC method from several random starting points.
28 CC ref. -- Goldstein and Price, Math.Comp. 25, 569 (1971)
29 CC
30 *
31 * $Id: d506cm.inc,v 1.1.1.1 1996/03/07 14:31:32 mclareni Exp $
32 *
33 * $Log: d506cm.inc,v $
34 * Revision 1.1.1.1 1996/03/07 14:31:32 mclareni
35 * Minuit
36 *
37 *
38 *
39 *
40 * d506cm.inc
41 *
42  parameter(mne=100 , mni=50)
43  parameter(mnihl=mni*(mni+1)/2)
44  CHARACTER*10 cpnam
45  COMMON
46  1/mn7nam/ cpnam(mne)
47  2/mn7ext/ u(mne) ,alim(mne) ,blim(mne)
48  3/mn7err/ erp(mni) ,ern(mni) ,werr(mni) ,globcc(mni)
49  4/mn7inx/ nvarl(mne) ,niofex(mne),nexofi(mni)
50  5/mn7int/ x(mni) ,xt(mni) ,dirin(mni)
51  6/mn7fx2/ xs(mni) ,xts(mni) ,dirins(mni)
52  7/mn7der/ grd(mni) ,g2(mni) ,gstep(mni) ,gin(mne) ,dgrd(mni)
53  8/mn7fx3/ grds(mni) ,g2s(mni) ,gsteps(mni)
54  9/mn7fx1/ ipfix(mni) ,npfix
55  a/mn7var/ vhmat(mnihl)
56  b/mn7vat/ vthmat(mnihl)
57  c/mn7sim/ p(mni,mni+1),pstar(mni),pstst(mni) ,pbar(mni),prho(mni)
58 C
59  parameter(maxdbg=10, maxstk=10, maxcwd=20, maxp=30, maxcpt=101)
60  parameter(zero=0.0, one=1.0, half=0.5)
61  COMMON
62  d/mn7npr/ maxint ,npar ,maxext ,nu
63  e/mn7iou/ isysrd ,isyswr ,isyssa ,npagwd ,npagln ,newpag
64  e/mn7io2/ istkrd(maxstk) ,nstkrd ,istkwr(maxstk) ,nstkwr
65  f/mn7tit/ cfrom ,cstatu ,ctitl ,cword ,cundef ,cvrsn ,covmes
66  g/mn7flg/ isw(7) ,idbg(0:maxdbg) ,nblock ,icomnd
67  h/mn7min/ amin ,up ,edm ,fval3 ,epsi ,apsi ,dcovar
68  i/mn7cnv/ nfcn ,nfcnmx ,nfcnlc ,nfcnfr ,itaur,istrat,nwrmes(2)
69  j/mn7arg/ word7(maxp)
70  k/mn7log/ lwarn ,lrepor ,limset ,lnolim ,lnewmn ,lphead
71  l/mn7cns/ epsmac ,epsma2 ,vlimlo ,vlimhi ,undefi ,bigedm,updflt
72  m/mn7rpt/ xpt(maxcpt) ,ypt(maxcpt)
73  n/mn7cpt/ chpt(maxcpt)
74  o/mn7xcr/ xmidcr ,ymidcr ,xdircr ,ydircr ,ke1cr ,ke2cr
75  CHARACTER ctitl*50, cword*(maxcwd), cundef*10, cfrom*8,
76  + cvrsn*6, covmes(0:3)*22, cstatu*10, chpt*1
77  LOGICAL lwarn, lrepor, limset, lnolim, lnewmn, lphead
78  EXTERNAL fcn,futil
79  dimension dsav(mni), y(mni+1)
80  parameter(alpha=1.,beta=0.5,gamma=2.0)
81  DATA rnum/0./
82  IF (npar .LE. 0) RETURN
83  IF (amin .EQ. undefi) CALL mnamin(fcn,futil)
84  cstatu = 'UNCHANGED '
85  itaur = 1
86  epsi = 0.1*up
87  npfn=nfcn
88  nloop = word7(2)
89  IF (nloop .LE. 0) nloop = npar + 4
90  nparx = npar
91  nparp1=npar+1
92  wg = 1.0/float(npar)
93  sigsav = edm
94  apsi = amin
95  DO 2 i= 1, npar
96  xt(i) = x(i)
97  dsav(i) = werr(i)
98  DO 2 j = 1, i
99  ndex = i*(i-1)/2 + j
100  p(i,j) = vhmat(ndex)
101  2 p(j,i) = p(i,j)
102  CALL mnvert(p,maxint,maxint,npar,ifail)
103  IF (ifail .GE. 1) GO TO 280
104 C Save inverted matrix in VT
105  DO 12 i= 1, npar
106  ndex = i*(i-1)/2
107  DO 12 j= 1, i
108  ndex = ndex + 1
109  12 vthmat(ndex) = p(i,j)
110  loop = 0
111 C
112  20 CONTINUE
113  DO 25 i= 1, npar
114  dirin(i) = 2.0*dsav(i)
115  CALL mnrn15(rnum,iseed)
116  25 x(i) = xt(i) + 2.0*dirin(i)*(rnum-0.5)
117  loop = loop + 1
118  reg = 2.0
119  IF (isw(5) .GE. 0) WRITE (isyswr, 1040) loop
120  30 CALL mncalf(fcn,x,ycalf,futil)
121  amin = ycalf
122 C . . . . set up random simplex
123  jl = nparp1
124  jh = nparp1
125  y(nparp1) = amin
126  amax = amin
127  DO 45 i= 1, npar
128  xi = x(i)
129  CALL mnrn15(rnum,iseed)
130  x(i) = xi - dirin(i) *(rnum-0.5)
131  CALL mncalf(fcn,x,ycalf,futil)
132  y(i) = ycalf
133  IF (y(i) .LT. amin) THEN
134  amin = y(i)
135  jl = i
136  ELSE IF (y(i) .GT. amax) THEN
137  amax = y(i)
138  jh = i
139  ENDIF
140  DO 40 j= 1, npar
141  40 p(j,i) = x(j)
142  p(i,nparp1) = xi
143  x(i) = xi
144  45 CONTINUE
145 C
146  edm = amin
147  sig2 = edm
148 C . . . . . . . start main loop
149  50 CONTINUE
150  IF (amin .LT. zero) GO TO 95
151  IF (isw(2) .LE. 2) GO TO 280
152  ep = 0.1*amin
153  IF (sig2 .LT. ep .AND. edm.LT.ep ) GO TO 100
154  sig2 = edm
155  IF ((nfcn-npfn) .GT. nfcnmx) GO TO 300
156 C calculate new point * by reflection
157  DO 60 i= 1, npar
158  pb = 0.
159  DO 59 j= 1, nparp1
160  59 pb = pb + wg * p(i,j)
161  pbar(i) = pb - wg * p(i,jh)
162  60 pstar(i)=(1.+alpha)*pbar(i)-alpha*p(i,jh)
163  CALL mncalf(fcn,pstar,ycalf,futil)
164  ystar = ycalf
165  IF(ystar.GE.amin) GO TO 70
166 C point * better than jl, calculate new point **
167  DO 61 i=1,npar
168  61 pstst(i)=gamma*pstar(i)+(1.-gamma)*pbar(i)
169  CALL mncalf(fcn,pstst,ycalf,futil)
170  ystst = ycalf
171  66 IF (ystst .LT. y(jl)) GO TO 67
172  CALL mnrazz(ystar,pstar,y,jh,jl)
173  GO TO 50
174  67 CALL mnrazz(ystst,pstst,y,jh,jl)
175  GO TO 50
176 C point * is not as good as jl
177  70 IF (ystar .GE. y(jh)) GO TO 73
178  jhold = jh
179  CALL mnrazz(ystar,pstar,y,jh,jl)
180  IF (jhold .NE. jh) GO TO 50
181 C calculate new point **
182  73 DO 74 i=1,npar
183  74 pstst(i)=beta*p(i,jh)+(1.-beta)*pbar(i)
184  CALL mncalf(fcn,pstst,ycalf,futil)
185  ystst = ycalf
186  IF(ystst.GT.y(jh)) GO TO 30
187 C point ** is better than jh
188  IF (ystst .LT. amin) GO TO 67
189  CALL mnrazz(ystst,pstst,y,jh,jl)
190  GO TO 50
191 C . . . . . . end main loop
192  95 IF (isw(5) .GE. 0) WRITE (isyswr,1000)
193  reg = 0.1
194 C . . . . . ask if point is new
195  100 CALL mninex(x)
196  CALL fcn(nparx,gin,amin,u,4,futil)
197  nfcn = nfcn + 1
198  DO 120 i= 1, npar
199  dirin(i) = reg*dsav(i)
200  IF (abs(x(i)-xt(i)) .GT. dirin(i)) GO TO 150
201  120 CONTINUE
202  GO TO 230
203  150 nfcnmx = nfcnmx + npfn - nfcn
204  npfn = nfcn
205  CALL mnsimp(fcn,futil)
206  IF (amin .GE. apsi) GO TO 325
207  DO 220 i= 1, npar
208  dirin(i) = 0.1 *dsav(i)
209  IF (abs(x(i)-xt(i)) .GT. dirin(i)) GO TO 250
210  220 CONTINUE
211  230 IF (amin .LT. apsi) GO TO 350
212  GO TO 325
213 C . . . . . . truly new minimum
214  250 lnewmn = .true.
215  IF (isw(2) .GE. 1) THEN
216  isw(2) = 1
217  dcovar = max(dcovar,half)
218  ELSE
219  dcovar = 1.
220  ENDIF
221  itaur = 0
222  nfcnmx = nfcnmx + npfn - nfcn
223  cstatu = 'NEW MINIMU'
224  IF (isw(5) .GE. 0) WRITE (isyswr,1030)
225  RETURN
226 C . . . return to previous region
227  280 IF (isw(5) .GT. 0) WRITE (isyswr,1020)
228  GO TO 325
229  300 isw(1) = 1
230  325 DO 330 i= 1, npar
231  dirin(i) = 0.01*dsav(i)
232  330 x(i) = xt(i)
233  amin = apsi
234  edm = sigsav
235  350 CALL mninex(x)
236  IF (isw(5) .GT. 0) WRITE (isyswr,1010)
237  cstatu= 'UNCHANGED '
238  CALL mnrset(0)
239  IF (isw(2) .LT. 2) GO TO 380
240  IF (loop .LT. nloop .AND. isw(1) .LT. 1) GO TO 20
241  380 CALL mnprin (5,amin)
242  itaur = 0
243  RETURN
244  1000 FORMAT (54h an improvement on the previous minimum has been found)
245  1010 FORMAT (51h improve has returned to region of original minimum)
246  1020 FORMAT (/44h covariance matrix was not positive-definite)
247  1030 FORMAT (/38h improve has found a truly new minimum/1h ,37(1h*)/)
248  1040 FORMAT (/18h start attempt no.,i2, 20h to find new minimum)
parameter
integer npitbl real *nx parameter(n=101, npitbl=46, nx=n-1) real *8 uconst
z
nodes z
Definition: cstdatmos0.fNew2.h:93
up
block data include Zlatfit h c fitting region data data data data data d0 data data d0 data data h g *is for param c g data up(2, 1)/7.0d0/
e
dE dx *! Nuc Int sampling table e
Definition: cblkMuInt.h:130
i
nodes i
Definition: cstdatmos0.fNew2.h:93
mnrset
subroutine mnrset(IOPT)
Definition: mnrset.f:10
true
block data cblkElemag data *AnihiE ! Eposi< 1 TeV, anihilation considered *X0/365.667/, ! radiation length of air in kg/m2 *Ecrit/81.e-3/, ! critical energy of air in GeV *MaxComptonE/1./, ! compton is considered below 1 GeV *MaxPhotoE/1.e-3/, ! above this, PhotoElectric effect neg. *MinPhotoProdE/153.e-3/, ! below 153 MeV, no gp --> hadrons ! scattering const not MeV *Knockon true
Definition: cblkElemag.h:7
mnvert
subroutine mnvert(A, L, M, N, IFAIL)
Definition: mnvert.f:25
nrl_atmos::h
real(4), dimension(:), allocatable, save h
Definition: cNRLAtmos.f:28
mnprin
subroutine mnprin(INKODE, FVAL)
Definition: mnprin.f:10
xs
! constants thru Cosmos real ! if multiplied to deg radian Torad ! light velocity m sec ! infinty ! kg m2 *Togpcm2 g cm2 ! g cm2 *Tokgpm2 kg m2 ! cm *Tom m ! m *Tocm cm ! g cm3 *Tokgpm3 kg m3 ! kg m3 *Togpcm3 g cm3 ! sec *Tonsec nsec ! Tesla m ! Avogadro *A2deninv ! mfp *n * xs
Definition: Zglobalc.h:18
mnamin
subroutine mnamin(FCN, FUTIL)
Definition: mnamin.f:10
p
********************block data cblkHeavy ********************integer j data *HeavyG2symbol p
Definition: cblkHeavy.h:7
j
********************block data cblkHeavy ********************integer j data *HeavyG2symbol *data *HeavyG2code kiron data j
Definition: cblkHeavy.h:36
to
integer maxbin nregion c minsave drx2 ! drawing region real maxdep integer maxpos integer to
Definition: Zfit.h:15
mncalf
subroutine mncalf(FCN, PVEC, YCALF, FUTIL)
Definition: mncalf.f:10
o
struct ob o[NpMax]
Definition: Zprivate.h:34
mninex
subroutine mninex(PINT)
Definition: mninex.f:10
y
latitude latitude this system is used *****************************************************************! type coord sequence union map real y
Definition: Zcoord.h:25
d
dE dx *! Nuc Int sampling table d
Definition: cblkMuInt.h:130
mnsimp
subroutine mnsimp(FCN, FUTIL)
Definition: mnsimp.f:25
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
g
dE dx *! Nuc Int sampling table g
Definition: cblkMuInt.h:130
mnrazz
subroutine mnrazz(YNEW, PNEW, Y, JH, JL)
Definition: mnrazz.f:10
nrl_atmos::b
real(4), save b
Definition: cNRLAtmos.f:21
maxp
!onst int maxp
Definition: Zprivate.h:3
mnrn15
subroutine mnrn15(VAL, INSEED)
Definition: mnrn15.f:10
n
integer n
Definition: Zcinippxc.h:1
softenpik::half
integer, parameter half
Definition: csoftenPiK.f:108
x
! 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 * x
Definition: Zptcl.h:21
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: