Query reagrding Simulating flow over staggered cylinders.

[Madhav-46]

Hello Dr. Laizet and the Incompact3d Development Team,

I am trying to simulate flow over an array of 8 cylinders arranged in a staggered pattern. For this, I tried two approaches:

• Importing a .STL file using the Python script provided by Incompact3d for the sphere example case. However, the importing process and generation of initial conditions are taking a lot of time in my case.
• Modifying the cylinder example case to introduce additional cylinders. Unfortunately, I was unable to implement the required modifications for the case setup.

I would be grateful for any insights or documentation from the Incompact3d Development Team that could assist me in this regard.

Best regards,
Madhav.

[slaizet]

Hi Madhav,
For your setup for which there is an analytical expression for each cylinder, I would suggest to not use a STL file. I would modify the cylinder example case, starting with iibm=1 (no reconstruction), then I would move to iibm=2/3 but by adjusting the number of immersed objects allowed in each direction (nobjmax).
Best,
Sylvain

[Madhav-46]

Hello,
Thank you for your response. Could you please elaborate more on how and what files I shall modify for my case setup. Earlier, I attempted to modify the example case as follows:

! -- mode: f90 --

!===================
&BasicParam
!===================

! Flow type (1=Lock-exchange, 2=TGV, 3=Channel, 4=Periodic hill, 5=Cylinder, 6=dbg-schemes)
itype = 5

! Domain decomposition
p_row=0 ! Row partition
p_col=0 ! Column partition

! Mesh
nx=513 ! X-direction nodes
ny=100 ! Y-direction nodes
nz=65 ! Z-direction nodes
istret = 0 ! y mesh refinement (0:no, 1:center, 2:both sides, 3:bottom)
beta = 0.259065151 ! Refinement parameter (beta)

! Domain
xlx = 50 ! Lx (Size of the box in x-direction)
yly = 3.75 ! Ly (Size of the box in y-direction)
zlz = 3.0 ! Lz (Size of the box in z-direction)

! Boundary conditions
nclx1 = 2
nclxn = 2
ncly1 = 0
nclyn = 0
nclz1 = 2
nclzn = 2

! Flow parameters
iin = 1 ! Inflow conditions (1: classic, 2: turbinit)
re = 3900. ! nu=1/re (Kinematic Viscosity)
u1 = 1 ! u1 (max velocity) (for inflow condition)
u2 = 1 ! u2 (min velocity) (for inflow condition)
init_noise = 0.125 ! Turbulence intensity (1=100%) !! Initial condition
inflow_noise = 0.0 ! Turbulence intensity (1=100%) !! Inflow condition

! Time stepping
dt = 0.000005 ! Time step 4.79599e-05
ifirst = 1 ! First iteration
ilast = 400000 ! Last iteration

! Enable modelling tools
ilesmod=1 ! if 0 then DNS
iscalar=0 ! If iscalar=0 (no scalar), if iscalar=1 (scalar)
iibm=2 ! Flag for immersed boundary method

/End

!====================
&NumOptions
!====================

! Spatial derivatives
ifirstder = 4 ! (1->2nd central, 2->4th central, 3->4th compact, 4-> 6th compact)
isecondder = 5 ! (1->2nd central, 2->4th central, 3->4th compact, 4-> 6th compact, 5->hyperviscous 6th)
ipinter = 3 ! interpolation scheme (1: classic, 2: optimized, 3: optimized agressive)

! Time scheme
itimescheme = 3 ! Time integration scheme (1->Euler,2->AB2, 3->AB3, 4->AB4,5->RK3,6->RK4)

/End

!=================
&InOutParam
!=================

! Basic I/O
irestart = 0 ! Read initial flow field ?
icheckpoint = 2000 ! Frequency for writing backup file
ioutput = 2000 ! Frequency for visualization
ilist = 25 ! Frequency for the output to screen
nvisu = 1 ! Size for visualisation collection

/End

!=================
&Statistics
!=================

nstat = 1 ! Size arrays for statistic collection
initstat = 10000 ! Time steps after which statistics are collected

/End

!########################
! OPTIONAL PARAMETERS
!#######################

!================
&ScalarParam
!================

