GOLDAK HEAT SOURCE

classic Classic list List threaded Threaded
1 message Options
Reply | Threaded
Open this post in threaded view
|

GOLDAK HEAT SOURCE

Ahmad Abbasi
HI GUYSI HAVE BEEN TRYING TO SIMULATE GOLDAK DOUBLE ELLIPSOIDAL HEAT SOURCE BUT MY SUBROUTINE WORKS ONLY WITH A PART WITH SPECIAL PARTITION. I WILL ATTACH TO MODEL. MY SUB WORKS WITH OK.CAE CORRECTLY BUT IT DOES'NT WORK WITH WRONG.CAE . CAN YOU HELP ME FIGURE THIS OUT? BEST REGARDS  
  ----------

       SUBROUTINE DFLUX(FLUX,SOL,KSTEP,KINC,TIME,NOEL,NPT,COORDS,
     1 JLTYP,TEMP,PRESS,SNAME)
C
      INCLUDE 'ABA_PARAM.INC'
C
      DIMENSION FLUX(2), TIME(2), COORDS(3)
      CHARACTER*80 SNAME
     
      DOUBLE PRECISION  FF, FR, QR, QF, HF, HR, SF, SR, Q, T
      REAL              X, Y, Z, X0, Y0, Z0, AF, AR, B, C, V
     1 PI, ETA, CURRENT, VOLTAGE
      !-----------------------------------------------------------------------
      V = 1.5
      T = TIME(1)
      X0 = 0
      Z0 = 2
      Y0 = V * T
      !-------------------------------------------------------------------------
      PI = 3.14159265359
      CURRENT = 120.0
      VOLTAGE = 28.2
      ETA = 0.85
      AF = 0.5
      AR = 1.5
      B = 0.5
      C = 0.5
      X = COORDS(1) - X0
      Y = COORDS(2) - Y0
      Z = COORDS(3) - Z0
      !-------------------------------------------------------------------------
      Q = CURRENT * VOLTAGE
      FF = (2 * AF) / (AF + AR)
      FR = (2 * AR) / (AF + AR)
      HF = (6*SQRT(3)*FF*ETA*Q) / (AF*B*C*PI*SQRT(PI))
      HR = (6*SQRT(3)*FR*ETA*Q) / (AR*B*C*PI*SQRT(PI))
      SF = EXP(-(3*X**2/AF**2)-(3*Y**2/B**2)-(3*Z**2/C**2))
      SR = EXP(-(3*X**2/AR**2)-(3*Y**2/B**2)-(3*Z**2/C**2))
      QF = HF * SF
      QR = HR * SR
      !------------------------------------------------------------------------
      IF (Y0 - AR <= COORDS(2) <= Y0) THEN
          FLUX(1) = QR
      ELSEIF (Y0 <= COORDS(2) <= Y0+AF ) THEN
          FLUX(1) = QF
      ELSE
          FLUX(1) = 0
      ENDIF
      RETURN
      END
  ----------

       SUBROUTINE DFLUX(FLUX,SOL,KSTEP,KINC,TIME,NOEL,NPT,COORDS,
     1 JLTYP,TEMP,PRESS,SNAME)
C
      INCLUDE 'ABA_PARAM.INC'
C
      DIMENSION FLUX(2), TIME(2), COORDS(3)
      CHARACTER*80 SNAME
     
      DOUBLE PRECISION  FF, FR, QR, QF, HF, HR, SF, SR, Q, T
      REAL              X, Y, Z, X0, Y0, Z0, AF, AR, B, C, V
     1 PI, ETA, CURRENT, VOLTAGE
C-----------------------------------------------------------------------
      V = 0.00333
      T = TIME(1)
      X0 = 0
      Z0 = V * T
      Y0 = 0
C-------------------------------------------------------------------------
      PI = 3.14159265359
      CURRENT = 120.0
      VOLTAGE = 28.2
      ETA = 0.85
      AF = 4.35E-3
      AR = 13.05E-3
      B = 4.35E-3
      C = 4.35E-3
      X = COORDS(1) - X0
      Y = COORDS(2) - Y0
      Z = COORDS(3) - Z0
C-------------------------------------------------------------------------
      Q = CURRENT * VOLTAGE
      FF = (2 * AF) / (AF + AR)
      FR = (2 * AR) / (AF + AR)
      HF = (6*SQRT(3)*FF*ETA*Q) / (AF*B*C*PI*SQRT(PI))
      HR = (6*SQRT(3)*FR*ETA*Q) / (AR*B*C*PI*SQRT(PI))
      SF = EXP(-(3*X**2/AF**2)-(3*Y**2/B**2)-(3*Z**2/C**2))
      SR = EXP(-(3*X**2/AR**2)-(3*Y**2/B**2)-(3*Z**2/C**2))
      QF = HF * SF
      QR = HR * SR
C-------------------------------------------------------------------------
      IF (Z0-AR <= COORDS(3) <= Z0) THEN
          FLUX(1) = QR
      ELSEIF (Z0 <= COORDS(3) <= Z0+AF ) THEN
          FLUX(1) = QF
      ELSE
          FLUX(1) = 0
      ENDIF  
      RETURN
      END


[Non-text portions of this message have been removed]