MATHEMATICAL MODELING

1. MATHEMATICAL MODELING

2. Example 7: A Three Bar Linkage Problem Consider the three-bar-linkage mechanism shown below. For a constant rotation rate ω of link L1, determine and plot the angular velocities and accelerations of links L2 and L3 for one cycle of rotation of L1. Choose L1, L2 and L3 as 0.35m, 1m and 1m respectively. Also choose ‘a’ and ‘b’ as 0.6m and 0.4m respectively. The angular velocity, ω of link L1 is chosen to be 3 rad/sec.

3. Example 7

4. Solution • From the above figure, the geometric relations can be derived as • First we have to solve for and for values of ranging from 0 to 2π given the values of ‘a’ and ‘b’. The value of alpha should be given in radians.

5. To find out the angular velocities, differentiate the above two equations (geometric relations) with respect to time and put them into matrix form. We get the following equation.

6. To calculate the angular accelerations, the geometric relations given by the second equation set must be differentiated to obtain the following relation,

7. Computational Model • The MATLAB code for solving the above equation is given below. p='cos(b)+cos(c)=0.6-0.35*cos(a)'; q='sin(b)-sin(c)=0.4-0.35*sin(a)'; [b,c]=solve(p,q) i=1; for a=0:0.1:6.3

8. Computational Model beta(:,i)=subs(b) gamma(:,i)=subs(c) i=i+1; end alpha=0:0.1:6.3 figure(1) plot(alpha,beta(1,:));grid xlabel('alpha')

9. Computational Model ylabel('beta') title('beta VS alpha') figure(2) plot(alpha,gamma(1,:));grid xlabel('alpha') ylabel('gama') title('gama VS alpha')

10. Computational Model • The MATLAB code to solve this equation to get the angular velocities and then plot them with respect to is given below. L1 = 0.35; L2 = 1; L3 = 1; omega = 3;

11. Computational Model for i=1:size(alpha,2) A=[-L2*sin(beta(1,i)) –L3*sin(gamma(1,i)); L2*cos(beta(1,i)) -L3*cos(gamma(1,i))]; B=[omega*L1*sin(alpha(i));-omega*L1*cos(alpha(i))]; C=inv(A)*B; betadot(1,i)= C(1); gammadot(1,i)= C(2); end

12. Computational Model figure(3) plot(alpha,betadot(1,:));grid grid xlabel('alpha') ylabel('beta dot') title('betadot VS aplha')

13. Computational Model figure(4) plot(alpha,gammadot(1,:));grid grid xlabel('alpha') ylabel('gamma dot') title('gammadot VS aplha') In the above code ‘A’ contains the coefficient matrix and ‘B’ contains the RHS matrix.

14. Computational Model • To calculate the angular accelerations, the MATLAB code for solving this equation is given below. for i=1:size(alpha,2) A1=[-L2*sin(beta(1,i)) -L3*sin(gamma(1,i)); L2*cos(beta(1,i)) -L3*cos(gamma(1,i))]; B1=[omega*omega*L1*cos(alpha(i))+ L2*(betadot(1,i)^2)*cos(beta(1,i))+ L3*(gammadot(1,i)^2)*cos(gamma(1,i));...

15. Computational Model omega*omega*L1*sin(alpha(i))+(betadot(1,i)^2)*L2*sin(beta(1,i))-(gammadot(1,i)^2)*L3*sin(gamma(1,i))]; C2=inv(A1)*B1; beta_accl(1,i)=C2(1); gamma_accl(1,i)=C2(2); end

16. Computational Model figure(5);plot(alpha,beta_accl(1,:));grid xlabel('alpha') ylabel('beta acceleration') title('beta accl VS alpha') figure(6);plot(alpha,gamma_accl(1,:));grid xlabel('alpha') ylabel('gamma_accl') title('gamma_accl VS alpha') • In the above code ‘A1’ contains the coefficient matrix and ‘B1’ contains the RHS matrix.

17. Plots

18. Plots

19. Plots

20. Plots

21. Plots

22. Plots

23. Conclusions It can be seen from the plots that, curves obtained for and angular velocities and angular accelerations against nn are having exactly opposite shapes.

24. Example 8-Tumor Growth • Cancer, the name given to a group of diseases in which abnormal cells grow and reproduce uncontrollably, is one of the major causes of death in the world. • Under normal conditions, cells in our body grow, divide, die and replace themselves in an orderly, controlled manner. • However, if the process gets out of control, cells can grow too rapidly without any order and develop into a lump which is called a tumor .

25. These tumors can be benign or malignant, and the latter is cancerous. • A malignant tumor essentially consists of cancer cells that are able to spread beyond the original site, invade neighbouring tissues, and spread to other parts of the body in a process called metastasis.

