97 lines
2.2 KiB
Fortran
97 lines
2.2 KiB
Fortran
subroutine pressure
|
|
|
|
use m_parameters
|
|
use m_fields
|
|
use m_work
|
|
use x_fftw
|
|
|
|
implicit none
|
|
|
|
integer :: i, j, k, n
|
|
real*8 :: div1, div2, lapl1, lapl2, p1, p2
|
|
|
|
do k = 1,nz
|
|
do j = 1,ny
|
|
do i = 1,nx+1,2
|
|
|
|
! getting the divergence (i*k*\hat(u))
|
|
! remember that in the Fourier space indicies go as (ix,iz,iy)
|
|
|
|
! getting divergence
|
|
div2 = akx(i )*fields(i ,j,k,1) + aky(k)*fields(i ,j,k,2) + akz(j)*fields(i ,j,k,3)
|
|
div1 = - ( akx(i+1)*fields(i+1,j,k,1) + aky(k)*fields(i+1,j,k,2) + akz(j)*fields(i+1,j,k,3) )
|
|
|
|
! inverce laplace operator
|
|
lapl1 = akx(i )**2 + aky(k)**2 + akz(j)**2
|
|
lapl2 = akx(i+1)**2 + aky(k)**2 + akz(j)**2
|
|
|
|
if (lapl1.eq.0.d0) lapl1 = 9e20
|
|
if (lapl2.eq.0.d0) lapl2 = 9e20
|
|
|
|
! calculating pressure
|
|
p1 = - div1 / lapl1
|
|
p2 = - div2 / lapl2
|
|
|
|
! Taking derivatives of the pressure and subtracting from the corresponding velocities
|
|
fields(i ,j,k,1) = fields(i ,j,k,1) + p2 * akx(i+1)
|
|
fields(i+1,j,k,1) = fields(i+1,j,k,1) - p1 * akx(i )
|
|
|
|
fields(i ,j,k,2) = fields(i ,j,k,2) + p2 * aky(k)
|
|
fields(i+1,j,k,2) = fields(i+1,j,k,2) - p1 * aky(k)
|
|
|
|
fields(i ,j,k,3) = fields(i ,j,k,3) + p2 * akz(j)
|
|
fields(i+1,j,k,3) = fields(i+1,j,k,3) - p1 * akz(j)
|
|
|
|
end do
|
|
end do
|
|
end do
|
|
|
|
|
|
return
|
|
end subroutine pressure
|
|
|
|
|
|
|
|
subroutine divergence
|
|
|
|
use m_openmpi
|
|
use m_parameters
|
|
use m_fields
|
|
use m_work
|
|
use x_fftw
|
|
|
|
implicit none
|
|
|
|
integer :: i,j,k
|
|
real*8 :: dmin, dmax, d1
|
|
|
|
|
|
wrk(:,:,:,1:3) = fields(:,:,:,1:3)
|
|
|
|
call x_derivative(1,'x',4)
|
|
call x_derivative(2,'y',5)
|
|
call x_derivative(3,'z',6)
|
|
|
|
call xFFT3d(-1,4)
|
|
call xFFT3d(-1,5)
|
|
call xFFT3d(-1,6)
|
|
|
|
wrk(:,:,:,0) = wrk(:,:,:,4) + wrk(:,:,:,5) + wrk(:,:,:,6)
|
|
|
|
d1 = minval(wrk(1:nx,:,:,0))
|
|
call MPI_REDUCE(d1,dmin,1,MPI_REAL8,MPI_MIN,0,MPI_COMM_TASK,mpi_err)
|
|
d1 = maxval(wrk(1:nx,:,:,0))
|
|
call MPI_REDUCE(d1,dmax,1,MPI_REAL8,MPI_MAX,0,MPI_COMM_TASK,mpi_err)
|
|
|
|
if (myid.eq.0) then
|
|
write(out,*)
|
|
write(out,'(a30,3x,2e20.10)') 'divergence: min, max',dmin,dmax
|
|
write(out,*)
|
|
call flush(out)
|
|
end if
|
|
|
|
return
|
|
end subroutine divergence
|
|
|
|
|
|
|