[Nek5000-users] Probability Density Function

[nek5000-users at lists.mcs.anl.gov]

Dear Neks, 

I want to calculate the probability density function of
a scalar field ff in the attached subroutine.
Since the element mesh is
nonuniform I have to assign a weight "wght" with each grid point of
ff.
The weight is given by the volume of the associated grid cell.
Coordinates of x,y,z of all vertices
(primary and secondary nodes) are
stored in xm1,ym1,zm1 as far as I know. 
This weight turns out to have
in parts negative values which is unphysical. What went wrong? 


subroutine pdf_calc(ff,step,i_offset,i_name)
 include 'SIZE'
 include
'TOTAL' 

 parameter(npdf=501)

 real ff(lx1,ly1,lz1,lelt)
 real
pdf(npdf)
 real work(npdf)
 real val, vol, offset, wght

 integer
e,eg,ex,ey,ez,f

!-----Set arrays to zero
 call rzero(pdf,npdf)
 call
rzero(work,npdf)

!-----Offset
 if(i_offset==0) offset=0.0

if(i_offset==1) offset=int(npdf/2)*step

!-----Pick face 5
 f = 5


vol=atan(1.0)

 do e=1,nelv
 do k=1,nz1
 do j=1,ny1
 do i=1,nx1
 wght=
(zm1(i,j,k+1,e)-zm1(i,j,k,e)) * area(i,j,f,e)  mpi_allreduce)
 call
gop(pdf,work,'+
',npdf)

!------------------------------------------------------------

if(nid.eq.0)then


if(i_name.eq.1)OPEN(10,file="pdf_uzte.dat",position="append")

if(i_name.eq.2)OPEN(10,file="pdf_epst.dat",position="append")

if(i_name.eq.3)OPEN(10,file="pdf_epsv.dat",position="append")

if(i_name.eq.4)OPEN(10,file="pdf_temp.dat",position="append")

if(i_name.eq.5)OPEN(10,file="pdf_dtdz.dat",position="append")

 do
ipdf=1,npdf
 write(10,*) (ipdf-1)*step-offset, pdf(ipdf)
 enddo

CLOSE(10)

 endif

 return
 end

Thanks in advance and best wishes,
Joerg. 

----------------------------------------------------
Joerg
Schumacher 

Heisenberg Professor for Theoretical Fluid
Mechanics
Institute of Thermodynamics and Fluid Mechanics
Department of
Mechanical Engineering
Ilmenau University of Technology
P.O. Box 100565

D-98684 Ilmenau 
Germany 

E-mail:
joerg.schumacher at tu-ilmenau.de
http://www.tu-ilmenau.de/tsm
Phone:
+49-3677-69-2428
Fax: +49-3677-69-2411
 
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/nek5000-users/attachments/20121217/8e71d170/attachment.html>