[Nek5000-users] Temperature gradient at a point
nek5000-users at lists.mcs.anl.gov
nek5000-users at lists.mcs.anl.gov
Mon Aug 23 16:31:37 CDT 2010
Hi Pradeep,
I can show you how to do this off-list, there is a way to use the uservp (variable properties) routine in the usr file. Although I am not sure
how this effects using the IFCHAR (if you are). I seem to recall Markus having an issue with PN/PN versus PN/PN-2 when using
variable properties, and the PN/PN-2 routine seemed to correct the issue he was having.
- Michael
----- Original Message -----
From: nek5000-users at lists.mcs.anl.gov
To: nek5000-users at lists.mcs.anl.gov
Sent: Monday, August 23, 2010 4:13:48 PM GMT -06:00 US/Canada Central
Subject: Re: [Nek5000-users] Temperature gradient at a point
Hi Paul,
Thanks for the detailed reply. The reason I'm not solving it as a conjugate heat transfer problem, is that thermal conductivity is a function of temperature based on some curve fit equations, and I am not sure how to implement that.
Thanks,
Pradeep
On Mon, Aug 23, 2010 at 4:04 PM, < nek5000-users at lists.mcs.anl.gov > wrote:
Hi Pradeep,
Why not just solve the conjugate heat transfer problem directly
using fluid + solid elements in nek?
Also, nek supports full Robin boundary conditions if you wish
to do a Newton law of cooling: k*dT/dn . n_hat = h*(T-Tinf), where Tinf is the external temperature and h is the heat transfer coefficient, both of which can be functions of time and space.
Regarding gradm1, you would call it from userchk, and store
the output in arrays in a common block, e.g., as below.
Paul
subroutine userchk
:
common /mygrad/ tx(lx1,ly1,lz1,lelt)
$ , ty(lx1,ly1,lz1,lelt)
$ , tz(lx1,ly1,lz1,lelt)
call gradm1(tx,ty,tz,t)
:
:
subroutine userbc (ix,iy,iz,iside,eg)
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
common /mygrad/ tx(lx1,ly1,lz1,lelt)
$ , ty(lx1,ly1,lz1,lelt)
$ , tz(lx1,ly1,lz1,lelt)
integer e,eg
e = gllel(eg) ! global element number to processor-local el. #
gtx=tx(ix,iy,iz,e)
gty=ty(ix,iy,iz,e)
gtz=tz(ix,iy,iz,e)
On Mon, 23 Aug 2010, nek5000-users at lists.mcs.anl.gov wrote:
Hi Paul,
I am basically trying to solve a conjugate heat transfer problem in an
iterative manner, for flow over an infinitely long cylinder (2D).
I need to use the heat transfer at the boundary, to calculate the new
temperature at the boundary for the next time step. The
temperature for the next time step is solved for using this heat flux, by a
function in the usr file using an FEM algorithm for the solid part
(cylinder). The bc type I am using is Temperature - fortran function.
Regards,
Pradeep
On Mon, Aug 23, 2010 at 2:50 PM, < nek5000-users at lists.mcs.anl.gov > wrote:
Pradeep,
if you give me some idea of the nature of your bc, I can
perhaps help --- there are a large number of bc types already
supported inside nek
Paul
On Mon, 23 Aug 2010, nek5000-users at lists.mcs.anl.gov wrote:
Hi,
I wanted to know if there was a way to find the temperature gradient at a
point. I need that information in the userbc function.
I tried using gradm1(), but I am not sure how to get the value at a given
point.
Thanks,
Pradeep
_______________________________________________
Nek5000-users mailing list
Nek5000-users at lists.mcs.anl.gov
https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
--
Pradeep C. Rao
Graduate Research Assistant for FT2L ( http://www1.mengr.tamu.edu/FT2L/ )
Department of Mechanical Engineering
Texas A&M University
College Station, TX 77843-3123
428 Engineering Physics Building
(713) 210-9769
uuuu
c-----------------------------------------------------------------------
C
C USER SPECIFIED ROUTINES:
C
C - boundary conditions
C - initial conditions
C - variable properties
C - local acceleration for fluid (a)
C - forcing function for passive scalar (q)
C - general purpose routine for checking errors etc.
C
c-----------------------------------------------------------------------
subroutine uservp (ix,iy,iz,eg)
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
integer e,f,eg
c e = gllel(eg)
udiff =0.
utrans=0.
return
end
c-----------------------------------------------------------------------
subroutine userf (ix,iy,iz,eg)
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
integer e,f,eg
c e = gllel(eg)
c Note: this is an acceleration term, NOT a force!
c Thus, ffx will subsequently be multiplied by rho(x,t).
ffx = 0.0
ffy = 0.0
ffz = 0.0
return
end
c-----------------------------------------------------------------------
subroutine userq (ix,iy,iz,eg)
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
integer e,f,eg
c e = gllel(eg)
qvol = 0.0
return
end
c-----------------------------------------------------------------------
subroutine userchk
include 'SIZE'
include 'TOTAL'
return
end
c-----------------------------------------------------------------------
subroutine userbc (ix,iy,iz,iside,ieg)
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
ux=0.0
uy=0.0
uz=0.0
temp=0.0
return
end
c-----------------------------------------------------------------------
subroutine useric (ix,iy,iz,ieg)
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
ux=0.0
uy=0.0
uz=0.0
temp=0
return
end
c-----------------------------------------------------------------------
subroutine usrdat
include 'SIZE'
include 'TOTAL'
c
return
end
c-----------------------------------------------------------------------
subroutine usrdat2
include 'SIZE'
include 'TOTAL'
param(66) = 4. ! These give the std nek binary i/o and are
param(67) = 4. ! good default values
return
end
c-----------------------------------------------------------------------
subroutine usrdat3
include 'SIZE'
include 'TOTAL'
c
return
end
c-----------------------------------------------------------------------
_______________________________________________
Nek5000-users mailing list
Nek5000-users at lists.mcs.anl.gov
https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
--
Pradeep C. Rao
Graduate Research Assistant for FT2L ( http://www1.mengr.tamu.edu/FT2L/ )
Department of Mechanical Engineering
Texas A&M University
College Station, TX 77843-3123
428 Engineering Physics Building
(713) 210-9769
_______________________________________________ Nek5000-users mailing list Nek5000-users at lists.mcs.anl.gov https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/nek5000-users/attachments/20100823/ac2b5f82/attachment.html>
More information about the Nek5000-users
mailing list