93 allocate(
aa(zsize(1),zsize(2),zsize(3)+2*
ghst) )
94 allocate(
bb(zsize(1),zsize(2),zsize(3)+2*
ghst) )
95 allocate(
cc(zsize(1),zsize(2),zsize(3)+2*
ghst) )
97 allocate(
bb_u(zsize(1),zsize(2),zsize(3)+2*
ghst) )
100 allocate(
bb_w(zsize(1),zsize(2),zsize(3)+2*
ghst) )
103 allocate(
bb_p(zsize(1),zsize(2),zsize(3)+2*
ghst) )
106 allocate(
bb_t(zsize(1),zsize(2),zsize(3)+2*
ghst) )
111 allocate(
c_src(zsize(1),zsize(2),zsize(3)) )
112 allocate(
c_var(zsize(1),zsize(2),zsize(3)) )
114 allocate(
c_tmp1(xsize(1),xsize(2),xsize(3)) )
115 allocate(
c_tmp2(ysize(1),ysize(2),ysize(3)) )
120 allocate(
c_u(zsize(1),zsize(2),zsize(3)) )
121 allocate(
c_v(zsize(1),zsize(2),zsize(3)) )
122 allocate(
c_w(zsize(1),zsize(2),zsize(3)) )
123 allocate(
c_p(zsize(1),zsize(2),zsize(3)) )
124 allocate(
c_phi(zsize(1),zsize(2),zsize(3)) )
125 allocate(
c_t(zsize(1),zsize(2),zsize(3)) )
127 allocate(
c_qu(zsize(1),zsize(2),zsize(3)) )
128 allocate(
c_qv(zsize(1),zsize(2),zsize(3)) )
129 allocate(
c_qw(zsize(1),zsize(2),zsize(3)) )
130 allocate(
c_ff(zsize(1),zsize(2),zsize(3)) )
131 allocate(
c_gg(zsize(1),zsize(2),zsize(3)) )
132 allocate(
c_hh(zsize(1),zsize(2),zsize(3)) )
133 allocate(
c_ss(zsize(1),zsize(2),zsize(3)) )
177 bb(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
180 bb_u(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
183 bb_w(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
186 bb_p(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
189 bb_t(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
206 c_src(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
207 c_var(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
208 c_tmp1(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
209 c_tmp2(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
213 c_u(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
214 c_v(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
215 c_w(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
216 c_t(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
217 c_p(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
218 c_phi(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
220 c_qu(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
221 c_qv(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
222 c_qw(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
223 c_ff(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
224 c_gg(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
225 c_hh(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
226 c_ss(:,:,:)=cmplx(0.d0,0.d0,kind=mytype)
250 print *,
'restart error!' 266 real(mytype),
dimension(:,:,:),
intent(out) :: u, v, w, t, p
267 real(mytype),
dimension(:,:,:),
intent(out) :: ff,gg,hh,ss
269 if (
myid .eq. 0)
then 270 print *,
'Reading Initial Fields.' 274 'results/init_u.dat')
277 'results/init_v.dat')
280 'results/init_w.dat')
283 'results/init_p.dat')
286 'results/init_t.dat')
296 'results/init_ff.dat')
299 'results/init_gg.dat')
302 'results/init_hh.dat')
305 'results/init_ss.dat')
310 if (
myid .eq. 0)
then 311 print *,
'Reading Initial Fields Finished!' 325 real(mytype),
dimension(:,:,:),
intent(inout) :: u, v, w, t, p
330 if (
myid .eq. 0)
then 331 print *,
'Adding some noise to the initial fields' 340 call random_number(tmp)
342 call random_number(tmp)
344 call random_number(tmp)
346 call random_number(tmp)
348 call random_number(tmp)
363 (u,v,w,t,uhx,vhy,mean_u,mean_v,mean_w,mean_t, mean_uu,mean_vv, &
364 mean_ww,mean_uw,mean_vw,mean_tt,mean_tw)
368 real(mytype),
dimension(:,:,:),
intent(in) :: u,v,w,t,uhx, vhy
369 real(mytype),
dimension(:,:,:),
intent(out):: &
370 mean_u,mean_v,mean_w,mean_t, mean_uu,mean_vv,mean_ww,mean_uw, mean_vw, &
380 mean_uu(i,j,k)=mean_uu(i,j,k)+(u(i,j,k)+
u_mrf)**2.d0
381 mean_vv(i,j,k)=mean_vv(i,j,k)+v(i,j,k)**2.d0
382 mean_ww(i,j,k)=mean_ww(i,j,k)+w(i,j,k)**2.d0
383 mean_uw(i,j,k)=mean_uw(i,j,k)+0.5d0*( w(i,j,k)+w(i,j,k-1) ) * &
385 mean_vw(i,j,k)=mean_vw(i,j,k)+0.5d0*( w(i,j,k)+w(i,j,k-1) ) * vhy(i,j,k)
387 mean_tt(i,j,k)=mean_tt(i,j,k)+t(i,j,k)**2.d0
388 mean_tw(i,j,k)=mean_tw(i,j,k)+t(i,j,k) * 0.5d0 * &
389 ( w(i,j,k)+w(i,j,k-1) )
392 mean_u(i,j,k)=mean_u(i,j,k) + u(i,j,k)+
u_mrf 393 mean_v(i,j,k)=mean_v(i,j,k) + v(i,j,k)
394 mean_w(i,j,k)=mean_w(i,j,k) + w(i,j,k)
395 mean_t(i,j,k)=mean_t(i,j,k) + t(i,j,k)
410 subroutine get_interp_fields_uvwhxyz(u,v,w,uhx,uhy,uhz,vhx,vhy,vhz,whx,why)
414 real(mytype),
dimension(:,:,:),
intent(in) :: u,v,w
415 real(mytype),
dimension(:,:,:),
intent(out):: uhx,uhy,uhz,vhx,vhy,vhz,&
427 uhx(i,j,k)=0.5d0*( u(i,j,k)+u(i-1,j,k) )
428 uhy(i,j,k)=0.5d0*( u(i,j+1,k)+u(i,j,k) )
429 vhx(i,j,k)=0.5d0*( v(i+1,j,k)+v(i,j,k) )
430 vhy(i,j,k)=0.5d0*( v(i,j,k)+v(i,j-1,k) )
431 whx(i,j,k)=0.5d0*( w(i+1,j,k)+w(i,j,k) )
432 why(i,j,k)=0.5d0*( w(i,j+1,k)+w(i,j,k) )
434 elseif (
cds .eq. 2)
then 436 uhx(i,j,k)=( -u(i+1,j,k)+9.d0*u(i,j,k)+ &
437 9.d0*u(i-1,j,k)-u(i-2,j,k) ) / 16.d0
439 uhy(i,j,k)=( -u(i,j+2,k)+9.d0*u(i,j+1,k)+ &
440 9.d0*u(i,j,k)-u(i,j-1,k) ) / 16.d0
442 vhx(i,j,k)=( -v(i+2,j,k)+9.d0*v(i+1,j,k)+ &
443 9.d0*v(i,j,k)-v(i-1,j,k) ) / 16.d0
445 vhy(i,j,k)=( -v(i,j+1,k)+9.d0*v(i,j,k)+ &
446 9.d0*v(i,j-1,k)-v(i,j-2,k) ) / 16.d0
448 whx(i,j,k)=( -w(i+2,j,k)+9.d0*w(i+1,j,k)+ &
449 9.d0*w(i,j,k)-w(i-1,j,k) ) / 16.d0
451 why(i,j,k)=( -w(i,j+2,k)+9.d0*w(i,j+1,k)+ &
452 9.d0*w(i,j,k)-w(i,j-1,k) ) / 16.d0
457 uhz(i,j,k)=u(i,j,k+1)*
l_t(k)+u(i,j,k)*(1.d0-
l_t(k))
459 vhz(i,j,k)=v(i,j,k+1)*
l_t(k)+v(i,j,k)*(1.d0-
l_t(k))
477 real(mytype),
dimension(:,:,:),
intent(in) :: t
478 real(mytype),
dimension(:,:,:),
intent(out):: thx,thy,thz
489 thx(i,j,k)=0.5d0*( t(i+1,j,k)+t(i,j,k) )
490 thy(i,j,k)=0.5d0*( t(i,j+1,k)+t(i,j,k) )
492 elseif (
cds .eq. 2)
then 494 thx(i,j,k)=( -t(i+2,j,k)+9.d0*t(i+1,j,k)+ &
495 9.d0*t(i,j,k)-t(i-1,j,k) ) / 16.d0
497 thy(i,j,k)=( -t(i,j+2,k)+9.d0*t(i,j+1,k)+ &
498 9.d0*t(i,j,k)-t(i,j-1,k) ) / 16.d0
503 thz(i,j,k)=t(i,j,k+1)*
l_t(k)+t(i,j,k)*(1.d0-
l_t(k))
521 use,
intrinsic :: iso_c_binding
528 include
'include/fftw3.f03' 531 integer,
intent(in) :: direction
533 complex(mytype),
dimension(:,:,:),
intent(inout) :: c_inout3
539 if (direction .eq. 1)
then 560 if (l .eq. -
nx/2)
then 561 c_inout3(i,j,k)=cmplx(0.d0,0.d0,kind=mytype)
564 if (m .eq. -
ny/2)
then 565 c_inout3(i,j,k)=cmplx(0.d0,0.d0,kind=mytype)
577 call transpose_z_to_y(c_inout3,
c_tmp2)
581 if (direction .eq. -1)
then 596 elseif (direction .eq. 1)
then 610 print *,
"direction is either -1 (forward) or 1 (backward)" 620 if (direction .eq. -1)
then 635 elseif (direction .eq. 1)
then 649 print *,
"direction is either -1 (forward) or 1 (backward)" 657 call transpose_y_to_z(
c_tmp2,c_inout3)
660 if (direction .eq. -1)
then 681 if (l .eq. -
nx/2)
then 682 c_inout3(i,j,k)=cmplx(0.d0,0.d0,kind=mytype)
685 if (m .eq. -
ny/2)
then 686 c_inout3(i,j,k)=cmplx(0.d0,0.d0,kind=mytype)
711 real(mytype),
dimension(:,:,:),
intent(inout) :: varin
713 real(mytype) :: zd,sigz,mean_local,mean_all
724 mean_local=mean_local+varin(i,j,k)
727 mean_local=mean_local/zsize(1)/zsize(2)
731 call mpi_allreduce(mean_local,mean_all,1,real_type,mpi_sum, &
733 mean_all=mean_all/
nprc 737 sigz=
cdamp*0.5d0*( 1.d0-cos(
pi*zd) )
740 varin(i,j,k)=varin(i,j,k)-sigz*( varin(i,j,k)-mean_all )
real(mytype), dimension(:,:,:), allocatable, save mean_uw
real(mytype), dimension(:), allocatable, save aa_w
integer, save, protected myid
real(mytype), dimension(:), allocatable, save aa_t
complex(mytype), dimension(:,:,:), allocatable, save c_qu
integer, save, protected nprc
integer, save, protected isnoise
integer, save, protected ny
real(mytype), dimension(:,:,:), allocatable, save vhy
complex(mytype), dimension(:,:,:), allocatable, save bb_u
real(mytype), dimension(:,:,:), allocatable, save thx
real(mytype), dimension(:,:,:), allocatable, save ss
complex(mytype), dimension(:,:,:), allocatable, save c_u
complex(mytype), dimension(:,:,:), allocatable, save c_qv
real(mytype), dimension(:,:,:), allocatable, save uhz
complex(mytype), dimension(:,:,:), allocatable, save c_src
real(mytype), parameter pi
real(mytype), dimension(:), allocatable, save aa_p
integer, save, protected dts
complex(mytype), dimension(:,:,:), allocatable, save c_ss
real(mytype), dimension(:,:,:), allocatable, save mean_w
complex(mytype), dimension(:,:,:), allocatable, save c_hh
real(mytype), dimension(:,:,:), allocatable, save ff0
integer, save, protected jend3
real(mytype), dimension(:,:,:), allocatable, save mean_u
real(mytype), dimension(:), allocatable, save cc_p
complex(mytype), dimension(:,:,:), allocatable, save c_p
real(mytype), dimension(:), allocatable, save cc_w
real(mytype), dimension(:,:,:), allocatable, save ss0
real(mytype), dimension(:,:,:), allocatable, save mean_v
real(mytype), dimension(:,:,:), allocatable, save cc
real(mytype), dimension(:,:,:), allocatable, save p
real(mytype), dimension(:,:,:), allocatable, save thz
real(mytype), dimension(:,:,:), allocatable, save mean_vw
integer, save, protected nia3
type(c_ptr), save, protected fft_plan1
real(mytype), dimension(:,:,:), allocatable, save aa
complex(mytype), dimension(:,:,:), allocatable, save c_tmp2
complex(mytype), dimension(:,:,:), allocatable, save c_phi
complex(mytype), dimension(:,:,:), allocatable, save bb_t
real(mytype), dimension(:), allocatable, save cc_u
real(mytype), dimension(:,:,:), allocatable, save mean_tt
real(mytype), save, protected cdamp
complex(mytype), dimension(:,:,:), allocatable, save bb_w
complex(mytype), dimension(:,:,:), allocatable, save c_gg
real(mytype), dimension(:,:,:), allocatable, save gg
real(mytype), dimension(:,:,:), allocatable, save vhx
integer, save, protected istr3
real(mytype), dimension(:,:,:), allocatable, save hh0
complex(mytype), dimension(:,:,:), allocatable, save c_tmp1
real(mytype), dimension(:,:,:), allocatable, save mean_ww
integer, save, protected jstr3
real(mytype), dimension(:,:,:), allocatable, save mean_t
real(mytype), dimension(:,:,:), allocatable, save w
real(mytype), dimension(:,:,:), allocatable, save u
type(c_ptr), save, protected ifft_plan1
integer, save, protected i_offset
integer, save, protected nja3
type(c_ptr), save, protected ifft_plan2
integer, save, protected j_offset
real(mytype), save, protected u_mrf
real(mytype), dimension(:,:,:), allocatable, save hh
complex(mytype), dimension(:,:,:), allocatable, save c_v
complex(mytype), dimension(:,:,:), allocatable, save c_ff
complex(mytype), dimension(:,:,:), allocatable, save c_qw
integer, save, protected cds
real(mytype), dimension(:,:,:), allocatable, save mean_vv
integer, save, protected restart
real(mytype), dimension(:), allocatable, save aa_u
complex(mytype), dimension(:,:,:), allocatable, save c_w
real(mytype), dimension(:,:,:), allocatable, save t
integer, save, protected nzdamp
real(mytype), save, protected noise_mag
real(mytype), dimension(:,:,:), allocatable, save qw
real(mytype), dimension(:,:,:), allocatable, save whx
real(mytype), dimension(:,:,:), allocatable, save uhx
complex(mytype), dimension(:,:,:), allocatable, save c_t
real(mytype), dimension(:,:,:), allocatable, save why
integer, save, protected nx
real(mytype), dimension(:,:,:), allocatable, save mean_tw
real(mytype), dimension(:), allocatable, save cc_t
real(mytype), dimension(:,:,:), allocatable, save ff
real(mytype), dimension(:,:,:), allocatable, save mean_uu
real(mytype), dimension(:,:,:), allocatable, save v
type(c_ptr), save, protected fft_plan2
real(mytype), dimension(:,:,:), allocatable, save gg0
real(mytype), dimension(:,:,:), allocatable, save qu
integer, save, protected iend3
real(mytype), dimension(:), allocatable, save, protected l_t
real(mytype), dimension(:,:,:), allocatable, save phi
complex(mytype), dimension(:,:,:), allocatable, save bb
real(mytype), dimension(:,:,:), allocatable, save thy
real(mytype), dimension(:,:,:), allocatable, save qv
complex(mytype), dimension(:,:,:), allocatable, save c_var
real(mytype), dimension(:,:,:), allocatable, save uhy
complex(mytype), dimension(:,:,:), allocatable, save bb_p
integer, save, protected kstr3
integer, save, protected nka3
real(mytype), dimension(:,:,:), allocatable, save vhz
integer, save, protected kend3