diff --git a/LES_Filtering_eig33.f90 b/LES_Filtering_eig33.f90 index 1d80e1e..15a7e65 100644 --- a/LES_Filtering_eig33.f90 +++ b/LES_Filtering_eig33.f90 @@ -20,9 +20,9 @@ SUBROUTINE EIG33(A,eig) COMPLEX(KIND=8) :: eig(3) CALL MMDOT(A,A,B,3) - !---------------------------------------------------------------------! - ! Trace function of 3 X 3 matrix ! - !----------------------------------------------------------------------! + !--------------------------------------------------------------------! + ! Trace function of 3 X 3 matrix ! + !--------------------------------------------------------------------! a0 = det(A) a1 = 0.5*(tr(B) - tr(A)**2) a2 = tr(A) diff --git a/LES_Filtering_module.f90 b/LES_Filtering_module.f90 index 199ae06..f09f88d 100644 --- a/LES_Filtering_module.f90 +++ b/LES_Filtering_module.f90 @@ -11,11 +11,12 @@ MODULE LES_FILTERING_module IMPLICIT NONE - INTEGER :: Nx, Ny, Nz, Nx_fil, Nz_fil, VS_CASE, FILTER_OX + INTEGER :: Nx, Ny, Nz, Nx_fil, Nz_fil, VS_CASE, FILTER_OX, VS_ONLY, & + Y_ORDER REAL(KIND=8) :: Del,dx,dz,FW,pi,tol CHARACTER(LEN=65) :: file_name, dir_name, path_name - REAL(KIND=8),DIMENSION(:),ALLOCATABLE :: X,Y,Z,dy + REAL(KIND=8),DIMENSION(:),ALLOCATABLE :: X,Y,Z,dy,YP REAL(KIND=8),DIMENSION(:,:,:),ALLOCATABLE :: U,V,W,U_Fil,V_Fil,W_Fil, & U_Fil_2,V_Fil_2,W_Fil_2, & VS,Cs diff --git a/LES_Filtering_output.f90 b/LES_Filtering_output.f90 index 70572b3..d40937c 100644 --- a/LES_Filtering_output.f90 +++ b/LES_Filtering_output.f90 @@ -12,47 +12,49 @@ SUBROUTINE OUTPUT ONLY : J_det USE LES_FILTERING_module, & - ONLY : Nx, Ny, Nz, file_name, dir_name, path_name, VS_CASE + ONLY : Nx, Ny, Nz, file_name, dir_name, path_name, VS_CASE, & + VS_ONLY USE LES_FILTERING_module, & ONLY : X, Y, Z, U, V, W, U_Fil, V_Fil, W_Fil, dy, Resi_T, & - S_T, S_T_Fil, NU_R, Cs, VS + S_T, S_T_Fil, NU_R, Cs, VS, YP IMPLICIT NONE INTEGER :: i,j,k,it,J_loc,v_i,v_j REAL(KIND=8) :: time_sta, time_end, U_ave, V_ave, W_ave, Cs_ave - REAL(KIND=8) :: YP(1:3), RESI_ave(1:3,1:3), NU_R_ave(1:3,1:3) & + REAL(KIND=8) :: RESI_ave(1:3,1:3), NU_R_ave(1:3,1:3) & ,S_T_ave(1:3,1:3,1:2) WRITE(*,*) '----------------------------------------------------' WRITE(*,*) ' WRITING PROCESS STARTED ' CALL CPU_TIME(time_sta) + IF (VS_ONLY - 2 < 0) THEN + !----------------------------------------------------------------! - ! Outputs for U_Fil(Filtered velocity) + ! Outputs for U_Fil(Filtered velocity) ! !----------------------------------------------------------------! dir_name = 'RESULT/U' - ! file_name = '/U_filtered.plt' - ! path_name = TRIM(dir_name)//TRIM(file_name) - ! OPEN(100,FILE=path_name,FORM='FORMATTED',POSITION='APPEND') - ! WRITE(100,*) 'VARIABLES = X,Y,Z,U_fil,V_Fil,W_Fil' - ! WRITE(100,"(3(A,I3,2X))")' ZONE I = ',Nx,' J = ',Ny, ' K = ', Nz - ! DO k = 1,Nz - ! DO j = 1,Ny - ! DO i = 1,Nx - ! WRITE(100,"(6F15.9)") X(i),Y(j),Z(k), & - ! U_fil(i,j,k),V_Fil(i,j,k),W_Fil(i,j,k) - ! END DO - ! END DO - ! END DO - ! CLOSE(100) + file_name = '/U_filtered.plt' + path_name = TRIM(dir_name)//TRIM(file_name) + + OPEN(100,FILE=path_name,FORM='FORMATTED',POSITION='APPEND') + WRITE(100,*) 'VARIABLES = X,Y,Z,U_fil,V_Fil,W_Fil' + WRITE(100,"(3(A,I3,2X))")' ZONE I = ',Nx,' J = ',Ny, ' K = ', Nz + DO k = 1,Nz + DO j = 1,Ny + DO i = 1,Nx + WRITE(100,"(6F15.9)") X(i),Y(j),Z(k), & + U_fil(i,j,k),V_Fil(i,j,k),W_Fil(i,j,k) + END DO + END DO + END DO + CLOSE(100) !----------------------------------------------------------------! ! Outputs for instaneous U at Y+ = 5,30,200 ! !----------------------------------------------------------------! - YP = [5,30,200] - DO it = 1,3 WRITE(file_name,"(I3.3,A)")INT(YP(it)),'.U_slice.plt' @@ -106,8 +108,6 @@ SUBROUTINE OUTPUT !----------------------------------------------------------------! dir_name = 'RESULT/RESIDUAL_STRESS' - YP = [5,30,200] - DO it = 1,3 WRITE(file_name,"(I3.3,A)")INT(YP(it)),'.Resi_slice.plt' @@ -161,8 +161,6 @@ SUBROUTINE OUTPUT !----------------------------------------------------------------! dir_name = 'RESULT/STRAIN_RATE' - YP = [5,30,200] - DO it = 1,3 WRITE(file_name,"(I3.3,A)")INT(YP(it)),'.S_T_slice.plt' @@ -221,13 +219,13 @@ SUBROUTINE OUTPUT END DO CLOSE(100) + WRITE(*,*) "S_T_AV COMPLETED" + !----------------------------------------------------------------! ! Outputs for eddy-viscosity at Y+ = 5,30,200 ! !----------------------------------------------------------------! dir_name = 'RESULT/EDDY_VISCOSITY' - YP = [5,30,200] - DO it = 1,3 WRITE(file_name,"(I3.3,A)")INT(YP(it)),'.NU_R_slice.plt' @@ -282,15 +280,11 @@ SUBROUTINE OUTPUT END DO CLOSE(100) - CALL CPU_TIME(time_end) - !----------------------------------------------------------------! ! Outputs for Smagorinsky coefficient at Y+ = 5,30,200 ! !----------------------------------------------------------------! dir_name = 'RESULT/SMARGORINSKY_COEFFICIENT' - YP = [5,30,200] - DO it = 1,3 WRITE(file_name,"(I3.3,A)")INT(YP(it)),'.CS_slice.plt' @@ -332,9 +326,11 @@ SUBROUTINE OUTPUT END DO CLOSE(100) + END IF !----------------------------------------------------------------! ! Outputs for Vortical Structures ! !----------------------------------------------------------------! + IF ( mod(VS_ONLY,2) == 0 ) THEN dir_name = 'RESULT/VORTICAL_STRUCTURE' file_name = '/VORTICAL_STRUCTURE.plt' @@ -358,6 +354,8 @@ SUBROUTINE OUTPUT END DO CLOSE(100) + END IF + CALL CPU_TIME(time_end) WRITE(*,*) ' WRITING PROCESS IS COMPLETED ' @@ -365,6 +363,6 @@ SUBROUTINE OUTPUT WRITE(*,*) '----------------------------------------------------' WRITE(*,*) '' - DEALLOCATE(X,Y,Z,U,V,W,U_Fil,V_Fil,W_Fil,dy,S_T,S_T_Fil) + DEALLOCATE(X,Y,Z,U,V,W,U_Fil,V_Fil,W_Fil,dy,S_T,S_T_Fil,YP) END SUBROUTINE OUTPUT diff --git a/LES_Filtering_read.f90 b/LES_Filtering_read.f90 index 3ffe890..e70b1c4 100644 --- a/LES_Filtering_read.f90 +++ b/LES_Filtering_read.f90 @@ -15,7 +15,7 @@ SUBROUTINE READ_DNS USE LES_FILTERING_module, & ONLY : Nx, Ny, Nz, dx, dz, Del, FW, tol, Nx_fil, Nz_fil, & - file_name, dir_name, path_name + file_name, dir_name, path_name, Y_order USE LES_FILTERING_module, & ONLY : X, Y, Z, U, V, W, dy, S_T, O_T @@ -41,18 +41,35 @@ SUBROUTINE READ_DNS !--------------------------------------------------------------------! ! Main loop of reading DNS datas ! !--------------------------------------------------------------------! - DO k = 1,Nz - DO j = Ny,1,-1 - DO i = 1,Nx - READ(100,*) tmp_x, tmp_y, tmp_z, U(i,j,k), V(i,j,k), W(i,j,k) - - IF (j==1 .AND. k==1) X(i) = tmp_x - IF (i==1 .AND. k==1) Y(j) = tmp_y - IF (i==1 .AND. j==1) Z(k) = tmp_z + IF ( Y_ORDER == 0 ) THEN + + DO k = 1,Nz + DO j = 1,Ny + DO i = 1,Nx + READ(100,*) tmp_x, tmp_y, tmp_z, U(i,j,k), V(i,j,k), W(i,j,k) + + IF (j==1 .AND. k==1) X(i) = tmp_x + IF (i==1 .AND. k==1) Y(j) = tmp_y + IF (i==1 .AND. j==1) Z(k) = tmp_z + END DO + END DO + END DO + + ELSE + + DO k = 1,Nz + DO j = Ny,1,-1 + DO i = 1,Nx + READ(100,*) tmp_x, tmp_y, tmp_z, U(i,j,k), V(i,j,k), W(i,j,k) + + IF (j==1 .AND. k==1) X(i) = tmp_x + IF (i==1 .AND. k==1) Y(j) = tmp_y + IF (i==1 .AND. j==1) Z(k) = tmp_z + END DO END DO END DO - END DO + END IF CLOSE(100) !--------------------------------------------------------------------! diff --git a/LES_Filtering_setup.f90 b/LES_Filtering_setup.f90 index c85d908..c970a36 100644 --- a/LES_Filtering_setup.f90 +++ b/LES_Filtering_setup.f90 @@ -12,22 +12,24 @@ SUBROUTINE SETUP USE LES_FILTERING_module, & ONLY : Nx, Ny, Nz, dx, dz, FW, pi, tol, & - file_name, dir_name, path_name, VS_CASE, FILTER_OX + file_name, dir_name, path_name, VS_CASE, FILTER_OX, & + VS_ONLY, Y_ORDER USE LES_FILTERING_module, & ONLY : X, Y, Z, dy, U, V, W, U_Fil, V_Fil, W_Fil, & Resi_T, S_T, S_T_Fil, NU_R, O_T, O_T_Fil, VS, & - L_T, M_T, U_Fil_2, V_Fil_2, W_Fil_2, S_T_Fil_2, Cs + L_T, M_T, U_Fil_2, V_Fil_2, W_Fil_2, S_T_Fil_2, Cs, YP IMPLICIT NONE - INTEGER :: i,j,k + INTEGER :: i,j,k,N pi = atan(1.0)*4 !------------------------------------------------------------------! ! Make & Initialize Result folder ! !------------------------------------------------------------------! - file_name = 'instantaneous_velocity_field_re644.plt' + ! file_name = 'instantaneous_velocity_field_re644.plt' + file_name = 'INSU_XYZ.plt' dir_name = 'RESULT' CALL SYSTEM('mkdir '//TRIM(dir_name)) @@ -48,6 +50,16 @@ SUBROUTINE SETUP CALL SYSTEM('rm -rf ./'//TRIM(dir_name)//'/VORTICAL_STRUCTURE'//'/*.plt') CALL SYSTEM('rm -rf ./'//TRIM(dir_name)//'/SMARGORINSKY_COEFFICIENT'//'/*.plt') + !------------------------------------------------------------------! + ! Statistic type ! + ! ! + ! (a) All Statistic : 0 ! + ! (b) Excpet Vortical Structure : 1 ! + ! (c) Only Vortical Structure : 2 ! + ! ! + !------------------------------------------------------------------! + VS_ONLY = 2 + !------------------------------------------------------------------! ! Vortical Structure methods ! ! ! @@ -56,7 +68,7 @@ SUBROUTINE SETUP ! (c) Lambda_ci criteria : 3 ! ! ! !------------------------------------------------------------------! - VS_CASE = 1 + VS_CASE = 2 !------------------------------------------------------------------! ! The number of Filter ! @@ -66,14 +78,34 @@ SUBROUTINE SETUP ! (c) Second Filter : 2 ! ! ! !------------------------------------------------------------------! - FILTER_OX = 1 + FILTER_OX = 0 + + !------------------------------------------------------------------! + ! Y sorting order ! + ! ! + ! (a) Ascending Order : 0 ! + ! (b) Descending Order : 1 ! + ! ! + !------------------------------------------------------------------! + Y_ORDER = 0 + + !------------------------------------------------------------------! + ! Statistic Slice Point ! + ! ! + ! N : Total Number of Slice Point ! + ! YP : Slice Point Position ! + ! ! + !------------------------------------------------------------------! + N = 3 + ALLOCATE( YP(1:N) ) + YP = [5,30,200] !------------------------------------------------------------------! ! Constants for LES filtering ! !------------------------------------------------------------------! - Nx = 288 - Ny = 257 - Nz = 288 + Nx = 128 + Ny = 191 + Nz = 159 FW = 4 ! Filter width constant tol = 1e-8 ! Tolerance for the number of nodes in x,z directions diff --git a/Makefile b/Makefile index 79044a3..2def6cb 100644 --- a/Makefile +++ b/Makefile @@ -1,5 +1,5 @@ F90=ifort -FCFLAGS=-O3 -qopenmp -I${MKL_HOME}/include -L${MKL_HOME}/lib -lmkl_intel_lp64 \ +FCFLAGS=-O3 -qopenmp#-I${MKL_HOME}/include/intel64/ilp64 -L${MKL_HOME}/lib -lmkl_intel_lp64 \ -lmkl_intel_thread -lmkl_core -lpthread TARGET= LES_Filtering diff --git a/Vortical_Structure.f90 b/Vortical_Structure.f90 index 4cf2b9e..441dc1b 100644 --- a/Vortical_Structure.f90 +++ b/Vortical_Structure.f90 @@ -87,7 +87,7 @@ SUBROUTINE VORTICAL_STRUCTURE ! LWORK = MIN( LWMAX, INT( WORK( 1 ) ) ) ! CALL DSYEV('N','U',3,D_T_tmp,3,EIG_R,WORK,LWORK,INFO) ! VS(i,j,k) = EIG_R(2) - + ! WRITE(*,"(3I,6F15.9)") i,j,k,REAL(eig(1:3)),AIMAG(eig(1:3)) CALL EIG33(D_T_tmp,eig) VS(i,j,k) = REAL(eig(2))