]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ISAJET/isasusy/sugra.F
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Adding MUON HLT code to the repository.
[u/mrichter/AliRoot.git] / ISAJET / isasusy / sugra.F
CommitLineData
0795afa3 1#include "isajet/pilot.h"
2C--------------------------------------------------------------------
3 SUBROUTINE SUGRA(M0,MHF,A0,TANB,SGNMU,MT,IMODEL)
4C--------------------------------------------------------------------
5C
6C Calculate supergravity spectra for ISAJET using as inputs
7C M0 = M_0 = common scalar mass at GUT scale
8C MHF = M_(1/2) = common gaugino mass at GUT scale
9C A0 = A_0 = trilinear soft breaking parameter at GUT scale
10C TANB = tan(beta) = ratio of vacuum expectation values v_1/v_2
11C SGNMU = sgn(mu) = +-1 = sign of Higgsino mass term
12C MT = M_t = mass of t quark
13C M0 = Lambda = ratio of vevs <F>/<S>
14C MHF = M_Mes = messenger scale
15C A0 = n_5 = number of messenger fields
16C IMODEL = 1 for SUGRA model
17C = 2 for GMSB model
18C = 7 for AMSB model
19C
20C Uses Runge-Kutta method to integrate RGE's from M_Z to M_GUT
21C and back, putting in correct thresholds. For the first iteration
22C only the first 6 couplings are included and a common threshold
23C is used.
24C
25C See /SUGMG/ for definitions of couplings and masses.
26C
27#if defined(CERNLIB_IMPNONE)
28 IMPLICIT NONE
29#endif
30#include "isajet/sslun.inc"
31#include "isajet/sspar.inc"
32#include "isajet/sssm.inc"
33#include "isajet/sugxin.inc"
34#include "isajet/sugmg.inc"
35#include "isajet/sugpas.inc"
36#include "isajet/sugnu.inc"
37#include "isajet/ssinf.inc"
38 REAL GY(7),W1(21),G(29),W2(87)
39 REAL G0(29)
40 COMPLEX*16 SSB0,SSB1
41 DOUBLE PRECISION DDILOG,XLM
42 INTEGER IG(29)
43 EXTERNAL SURG06,SURG26
44 REAL M0,MHF,A0,TANB,SGNMU,MT,XLAMGM,XMESGM,XN5GM
45 INTEGER NSTEP
46 REAL M2,SUALFE,SUALFS,Q,T,A1I,AGUT,A3I,A2I,MTMT,ASMT,DT,
47 $TGUT,TZ,GGUT,SIG2,SIG1,MH1S,MH2S,AGUTI,
48 $MUS,MBMZ,MB,MTAU,MZ,MW,SR2,PI,ALEM,MTAMZ,
49 $MTAMB,MTAMTA,MBMB,ASMB,BETA,COTB,SINB,COS2B,COSB,XC,
50 $MSN,MG,MT1,MT2,MB1,MB2,MW1,MW2,AMU,BTHAT,BBHAT,BLHAT,AM2
51 INTEGER II,I,J,IMODEL
52 REAL G0SAVE(26),DELG0,DELLIM,THRF,THRG,DY,QOLD
53 INTEGER MXITER,NSTEP0
54 COMPLEX*16 ZZZ
55 REAL*8 REAL8
56C
57 DATA MZ/91.187/,MTAU/1.777/,MB/4.9/,ALEM/.0078186/
58C This choice is a compromise between precision and speed:
59 DATA MXITER/20/,NSTEP0/200/,DELLIM/2.E-2/
60C
61C Define REAL(COMPLEX*16) for g77. This might need to be
62C changed for 64-bit machines?
63C
64 REAL8(ZZZ)=DREAL(ZZZ)
65C
66C Save input parameters
67C
68 XSUGIN(1)=M0
69 XSUGIN(2)=MHF
70 XSUGIN(3)=A0
71 XSUGIN(4)=TANB
72 XSUGIN(5)=SGNMU
73 XSUGIN(6)=MT
74 XLAMGM=M0
75 XMESGM=MHF
76 XN5GM=A0
77 XGMIN(1)=XLAMGM
78 XGMIN(2)=XMESGM
79 XGMIN(3)=XN5GM
80 XGMIN(4)=TANB
81 XGMIN(5)=SGNMU
82 XGMIN(6)=MT
83 IF (XGMIN(12).EQ.0.) XGMIN(12)=XN5GM
84 IF (XGMIN(13).EQ.0.) XGMIN(13)=XN5GM
85 IF (XGMIN(14).EQ.0.) XGMIN(14)=XN5GM
86C
87C Compute gauge mediated threshold functions
88C
89 IF (IMODEL.EQ.2) THEN
90 XLM=XLAMGM/XMESGM
91 THRF=((1.D0+XLM)*(LOG(1.D0+XLM)-2*DDILOG(XLM/(1.D0+XLM))+
92 , .5*DDILOG(2*XLM/(1.D0+XLM)))+
93 , (1.D0-XLM)*(LOG(1.D0-XLM)-2*DDILOG(-XLM/(1.D0-XLM))+
94 , .5*DDILOG(-2*XLM/(1.D0-XLM))))/XLM**2
95 THRG=((1.D0+XLM)*LOG(1.D0+XLM)+(1.D0-XLM)*LOG(1.D0-XLM))/XLM**2
96 END IF
97C
98C Initialize standard model parameters in /SSSM/:
99C
100 AMUP=0.0056
101 AMDN=0.0099
102 AMST=0.199
103 AMCH=1.35
104 AMBT=5.0
105 AMTP=MT
106 AMT=MT
107 AME=0.511E-3
108 AMMU=0.105
109 AMTAU=1.777
110 AMZ=91.17
111 GAMW=2.12
112 GAMZ=2.487
113 ALFAEM=1./128.
114 SN2THW=0.232
115 ALFA2=ALFAEM/SN2THW
116 ALQCD4=0.177
117 ALFA3=0.118
118C
119 NOGOOD=0
120 ITACHY=0
121 PI=4.*ATAN(1.)
122 SR2=SQRT(2.)
123 XW=.2324-1.03E-7*(MT**2-138.**2)
124 MW=MZ*SQRT(1.-XW)
125 AMW=MW
126 A1MZ=5*ALEM/3./(1.-XW)
127 A2MZ=ALEM/XW
128 G2=SQRT(4*PI*A2MZ)
129 GP=SQRT(3./5.*A1MZ*4.*PI)
130 XTANB=TANB
131 COTB=1./TANB
132 BETA=ATAN(TANB)
133 SINB=SIN(BETA)
134 COSB=COS(BETA)
135 SIN2B=SIN(2*BETA)
136 COS2B=COS(2*BETA)
137 IF (IMODEL.EQ.1) THEN
138 MSUSY=SQRT(M0**2+4*MHF**2)
139 ELSE IF (IMODEL.EQ.2) THEN
140 MSUSY=XLAMGM/100.
141 ELSE IF (IMODEL.EQ.7) THEN
142 MSUSY=SQRT(M0**2+(.01*MHF)**2)
143 END IF
144C USE PIERCE PRESCRIPTION FOR MAGNITUDE OF VEV
145C VEV=SR2*(248.6+0.9*LOG(MSUSY/AMZ)
146C V=SQRT(VEV**2/(1.+COTB))
147C PREVIOUS PRESCRIPTION
148 V=SQRT(2*MW**2/G2**2/(1.+COTB**2))
149 VP=V/TANB
150 VEV=SQRT(V**2+VP**2)
151C
152C Compute m(tau), m(b) at z scale using qcd, qed
153C
154 MTAMTA=MTAU*(1.-SUALFE(MTAU**2)/PI)
155 MTAMB=MTAMTA*(SUALFE(MB**2)/SUALFE(MTAU**2))**(-27./76.)
156 MTAMZ=MTAMB*(SUALFE(MZ**2)/SUALFE(MB**2))**(-27./80.)
157 FTAMZ=MTAMZ/COSB/VEV
158 ASMB=SUALFS(MB**2,.36,MT,3)
159 MBMB=MB*(1.-4*ASMB/3./PI)
160 ASMZ=SUALFS(MZ**2,.36,MT,3)
161 MBMZ=MBMB*(ASMZ/ASMB)**(12./23.)*
162 $ (SUALFE(MZ**2)/SUALFE(MB**2))**(-3./80.)
163 FBMZ=MBMZ/COSB/VEV
164 ASMT=SUALFS(MT**2,.36,MT,3)
165 MTMT=MT/(1.+4*ASMT/3./PI+(16.11-1.04*(5.-6.63/MT))*(ASMT/PI)**2)
166 FTMT=MTMT/SINB/VEV
167 FNMZ=SQRT(XNRIN(2)*XNRIN(1)/(SINB*VEV)**2)
168 AMNRMJ=XNRIN(2)
169C
170C Run the 3 gauge and 3 Yukawa's up to find M_GUT ,A_GUT and
171C Yukawa_GUT
172C
173C
174 NSTEP=NSTEP0
175 GY(1)=SQRT(4*PI*A1MZ)
176 GY(2)=SQRT(4*PI*A2MZ)
177 GY(3)=SQRT(4*PI*ALFA3)
178 GY(4)=FTAMZ
179 GY(5)=FBMZ
180 GY(6)=0.
181 GY(7)=0.
182 IF (IMODEL.EQ.1.OR.IMODEL.EQ.7) THEN
183 IF (XSUGIN(7).EQ.0.) THEN
184 MGUT=1.E19
185 ELSE
186 MGUT=XSUGIN(7)
187 END IF
188 ELSE IF (IMODEL.EQ.2) THEN
189 MGUT=XMESGM
190 END IF
191 TZ=LOG(MZ/MGUT)
192 TGUT=0.
193 DT=(TGUT-TZ)/FLOAT(NSTEP)
194 DO 200 II=1,NSTEP
195 T=TZ+(TGUT-TZ)*FLOAT(II-1)/FLOAT(NSTEP)
196 Q=MGUT*EXP(T)
197 IF (Q.GT.MT.AND.GY(6).EQ.0.) GY(6)=FTMT
198 IF (Q.GT.XNRIN(2).AND.GY(7).EQ.0.) GY(7)=FNMZ
199 CALL RKSTP(7,DT,T,GY,SURG06,W1)
200 A1I=4*PI/GY(1)**2
201 A2I=4*PI/GY(2)**2
202 A3I=4*PI/GY(3)**2
203 IF (GY(5).GT.10..OR.GY(6).GT.10..OR.GY(7).GT.10.) THEN
204 NOGOOD=4
205 GO TO 100
206 END IF
207 IF (A1I.LT.A2I.AND.XSUGIN(7).EQ.0.) GO TO 10
208200 CONTINUE
209 IF (MGUT.EQ.1.E19) THEN
210 WRITE(LOUT,*) 'SUGRA: NO UNIFICATION FOUND'
211 GO TO 100
212 END IF
21310 IF (XSUGIN(7).EQ.0.) THEN
214 MGUT=Q
215 ELSE
216 MGUT=XSUGIN(7)
217 END IF
218 AGUT=(GY(1)**2/4./PI+GY(2)**2/4./PI)/2.
219 GGUT=SQRT(4*PI*AGUT)
220 AGUTI=1./AGUT
221 FTAGUT=GY(4)
222 FBGUT=GY(5)
223 FTGUT=GY(6)
224 IF (XNRIN(1).EQ.0..AND.XNRIN(2).LT.1.E19) THEN
225C UNIFY FN-FT
226 FNGUT=GY(6)
227 ELSE
228 FNGUT=GY(7)
229 END IF
230C
231C Define parameters at GUT scale
232C
233 DO 210 J=1,3
234 IF (IMODEL.EQ.1) THEN
235 G(J)=GY(J)
236 G(J+6)=MHF
237 G(J+9)=A0
238 ELSE IF (IMODEL.EQ.2) THEN
239 G(J)=GY(J)
240 G(J+6)=XGMIN(11+J)*XGMIN(8)*THRG*(GY(J)/4./PI)**2*XLAMGM
241 G(J+9)=0.
242 END IF
243210 CONTINUE
244C OVERWRITE ALFA_3 UNIFICATION TO GET ALFA_3(MZ) RIGHT
245 IF (IMODEL.EQ.1.AND.IAL3UN.NE.0) G(3)=GGUT
246 G(4)=FTAGUT
247 G(5)=FBGUT
248 G(6)=FTGUT
249C IF NR MAJORANA MASS EXISTS, SET EXTRA NR RGE PARAMETERS
250 IF (XNRIN(2).LT.1.E19) THEN
251 G(27)=FNGUT
252 G(28)=XNRIN(4)**2
253 G(29)=XNRIN(3)
254 ELSE
255 G(27)=0.
256 G(28)=0.
257 G(29)=0.
258 END IF
259 IF (IMODEL.EQ.1) THEN
260 DO 220 J=13,24
261 G(J)=M0**2
262220 CONTINUE
263C Set possible non-universal boundary conditions
264 DO 230 J=1,6
265 IF (XNUSUG(J).LT.1.E19) THEN
266 G(J+6)=XNUSUG(J)
267 END IF
268230 CONTINUE
269 DO 231 J=7,18
270 IF (XNUSUG(J).LT.1.E19) THEN
271 G(J+6)=XNUSUG(J)**2
272 END IF
273231 CONTINUE
274 ELSE IF (IMODEL.EQ.2) THEN
275 XC=2*THRF*XLAMGM**2
276 DY=SQRT(3./5.)*GY(1)*XGMIN(11)
277 G(13)=XC*(.75*XGMIN(13)*(GY(2)/4./PI)**4+.6*.25*
278 ,XGMIN(12)*(GY(1)/4./PI)**4)+XGMIN(9)-DY
279 G(14)=XC*(.75*XGMIN(13)*(GY(2)/4./PI)**4+.6*.25*
280 ,XGMIN(12)*(GY(1)/4./PI)**4)+XGMIN(10)+DY
281 G(15)=XC*(.6*XGMIN(12)*(GY(1)/4./PI)**4)+2*DY
282 G(16)=XC*(.75*XGMIN(13)*(GY(2)/4./PI)**4+.6*.25*
283 ,XGMIN(12)*(GY(1)/4./PI)**4)-DY
284 G(17)=XC*(4*XGMIN(14)*(GY(3)/4./PI)**4/3.+.6*XGMIN(12)*
285 ,(GY(1)/4./PI)**4/9.)+2*DY/3.
286 G(18)=XC*(4*XGMIN(14)*(GY(3)/4./PI)**4/3.+.6*4*XGMIN(12)*
287 ,(GY(1)/4./PI)**4/9.)-4*DY/3.
288 G(19)=XC*(4*XGMIN(14)*(GY(3)/4./PI)**4/3.+.75*XGMIN(13)*
289 ,(GY(2)/4./PI)**4+.6*XGMIN(12)*(GY(1)/4./PI)**4/36.)+DY/3.
290 G(20)=G(15)
291 G(21)=G(16)
292 G(22)=G(17)
293 G(23)=G(18)
294 G(24)=G(19)
295 ELSE IF (IMODEL.EQ.7) THEN
296 G(1)=GY(1)
297 G(2)=GY(2)
298 G(3)=GY(3)
299 BLHAT=G(4)*(-9*G(1)**2/5.-3*G(2)**2+3*G(5)**2+4*G(4)**2)
300 BBHAT=G(5)*(-7*G(1)**2/15.-3*G(2)**2-16*G(3)**2/3.+
301 , G(6)**2+6*G(5)**2+G(4)**2)
302 BTHAT=G(6)*(-13*G(1)**2/15.-3*G(2)**2-16*G(3)**2/3.+
303 , 6*G(6)**2+G(5)**2)
304 G(7)=-33*MHF*G(1)**2/5./16./PI**2
305 G(8)=-MHF*G(2)**2/16./PI**2
306 G(9)=3*MHF*G(3)**2/16./PI**2
307 G(10)=BLHAT*MHF/G(4)/16./PI**2
308 G(11)=BBHAT*MHF/G(5)/16./PI**2
309 G(12)=BTHAT*MHF/G(6)/16./PI**2
310 G(13)=(-99*G(1)**4/50.-3*G(2)**4/2.+3*G(5)*BBHAT+G(4)*BLHAT)*
311 , MHF**2/(16*PI**2)**2
312 G(14)=(-99*G(1)**4/50.-3*G(2)**4/2.+3*G(6)*BTHAT)*
313 , MHF**2/(16*PI**2)**2
314 G(15)=(-198*G(1)**4/25.)*MHF**2/(16*PI**2)**2
315 G(16)=(-99*G(1)**4/50.-3*G(2)**4/2.)*MHF**2/(16*PI**2)**2
316 G(17)=(-22*G(1)**4/25.+8*G(3)**4)*MHF**2/(16*PI**2)**2
317 G(18)=(-88*G(1)**4/25.+8*G(3)**4)*MHF**2/(16*PI**2)**2
318 G(19)=(-11*G(1)**4/50.-3*G(2)**4/2.+8*G(3)**4)*
319 , MHF**2/(16*PI**2)**2
320 G(20)=(-198*G(1)**4/25.+2*G(4)*BLHAT)*MHF**2/(16*PI**2)**2
321 G(21)=(-99*G(1)**4/50.-3*G(2)**4/2.+G(4)*BLHAT)*
322 , MHF**2/(16*PI**2)**2
323 G(22)=(-22*G(1)**4/25.+8*G(3)**4+2*G(5)*BBHAT)*
324 ,MHF**2/(16*PI**2)**2
325 G(23)=(-88*G(1)**4/25.+8*G(3)**4+2*G(6)*BTHAT)*
326 ,MHF**2/(16*PI**2)**2
327 G(24)=(-11*G(1)**4/50.-3*G(2)**4/2.+8*G(3)**4+G(5)*BBHAT+
328 , G(6)*BTHAT)*MHF**2/(16*PI**2)**2
329 DO 234 I=13,24
330234 G(I)=G(I)+M0**2
331 END IF
332 G(25)=0.
333 G(26)=0.
334 DO 235 I=1,29
335 IG(I)=0
336235 CONTINUE
337C Check for tachyonic sleptons at GUT scale
338 IF (G(15).LT.0..OR.G(16).LT.0.) THEN
339 ITACHY=1
340 END IF
341C
342C Initialize thresholds
343C
344 MSS(1)=MSUSY
345 MSS(2)=MSUSY
346 MSS(17)=MSUSY
347 MSS(27)=MSUSY
348 MSS(31)=MSUSY
349 MU=MSUSY
350C
351C Evolve parameters from mgut to mz
352C
353 TZ=LOG(MZ/MGUT)
354 TGUT=0.
355 DT=(TZ-TGUT)/FLOAT(NSTEP)
356C Freeze Higgs parameters at HIGFRZ = Drees' value
357C AMTLSS, AMTRSS initialized to 0 for later use in HIGFRZ
358 IF (IMODEL.EQ.1) THEN
359 HIGFRZ=SQRT(M0**2+3*MHF**2)
360 ELSE IF (IMODEL.EQ.2) THEN
361 HIGFRZ=MSUSY
362 ELSE IF (IMODEL.EQ.7) THEN
363 HIGFRZ=SQRT(M0**2+(.01*MHF)**2)
364 END IF
365 AMTLSS=0
366 AMTRSS=0
367 DO 240 II=1,NSTEP+2
368 T=TGUT+(TZ-TGUT)*FLOAT(II-1)/FLOAT(NSTEP)
369 QOLD=Q
370 Q=MGUT*EXP(T)
371 CALL RKSTP(29,DT,T,G,SURG26,W2)
372 IF (Q.LT.AMNRMJ.AND.QOLD.GE.AMNRMJ.AND.FNMZ.EQ.0.) THEN
373 FNMZ=G(27)
374 END IF
375 IF (Q.LT.AMNRMJ) THEN
376 G(27)=0.
377 G(28)=0.
378 G(29)=0.
379 END IF
380 CALL SUGFRZ(Q,G,G0,IG)
381 IF (NOGOOD.NE.0) GO TO 100
382 IF (Q.LT.MZ) GO TO 20
383240 CONTINUE
38420 CONTINUE
385 ASMZ=G0(3)**2/4./PI
386C Electroweak breaking constraints; tree level
387 MUS=(G0(13)-G0(14)*TANB**2)/(TANB**2-1.)-MZ**2/2.
388 IF (MUS.LT.0.) THEN
389 NOGOOD=2
390 GO TO 100
391 END IF
392 MU=SQRT(MUS)*SIGN(1.,SGNMU)
393 B=(G0(13)+G0(14)+2*MUS)*SIN2B/MU/2.
394C Compute tree level masses
395 CALL SUGMAS(G0,0,IMODEL)
396 IF (NOGOOD.NE.0) GO TO 100
397C Compute effective potential corrections
398 CALL SUGEFF(G0,SIG1,SIG2)
399 MH1S=G0(13)+SIG1
400 MH2S=G0(14)+SIG2
401 MUS=(MH1S-MH2S*TANB**2)/(TANB**2-1.)-MZ**2/2.
402 IF (MUS.LT.0.) THEN
403 NOGOOD=2
404 GO TO 100
405 END IF
406 MU=SQRT(MUS)*SIGN(1.,SGNMU)
407 B=(MH1S+MH2S+2*MUS)*SIN2B/MU/2.
408C
409C Recompute weak scale Yukawa couplings including SUSY loops
410C Follow formulae of Pierce et al. NPB491, 3 (1997)
411C
412 M2=G0(8)
413 AM2=ABS(M2)
414 MSN=MSS(16)
415 MG=MSS(1)
416 MT1=MSS(12)
417 MT2=MSS(13)
418 MB1=MSS(10)
419 MB2=MSS(11)
420 MW1=ABS(MSS(27))
421 MW2=ABS(MSS(28))
422 AMU=ABS(MU)
423 XLAM=LOG(MT**2)
424C Be careful in using our convention vs Pierce et al.
425 IF (ABS(COS(THETAT)).LT..707107) THEN
426 MTMT=MT/(1.+4*ASMT/3./PI+(16.11-1.04*(5.-6.63/MT))*(ASMT/PI)**2
427 $ -ASMT/3./PI*(REAL8(SSB1(MT**2,MG,MT1))+
428 $ REAL8(SSB1(MT**2,MG,MT2))-SIN(2*THETAT)*MG/MT*
429 $ (REAL8(SSB0(MT**2,MG,MT1))-REAL8(SSB0(MT**2,MG,MT2)))))
430 ELSE
431 MTMT=MT/(1.+4*ASMT/3./PI+(16.11-1.04*(5.-6.63/MT))*(ASMT/PI)**2
432 $ -ASMT/3./PI*(REAL8(SSB1(MT**2,MG,MT2))+
433 $ REAL8(SSB1(MT**2,MG,MT1))+SIN(2*THETAT)*MG/MT*
434 $ (REAL8(SSB0(MT**2,MG,MT2))-REAL8(SSB0(MT**2,MG,MT1)))))
435 END IF
436 FTMT=MTMT/SINB/VEV
437 XLAM=LOG(MZ**2)
438 IF (ABS(COS(THETAB)).LT..707107) THEN
439 MBMZ=MBMZ*(1.+ASMZ/3./PI*(REAL8(SSB1(MZ**2,MG,MB1))+
440 $ REAL8(SSB1(MZ**2,MG,MB2))-SIN(2*THETAB)*MG/MB*
441 $ (REAL8(SSB0(MZ**2,MG,MB1))-REAL8(SSB0(MZ**2,MG,MB2))))
442 $ -FTMT**2*MU*(-AAT*TANB+MU)/16./PI**2/(MT1**2-MT2**2)*
443 $ (REAL8(SSB0(MZ**2,AMU,MT1))-REAL8(SSB0(MZ**2,AMU,MT2)))+
444 $ G2**2*MU*M2*TANB/16./PI**2/(AMU**2-M2**2)*
445 $ (SIN(THETAT)**2*(REAL8(SSB0(MZ**2,AM2,MT1))-
446 $ REAL8(SSB0(MZ**2,AMU,MT1)))+
447 $ COS(THETAT)**2*(REAL8(SSB0(MZ**2,AM2,MT2))-
448 $ REAL8(SSB0(MZ**2,AMU,MT2)))))
449 ELSE
450 MBMZ=MBMZ*(1.+ASMZ/3./PI*(REAL8(SSB1(MZ**2,MG,MB2))+
451 $ REAL8(SSB1(MZ**2,MG,MB1))+SIN(2*THETAB)*MG/MB*
452 $ (REAL8(SSB0(MZ**2,MG,MB2))-REAL8(SSB0(MZ**2,MG,MB1))))
453 $ -FTMT**2*MU*(-AAT*TANB+MU)/16./PI**2/(MT2**2-MT1**2)*
454 $ (REAL8(SSB0(MZ**2,AMU,MT2))-REAL8(SSB0(MZ**2,AMU,MT1)))+
455 $ G2**2*MU*M2*TANB/16./PI**2/(AMU**2-M2**2)*
456 $ (COS(THETAT)**2*(REAL8(SSB0(MZ**2,AM2,MT2))-
457 $ REAL8(SSB0(MZ**2,AMU,MT2)))+
458 $ SIN(THETAT)**2*(REAL8(SSB0(MZ**2,AM2,MT1))-
459 $ REAL8(SSB0(MZ**2,AMU,MT1)))))
460 END IF
461 FBMZ=MBMZ/COSB/VEV
462 MTAMZ=MTAMZ*(1.+G2**2*MU*M2*TANB/16./PI**2/(MUS-M2**2)*
463 $(REAL8(SSB0(MZ**2,AM2,MSN))-REAL8(SSB0(MZ**2,AMU,MSN))))
464 FTAMZ=MTAMZ/COSB/VEV
465C
466C Iterate entire process, increasing NSTEP each time
467C This time, freeze out parameters at sqrt(t_l t_r)
468C
469 HIGFRZ=MAX(AMZ,(G0(23)*G0(24))**0.25)
470 DO 300 I=1,MXITER
471 DO 310 J=1,26
472310 G0SAVE(J)=G0(J)
473 NSTEP=1.2*NSTEP
474 CALL SUGRGE(M0,MHF,A0,TANB,SGNMU,MT,G,G0,IG,W2,NSTEP,IMODEL)
475 IF(NOGOOD.NE.0) GO TO 100
476 DELG0=0.
477 DO 320 J=1,24
478320 DELG0=MAX(DELG0,ABS((G0(J)-G0SAVE(J))/G0(J)))
479 IF(DELG0.LT.DELLIM) GO TO 400
480300 CONTINUE
481 WRITE(LOUT,1000) MXITER
4821000 FORMAT(/' SUGRA WARNING: NO CONVERGENCE IN',I4,' ITERATIONS')
483C
484C Save results
485C
486400 DO 410 I=1,26
487 GSS(I)=G0(I)
488410 CONTINUE
489 MGUTSS=MGUT
490 AGUTSS=AGUT
491 GGUTSS=GGUT
492C
493C Fill XISAIN common block
494C
495 XISAIN(1)=MSS(1)
496 XISAIN(2)=MU
497 XISAIN(3)=MSS(31)
498 XISAIN(4)=TANB
499 XISAIN(5)=SQRT(G0(19))
500 XISAIN(6)=SQRT(G0(17))
501 XISAIN(7)=SQRT(G0(18))
502 XISAIN(8)=SQRT(G0(16))
503 XISAIN(9)=SQRT(G0(15))
504 XISAIN(10)=XISAIN(5)
505 XISAIN(11)=XISAIN(6)
506 XISAIN(12)=XISAIN(7)
507 XISAIN(13)=XISAIN(8)
508 XISAIN(14)=XISAIN(9)
509 XISAIN(15)=SQRT(G0(24))
510 XISAIN(16)=SQRT(G0(22))
511 XISAIN(17)=SQRT(G0(23))
512 XISAIN(18)=SQRT(G0(21))
513 XISAIN(19)=SQRT(G0(20))
514 XISAIN(20)=G0(12)
515 XISAIN(21)=G0(11)
516 XISAIN(22)=G0(10)
517 XISAIN(23)=G0(7)
518 XISAIN(24)=G0(8)
519 M2=G0(8)
520100 RETURN
521 END