📄 exm_10.bas
字号:
Attribute VB_Name = "Exm_10"
Option Explicit
Dim i%, j%, FLOWIN#, FL#, AFL#, VMIN#, FACTOR#, REL#, AMT#, GAMP#, GAMM#
Dim DUDX#, DVDY#, DUDY#, DVDX#, DISS#
Dim AMU#, CMU#, C1#, C2#, PRT#, PRK#, PRD#, PR#, PRPRT#, PFN#, CMU4#
Dim XPLUS(22) As Double
Sub GRID()
For i = 1 To 6
LSolve(i) = True
LPrint(i) = True
Next i
LPrint(7) = True
LPrint(11) = True
RELAX(1) = 0.5
RELAX(2) = 0.5
RELAX(5) = 0.4
RELAX(6) = 0.4
RELAX(11) = 0.8
RELAX(13) = 0.5
Title(1) = ".VEL U."
Title(2) = ".VEL V."
Title(3) = ".STR FN."
Title(4) = ". TEMP ."
Title(5) = "KIN ENE"
Title(6) = ".DISIPA."
Title(7) = "TURB VI"
Title(11) = "PRESSURE"
Title(12) = "DENSITY"
LAST = 50
XL = 1#
YL = 4#
L1 = 7
M1 = 9
Call UGRID
End Sub
Sub START()
For j = 1 To M1
For i = 1 To L1
F(i, j, 1) = 0#
F(i, j, 2) = 10#
F(1, j, 2) = 0#
F(i, 2, 2) = 10#
If (i > 4) Then F(i, 2, 2) = 100#
F(i, j, 4) = 100#
F(1, j, 4) = 0#
If (i > 4) Then F(i, 1, 4) = 400#
F(i, j, 5) = 0.005 * F(i, 2, 2) ^ 2
F(i, j, 6) = 0.1 * F(i, j, 5) ^ 2
Next i
Next j
AMU = 0.000001
CMU = 0.09
C1 = 1.44
C2 = 1.92
PRT = 0.9
PRK = 1#
PRD = 1.3
PR = 0.7
PRPRT = PR / PRT
PFN = 9# * (PRPRT - 1#) / PRPRT ^ 0.25
CMU4 = CMU ^ 0.25
End Sub
Sub DENSE()
End Sub
Sub BOUND()
If (ITER = 0) Then
FLOWIN = 0#
For i = 2 To L2
FLOWIN = FLOWIN + RHO(i, 1) * F(i, 2, 2) * XCV(i)
Next i
End If
FL = 0#
AFL = 0#
VMIN = 0#
For i = 2 To L2
If (F(i, M2, 2) < 0#) Then VMIN = Max(VMIN, -1 * F(i, M2, 2))
AFL = AFL + RHO(i, M1) * XCV(i)
FL = FL + RHO(i, M1) * F(i, M2, 2) * XCV(i)
Next i
FACTOR = FLOWIN / (FL + AFL * VMIN)
For i = 2 To L2
F(i, M1, 2) = (F(i, M2, 2) + VMIN) * FACTOR
Next i
For j = 2 To M2
F(L1, j, 5) = F(L2, j, 5)
F(L1, j, 6) = F(L2, j, 6)
Next j
End Sub
Sub OUTPUT()
If (ITER = 0) Then
'WRITE(8,401)
Print #1,
'401 FORMAT(1X,' ITER',6X,'SMAX',6X,'SSUM',5X,'V(6,6)', 4X,'T(5,6)',4X,'KE(5,6)')
Print #1,
Print #1, " ITER", Space(6), "SMAX", Space(6), "SSUM", Space(5), "V(6,6)", Space(4), "T(5,6)", Space(4), "KE(5,6)"
End If
'400 WRITE(8,403) ITER,SMAX,SSUM,V(6,6),T(5,6),AKE(5,6)
Print #1,
Print #1, Format(ITER, "0# "); Format(SMAX, "#.000E+###"); Space(4); Format(SSUM, "#.000E+###"); _
Space(4); Format(F(6, 6, 2), "#.000E+###"); Space(4); Format(F(5, 6, 4), "#.000E+###"); _
Space(4); Format(F(5, 6, 5), "#.000E+###");
'403 FORMAT(1X,I6,1P5E11.3)
If (ITER = LAST) Then
Call ToPrint
End If
End Sub
Sub GAMSOR()
If NF = 3 Then
Exit Sub
End If
If NF = 1 Then
REL = 1# - RELAX(NGAM)
For j = 1 To M1
For i = 1 To L1
AMT = CMU * RHO(i, j) * F(i, j, 5) ^ 2 / (F(i, j, 6) + 1E-30)
If (ITER = 0) Then
F(i, j, 7) = AMT
End If
F(i, j, 7) = RELAX(NGAM) * AMT + REL * F(i, j, 7)
Next i
Next j
End If
FACTOR = 1#
If (NF = 4) Then FACTOR = 1# / PRT
If (NF = 5) Then FACTOR = 1# / PRK
If (NF = 6) Then FACTOR = 1# / PRD
For j = 1 To M1
For i = 1 To L1
GAM(i, j) = F(i, j, 7) * FACTOR
If (NF <> 1) Then
GAM(L1, j) = 0#
End If
GAM(i, M1) = 0#
Next i
Next j
For j = 2 To M2
' GOTO(531,532,531,533,531,531),NF
'531 GAM(1, j) = 0#
' GoTo 530
'532 XPLUS(j) = RHO(2, j) * SQRT(AKE(2, j)) * CMU4 * XDIF(2) / AMU
' GAM(1, j) = AMU
' IF(XPLUS(J).GT.11.5) GAM(1,J)=AMU*XPLUS(J)/(ALOG(9.*XPLUS(J))*2.5)
' GoTo 530
'533 GAM(1, j) = AMU / PR
' IF(XPLUS(J).GT.11.5) GAM(1,J)=AMU/PRT*XPLUS(J)
' & /(2.5*ALOG(9.*XPLUS(J))+PFN)
Select Case NF
Case 1, 3, 5, 6
GAM(1, j) = 0#
Case 2
XPLUS(j) = RHO(2, j) * Sqr(F(2, j, 5)) * CMU4 * XDIF(2) / AMU
GAM(1, j) = AMU
If (XPLUS(j) > 11.5) Then
GAM(1, j) = AMU * XPLUS(j) / (Log(9# * XPLUS(j)) * 2.5)
End If
Case 4
GAM(1, j) = AMU / PR
If (XPLUS(j) > 11.5) Then
GAM(1, j) = AMU / PRT * XPLUS(j) / (2.5 * Log(9# * XPLUS(j)) + PFN)
End If
End Select
Next j
If (NF = 1) Then
For j = 2 To M2
For i = 3 To L2
CON(i, j) = (GAM(i, j) * (F(i + 1, j, 1) - F(i, j, 1)) / XCV(i) _
- GAM(i - 1, j) * (F(i, j, 1) - F(i - 1, j, 1)) / XCV(i - 1)) / XDIF(i)
GAMP = GAM(i, j + 1) * GAM(i - 1, j + 1) / (GAM(i, j + 1) + GAM(i - 1, j + 1) + 1E-30)
GAMP = GAMP + GAM(i, j) * GAM(i - 1, j) / (GAM(i, j) + GAM(i - 1, j) + 1E-30)
GAMM = GAM(i, j - 1) * GAM(i - 1, j - 1) / (GAM(i, j - 1) + GAM(i - 1, j - 1) + 1E-30)
GAMM = GAMM + GAM(i, j) * GAM(i - 1, j) / (GAM(i, j) + GAM(i - 1, j) + 1E-30)
CON(i, j) = CON(i, j) + (GAMP * (F(i, j + 1, 2) - F(i - 1, j + 1, 2)) _
- GAMM * (F(i, j, 2) - F(i - 1, j, 2))) / (YCV(j) * XDIF(i))
AP(i, j) = 0#
Next i
Next j
Exit Sub
End If
If (NF = 2) Then
For j = 3 To M2
For i = 2 To L2
CON(i, j) = (GAM(i, j) * (F(i, j + 1, 2) - F(i, j, 2)) / YCV(j) - _
GAM(i, j - 1) * (F(i, j, 2) - F(i, j - 1, 2)) / YCV(j - 1)) / (YDIF(j))
GAMP = GAM(i + 1, j) * GAM(i + 1, j - 1) / (GAM(i + 1, j) + GAM(i + 1, j - 1) + 1E-30)
GAMP = GAMP + GAM(i, j) * GAM(i, j - 1) / (GAM(i, j) + GAM(i, j - 1) + 1E-30)
GAMM = GAM(i - 1, j) * GAM(i - 1, j - 1) / (GAM(i - 1, j) + GAM(i - 1, j - 1) + 1E-30)
GAMM = GAMM + GAM(i, j) * GAM(i, j - 1) / (GAM(i, j) + GAM(i, j - 1) + 1E-30)
CON(i, j) = CON(i, j) + (GAMP * (F(i + 1, j, 1) - F(i + 1, j - 1, 1)) _
- GAMM * (F(i, j, 1) - F(i, j - 1, 1))) / (XCV(i) * YDIF(j))
AP(i, j) = 0#
Next i
Next j
Exit Sub
End If
If (NF = 4) Then
For j = 2 To M2
For i = 2 To L2
CON(i, j) = 0#
AP(i, j) = 0#
Next i
Next j
Exit Sub
End If
If (NF = 5) Then
For j = 2 To M2
For i = 2 To L2
DUDX = (F(i + 1, j, 1) - F(i, j, 1)) / XCV(i)
DVDY = (F(i, j + 1, 2) - F(i, j, 2)) / YCV(j)
' IF(J.EQ.2) GOTO 581
' If (j.EQ.M2) Then GoTo 582
' DUDY = (0.5 * (U(i, j + 1) - U(i, j - 1)) + 0.5 * (U(i + 1, j + 1) - U(i + 1, j - 1)))
' & /(YDIF(J)+YDIF(J+1))
' GoTo 583
'581 DUDY = (0.5 * (U(i, j + 1) - U(i, j)) + 0.5 * (U(i + 1, j + 1) - U(i + 1, j))) / YDIF(j + 1)
' GoTo 583
'582 DUDY = (0.5 * (U(i, j) - U(i, j - 1)) + 0.5 * (U(i + 1, j) - U(i + 1, j - 1)))
' & /YDIF(J)
If j = 2 Then
DUDY = (0.5 * (F(i, j + 1, 1) - F(i, j, 1)) + 0.5 * (F(i + 1, j + 1, 1) - F(i + 1, j, 1))) / YDIF(j + 1)
ElseIf j = M2 Then
DUDY = (0.5 * (F(i, j, 1) - F(i, j - 1, 1)) + 0.5 * (F(i + 1, j, 1) - F(i + 1, j - 1, 1))) / YDIF(j)
Else
DUDY = (0.5 * (F(i, j + 1, 1) - F(i, j - 1, 1)) + 0.5 * (F(i + 1, j + 1, 1) - F(i + 1, j - 1, 1))) / (YDIF(j) + YDIF(j + 1))
End If
' 583 IF(I.EQ.2) GOTO 584
' If (i.EQ.L2) Then GoTo 585
' DVDX=(0.5*(V(I+1,J)-V(I-1,J))+0.5*(V(I+1,J+1)
' & -V(I-1,J+1)))/(XDIF(I)+XDIF(I+1))
' GoTo 586
' 584 DVDX=(0.5*(V(I+1,J)-V(I,J))+0.5*(V(I+1,J+1)
' & -V(I,J+1)))/XDIF(I+1)
' GoTo 586
' 585 DVDX=(0.5*(V(I,J)-V(I-1,J))+0.5*(V(I,J+1)
' & -V(I-1,J+1)))/XDIF(I)
If i = 2 Then
DVDX = (0.5 * (F(i + 1, j, 2) - F(i, j, 2)) + 0.5 * (F(i + 1, j + 1, 2) - F(i, j + 1, 2))) / XDIF(i + 1)
ElseIf i = L2 Then
DVDX = (0.5 * (F(i, j, 2) - F(i - 1, j, 2)) + 0.5 * (F(i, j + 1, 2) - F(i - 1, j + 1, 2))) / XDIF(i)
Else
DVDX = (0.5 * (F(i + 1, j, 2) - F(i - 1, j, 2)) + 0.5 * (F(i + 1, j + 1, 2) - F(i - 1, j + 1, 2))) / (XDIF(i) + XDIF(i + 1))
End If
F(i, j, 8) = 2# * (DUDX ^ 2 + DVDY ^ 2) + (DUDY + DVDX) ^ 2
CON(i, j) = F(i, j, 8) * F(i, j, 7)
AP(i, j) = -RHO(i, j) * F(i, j, 6) / (F(i, j, 5) + 1E-30)
Next
Next
Exit Sub
End If
For j = 2 To M2
For i = 2 To L2
CON(i, j) = C1 * F(i, j, 8) * CMU * RHO(i, j) * F(i, j, 5)
AP(i, j) = -C2 * RHO(i, j) * F(i, j, 6) / (F(i, j, 5) + 1E-30)
Next i
Next j
For j = 2 To M2
DISS = CMU * F(2, j, 5) ^ 1.5 / (0.4 * CMU4 * XDIF(2))
CON(2, j) = 1E+30 * DISS
AP(2, j) = -1E+30
Next j
End Sub
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -