Viscous Flow with Julia

In this lab, you will use the Julia package ViscousFlow.jl to conduct numerical simulations of fluid flows at moderate Reynolds numbers.

You will submit your work as Jupyter Notebooks as well as code files. However, don't worry about how to use Jupyter Notebooks for now -- once you've completed the coding tasks, use the tutorial here (https://www.jousefmurad.com/coding/install-julia-jupyter-notebook/) to learn how to set up a Jupyter notebook.

Setting up Julia: environments and the ViscousFlow package.

Download and install Julia, using the instructions specific to your operating system. https://julialang.org/

  • Run the program 'Julia', which should bring up the 'REPL' interface ('Read - Evaluate - Print - Loop'), which is essentially a terminal instance that is running julia. On Unix-like systems, you can do this by entering julia on the command line.

  • Create a new folder named lab4 or similar.

  • Inside julia, type cd("path to lab4"), which should change the present working directory to that folder.

  • At any time, you can check Julia's present working directory by running the function pwd() with no arguments.

  • Activate an environment at the folder you created. Type the character ]. This brings you to the 'package manager' of Julia, and the prompt message changes. Then use the command activate . where "." is shorthand for "the present working directory". You can also activate an environment located in a folder that's not your present working directory, but let's avoid that for now.

  • Type the character ]. This brings you to the 'package manager' of Julia, and the prompt message changes. Your prompt message should now show the name of your environment in parentheses. Here, you can add packages from the extensive library of 'official' packages located at https://github.com/JuliaRegistries/General. You can also add custom packages, your own packages, packages downloaded from other websites, etc.

  • Use add ViscousFlow and add Plots to add these two packages to your current environment. You only need to do this once.

  • The following notation is an alternative way to add packages.

Pkg.add("ViscousFlow")
Pkg.add("Plots")
1.1s
  • To check what packages are installed in your currently-activated environment, enter status from inside the package manager (i.e., after you have pressed ] and changed your prompt message.

  • If you close Julia now, reopen it, and activate this package again, it will still contain these packages. So you only need to add packages once inside any particular environment.

  • Inside Julia, type using ViscousFlow. Check that you indeed have this package installed by doing something that will only work if it's installed correctly: Enter ? to get into help mode, type setup_grid and press enter. This should show you the help message for the built-in function setup_grid. Verify that tab-completion works as well.

Note

There are two alternative ways of running code in Julia.

  • One way is to use the REPL to execute code line-by-line.

  • The other is to write your code in a separate file; the convention is to name the file, say, script1.jl, and then from within the Julia instance, run include("script1.jl").

A simple script making use of the ViscousFlow package for flow past cylinder at Re = 40

Start by typing out the following code line-by-line into the REPL, making sure you understand what each line is doing. Later, you should copy-paste the text blocks into your own script file and run using the include command. Turn in this completed julia file as "Cylinder_Re40.jl"

using ViscousFlow
using Plots
# initialize a 'dictionary' that will be used for parameters such as the Reynolds number
simulation_parameters = Dict() 
simulation_parameters["Re"] = 40
simulation_parameters["grid Re"] = 1.5
simulation_parameters["freestream speed"] = 1.0
simulation_parameters["freestream angle"] = 0
# Create tuples that contain the x and y limits of your domain.
xdomain = (-2.0,10.0)
ydomain = (-3.0,3.0)
# Using the constructor function setup_grid, create an object of type 'PhysicalGrid' and assign it to the variable "g".
grid1 = setup_grid(xdomain,ydomain,simulation_parameters)

Task: What information is contained in grid?

Using the help command in Julia (type ? followed by the name of the item you want to know about), explore the object that you have just created, grid. What information does it contain? You should start by simply entering grid as a command in Julia. This will tell you what 'type' this object is. Then, get help on this 'type'.

