matlab question

Drexel University, College of Engineering 2015-2016 Academic Year

1

Drexel University Office of the Dean of the College of Engineering

ENGR 232 – Dynamic Engineering Systems

Week 1 Laboratory Assignment

In this lab we will investigate the use of direction fields to evaluate solutions to first-order differential equations and introduce many important Matlab commands for working exactly with integrals and differential equations.

For this lab, we will study these qualitative and exact analysis tools via the example of population growth using the Gompertz equation as a model. This equation is used to model the growth of tumor cells, and gives an expression of the dependence of tumor size on time.

The differential equation governing the growth is:

𝑑𝑉 𝑑𝑡

= 𝛽𝑉 ln 𝑘 𝑉

Here 𝛽 is the growth rate and 𝑘 is a reflection of the carrying capacity of the tumor. It’s growth is limited by the supply of nutrients.

Let 𝛽 = 2; 𝑘 = 100; and V range from 0 to 140

Part A. Open Matlab, and review basic techniques like symbolic integration and solution of differential equations using dsolve.

It is always a good idea to initialize your workspace when you start a new problem. This especially helps clear out any variables that may have been used in previous problems. Just to review how the clear command works, enter the following two values for the growth rate and carrying capacity of the tumor.

B = 2; K = 100;

If you type in either B or K, Matlab returns the value you just entered. OK, let’s clear these values. Enter clear, then see if Matlab still knows B and K. Next type in clc. What does it do? Explore the help facility by entering: help clear and help clc.

Solving differential equations usually involves finding integrals. Matlab’s int command can help! First, get help on the int command. help sym/int

Note: You must always declare any symbolic variables before attempting a symbolic integration.

Here is how to find for example, the integral of x2.

syms x int(x^2, x)

Question 1: Find this integral and display the answer in the box below. Find: ,- .

/0-1 𝑑𝑥

Note that by default, Matlab does not display the constant of integration c. You should include that in your answer.

Also include the code you used in the answer box below.

Alert: If you don’t declare x to be a symbolic value you will get an error.

1. The integral is:

Drexel University, College of Engineering 2015-2016 Academic Year

2

Matlab can also solve many differential equations exactly using the dsolve command. First clear your variables and command window, then get help on the dsolve command. Spend a few minutes reviewing that help info. There is lots of good stuff!

Try this simple example to solve the second-order differential equation: 31- 341

= −4𝑥

You may recognize this as simple harmonic motion and should expect two constants of integration.

syms x(t) eqn = diff(x,t,t) == -4*x dsolve(eqn) Now find the solution with the initial conditions x(0) = 1, dx/dt(0) = 0. Dx = diff(x) dsolve(eqn, x(0)==1, Dx(0)==0) Question 2: Find the exact solution to this new differential equation with the given initial conditions. Place your code and answers in the numbered box below. You can save space by placing multiple commands on the same line separated by a ;

31- 341

= −9𝑥 + 𝑡 with the initial conditions x(0) = 1, dx/dt(0) = 2.

Question 3: The Gompertz differential equation at the top of the first page is a very difficult non-linear differential equation. But let’s see if Matlab can solve it exactly. No harm in at least trying. Clear all your variables. Now enter the values for the growth rate and carrying capacity of the tumor once again.

B = 2; K = 100; The volume V of the tumor is a function of time – it’s growing!

syms V(t); Complete the next line to enter the Gompertz DE into Matlab. You can use B and K in your code. Matlab uses log for the natural log. Don’t forget to use == to denote equality in the equation.

𝑑𝑉 𝑑𝑡

= 𝛽𝑉 ln 𝑘 𝑉

eqn(t) = diff(V(t), t) … ß add more to complete this line. Now solve the Gompertz exactly using dsolve and record your command and the specific solution starting with a volume of 1 at time 0. Record your code and the answer in the box below.

Note the use of == so you can give Matlab actual equations!

2.

3.

Drexel University, College of Engineering 2015-2016 Academic Year