numscalar = 0 ! How many scalars?
Sc = 0.2 ! Schmidt number

/End

!================
&LESModel
!================

jles = 1 ! LES Model (1: Phys Smag, 2: Phys WALE, 3: Phys dyn. Smag, 4: iSVV, 5: dyn SEV)
smagcst = 0.14 ! Smagorinsky constant
SmagWallDamp = 1 ! 1: Mason and Thomson Damping function, otherwise OFF
walecst = 0.5 ! WALES Model Coefficient
iconserv = 0 ! Formulation SGS divergence (0: non conservative, 1: conservative)

/End

!================
&ibmstuff
!================

cex(1)=17.5 !when simulating a cylinder x coordinate of the center
cey(1)=0.625 !when simulating a cylinder y coordinate of the center
ra(1)=0.5 !when simulating a cylinder, radius

cex(2)=17.5 !when simulating a cylinder x coordinate of the center
cey(2)=3.125 !when simulating a cylinder y coordinate of the center
ra(2)=0.5 !when simulating a cylinder, radius

cex(3)=20
cey(3)=1.875
ra(3)=0.5

cex(4)=22.5
cey(4)=0.625
ra(4)=0.5

cex(5)=22.5
cey(5)=3.125
ra(5)=0.5

cex(6)=25
cey(6)=1.875
ra(6)=0.5

cex(7)=27.5
cey(7)=0.625
ra(7)=0.5

cex(8)=27.5
cey(8)=3.125
ra(8)=0.5

nraf=10 !level of refinement for iibm==2 to find the surface of the immersed object
nobjmax=8 !number of immersed objects (DO NOT USE ZERO OBJECTS)
iforces=1
nvol=8 !Number of volumes for computing force balance

/End

&ForceCVs

xld(1) = 16.5
xrd(1) = 18.5
yld(1) = 0.1
yud(1) = 1.15

xld(2) = 16.5
xrd(2) = 18.5
yld(2) = 2.6
yud(2) = 3.65

xld(3) = 19.0
xrd(3) = 21.0
yld(3) = 1.35
yud(3) = 2.4

xld(4) = 21.5
xrd(4) = 23.5
yld(4) = 0.1
yud(4) = 1.15

xld(5) = 21.5
xrd(5) = 23.5
yld(5) = 2.6
yud(5) = 3.65

xld(6) = 24.0
xrd(6) = 26.0
yld(6) = 1.35
yud(6) = 2.4

xld(7) = 26.5
xrd(7) = 28.5
yld(7) = 0.1
yud(7) = 1.15

xld(8) = 26.5
xrd(8) = 28.5
yld(8) = 2.6
yud(8) = 3.65

/End

&CASE
/End

Do I need to modify the Case-Cylinder-wake.f90 file or create a new generic .f90 file in the src folder? I was unable to implement the required modifications in the Case-Cylinder-wake.f90 file. Could you guide me on how to proceed further? Also, after modifying or creating the .f90 file, what steps are required to compile this new solver?

I am attaching an image of the domain for your reference.

Best regards,
Madhav

[slaizet]

Hi,
You made a lot of modifications like the BC in the Z direction, the Reynolds number which is too high for that resolution, the LES mode that is not set up properly, etc. My recommendation would be to go step by step, starting from input_DNS300_LR.i3d. You should be able to add several cylinders by just hard coding new cylinder in subroutine geomcomplex_cyl and using iibm=1 as suggested in my previous email. Then you can add the flow reconstruction by defining nobjmax as the maximum number of fluid/solid and/or solid/fluid interfaces.

Below is an example with two cylinders, obtained with minimal changes.

I just added the loop below in subroutine geomcomplex_cyl

do k=nzi,nzf
       do j=nyi,nyf
          ym=yp(j)
          do i=nxi,nxf
             xm=real(i-1,mytype)*dx
             r=sqrt((xm-cexx+2.)**two+(ym-ceyy+2.)**two)
             if (r-ra.gt.zeromach) then
                cycle
             endif
             epsi(i,j,k)=remp
          enddo
       enddo
    enddo

and nobjmax=4 in the input file.

Best,
Sylvain