giant.f

Go to the documentation of this file.

      SUBROUTINE giant(IPAIR,I,W,Q,WSCALE,QSCALE,XN,QL)
*
*
*       Structure constants of giant star.
*       ----------------------------------
*
*       Theory of Rosemary Mardling, Ap. J. XX, YYY, 1995.
*       @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
*
      include 'common6.h'
      common/modes/  eb0(ntmax),zj0(ntmax),ecrit(ntmax),ar(ntmax),
     &               br(ntmax),eosc(4,ntmax),edec(ntmax),tosc(ntmax),
     &               rp(ntmax),es(ntmax),cm(2,ntmax),iosc(ntmax),
     &               namec(ntmax)
      REAL*8  ww(6),qq(6),w(2),q(2),wscale(2),qscale(2),sw(2)
      DATA  ww  /2.119,3.113,8.175,3.742,4.953,9.413/
      DATA  qq  /0.4909,0.4219,0.2372,0.4677,0.3560,0.1519/
      DATA  a0,a1,a2,a3  /0.944525,-0.392030,6.01655d-02,-3.34790d-03/
      DATA  b0,b1,b2  /-0.3789,1.481,-0.1018/
      DATA  c0,c1,c2,c3  /1.452104,3.923872,-11.88722,13.46106/
      DATA  e0,e1,e2,e3  /1.934977,2.214222,-4.855796,4.025394/
*
*
*       Determine the appropriate chaos index (adopt NCHAOS+1 if not found).
      ic = 0
      j = n + ipair
      DO 1 k = 1,nchaos
          IF (namec(k).EQ.name(j)) ic = k
    1 CONTINUE
*
*       Decide between first & second binary component for saving core mass.
      i1 = 2*ipair - 1
      i2 = i1 + 1
      l = 1
*       Note possibility I = I2 on first call.
      IF (radius(i1).GE.radius(i2)) THEN
          IF (i.EQ.i2) l = 2
      END IF
*
*       Set typical core mass of 0.3/0.5 Msun for general binary in next bin.
      IF (ic.EQ.0) THEN
          ic = nchaos + 1
          cm(l,ic) = (0.3 + 0.1*float(kstar(i) - 3))/zmbar
*       Include rare case of hyperbolic encounter as the first event.
      ELSE IF (cm(l,ic).LE.0.0d0) THEN
          cm(l,ic) = (0.3 + 0.1*float(kstar(i) - 3))/zmbar
      END IF
      IF(kstar(i).EQ.9) cm(l,ic) = 0.5*body(i)
*
*       Form ratio of core mass and mass.
      sig = cm(l,ic)/body(i)
*
*       Include safety check on mass ratio just in case.
      IF (sig.GT.0.9.OR.sig.LT.0.0) THEN
*         WRITE (6,5)  IC, KSTAR(I), NAME(I), BODY0(I)*ZMBAR,
*    &                 BODY(I)*ZMBAR, CM(L,IC)*ZMBAR
*   5     FORMAT (' WARNING!    GIANT    IC K* NM M0 M MC ',
*    &                                   2I4,I6,3F7.2)
          sig = 0.9
          cm(l,ic) = 0.9*body(i)
      END IF
*
*       Define mass, core mass, radius, envelope mass and luminosity in S.U.
      zm = body(i)*zmbar
      zmc = cm(l,ic)*zmbar
      rsi = radius(i)*su
      zme = zm - zmc
*       Obtain L^{1/3} from GB L(Mc) relation for Pop II M = 0.8Msun.
      zl = 1.98d+05*zmc**6
*       Evaluate damping constant from Zahn's theory (R.M. 13/5/97).
*       FAC = (GM)^{1/2}*M_{env}^{1/3)/(R^{5/6}*L^{1/3}) = 8.48D+03 for S.U.
      ql = 8.48d+03*sqrt(zm)*(zme/zl)**0.33/rsi**0.833
*
*       Set effective frequencies, overlap integrals and structure constants.
      DO 10 k = 1,2
          k1 = 3*k - 2
      IF (k.EQ.1) THEN
              sw(k) = ((c3*sig + c2)*sig + c1)*sig + c0
          ELSE
              sw(k) = ((e3*sig + e2)*sig + e1)*sig + e0
      END IF
      w(k) = sw(k)**2
          q(k) = ((a3*sw(k) + a2)*sw(k) + a1)*sw(k) + a0
          wscale(k) = sqrt(w(k)/ww(k1))
          qscale(k) = (q(k)/qq(k1)/wscale(k))**2
          qscale(k) = max(qscale(k),0.0001d0)
   10 CONTINUE
*
*       Evaluate new polytropic index.
      xn = (b2*sw(1) + b1)*sw(1) + b0
*     WRITE (24,20)  IC, IPAIR, KSTAR(J), CM(L,IC)/BODY(I), RSI, XN, QL
*  20 FORMAT (' GIANT:    IC KS K* MC/M R* n Q ',3I4,3F6.2,F7.1)
*     CALL FLUSH(24)
*
*       Include warning if n > 5 (note limit QSCALE >= 0.0001).
*     QSM = MIN(QSCALE(1),QSCALE(2))
*     IF (XN.GE.5.0.OR.QSM.LT.0.00011) THEN
*         WRITE (6,30)  IC, IPAIR, KSTAR(J), CM(L,IC)/BODY(I),
*    &                  RADIUS(I)*SU, QSM, XN, QL
*  30     FORMAT (' GIANT:    WARNING!    IC KS K* MC/M R/R0 QSM n QL ',
*    &                                    3I4,F6.2,F7.1,F8.4,F6.2,F7.1)
*     END IF

      RETURN
*
      END