3

Once you get your answer, you can confirm it approaches a limiting volume of 100 using: eval( subs(ans, t, 100) ) Part B: Qualitative Analysis The Gompertz differential equation is notoriously difficult to solve, but Matlab quickly dispatched with the challenge. That’s really impressive raw power! Now let’s explore the solutions graphically to get a better intuition about the Gompertz model. We will now use the direction field tool to examine the solution. Download the Java JAR named dfield available in the Software folder. Use only the Java version for this lab.

Plot the direction field for the Gompertz differential equation using the Java program named dfield. Set time from 0 to 10, and V from 0 to 140. Let 𝛽 = 2 and 𝑘 = 100 as before.

𝑑𝑉 𝑑𝑡

= 𝛽𝑉 ln 𝑘 𝑉

The graph below shows two solutions starting at V(0) = 50 and V(0) = 150 obtained using the keyboard input option available under the Solution menu of the JAR. The curve approaching from above may correspond to a large tumor, which has lost some of its former blood supply due to necrosis and is now shrinking.

Drexel University, College of Engineering 2015-2016 Academic Year

4

Question 4: Add 7 more solutions so that the volume at times 1, 2, 3, 4, 5, 6 and 7 is 50. These tumors are growing. Add solution curves so that the volume at times 1, 2, 3, 4, 5, 6 and 7 is 150. These tumors are shrinking. Note some of these curves may need to be recalculated using the “Recalculate All Solutions Using Current Settings” under the Solution menu. Paste your image here.

Solution curves for Question 4 with B = 2 Question 5: Estimate how long it takes a tumor of size V = 20 to double in size to V = 40. Zoom in! Use B = 2.

Question 6: The critical values for this equation are V = K = 100 and V = 0. Classify each as stable or unstable in the box below.

(Don’t worry about semi-stable. The DE is not even defined for negative volumes.) Question 7: Simply change B (beta) from 2 to 4. Obtain a new direction field, and include multiple solutions curves passing through V = 50 and 150 as before. Paste your plot below. Solution curves for Question 7 with B = 4

5. With B = 2, it takes a tumor starting at V = 20 about this long to double to V = 40: ANS: _ _ _ _ time units

The critical value V = K = 100 is ... The critical value V = 0 is …

Drexel University, College of Engineering 2015-2016 Academic Year

5

Question 8: Estimate how long it takes a tumor of size V = 20 to double in size to V = 40. Zoom in! Now B = 4. These tumors are growing much faster!

Question 9: Fastest Growth The tumors are growing the fastest at the maximum for:

𝑑𝑉 𝑑𝑡

= 𝛽𝑉 ln 𝑘 𝑉

Let’s find the maximum growth rate using Matlab’s diff command. Start by clearing all variables in the Matlab window. Then declare B, K and V(t) to be symbols. syms V(t) B K Enter the Gompertz equation again, but this time B and K are also symbolic. eqn1 = diff(V,t) ==B*V*log(K/V); To find the maximum, we need to set the double derivative to zero. (inflection point). eqn2 = diff(eqn1, t) eqn2(t) = diff(V(t), t, t) == B*log(K/V(t))*diff(V(t), t) - B*diff(V(t), t) Simplifying by hand we see the inflection point is where: log(K/V(t)) = 1 These tumors are growing fastest at the inflection point: (answer is a simple fraction) Place your answer in the box below. Question 10: Autonomous Property The Gompertz equation does not depend on time explicitly. Such equations are called autonomous. All the solution curves for growing tumors (below the value K) can be obtained from just one by shifting it in the … direction. Enter horizontal or vertical below. Be sure all ten questions are answered, then submit your Lab1 file as a single PDF. We only need to see your name and the answers labelled 1 – 10, You can delete the explanatory material when you submit if you wish. Submission must be a single PDF file!

8. With B = 4, it takes a tumor starting at V = 20 about this long to double to V = 40: ANS: _ _ _ _ time units

9. V = _ _ _

10.