MECH4001辅导、Java，c++程序辅导

日期：2021-12-17 07:47

Project 1
Unit Number: MECH4001
Unit Name: Computational Fluid Dynamics
Please read the first page of this document carefully before you doing your project.
Instruction
? Students will develop CFD programs to solve simple problems in this project.
? There are two parts in this project, and each part worth 50% of the marks of this project. Each
student needs to complete the work assigned in each part and submit a written report.
? Each student needs to submit the following TWO documents through vUWS:
(1) The report without the code including both Part 1 and Part 2 into assignment box Project
1 report submission on vUWS.
(2) Copy the program code of both parts I and II into a word document and submit it to
assignment box Project 1 code submission on vUWS.
Note: submitting only one document will lead to a zero mark.
? The reports must be in the format of Microsoft Word file and submitted through vUWS.
Reports in hard copy or other file types will not be accepted.
Note:
Writing
The report will be assessed based on the accuracy of the results, discussion and writing. Parameters in equations need to
be defined and clearly explained if necessary. Figures should have clear captions and legends. It is not acceptable to only
list all the formulae and equations without any explanation. You need to type formulae in your report. Copying formulae
as pictures form books and lecture notes are not acceptable. You can use Excel to draw diagrams or tables based on the
calculated results.
Computer program code
You only submit one set of computer program code in each part. If you need to do many cases of simulation, in your
report, explain clearly what parameters in the code you need to change or modify for doing different cases, and highlight
the parameters in the submitted computer code.
Plagiarism
Plagiarism is unacceptable. If plagiarism is suspected to be involved in a report, it will be reported to the school for
further investigation.
Late submission
- Submitting the report late due to your own internet connection problem or computer problem will be treated as late
submission. To avoid delay, you need to submit the report before the due date.
Page?2?of?6
Part one: Simple initial value problem
Assignment (please read the problem description before doing the assignment)
(1) Find out the equation of motion of the cylinder. The derivation procedure of the equation
must be included in the report.
(2) Develop a numerical method for predicting the vibration of the cylinder. The procedure of
deriving the FDM formula must be included in the report.
(3) Develop a MATLAB program for solving the problem. The MATLAB program code must be
submitted (see the instruction on the first page of this document).
(4) Keeping the damping coefficient of c = 100 N?m/s, conduct numerical simulations for KC=2
to 20 with an interval of 2 and find out the variation of the power with the KC number.
(5) Keeping the KC number to be constant of KC=10, conduct numerical simulations for c=0 to
2000 N?m/s with an interval of 50 and find the variation of the power with c.
(6) In the report, show the time histories of the vibration for all the calculated cases. Include
them in the report as Appendix.
(7) Discuss how KC and c affects the power generation.
Problem description
An offshore wave energy extraction device includes a cylinder elastically mounted on a spring
and the electricity generator. The cylinder vibrates only in a direction parallel to the flow and drives
an electricity generator. The electricity generator extracts the energy from the cylinder in the same
way as a damper. It is modelled as a damper with a damping constant c when the vibration is studied.
The flow velocity is an oscillatory flow with two frequencies
where Um is the amplitude of the oscillatory flow velocity, ω=2π/T is the angular frequency, T is the
period of the oscillatory flow and t is the time.
The KC number is defined as
D
U T KC m ?
The electricity generator can modelled by a damper that extract energy from the motion of the
cylinder with a constant damping constant of c. The structural damping is considered to be negligibly
small. The instantaneous power extracted from the motion by the electricity generator (damper) can
be calculated by
?? ? ????
k
c
Oscillatory?flow
Page?3?of?6
where P is the instantaneous power of the system, F is the hydrodynamic force on the cylinder, V is
the vibration speed of the cylinder. The time averaged power can be calculated by
??? ? ? ????d?? ????
?
The power must be calculated by stable numerical solution. The hydrodynamic force on the cylinder
can be predicted by the Morison equation
inertia drag
where Vr=u-V is the velocity of the flow relative to the cylinder, CA and CD and inertia and drag
coefficients, respectively, md is the displaced fluid mass and Ap is the projected area. The calculating
parameters are listed in the table below.
Table 1 Parameter table
Note: In the table, the number j stands for the last digit of your student ID. For example, if your
student ID is 2345673, 10+j = 13.
Amplitude of the oscillatory flow velocity, Um
10
1 j ? m/s
Steady flow velocity U0 0.1×Um m/s
Diameter of the cylinder, D 10+j cm
Length of the cylinder, L 1 m
Mass of the cylinder, m 50 kg
Density of the water, ρ 1024 kg/m3
Stiffness of the spring (total), K 200 N/m
Damping coefficient (total), c See the assignment
Added mass coefficient (assume to be constant) 1
Drag coefficient (assume to be constant) 1.8
Note: Stable, periodic vibration can be achieved only if the simulated time is sufficiently long.
? ?
Page?4?of?6
Part 2: Hydrodynamic force on a square prism
The figure below shows a water flow past a square/rectangular prism with a boundary length of Bx
and By in the x- and y-directions, respectively. The flow will be simulated by solving the NavierStokes
Equations. Water density is 1000 kg/m3
and viscosity is 1.3×10-6 m2
/s. The Navier-Stokes
equations and the numerical method for solving the Navier-Stokes equations have been explained in
Lecture 7. The drag force and lift force on the prism will be analysed.
The drag force FD has two components: the drag force from the pressure FDP and the drag from the
shear stress FDS. Similarly, the lift force FL also has two components: FLP and FLS.
An example MATLAB program has been given in Lecture 7 (Chapter 7).
Assignment
The velocity of the incoming flow is U=(10+j)/100 m/s, where j is the last digit of your student
number. For example if your student ID is 12345, the velocity U=(10+5)/100=0.15 m/s.
1. Clearly specify the boundary conditions of the problem.
2. You can choose ?x=?y=0.0001 m. Using CFL condition, find out the maximum computational
time step for stability. Discuss how accuracy this maximum time step is. Choose a time step
satisfy the CFL condition and workout the Courant number of your case.
3. Find out the time histories of the drag and lift force for square prisms with Bx=By=0.01 m, 0.012
m, 0.014 m, 0.016 m and 0.02 m. You need to show the plot of the time histories of the forces for
all the cases in the report.
4. Calculate the mean drag coefficient (averaged over time), amplitude of the drag coefficient
amplitude of the lift coefficient of the prism for all the cases. Discuss how the Reynolds number
affects the force coefficients (draw coefficient vs Re diagrams). The drag and lift coefficient is
defined as 2
5. Calculate the Strouhal number based on the numerical results and discuss how the Reynolds
number affects the Strouhal number. The Strouhal number is St=fB/U, where f is the vibration
frequency of the lift force and U is the incoming flow velocity.
6. Show the velocity vectors for one of the cases you have calculated and make comments on the
vortex shedding flow.
7. Repeat steps 2 to 4 using Rectangular prism 1 with (Bx,By)=(0.004 m, 0.01 m) and Rectangular
prism 2 with (Bx,By)=(0.04 m, 0.01 m). Discuss the effect of the prism size in the x-direction on
the force coefficients and the Strouhal number. To facilitate the comparison, you can draw
variation of drag coefficient with Re for all the three cases (square prism, rectangular prism 1 and
2) in the same diagram, and do the same for lift coefficient and Strouhal number.
Note: In this part, your ability of developing a method for calculating the forces in the code is
examined.

Page?6?of?6
A method for calculating forces:
Shear stress on the wall is ? ? μ
??
?? on the top and bottom surface and ? ? ??
??
?? on the left and right
boundary, where μ is the dynamic viscosity (=?ν). The method for calculating the shear stress
numerically is given in the figure below (the sign (direction) of the shear stress must be correct).
On each cell on the four wall boundaries, the forces Fxe and Fye in the x- and y-directions can be
calculated (see the picture below). The directions (signs) of the forces must be correct when
calculating Fex and Fey.
The force on the prism in the x- and y-directions are
The part of the program code for calculating Fx and Fy automatically needs to be developed.
A schematic figure demonstrating a numerical method for calculating forces of a rectangular prism
?x
?y
Shear stress ?
Pressure p
Fxe=0.5(pi1,j+pi1,j+1)?y
i1
j
j+1
Fye=0.5(?i1,j+?i1,j+1)?y
j1
i i+1
Fye=0.5(pi,j1+pi+1,j1)?x Fxe=0.5(?i,j1+?i+1,j1)?x
?=μu/?y, u is the velocity next to the wall
u
?=μv/?x,
v is the velocity
next to the wall
v
Wall
Boundaries