\(\require{physics}\)

ポアズイユ流れ

package.cpath = "../fantom_libraries/Release/?.dll;" .. package.cpath
require("fantom")
require("writer")
-- basic parameters
case_name              = "plane_poiseuille"
total_count            = 1000
save_interval          = 10
dx                     = 1.0
dy                     = 1.0
dt                     = 1000.0
fantom.Parameter_THETA = 0.5
num_threads            = 1
-- create vertical grid info
dz    = fantom.VectorDouble(50)
for i = 0, dz:size()-1 do
  dz[i]=0.02 -- meter
end
-- creating geometry and refinement files (one vertical column)
fg  = io.open("geom.txt", "w")
fxr = io.open("xrfm.txt", "w")
fyr = io.open("yrfm.txt", "w")
fg:write ("0.0\n")
fxr:write("1\n")
fyr:write("1\n")
fg:close()
fxr:close()
fyr:close()
-- creating objects
water_body   = fantom.WaterBody(dx, dy, dz, "geom.txt", "xrfm.txt", "yrfm.txt", false, true) -- periodic in west-east
diffusion    = fantom.VerticalDiffusionHorizontalVelocityNonslip()
txt_writer   = writer.TextWriter(water_body, case_name)
-- object settings
water_body:set_num_threads(num_threads) 
water_body:set_var("u",  0.0)
water_body:set_var("v",  0.0)
water_body:set_var("w",  0.0)
-- forcing function
function pressure_gradient(water_body)
  dpdx = 1.0
  rho  = 1000.0
  for i=0,water_body:vertical_face_columns():size()-1 do
    local column = water_body:vertical_face_columns()[i]
    if column:direction()==fantom.Parameter.WEST_EAST then
      for k=0, column:original_components():size()-1 do
        column:get_face(k):velocity():add_increment(1.0/rho*dpdx) -- 1/rho*dp/dx
      end
    end
  end
end

-- time marching
for counter = 0, total_count do
  water_body:construct_box_variables()  
  print("----- count = ", counter)
  if counter % save_interval == 0 then
    txt_writer:write(dt*counter, counter)
  end
  pressure_gradient(water_body)
  water_body :update_explicit_terms(dt)
  diffusion:apply(water_body, dt)
end