Commit 2a0b657c authored by ylan's avatar ylan
Browse files

[peb1568] add new mesh for 1568 pebbles case, run with NekRS

parent 73fd1dd5
This diff is collapsed.
## 1568 (non-touching) pebbles in a cylinder (Last updated: Mar. 02, 2020)
Peb1568 N, Version 2.0 test 1, E524386
This is a better mesh using improve technique developped within last few months.
Please note that the improvement has not finished.
- New features:
- the resolution is refined, close to the octant refinement, but the elements is only 4x
- remove/merge small facet/edge redistribute the points
- no negative Jacobain in NekRS multigrids
- better convergence in pressure solver, so NekRS can start from constant initial flow.
- Geometry and BC:
- Sphere radius = 1.0 (Wall)
- Cylinder radius = 13.858 (SYM), (convert SYM to Wall in usr file)
- Zbot of the cylinder = 17.691 (Outflow)
- Ztop of the cylinder = -14.673 (Inflow)
- We use dimless scale to run Re5000 with vz=1.0
- Case files:
- xn3b.re2, a curved (projected at lx1=3) hex20 re2 mesh, converted by a hex20 rea: xn3.rea+reatore2, R_peb=1.0
```
xyz min -13.858 -13.858 -14.673
xyz max 13.858 13.858 17.691
mesh metrics:
GLL grid spacing min/max : 8.26E-05 3.02E-01
scaled Jacobian min/max/avg: 5.87E-03 9.91E-01 3.71E-01
aspect ratio min/max/avg: 1.12E+00 4.10E+02 1.26E+01
```
module load cmake gcc cuda
. ~/.bashrc
export NEKRS_HOME=$HOME/.local/nekrs_repo0204
PATH=${NEKRS_HOME}/bin:${PATH}
export OCCA_CACHE_DIR=/gpfs/alpine/proj-shared/csc262/ylan/cache_repo0204
which nekrs
echo $OCCA_CACHE_DIR
// Boundary conditions
/* wall 1, inflow 2, outflow 3, x-slip 4, y-slip 5, z-slip 6 */
void insVelocityDirichletConditions3D(bcData *bc)
{
bc->uP = 0.0;
bc->vP = 0.0;
bc->wP = 1.0;
}
void insVelocityNeumannConditions3D(bcData *bc)
{
}
// Stabilized outflow (Dong et al)
void insPressureDirichletConditions3D(bcData *bc)
{
const dfloat iU0delta = 10.0;
const dfloat un = bc->uM*bc->nx + bc->vM*bc->ny + bc->wM*bc->nz;
const dfloat s0 = 0.5 * (1.0 - tanh(un*iU0delta));
bc->pP = -0.5 * (bc->uM*bc->uM + bc->vM*bc->vM + bc->wM*bc->wM) * s0;
}
[OCCA]
backend = CUDA
deviceNumber = 0 #LOCAL-RANK
[GENERAL]
polynomialOrder = 7
startFrom = r17.fld
stopAt = endTime
endTime = 200
dealiasing = yes
dt = 1.e-4 #1.0e-3
timeStepper = tombo2
extrapolation = subCycling
subCyclingSteps = 2
targetCFL = 3.8
writeControl = TIMESTEP
writeInterval = 1000
filtering = hpfrt
filterWeight = 250
filterCutoffRatio = 0.9
[PRESSURE]
preconditioner = semg_amg
#amgSolver = parAlmond
galerkincoarseoperator = yes
residualTol = 1e-04
residualProj = yes
[VELOCITY]
boundaryTypeMap = inlet, outlet, wall
density = 1.0
viscosity = -5000.0 #0.25 0.000405
residualTol = 1e-06
#[TEMPERATURE]
#solver=none
#boundaryTypeMap = t, O, I, t
#rhoCp = 10
#conductivity = 5000
#residualTol = 1e-08
****** PARAMETERS *****
2.60999990 NEKTON VERSION
3 DIMENSIONAL RUN
103 PARAMETERS FOLLOW
1.00000 p001 DENSITY
-1000.000 p002 VISCOS
0.00000 p003
0.00000 p004
0.00000 p005
0.00000 p006
1.00000 p007 RHOCP
-1000.00000 p008 CONDUCT
0.00000 p009
0.00000 p010 FINTIME
2000.00 p011 NSTEPS
-4.00000E-04 p012 DT
0.00000 p013 IOCOMM
0.00000 p014 IOTIME
0.00000 p015 IOSTEP
0.00000 p016 p16 PSSOLVER
0.00000 p017
0.00000 p018
0.00000 p019
0.00000 p020
0.100000E-05 p021 DIVERGENCE
0.100000E-07 p022 HELMHOLTZ
0.00000 p023 NPSCAL
0.100000E-09 p024 TOLREL
0.100000E-09 p025 TOLABS
0.250000 p026 COURANT
2.00000 p027 TORDER
0.00000 p028 TORDER: mesh velocity (0: p28=p27)
0.00000 p029 = magnetic visc if > 0, = -1/Rm if < 0
0.00000 p030 > 0 properties set in uservp()
0.00000 p031 NPERT: #perturbation modes
0.00000 p032 #BCs in re2 file, if > 0
0.00000 p033
0.00000 p034
0.00000 p035
0.00000 p036
0.00000 p037
0.00000 p038
0.00000 p039
0.00000 p040
0.00000 p041
0.00000 p042 p42 0=gmres/1=pcg
0.00000 p043 p43 0=semg/1=schwarz
0.00000 p044 p44 0=E-based/1=A-based prec.
0.00000 p045 p45 Relaxation factor for DTFS
0.00000 p046
0.00000 p047 viscosity for mesh elasticity solver
0.00000 p048
0.00000 p049
0.00000 p050
0.00000 p051
0.00000 p052 p52 IOHIS
0.00000 p053
0.00000 p054 fixed flow rate dir: |p54|=1,2,3=x,y,z
0.00000 p055 vol.flow rate (p54>0) or Ubar (p54<0)
0.00000 p056
0.00000 p057
0.00000 p058
0.00000 p059
0.00000 p060 !=0 --> init. velocity to small nonzero
0.00000 p061
0.00000 p062 >0 --> force byte_swap for output
8 0.00000 p063 =8 --> force 8-byte output
0.00000 p064 =1 --> perturbation restart
0.00000 p065 #iofiles (eg, 0 or 64); <0 --> sep. dirs
6.00000 p066 p66 write fmt:ONLY postx uses rea v
6.00000 p067
0.00000 p068 iastep: freq for avg_all (0=iostep)
0.00000 p069
0.00000 p070
0.00000 p071
0.00000 p072
0.00000 p073
0.00000 p074 p74 verbose Helmholtz
0.00000 p075
0.00000 p076
0.00000 p077
0.00000 p078
0.00000 p079
0.00000 p080
0.00000 p081
0.00000 p082
0.00000 p083
0.00000 p084 !=0 --> sets initial timestep if p12>0
0.00000 p085 dt ratio if p84 !=0, for timesteps>0
0.00000 p086 p86 reserved
0.00000 p087
0.00000 p088
0.00000 p089 p89 reserved
0.00000 p090
0.00000 p091
0.00000 p092
20.00000 p093 Number of previous pressure solns saved
5.00000 p094 start projecting velocity after p94 step
5.00000 p095 start projecting pressure after p95 step
0.00000 p096
0.00000 p097
0.00000 p098
3.00000 p099 dealiasing: <0 off
0.00000 p100
0.00000 p101 Number of additional modes to filter
0.00000 p102
0.500000E-01 p103 weight of stabilizing filter (.01)
4 Lines of passive scalar data follows2 CONDUCT; 2RHOCP
1.00000 1.00000 1.00000 1.00000 1.00000
1.00000 1.00000 1.00000 1.00000
1.00000 1.00000 1.00000 1.00000 1.00000
1.00000 1.00000 1.00000 1.00000
13 LOGICAL SWITCHES FOLLOW
T IFFLOW
T IFHEAT
T IFTRAN
T T T T T F F F F F F IFNAV & IFADVC (convection in P.S. fields)
F F F F F F T T T T T T IFTMSH (IF mesh for this field is T mesh)
F IFAXIS
T IFSTRS
F IFSPLIT
F IFMGRID
F IFMODEL
F IFKEPS
F IFMVBD
F IFCHAR
8.00000 8.00000 -4.00000 -4.00000 XFAC,YFAC,XZERO,YZERO
**MESH DATA** 6 lines are X,Y,Z;X,Y,Z. Columns corners 1-4;5-8
-524386 3 524386 NELT,NDIM,NELV
0 PRESOLVE/RESTART OPTIONS *****
7 INITIAL CONDITIONS *****
C Default
C Default
C Default
C Default
C Default
C Default
C Default
***** DRIVE FORCE DATA ***** BODY FORCE, FLOW, Q
4 Lines of Drive force data follow
C
C
C
C
***** Variable Property Data ***** Overrrides Parameter data.
1 Lines follow.
0 PACKETS OF DATA FOLLOW
***** HISTORY AND INTEGRAL DATA *****
0 POINTS. Hcode, I,J,H,IEL
***** OUTPUT FIELD SPECIFICATION *****
6 SPECIFICATIONS FOLLOW
T COORDINATES
T VELOCITY
T PRESSURE
T TEMPERATURE
F TEMPERATURE GRADIENT
0 PASSIVE SCALARS
***** OBJECT SPECIFICATION *****
0 Surface Objects
0 Volume Objects
0 Edge Objects
0 Point Objects
//
// nekRS User Defined File
//
#include <math.h>
#include "udf.hpp"
/* UDF Functions */
void UDF_LoadKernels(ins_t *ins)
{
}
void UDF_Setup(ins_t *ins)
{
// get IC from nek
if (!ins->readRestartFile) nek_copyTo(ins, ins->startTime);
}
void UDF_ExecuteStep(ins_t *ins, dfloat time, int tstep)
{
if (ins->isOutputStep) nek_userchk();
}
c-----------------------------------------------------------------------
subroutine uservp(i,j,k,eg) ! set variable properties
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
integer i,j,k,e,eg
c e = gllel(eg)
udiff = 0.0
utrans = 0.0
return
end
c-----------------------------------------------------------------------
subroutine userf(i,j,k,eg) ! set acceleration term
c
c Note: this is an acceleration term, NOT a force!
c Thus, ffx will subsequently be multiplied by rho(x,t).
c
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
integer i,j,k,e,eg
c e = gllel(eg)
ffx = 0.0
ffy = 0.0
ffz = 0.0
return
end
c-----------------------------------------------------------------------
subroutine userq(i,j,k,eg) ! set source term
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
integer i,j,k,e,eg
c e = gllel(eg)
qvol = 0.0
return
end
c-----------------------------------------------------------------------
subroutine userbc(i,j,k,f,eg) ! set up boundary conditions
c NOTE ::: This routine MAY NOT be called by every process
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
integer i,j,k,f,e,eg
e=gllel(eg)
ux = 0.0
uy = 0.0
uz = 1.0
flux = 0.0
temp = 0.0
if (cbc(f,e,1).eq.'W '.AND.cbc(f,e,2).eq.'t ') temp=1.0
return
end
c-----------------------------------------------------------------------
subroutine useric(i,j,k,eg) ! set up initial conditions
include 'SIZE'
include 'TOTAL'
include 'NEKUSE'
integer i,j,k,e,eg
ux = 0.0
uy = 0.0
uz = 1.0
temp = 0.0
return
end
c-----------------------------------------------------------------------
subroutine userchk
include 'SIZE'
include 'TOTAL'
c common /myoutflow/ d(lx1,ly1,lz1,lelt),m1(lx1*ly1*lz1,lelt)
c real m1
c rq = 2.
c uin = 0.
c call turb_outflow(d,m1,rq,uin)
n=lx1*ly1*lz1*nelt
if (istep.eq.0) then
c time=0.0
c
c Rfld = 0.9 ! radii of peb from restart file
c Rtar = 1.5 ! target radii matching Re/density
c scale = 1.5/0.9
c do i=1,n
c xm1(i,1,1,1)=xm1(i,1,1,1)*scale
c ym1(i,1,1,1)=ym1(i,1,1,1)*scale
c zm1(i,1,1,1)=zm1(i,1,1,1)*scale
c enddo
c call fix_geom
c call geom_reset(1)
xxmax = glmax(xm1,n)
yymax = glmax(ym1,n)
zzmax = glmax(zm1,n)
xxmin = glmin(xm1,n)
yymin = glmin(ym1,n)
zzmin = glmin(zm1,n)
if (nio.eq.0) then
write(6,18) xxmin,yymin,zzmin
write(6,19) xxmax,yymax,zzmax
18 format(' xyz min ',5g13.5)
19 format(' xyz max ',5g13.5)
endif
ifxyo=.true.
c call outpost(vx,vy,vz,t,pr,' ')
endif
c if (mod(istep,500).eq.0) then
umin=glmin(vx,n)
umax=glmax(vx,n)
vmin=glmin(vy,n)
vmax=glmax(vy,n)
wmin=glmin(vz,n)
wmax=glmax(vz,n)
tmin=glmin(t ,n)
tmax=glmax(t ,n)
if (nio.eq.0) write(6,1)
$ istep,time,umin,umax,vmin,vmax,wmin,wmax,tmin,tmax
1 format(i9,1p9e11.3,' tmax')
c
c call compute_cfl(cfl,vx,vy,vz,dt)
c call outpost(cflf,vy,vz,pr,t,'cfl')
c endif
return
end
c-----------------------------------------------------------------------
subroutine userqtl ! Set thermal divergence
call userqtl_scig
return
end
c-----------------------------------------------------------------------
subroutine usrdat ! This routine to modify element vertices
include 'SIZE'
include 'TOTAL'
return
end
c-----------------------------------------------------------------------
subroutine usrdat3
include 'SIZE'
include 'TOTAL'
return
end
c-----------------------------------------------------------------------
subroutine usrdat2
include 'SIZE'
include 'TOTAL'
integer e,f
c ifheat = .true.
nface = 2*ldim
do e=1,nelv
do f=1,nface
if(cbc(f,e,1).eq.'E '.or.cbc(f,e,1).eq.' ')then
cbc(f,e,1)='E '
cbc(f,e,2)='E '
elseif(cbc(f,e,1).eq.'v ')then ! z_bottom
cbc(f,e,2)='t '
elseif(cbc(f,e,1).eq.'O ')then ! z_top
cbc(f,e,2)='O '
elseif(cbc(f,e,1).eq.'SYM')then ! side of cylinder
cbc(f,e,1)='W '
cbc(f,e,2)='I '
elseif(cbc(f,e,1).eq.'W ')then ! pebbles
cbc(f,e,2)='t '
else
write(*,*)'bc incorrect',e,f,cbc(f,e,1)
endif
enddo
enddo
do iel=1,nelt
do ifc=1,2*ndim
if (cbc(ifc,iel,1) .eq. 'v ') boundaryID(ifc,iel) = 1
if (cbc(ifc,iel,1) .eq. 'O ') boundaryID(ifc,iel) = 2
if (cbc(ifc,iel,1) .eq. 'W ') boundaryID(ifc,iel) = 3
enddo
enddo
return
end
c-----------------------------------------------------------------------
#!/bin/bash
bin=nekrs
project_id=nfi114
project_id=csc262
time="00:30"
#source $NEKRS_INSTALL_DIR/bin/nekrs.bashrc
export OCCA_VERBOSE
case=$1
#if [ $# -ne 3 ]; then
# echo "usage: $0 <setup file> <number of nodes> <GPUs per node> <hh:mm>"
# echo "usage: $0 <setup file> <number of nodes> <GPUs per node>"
# exit 0
#fi
if [ $# -eq 4 ]; then
time=$4
# Bin Min Nodes Max Nodes Max Walltime (Hours) Aging Boost (Days)
# 1 2,765 4,608 24.0 15
# 2 922 2,764 24.0 10
# 3 92 921 12.0 0
# 4 46 91 6.0 0
# 5 1 45 2.0 0
fi
if [ ! -f $case.par ]; then
echo "Cannot find" $case.par
exit 1
fi
mkdir -p $OCCA_CACHE_DIR 2>/dev/null
nodes=$2
gpu_per_node=$3
let nn=$nodes*$gpu_per_node
jobname="NekRS_"$case
cmd="bsub -nnodes $nodes -W $time -P $project_id -J $jobname jsrun -n$nn -r$gpu_per_node -a1 -c1 -g1 -b packed:1 -d packed $bin --setup $case"
echo $cmd
echo $cmd > qq
$cmd
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment