ادامه

ادامه برنامه های 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