<p>Dear Neks,</p>
<p>I want to calculate the probability density function of a scalar field ff in the attached subroutine.<br />Since the element mesh is nonuniform I have to assign a weight "wght" with each grid point of ff.<br />The weight is given by the volume of the associated grid cell. Coordinates of x,y,z of all vertices<br />(primary and secondary nodes) are stored in xm1,ym1,zm1 as far as I know.  <br />This weight turns out to have in parts negative values which is unphysical. What went wrong?</p>
<p> </p>
<p>      subroutine pdf_calc(ff,step,i_offset,i_name)<br />      include 'SIZE'<br />      include 'TOTAL'</p>
<p>      parameter(npdf=501)<br /><br />      real ff(lx1,ly1,lz1,lelt)<br />      real pdf(npdf)<br />      real work(npdf)<br />      real val, vol, offset, wght<br /><br />      integer e,eg,ex,ey,ez,f<br /><br />!-----Set arrays to zero<br />      call rzero(pdf,npdf)<br />      call rzero(work,npdf)<br /><br />!-----Offset<br />      if(i_offset==0) offset=0.0<br />      if(i_offset==1) offset=int(npdf/2)*step<br /><br />!-----Pick face 5<br />      f = 5<br /><br />      vol=atan(1.0)<br /><br />      do e=1,nelv<br />        do k=1,nz1<br />          do j=1,ny1<br />            do i=1,nx1<br />               wght= (zm1(i,j,k+1,e)-zm1(i,j,k,e)) * area(i,j,f,e)          <-----This seems to be not correct and result to negative values<br />               val= ff(i,j,k,e)<br /><br />               do ipdf=1,npdf<br />                dm1=(ipdf -1 )*step-offset<br />                dm2= ipdf     *step-offset<br /><br />                if((val.ge.dm1).and.(val.lt.dm2))then<br />                  pdf(ipdf)=pdf(ipdf)+wght<br />                endif<br /><br />               enddo<br />            enddo<br />          enddo<br />        enddo<br />      enddo<br /><br />      do i=1,npdf<br />        pdf(i)=pdf(i)/vol<br />      enddo<br /><br /><br />!-----Gather over all processes (-> mpi_allreduce)<br />      call gop(pdf,work,'+  ',npdf)<br /><br />!------------------------------------------------------------<br />      if(nid.eq.0)then<br /><br />        if(i_name.eq.1)OPEN(10,file="pdf_uzte.dat",position="append")<br />        if(i_name.eq.2)OPEN(10,file="pdf_epst.dat",position="append")<br />        if(i_name.eq.3)OPEN(10,file="pdf_epsv.dat",position="append")<br />        if(i_name.eq.4)OPEN(10,file="pdf_temp.dat",position="append")<br />        if(i_name.eq.5)OPEN(10,file="pdf_dtdz.dat",position="append")<br /><br />        do ipdf=1,npdf<br />          write(10,*) (ipdf-1)*step-offset, pdf(ipdf)<br />        enddo<br />        CLOSE(10)<br /><br />      endif<br /><br />      return<br />      end<br /><br /></p>
<p>Thanks in advance and best wishes, Joerg.</p>
<div>
<pre>----------------------------------------------------
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@tu-ilmenau.de
http://www.tu-ilmenau.de/tsm
Phone: +49-3677-69-2428
Fax:   +49-3677-69-2411</pre>
</div>