[Nek5000-users] Robin boundary condition on velocity

nek5000-users at lists.mcs.anl.gov nek5000-users at lists.mcs.anl.gov
Tue Apr 25 03:17:27 CDT 2017


Hi Francesco,

BCNEUSC computes the contribution of Robin boundary condition to the 
diagonal components of the matrix (if ITYPE=-1), and to the RHS (if 
ITYPE=1). You can read the nek5000 manual about Newton cooling BC, and 
you will understand immediately what is being done.

In your specific case, you do not have an RHS contribution. You can see 
that this subroutine BCNEUSC computes matrix "S" which is added to "H2" 
(see subroutine CDSCAL). In your case, you need to do the same thing for 
Navier-Stokes. If you have Robin BC for vx and other BC for vy, you need 
to define H2 for each direction.

Note that this works only when IFSTRS=false your rea file. Hope this helps!

Sudhakar.


On 04/24/2017 05:20 PM, nek5000-users at lists.mcs.anl.gov wrote:
> Hi Sudhakar,
> thank you for your reply.
>
> It seems to me that subroutine BCNEUSC is in charge of reading the 
> "HC" value for all element faces withstanding  the "c  " boundary 
> condition.
> Nevertheless there is no info about how to Robin boundary condition is 
> effectively implemented within the solver.
>
> May I ask you to give me some more suggestions, in order to implement 
> such BC on a velocity component?
>
> Best regards,
>
> Francesco
>
> 2017-04-22 12:36 GMT+02:00 <nek5000-users-request at lists.mcs.anl.gov 
> <mailto:nek5000-users-request at lists.mcs.anl.gov>>:
>
>     Send Nek5000-users mailing list submissions to
>     nek5000-users at lists.mcs.anl.gov
>     <mailto:nek5000-users at lists.mcs.anl.gov>
>
>     To subscribe or unsubscribe via the World Wide Web, visit
>     https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>     <https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users>
>     or, via email, send a message with subject or body 'help' to
>     nek5000-users-request at lists.mcs.anl.gov
>     <mailto:nek5000-users-request at lists.mcs.anl.gov>
>
>     You can reach the person managing the list at
>     nek5000-users-owner at lists.mcs.anl.gov
>     <mailto:nek5000-users-owner at lists.mcs.anl.gov>
>
>     When replying, please edit your Subject line so it is more specific
>     than "Re: Contents of Nek5000-users digest..."
>
>
>     Today's Topics:
>
>        1. Robin boundary condition on velocity
>           (nek5000-users at lists.mcs.anl.gov
>     <mailto:nek5000-users at lists.mcs.anl.gov>)
>        2. Re: Robin boundary condition on velocity
>           (nek5000-users at lists.mcs.anl.gov
>     <mailto:nek5000-users at lists.mcs.anl.gov>)
>
>
>     ----------------------------------------------------------------------
>
>     Message: 1
>     Date: Fri, 21 Apr 2017 20:35:23 +0200
>     From: nek5000-users at lists.mcs.anl.gov
>     <mailto:nek5000-users at lists.mcs.anl.gov>
>     To: nek5000-users at lists.mcs.anl.gov
>     <mailto:nek5000-users at lists.mcs.anl.gov>
>     Subject: [Nek5000-users] Robin boundary condition on velocity
>     Message-ID:
>            
>     <mailman.5414.1492799726.2967.nek5000-users at lists.mcs.anl.gov
>     <mailto:mailman.5414.1492799726.2967.nek5000-users at lists.mcs.anl.gov>>
>     Content-Type: text/plain; charset="utf-8"
>
>     Hi Neks
>     I'm trying to implement a slip-length (Robin) as boundary
>     condition for a
>     plane channel flow, namely:
>     (u+du/dy)_{wall}=0
>     where x and y are the streamwise and wall normal directions
>     respectively.
>
>     Following the latter post:
>     http://lists.mcs.anl.gov/pipermail/nek5000-users/2011-April/001313.html
>     <http://lists.mcs.anl.gov/pipermail/nek5000-users/2011-April/001313.html>
>
>     here there's what I came up with:
>
>     c-----------------------------------------------------------------------
>           subroutine userchk
>           include 'SIZE'
>           include 'TOTAL'
>           common /my_DEBUG/ rmodulate
>           common /my_ROBIN/ ROBIN_mask(lx1,ly1,lz1,lelt)
>          $                ,      work1(lx1,ly1,lz1,lelt)
>          $                ,      work2(lx1,ly1,lz1,lelt)
>          $                ,     graduy(lx1,ly1,lz1,lelt)
>           n = nx1*ny1*nz1*nelv
>     !!!Build Up MASKS!
>           if (istep.eq.0) then
>             do i=1,n
>               x=xm1(i,1,1,1)
>               y=ym1(i,1,1,1)
>               if(x.gt.3.and.x.lt.4)then
>                 ROBIN_mask(i,1,1,1)=1.0
>               else
>                 ROBIN_mask(i,1,1,1)=0.0
>               endif
>             enddo
>           endif
>     !Compute velocity gradient
>           call rzero(graduy,n)
>           call gradm1(work1,graduy,work2,vx)
>           call dsavg(graduy)
>     !      call outpost(work1,graduy,work2,pr,t,'GRD')
>     !Da aggiungere l'eventualit? del 3D...
>            rmodulate = 0.5*(tanh(0.25*time-10)+1)
>           if (nid.eq.0) write(6,*)"MODULATE", rmodulate,time
>           return
>           end
>     c-----------------------------------------------------------------------
>     c-----------------------------------------------------------------------
>           subroutine userbc (ix,iy,iz,iside,eg)
>     c     NOTE ::: This subroutine MAY NOT be called by every process
>           include 'SIZE'
>           include 'TOTAL'
>           include 'NEKUSE'
>           integer e, eg
>     !Imposing Fancy Boundary Conditions...
>     ![Nek5000-users] Robin boundary conditions for the velocity
>           common /my_ROBIN/ ROBIN_mask(lx1,ly1,lz1,lelt)
>          $                ,      work1(lx1,ly1,lz1,lelt)
>          $                ,      work2(lx1,ly1,lz1,lelt)
>          $                ,     graduy(lx1,ly1,lz1,lelt)
>           common /my_DEBUG/ rmodulate
>           e=gllel(eg) ! global element number to processor-local el. #
>           ux=1*graduy(ix,iy,iz,e)*rmodulate
>           uy=0.0
>           uz=0.0
>           return
>           end
>     c-----------------------------------------------------------------------
>           subroutine useric (ix,iy,iz,ieg)
>           include 'SIZE'
>           include 'TOTAL'
>           include 'NEKUSE'
>           ux=1.0*y**2
>           uy=0.0
>           uz=0.0
>           temp=0
>           return
>           end
>     c-----------------------------------------------------------------------
>
>     The code becomes unstable in the gradient computation due to some
>     point-to-point oscillations.
>     As you can see I've also tried to use a ramp function in order to
>     make the
>     onset of Robin condition more gentle,
>     but at a given "slip-length amplitude" oscillations deteriorates the
>     solution.
>     I use ROBIN_mask in order to chose where to have such BC, but it
>     seems that
>     I'm still missing something in computing the velocity gradient GRD.
>
>     By the way, this is my box.BOX file:
>
>     base.rea
>     2                      spatial dimension  ( < 0 --> generate
>     .rea/.re2 pair)
>     1                      number of fields
>     #=======================================================================
>     #
>     #    Example of .box file for channel flow
>     #
>     #    If nelx (y or z) < 0, then genbox automatically generates the
>     #                          grid spacing in the x (y or z) direction
>     #                          with a geometric ratio given by "ratio".
>     #                          ( ratio=1 implies uniform spacing )
>     #
>     #    Note that the character bcs _must_ have 3 spaces.
>     #
>     #=======================================================================
>     #
>     Box
>     -32  -16 nelx,nely,nelz for Box
>     0  6.2832 1.  x0,x1,gain
>     -1  1        1.  y0,y1,gain
>     P  ,P  ,v  ,v  ,   ,                                bc's  (3 chars
>     each!)
>
>
>     I was wondering if you might be able to give me some advice.
>
>     Best regards,
>
>     Francesco
>     -------------- next part --------------
>     An HTML attachment was scrubbed...
>     URL:
>     <http://lists.mcs.anl.gov/pipermail/nek5000-users/attachments/20170421/fcdeceb5/attachment-0001.html
>     <http://lists.mcs.anl.gov/pipermail/nek5000-users/attachments/20170421/fcdeceb5/attachment-0001.html>>
>
>     ------------------------------
>
>     Message: 2
>     Date: Sat, 22 Apr 2017 12:35:45 +0200
>     From: nek5000-users at lists.mcs.anl.gov
>     <mailto:nek5000-users at lists.mcs.anl.gov>
>     To: nek5000-users at lists.mcs.anl.gov
>     <mailto:nek5000-users at lists.mcs.anl.gov>
>     Subject: Re: [Nek5000-users] Robin boundary condition on velocity
>     Message-ID:
>            
>     <mailman.5458.1492857360.2967.nek5000-users at lists.mcs.anl.gov
>     <mailto:mailman.5458.1492857360.2967.nek5000-users at lists.mcs.anl.gov>>
>     Content-Type: text/plain; charset="utf-8"
>
>     Hi Francesco,
>
>     Robin boundary condition is implemented already for the passive
>     scalar, in the form of Newton cooling boundary condition. You can
>     see how this is implemented in subroutine BCNEUSC. In your case
>     ITYPE=-1 and CB=?c  ' in the subroutine, and you need to only set
>     ?HC? accordingly. We tested this and even our turbulent
>     simulations are stable with Robin BC.
>
>     Best regards,
>     Sudhakar.
>
>     > On 21 Apr 2017, at 20:35, nek5000-users at lists.mcs.anl.gov
>     <mailto:nek5000-users at lists.mcs.anl.gov> wrote:
>     >
>     > Hi Neks
>     > I'm trying to implement a slip-length (Robin) as boundary
>     condition for a plane channel flow, namely:
>     > (u+du/dy)_{wall}=0
>     > where x and y are the streamwise and wall normal directions
>     respectively.
>     >
>     > Following the latter post:
>     >
>     http://lists.mcs.anl.gov/pipermail/nek5000-users/2011-April/001313.html
>     <http://lists.mcs.anl.gov/pipermail/nek5000-users/2011-April/001313.html>
>     <http://lists.mcs.anl.gov/pipermail/nek5000-users/2011-April/001313.html
>     <http://lists.mcs.anl.gov/pipermail/nek5000-users/2011-April/001313.html>>
>     >
>     > here there's what I came up with:
>     >
>     >
>     c-----------------------------------------------------------------------
>     >       subroutine userchk
>     >       include 'SIZE'
>     >       include 'TOTAL'
>     >       common /my_DEBUG/ rmodulate
>     >       common /my_ROBIN/ ROBIN_mask(lx1,ly1,lz1,lelt)
>     >      $                ,      work1(lx1,ly1,lz1,lelt)
>     >      $                ,      work2(lx1,ly1,lz1,lelt)
>     >      $                ,     graduy(lx1,ly1,lz1,lelt)
>     >       n = nx1*ny1*nz1*nelv
>     > !!!Build Up MASKS!
>     >       if (istep.eq.0) then
>     >         do i=1,n
>     >           x=xm1(i,1,1,1)
>     >           y=ym1(i,1,1,1)
>     >           if(x.gt.3.and.x.lt.4)then
>     >             ROBIN_mask(i,1,1,1)=1.0
>     >           else
>     >             ROBIN_mask(i,1,1,1)=0.0
>     >           endif
>     >         enddo
>     >       endif
>     > !Compute velocity gradient
>     >       call rzero(graduy,n)
>     >       call gradm1(work1,graduy,work2,vx)
>     >       call dsavg(graduy)
>     > !      call outpost(work1,graduy,work2,pr,t,'GRD')
>     > !Da aggiungere l'eventualit? del 3D...
>     >        rmodulate = 0.5*(tanh(0.25*time-10)+1)
>     >       if (nid.eq.0) write(6,*)"MODULATE", rmodulate,time
>     >       return
>     >       end
>     >
>     c-----------------------------------------------------------------------
>     >
>     c-----------------------------------------------------------------------
>     >       subroutine userbc (ix,iy,iz,iside,eg)
>     > c     NOTE ::: This subroutine MAY NOT be called by every process
>     >       include 'SIZE'
>     >       include 'TOTAL'
>     >       include 'NEKUSE'
>     >       integer e, eg
>     > !Imposing Fancy Boundary Conditions...
>     > ![Nek5000-users] Robin boundary conditions for the velocity
>     >       common /my_ROBIN/ ROBIN_mask(lx1,ly1,lz1,lelt)
>     >      $                ,      work1(lx1,ly1,lz1,lelt)
>     >      $                ,      work2(lx1,ly1,lz1,lelt)
>     >      $                ,     graduy(lx1,ly1,lz1,lelt)
>     >       common /my_DEBUG/ rmodulate
>     >       e=gllel(eg) ! global element number to processor-local el. #
>     >       ux=1*graduy(ix,iy,iz,e)*rmodulate
>     >       uy=0.0
>     >       uz=0.0
>     >       return
>     >       end
>     >
>     c-----------------------------------------------------------------------
>     >       subroutine useric (ix,iy,iz,ieg)
>     >       include 'SIZE'
>     >       include 'TOTAL'
>     >       include 'NEKUSE'
>     >       ux=1.0*y**2
>     >       uy=0.0
>     >       uz=0.0
>     >       temp=0
>     >       return
>     >       end
>     >
>     c-----------------------------------------------------------------------
>     >
>     > The code becomes unstable in the gradient computation due to
>     some point-to-point oscillations.
>     > As you can see I've also tried to use a ramp function in order
>     to make the onset of Robin condition more gentle,
>     > but at a given "slip-length amplitude" oscillations deteriorates
>     the solution.
>     > I use ROBIN_mask in order to chose where to have such BC, but it
>     seems that I'm still missing something in computing the velocity
>     gradient GRD.
>     >
>     > By the way, this is my box.BOX file:
>     >
>     > base.rea
>     > 2                      spatial dimension  ( < 0 --> generate
>     .rea/.re2 pair)
>     > 1                      number of fields
>     >
>     #=======================================================================
>     > #
>     > #    Example of .box file for channel flow
>     > #
>     > #    If nelx (y or z) < 0, then genbox automatically generates the
>     > #                          grid spacing in the x (y or z) direction
>     > #                          with a geometric ratio given by "ratio".
>     > #                          ( ratio=1 implies uniform spacing )
>     > #
>     > #    Note that the character bcs _must_ have 3 spaces.
>     > #
>     >
>     #=======================================================================
>     > #
>     > Box
>     > -32  -16 nelx,nely,nelz for Box
>     > 0  6.2832 1.  x0,x1,gain
>     > -1  1        1.  y0,y1,gain
>     > P  ,P  ,v  ,v  ,   , bc's  (3 chars each!)
>     >
>     >
>     > I was wondering if you might be able to give me some advice.
>     >
>     > Best regards,
>     >
>     > Francesco
>     > _______________________________________________
>     > Nek5000-users mailing list
>     > Nek5000-users at lists.mcs.anl.gov
>     <mailto:Nek5000-users at lists.mcs.anl.gov>
>     > https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>     <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/20170422/42c96cbe/attachment.html
>     <http://lists.mcs.anl.gov/pipermail/nek5000-users/attachments/20170422/42c96cbe/attachment.html>>
>
>     ------------------------------
>
>     _______________________________________________
>     Nek5000-users mailing list
>     Nek5000-users at lists.mcs.anl.gov
>     <mailto:Nek5000-users at lists.mcs.anl.gov>
>     https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>     <https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users>
>
>
>     End of Nek5000-users Digest, Vol 98, Issue 10
>     *********************************************
>
>
>
>
> _______________________________________________
> 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/20170425/cd652eae/attachment-0001.html>


More information about the Nek5000-users mailing list