-- Musubi configuration file. 
-- This is a LUA script.
require 'seeder'

if qVal then
  bc_cyl = 'wall_linearInterpolation'
else
  bc_cyl = 'wall'
end 
bc_wall = 'wall'

--flow parameters
Re = 20
rho0_p = 1.0
rho0_L = 1.0
cs2 = 1.0/3.0
nu_phy = 1e-3 --m^2/s
u_mean_phy = Re*nu_phy/Dia --s

u_max_phy = 3.0*u_mean_phy/2.0

scaling = 'diffusive'
scaling = 'acoustic'
if scaling == 'acoustic' then
  u_max_L = 0.025
  dt = u_max_L * dx / u_max_phy
  nu_L = nu_phy * dt / dx^2.
  omega = 1.0 / (3.0 * nu_L + 0.5 )
elseif scaling == 'diffusive' then
  omega = 1.85
  nu_L = (1.0/omega-0.5)/3.0
  dt = nu_L*dx^2/nu_phy
  u_max_L = u_max_phy*dt/dx
else
  print 'ERROR: scaling not supported'
  os.exit()
end

u_mean_L = 2.0*u_max_L/3.0
cs = math.sqrt(1.0/3.0)
machNr = u_max_L / cs
if machNr > 0.1 then
  print ('ERROR: Mach number is above LBM limit. Ma =',machNr)
  os.exit()
end 

ref_press = 1.0
press_drop = 8*u_max_phy*rho0_p*nu_phy*length/height^2
press_grad = press_drop/length

function u_inflow(x,y,z) 
  --return (0.5/(rho0_p*nu_phy))*press_grad*y*(height-y)
  return 4.0*u_max_phy*y*(height-y)/height^2.0
end

-- draf and lift coeff factors
cD = 2 / (rho0_p * u_mean_phy * u_mean_phy * Dia * dx)
cL = cD

-- Simulation name
simulation_name = 'channel'
mesh = 'mesh/' -- Mesh information
printRuntimeInfo = false

interpolation_method = 'linear'

-- Time step settigs
if Re == 20 then
  tmax_p = 20 --s
else
  tmax_p = 20 --s
end
tmax  =  tmax_p/dt    -- total iteration number
interval_p = 1 -- s
sim_control = {
  time_control = {
    max = {sim=tmax_p},
    interval = interval_p,
  },
  abort_criteria = {
    stop_file = 'stop',
    steady_state = false
  }
}  

--debug = {debugMode = true, debugFiles=true, debugMesh='debug/'}
-- logging
logging = {level=5}

-- needed to dump variable in physical unit
physics = { dt = dt, rho0 = rho0_p }            

fluid = { 
         omega = omega, rho0 = rho0_p }

-- Initial condition 
initial_condition = { pressure = ref_press, 
                      velocityX = 0.0,
                      velocityY = 0.0,
                      velocityZ = 0.0, }

identify = {label='2D',layout='d3q19', relaxation='bgk',
            kind = 'lbm_incomp'
            }

-- Boundary conditions
boundary_condition = {  
{ label = 'west', 
  kind = 'inlet_ubb',
 -- kind = 'wall',
 -- kind = 'outlet_expol',
 -- pressure=ref_press,
  velocityX = { predefined = 'combined',
                temporal = {predefined='smooth', min_factor = 0, 
                            max_factor=1.0, from_time=0, to_time=1}, 
                spatial = u_inflow
	      },
   velocityY = 0.0, velocityZ = 0.0},
{ label = 'east',
  --kind = 'wall',
   kind = 'outlet_expol',
   pressure = ref_press}, 
{ label = 'north', 
  kind = bc_wall ,
   fac = 1.0 },
{ label = 'south', 
  kind = bc_wall},
{ label = 'obst', 
  kind = bc_cyl},
 }

-- lua variables
variable = {
  { name = 'coeff_fac', 
    ncomponents = 3, 
    vartype = 'st_fun',
    st_fun = {cD, cL, 0} 
  },
  { name  = 'coeff',
    ncomponents = 3, 
    vartype = 'operation',
    operation = { 
      kind = 'multiplication',
      input_varname = { 'bnd_force_phy', 'coeff_fac'}
    }  
  }
}

-- Tracking              
tracking = {
{
  label = 'probe', 
  folder = 'tracking/',
  variable = {'pressure_phy','velocity_phy'}, 
  shape = {kind = 'canoND', object = {origin ={length*0.5,height*0.5,zpos} } },
  time_control = {min = 0, max = tmax_p, interval = {iter=10}},
  format = 'ascii'      
},
-- pressure drop at two points over time: P(0.15,0.2,t) and P(0.25,0.2,t)
{
  label = 'probe_pressure_1',
  folder = 'tracking/',
  variable = {'pressure_phy'}, 
  shape = {
          {kind = 'canoND', object = {origin ={offset_cyl+0.15-dx,0.2,zpos} }},
          },
  time_control = {min = 0, max = tmax_p, interval = 10*dt},
  format = 'ascii'      -- over time
},
{
  label = 'probe_pressure_2',
  folder = 'tracking/',
  variable = {'pressure_phy'}, 
  shape = {
          {kind = 'canoND', object = {origin ={offset_cyl+0.25+dx,0.2,zpos} }}
          },
  time_control = {min = 0, max = tmax_p, interval = 10*dt},
  format = 'ascii'      -- over time
},
{
  label = 'centerLine', 
  folder = 'tracking/',
  variable = {'pressure_phy','velocity_phy'}, 
  shape = {kind = 'canoND', object = {origin ={0.0,0.2,zpos},
                                     vec = {length,0.0,0.0},
                                     segments = nLength+2} },
  time_control = { min = tmax_p, max = tmax_p, interval = tmax_p},
  format = 'asciiSpatial'      
},
{ label = 'convergence',
  variable = {'pressure_phy'}, 
  shape = {
          {kind = 'canoND', object = {origin ={offset_cyl+0.15-dx,0.2,zpos} }},
          {kind = 'canoND', object = {origin ={offset_cyl+0.25+dx,0.2,zpos} }}
          },
  time_control = {min = 1, max = tmax_p, interval = {iter=10}},
  format='convergence',
  reduction = {'average'},
  convergence = {norm='average', nvals = 50, absolute = true,
  condition = { threshold = 1.e-10, operator = '<=' }}
},
{ label = 'coeff',
  folder = 'tracking/',
  variable = {'coeff'}, 
  shape = {kind = 'boundary', boundary = {'obst'}},
  time_control = {min = 1, max = tmax_p, interval = {iter=100}},
  format='ascii',
  reduction = {'sum'}
},
{ label = 'hvs',
  folder = 'tracking/',
  variable = {'pressure_phy','velocity_phy','process','shear_stress'}, 
  shape = {kind = 'all'},
  time_control = {min = tmax_p, max = tmax_p, interval = interval_p},
  format='harvester'
}

}

estart = {
      ead = 'restart/'..'/channel2D_lastHeader.lua',
      write = 'restart/',
      time_control = { min = tmax_p, max = tmax_p, interval = tmax_p}
 }
