[Nek5000-users] History points
nek5000-users at lists.mcs.anl.gov
nek5000-users at lists.mcs.anl.gov
Thu May 6 08:49:41 CDT 2010
Hi Guys,
yes the history point module needs a major overhaul. I don't like the way it works!
In fact we could do a much better job using our interpolation tool.
I just coded something up (without testing the code snipplet) very quicky but I think I should work out of the box.
How does it work:
- create a 'hpts.in' file
line 1: number of points
line 2: X,Y,Z coordinates of a point 1
....
line n+1: X,Y,Z coordinates of point n
- set lpart >= n in the SIZE file
- call hpts() in userchk()
The results will be written to 'hpts.out' .
Enjoy
-Stefan
----------------------------------------------------------------------------------------------------------------------
subroutine hpts
c
c dump velocity, temperature and ps-scalars for a given list
c of points (read from file hpts.in) into a file (hpts.out)
c
INCLUDE 'SIZE'
INCLUDE 'TOTAL'
parameter(nmax=lpart,nfldmax=ldim+ldimt)
parameter(mi=4,mr=1+2*ldim+nfldmax)
real rTL(mr,nmax)
integer iTL(mi,nmax)
common /itlcb/ iTL
common /rtlcb/ rTL
common /outtmp / wrk(lx1,ly1,lz1,lelt,nfldmax)
character*80 fname
integer icalld,npoints,icount
save icalld,npoints,icount
data icalld /0/
data npoints /0/
data icount /0/
nxyz = nx1*ny1*nz1
ntot = nxyz*nelt
if(icalld.eq.0) then
icalld = 1
npoints = 0
write(6,*) 'reading hpts.in'
open(50,file='hpts.in',status='old')
read(50,*) npoints
write(6,*) 'found ', npoints, ' points'
do i = 1,npoints
read(50,*) (rTL(1+j,i),j=1,ndim)
enddo
close(50)
bb_tol = 1e-8 ! bounding-box test tolerance
call intpts_setup(bb_tol)
open(50,file='hpts.out',status='new')
write(50,'(A)') '# time vx vy [vz] T PS1 PS2 ...'
endif
nflds = nfield + ndim-1 ! number of fields you want to interpolate
! pack working array
call copy(wrk(1,1,1,1,1),vx,ntot)
call copy(wrk(1,1,1,1,2),vy,ntot)
if(if3d) call copy(wrk(1,1,1,1,2),vz,ntot)
do i = 1,nfield-1
call copy(wrk(1,1,1,1,ndim+i),T(1,1,1,1,i),ntot)
enddo
! interpolate
call intpts(wrk,nflds,iTL,mi,rTL,mr,npoints,nmax)
! write interpolation results to file
if(nid.eq.0) then
do ip = 1,npoints
write(50,'(1p20E15.7)') time,
& (rTL(1+2*ndim+ifld,ip), ifld=1,nflds)
enddo
endif
return
end
On May 6, 2010, at 3:24 PM, nek5000-users at lists.mcs.anl.gov wrote:
> P.S.
>
> Yes, Frank, the exact coordinates of the history points are written only when the code rans.
>
> But if for history points you choose the collocation points on the boundary of a spectral element -- those that have I,J,K equal to either 1 or lx1 -- than you can find out their coordinates from .rea ASCII file after the line
>
> *** MESH DATA ***
>
> where coordinates of the element edges are written.
>
> Aleks
>
>
>
>
>
> ----- Original Message -----
> From: nek5000-users at lists.mcs.anl.gov
> To: nek5000-users at lists.mcs.anl.gov
> Sent: Wednesday, May 5, 2010 4:45:11 PM GMT -06:00 US/Canada Central
> Subject: Re: [Nek5000-users] History points
>
> Hi Aleks,
>
> My question is what is the logfile? Perhaps you mean the output that
> NEK writes to the screen? If so, then one would know the location of
> the point at which data is being dumped only after entering the i,j,k
> and element # of the point in the *.rea file....
>
> Cheers,
> Frank
>
> On Wed, 2010-05-05 at 23:37 +0200, nek5000-users at lists.mcs.anl.gov
> wrote:
>> Thanks Aleks.
>> Concerning grep, it is a linux command that searches for the string
>> History in the logfile.
>> _______________________________________________
>> Nek5000-users mailing list
>> Nek5000-users at lists.mcs.anl.gov
>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
> --
> Frank Herbert Muldoon, Ph.D. Mechanical Engineering
> Technische Universität Wien (Technical University of Vienna)
> Inst. f. Strömungsmechanik und Wärmeübertragung (Institute of Fluid
> Mechanics and Heat Transfer)
> Resselgasse 3
> 1040 Wien
> Tel: +4315880132232
> Fax: +4315880132299
> Cell:+436765203470
> fmuldoo (skype)
> http://tetra.fluid.tuwien.ac.at/fmuldoo/public_html/webpage/frank-muldoon.html
>
> _______________________________________________
> Nek5000-users mailing list
> Nek5000-users at lists.mcs.anl.gov
> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
> _______________________________________________
> Nek5000-users mailing list
> Nek5000-users at lists.mcs.anl.gov
> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
More information about the Nek5000-users
mailing list