88 subroutine drag_torque_zone(dgtq, tstep, zone, center, s11, s22, s33, s12, s13, s23,&
90 integer,
intent(in) :: tstep
92 type(
coef_t),
intent(inout) :: coef
93 real(kind=
rp),
intent(inout) :: s11(coef%Xh%lx,coef%Xh%lx,coef%Xh%lz,coef%msh%nelv)
94 real(kind=
rp),
intent(inout) :: s22(coef%Xh%lx,coef%Xh%lx,coef%Xh%lz,coef%msh%nelv)
95 real(kind=
rp),
intent(inout) :: s33(coef%Xh%lx,coef%Xh%lx,coef%Xh%lz,coef%msh%nelv)
96 real(kind=
rp),
intent(inout) :: s12(coef%Xh%lx,coef%Xh%lx,coef%Xh%lz,coef%msh%nelv)
97 real(kind=
rp),
intent(inout) :: s13(coef%Xh%lx,coef%Xh%lx,coef%Xh%lz,coef%msh%nelv)
98 real(kind=
rp),
intent(inout) :: s23(coef%Xh%lx,coef%Xh%lx,coef%Xh%lz,coef%msh%nelv)
99 type(
field_t),
intent(inout) :: p
100 real(kind=
rp),
intent(in) :: visc, center(3)
101 real(kind=
rp) :: dgtq(3,4)
102 real(kind=
rp) :: dragpx = 0.0_rp
103 real(kind=
rp) :: dragpy = 0.0_rp
104 real(kind=
rp) :: dragpz = 0.0_rp
105 real(kind=
rp) :: dragvx = 0.0_rp
106 real(kind=
rp) :: dragvy = 0.0_rp
107 real(kind=
rp) :: dragvz = 0.0_rp
108 real(kind=
rp) :: torqpx = 0.0_rp
109 real(kind=
rp) :: torqpy = 0.0_rp
110 real(kind=
rp) :: torqpz = 0.0_rp
111 real(kind=
rp) :: torqvx = 0.0_rp
112 real(kind=
rp) :: torqvy = 0.0_rp
113 real(kind=
rp) :: torqvz = 0.0_rp
114 real(kind=
rp) :: dragx, dragy, dragz
115 real(kind=
rp) :: torqx, torqy, torqz
116 integer :: ie, ifc, mem, ierr
131 ie = zone%facet_el(mem)%x(2)
132 ifc = zone%facet_el(mem)%x(1)
135 s11, s22, s33, s12, s13, s23,&
136 p%x,visc,ifc,ie, coef, coef%Xh)
138 dragpx = dragpx + dgtq(1,1)
139 dragpy = dragpy + dgtq(2,1)
140 dragpz = dragpz + dgtq(3,1)
142 dragvx = dragvx + dgtq(1,2)
143 dragvy = dragvy + dgtq(2,2)
144 dragvz = dragvz + dgtq(3,2)
146 torqpx = torqpx + dgtq(1,3)
147 torqpy = torqpy + dgtq(2,3)
148 torqpz = torqpz + dgtq(3,3)
150 torqvx = torqvx + dgtq(1,4)
151 torqvy = torqvy + dgtq(2,4)
152 torqvz = torqvz + dgtq(3,4)
157 call mpi_allreduce(mpi_in_place,dragpx, 1, &
159 call mpi_allreduce(mpi_in_place,dragpy, 1, &
161 call mpi_allreduce(mpi_in_place,dragpz, 1, &
163 call mpi_allreduce(mpi_in_place,dragvx, 1, &
165 call mpi_allreduce(mpi_in_place,dragvy, 1, &
167 call mpi_allreduce(mpi_in_place,dragvz, 1, &
170 call mpi_allreduce(mpi_in_place,torqpx, 1, &
172 call mpi_allreduce(mpi_in_place,torqpy, 1, &
174 call mpi_allreduce(mpi_in_place,torqpz, 1, &
176 call mpi_allreduce(mpi_in_place,torqvx, 1, &
178 call mpi_allreduce(mpi_in_place,torqvy, 1, &
180 call mpi_allreduce(mpi_in_place,torqvz, 1, &
213 s11, s22, s33, s12, s13, s23,&
214 pm1,visc,f,e, coef, Xh)
215 type(
coef_t),
intent(in) :: coef
216 type(
space_t),
intent(in) :: xh
217 real(kind=
rp),
intent(out) :: dgtq(3,4)
218 real(kind=
rp),
intent(in) :: center(3)
219 real(kind=
rp),
intent(in) :: xm0(xh%lx,xh%ly,xh%lz,coef%msh%nelv)
220 real(kind=
rp),
intent(in) :: ym0(xh%lx,xh%ly,xh%lz,coef%msh%nelv)
221 real(kind=
rp),
intent(in) :: zm0(xh%lx,xh%ly,xh%lz,coef%msh%nelv)
222 real(kind=
rp),
intent(in) :: s11(xh%lx,xh%ly,xh%lz,coef%msh%nelv)
223 real(kind=
rp),
intent(in) :: s22(xh%lx,xh%ly,xh%lz,coef%msh%nelv)
224 real(kind=
rp),
intent(in) :: s33(xh%lx,xh%ly,xh%lz,coef%msh%nelv)
225 real(kind=
rp),
intent(in) :: s12(xh%lx,xh%ly,xh%lz,coef%msh%nelv)
226 real(kind=
rp),
intent(in) :: s13(xh%lx,xh%ly,xh%lz,coef%msh%nelv)
227 real(kind=
rp),
intent(in) :: s23(xh%lx,xh%ly,xh%lz,coef%msh%nelv)
228 real(kind=
rp),
intent(in) :: pm1(xh%lx,xh%ly,xh%lz,coef%msh%nelv)
229 real(kind=
rp),
intent(in) :: visc
230 integer,
intent(in) :: f,e
231 integer :: pf,l, k, i, j1, j2
232 real(kind=
rp) :: n1,n2,n3, j, a, r1, r2, r3, v, dgtq_i(3,4)
233 integer :: skpdat(6,6), nx, ny, nz
240 real(kind=
rp) :: s11_, s21_, s31_, s12_, s22_, s32_, s13_, s23_, s33_
247 skpdat(2,1)=nx*(ny-1)+1
250 skpdat(5,1)=ny*(nz-1)+1
253 skpdat(1,2)=1 + (nx-1)
254 skpdat(2,2)=nx*(ny-1)+1 + (nx-1)
257 skpdat(5,2)=ny*(nz-1)+1
264 skpdat(5,3)=ny*(nz-1)+1
267 skpdat(1,4)=1 + nx*(ny-1)
268 skpdat(2,4)=nx + nx*(ny-1)
271 skpdat(5,4)=ny*(nz-1)+1
281 skpdat(1,6)=1 + nx*ny*(nz-1)
282 skpdat(2,6)=nx + nx*ny*(nz-1)
290 jskip1 = skpdat(3,pf)
293 jskip2 = skpdat(6,pf)
300 n1 = coef%nx(i,1,f,e)*coef%area(i,1,f,e)
301 n2 = coef%ny(i,1,f,e)*coef%area(i,1,f,e)
302 n3 = coef%nz(i,1,f,e)*coef%area(i,1,f,e)
303 a = a + coef%area(i,1,f,e)
305 s11_ = s11(j1,j2,1,e)
306 s12_ = s12(j1,j2,1,e)
307 s22_ = s22(j1,j2,1,e)
308 s13_ = s13(j1,j2,1,e)
309 s23_ = s23(j1,j2,1,e)
310 s33_ = s33(j1,j2,1,e)
311 call drag_torque_pt(dgtq_i,xm0(j1,j2,1,e), ym0(j1,j2,1,e),zm0(j1,j2,1,e), center,&
312 s11_, s22_, s33_, s12_, s13_, s23_,&
313 pm1(j1,j2,1,e), n1, n2, n3, v)
331 subroutine drag_torque_pt(dgtq,x,y,z, center, s11, s22, s33, s12, s13, s23,&
333 real(kind=
rp),
intent(inout) :: dgtq(3,4)
334 real(kind=
rp),
intent(in) :: x
335 real(kind=
rp),
intent(in) :: y
336 real(kind=
rp),
intent(in) :: z
337 real(kind=
rp),
intent(in) :: p
338 real(kind=
rp),
intent(in) :: v
339 real(kind=
rp),
intent(in) :: n1, n2, n3, center(3)
340 real(kind=
rp),
intent(in) :: s11, s12, s22, s13, s23, s33
341 real(kind=
rp) :: s21, s31, s32, r1, r2, r3
351 dgtq(1,2) = -v*(s11*n1 + s12*n2 + s13*n3)
352 dgtq(2,2) = -v*(s21*n1 + s22*n2 + s23*n3)
353 dgtq(3,2) = -v*(s31*n1 + s32*n2 + s33*n3)
358 dgtq(1,3) = (r2*dgtq(3,1)-r3*dgtq(2,1))
359 dgtq(2,3) = (r3*dgtq(1,1)-r1*dgtq(3,1))
360 dgtq(3,3) = (r1*dgtq(2,1)-r2*dgtq(1,1))
362 dgtq(1,4) = (r2*dgtq(3,2)-r3*dgtq(2,2))
363 dgtq(2,4) = (r3*dgtq(1,2)-r1*dgtq(3,2))
364 dgtq(3,4) = (r1*dgtq(2,2)-r2*dgtq(1,2))
type(mpi_comm) neko_comm
MPI communicator.
type(mpi_datatype) mpi_real_precision
MPI type for working precision of REAL types.
subroutine, public drag_torque_zone(dgtq, tstep, zone, center, s11, s22, s33, s12, s13, s23, p, coef, visc)
Some functions to calculate the lift/drag and torque Calculation can be done on a zone,...
subroutine, public drag_torque_facet(dgtq, xm0, ym0, zm0, center, s11, s22, s33, s12, s13, s23, pm1, visc, f, e, coef, Xh)
Calculate drag and torque over a facet.
subroutine, public drag_torque_pt(dgtq, x, y, z, center, s11, s22, s33, s12, s13, s23, p, n1, n2, n3, v)
Calculate drag and torque from one point.
Defines a zone as a subset of facets in a mesh.
subroutine, public rzero(a, n)
Zero a real vector.
integer, parameter, public rp
Global precision used in computations.
Defines a function space.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
The function space for the SEM solution fields.