PH36010 Numerical Methods

1. PH36010 Numerical Methods Solving Differential Equations using MATHCAD

2. Solving ODEs numerically • Produce numeric solution to system of ODEs. • Must have initial conditions • Use one of several different solvers • Produces matrix of solutions

3. Steps to solving ODEs • Scale equations, parameters & initial conditions to remove units • Manipulate equations to give vector of derivatives • Give vector of initial conditions • Call solver • Plot results

4. Solution of First Order ODE • Radioactive decay, Newton’s law of cooling etc • A is amount of material, temperature difference, etc • k is rate constant

5. Forming the derivative vector • “The derivative of A with respect to t is –k times A”

6. Initial Conditions • First order ODE • 1 member of D(t,A) vector • 1 member of ic vector Although only 1 member, still needs to be a vector

7. Solution parameters • TStart and TFinish times • Number of points, N

8. Create Solution Matrix • Use rkfixed() function • Return matrix with 1 column per DE + 1 for independent variable • Use column operator to strip off columns

9. Plot Results

10. What can go wrong • Found number greater than 10^307 • Use more points • Reduce TFinish • Can’t have anything with dimensions here • Strip units from system before solution

11. More complex first order system • Double decay • A => B => C K1 K2

12. Forming the derivative vector

13. Initial Conditions • System of 2 DEs • 2 members in D(t,F) vector • 2 members in ic vector • F0 refers to A => ic0 • F1 refers to B => ic1 Initially 100% of A, 0% of B

14. Form Solution of Double Decay system • Solution parameters as before • Call rkfixed() as before to create matrix

15. Plot results of double decay

16. Second Order SystemDamped SHM • Rewrite as system of first order equations • Solve as before

17. Second Order systemLCR Circuit Kirchoff’s Voltage Law =>

18. Second order SystemForming the equations Use dq/dt=i and divide by Lcoil To get homogeneous equation in q

19. Second Order SystemWriting standard form #1 Rewrite again, getting rid of d/dt

20. Second Order SystemWriting standard form #2 Use symbolic solver to get q2

21. Second Order SystemFill in D vector • Have expressions for q1 & q2 • Now fill in D(t,q) vector • Replace • q0  q0 ,q1  q1 ,q2  q2

22. Second Order SystemDefine Initial Conditions • Start with 1V across capacitor & no current flowing

23. Second Order SystemCreate Solution • Examine over 0.1s • Use 1000 points

24. Second Order SystemPlot Solution

25. Driven Systems • So far all systems have been ‘Relaxation to steady state’ • Can also model systems driven by ‘Forcing Function’

26. Forcing Function for LCR circuit • Enables us to see resonance • Put voltage source in loop

27. Solution for Forced Oscillation • Use symbolic solver to solve for q2 as before • Compare with undriven case…

28. Forcing function • Drive with sinusoidal waveform • Substitute to give…

29. Derivative Vector for Forced SHM • No initial charge or current

30. Form solution • Use rkfixed as before…

31. Extract values & plot

32. Phase Plot • Plot Current (q1) vs Charge (q0)