The object grid is of type <---->. It contains the following four attributes:

  1. attribute 1

  2. attribute 2

  3. attribute 3

  4. attribute 4

  5. attribute 5

Note that if an object named A is of type thistype, and objects of type thistype have attributes including B and C, then you can access the attributes of object A using the notation A.B and A.C

Setting up a solid body immersed in the fluid

For now, we will use a cylinder of radius 1.

Δs = surface_point_spacing(grid1,simulation_parameters)
body1 = Circle(0.5,Δs)
location = [1.0,0.0] 
α = 30π/180
X = MotionTransform(location,α)
joint = Joint(X)

Task: Explain the preceding lines of code

Explain what each line of code does here in the form of comments. You may extensively borrow information from the help pages for this task. You are encouraged to use ? <function name> etc to find these help pages.

Set up 'Rigid Body Motion' (but with no motion for now)

The following code is some package-specific machinery that is not important to understand at this stage, but must still be included and run together with

motion = RigidBodyMotion(joint,body1)
x = init_motion_state(body1,motion)
update_body!(body1,x,motion)

Task: Plot a picture of the body in the domain

Use code similar to the following

plot(body1,xlim=xdomain,ylim=ydomain)
savefig("domain.png")

Set up the simulation

sys = viscousflow_system(grid1,body1,phys_params=simulation_parameters,motions=motion)
u0 = init_sol(sys)
tspan = (0.0,20.0)
integrator = init(u0,tspan,sys)

The object 'integrator' now needs to be 'stepped forward' in time. The program will find the best Δt to use for this; you just need to specify how far you would like to advance the problem.

step!(integrator,1.0)

Task: Plot the velocity, streamfunction, and vorticity after one time unit has passed.

The ViscousFlow package can be used together with Julia's plotting package (https://docs.juliaplots.org/stable/). Extract the velocity, vorticity and streamfunction from integrator using code such as the following:

uv1 = velocity(integrator)
plot(uv1,sys,savefig="velocity1.png")
ω1  = vorticity(integrator)
ψ1  = streamfunction(integrator)
p1  = pressure(integrator)

and then plot these quantities. You will produce four images; save each as a separate file with an appropriate naming convention.

Task: Write code that advances the solution forward until time t = 10

and make a plot of the streamlines and the vorticity contours after 10 time units have passed.

Task: Determine how long before no significant change happens

Continue to step forward the solution in time until you determine that no further change happens. Record the time at which "steady state" is reached, if such a state is possible.

Flow past a square

Now that you have conducted a simulation for the flow past a circle, conduct a study of the flow past a square at the same Reynolds number. You will need to make use of the function Rectangle instead of Circle.

Merger of two vortices

Create a numerical simulation of the process of two Gaussian vortices located at the locations (+2,+2) and (-2,-2). Each vortex should be a Gaussian with radius 0.25 and unit strength. You should refer to the tutorial/example notebook at https://github.com/JuliaIBPM/ViscousFlow.jl/blob/master/examples/1.-Basic-viscous-flow.ipynb to get started.

You will have to decide how long to run the simulation. Your goal is march it forward in time until the two vortices have, for all intents and purposes, merged with each other.

Generate three plots of the vorticity, one showing the initial state, one showing an intermediate state, and one showing the final state. Also generate three plots of the streamfunction at each of these states.

Note that your domain should be large enough that all the interesting things happen inside the domain.

Turn in the code for this task as "vortex_merger.jl".

Merger of Elliptical Vortices

Repeat the above task but for two elliptical vortices located a similar distance apart. Document your findings in a Jupyter Notebook.

Summary of deliverables.

There are four high-level tasks:

  • Flow past circular cross-section

  • Flow past square cross-section

  • Merger of circular vortices

  • Merger of elliptical vortices

Submit two Jupyter notebooks (one for flow with a body and one for flow without a body) and four .jl files, one for each high-level task.

When I run your .jl file, I should be able to reproduce your figures.

Runtimes (1)