26. The early stages of tumor growth are difficult to be studied clinically as the size of the tumor is too small. • However, experiments have been carried out to study early growth of tumor in vitro, using the multicellular spheroid approach.

27. Three-layer structure of a tumor spheroid As depicted, a typical multicell spheroid consists of an outer shell of proliferating cells, an inner layer of quiescent cells which are dormant but viable, and a central core of necrotic material.

28. Tumor growth: Sherratt-Chaplain Mathematical model • In this model, the cell densities for the proliferating, quiescent and necrotic cells are denoted by p(x, t), q(x, t) and n(x, t) respectively, where t represents time and x is the spatial coordinate. • Necrotic cells are dead and hence non-motile, while proliferating and quiescent cells can move.

29. However, in a close-packed environment such as a tumor spheroid, movement of cells can be restricted. • Overall viable cell flux is given by . • Assume that the two cell populations have equal motility.

30. The movement terms of the proliferating and quiescent cells will be given as respectively. • Suppose proliferating cells grow at a rate that is limited by crowding of the total cell population, and that they become quiescent at a rate that depends on the concentration c(x, t) of some nutrients.

31. Suppose quiescent cells become necrotic also at a rate that depends on c(x, t). • Then, we can write down the following set of model equations.

32. The value of 1 corresponds to a completely close packed population. • g(0) is set to 1 to provide a suitable initial condition. • The functions f and h are decreasing functions, with both tending to zero as c tends to +1, while g is an increasing function.

33. Since the rate of cells becoming quiescent is normally higher than the rate of cells entering necrosis, we expect f(c) > h(c). • Given these conditions, one can specify suitable or appropriate functional forms for f, g and h. where α and γ are dimensionless parameters.

34. The system of equations may be solved numerically using the finite difference formulations for the equations. • Using forward differencing for time and central differencing for space, we obtain the following set of finite difference equations.

36. In the above finite difference equations, Δt and Δx refer to the time steps and space intervals respectively, and the superscript for the dependent variables indicate the time level and space position respectively. • For instance, pji would refer to the density of the proliferating cells at the jth time level and ith space position.

37. Computational Model • we choose f(c) = (1−tanh(4c−2))/2 and h(c) = f(c)/2, g(c) =βeβc with β= 0.5.

38. Initial and boundary conditions • Let α= 0.8 and γ = 10 (chosen arbitrarily). • Assume that at t = 0, q(x, 0) = 0, n(x, 0) = 0, c = 1, and the proliferating cell density decreases exponentially as x increases. • Thus, we let p(x, 0) = e0.1x.

39. As for the boundary conditions, it is usual to assume zero-flux at the boundaries. • In terms of implementation, it is not possible to let x .

40. We therefore need to choose a value of x that is sufficiently large. • After experimenting with various values, it is found that x = 210 is appropriate and good enough to observe patterns of the evolution of cell densities.

41. Model Calibration Comparison of cell counts between calibrated model and the experimental data.

42. The model is calibrated by choosing an appropriate scaling factor to transform cell densities to match the cell counts in the experimental data. • The reference point in this calibration is chosen to be t = 1, at which time the experimental total cell count is about 7015. • The total cell density from the model is around 15.809; hence a scaling factor of about 443.7249 is used. • With this value, the total live cells (that is, proliferating and quiescent) and dead cells (necrotic) as predicted by the model can be compared against the experimental data.

43. Results and Discussion

44. It is observed that an advancing pulse of proliferating cells is accompanied by a corresponding band of quiescent cells as time increases. • The necrotic cells seem to be concentrated at the core initially (t = 2 to t = 8). • As time passes, the necrotic cells continue to develop at the core, but begin to disperse towards the outer edge of the tumor spheroid.

45. Although necrotic cells are building up with time, no tumor regression is detected. • This is expected as the model does not restrict or limit nutrient flow sufficiently, and no therapy or control measure is introduced. • Both the current model and the experiment results report no limiting spheroid volume.

46. Visualizing Tumor Growth • Snapshots of simulated tumor growth at t = 2, 4, . . . , 16 units, with blue, red and black colored dots representing proliferating, quiescent and necrotic cells respectively.

47. Conclusion • The model had produced results which compare fairly well with the experimental results. • This model provides qualitative results which can be further analyzed and studied. • The use of MATLAB in this discussion illustrates the important role of technology in research in mathematical modeling. • Not only does MATLAB help in providing a computing platform for implementing the numerical scheme efficiently, it also serves as a useful tool for generating visual images resulting from the model.

48. Thank You!