Main Cavity.f90 - Main modul
From CFD-Wiki
(Difference between revisions)
(changed core code with fixed bugs) |
|||
| Line 4: | Line 4: | ||
! Main modul | ! Main modul | ||
!Copyright (C) 2010 Michail Kiričkov | !Copyright (C) 2010 Michail Kiričkov | ||
| - | + | ! 2016 - I fixed some bugs | |
!This program is free software; you can redistribute it and/or | !This program is free software; you can redistribute it and/or | ||
!modify it under the terms of the GNU General Public License | !modify it under the terms of the GNU General Public License | ||
| Line 21: | Line 21: | ||
!----------------------------------------------------------------------- | !----------------------------------------------------------------------- | ||
Program Main | Program Main | ||
| - | |||
include 'icomm_1.f90' | include 'icomm_1.f90' | ||
| - | |||
Dimension F_out(nx,ny),x_exct(30), v_exct(30) | Dimension F_out(nx,ny),x_exct(30), v_exct(30) | ||
Dimension y_exct(30), u_exct(30) | Dimension y_exct(30), u_exct(30) | ||
Character Filename*10 | Character Filename*10 | ||
| + | real qqq | ||
| + | iii = 150 | ||
Call Ini | Call Ini | ||
Call Grid_rectangular | Call Grid_rectangular | ||
| - | |||
Call Geom | Call Geom | ||
Call Init_all_cavity | Call Init_all_cavity | ||
| - | + | write(*,*) Yc(1,N_centr) | |
| + | read(*,*) qqq | ||
| + | x_exct(1)=0.0000 ; v_exct(1)= 0.00000 | ||
| + | x_exct(2)=0.0625 ; v_exct(2)= 0.27485 | ||
| + | x_exct(3)=0.0703 ; v_exct(3)= 0.29012 | ||
| + | x_exct(4)=0.0781 ; v_exct(4)= 0.30353 | ||
| + | x_exct(5)=0.0938 ; v_exct(5)= 0.32627 | ||
| + | x_exct(6)=0.1563 ; v_exct(6)= 0.37095 | ||
| + | x_exct(7)=0.2266 ; v_exct(7)= 0.33075 | ||
| + | x_exct(8)=0.2344 ; v_exct(8)= 0.32235 | ||
| + | x_exct(9)=0.5000 ; v_exct(9)= 0.02526 | ||
| + | x_exct(10)=0.8047 ; v_exct(10)=-0.31966 | ||
| + | x_exct(11)=0.8594 ; v_exct(11)=-0.42665 | ||
| + | x_exct(12)=0.9063 ; v_exct(12)=-0.51550 | ||
| + | x_exct(13)=0.9453 ; v_exct(13)=-0.39188 | ||
| + | x_exct(14)=0.9531 ; v_exct(14)=-0.33714 | ||
| + | x_exct(15)=0.9609 ; v_exct(15)=-0.27669 | ||
| + | x_exct(16)=0.9688 ; v_exct(16)=-0.21388 | ||
| + | x_exct(17)=1.0000 ; v_exct(17)=-0.00000 | ||
| + | |||
| + | y_exct(1)=0.00000 ; u_exct(1)=-0.00000 | ||
| + | y_exct(2)=0.0547 ; u_exct(2)=-0.18109 | ||
| + | y_exct(3)=0.0625 ; u_exct(3)=-0.20196 | ||
| + | y_exct(4)=0.0703 ; u_exct(4)=-0.22220 | ||
| + | y_exct(5)=0.1016 ; u_exct(5)=-0.29730 | ||
| + | y_exct(6)=0.1719 ; u_exct(6)=-0.38289 | ||
| + | y_exct(7)=0.2813 ; u_exct(7)=-0.27805 | ||
| + | y_exct(8)=0.4531 ; u_exct(8)=-0.10648 | ||
| + | y_exct(9)=0.5000 ; u_exct(9)=-0.06080 | ||
| + | y_exct(10)=0.6172 ; u_exct(10)=0.05702 | ||
| + | y_exct(11)=0.7344 ; u_exct(11)=0.18719 | ||
| + | y_exct(12)=0.8516 ; u_exct(12)=0.33304 | ||
| + | y_exct(13)=0.9531 ; u_exct(13)=0.46604 | ||
| + | y_exct(14)=0.9609 ; u_exct(14)=0.51117 | ||
| + | y_exct(15)=0.9688 ; u_exct(15)=0.57492 | ||
| + | y_exct(16)=0.9766 ; u_exct(16)=0.65928 | ||
| + | y_exct(17)=1.0000 ; u_exct(17)=1.00000 | ||
| + | |||
| + | open(55,file='Hystory_UVV.dat') | ||
| + | !D_max_V = max_V - 0.3709 | ||
| + | !D_min_V = min_V - 0.5155 | ||
| + | !D_max_U = max_U - 0.38290 | ||
| + | WRITE(55,*)'VARIABLES = "Iter", "Umin" , "Vmin" , "Vmax" , "Vmax_c", & | ||
| + | "Vmin_c" , "Umin_c" ' | ||
| + | WRITE(55,*)' ZONE F=POINT, I=', iii | ||
| + | open(56,file='Hystory_UVV_proc.dat') | ||
| + | WRITE(56,*)'VARIABLES = "Iter", "D_Umin_pr" , "D_Vmin_pr" , "D_Vmax_pr" ' | ||
| + | WRITE(56,*)' ZONE F=POINT, I=', iii | ||
! ------------------------------------------------ | ! ------------------------------------------------ | ||
! Open(10,file='binary_dat_all.dat') | ! Open(10,file='binary_dat_all.dat') | ||
| - | + | ! read(10,*)F(:,:,1:5) | |
| - | ! | + | |
| - | + | ||
! close(10) | ! close(10) | ||
| - | |||
! ------------------------------------------------ | ! ------------------------------------------------ | ||
Niter_last = 0 | Niter_last = 0 | ||
Niter_print = 0 | Niter_print = 0 | ||
| - | + | Do 100 N_iter=1,iii | |
| - | Do 100 | + | N_golbal_iter = N_iter |
| - | + | write(*,*) '------------------------',N_golbal_iter,'------------------------' | |
| - | + | ||
| - | write(*,*) '------------------------', | + | |
| - | + | ||
Call Solve_UV | Call Solve_UV | ||
Call Solve_Pressure_Correction | Call Solve_Pressure_Correction | ||
| - | + | !-------------------------------------------------------------------------- | |
| - | + | DO 7782 ij=1,NYmaxC | |
| + | if(Xc(ij,N_centr) < 0.5) then | ||
| + | max_V = MAX(max_V,F(ij,N_centr,2)) | ||
| + | end if | ||
| + | if(Xc(ij,N_centr) > 0.5) then | ||
| + | min_V = MAX(min_V,ABS(F(ij,N_centr,2))) | ||
| + | end if | ||
| + | if(Yc(N_centr,ij) < 0.5) then | ||
| + | max_U = MAX(max_U,ABS(F(N_centr,ij,1))) | ||
| + | end if | ||
| + | 7782 continue | ||
| + | !------------------------------------------------------------ | ||
| + | write(55,*) N_golbal_iter, -1.*max_U, -1.*min_V, max_V & | ||
| + | , ' 0.3709 -0.5155 -0.38290' | ||
| + | D_max_V = max_V - 0.3709 | ||
| + | D_min_V = min_V - 0.5155 | ||
| + | D_max_U = max_U - 0.38290 | ||
| + | D_max_V_pr = ( D_max_V / 0.3709 ) * 100. | ||
| + | D_min_V_pr = ( D_min_V / 0.5155 ) * 100. | ||
| + | D_max_U_pr = ( D_max_U / 0.38290) * 100. | ||
| + | write(56,*) N_golbal_iter, D_max_U_pr, D_min_V_pr, D_max_V_pr | ||
!-------------------------------------------------------------------------- | !-------------------------------------------------------------------------- | ||
| - | + | if((N_iter-Niter_print).eq.500)then | |
| - | + | Niter_print = N_iter | |
| - | + | Open(48,file='binary_dat_all.dat') | |
| - | + | write(48,*)F(:,:,1:5) | |
| - | open ( | + | close(48) |
| - | + | open (47,file='Domain_all.dat') | |
| - | + | WRITE(47,*)'VARIABLES = "Xp", "Yp" , "Up" , "Vp" , "Pp" ' | |
| - | + | WRITE (47,*)' ZONE I=' ,NXmaxC, ', J=', NYmaxC, ', F=POINT' | |
| - | + | ||
| - | + | ||
DO 88 J=1, NYmaxC | DO 88 J=1, NYmaxC | ||
DO 88 I=1, NXmaxC | DO 88 I=1, NXmaxC | ||
| - | + | WRITE (47,*) Xc(I,J), Yc(I,J) , F(i,j,1) , F(i,j,2) , F(i,j,4) | |
| - | + | ||
| - | + | ||
88 continue | 88 continue | ||
| - | + | close(47) | |
| - | close( | + | !------------------------------------------------------ |
| - | !------------------------------------------------------ | + | end if |
| + | !************************************************************************* | ||
| + | 100 continue | ||
open (23,file='Profiles_V.dat') | open (23,file='Profiles_V.dat') | ||
| - | + | WRITE(23,*)'VARIABLES = "Xcalc", "Vcalc" ' | |
| - | + | WRITE (23,*)' ZONE F=POINT, I=', NYmaxC | |
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
DO 881 i=1,NYmaxC | DO 881 i=1,NYmaxC | ||
| - | + | WRITE (23,*) Xc(i,N_centr), F(i,N_centr,2)*0.5 + 0.5 | |
| - | + | 881 continue | |
| - | + | WRITE (23,*)' ZONE F=POINT, I=', 17 | |
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
DO 882 i=1,17 | DO 882 i=1,17 | ||
| - | + | WRITE (23,*) x_exct(i), v_exct(i)*0.5 + 0.5 | |
| - | + | ||
| - | + | ||
882 continue | 882 continue | ||
| - | |||
| - | |||
close(23) | close(23) | ||
!-------------------------------------------------------- | !-------------------------------------------------------- | ||
| - | + | open (24,file='Profiles_U.dat') | |
| - | + | WRITE(24,*)'VARIABLES = "Ycalc", "Ucalc" ' | |
| - | + | WRITE (24,*)' ZONE F=POINT, I=', NYmaxC | |
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
DO 871 j=1,NYmaxC | DO 871 j=1,NYmaxC | ||
| - | + | WRITE (24,*) Yc(N_centr,j), F(N_centr,j,1)*0.5 + 0.5 | |
| - | + | 871 continue | |
| - | + | WRITE (24,*)' ZONE F=POINT, I=', 17 | |
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
DO 872 j=1,17 | DO 872 j=1,17 | ||
| - | + | WRITE (24,*) y_exct(j), u_exct(j)*0.5 + 0.5 | |
| - | + | 872 continue | |
| - | + | close(24) | |
| - | + | close(55) ! hystory-1 file closed | |
| - | + | close(56) ! hystory-2 file closed | |
| - | + | !********************************************************************** | |
| - | + | ||
| - | + | ||
| - | !- | + | |
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
Open(10,file='binary_dat_all.dat') | Open(10,file='binary_dat_all.dat') | ||
| - | |||
write(10,*)F(:,:,1:5) | write(10,*)F(:,:,1:5) | ||
close(10) | close(10) | ||
| - | |||
| - | |||
| - | |||
| - | |||
| - | |||
| - | |||
| - | |||
| - | |||
| - | |||
| - | |||
| - | |||
| - | |||
| - | |||
!---------------------------------------------------------------- | !---------------------------------------------------------------- | ||
NImax = NXmaxC | NImax = NXmaxC | ||
NJmax = NYmaxC | NJmax = NYmaxC | ||
| - | |||
F_out = F(:,:,1) | F_out = F(:,:,1) | ||
| - | |||
! 1234567890 | ! 1234567890 | ||
Filename ='1_U_s.txt' | Filename ='1_U_s.txt' | ||
| - | + | ! Call Out_array(F_out,NImax,NJmax,Filename) | |
| - | + | ||
!------------------------------------------------------------------- | !------------------------------------------------------------------- | ||
| - | |||
NImax = NXmaxC | NImax = NXmaxC | ||
NJmax = NYmaxC | NJmax = NYmaxC | ||
| - | |||
F_out = F(:,:,2) | F_out = F(:,:,2) | ||
| - | |||
! 1234567890 | ! 1234567890 | ||
Filename ='1_V_s.txt' | Filename ='1_V_s.txt' | ||
| - | + | ! Call Out_array(F_out,NImax,NJmax,Filename) | |
| - | + | ||
| - | + | ||
!---------------------------------------------------------------- | !---------------------------------------------------------------- | ||
NImax = NXmaxC | NImax = NXmaxC | ||
NJmax = NYmaxC | NJmax = NYmaxC | ||
| - | |||
F_out = F(:,:,5) | F_out = F(:,:,5) | ||
| - | |||
! 1234567890 | ! 1234567890 | ||
Filename ='1_T_s.txt' | Filename ='1_T_s.txt' | ||
| - | |||
! Call Out_array(F_out,NImax,NJmax,Filename) | ! Call Out_array(F_out,NImax,NJmax,Filename) | ||
!------------------------------------------------------------------- | !------------------------------------------------------------------- | ||
| - | |||
| - | |||
open (23,file='Domain_all.dat') | open (23,file='Domain_all.dat') | ||
| - | |||
WRITE(23,*)'VARIABLES = "Xp", "Yp" , "Up" , "Vp" , "Pp" ' | WRITE(23,*)'VARIABLES = "Xp", "Yp" , "Up" , "Vp" , "Pp" ' | ||
| - | + | WRITE (23,*)' ZONE I=' ,NXmaxC, ', J=', NYmaxC, ', F=POINT' | |
| - | + | ||
| - | + | ||
DO 4 J=1, NYmaxC | DO 4 J=1, NYmaxC | ||
DO 4 I=1, NXmaxC | DO 4 I=1, NXmaxC | ||
| - | + | WRITE (23,*) Xc(I,J), Yc(I,J) , F(i,j,1) , F(i,j,2) , F(i,j,4) | |
| - | + | 4 continue | |
| - | + | ||
| - | + | ||
| - | + | ||
close(23) | close(23) | ||
!-------------------------------------------------------------------------- | !-------------------------------------------------------------------------- | ||
| - | |||
WRITE(*,*) 'PRIVET' | WRITE(*,*) 'PRIVET' | ||
STOP | STOP | ||
END | END | ||
| - | |||
</pre> | </pre> | ||
Revision as of 14:29, 19 May 2016
!Sample program for solving Lid-Driven cavity test using SIMPLE-algorithm
! Main modul
!Copyright (C) 2010 Michail Kiričkov
! 2016 - I fixed some bugs
!This program is free software; you can redistribute it and/or
!modify it under the terms of the GNU General Public License
!as published by the Free Software Foundation; either version 2
!of the License, or (at your option) any later version.
!This program is distributed in the hope that it will be useful,
!but WITHOUT ANY WARRANTY; without even the implied warranty of
!MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
!GNU General Public License for more details.
!You should have received a copy of the GNU General Public License
!along with this program; if not, write to the Free Software
!Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
!-----------------------------------------------------------------------
Program Main
include 'icomm_1.f90'
Dimension F_out(nx,ny),x_exct(30), v_exct(30)
Dimension y_exct(30), u_exct(30)
Character Filename*10
real qqq
iii = 150
Call Ini
Call Grid_rectangular
Call Geom
Call Init_all_cavity
write(*,*) Yc(1,N_centr)
read(*,*) qqq
x_exct(1)=0.0000 ; v_exct(1)= 0.00000
x_exct(2)=0.0625 ; v_exct(2)= 0.27485
x_exct(3)=0.0703 ; v_exct(3)= 0.29012
x_exct(4)=0.0781 ; v_exct(4)= 0.30353
x_exct(5)=0.0938 ; v_exct(5)= 0.32627
x_exct(6)=0.1563 ; v_exct(6)= 0.37095
x_exct(7)=0.2266 ; v_exct(7)= 0.33075
x_exct(8)=0.2344 ; v_exct(8)= 0.32235
x_exct(9)=0.5000 ; v_exct(9)= 0.02526
x_exct(10)=0.8047 ; v_exct(10)=-0.31966
x_exct(11)=0.8594 ; v_exct(11)=-0.42665
x_exct(12)=0.9063 ; v_exct(12)=-0.51550
x_exct(13)=0.9453 ; v_exct(13)=-0.39188
x_exct(14)=0.9531 ; v_exct(14)=-0.33714
x_exct(15)=0.9609 ; v_exct(15)=-0.27669
x_exct(16)=0.9688 ; v_exct(16)=-0.21388
x_exct(17)=1.0000 ; v_exct(17)=-0.00000
y_exct(1)=0.00000 ; u_exct(1)=-0.00000
y_exct(2)=0.0547 ; u_exct(2)=-0.18109
y_exct(3)=0.0625 ; u_exct(3)=-0.20196
y_exct(4)=0.0703 ; u_exct(4)=-0.22220
y_exct(5)=0.1016 ; u_exct(5)=-0.29730
y_exct(6)=0.1719 ; u_exct(6)=-0.38289
y_exct(7)=0.2813 ; u_exct(7)=-0.27805
y_exct(8)=0.4531 ; u_exct(8)=-0.10648
y_exct(9)=0.5000 ; u_exct(9)=-0.06080
y_exct(10)=0.6172 ; u_exct(10)=0.05702
y_exct(11)=0.7344 ; u_exct(11)=0.18719
y_exct(12)=0.8516 ; u_exct(12)=0.33304
y_exct(13)=0.9531 ; u_exct(13)=0.46604
y_exct(14)=0.9609 ; u_exct(14)=0.51117
y_exct(15)=0.9688 ; u_exct(15)=0.57492
y_exct(16)=0.9766 ; u_exct(16)=0.65928
y_exct(17)=1.0000 ; u_exct(17)=1.00000
open(55,file='Hystory_UVV.dat')
!D_max_V = max_V - 0.3709
!D_min_V = min_V - 0.5155
!D_max_U = max_U - 0.38290
WRITE(55,*)'VARIABLES = "Iter", "Umin" , "Vmin" , "Vmax" , "Vmax_c", &
"Vmin_c" , "Umin_c" '
WRITE(55,*)' ZONE F=POINT, I=', iii
open(56,file='Hystory_UVV_proc.dat')
WRITE(56,*)'VARIABLES = "Iter", "D_Umin_pr" , "D_Vmin_pr" , "D_Vmax_pr" '
WRITE(56,*)' ZONE F=POINT, I=', iii
! ------------------------------------------------
! Open(10,file='binary_dat_all.dat')
! read(10,*)F(:,:,1:5)
! close(10)
! ------------------------------------------------
Niter_last = 0
Niter_print = 0
Do 100 N_iter=1,iii
N_golbal_iter = N_iter
write(*,*) '------------------------',N_golbal_iter,'------------------------'
Call Solve_UV
Call Solve_Pressure_Correction
!--------------------------------------------------------------------------
DO 7782 ij=1,NYmaxC
if(Xc(ij,N_centr) < 0.5) then
max_V = MAX(max_V,F(ij,N_centr,2))
end if
if(Xc(ij,N_centr) > 0.5) then
min_V = MAX(min_V,ABS(F(ij,N_centr,2)))
end if
if(Yc(N_centr,ij) < 0.5) then
max_U = MAX(max_U,ABS(F(N_centr,ij,1)))
end if
7782 continue
!------------------------------------------------------------
write(55,*) N_golbal_iter, -1.*max_U, -1.*min_V, max_V &
, ' 0.3709 -0.5155 -0.38290'
D_max_V = max_V - 0.3709
D_min_V = min_V - 0.5155
D_max_U = max_U - 0.38290
D_max_V_pr = ( D_max_V / 0.3709 ) * 100.
D_min_V_pr = ( D_min_V / 0.5155 ) * 100.
D_max_U_pr = ( D_max_U / 0.38290) * 100.
write(56,*) N_golbal_iter, D_max_U_pr, D_min_V_pr, D_max_V_pr
!--------------------------------------------------------------------------
if((N_iter-Niter_print).eq.500)then
Niter_print = N_iter
Open(48,file='binary_dat_all.dat')
write(48,*)F(:,:,1:5)
close(48)
open (47,file='Domain_all.dat')
WRITE(47,*)'VARIABLES = "Xp", "Yp" , "Up" , "Vp" , "Pp" '
WRITE (47,*)' ZONE I=' ,NXmaxC, ', J=', NYmaxC, ', F=POINT'
DO 88 J=1, NYmaxC
DO 88 I=1, NXmaxC
WRITE (47,*) Xc(I,J), Yc(I,J) , F(i,j,1) , F(i,j,2) , F(i,j,4)
88 continue
close(47)
!------------------------------------------------------
end if
!*************************************************************************
100 continue
open (23,file='Profiles_V.dat')
WRITE(23,*)'VARIABLES = "Xcalc", "Vcalc" '
WRITE (23,*)' ZONE F=POINT, I=', NYmaxC
DO 881 i=1,NYmaxC
WRITE (23,*) Xc(i,N_centr), F(i,N_centr,2)*0.5 + 0.5
881 continue
WRITE (23,*)' ZONE F=POINT, I=', 17
DO 882 i=1,17
WRITE (23,*) x_exct(i), v_exct(i)*0.5 + 0.5
882 continue
close(23)
!--------------------------------------------------------
open (24,file='Profiles_U.dat')
WRITE(24,*)'VARIABLES = "Ycalc", "Ucalc" '
WRITE (24,*)' ZONE F=POINT, I=', NYmaxC
DO 871 j=1,NYmaxC
WRITE (24,*) Yc(N_centr,j), F(N_centr,j,1)*0.5 + 0.5
871 continue
WRITE (24,*)' ZONE F=POINT, I=', 17
DO 872 j=1,17
WRITE (24,*) y_exct(j), u_exct(j)*0.5 + 0.5
872 continue
close(24)
close(55) ! hystory-1 file closed
close(56) ! hystory-2 file closed
!**********************************************************************
Open(10,file='binary_dat_all.dat')
write(10,*)F(:,:,1:5)
close(10)
!----------------------------------------------------------------
NImax = NXmaxC
NJmax = NYmaxC
F_out = F(:,:,1)
! 1234567890
Filename ='1_U_s.txt'
! Call Out_array(F_out,NImax,NJmax,Filename)
!-------------------------------------------------------------------
NImax = NXmaxC
NJmax = NYmaxC
F_out = F(:,:,2)
! 1234567890
Filename ='1_V_s.txt'
! Call Out_array(F_out,NImax,NJmax,Filename)
!----------------------------------------------------------------
NImax = NXmaxC
NJmax = NYmaxC
F_out = F(:,:,5)
! 1234567890
Filename ='1_T_s.txt'
! Call Out_array(F_out,NImax,NJmax,Filename)
!-------------------------------------------------------------------
open (23,file='Domain_all.dat')
WRITE(23,*)'VARIABLES = "Xp", "Yp" , "Up" , "Vp" , "Pp" '
WRITE (23,*)' ZONE I=' ,NXmaxC, ', J=', NYmaxC, ', F=POINT'
DO 4 J=1, NYmaxC
DO 4 I=1, NXmaxC
WRITE (23,*) Xc(I,J), Yc(I,J) , F(i,j,1) , F(i,j,2) , F(i,j,4)
4 continue
close(23)
!--------------------------------------------------------------------------
WRITE(*,*) 'PRIVET'
STOP
END
