Horizontal average/integral of a variable for domain with topography

Barbara Zemskova

unread,

Nov 29, 2019, 4:13:39 PM11/29/19

to Nek5000

Hello all,

I have a 3D domain with bottom hill topography, and I am trying to extract 1D profile of a horizontal integral (over xy-plane) as a function of z for some variables (dissipation, for example). I have been successful at using planar_avg routine, but that outputs a 3D array, not 1D.

I looked up previous examples that have been given in forums, but the results that I am getting don't make sense - I tried to get the average z value, i.e. Z(z) as well, and the output that I get makes no sense given my domain. Below is part of the code that I am using in userchk().

I would appreciate if someone could point out where my problems are with computing these averages.

Thank you!

Barbara

common /fbar/ fbar(lz1,lelz),wght(lz1,lelz)

common /depth/ dpth(lz1,lelz),work2(lz1,lelz)

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

ntot = nx1*ny1*nz1*nelv

nelx = 32

nely = 32

nelz = 64

f = 5

nelxy = nelgv/nelz

call rzero(fbar,nz1*nelz)

call rzero(wght,nz1*nelz)

call rzero(dpth,nz1*nelz)

do e=1,nelv

eg = lglel(e)

call get_exyz(ex,ey,ez,eg,nelxy,1,nelz)

do k=1,nz1

do in=1,nx1*ny1

fbar(k,ez) = fbar(k,ez)+area(in,1,f,e)*dke2(in,1,k,e)

wght(k,ez) = wght(k,ez)+area(in,1,f,e)

dpth(k,ez) = dpth(k,ez)+zm1(in,1,k,e)

enddo

call gop(fbar,work2,'+ ',nz1*nelz)

call gop(wght,work2,'+ ',nz1*nelz)

call gop(dpth,work2,'+ ',nz1*nelz)

do in=1,nz1*nelz

dpth(in,1)=dpth(in,1)/wght(in,1)

enddo

do k=1,nelz

do in=1,nz1

write(1,*) dpth(in,k),wght(in,k),fbar(in,k)

enddo

Fischer, Paul

unread,

Nov 30, 2019, 9:26:21 AM11/30/19

to Barbara Zemskova, Nek5000

Hi Barbara,

It looks like your dpth is not multiplied by area() ?

dpth(k,ez) = dpth(k,ez)+zm1(in,1,k,e)*area(in,1,f,e)

and your fbar is not divided by wght ?

Best,

Paul

From: nek...@googlegroups.com <nek...@googlegroups.com> on behalf of Barbara Zemskova <b.e.ze...@gmail.com>
Sent: Friday, November 29, 2019 3:13 PM
To: Nek5000 <nek...@googlegroups.com>
Subject: [nek5000] Horizontal average/integral of a variable for domain with topography

--
You received this message because you are subscribed to the Google Groups "Nek5000" group.
To unsubscribe from this group and stop receiving emails from it, send an email to nek5000+u...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/nek5000/b6fefe19-3945-45b6-9efa-5e90d9cbf3e6%40googlegroups.com.

Barbara Zemskova

unread,

Dec 3, 2019, 12:21:48 PM12/3/19

to Nek5000

Hi Paul,

Thank you very much for looking into the issue. I tried incorporating your suggestions (please see code below), however, I am still getting non-sensical results.

I am attaching a file with my output as well. The domain size is (x,y,z) = (0.5,0.5,1), so dpth should be between [0,1] and the horizontal area should be around 0.25 or smaller (at depths with topography). Velocity (vx) ranges from -0.25 to 0.25. However, those are not at all the numbers that I am getting.

I would appreciate any input regarding this matter. Thank you for your help again.

-Barbara

Code:

call rzero(fbar,nz1*nelz)

call rzero(wght,nz1*nelz)

call rzero(dpth,nz1*nelz)

do e=1,nelv

eg = lglel(e)

call get_exyz(ex,ey,ez,eg,nelxy,1,nelz)

do k=1,nz1

do in=1,nx1*ny1

fbar(k,ez) = fbar(k,ez)+area(in,1,5,e)*vx(in,1,k,e)

wght(k,ez) = wght(k,ez)+area(in,1,5,e)

dpth(k,ez) = dpth(k,ez)+zm1(in,1,k,e)*area(in,1,5,e)

enddo

call gop(fbar,work2,'+ ',nz1*nelz)

call gop(wght,work2,'+ ',nz1*nelz)

call gop(dpth,work2,'+ ',nz1*nelz)

do in=1,nz1*nelz

dpth(in,1)=dpth(in,1)/wght(in,1)

fbar(in,1)=fbar(in,1)/wght(in,1)

enddo

do k=1,nelz

do in=1,nz1

write(1,*) dpth(in,k),wght(in,k),fbar(in,k)

enddo

To unsubscribe from this group and stop receiving emails from it, send an email to nek5000+unsubscribe@googlegroups.com.

MeanUProfile

Reply all

Reply to author

Forward