Introduction
The usvpst user subroutine allows one to output additional variables on a post-processing mesh. Several "automatic" post-processing meshes may be defined :
- The volume mesh (
ipart=-1
)
- The boundary mesh (
ipart=-2
)
- SYRTHES coupling surface (
ipart
< -2)
- Cooling tower exchange zone meshes (
ipart
< -2) if ichrze
= 1
Additional meshes (cells or faces) may also be defined through the GUI or using the cs_user_postprocess_meshes function from the cs_user_postprocess.c file.
The examples of post-processing given below are using the meshes defined here.
Output on the volume mesh (ipart = -1)
Output of a combination of moments
A combination of moments can also be post-processed using the usvpst subroutine.
imom1 = 1
imom2 = 2
allocate(scel(ncelps))
do iloc = 1, ncelps
iel = lstcel(iloc)
scel(iloc) = cmom_2(iel) - (cmom_1(iel))**2
enddo
idimt = 1
ientla = .true.
ivarpr = .false.
deallocate(scel)
endif
Output on the boundary mesh (ipart = -2)
Variables can be post-processed on the boundary mesh even if they are orignally define at cell centers. The following code block illustrates the post-processing of the density on the boundary mesh.
else if (ipart .eq. -2) then
idimt = 1
ientla = .true.
ivarpr = .true.
call
post_write_var(ipart,
'Density at boundary', idimt, ientla, ivarpr, &
Output on user meshes 1 or 2
User meshes appearing in the examples are defined here.
Output of an interpolated velocity on user meshes
An interpolated velocity is computed on both interior and boundary faces using a simple linear interpolation on user meshes 1 or 2.
else if (ipart.eq.1 .or. ipart.eq.2) then
endif
allocate(vfac(3,nfacps), vfbr(3,nfbrps))
do iloc = 1, nfacps
ifac = lstfac(iloc)
vfac(1,iloc) = pnd * cvar_vel(1,ii) + (1.d0 - pnd) * cvar_vel(1,jj)
vfac(2,iloc) = pnd * cvar_vel(2,ii) + (1.d0 - pnd) * cvar_vel(2,jj)
vfac(3,iloc) = pnd * cvar_vel(3,ii) + (1.d0 - pnd) * cvar_vel(3,jj)
enddo
do iloc = 1, nfbrps
ifac = lstfbr(iloc)
vfbr(1,iloc) = cvar_vel(1,ii)
vfbr(2,iloc) = cvar_vel(2,ii)
vfbr(3,iloc) = cvar_vel(3,ii)
enddo
idimt = 3
ientla = .true.
ivarpr = .false.
call
post_write_var(ipart,
'Interpolated velocity', idimt, ientla, ivarpr, &
deallocate(vfac, vfbr)
Output of the pressure on user meshes
Similarly, the pressure is computed on both interior and boundary faces and then post-processed on user meshes 1 or 2.
endif
allocate(sfac(nfacps), sfbr(nfbrps))
do iloc = 1, nfacps
ifac = lstfac(iloc)
sfac(iloc) = pnd * cvar_var(ii) &
+ (1.d0 - pnd) * cvar_var(jj)
enddo
do iloc = 1, nfbrps
ifac = lstfbr(iloc)
sfbr(iloc) = cvar_var(ii)
enddo
idimt = 1
ientla = .true.
ivarpr = .false.
call
post_write_var(ipart,
'Interpolated pressure', idimt, ientla, ivarpr, &
deallocate(sfac, sfbr)
Output of the centers of gravity in different ways
The examples below illustrate how to output a same variable in different ways (interlaced or not, using an indirection or not).
Output of the centers of gravity, interleaved
if (intpst.eq.1) then
allocate(vfac(3,nfacps), vfbr(3,nfbrps))
do iloc = 1, nfacps
ifac = lstfac(iloc)
vfac(1,iloc) =
cdgfac(1, ifac)
vfac(2,iloc) =
cdgfac(2, ifac)
vfac(3,iloc) =
cdgfac(3, ifac)
enddo
do iloc = 1, nfbrps
ifac = lstfbr(iloc)
vfbr(1, iloc) =
cdgfbo(1, ifac)
vfbr(2, iloc) =
cdgfbo(2, ifac)
vfbr(3, iloc) =
cdgfbo(3, ifac)
enddo
ntindp = -1
idimt = 3
ientla = .true.
ivarpr = .false.
ientla, ivarpr, &
ntindp,
ttcabs, rvoid, vfac, vfbr)
deallocate(vfac, vfbr)
endif
Output of the centers of gravity, non-interleaved, time-dependant
if (intpst.eq.1) then
allocate(vfac(nfacps, 3), vfbr(nfbrps, 3))
do iloc = 1, nfacps
ifac = lstfac(iloc)
vfac(iloc,1) =
cdgfac(1, ifac)
vfac(iloc,2) =
cdgfac(2, ifac)
vfac(iloc,3) =
cdgfac(3, ifac)
enddo
do iloc = 1, nfbrps
ifac = lstfbr(iloc)
vfbr(iloc,1) =
cdgfbo(1, ifac)
vfbr(iloc,2) =
cdgfbo(2, ifac)
vfbr(iloc,3) =
cdgfbo(3, ifac)
enddo
ntindp = -1
idimt = 3
ientla = .false.
ivarpr = .false.
ientla, ivarpr, &
ntindp,
ttcabs, rvoid, vfac, vfbr)
deallocate(vfac, vfbr)
endif
Output of the centers of gravity, with indirection (parent-based)
if (intpst.eq.1) then
ntindp = -1
idimt = 3
ientla = .true.
ivarpr = .true.
call
post_write_var(ipart,
'face cog (parent)', idimt, ientla, ivarpr, &
endif
Output on user meshes 3 or 4
Output of an interpolated velocity on user meshes
else if (ipart.ge.3 .and. ipart.le.4) then
endif
allocate(vfac(3,nfacps), vfbr(3,nfbrps))
do iloc = 1, nfacps
ifac = lstfac(iloc)
vfac(1,iloc) = pnd * cvar_vel(1,ii) &
+ (1.d0 - pnd) * cvar_vel(1,jj)
vfac(2,iloc) = pnd * cvar_vel(2,ii) &
+ (1.d0 - pnd) * cvar_vel(2,jj)
vfac(3,iloc) = pnd * cvar_vel(3,ii) &
+ (1.d0 - pnd) * cvar_vel(3,jj)
enddo
do iloc = 1, nfbrps
ifac = lstfbr(iloc)
vfbr(1,iloc) = cvar_vel(1,ii)
vfbr(2,iloc) = cvar_vel(2,ii)
vfbr(3,iloc) = cvar_vel(3,ii)
enddo
idimt = 3
ientla = .true.
ivarpr = .false.
deallocate(vfac, vfbr)
Output of the pressure on user meshes
endif
allocate(sfac(nfacps), sfbr(nfbrps))
do iloc = 1, nfacps
ifac = lstfac(iloc)
sfac(iloc) = pnd * cvar_var(ii) &
+ (1.d0 - pnd) * cvar_var(jj)
enddo
do iloc = 1, nfbrps
ifac = lstfbr(iloc)
sfbr(iloc) = cvar_var(ii)
enddo
idimt = 1
ientla = .true.
ivarpr = .false.
deallocate(sfac, sfbr)
endif