programmer's documentation
Output additional variables on a post-processing mesh (cs_user_postprocess_var.f90)

Introduction

The usvpst user subroutine allows one to output additional variables on a post-processing mesh. Several "automatic" post-processing meshes may be defined :

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 the turbulent kinetic energy for the Rij-Epsilon model on the volume mesh

One can define, compute and post-process the turbulent kinetic energy for the Rij-Epsilon as shown in the following example.

!===============================================================================
! Examples of volume variables on the main volume mesh (ipart = -1)
!===============================================================================
if (ipart .eq. -1) then
! Output of k=1/2(R11+R22+R33) for the Rij-epsilon model
! ------------------------------------------------------
if (itytur .eq. 3) then
allocate(scel(ncelps))
call field_get_val_s(ivarfl(ir11), cvar_r11)
call field_get_val_s(ivarfl(ir22), cvar_r22)
call field_get_val_s(ivarfl(ir33), cvar_r33)
do iloc = 1, ncelps
iel = lstcel(iloc)
scel(iloc) = 0.5d0*( cvar_r11(iel) &
+ cvar_r22(iel) &
+ cvar_r33(iel))
enddo
idimt = 1 ! 1: scalar, 3: vector, 6/9: symm/non-symm tensor
ientla = .true. ! dimension 1 here, so no effect
ivarpr = .false. ! defined on the work array, not on the parent
! Output values; as we have no face values, we can pass a
! trivial array rvoid for those.
call post_write_var(ipart, 'Turb energy', idimt, ientla, ivarpr, &
ntcabs, ttcabs, scel, rvoid, rvoid)
deallocate(scel)
endif

Output of a combination of moments

A combination of moments can also be post-processed using the usvpst subroutine.

! Output of a combination of moments
! ----------------------------------
! We assume in this example that we have 2 temporal means (moments):
! <u> for imom=1
! <uu> for imom=2
! We seek to plot <u'u'>=<uu>-<U>**2
if (cs_time_moment_n_moments() .ge. 2) then
! Moment numbers:
imom1 = 1
imom2 = 2
! Temporal accumulation for moments
call field_get_val_s(time_moment_field_id(imom1), cmom_1)
call field_get_val_s(time_moment_field_id(imom2), cmom_2)
! To improve this example's readability, we assume moments imom1 and imom2
! have been computed on the same time window.
allocate(scel(ncelps))
do iloc = 1, ncelps
iel = lstcel(iloc)
scel(iloc) = cmom_2(iel) - (cmom_1(iel))**2
enddo
idimt = 1 ! 1: scalar, 3: vector, 6/9: symm/non-symm tensor
ientla = .true. ! dimension 1 here, so no effect
ivarpr = .false. ! defined on the work array, not on the parent
! Output values; as we have no face values, we can pass a
! trivial array for those.
call post_write_var(ipart, '<upup>', idimt, ientla, ivarpr, &
ntcabs, ttcabs, scel, rvoid, rvoid)
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.

!===============================================================================
! Examples of volume variables on the boundary mesh (ipart = -2)
!===============================================================================
else if (ipart .eq. -2) then
! Output of the density at the boundary
! -------------------------------------
idimt = 1 ! 1: scalar, 3: vector, 6/9: symm/non-symm tensor
ientla = .true. ! dimension 1 here, so no effect
ivarpr = .true. ! we use the bfpro_rom array defined on the parent mesh
! Output values; as we have no cell or interior face values, we can pass a
! trivial array for those.
call field_get_val_s(ibrom, bfpro_rom)
call post_write_var(ipart, 'Density at boundary', idimt, ientla, ivarpr, &
ntcabs, ttcabs, rvoid, rvoid, bfpro_rom)

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.

!===============================================================================
! Examples of volume variables on user meshes 1 or 2
!===============================================================================
else if (ipart.eq.1 .or. ipart.eq.2) then
! Map field arrays
call field_get_val_v(ivarfl(iu), cvar_vel)
! Output of the velocity
! ----------------------
! Compute variable values on interior faces.
! In this example, we use a simple linear interpolation.
! For parallel calculations, if neighbors are used, they must be synchronized
! first. This also applies for periodicity.
if (irangp.ge.0.or.iperio.eq.1) then
call synvin(cvar_vel)
endif
allocate(vfac(3,nfacps), vfbr(3,nfbrps))
do iloc = 1, nfacps
ifac = lstfac(iloc)
ii = ifacel(1, ifac)
jj = ifacel(2, ifac)
pnd = pond(ifac)
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
! Compute variable values on boundary faces.
! In this example, we use a simple copy of the adjacent cell value.
do iloc = 1, nfbrps
ifac = lstfbr(iloc)
ii = ifabor(ifac)
vfbr(1,iloc) = cvar_vel(1,ii)
vfbr(2,iloc) = cvar_vel(2,ii)
vfbr(3,iloc) = cvar_vel(3,ii)
enddo
idimt = 3 ! 1: scalar, 3: vector, 6/9: symm/non-symm tensor
ientla = .true. ! interleaved
ivarpr = .false. ! defined on the work array, not on the parent
! Output values; as we have no cell values, we can pass a
! trivial array for those.
call post_write_var(ipart, 'Interpolated velocity', idimt, ientla, ivarpr, &
ntcabs, ttcabs, rvoid, vfac, vfbr)
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.

! Output of the pressure
! ----------------------
! Variable number
ivar = ipr
call field_get_val_s(ivarfl(ivar), cvar_var)
! Compute variable values on interior faces.
! In this example, we use a simple linear interpolation.
! For parallel calculations, if neighbors are used, they must be synchronized
! first. This also applies for periodicity.
if (irangp.ge.0.or.iperio.eq.1) then
call synsca(cvar_var)
endif
allocate(sfac(nfacps), sfbr(nfbrps))
do iloc = 1, nfacps
ifac = lstfac(iloc)
ii = ifacel(1, ifac)
jj = ifacel(2, ifac)
pnd = pond(ifac)
sfac(iloc) = pnd * cvar_var(ii) &
+ (1.d0 - pnd) * cvar_var(jj)
enddo
! Compute variable values on boundary faces.
! In this example, we use a simple copy of the adjacent cell value.
do iloc = 1, nfbrps
ifac = lstfbr(iloc)
ii = ifabor(ifac)
sfbr(iloc) = cvar_var(ii)
enddo
idimt = 1 ! 1: scalar, 3: vector, 6/9: symm/non-symm tensor
ientla = .true. ! dimension 1 here, so no effect
ivarpr = .false. ! defined on the work array, not on the parent
! Output values; as we have no cell values, we can pass a
! trivial array for those.
call post_write_var(ipart, 'Interpolated pressure', idimt, ientla, ivarpr, &
ntcabs, ttcabs, rvoid, sfac, sfbr)
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

! Output of the centers of gravity, interlaced
! --------------------------------
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
! Compute variable values on boundary faces
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
! We assign a negative time step and output this variable once only
! to avoid duplicating it at each output time (assuming a fixed mesh).
ntindp = -1
idimt = 3 ! 1: scalar, 3: vector, 6/9: symm/non-symm tensor
ientla = .true. ! interleaved
ivarpr = .false. ! defined on the work array, not on the parent
! Output values; as we have no cell values, we can pass a
! trivial array for those.
call post_write_var(ipart, 'face cog (interlaced)', idimt, &
ientla, ivarpr, &
ntindp, ttcabs, rvoid, vfac, vfbr)
deallocate(vfac, vfbr)
endif

Output of the centers of gravity, non-interleaved, time-dependant

! Output of the centers of gravity, non-interlaced, time independent
! --------------------------------
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
! Compute variable values on boundary faces
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
! We assign a negative time step and output this variable once only
! to avoid duplicating it at each output time (assuming a fixed mesh).
ntindp = -1
idimt = 3 ! 1: scalar, 3: vector, 6/9: symm/non-symm tensor
ientla = .false. ! not interleaved
ivarpr = .false. ! defined on the work array, not on the parent
! Output values; as we have no cell values, we can pass a
! trivial array for those.
call post_write_var(ipart, 'face cog (non interlaced)', idimt, &
ientla, ivarpr, &
ntindp, ttcabs, rvoid, vfac, vfbr)
deallocate(vfac, vfbr)
endif

