[Nek5000-users] multiple runs and poisson for velocity

nek5000-users at lists.mcs.anl.gov nek5000-users at lists.mcs.anl.gov
Fri Apr 4 12:40:07 CDT 2014


Hi Giuseppe,

Regarding the first question, you need to change not the param(2) but copy the value of dynamic viscosity (which 1/Re for non-dimentional case with rho=1) into an array

VDIFF (1,1,1,1,1)

as in vprops() of subs1.f


...
      DO 1000 IEL=1,NEL
...
         ELSE IF(MATYPE(IGRP,IFIELD).EQ.0)THEN
C
C           Default constant property
C
            CDIFF  = CPFLD(IFIELD,1)
...
            CALL CFILL(VDIFF (1,1,1,IEL,IFIELD),CDIFF,NXYZ1)


where cpfld(1,1) array is filled in rdparam() of connect2.f during the ininitialization:

       CPFLD(1,1)=PARAM(2)

Aleks







________________________________________
From: nek5000-users-bounces at lists.mcs.anl.gov [nek5000-users-bounces at lists.mcs.anl.gov] on behalf of nek5000-users at lists.mcs.anl.gov [nek5000-users at lists.mcs.anl.gov]
Sent: Friday, April 04, 2014 7:09 AM
To: nek5000-users at lists.mcs.anl.gov
Subject: [Nek5000-users] multiple runs and poisson for velocity

Dear Neks,
I have two issues that I would like to discuss with you.
Firstly, I would like to perform many steady state simulations with different values of the viscosity, each time using the previous solution as starting condition.
For this reason, I have built the following subroutine:

c-----------------------------------------------------------------------
      subroutine userchk
      include 'SIZE'
      include 'TOTAL'


      do 120 ire=1,10
      param(2) = -1*ire
      call setprop
      igeom=1
      call fluid(igeom)
120 continue

      ifxyo = .true.
      call outpost(vx,vy,vz,pr,t,'ste')
      call exitt

      return
      end
c-----------------------------------------------------------------------

but this does not work as I expected, could somebody suggest me what I should do instead?

Secondly, I would like to solve (in postprocessing mode) a Poisson problem for the velocity components, in the form -lapl(vx) = fx and -lapl(vy)=fy, with fx and fy known.
Is it possible to do this within NEK? I think this could be done by adapting the Uzawa pressure solver in navier1.f, but I am not sure if this could work.
Best regards,
Giuseppe
_______________________________________________
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