seperate fsd and masked fsd

This commit is contained in:
ignis 2019-06-04 01:58:33 +09:00
parent b72c633b5f
commit 71504d516c
2 changed files with 129 additions and 108 deletions

View file

@ -4,13 +4,14 @@ c = 1.0 - y
wrate = rxn_rate(c)
mask = $threshold_min_max(c)
fsd = sqrt (sqr(ddx(c)) + sqr(ddy(c)) + sqr(ddz(c))) * mask
fsd_orig = sqrt (sqr(ddx(c)) + sqr(ddy(c)) + sqr(ddz(c))) * mask
fsd = fsd_orig * mask
sd = ((d2dx(c) + d2dy(c) + d2dz(c)) * $rod + wrate) / fsd
sd = ((d2dx(c) + d2dy(c) + d2dz(c)) * $rod + wrate) / fsd_orig
nx = - ddx(c) / fsd
ny = - ddy(c) / fsd
nz = - ddz(c) / fsd
nx = - ddx(c) / fsd_orig
ny = - ddy(c) / fsd_orig
nz = - ddz(c) / fsd_orig
vn = u * nx + v * ny + w * nz

View file

@ -261,25 +261,7 @@ end do
end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_u(i) = fsd_avg_u(i) + u(i,j,k) * xyzbuffer7(i,j,k)
end do
end do
end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_w(i) = fsd_avg_w(i) + w(i,j,k) * xyzbuffer7(i,j,k)
end do
end do
end do
! ( ( - ddx_c(i,j,k) ) / fsd(i,j,k) )
! ( ( - ddx_c(i,j,k) ) / fsd_orig(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
@ -289,7 +271,7 @@ end do
end do
! ( ( - ddz_c(i,j,k) ) / fsd(i,j,k) )
! ( ( - ddz_c(i,j,k) ) / fsd_orig(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
@ -319,10 +301,11 @@ end do
call d2dz ( xyzbuffer4, xyzbuffer11 )
call ddz ( xyzbuffer11, xyzbuffer5 )
! ( fsd_orig(i,j,k) * mask(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_nx(i) = fsd_avg_nx(i) + xyzbuffer6(i,j,k) * xyzbuffer7(i,j,k)
xyzbuffer9(i,j,k) = ( xyzbuffer7(i,j,k) * xyzbuffer2(i,j,k) )
end do
end do
end do
@ -331,7 +314,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_nz(i) = fsd_avg_nz(i) + xyzbuffer5(i,j,k) * xyzbuffer7(i,j,k)
fsd_avg_nx(i) = fsd_avg_nx(i) + xyzbuffer6(i,j,k) * xyzbuffer9(i,j,k)
end do
end do
end do
@ -340,17 +323,35 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_v(i) = fsd_avg_v(i) + v(i,j,k) * xyzbuffer7(i,j,k)
fsd_avg_nz(i) = fsd_avg_nz(i) + xyzbuffer5(i,j,k) * xyzbuffer9(i,j,k)
end do
end do
end do
! ( ( - ddy_c(i,j,k) ) / fsd(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer9(i,j,k) = ( ( - xyzbuffer8(i,j,k) ) / xyzbuffer7(i,j,k) )
fsd_avg_u(i) = fsd_avg_u(i) + u(i,j,k) * xyzbuffer9(i,j,k)
end do
end do
end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_w(i) = fsd_avg_w(i) + w(i,j,k) * xyzbuffer9(i,j,k)
end do
end do
end do
! ( ( - ddy_c(i,j,k) ) / fsd_orig(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer3(i,j,k) = ( ( - xyzbuffer8(i,j,k) ) / xyzbuffer7(i,j,k) )
end do
end do
end do
@ -360,7 +361,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer8(i,j,k) = ( ( ( u(i,j,k) * xyzbuffer6(i,j,k) ) + ( v(i,j,k) * xyzbuffer9(i,j,k) ) ) + ( w(i,j,k) * xyzbuffer5(i,j,k) ) )
xyzbuffer2(i,j,k) = ( ( ( u(i,j,k) * xyzbuffer6(i,j,k) ) + ( v(i,j,k) * xyzbuffer3(i,j,k) ) ) + ( w(i,j,k) * xyzbuffer5(i,j,k) ) )
end do
end do
end do
@ -374,12 +375,50 @@ end do
end do
end do
call ddx ( xyzbuffer3, xyzbuffer6 )
call ddx ( xyzbuffer5, xyzbuffer6 )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_ny(i) = fsd_avg_ny(i) + xyzbuffer9(i,j,k) * xyzbuffer7(i,j,k)
fsd_avg_ny(i) = fsd_avg_ny(i) + xyzbuffer3(i,j,k) * xyzbuffer9(i,j,k)
end do
end do
end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_v(i) = fsd_avg_v(i) + v(i,j,k) * xyzbuffer9(i,j,k)
end do
end do
end do
! ( ( ( ( ( d2dx_c(i,j,k) + d2dy_c(i,j,k) ) + d2dz_c(i,j,k) ) * rod ) + wrate(i,j,k) ) / fsd_orig(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer8(i,j,k) = ( ( ( ( ( xyzbuffer0(i,j,k) + xyzbuffer10(i,j,k) ) + xyzbuffer4(i,j,k) ) * rod ) + xyzbuffer1(i,j,k) ) / xyzbuffer7(i,j,k) )
end do
end do
end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
avg_fsd(i) = avg_fsd(i) + xyzbuffer9(i,j,k)
end do
end do
end do
call ddy ( xyzbuffer1, xyzbuffer3 )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_sd(i) = fsd_avg_sd(i) + xyzbuffer8(i,j,k) * xyzbuffer9(i,j,k)
end do
end do
end do
@ -389,7 +428,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer2(i,j,k) = ( ( xyzbuffer0(i,j,k) + xyzbuffer10(i,j,k) ) + xyzbuffer4(i,j,k) )
xyzbuffer7(i,j,k) = ( ( xyzbuffer0(i,j,k) + xyzbuffer10(i,j,k) ) + xyzbuffer4(i,j,k) )
end do
end do
end do
@ -398,36 +437,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
avg_fsd(i) = avg_fsd(i) + xyzbuffer7(i,j,k)
end do
end do
end do
call ddy ( xyzbuffer5, xyzbuffer9 )
! ( ( ( ( ( d2dx_c(i,j,k) + d2dy_c(i,j,k) ) + d2dz_c(i,j,k) ) * rod ) + wrate(i,j,k) ) / fsd(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer9(i,j,k) = ( ( ( ( ( xyzbuffer0(i,j,k) + xyzbuffer10(i,j,k) ) + xyzbuffer4(i,j,k) ) * rod ) + xyzbuffer1(i,j,k) ) / xyzbuffer7(i,j,k) )
end do
end do
end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
avg_lapc(i) = avg_lapc(i) + xyzbuffer2(i,j,k)
end do
end do
end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_sd(i) = fsd_avg_sd(i) + xyzbuffer9(i,j,k) * xyzbuffer7(i,j,k)
avg_lapc(i) = avg_lapc(i) + xyzbuffer7(i,j,k)
end do
end do
end do
@ -437,7 +447,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer2(i,j,k) = ( ( xyzbuffer8(i,j,k) + xyzbuffer9(i,j,k) ) * xyzbuffer6(i,j,k) )
xyzbuffer4(i,j,k) = ( ( xyzbuffer2(i,j,k) + xyzbuffer8(i,j,k) ) * xyzbuffer6(i,j,k) )
end do
end do
end do
@ -447,7 +457,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer1(i,j,k) = ( ( xyzbuffer3(i,j,k) + xyzbuffer5(i,j,k) ) + xyzbuffer11(i,j,k) )
xyzbuffer10(i,j,k) = ( ( xyzbuffer5(i,j,k) + xyzbuffer1(i,j,k) ) + xyzbuffer11(i,j,k) )
end do
end do
end do
@ -456,7 +466,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_divn(i) = fsd_avg_divn(i) + xyzbuffer1(i,j,k) * xyzbuffer7(i,j,k)
fsd_avg_divn(i) = fsd_avg_divn(i) + xyzbuffer10(i,j,k) * xyzbuffer9(i,j,k)
end do
end do
end do
@ -466,7 +476,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer0(i,j,k) = ( dabs ( xyzbuffer1(i,j,k) ) )
xyzbuffer8(i,j,k) = ( dabs ( xyzbuffer10(i,j,k) ) )
end do
end do
end do
@ -475,7 +485,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_t1(i) = fsd_avg_t1(i) + xyzbuffer2(i,j,k) * xyzbuffer7(i,j,k)
fsd_avg_t1(i) = fsd_avg_t1(i) + xyzbuffer4(i,j,k) * xyzbuffer9(i,j,k)
end do
end do
end do
@ -484,7 +494,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_absk(i) = fsd_avg_absk(i) + xyzbuffer0(i,j,k) * xyzbuffer7(i,j,k)
fsd_avg_absk(i) = fsd_avg_absk(i) + xyzbuffer8(i,j,k) * xyzbuffer9(i,j,k)
end do
end do
end do
@ -664,21 +674,43 @@ end do
end do
! ( ( - ddx_c(i,j,k) ) / fsd(i,j,k) )
! ( ( - ddy_c(i,j,k) ) / fsd_orig(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer5(i,j,k) = ( ( - xyzbuffer8(i,j,k) ) / xyzbuffer6(i,j,k) )
xyzbuffer5(i,j,k) = ( ( - xyzbuffer3(i,j,k) ) / xyzbuffer6(i,j,k) )
end do
end do
end do
call d2dz ( xyzbuffer4, xyzbuffer10 )
call ddy ( xyzbuffer10, xyzbuffer5 )
! ( fsd_orig(i,j,k) * mask(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer3(i,j,k) = ( xyzbuffer6(i,j,k) * xyzbuffer1(i,j,k) )
end do
end do
end do
! ( ( - ddz_c(i,j,k) ) / fsd(i,j,k) )
! ( ( - ddx_c(i,j,k) ) / fsd_orig(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer4(i,j,k) = ( ( - xyzbuffer2(i,j,k) ) / xyzbuffer6(i,j,k) )
xyzbuffer1(i,j,k) = ( ( - xyzbuffer8(i,j,k) ) / xyzbuffer6(i,j,k) )
end do
end do
end do
! ( ( ( ( ( d2dx_c(i,j,k) + d2dy_c(i,j,k) ) + d2dz_c(i,j,k) ) * rod ) + wrate(i,j,k) ) / fsd_orig(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer8(i,j,k) = ( ( ( ( ( xyzbuffer7(i,j,k) + xyzbuffer9(i,j,k) ) + xyzbuffer4(i,j,k) ) * rod ) + xyzbuffer0(i,j,k) ) / xyzbuffer6(i,j,k) )
end do
end do
end do
@ -688,28 +720,27 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer8(i,j,k) = (((xyzbuffer5(i,j,k) - fsd_avg_nx(i)))*((xyzbuffer5(i,j,k) - fsd_avg_nx(i))))
xyzbuffer0(i,j,k) = (((xyzbuffer1(i,j,k) - fsd_avg_nx(i)))*((xyzbuffer1(i,j,k) - fsd_avg_nx(i))))
end do
end do
end do
call d2dz ( xyzbuffer2, xyzbuffer10 )
call ddz ( xyzbuffer1, xyzbuffer4 )
call ddx ( xyzbuffer7, xyzbuffer1 )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_t6(i) = fsd_avg_t6(i) + xyzbuffer8(i,j,k) * xyzbuffer6(i,j,k)
fsd_avg_t6(i) = fsd_avg_t6(i) + xyzbuffer0(i,j,k) * xyzbuffer3(i,j,k)
end do
end do
end do
! ( ( - ddy_c(i,j,k) ) / fsd(i,j,k) )
! ( ( - ddz_c(i,j,k) ) / fsd_orig(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer10(i,j,k) = ( ( - xyzbuffer3(i,j,k) ) / xyzbuffer6(i,j,k) )
xyzbuffer4(i,j,k) = ( ( - xyzbuffer2(i,j,k) ) / xyzbuffer6(i,j,k) )
end do
end do
end do
@ -719,29 +750,18 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer8(i,j,k) = ( ( ( u(i,j,k) * xyzbuffer5(i,j,k) ) + ( v(i,j,k) * xyzbuffer10(i,j,k) ) ) + ( w(i,j,k) * xyzbuffer4(i,j,k) ) )
xyzbuffer9(i,j,k) = ( ( ( u(i,j,k) * xyzbuffer1(i,j,k) ) + ( v(i,j,k) * xyzbuffer5(i,j,k) ) ) + ( w(i,j,k) * xyzbuffer4(i,j,k) ) )
end do
end do
end do
call ddx ( xyzbuffer3, xyzbuffer5 )
! ( ( ( ( ( d2dx_c(i,j,k) + d2dy_c(i,j,k) ) + d2dz_c(i,j,k) ) * rod ) + wrate(i,j,k) ) / fsd(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer4(i,j,k) = ( ( ( ( ( xyzbuffer7(i,j,k) + xyzbuffer9(i,j,k) ) + xyzbuffer2(i,j,k) ) * rod ) + xyzbuffer0(i,j,k) ) / xyzbuffer6(i,j,k) )
end do
end do
end do
call ddy ( xyzbuffer9, xyzbuffer10 )
call ddz ( xyzbuffer2, xyzbuffer4 )
! ( ( vn(i,j,k) + sd(i,j,k) ) * nx(i,j,k) )
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer2(i,j,k) = ( ( xyzbuffer8(i,j,k) + xyzbuffer4(i,j,k) ) * xyzbuffer5(i,j,k) )
xyzbuffer0(i,j,k) = ( ( xyzbuffer9(i,j,k) + xyzbuffer8(i,j,k) ) * xyzbuffer1(i,j,k) )
end do
end do
end do
@ -751,7 +771,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer0(i,j,k) = (((xyzbuffer2(i,j,k) - fsd_avg_t1(i)))*((xyzbuffer2(i,j,k) - fsd_avg_t1(i))))
xyzbuffer6(i,j,k) = (((xyzbuffer0(i,j,k) - fsd_avg_t1(i)))*((xyzbuffer0(i,j,k) - fsd_avg_t1(i))))
end do
end do
end do
@ -761,7 +781,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer7(i,j,k) = ( ( xyzbuffer3(i,j,k) + xyzbuffer9(i,j,k) ) + xyzbuffer1(i,j,k) )
xyzbuffer5(i,j,k) = ( ( xyzbuffer7(i,j,k) + xyzbuffer10(i,j,k) ) + xyzbuffer2(i,j,k) )
end do
end do
end do
@ -770,7 +790,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_t5(i) = fsd_avg_t5(i) + xyzbuffer0(i,j,k) * xyzbuffer6(i,j,k)
fsd_avg_t5(i) = fsd_avg_t5(i) + xyzbuffer6(i,j,k) * xyzbuffer3(i,j,k)
end do
end do
end do
@ -780,7 +800,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer4(i,j,k) = ( (xyzbuffer2(i,j,k) - fsd_avg_t1(i)) * (xyzbuffer5(i,j,k) - fsd_avg_nx(i)) )
xyzbuffer4(i,j,k) = ( (xyzbuffer0(i,j,k) - fsd_avg_t1(i)) * (xyzbuffer1(i,j,k) - fsd_avg_nx(i)) )
end do
end do
end do
@ -790,7 +810,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer10(i,j,k) = (((xyzbuffer7(i,j,k) - fsd_avg_divn(i)))*((xyzbuffer7(i,j,k) - fsd_avg_divn(i))))
xyzbuffer10(i,j,k) = (((xyzbuffer5(i,j,k) - fsd_avg_divn(i)))*((xyzbuffer5(i,j,k) - fsd_avg_divn(i))))
end do
end do
end do
@ -800,7 +820,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer9(i,j,k) = ( dabs ( xyzbuffer7(i,j,k) ) )
xyzbuffer9(i,j,k) = ( dabs ( xyzbuffer5(i,j,k) ) )
end do
end do
end do
@ -809,7 +829,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_t7(i) = fsd_avg_t7(i) + xyzbuffer10(i,j,k) * xyzbuffer6(i,j,k)
fsd_avg_t7(i) = fsd_avg_t7(i) + xyzbuffer10(i,j,k) * xyzbuffer3(i,j,k)
end do
end do
end do
@ -819,7 +839,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer8(i,j,k) = ( (xyzbuffer2(i,j,k) - fsd_avg_t1(i)) * (xyzbuffer7(i,j,k) - fsd_avg_divn(i)) )
xyzbuffer8(i,j,k) = ( (xyzbuffer0(i,j,k) - fsd_avg_t1(i)) * (xyzbuffer5(i,j,k) - fsd_avg_divn(i)) )
end do
end do
end do
@ -829,7 +849,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer3(i,j,k) = (((xyzbuffer9(i,j,k) - fsd_avg_absk(i)))*((xyzbuffer9(i,j,k) - fsd_avg_absk(i))))
xyzbuffer2(i,j,k) = (((xyzbuffer9(i,j,k) - fsd_avg_absk(i)))*((xyzbuffer9(i,j,k) - fsd_avg_absk(i))))
end do
end do
end do
@ -838,7 +858,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_t2(i) = fsd_avg_t2(i) + xyzbuffer4(i,j,k) * xyzbuffer6(i,j,k)
fsd_avg_t2(i) = fsd_avg_t2(i) + xyzbuffer4(i,j,k) * xyzbuffer3(i,j,k)
end do
end do
end do
@ -847,7 +867,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_t8(i) = fsd_avg_t8(i) + xyzbuffer3(i,j,k) * xyzbuffer6(i,j,k)
fsd_avg_t8(i) = fsd_avg_t8(i) + xyzbuffer2(i,j,k) * xyzbuffer3(i,j,k)
end do
end do
end do
@ -856,7 +876,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_t3(i) = fsd_avg_t3(i) + xyzbuffer8(i,j,k) * xyzbuffer6(i,j,k)
fsd_avg_t3(i) = fsd_avg_t3(i) + xyzbuffer8(i,j,k) * xyzbuffer3(i,j,k)
end do
end do
end do
@ -866,7 +886,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
xyzbuffer1(i,j,k) = ( (xyzbuffer2(i,j,k) - fsd_avg_t1(i)) * (xyzbuffer9(i,j,k) - fsd_avg_absk(i)) )
xyzbuffer1(i,j,k) = ( (xyzbuffer0(i,j,k) - fsd_avg_t1(i)) * (xyzbuffer9(i,j,k) - fsd_avg_absk(i)) )
end do
end do
end do
@ -875,7 +895,7 @@ end do
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
fsd_avg_t4(i) = fsd_avg_t4(i) + xyzbuffer1(i,j,k) * xyzbuffer6(i,j,k)
fsd_avg_t4(i) = fsd_avg_t4(i) + xyzbuffer1(i,j,k) * xyzbuffer3(i,j,k)
end do
end do
end do