Output of the centers of gravity, with indirection (parent-based)

! Output of the centers of gravity, with indirection (parent-based)
! --------------------------------
if (intpst.eq.1) then
! We assign a negative time step and output this variable once only
! to avoid duplicating it at each output time (assuming a fixed mesh).
ntindp = -1
idimt = 3 ! 1: scalar, 3: vector, 6/9: symm/non-symm tensor
ientla = .true. ! interleaved
ivarpr = .true. ! defined on the parent
! Output values; as we have no cell values, we can pass a
! trivial array for those.
call post_write_var(ipart, 'face cog (parent)', idimt, ientla, ivarpr, &
ntindp, ttcabs, rvoid, cdgfac, cdgfbo)
endif

Output on user meshes 3 or 4

Output of an interpolated velocity on user meshes

!===============================================================================
! Examples of volume variables on user meshes 3 or 4
!===============================================================================
else if (ipart.ge.3 .and. ipart.le.4) then
! Map field arrays
call field_get_val_v(ivarfl(iu), cvar_vel)
! Output of the velocity
! ----------------------
! Compute variable values on interior faces.
! In this example, we use a simple linear interpolation.
! For parallel calculations, if neighbors are used, they must be synchronized
! first. This also applies for periodicity.
if (irangp.ge.0.or.iperio.eq.1) then
call synvin(cvar_vel)
endif
allocate(vfac(3,nfacps), vfbr(3,nfbrps))
do iloc = 1, nfacps
ifac = lstfac(iloc)
ii = ifacel(1, ifac)
jj = ifacel(2, ifac)
pnd = pond(ifac)
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
! Compute variable values on boundary faces.
! In this example, we use a simple copy of the adjacent cell value.
do iloc = 1, nfbrps
ifac = lstfbr(iloc)
ii = ifabor(ifac)
vfbr(1,iloc) = cvar_vel(1,ii)
vfbr(2,iloc) = cvar_vel(2,ii)
vfbr(3,iloc) = cvar_vel(3,ii)
enddo
idimt = 3 ! 1: scalar, 3: vector, 6/9: symm/non-symm tensor
ientla = .true. ! interleaved
ivarpr = .false. ! defined on the work array
! Output values; as we have no cell values, we can pass a
! trivial array for those.
call post_write_var(ipart, 'Velocity', idimt, ientla, ivarpr, &
ntcabs, ttcabs, rvoid, vfac, vfbr)
deallocate(vfac, vfbr)

Output of the pressure on user meshes

! Output of the pressure
! ----------------------
! Variable number
ivar = ipr
call field_get_val_s(ivarfl(ivar), cvar_var)
! Compute variable values on interior faces.
! In this example, we use a simple linear interpolation.
! For parallel calculations, if neighbors are used, they must be synchronized
! first. This also applies for periodicity.
if (irangp.ge.0.or.iperio.eq.1) then
call synsca(cvar_var)
endif
allocate(sfac(nfacps), sfbr(nfbrps))
do iloc = 1, nfacps
ifac = lstfac(iloc)
ii = ifacel(1, ifac)
jj = ifacel(2, ifac)
pnd = pond(ifac)
sfac(iloc) = pnd * cvar_var(ii) &
+ (1.d0 - pnd) * cvar_var(jj)
enddo
! Compute variable values on boundary faces.
! In this example, we use a simple copy of the adjacent cell value.
do iloc = 1, nfbrps
ifac = lstfbr(iloc)
ii = ifabor(ifac)
sfbr(iloc) = cvar_var(ii)
enddo
idimt = 1 ! 1: scalar, 3: vector, 6/9: symm/non-symm tensor
ientla = .true. ! interleaved
ivarpr = .false. ! defined on the work array
! Output values; as we have no cell values, we can pass a
! trivial array for those.
call post_write_var(ipart, 'Pressure', idimt, ientla, ivarpr, &
ntcabs, ttcabs, rvoid, sfac, sfbr)
deallocate(sfac, sfbr)
endif ! end of test on post-processing mesh number