ادامه
ادامه برنامه های CFD
Rectangle heating
Dim Lx, Ly, Delta_x, Delta_y, x(1000, 1000), Y(1000, 1000), Beta, Delta_M, T1, T2, D_T As Double
Dim Temp_A, Temp_B, Temp_C, Temp_INI, T(1000, 1000), U(1000, 1000), Epsilon, W As Double
Dim A(1000), B(1000), C(1000), U2(1000), R(1000), Wopt, a1
Dim Opt, Imax, Jmax, i, j As Integer
Dim exportfile
Private Sub Command1_Click()
L1:
T1 = Timer
If Opt_Jacobi = True Then Opt = 1
If Opt_G_P_By_P = True Then Opt = 2
If Opt_P_S_O_R = True Then Opt = 3
If Opt_L_by_L_GS = True Then Opt = 4
Imax = Val(txt_Imax.Text)
Jmax = Val(txt_Jmax.Text)
Lx = Val(txt_Lx.Text)
Ly = Val(txt_Ly.Text)
Delta_x = Lx / (Imax - 1)
Delta_y = Ly / (Jmax - 1)
Beta = Delta_x / Delta_y
Epsilon = Val(txt_Epsilon.Text)
Temp_A = Val(Txt_Temp_A.Text)
Temp_B = Val(Txt_temp_B.Text)
Temp_C = Val(Txt_Temp_C.Text)
Temp_INI = Val(Txt_INI_Temp.Text)
' Grid Generation
For j = 1 To Jmax
For i = 1 To Imax
x(i, j) = (i * Delta_x) - Delta_x
Next i
Next j
For i = 1 To Imax
For j = 1 To Jmax
Y(i, j) = (j * Delta_y) - Delta_y
Next j
Next i
' End Of Grid Generation
'----------------------------------------------------------------------------
For j = 1 To Jmax
T(1, j) = Temp_A
U(1, j) = Temp_A
Next j
For i = 2 To Imax
T(i, 1) = (((Temp_A - Temp_C) / (1 - Imax)) * (i - Imax)) + Temp_C
U(i, 1) = (((Temp_A - Temp_C) / (1 - Imax)) * (i - Imax)) + Temp_C
Next i
For i = 2 To Imax
T(i, Jmax) = (((Temp_A - Temp_B) / (1 - Imax)) * (i - Imax)) + Temp_B
U(i, Jmax) = (((Temp_A - Temp_B) / (1 - Imax)) * (i - Imax)) + Temp_B
Next i
'----------------------------------------------------------------------------
For i = 2 To Imax - 1
For j = 2 To Jmax - 1
T(i, j) = Temp_INI
U(i, j) = Temp_INI
Next j
Next i
'----------------------------------------------------------------------------
' Jacobi Method
If Opt = 1 Then
exportfile = "Solve With Jacobi Meth.dat"
Open exportfile For Output As 1
Print #1, " Solved By Jacobi Metod "
Print #1, "variables= x , y , T"
Print #1, "zone i=", Jmax, "j=", Imax
nj = 0
J1:
nj = nj + 1
For j = 2 To Jmax - 1
For i = 2 To Imax
T(i, j) = (1 / (2 * (1 + Beta ^ 2))) * (U(i + 1, j) + U(i - 1, j) + (Beta ^ 2) * (U(i, j + 1) + U(i, j - 1)))
If i = Imax Then
T(i, j) = T(i - 1, j)
U(i, j) = U(i - 1, j)
End If
Next i
Next j
For i = 2 To Imax - 1
For j = 2 To Jmax - 1
Delta_M = Abs(U(i, j) - T(i, j))
If Delta_M > Epsilon Then
For m = 2 To Imax - 1
For N = 2 To Jmax - 1
U(m, N) = T(m, N)
Next N
Next m
GoTo J1
End If
Next j
Next i
T2 = Timer
D_T = T2 - T1
Text1.Text = Str(nj)
Text2.Text = Str(D_T) + " s"
End If
' End Jacobi Method
'----------------------------------------------------------------------------
' Gauss-Sidel(Point by point) Method
If Opt = 2 Then
exportfile = "S With GSPBP Meth.dat"
Open exportfile For Output As 1
Print #1, " Solved By Gauss-Sidel(Point by point) Method"
Print #1, "variables= x , y , T"
Print #1, "zone i=", Jmax, "j=", Imax
ngp = 0
Gp1:
ngp = ngp + 1
For j = 2 To Jmax - 1
For i = 2 To Imax
T(i, j) = (1 / (2 * (1 + Beta ^ 2))) * (T(i + 1, j) + T(i - 1, j) + (Beta ^ 2) * (T(i, j + 1) + T(i, j - 1)))
If i = Imax Then
T(i, j) = T(i - 1, j)
End If
Next i
Next j
For i = 2 To Imax - 1
For j = 2 To Jmax - 1
Delta_M = Abs(U(i, j) - T(i, j))
If Delta_M > Epsilon Then
For m = 2 To Imax - 1
For N = 2 To Jmax - 1
U(m, N) = T(m, N)
Next N
Next m
GoTo Gp1
End If
Next j
Next i
T2 = Timer
D_T = T2 - T1
Text3.Text = ngp
Text4.Text = Str(D_T) + " s"
End If
' End Gauss-Sidel(Point by point) Method
'----------------------------------------------------------------------------
' The P.S.O.R Method
If Opt = 3 Then
a1 = ((Cos(4 * Atn(1) / (Imax - 1)) + (Beta ^ 2) * Cos(4 * Atn(1) / (Jmax - 1))) / (1 + Beta ^ 2)) ^ 2
Wopt = (2 - 2 * (1 - a1) ^ 0.5) / a1
W = InputBox("Inter The Under(Over) Relaxation Value(0<2)", "CFD 6", Wopt)
T1 = Timer
exportfile = "Solve With PSOR Meth.dat"
Open exportfile For Output As 1
Print #1, " Solved By P.S.O.R Method"
Print #1, "variables= x , y , T"
Print #1, "zone i=", Jmax, "j=", Imax
np = 0
p1:
nj = nj + 1
For j = 2 To Jmax - 1
For i = 2 To Imax
T(i, j) = ((1 - W) * U(i, j)) + (W / (2 * (1 + Beta ^ 2))) * (U(i + 1, j) + T(i - 1, j) + (Beta ^ 2) * (U(i, j + 1) + T(i, j - 1)))
If i = Imax Then
T(i, j) = T(i - 1, j)
U(i, j) = U(i - 1, j)
End If
Next i
Next j
For i = 2 To Imax - 1
For j = 2 To Jmax - 1
Delta_M = Abs(U(i, j) - T(i, j))
If Delta_M > 1E+308 Then
dd = MsgBox("Divergence Detected Please Enter A New Value For W", vbOKOnly + vbExclamation, "CFD 6")
Close #1
GoTo L1
