%23%20%2F%2F%2F%20script%0A%23%20requires-python%20%3D%20%22%3E%3D3.11%22%0A%23%20dependencies%20%3D%20%5B%22marimo%3E%3D0.23.9%22%2C%20%22numpy%3E%3D2.0%22%2C%20%22matplotlib%3E%3D3.9%22%2C%20%22scipy%3E%3D1.14%22%5D%0A%23%20%2F%2F%2F%0A%0Aimport%20marimo%0A%0A__generated_with%20%3D%20%220.25.2%22%0Aapp%20%3D%20marimo.App(width%3D%22medium%22)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%20Introduction%20to%20numerical%20solution%20of%20Ordinary%20Differential%20Equation%0A%0A%20%20%20%20%3E%20Marcos%20Duarte%2C%20Renato%20Naville%20Watanabe%3Cbr%3E%0A%20%20%20%20%3E%20%5BLaboratory%20of%20Biomechanics%20and%20Motor%20Control%5D(https%3A%2F%2Fbmclab.pesquisa.ufabc.edu.br%2F)%3Cbr%3E%0A%20%20%20%20%3E%20Federal%20University%20of%20ABC%2C%20Brazil%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20How%20to%20use%20this%20tutorial%0A%0A%20%20%20%20This%20notebook%20is%20written%20as%20a%20guided%20study%20session.%20You%20will%20read%20a%20short%20explanation%2C%20predict%20what%20should%20happen%2C%20run%20a%20small%20computation%2C%20and%20then%20use%20the%20result%20to%20answer%20the%20next%20question.%20Keep%20a%20scratchpad%20nearby%20and%20write%20your%20answers%20before%20you%20reveal%20the%20plots.%0A%0A%20%20%20%20You%20do%20not%20need%20to%20finish%20every%20challenge%20on%20the%20first%20pass.%20Start%20by%20running%20each%20cell%20in%20order.%20Then%20return%20to%20the%20challenges%20and%20change%20one%20parameter%20at%20a%20time.%20The%20goal%20is%20not%20only%20to%20make%20the%20code%20run%3B%20the%20goal%20is%20to%20connect%20the%20equations%2C%20the%20numerical%20method%2C%20and%20the%20motion%20you%20see%20in%20the%20plots.%0A%0A%20%20%20%20**Challenge%200.**%20Before%20you%20begin%2C%20write%20down%20one%20mechanical%20system%20you%20already%20know%20that%20changes%20over%20time%3A%20a%20bouncing%20ball%2C%20a%20pendulum%2C%20a%20person%20walking%2C%20a%20mass%20on%20a%20spring%2C%20or%20something%20else.%20What%20variable%20would%20you%20measure%20to%20describe%20its%20state%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Python%20setup%0A%0A%20%20%20%20NumPy%20stores%20the%20numerical%20states%2C%20Matplotlib%20draws%20the%20trajectories%2C%20and%20SciPy%20provides%20%60solve_ivp%60.%20All%20inputs%20are%20stated%20model%20parameters%3B%20no%20external%20dataset%20is%20needed.%20Run%20the%20examples%20and%20compare%20step%20sizes%20before%20trusting%20a%20trajectory.%20Expand%20exercise%20solutions%20only%20after%20trying%20the%20update%20yourself.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20An%20unassuming%20derivation%20of%20a%20numerical%20solution%20to%20an%20Ordinary%20Differential%20Equation%20(ODE)%0A%0A%20%20%20%20Imagine%20you%20are%20in%20a%20car.%20You%20want%20to%20know%20where%20the%20car%20is%20at%20each%20moment%2C%20but%20you%20do%20not%20have%20a%20formula%20for%20its%20position.%20What%20you%20do%20have%20is%20simpler%3A%20the%20initial%20position%20of%20the%20car%20and%20measurements%20of%20its%20speed%20from%20time%20to%20time.%0A%0A%20%20%20%20That%20is%20already%20enough%20to%20make%20progress.%0A%0A%20%20%20%20Velocity%20tells%20us%20how%20fast%20position%20is%20changing%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20v(t)%20%3D%20%5Cfrac%7B%5Cmathrm%7Bd%7Dx(t)%7D%7B%5Cmathrm%7Bd%7Dt%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20The%20notation%20%5C(x(t)%5C)%20and%20%5C(v(t)%5C)%20is%20just%20a%20reminder%20that%20both%20position%20and%20velocity%20can%20change%20with%20time.%0A%0A%20%20%20%20If%20we%20look%20at%20two%20nearby%20instants%2C%20we%20can%20approximate%20the%20velocity%20by%20asking%20how%20much%20the%20position%20changed%20during%20that%20time%20interval%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20v_i%20%5Capprox%20%5Cfrac%7Bx_%7Bi%2B1%7D%20-%20x_i%7D%7Bt_%7Bi%2B1%7D%20-%20t_i%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Here%2C%20%5C(i%5C)%20labels%20one%20measured%20instant%20and%20%5C(i%2B1%5C)%20labels%20the%20next%20one.%20To%20keep%20the%20first%20idea%20simple%2C%20suppose%20the%20measurements%20are%20equally%20spaced%20in%20time.%20Then%0A%0A%20%20%20%20%24%24%0A%20%20%20%20t_%7Bi%2B1%7D%20-%20t_i%20%3D%20%5CDelta%20t%2C%0A%20%20%20%20%24%24%0A%0A%20%20%20%20and%20the%20approximation%20becomes%0A%0A%20%20%20%20%24%24%0A%20%20%20%20v_i%20%5Capprox%20%5Cfrac%7Bx_%7Bi%2B1%7D%20-%20x_i%7D%7B%5CDelta%20t%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Now%20comes%20the%20useful%20trick.%20Rearrange%20this%20expression%20to%20isolate%20the%20next%20position%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20x_%7Bi%2B1%7D%20%5Capprox%20x_i%20%2B%20v_i%20%5CDelta%20t.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20This%20says%20something%20very%20practical%3A%20if%20you%20know%20where%20the%20car%20is%20now%2C%20and%20you%20know%20its%20current%20speed%2C%20you%20can%20estimate%20where%20it%20will%20be%20after%20one%20small%20time%20step.%0A%0A%20%20%20%20For%20example%2C%20suppose%0A%0A%20%20%20%20%24%24%0A%20%20%20%20x_0%20%3D%20100%5C%20%5Cmathrm%7Bm%7D%2C%20%5Cqquad%20v_0%20%3D%2020%5C%20%5Cmathrm%7Bm%2Fs%7D%2C%20%5Cqquad%20%5CDelta%20t%20%3D%2010%5C%20%5Cmathrm%7Bs%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Then%0A%0A%20%20%20%20%24%24%0A%20%20%20%20x_1%20%5Capprox%20100%20%2B%2020%20%5Ccdot%2010%20%3D%20300%5C%20%5Cmathrm%7Bm%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20The%20index%20%5C(1%5C)%20means%20the%20next%20measured%20instant%2C%20which%20in%20this%20example%20is%20%5C(t%3D10%5C)%20s.%0A%0A%20%20%20%20Now%20suppose%20that%20at%20%5C(t%3D10%5C)%20s%20the%20measured%20speed%20is%0A%0A%20%20%20%20%24%24%0A%20%20%20%20v_1%20%3D%2025%5C%20%5Cmathrm%7Bm%2Fs%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Because%20we%20already%20estimated%20%5C(x_1%5C)%2C%20we%20can%20repeat%20the%20same%20step%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20x_2%20%5Capprox%20x_1%20%2B%20v_1%20%5CDelta%20t%2C%0A%20%20%20%20%24%24%0A%0A%20%20%20%20so%0A%0A%20%20%20%20%24%24%0A%20%20%20%20x_2%20%5Capprox%20300%20%2B%2025%20%5Ccdot%2010%20%3D%20550%5C%20%5Cmathrm%7Bm%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Another%20way%20to%20see%20the%20same%20calculation%20is%20to%20keep%20a%20running%20table.%20The%20position%20in%20each%20row%20is%20the%20position%20after%20applying%20the%20speed%20from%20the%20previous%20interval%3A%0A%0A%20%20%20%20%7C%20%5C(i%5C)%20%7C%20Time%20%5C(t_i%5C)%20%5Bs%5D%20%7C%20Speed%20used%20over%20previous%20interval%20%5Bm%2Fs%5D%20%7C%20Position%20%5C(x_i%5C)%20%5Bm%5D%20%7C%0A%20%20%20%20%7C---%3A%7C---%3A%7C---%3A%7C---%3A%7C%0A%20%20%20%20%7C%200%20%7C%200%20%7C%20--%20%7C%20100%20%7C%0A%20%20%20%20%7C%201%20%7C%2010%20%7C%2020%20%7C%20300%20%7C%0A%20%20%20%20%7C%202%20%7C%2020%20%7C%2025%20%7C%20550%20%7C%0A%20%20%20%20%7C%203%20%7C%2030%20%7C%2030%20%7C%20%3F%20%7C%0A%0A%20%20%20%20The%20first%20row%20is%20the%20initial%20condition.%20The%20second%20row%20uses%20%5C(v_0%3D20%5C%20%5Cmathrm%7Bm%2Fs%7D%5C)%2C%20and%20the%20third%20row%20uses%20%5C(v_1%3D25%5C%20%5Cmathrm%7Bm%2Fs%7D%5C).%20Calculate%20the%20new%20position%20at%20%5C(t%3D30%5C%20%5Cmathrm%7Bs%7D%5C).%0A%0A%20%20%20%20And%20there%20is%20the%20pattern%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20x_%7Bi%2B1%7D%20%5Capprox%20x_i%20%2B%20v_i%20%5CDelta%20t.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20To%20estimate%20the%20position%20farther%20into%20the%20future%2C%20we%20do%20not%20jump%20there%20all%20at%20once.%20We%20move%20step%20by%20step%2C%20always%20using%20the%20current%20position%20and%20the%20current%20velocity%20to%20estimate%20the%20next%20position.%0A%0A%20%20%20%20This%20is%20the%20basic%20idea%20behind%20Euler's%20method%20for%20numerically%20solving%20a%20first-order%20ordinary%20differential%20equation%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cfrac%7B%5Cmathrm%7Bd%7Dx(t)%7D%7B%5Cmathrm%7Bd%7Dt%7D%20%3D%20v(t).%0A%20%20%20%20%24%24%0A%0A%20%20%20%20In%20this%20example%3A%0A%0A%20%20%20%20-%20%5C(x_0%5C)%20is%20the%20initial%20condition%3B%0A%20%20%20%20-%20%5C(v(t)%5C)%20gives%20the%20rate%20of%20change%20of%20position%3B%0A%20%20%20%20-%20%5C(%5CDelta%20t%5C)%20is%20the%20time%20step%3B%0A%20%20%20%20-%20%5C(x_%7Bi%2B1%7D%20%5Capprox%20x_i%20%2B%20v_i%5CDelta%20t%5C)%20is%20the%20Euler%20update.%0A%0A%20%20%20%20Notice%20what%20we%20obtained%3A%20a%20list%20of%20estimated%20positions%20at%20selected%20times.%20We%20did%20not%20obtain%20a%20single%20formula%20for%20%5C(x(t)%5C).%20That%20is%20what%20makes%20this%20a%20numerical%20solution.%20Because%20we%20started%20with%20speed%20and%20estimated%20position%20by%20summing%20increments%2C%20this%20process%20is%20numerical%20integration.%0A%20%20%20%20Because%20the%20calculation%20starts%20from%20a%20known%20initial%20value%2C%20this%20kind%20of%20problem%20is%20called%20an%20initial%20value%20problem%2C%20or%20IVP.%0A%0A%20%20%20%20**The%20case%20when%20only%20the%20acceleration%20is%20known**%0A%20%20%20%20Now%20imagine%20that%20we%20only%20have%20access%20to%20the%20car's%20acceleration%2C%20for%20example%20from%20a%20smartphone%20accelerometer%2C%20and%20that%20we%20know%20the%20initial%20position%20and%20speed.%20Can%20we%20still%20find%20the%20position%20at%20a%20later%20time%3F%20Let's%20try%20a%20similar%20approach%20for%20this%20more%20challenging%20problem.%0A%0A%20%20%20%20Acceleration%20tells%20us%20how%20fast%20velocity%20is%20changing%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20a(t)%20%3D%20%5Cfrac%7B%5Cmathrm%7Bd%7Dv(t)%7D%7B%5Cmathrm%7Bd%7Dt%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20If%20the%20acceleration%20is%20measured%20at%20equally%20spaced%20instants%2C%20we%20can%20use%20the%20same%20finite-difference%20idea%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20a_i%20%5Capprox%20%5Cfrac%7Bv_%7Bi%2B1%7D%20-%20v_i%7D%7B%5CDelta%20t%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Rearranging%20gives%20an%20update%20for%20the%20next%20speed%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20v_%7Bi%2B1%7D%20%5Capprox%20v_i%20%2B%20a_i%5CDelta%20t.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Alongside%20that%2C%20we%20keep%20using%20the%20current%20speed%20to%20update%20position.%20So%20the%20state%20we%20carry%20from%20one%20row%20to%20the%20next%20is%20not%20just%20position%20anymore%3B%20it%20is%20the%20pair%20%5C((x_i%2C%20v_i)%5C).%20In%20the%20simplest%20Euler%20version%2C%20the%20next%20state%20is%20defined%20by%20these%20two%20coupled%20equations%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cleft%5C%7B%0A%20%20%20%20%5Cbegin%7Barray%7D%7Bl%7D%0A%20%20%20%20x_%7Bi%2B1%7D%20%5Capprox%20x_i%20%2B%20v_i%5CDelta%20t%2C%5C%5C%0A%20%20%20%20v_%7Bi%2B1%7D%20%5Capprox%20v_i%20%2B%20a_i%5CDelta%20t.%0A%20%20%20%20%5Cend%7Barray%7D%0A%20%20%20%20%5Cright.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20They%20are%20coupled%20because%20the%20position%20update%20needs%20the%20current%20speed%2C%20and%20the%20speed%20update%20needs%20the%20current%20acceleration.%20This%20idea%20will%20come%20back%20later%20when%20we%20rewrite%20higher-order%20ODEs%20as%20systems%20of%20first-order%20equations.%0A%0A%20%20%20%20Suppose%20the%20car%20starts%20at%20%5C(x_0%3D100%5C%20%5Cmathrm%7Bm%7D%5C)%2C%20starts%20with%20%5C(v_0%3D20%5C%20%5Cmathrm%7Bm%2Fs%7D%5C)%2C%20and%20%5C(%5CDelta%20t%3D10%5C%20%5Cmathrm%7Bs%7D%5C).%20If%20the%20first%20three%20acceleration%20measurements%20are%20%5C(a_0%3D0.5%5C%20%5Cmathrm%7Bm%2Fs%5E2%7D%5C)%2C%20%5C(a_1%3D0.5%5C%20%5Cmathrm%7Bm%2Fs%5E2%7D%5C)%2C%20and%20%5C(a_2%3D0.5%5C%20%5Cmathrm%7Bm%2Fs%5E2%7D%5C)%2C%20the%20calculation%20looks%20like%20this%3A%0A%0A%20%20%20%20%7C%20%5C(i%5C)%20%7C%20Time%20%5C(t_i%5C)%20%5Bs%5D%20%7C%20Acceleration%20used%20over%20next%20interval%20%5Bm%2Fs%5E2%5D%20%7C%20Speed%20%5C(v_i%5C)%20%5Bm%2Fs%5D%20%7C%20Position%20%5C(x_i%5C)%20%5Bm%5D%20%7C%0A%20%20%20%20%7C---%3A%7C---%3A%7C---%3A%7C---%3A%7C---%3A%7C%0A%20%20%20%20%7C%200%20%7C%200%20%7C%200.5%20%7C%2020%20%7C%20100%20%7C%0A%20%20%20%20%7C%201%20%7C%2010%20%7C%200.5%20%7C%2025%20%7C%20300%20%7C%0A%20%20%20%20%7C%202%20%7C%2020%20%7C%200.5%20%7C%2030%20%7C%20550%20%7C%0A%20%20%20%20%7C%203%20%7C%2030%20%7C%20--%20%7C%2035%20%7C%20%3F%20%7C%0A%0A%20%20%20%20Use%20the%20row%20at%20%5C(t%3D20%5C%20%5Cmathrm%7Bs%7D%5C)%20to%20calculate%20the%20new%20position%20at%20%5C(t%3D30%5C%20%5Cmathrm%7Bs%7D%5C).%20Here%2C%20acceleration%20is%20first%20accumulated%20into%20speed%2C%20and%20speed%20is%20then%20accumulated%20into%20position.%20This%20is%20a%20double%20numerical%20integration%3A%20acceleration%20to%20velocity%2C%20then%20velocity%20to%20position.%0A%0A%20%20%20%20Next%20we%20will%20switch%20to%20a%20more%20formal%20description%20of%20ordinary%20differential%20equations%20and%20work%20through%20more%20examples.%20If%20a%20problem%20starts%20to%20feel%20abstract%2C%20come%20back%20to%20the%20car%3A%20ask%20which%20quantities%20play%20the%20role%20of%20position%20and%20velocity%2C%20and%20which%20equation%20is%20giving%20you%20the%20acceleration-like%20rate%20of%20change.%20Once%20those%20pieces%20are%20named%2C%20the%20numerical%20task%20is%20the%20same%20kind%20of%20step-by-step%20update.%0A%0A%20%20%20%20**Challenge%201.**%20Write%20a%20pseudocode%20to%20implement%20the%20Euler's%20method%20for%20numerically%20solving%20a%20first-order%20ordinary%20differential%20equation.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%201.**%20Complete%20the%20car%20table%20at%2030%20s.%20What%20would%20the%20exact%20position%20be%20under%20constant%20acceleration%200.5%20m%2Fs%C2%B2%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%201%0A%20%20%20%20Forward%20Euler%20gives%20%24550%2B30(10)%3D850%24%20m.%20The%20exact%20constant-acceleration%20result%20is%20%24100%2B20(30)%2B0.5(0.5)(30)%5E2%3D925%24%20m.%20The%20method%20used%20the%20speed%20at%20the%20start%20of%20each%20interval%20and%20therefore%20underestimated%20the%20distance%20during%20acceleration.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Ordinary%20Differential%20Equation%0A%0A%20%20%20%20An%20ODE%20is%20an%20equation%20that%20relates%20a%20function%20of%20one%20independent%20variable%20to%20one%20or%20more%20of%20its%20derivatives.%20In%20this%20tutorial%2C%20the%20independent%20variable%20will%20usually%20be%20time.%20Solving%20an%20ODE%20means%20finding%20a%20function%20whose%20derivatives%20satisfy%20that%20relationship.%0A%0A%20%20%20%20The%20order%20of%20an%20ODE%20is%20the%20order%20of%20its%20highest%20derivative.%20For%20example%2C%20a%20first-order%20ODE%20contains%20only%20first%20derivatives.%20A%20linear%20ODE%20is%20linear%20in%20the%20unknown%20function%20and%20its%20derivatives%3B%20many%20important%20linear%20ODEs%20have%20analytical%20solutions%2C%20but%20others%20are%20still%20solved%20numerically%20in%20practice.%20Nonlinear%20ODEs%20can%20sometimes%20be%20solved%20exactly%2C%20but%20numerical%20methods%20are%20often%20the%20practical%20approach.%0A%0A%20%20%20%20A%20partial%20differential%20equation%20contains%20partial%20derivatives%20of%20an%20unknown%20function%20of%20several%20independent%20variables%2C%20such%20as%20both%20position%20and%20time.%0A%0A%20%20%20%20Newton's%20second%20law%20gives%20a%20simple%20example%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20m%5Cfrac%7B%5Cmathrm%7Bd%7D%5E2%20%5Cmathbf%7Bx%7D%7D%7B%5Cmathrm%7Bd%7Dt%5E2%7D(t)%20%3D%20%5Cmathbf%7BF%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Here%2C%20%24%5Cmathbf%7Bx%7D(t)%24%20is%20the%20position%20as%20a%20function%20of%20time%20and%20%24t%24%20is%20the%20independent%20variable.%20The%20force%20%24%5Cmathbf%7BF%7D%24%20may%20be%20constant%2C%20such%20as%20gravitational%20force%2C%20or%20it%20may%20depend%20on%20position%2C%20velocity%2C%20time%2C%20or%20other%20quantities.%20If%20%24%5Cmathbf%7BF%7D%24%20is%20constant%20or%20a%20linear%20function%20of%20%24%5Cmathbf%7Bx%7D%24%2C%20this%20equation%20is%20a%20second-order%20linear%20ODE.%0A%0A%20%20%20%20**Guiding%20questions%201.**%0A%20%20%20%201.%20Which%20quantity%20is%20changing%20with%20time%20in%20Newton's%20second%20law%3F%0A%20%20%20%202.%20Which%20derivative%20appears%20in%20the%20equation%3F%0A%20%20%20%203.%20What%20extra%20information%20would%20you%20need%20before%20predicting%20a%20unique%20motion%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20First-order%20ODE%0A%0A%20%20%20%20A%20first-order%20ODE%20has%20the%20general%20form%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cfrac%7B%5Cmathrm%7Bd%7D%20y%7D%7B%5Cmathrm%7Bd%7D%20x%7D%20%3D%20f(x%2C%20y).%0A%20%20%20%20%24%24%0A%0A%20%20%20%20The%20function%20%24f(x%2C%20y)%24%20gives%20the%20derivative%20of%20%24y%24%20at%20a%20given%20pair%20%24(x%2C%20y)%24.%20When%20%24f(x%2C%20y)%24%20is%20linear%20with%20respect%20to%20%24y%24%2C%20the%20equation%20is%20a%20first-order%20linear%20ODE%2C%20which%20can%20be%20written%20as%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cfrac%7B%5Cmathrm%7Bd%7D%20y%7D%7B%5Cmathrm%7Bd%7D%20x%7D%20%2B%20p(x)y%20%3D%20q(x)%0A%20%20%20%20%24%24%0A%0A%20%20%20%20where%20%5C(p(x)%5C)%20and%20%5C(q(x)%5C)%20are%20continuous%20functions%20of%20%5C(x%5C).%20In%20this%20form%2C%20the%20equation's%20linearity%20means%20the%20dependent%20variable%20%5C(y%5C)%20and%20its%20derivative%20%24%5Cmathrm%7Bd%7Dy%2F%5Cmathrm%7Bd%7Dx%24%20are%20raised%20to%20the%20first%20power%20and%20are%20not%20multiplied%20together.%0A%0A%20%20%20%20**Challenge%202.**%20Suppose%20%24y%24%20is%20position%20and%20%24x%24%20is%20time.%20In%20one%20sentence%2C%20explain%20why%20knowing%20only%20%24%5Cmathrm%7Bd%7Dy%2F%5Cmathrm%7Bd%7Dx%24%20at%20a%20single%20instant%20is%20not%20enough%20to%20reconstruct%20the%20full%20motion.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Initial%20Value%20Problems%20(IVP)%0A%0A%20%20%20%20When%20solving%20ODEs%20with%20a%20known%20set%20of%20initial%20conditions%2C%20we%20are%20solving%20an%20initial%20value%20problem%20(IVP).%20The%20solution%20often%20describes%20how%20a%20system%20changes%20over%20time%20from%20a%20specified%20initial%20state.%0A%0A%20%20%20%20A%20solution%20to%20an%20IVP%20is%20a%20function%20that%20satisfies%20both%20the%20differential%20equation%20and%20the%20initial%20condition.%20For%20a%20first-order%20problem%2C%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cdot%7By%7D(t)%20%3D%20f(t%2C%20y(t))%2C%20%5Cqquad%20y(t_0)%3Dy_0.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Higher-order%20IVPs%20can%20be%20rewritten%20as%20systems%20of%20first-order%20IVPs.%20For%20example%2C%20a%20second-order%20equation%20can%20be%20written%20as%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cddot%7By%7D(t)%20%3D%20f(t%2C%20y(t)%2C%20%5Cdot%7By%7D(t)).%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Define%20a%20state%20vector%20with%20the%20original%20variable%20and%20its%20derivative%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20x_1(t)%20%3D%20y(t)%2C%20%5Cqquad%20x_2(t)%20%3D%20%5Cdot%7By%7D(t).%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Then%20the%20same%20problem%20becomes%20two%20coupled%20first-order%20equations%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cdot%7Bx%7D_1(t)%3Dx_2(t)%2C%20%5Cqquad%20%5Cdot%7Bx%7D_2(t)%3Df(t%2C%20x_1(t)%2C%20x_2(t)).%0A%20%20%20%20%24%24%0A%0A%20%20%20%20In%20the%20same%20way%2C%20an%20ODE%20of%20order%20%24N%24%20can%20be%20represented%20as%20a%20system%20of%20%24N%24%20first-order%20ODEs.%20A%20simple%20pendulum%2C%20a%20projectile%20under%20gravity%2C%20and%20many%20biomechanical%20models%20are%20naturally%20written%20this%20way%20before%20numerical%20integration.%0A%0A%20%20%20%20**Guiding%20question%202.**%0A%20%20%20%201.%20What%20does%20the%20initial%20condition%20%24y(t_0)%3Dy_0%24%20tell%20you%3F%0A%0A%20%20%20%20**Challenge%203.**%20A%20mass-spring-damper%20model%20can%20be%20written%20as%20%24m%20%5Cddot%7Bx%7D%20%2B%20c%20%5Cdot%7Bx%7D%20%2B%20kx%20%3D%20F(t)%24.%20Before%20looking%20ahead%2C%20define%20two%20state%20variables%20that%20would%20convert%20this%20equation%20into%20two%20first-order%20ODEs.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%202.**%20Write%20the%20first-order%20state%20equations%20for%20%24m%5Cddot%20x%2Bc%5Cdot%20x%2Bkx%3DF(t)%24.%20Which%20initial%20values%20are%20needed%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%202%0A%20%20%20%20With%20state%20%24(x%2Cv)%24%2C%20%24%5Cdot%20x%3Dv%24%20and%20%24%5Cdot%20v%3D%5BF(t)-cv-kx%5D%2Fm%24.%20Specify%20both%20%24x(t_0)%24%20and%20%24v(t_0)%24%3B%20an%20acceleration%20law%20alone%20does%20not%20locate%20a%20unique%20trajectory.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Numerical%20methods%20for%20solving%20ODEs%0A%0A%20%20%20%20When%20an%20ODE%20is%20difficult%20or%20impossible%20to%20solve%20analytically%2C%20numerical%20methods%20approximate%20the%20solution%20at%20selected%20points.%20This%20process%20is%20also%20called%20numerical%20integration.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Euler%20method%0A%0A%20%20%20%20The%20Euler%20method%20is%20the%20simplest%20numerical%20method%20for%20solving%20an%20initial%20value%20problem.%20First%2C%20approximate%20the%20derivative%20of%20%24y%24%20with%20a%20finite%20difference%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cdot%7By%7D(t)%20%5Capprox%20%5Cfrac%7By(t%2Bh)-y(t)%7D%7Bh%7D%2C%0A%20%20%20%20%24%24%0A%0A%20%20%20%20where%20%24h%24%20is%20the%20step%20size.%20Rearranging%20gives%0A%0A%20%20%20%20%24%24%0A%20%20%20%20y(t%2Bh)%20%5Capprox%20y(t)%20%2B%20h%5Cdot%7By%7D(t).%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Replacing%20%24%5Cdot%7By%7D(t)%24%20by%20%24f(t%2C%20y(t))%24%20gives%0A%0A%20%20%20%20%24%24%0A%20%20%20%20y(t%2Bh)%20%5Capprox%20y(t)%20%2B%20hf(t%2C%20y(t)).%0A%20%20%20%20%24%24%0A%0A%20%20%20%20Starting%20from%20the%20known%20value%20%24y(t_0)%3Dy_0%24%2C%20this%20becomes%20the%20recursive%20update%0A%0A%20%20%20%20%24%24%0A%20%20%20%20y_%7Bn%2B1%7D%20%3D%20y_n%20%2B%20h%20f(t_n%2C%20y_n).%0A%20%20%20%20%24%24%0A%0A%0A%20%20%20%20**Guiding%20questions%203.**%0A%20%20%20%201.%20Why%20do%20you%20think%20numerical%20integration%20advances%20in%20small%20time%20steps%3F%0A%20%20%20%202.%20What%20might%20go%20wrong%20if%20the%20time%20step%20is%20too%20large%3F%0A%0A%20%20%20%20**Challenge%204.**%20Euler's%20method%20is%20simple%2C%20but%20it%20makes%20a%20local%20straight-line%20prediction.%20Sketch%20what%20you%20expect%20to%20happen%20when%20the%20real%20trajectory%20is%20strongly%20curved%20and%20the%20step%20size%20%24h%24%20is%20large.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Euler's%20method%20pseudocode%0A%0A%20%20%20%20%60%60%60bash%0A%20%20%20%20%23%20Initialization%0A%20%20%20%20Define%20the%20differential%20equation%20fdot(t%2C%20y)%0A%20%20%20%20Input%20initial%20time%20t0%2C%20initial%20value%20y0%0A%20%20%20%20Input%20step%20size%20h%2C%20and%20the%20number%20of%20steps%20n%0A%20%20%20%20Create%20array%20and%20store%20initial%20values%0A%0A%20%20%20%20%23%20Iteration%0A%20%20%20%20for%20i%20from%201%20to%20n%3A%0A%20%20%20%20%20%20%20%20%23%20Calculate%20the%20slope%20at%20the%20current%20point%0A%20%20%20%20%20%20%20%20slope%20%3D%20fdot(t0%2C%20y0)%0A%20%20%20%20%20%20%20%20%23%20Calculate%20the%20next%20time%20and%20y-value%0A%20%20%20%20%20%20%20%20ti%20%3D%20t0%20%2B%20h%0A%20%20%20%20%20%20%20%20yi%20%3D%20y0%20%2B%20h%20*%20slope%0A%20%20%20%20%20%20%20%20%23%20Store%20new%20values%20(ti%2C%20yi)%0A%20%20%20%20%20%20%20%20array%5Bi%5D%20%3D%20ti%2C%20yi%0A%20%20%20%20%20%20%20%20%23%20Update%20variables%20for%20the%20next%20iteration%0A%20%20%20%20%20%20%20%20t0%20%3D%20ti%0A%20%20%20%20%20%20%20%20y0%20%3D%20yi%0A%20%20%20%20%60%60%60%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Example%3A%20Exact%20and%20approximate%20solutions%20for%20an%20exponential%20growth%20ODE%0A%0A%20%20%20%20Consider%20%24%5Cdot%7Bf(t)%7D%3D.8*e%5E%7B2*t%7D%24%20with%20initial%20condition%20%24f_%7B0%7D%3D.4%24.%20In%20this%20case%2C%20the%20exact%20solution%20is%20known%3A%20%24f(t)%3D.4e%5E%7B2*t%7D%24.%0A%20%20%20%20Let's%20solve%20this%20ODE%20numerically%20and%20compare%20the%20solutions%20in%20the%20interval%20%24t%5B0%2C1%5D%24.%0A%0A%20%20%20%20First%2C%20let's%20import%20the%20Python%20libraries%20we%20will%20use%20and%20configure%20some%20settings%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20numpy%20as%20np%0A%0A%20%20%20%20import%20matplotlib%0A%20%20%20%20import%20matplotlib.pyplot%20as%20plt%0A%0A%20%20%20%20matplotlib.rcParams%5B%22lines.linewidth%22%5D%20%3D%203%0A%20%20%20%20matplotlib.rcParams%5B%22font.size%22%5D%20%3D%2013%0A%20%20%20%20matplotlib.rcParams%5B%22lines.markersize%22%5D%20%3D%205%0A%20%20%20%20matplotlib.rc(%22axes%22%2C%20grid%3DFalse%2C%20labelsize%3D14%2C%20titlesize%3D16%2C%20ymargin%3D0.05)%0A%20%20%20%20matplotlib.rc(%22legend%22%2C%20numpoints%3D1%2C%20fontsize%3D11)%0A%20%20%20%20return%20np%2C%20plt%0A%0A%0A%40app.cell%0Adef%20_(np%2C%20plt)%3A%0A%20%20%20%20def%20exp_growth(step%3D0.1)%3A%0A%20%20%20%20%20%20%20%20def%20fdot(t%2C%20y)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20return%200.8%20*%20np.exp(2%20*%20t)%0A%0A%20%20%20%20%20%20%20%20t0%2C%20y0%20%3D%200.0%2C%200.4%0A%20%20%20%20%20%20%20%20h%20%3D%20step%0A%20%20%20%20%20%20%20%20n%20%3D%20int(1%20%2F%20h%20%2B%201)%0A%20%20%20%20%20%20%20%20array%20%3D%20np.empty(shape%3D%5Bn%2C%202%5D)%0A%20%20%20%20%20%20%20%20array%5B0%2C%20%3A%5D%20%3D%20t0%2C%20y0%0A%0A%20%20%20%20%20%20%20%20for%20i%20in%20range(1%2C%20n)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20slope%20%3D%20fdot(t0%2C%20y0)%0A%20%20%20%20%20%20%20%20%20%20%20%20ti%20%3D%20t0%20%2B%20h%0A%20%20%20%20%20%20%20%20%20%20%20%20yi%20%3D%20y0%20%2B%20h%20*%20slope%0A%20%20%20%20%20%20%20%20%20%20%20%20array%5Bi%2C%20%3A%5D%20%3D%20ti%2C%20yi%0A%20%20%20%20%20%20%20%20%20%20%20%20t0%2C%20y0%20%3D%20ti%2C%20yi%0A%0A%20%20%20%20%20%20%20%20t%2C%20f_num%20%3D%20array.T%0A%20%20%20%20%20%20%20%20f_lit%20%3D%200.4%20*%20np.exp(2%20*%20t)%0A%20%20%20%20%20%20%20%20plt.figure()%0A%20%20%20%20%20%20%20%20plt.plot(t%2C%20f_num%2C%20%22bo%3A%22%2C%20label%3D%22Approximate%22)%0A%20%20%20%20%20%20%20%20plt.plot(t%2C%20f_lit%2C%20%22r%22%2C%20label%3D%22Exact%22)%0A%20%20%20%20%20%20%20%20plt.legend()%0A%20%20%20%20%20%20%20%20plt.show()%0A%0A%20%20%20%20exp_growth()%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Guiding%20questions%204**%0A%20%20%20%201.%20Within%20the%20iteration%20loop%2C%20which%20rows%20can%20be%20reordered%20and%20which%20cannot%3F%0A%20%20%20%202.%20Why%20pre-assign%20an%20array%3F%0A%0A%20%20%20%20**Challenge%205**%20Implement%20this%20solution%20without%20using%20a%20pre-assigned%20array.%0A%0A%20%20%20%20**Challenge%206**%20Modify%20the%20code%20to%20increase%20the%20step%20and%20compare%20the%20exact%20and%20approximate%20solutions.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%203.**%20Use%20one%20forward%20Euler%20step%20for%20%24%5Cdot%20y%3D2y%24%2C%20%24y(0)%3D1%24%2C%20and%20%24h%3D0.1%24.%20Compare%20with%20the%20exact%20value.%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%203%0A%20%20%20%20Euler%20gives%20%241%2B0.1(2)%3D1.2%24.%20The%20exact%20value%20is%20%24e%5E%7B0.2%7D%5Capprox1.2214%24.%20The%20tangent-line%20step%20underestimates%20this%20convex%20growth%20curve.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Explicit%20and%20Semi-Implicit%20Euler%20Methods%0A%0A%20%20%20%20For%20a%20state%20with%20position%20%24y%24%20and%20velocity%20%24v%24%2C%20the%20update%20order%20matters.%20The%20forward%2C%20or%20explicit%2C%20Euler%20method%20uses%20only%20the%20current%20state%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20y_%7Bi%2B1%7D%3Dy_i%2Bh%20v_i%2C%20%5Cqquad%20v_%7Bi%2B1%7D%3Dv_i%2Bh%20a(t_i%2C%20y_i%2C%20v_i).%0A%20%20%20%20%24%24%0A%0A%20%20%20%20A%20common%20alternative%20in%20mechanics%20is%20the%20semi-implicit%20Euler%20method%2C%20also%20called%20symplectic%20Euler%20or%20Euler-Cromer.%20It%20updates%20velocity%20first%20and%20then%20uses%20the%20new%20velocity%20to%20update%20position%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20v_%7Bi%2B1%7D%3Dv_i%2Bh%20a(t_i%2C%20y_i%2C%20v_i)%2C%20%5Cqquad%20y_%7Bi%2B1%7D%3Dy_i%2Bh%20v_%7Bi%2B1%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20This%20is%20not%20the%20fully%20implicit%20backward%20Euler%20method%2C%20which%20would%20evaluate%20the%20derivative%20at%20the%20future%20state%20and%20usually%20requires%20solving%20an%20algebraic%20equation%20at%20every%20step.%20The%20semi-implicit%20version%20is%20still%20simple%20to%20implement%2C%20and%20for%20many%20mechanical%20systems%20it%20behaves%20better%20than%20forward%20Euler.%0A%0A%20%20%20%20**Guiding%20questions%206.**%0A%20%20%20%201.%20In%20forward%20Euler%2C%20which%20value%20of%20velocity%20updates%20position%3F%0A%20%20%20%202.%20In%20semi-implicit%20Euler%2C%20which%20value%20of%20velocity%20updates%20position%3F%0A%20%20%20%203.%20Why%20might%20the%20second%20choice%20matter%20for%20oscillating%20systems%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%204.**%20With%20%24y_0%3D1%24%20m%2C%20%24v_0%3D2%24%20m%2Fs%2C%20%24a%3D-10%24%20m%2Fs%C2%B2%2C%20and%20%24h%3D0.1%24%20s%2C%20compare%20one%20explicit%20and%20one%20semi-implicit%20Euler%20step.%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%204%0A%20%20%20%20Both%20give%20%24v_1%3D1%24%20m%2Fs.%20Explicit%20Euler%20uses%20the%20old%20velocity%3A%20%24y_1%3D1.2%24%20m.%20Semi-implicit%20Euler%20uses%20the%20new%20velocity%3A%20%24y_1%3D1.1%24%20m.%20The%20exact%20constant-acceleration%20value%20is%201.15%20m%3B%20update%20order%20matters.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Other%20numerical%20methods%20for%20solving%20ODEs%0A%0A%20%20%20%20More%20accurate%20methods%20use%20additional%20evaluations%20of%20the%20derivative%20inside%20each%20interval%20%24%5Bt_n%2C%20t_%7Bn%2B1%7D%5D%24.%20A%20common%20family%20is%20the%20%5BRunge-Kutta%20methods%5D(https%3A%2F%2Fen.wikipedia.org%2Fwiki%2FRunge%25E2%2580%2593Kutta_methods).%0A%0A%20%20%20%20In%20SciPy%2C%20the%20current%20general-purpose%20interface%20for%20initial%20value%20problems%20is%20%5B%60scipy.integrate.solve_ivp%60%5D(https%3A%2F%2Fdocs.scipy.org%2Fdoc%2Fscipy%2Freference%2Fgenerated%2Fscipy.integrate.solve_ivp.html).%20SciPy%20recommends%20%60solve_ivp%60%20for%20new%20code.%20It%20provides%20explicit%20Runge-Kutta%20methods%20such%20as%20%60RK45%60%20and%20%60DOP853%60%2C%20implicit%20methods%20for%20stiff%20systems%20such%20as%20%60Radau%60%20and%20%60BDF%60%2C%20and%20%60LSODA%60%2C%20which%20automatically%20switches%20between%20nonstiff%20and%20stiff%20methods.%0A%0A%20%20%20%20**Challenge%207.**%20When%20you%20later%20compare%20Euler%20with%20%60solve_ivp%60%2C%20predict%20which%20method%20will%20be%20more%20accurate%20for%20the%20same%20output%20time%20step%20and%20explain%20why.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Examples%0A%0A%20%20%20%20Now%20you%20will%20work%20through%20examples.%20Each%20example%20starts%20with%20a%20model%2C%20rewrites%20it%20as%20an%20IVP%2C%20and%20then%20asks%20you%20to%20inspect%20or%20modify%20the%20numerical%20solution.%0A%0A%20%20%20%20First%2C%20import%20the%20required%20Python%20libraries%20and%20customize%20the%20plotting%20environment.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Example%3A%20Simple%20Pendulum%20under%20Gravity%0A%0A%20%20%20%20A%20simple%20pendulum%20is%20a%20compact%20example%20of%20an%20IVP.%20Let%20%24%5Ctheta%24%20be%20the%20angular%20position%20and%20%24%5Comega%3D%5Cdot%7B%5Ctheta%7D%24%20be%20the%20angular%20velocity.%20For%20a%20pendulum%20of%20length%20%24d%24%20under%20gravity%20%24g%24%2C%20the%20equation%20of%20motion%20is%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cddot%7B%5Ctheta%7D(t)%20%3D%20-%5Cfrac%7Bg%7D%7Bd%7D%5Csin(%5Ctheta(t)).%0A%20%20%20%20%24%24%0A%0A%20%20%20%20With%20initial%20conditions%20%24%5Ctheta(0)%3D45%5E%5Ccirc%24%20and%20%24%5Comega(0)%3D0%24%2C%20the%20second-order%20IVP%20can%20be%20written%20as%20two%20first-order%20equations%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cdot%7B%5Ctheta%7D%20%3D%20%5Comega%2C%20%5Cqquad%0A%20%20%20%20%5Cdot%7B%5Comega%7D%20%3D%20-%5Cfrac%7Bg%7D%7Bd%7D%5Csin(%5Ctheta).%0A%20%20%20%20%24%24%0A%0A%20%20%20%20We%20will%20solve%20this%20system%20with%20the%20semi-implicit%20Euler%20method.%0A%0A%20%20%20%20**Before%20you%20run%20the%20next%20cell**%2C%20predict%20whether%20the%20angular%20velocity%20will%20be%20largest%20near%20the%20top%20of%20the%20swing%20or%20near%20the%20bottom.%20Then%20compare%20your%20prediction%20with%20the%20plot.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.function%0Adef%20euler_method(acceleration%2C%20T%3D10%2C%20y0%3D(0.0%2C%200.0)%2C%20h%3D0.01%2C%20method%3D2)%3A%0A%20%20%20%20%22%22%22%0A%20%20%20%20First-order%20numerical%20approximation%20for%20two%20coupled%20first-order%20ODEs.%0A%0A%20%20%20%20The%20state%20is%20%5Bposition%2C%20velocity%5D.%20The%20acceleration%20function%20has%20the%0A%20%20%20%20signature%20acceleration%20(t%2C%20state)%20and%20returns%20dv%2Fdt.%0A%0A%20%20%20%20The%20Euler%20method%20can%20be%201%20(explicit)%2C%20or%202%20(semi-implicit).%0A%20%20%20%20%22%22%22%0A%20%20%20%20import%20numpy%20as%20np%0A%0A%20%20%20%20n_samples%20%3D%20int(np.ceil(T%20%2F%20h))%0A%20%20%20%20y%20%3D%20np.zeros((2%2C%20n_samples))%0A%20%20%20%20y%5B%3A%2C%200%5D%20%3D%20np.asarray(y0%2C%20dtype%3Dfloat)%0A%20%20%20%20t%20%3D%20h%20*%20np.arange(n_samples)%0A%0A%20%20%20%20for%20i%20in%20range(n_samples%20-%201)%3A%0A%20%20%20%20%20%20%20%20if%20method%20%3D%3D%201%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20y%5B0%2C%20i%20%2B%201%5D%20%3D%20y%5B0%2C%20i%5D%20%2B%20h%20*%20y%5B1%2C%20i%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20y%5B1%2C%20i%20%2B%201%5D%20%3D%20y%5B1%2C%20i%5D%20%2B%20h%20*%20acceleration(t%5Bi%5D%2C%20y%5B%3A%2C%20i%5D)%0A%20%20%20%20%20%20%20%20elif%20method%20%3D%3D%202%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20y%5B1%2C%20i%20%2B%201%5D%20%3D%20y%5B1%2C%20i%5D%20%2B%20h%20*%20acceleration(t%5Bi%5D%2C%20y%5B%3A%2C%20i%5D)%0A%20%20%20%20%20%20%20%20%20%20%20%20y%5B0%2C%20i%20%2B%201%5D%20%3D%20y%5B0%2C%20i%5D%20%2B%20h%20*%20y%5B1%2C%20i%20%2B%201%5D%0A%20%20%20%20%20%20%20%20else%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20raise%20ValueError(%22Valid%20options%20for%20method%20are%201%20or%202.%22)%0A%0A%20%20%20%20return%20t%2C%20y%0A%0A%0A%40app.function%0Adef%20pendulum_acceleration(t%2C%20state)%3A%0A%20%20%20%20%22%22%22Return%20angular%20acceleration%20for%20a%20simple%20pendulum.%22%22%22%0A%20%20%20%20import%20numpy%20as%20np%0A%0A%20%20%20%20del%20t%0A%20%20%20%20length%20%3D%200.5%0A%20%20%20%20gravity%20%3D%209.8%0A%20%20%20%20return%20-(gravity%20%2F%20length)%20*%20np.sin(state%5B0%5D)%0A%0A%0A%40app.function%0Adef%20plot_pendulum(t%2C%20y%2C%20labels)%3A%0A%20%20%20%20%22%22%22Plot%20angular%20position%20and%20velocity%20for%20the%20pendulum%20example.%22%22%22%0A%20%20%20%20import%20matplotlib.pyplot%20as%20plt%0A%0A%20%20%20%20_%2C%20ax1%20%3D%20plt.subplots(1%2C%201%2C%20figsize%3D(10%2C%204))%0A%20%20%20%20ax1.set_title(labels%5B0%5D)%0A%20%20%20%20ax1.plot(t%2C%20y%5B0%2C%20%3A%5D%2C%20%22b%22%2C%20label%3D%22Angular%20position%22)%0A%20%20%20%20ax1.set_xlabel(%22Time%20%5Bs%5D%22)%0A%20%20%20%20ax1.set_ylabel(labels%5B1%5D%2C%20color%3D%22b%22)%0A%20%20%20%20ax1.tick_params(%22y%22%2C%20colors%3D%22b%22)%0A%0A%20%20%20%20ax2%20%3D%20ax1.twinx()%0A%20%20%20%20ax2.plot(t%2C%20y%5B1%2C%20%3A%5D%2C%20%22r-.%22%2C%20label%3D%22Angular%20velocity%22)%0A%20%20%20%20ax2.set_ylabel(labels%5B2%5D%2C%20color%3D%22r%22)%0A%20%20%20%20ax2.tick_params(%22y%22%2C%20colors%3D%22r%22)%0A%20%20%20%20plt.tight_layout()%0A%20%20%20%20plt.show()%0A%0A%0A%40app.cell%0Adef%20_(np)%3A%0A%20%20%20%20T_pendulum%20%3D%2010%0A%20%20%20%20y0_pendulum%20%3D%20%5B45%20*%20np.pi%20%2F%20180%2C%200%5D%0A%20%20%20%20h_pendulum%20%3D%200.01%0A%20%20%20%20t_pendulum%2C%20theta%20%3D%20euler_method(%0A%20%20%20%20%20%20%20%20pendulum_acceleration%2C%0A%20%20%20%20%20%20%20%20T%3DT_pendulum%2C%0A%20%20%20%20%20%20%20%20y0%3Dy0_pendulum%2C%0A%20%20%20%20%20%20%20%20h%3Dh_pendulum%2C%0A%20%20%20%20%20%20%20%20method%3D2%2C%0A%20%20%20%20)%0A%20%20%20%20pendulum_labels%20%3D%20%5B%0A%20%20%20%20%20%20%20%20%22Trajectory%20of%20a%20simple%20pendulum%20under%20gravity%22%2C%0A%20%20%20%20%20%20%20%20%22Angular%20position%20%5Bdeg%5D%22%2C%0A%20%20%20%20%20%20%20%20%22Angular%20velocity%20%5Bdeg%2Fs%5D%22%2C%0A%20%20%20%20%5D%0A%20%20%20%20return%20pendulum_labels%2C%20t_pendulum%2C%20theta%0A%0A%0A%40app.cell%0Adef%20_(np%2C%20pendulum_labels%2C%20t_pendulum%2C%20theta)%3A%0A%20%20%20%20plot_pendulum(t_pendulum%2C%20np.rad2deg(theta)%2C%20pendulum_labels)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Challenge%208.**%20Change%20the%20initial%20angle%20to%20%2410%5E%5Ccirc%24%2C%20%2445%5E%5Ccirc%24%2C%20and%20%2490%5E%5Ccirc%24.%20Does%20the%20plot%20look%20like%20a%20perfect%20sine%20wave%20in%20all%20three%20cases%3F%20What%20does%20that%20tell%20you%20about%20the%20small-angle%20approximation%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%205.**%20For%20a%201%20m%20pendulum%20at%20%24%5Ctheta%3D0%24%2C%20where%20is%20speed%20largest%20during%20a%20swing%20released%20from%20rest%3F%20Why%20is%20angle%20in%20radians%20essential%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%205%0A%20%20%20%20At%20the%20bottom%2C%20gravitational%20potential%20energy%20has%20been%20converted%20into%20kinetic%20energy%2C%20so%20speed%20is%20largest.%20The%20derivatives%20and%20the%20sine%20in%20the%20equation%20use%20radians%3B%20substituting%20a%20degree%20value%20directly%20changes%20the%20model.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Motion%20under%20Constant%20Force%0A%0A%20%20%20%20Consider%20a%20football%20kicked%20upward%20from%20an%20initial%20height%20%24y_0%24%20with%20initial%20vertical%20velocity%20%24v_0%24.%20We%20want%20the%20equation%20of%20motion%20in%20the%20vertical%20direction.%0A%0A%20%20%20%20Neglecting%20air%20resistance%2C%20Newton's%20second%20law%20gives%0A%0A%20%20%20%20%24%24%0A%20%20%20%20m%5Cfrac%7B%5Cmathrm%7Bd%7D%5E2%20y%7D%7B%5Cmathrm%7Bd%7Dt%5E2%7D%20%3D%20-mg.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20We%20will%20use%20%24g%3D9.8%5C%2C%5Cmathrm%7Bm%2Fs%5E2%7D%24%2C%20%24y_0%3D1%5C%2C%5Cmathrm%7Bm%7D%24%20at%20%24t_0%3D0%24%2C%20and%20%24v_0%3D20%5C%2C%5Cmathrm%7Bm%2Fs%7D%24.%20The%20analytical%20solution%20is%0A%0A%20%20%20%20%24%24%0A%20%20%20%20y(t)%20%3D%20y_0%20%2B%20v_0%20t%20-%20%5Cfrac%7Bg%7D%7B2%7Dt%5E2.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20To%20solve%20the%20same%20problem%20numerically%2C%20rewrite%20the%20second-order%20ODE%20as%20two%20first-order%20ODEs%3A%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cdot%7By%7D%20%3D%20v%2C%20%5Cqquad%20%5Cdot%7Bv%7D%20%3D%20a.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20For%20constant%20gravitational%20acceleration%2C%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cleft%5C%7B%0A%20%20%20%20%5Cbegin%7Barray%7D%7Brl%7D%0A%20%20%20%20%5Cfrac%7B%5Cmathrm%7Bd%7D%20y%7D%7B%5Cmathrm%7Bd%7Dt%7D%20%26%3D%20v%2C%20%5Cquad%20y(t_0)%20%3D%20y_0%2C%5C%5C%0A%20%20%20%20%5Cfrac%7B%5Cmathrm%7Bd%7D%20v%7D%7B%5Cmathrm%7Bd%7Dt%7D%20%26%3D%20-g%2C%20%5Cquad%20v(t_0)%20%3D%20v_0.%0A%20%20%20%20%5Cend%7Barray%7D%0A%20%20%20%20%5Cright.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20**Guiding%20questions%205.**%0A%20%20%20%201.%20Which%20two%20variables%20define%20the%20state%20of%20the%20ball%3F%0A%20%20%20%202.%20Which%20equation%20updates%20position%3F%0A%20%20%20%203.%20Which%20equation%20updates%20velocity%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20This%20function%20calculates%20the%20vertical%20trajectory%20with%20the%20Euler%20method.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.function%0Adef%20ball_euler(t0%2C%20tend%2C%20y0%2C%20v0%2C%20h)%3A%0A%20%20%20%20%22%22%22Integrate%20vertical%20ball%20motion%20with%20the%20explicit%20Euler%20method.%22%22%22%0A%20%20%20%20import%20numpy%20as%20np%0A%0A%20%20%20%20t%20%3D%20%5Bt0%5D%0A%20%20%20%20y%20%3D%20%5By0%5D%0A%20%20%20%20v%20%3D%20%5Bv0%5D%0A%20%20%20%20a%20%3D%20-9.8%0A%20%20%20%20while%20t%5B-1%5D%20%3C%3D%20tend%20and%20y%5B-1%5D%20%3E%200%3A%0A%20%20%20%20%20%20%20%20y.append(y%5B-1%5D%20%2B%20h%20*%20v%5B-1%5D)%0A%20%20%20%20%20%20%20%20v.append(v%5B-1%5D%20%2B%20h%20*%20a)%0A%20%20%20%20%20%20%20%20t.append(t%5B-1%5D%20%2B%20h)%0A%20%20%20%20return%20np.array(t)%2C%20np.array(y)%2C%20np.array(v)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20Set%20the%20initial%20conditions.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20y0%20%3D%201%0A%20%20%20%20v0%20%3D%2020%0A%20%20%20%20return%20v0%2C%20y0%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20Now%20integrate%20the%20same%20motion%20with%20two%20step%20sizes.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(v0%2C%20y0)%3A%0A%20%20%20%20t100%2C%20y100%2C%20v100%20%3D%20ball_euler(0%2C%204%2C%20y0%2C%20v0%2C%200.1)%0A%20%20%20%20t10%2C%20y10%2C%20v10%20%3D%20ball_euler(0%2C%204%2C%20y0%2C%20v0%2C%200.01)%0A%20%20%20%20return%20t10%2C%20t100%2C%20v10%2C%20v100%2C%20y10%2C%20y100%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20Here%20are%20the%20numerical%20results%20compared%20with%20the%20analytical%20solution.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(plt%2C%20v0%2C%20y0)%3A%0A%20%20%20%20def%20plots(t100%2C%20y100%2C%20v100%2C%20t10%2C%20y10%2C%20v10%2C%20title)%3A%0A%20%20%20%20%20%20%20%20%22%22%22Plot%20numerical%20integration%20results%20against%20the%20analytical%20solution.%22%22%22%0A%20%20%20%20%20%20%20%20a%20%3D%20-9.8%0A%0A%20%20%20%20%20%20%20%20_%2C%20axs%20%3D%20plt.subplots(nrows%3D2%2C%20ncols%3D2%2C%20sharex%3DTrue%2C%20figsize%3D(10%2C%205))%0A%0A%20%20%20%20%20%20%20%20axs%5B0%2C%200%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t10%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20y0%20%2B%20v0%20*%20t10%20%2B%200.5%20*%20a%20*%20t10**2%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B0%2C%200%2C%201%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22Analytical%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20axs%5B0%2C%200%5D.plot(t100%2C%20y100%2C%20%22--%22%2C%20color%3D%5B0%2C%201%2C%200%2C%200.7%5D%2C%20label%3D%22h%20%3D%20100%20ms%22)%0A%20%20%20%20%20%20%20%20axs%5B0%2C%200%5D.plot(t10%2C%20y10%2C%20%22%3A%22%2C%20color%3D%5B1%2C%200%2C%200%2C%200.7%5D%2C%20label%3D%22h%20%3D%2010%20ms%22)%0A%0A%20%20%20%20%20%20%20%20axs%5B0%2C%201%5D.plot(t10%2C%20v0%20%2B%20a%20*%20t10%2C%20color%3D%5B0%2C%200%2C%201%2C%200.5%5D%2C%20label%3D%22Analytical%22)%0A%20%20%20%20%20%20%20%20axs%5B0%2C%201%5D.plot(t100%2C%20v100%2C%20%22--%22%2C%20color%3D%5B0%2C%201%2C%200%2C%200.7%5D%2C%20label%3D%22h%20%3D%20100%20ms%22)%0A%20%20%20%20%20%20%20%20axs%5B0%2C%201%5D.plot(t10%2C%20v10%2C%20%22%3A%22%2C%20color%3D%5B1%2C%200%2C%200%2C%200.7%5D%2C%20label%3D%22h%20%3D%2010%20ms%22)%0A%0A%20%20%20%20%20%20%20%20axs%5B1%2C%200%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t10%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20y0%20%2B%20v0%20*%20t10%20%2B%200.5%20*%20a%20*%20t10**2%20-%20(y0%20%2B%20v0%20*%20t10%20%2B%200.5%20*%20a%20*%20t10**2)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B0%2C%200%2C%201%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22Analytical%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%200%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t100%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20y100%20-%20(y0%20%2B%20v0%20*%20t100%20%2B%200.5%20*%20a%20*%20t100**2)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%22--%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B0%2C%201%2C%200%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22h%20%3D%20100%20ms%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%200%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t10%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20y10%20-%20(y0%20%2B%20v0%20*%20t10%20%2B%200.5%20*%20a%20*%20t10**2)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%22%3A%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B1%2C%200%2C%200%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22h%20%3D%2010%20ms%22%2C%0A%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20axs%5B1%2C%201%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t10%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20v0%20%2B%20a%20*%20t10%20-%20(v0%20%2B%20a%20*%20t10)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B0%2C%200%2C%201%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22Analytical%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%201%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t100%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20v100%20-%20(v0%20%2B%20a%20*%20t100)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%22--%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B0%2C%201%2C%200%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22h%20%3D%20100%20ms%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%201%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t10%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20v10%20-%20(v0%20%2B%20a%20*%20t10)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%22%3A%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B1%2C%200%2C%200%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22h%20%3D%2010%20ms%22%2C%0A%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20ylabel%20%3D%20%5B%22y%20%5Bm%5D%22%2C%20%22v%20%5Bm%2Fs%5D%22%2C%20%22y%20error%20%5Bm%5D%22%2C%20%22v%20error%20%5Bm%2Fs%5D%22%5D%0A%20%20%20%20%20%20%20%20axs%5B0%2C%200%5D.set_xlim(t10%5B0%5D%2C%20t10%5B-1%5D)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%200%5D.set_xlabel(%22Time%20%5Bs%5D%22)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%201%5D.set_xlabel(%22Time%20%5Bs%5D%22)%0A%20%20%20%20%20%20%20%20axs%5B0%2C%201%5D.legend()%0A%20%20%20%20%20%20%20%20axs%20%3D%20axs.flatten()%0A%20%20%20%20%20%20%20%20for%20i%2C%20ax%20in%20enumerate(axs)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20ax.set_ylabel(ylabel%5Bi%5D)%0A%20%20%20%20%20%20%20%20plt.suptitle(f%22Kinematics%20of%20a%20football%20-%20%7Btitle%7D%20method%22%2C%20y%3D1.02%2C%20fontsize%3D16)%0A%20%20%20%20%20%20%20%20plt.tight_layout()%0A%20%20%20%20%20%20%20%20plt.show()%0A%0A%20%20%20%20return%20(plots%2C)%0A%0A%0A%40app.cell%0Adef%20_(plots%2C%20t10%2C%20t100%2C%20v10%2C%20v100%2C%20y10%2C%20y100)%3A%0A%20%20%20%20plots(t100%2C%20y100%2C%20v100%2C%20t10%2C%20y10%2C%20v10%2C%20%22Euler%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Challenge%209.**%20Compare%20the%20error%20curves%20for%20%24h%3D100%24%20ms%20and%20%24h%3D10%24%20ms.%20Which%20error%20changes%20more%20when%20the%20step%20size%20decreases%3A%20position%20or%20velocity%3F%20Explain%20your%20answer%20from%20the%20update%20equations.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%60solve_ivp%60%20can%20call%20the%20LSODA%20method%20for%20the%20same%20problem.%20LSODA%20is%20useful%20as%20a%20robust%20general%20solver%20because%20it%20switches%20automatically%20between%20methods%20for%20nonstiff%20and%20stiff%20systems.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.function%0Adef%20ball_constant_force(t%2C%20yv)%3A%0A%20%20%20%20y%2C%20v%20%3D%20yv%0A%20%20%20%20a%20%3D%20-9.8%0A%20%20%20%20return%20%5Bv%2C%20a%5D%0A%0A%0A%40app.function%0Adef%20solve_trajectory(fun%2C%20t0%2C%20tend%2C%20yv0%2C%20h%2C%20method)%3A%0A%20%20%20%20import%20numpy%20as%20np%0A%20%20%20%20from%20scipy.integrate%20import%20solve_ivp%0A%0A%20%20%20%20t_eval%20%3D%20np.arange(t0%2C%20tend%20%2B%20h%20%2F%202%2C%20h)%0A%20%20%20%20solution%20%3D%20solve_ivp(%0A%20%20%20%20%20%20%20%20fun%3Dfun%2C%0A%20%20%20%20%20%20%20%20t_span%3D(t0%2C%20tend)%2C%0A%20%20%20%20%20%20%20%20y0%3Dyv0%2C%0A%20%20%20%20%20%20%20%20method%3Dmethod%2C%0A%20%20%20%20%20%20%20%20t_eval%3Dt_eval%2C%0A%20%20%20%20)%0A%20%20%20%20if%20not%20solution.success%3A%0A%20%20%20%20%20%20%20%20raise%20RuntimeError(solution.message)%0A%20%20%20%20return%20solution.t%2C%20solution.y%5B0%5D%2C%20solution.y%5B1%5D%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20t100_1%2C%20y100_1%2C%20v100_1%20%3D%20solve_trajectory(%0A%20%20%20%20%20%20%20%20ball_constant_force%2C%200%2C%204%2C%20%5B1%2C%2020%5D%2C%200.1%2C%20method%3D%22LSODA%22%0A%20%20%20%20)%0A%20%20%20%20t10_1%2C%20y10_1%2C%20v10_1%20%3D%20solve_trajectory(%0A%20%20%20%20%20%20%20%20ball_constant_force%2C%200%2C%204%2C%20%5B1%2C%2020%5D%2C%200.01%2C%20method%3D%22LSODA%22%0A%20%20%20%20)%0A%20%20%20%20return%20t100_1%2C%20t10_1%2C%20v100_1%2C%20v10_1%2C%20y100_1%2C%20y10_1%0A%0A%0A%40app.cell%0Adef%20_(plots%2C%20t100_1%2C%20t10_1%2C%20v100_1%2C%20v10_1%2C%20y100_1%2C%20y10_1)%3A%0A%20%20%20%20plots(t100_1%2C%20y100_1%2C%20v100_1%2C%20t10_1%2C%20y10_1%2C%20v10_1%2C%20%22LSODA%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%60solve_ivp%60%20also%20offers%20explicit%20Runge-Kutta%20methods.%20Its%20default%20method%2C%20%60RK45%60%2C%20is%20an%20explicit%20Dormand-Prince%20method%20of%20order%205(4)%2C%20similar%20in%20spirit%20to%20MATLAB's%20%60ode45%60.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20t100_2%2C%20y100_2%2C%20v100_2%20%3D%20solve_trajectory(%0A%20%20%20%20%20%20%20%20ball_constant_force%2C%200%2C%204%2C%20%5B1%2C%2020%5D%2C%200.1%2C%20method%3D%22RK45%22%0A%20%20%20%20)%0A%20%20%20%20t10_2%2C%20y10_2%2C%20v10_2%20%3D%20solve_trajectory(%0A%20%20%20%20%20%20%20%20ball_constant_force%2C%200%2C%204%2C%20%5B1%2C%2020%5D%2C%200.01%2C%20method%3D%22RK45%22%0A%20%20%20%20)%0A%20%20%20%20return%20t100_2%2C%20t10_2%2C%20v100_2%2C%20v10_2%2C%20y100_2%2C%20y10_2%0A%0A%0A%40app.cell%0Adef%20_(plots%2C%20t100_2%2C%20t10_2%2C%20v100_2%2C%20v10_2%2C%20y100_2%2C%20y10_2)%3A%0A%20%20%20%20plots(t100_2%2C%20y100_2%2C%20v100_2%2C%20t10_2%2C%20y10_2%2C%20v10_2%2C%20%22RK45%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Challenge%2010.**%20Compare%20the%20Euler%2C%20LSODA%2C%20and%20RK45%20plots.%20If%20you%20had%20to%20choose%20a%20method%20for%20a%20new%20biomechanical%20simulation%2C%20what%20evidence%20in%20these%20plots%20would%20influence%20your%20choice%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Motion%20under%20varying%20force%0A%0A%20%20%20%20Now%20include%20air%20resistance%20in%20the%20vertical%20trajectory%20of%20the%20football.%20Current%20IFAB%20Law%202%20specifies%20that%20a%20football%20must%20be%20spherical%2C%20have%20a%20circumference%20between%2068%20cm%20and%2070%20cm%2C%20and%20weigh%20between%20410%20g%20and%20450%20g%20at%20the%20start%20of%20the%20match.%20The%20representative%20values%20used%20here%2C%20circumference%20%240.69%5C%2C%5Cmathrm%7Bm%7D%24%20and%20mass%20%240.43%5C%2C%5Cmathrm%7Bkg%7D%24%2C%20are%20the%20midpoints%20of%20those%20ranges.%0A%0A%20%20%20%20Model%20the%20magnitude%20of%20the%20drag%20force%20as%0A%0A%20%20%20%20%24%24%0A%20%20%20%20F_d(v)%20%3D%20%5Cfrac%7B1%7D%7B2%7D%5Crho%20C_d%20A%20v%5E2%2C%0A%20%20%20%20%24%24%0A%0A%20%20%20%20where%20%24%5Crho%3D1.22%5C%2C%5Cmathrm%7Bkg%2Fm%5E3%7D%24%20is%20air%20density%2C%20%24A%3D0.0379%5C%2C%5Cmathrm%7Bm%5E2%7D%24%20is%20the%20ball%20cross-sectional%20area%2C%20and%20%24C_d%3D0.25%24%20is%20an%20approximate%20constant%20drag%20coefficient%20from%20Bray%20and%20Kerwin%20(2003).%20Drag%20acts%20opposite%20to%20velocity%2C%20so%20the%20vertical%20equation%20is%0A%0A%20%20%20%20%24%24%0A%20%20%20%20m%5Cfrac%7B%5Cmathrm%7Bd%7D%5E2%20y%7D%7B%5Cmathrm%7Bd%7Dt%5E2%7D%0A%20%20%20%20%3D%0A%20%20%20%20-mg%20-%20%5Cfrac%7B1%7D%7B2%7D%5Crho%20C_d%20A%20v%7Cv%7C.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20As%20a%20first-order%20system%2C%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cleft%5C%7B%0A%20%20%20%20%5Cbegin%7Barray%7D%7Brl%7D%0A%20%20%20%20%5Cfrac%7B%5Cmathrm%7Bd%7D%20y%7D%7B%5Cmathrm%7Bd%7Dt%7D%20%26%3D%20v%2C%20%5Cquad%20y(t_0)%20%3D%20y_0%2C%5C%5C%0A%20%20%20%20%5Cfrac%7B%5Cmathrm%7Bd%7D%20v%7D%7B%5Cmathrm%7Bd%7Dt%7D%20%26%3D%0A%20%20%20%20-g%20-%20%5Cfrac%7B1%7D%7B2m%7D%5Crho%20C_d%20A%20v%7Cv%7C%2C%20%5Cquad%20v(t_0)%20%3D%20v_0.%0A%20%20%20%20%5Cend%7Barray%7D%0A%20%20%20%20%5Cright.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20We%20will%20solve%20this%20nonlinear%20IVP%20numerically%20with%20%60solve_ivp%60.%0A%0A%20%20%20%20**Before%20you%20run%20the%20next%20cell**%2C%20predict%20what%20air%20resistance%20will%20do%20to%20the%20maximum%20height%20and%20time%20of%20flight.%20Will%20the%20ball%20spend%20more%20time%20rising%2C%20more%20time%20falling%2C%20or%20both%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.function%0Adef%20ball_with_drag(t%2C%20yv)%3A%0A%20%20%20%20import%20numpy%20as%20np%0A%0A%20%20%20%20y%2C%20v%20%3D%20yv%0A%20%20%20%20g%20%3D%209.8%0A%20%20%20%20m%20%3D%200.43%0A%20%20%20%20rho%20%3D%201.22%0A%20%20%20%20cd%20%3D%200.25%0A%20%20%20%20area%20%3D%200.0379%0A%20%20%20%20a%20%3D%20-g%20-%20rho%20*%20cd%20*%20area%20*%20v%20*%20np.abs(v)%20%2F%20(2%20*%20m)%0A%20%20%20%20return%20%5Bv%2C%20a%5D%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20t10_3%2C%20y10_3%2C%20v10_3%20%3D%20solve_trajectory(%0A%20%20%20%20%20%20%20%20ball_with_drag%2C%200%2C%204%2C%20%5B1%2C%2020%5D%2C%200.01%2C%20method%3D%22LSODA%22%0A%20%20%20%20)%0A%20%20%20%20return%20t10_3%2C%20v10_3%2C%20y10_3%0A%0A%0A%40app.cell%0Adef%20_(plt%2C%20v0%2C%20y0)%3A%0A%20%20%20%20def%20plots_1(t10%2C%20y10%2C%20v10)%3A%0A%20%20%20%20%20%20%20%20%22%22%22Plot%20vertical%20ball%20motion%20with%20and%20without%20air%20resistance.%22%22%22%0A%20%20%20%20%20%20%20%20a%20%3D%20-9.8%0A%20%20%20%20%20%20%20%20_%2C%20axs%20%3D%20plt.subplots(nrows%3D2%2C%20ncols%3D2%2C%20sharex%3DTrue%2C%20figsize%3D(10%2C%205))%0A%20%20%20%20%20%20%20%20axs%5B0%2C%200%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t10%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20y0%20%2B%20v0%20*%20t10%20%2B%200.5%20*%20a%20*%20t10**2%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B0%2C%200%2C%201%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22No%20resistance%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20axs%5B0%2C%200%5D.plot(t10%2C%20y10%2C%20%22-%22%2C%20color%3D%5B1%2C%200%2C%200%2C%200.7%5D%2C%20label%3D%22With%20resistance%22)%0A%20%20%20%20%20%20%20%20axs%5B0%2C%201%5D.plot(t10%2C%20v0%20%2B%20a%20*%20t10%2C%20color%3D%5B0%2C%200%2C%201%2C%200.7%5D%2C%20label%3D%22No%20resistance%22)%0A%20%20%20%20%20%20%20%20axs%5B0%2C%201%5D.plot(t10%2C%20v10%2C%20%22-%22%2C%20color%3D%5B1%2C%200%2C%200%2C%200.7%5D%2C%20label%3D%22With%20resistance%22)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%200%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t10%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20y0%20%2B%20v0%20*%20t10%20%2B%200.5%20*%20a%20*%20t10**2%20-%20(y0%20%2B%20v0%20*%20t10%20%2B%200.5%20*%20a%20*%20t10**2)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B0%2C%200%2C%201%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22No%20resistance%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%200%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t10%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20y10%20-%20(y0%20%2B%20v0%20*%20t10%20%2B%200.5%20*%20a%20*%20t10**2)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%22-%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B1%2C%200%2C%200%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22With%20resistance%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%201%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t10%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20v0%20%2B%20a%20*%20t10%20-%20(v0%20%2B%20a%20*%20t10)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B0%2C%200%2C%201%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22No%20resistance%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%201%5D.plot(%0A%20%20%20%20%20%20%20%20%20%20%20%20t10%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20v10%20-%20(v0%20%2B%20a%20*%20t10)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%22-%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%5B1%2C%200%2C%200%2C%200.7%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20label%3D%22With%20resistance%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20ylabel%20%3D%20%5B%22y%20%5Bm%5D%22%2C%20%22v%20%5Bm%2Fs%5D%22%2C%20%22y%20difference%20%5Bm%5D%22%2C%20%22v%20difference%20%5Bm%2Fs%5D%22%5D%0A%20%20%20%20%20%20%20%20axs%5B1%2C%200%5D.set_xlabel(%22Time%20%5Bs%5D%22)%0A%20%20%20%20%20%20%20%20axs%5B1%2C%201%5D.set_xlabel(%22Time%20%5Bs%5D%22)%0A%20%20%20%20%20%20%20%20axs%5B0%2C%201%5D.legend()%0A%20%20%20%20%20%20%20%20axs%20%3D%20axs.flatten()%0A%20%20%20%20%20%20%20%20for%20i%2C%20ax%20in%20enumerate(axs)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20ax.set_ylabel(ylabel%5Bi%5D)%0A%20%20%20%20%20%20%20%20plt.suptitle(%0A%20%20%20%20%20%20%20%20%20%20%20%20%22Kinematics%20of%20a%20football%20-%20effect%20of%20air%20resistance%22%2C%20y%3D1.02%2C%20fontsize%3D16%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20plt.tight_layout()%0A%20%20%20%20%20%20%20%20plt.show()%0A%0A%20%20%20%20return%20(plots_1%2C)%0A%0A%0A%40app.cell%0Adef%20_(plots_1%2C%20t10_3%2C%20v10_3%2C%20y10_3)%3A%0A%20%20%20%20plots_1(t10_3%2C%20y10_3%2C%20v10_3)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%206.**%20For%20downward%20velocity%20%24v%3C0%24%2C%20what%20is%20the%20sign%20of%20%24-k%20v%7Cv%7C%24%3F%20Why%20not%20use%20%24-kv%5E2%24%20in%20both%20directions%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%206%0A%20%20%20%20With%20%24v%3C0%24%2C%20%24v%7Cv%7C%3C0%24%2C%20so%20the%20drag%20term%20is%20positive%20(upward).%20A%20permanently%20negative%20%24-kv%5E2%24%20term%20would%20accelerate%20the%20falling%20body%20downward%20instead%20of%20resisting%20its%20motion.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Example%3A%20Mass-Spring-Damper%20System%0A%0A%20%20%20%20Now%20consider%20a%20system%20made%20of%20a%20mass%2C%20a%20spring%2C%20and%20a%20damper.%20This%20is%20a%20basic%20model%20for%20many%20mechanical%20and%20biomechanical%20systems%3A%20a%20shoe%20midsole%2C%20a%20vehicle%20suspension%2C%20a%20vibrating%20instrument%20string%2C%20or%20a%20simplified%20muscle-tendon%20element.%0A%0A%20%20%20%20Let%20%24x%24%20be%20displacement%20from%20equilibrium%20and%20%24v%3D%5Cdot%7Bx%7D%24%20be%20velocity.%20A%20linear%20mass-spring-damper%20model%20is%0A%0A%20%20%20%20%24%24%0A%20%20%20%20m%5Cddot%7Bx%7D%20%2B%20c%5Cdot%7Bx%7D%20%2B%20kx%20%3D%20F(t)%2C%0A%20%20%20%20%24%24%0A%0A%20%20%20%20where%20%24m%24%20is%20mass%2C%20%24c%24%20is%20damping%2C%20%24k%24%20is%20spring%20stiffness%2C%20and%20%24F(t)%24%20is%20an%20external%20force.%20As%20a%20first-order%20IVP%2C%0A%0A%20%20%20%20%24%24%0A%20%20%20%20%5Cdot%7Bx%7D%3Dv%2C%20%5Cqquad%0A%20%20%20%20%5Cdot%7Bv%7D%3D%5Cfrac%7BF(t)-cv-kx%7D%7Bm%7D.%0A%20%20%20%20%24%24%0A%0A%20%20%20%20**Challenge%2011.**%20Before%20running%20the%20simulation%2C%20decide%20what%20should%20happen%20when%20damping%20increases%3A%20should%20the%20oscillations%20grow%2C%20stay%20the%20same%2C%20or%20decay%20faster%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.function%0Adef%20mass_spring_damper_acceleration(%0A%20%20%20%20t%2C%0A%20%20%20%20state%2C%0A%20%20%20%20mass%3D1.0%2C%0A%20%20%20%20stiffness%3D25.0%2C%0A%20%20%20%20damping%3D1.0%2C%0A%20%20%20%20external_force%3D0.0%2C%0A)%3A%0A%20%20%20%20%22%22%22Return%20acceleration%20for%20a%20linear%20mass-spring-damper%20system.%22%22%22%0A%20%20%20%20del%20t%0A%20%20%20%20position%2C%20velocity%20%3D%20state%0A%20%20%20%20return%20(external_force%20-%20damping%20*%20velocity%20-%20stiffness%20*%20position)%20%2F%20mass%0A%0A%0A%40app.function%0Adef%20simulate_mass_spring_damper(%0A%20%20%20%20mass%3D1.0%2C%0A%20%20%20%20stiffness%3D25.0%2C%0A%20%20%20%20damping%3D1.0%2C%0A%20%20%20%20external_force%3D0.0%2C%0A%20%20%20%20initial_state%3D(0.1%2C%200.0)%2C%0A%20%20%20%20T%3D6.0%2C%0A%20%20%20%20h%3D0.01%2C%0A)%3A%0A%20%20%20%20def%20acceleration(t%2C%20state)%3A%0A%20%20%20%20%20%20%20%20return%20mass_spring_damper_acceleration(%0A%20%20%20%20%20%20%20%20%20%20%20%20t%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20state%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20mass%3Dmass%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20stiffness%3Dstiffness%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20damping%3Ddamping%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20external_force%3Dexternal_force%2C%0A%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20return%20euler_method(acceleration%2C%20T%3DT%2C%20y0%3Dinitial_state%2C%20h%3Dh%2C%20method%3D2)%0A%0A%0A%40app.function%0Adef%20plot_mass_spring_damper(time%2C%20state%2C%20title)%3A%0A%20%20%20%20import%20matplotlib.pyplot%20as%20plt%0A%0A%20%20%20%20_%2C%20ax1%20%3D%20plt.subplots(1%2C%201%2C%20figsize%3D(10%2C%204))%0A%20%20%20%20ax1.set_title(title)%0A%20%20%20%20ax1.plot(time%2C%20state%5B0%5D%2C%20%22b%22%2C%20label%3D%22Displacement%22)%0A%20%20%20%20ax1.set_xlabel(%22Time%20%5Bs%5D%22)%0A%20%20%20%20ax1.set_ylabel(%22Displacement%20%5Bm%5D%22%2C%20color%3D%22b%22)%0A%20%20%20%20ax1.tick_params(%22y%22%2C%20colors%3D%22b%22)%0A%0A%20%20%20%20ax2%20%3D%20ax1.twinx()%0A%20%20%20%20ax2.plot(time%2C%20state%5B1%5D%2C%20%22r-.%22%2C%20label%3D%22Velocity%22)%0A%20%20%20%20ax2.set_ylabel(%22Velocity%20%5Bm%2Fs%5D%22%2C%20color%3D%22r%22)%0A%20%20%20%20ax2.tick_params(%22y%22%2C%20colors%3D%22r%22)%0A%20%20%20%20plt.tight_layout()%0A%20%20%20%20plt.show()%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20mass_spring_time%2C%20mass_spring_state%20%3D%20simulate_mass_spring_damper(%0A%20%20%20%20%20%20%20%20mass%3D1.0%2C%0A%20%20%20%20%20%20%20%20stiffness%3D25.0%2C%0A%20%20%20%20%20%20%20%20damping%3D1.2%2C%0A%20%20%20%20%20%20%20%20initial_state%3D(0.1%2C%200.0)%2C%0A%20%20%20%20%20%20%20%20T%3D6.0%2C%0A%20%20%20%20%20%20%20%20h%3D0.01%2C%0A%20%20%20%20)%0A%20%20%20%20return%20mass_spring_state%2C%20mass_spring_time%0A%0A%0A%40app.cell%0Adef%20_(mass_spring_state%2C%20mass_spring_time)%3A%0A%20%20%20%20plot_mass_spring_damper(%0A%20%20%20%20%20%20%20%20mass_spring_time%2C%0A%20%20%20%20%20%20%20%20mass_spring_state%2C%0A%20%20%20%20%20%20%20%20%22Mass-spring-damper%20response%22%2C%0A%20%20%20%20)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Challenge%2012.**%20Try%20damping%20values%20%24c%3D0%24%2C%20%24c%3D1.2%24%2C%20and%20%24c%3D10%24%20(with%20%24m%3D1%24%20and%20%24k%3D25%24%2C%20so%20the%20natural%20frequency%20is%20%24%5Comega_n%3D5%5C%2C%5Cmathrm%7Brad%2Fs%7D%24).%20Classify%20each%20response%20as%20undamped%2C%20underdamped%2C%20or%20critically%20damped.%20Which%20case%20oscillates%20without%20ever%20settling%3F%20Which%20returns%20to%20equilibrium%20without%20oscillating%3F%20What%20changes%20if%20you%20double%20the%20stiffness%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%207.**%20Classify%20the%20damping%20cases%20%24c%3D0%24%2C%201.2%2C%20and%2010%20for%20%24m%3D1%24%2C%20%24k%3D25%24.%20Does%20arbitrarily%20large%20damping%20always%20return%20the%20system%20faster%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%207%0A%20%20%20%20The%20critical%20coefficient%20is%20%242%5Csqrt%7Bmk%7D%3D10%24.%20The%20cases%20are%20undamped%2C%20underdamped%2C%20and%20critically%20damped.%20Very%20large%20damping%20produces%20a%20slow%20overdamped%20response%2C%20so%20%E2%80%9Cmore%20damping%E2%80%9D%20is%20not%20synonymous%20with%20%E2%80%9Cfaster%20return%E2%80%9D.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Checkpoint%20questions%0A%0A%20%20%20%20Pause%20here%20before%20the%20final%20problems.%20These%20questions%20connect%20the%20examples%20in%20this%20notebook%20with%20mechanics%20examples%20used%20elsewhere%20in%20the%20course.%0A%0A%20%20%20%201.%20In%20the%20football-with-drag%20example%2C%20what%20would%20terminal%20velocity%20mean%3F%20Which%20term%20in%20the%20ODE%20prevents%20the%20downward%20speed%20from%20growing%20without%20bound%3F%0A%20%20%20%202.%20In%20a%20one-mass%20spring%20problem%2C%20why%20is%20the%20displacement%20usually%20measured%20from%20the%20equilibrium%20position%20instead%20of%20from%20the%20unstretched%20spring%20length%3F%0A%20%20%20%203.%20If%20you%20had%20a%20force-platform%20record%20from%20a%20vertical%20jump%2C%20which%20signal%20would%20you%20integrate%20to%20estimate%20velocity%3F%20What%20initial%20condition%20would%20you%20need%3F%0A%20%20%20%204.%20A%20two-mass%20spring-damper%20system%20needs%20more%20state%20variables%20than%20the%20one-mass%20example.%20What%20are%20the%20minimum%20state%20variables%20if%20the%20two%20masses%20move%20only%20along%20one%20horizontal%20line%3F%0A%20%20%20%205.%20In%20a%20simple%20muscle-tendon%20model%2C%20which%20elements%20would%20behave%20like%20springs%2C%20which%20element%20would%20behave%20like%20a%20damper%2C%20and%20which%20element%20would%20represent%20active%20force%20production%3F%0A%20%20%20%206.%20When%20would%20you%20trust%20a%20coarse%20Euler%20solution%20for%20teaching%20or%20estimation%2C%20and%20when%20would%20you%20switch%20to%20%60solve_ivp%60%20with%20tighter%20tolerances%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Checkpoint%20answers%0A%20%20%20%201.%20Terminal%20velocity%20is%20a%20downward%20steady%20speed%20where%20drag%20balances%20weight%3B%20the%20quadratic%20term%20limits%20acceleration.%0A%20%20%20%202.%20The%20static%20preload%20can%20then%20be%20absorbed%20into%20the%20equilibrium%20definition.%0A%20%20%20%203.%20Integrate%20%24(F_z-mg)%2Fm%24%2C%20starting%20from%20a%20known%20center-of-mass%20velocity.%0A%20%20%20%204.%20Two%20positions%20and%20two%20velocities.%0A%20%20%20%205.%20Tendon%20and%20passive%20elastic%20tissue%20are%20spring-like%3B%20viscous%20tissue%20is%20damper-like%3B%20the%20contractile%20element%20supplies%20active%20force.%0A%20%20%20%206.%20Use%20step%20refinement%20and%20comparison%20with%20a%20known%20solution%20or%20tighter%20reference.%20Switch%20when%20numerical%20error%20matters%20for%20the%20question%3B%20choose%20tolerances%20against%20a%20physical%20accuracy%20requirement.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Problems%0A%0A%20%20%20%20Work%20through%20these%20in%20order%20if%20you%20are%20studying%20independently.%20For%20each%20one%2C%20first%20write%20the%20IVP%20on%20paper%2C%20then%20modify%20or%20reuse%20code%20from%20this%20notebook.%0A%0A%20%20%20%201.%20**IVP%20audit.**%20Explain%20in%20your%20own%20words%20why%20an%20IVP%20needs%20both%20a%20differential%20equation%20and%20an%20initial%20condition.%20Then%20give%20one%20example%20where%20changing%20only%20the%20initial%20condition%20changes%20the%20entire%20motion.%0A%0A%20%20%20%202.%20**Pendulum%20amplitude.**%20For%20the%20pendulum%2C%20change%20only%20the%20initial%20angular%20velocity.%20What%20changes%20in%20the%20plot%20and%20what%20stays%20similar%3F%20Repeat%20with%20initial%20angles%20of%20%2410%5E%5Ccirc%24%2C%20%2445%5E%5Ccirc%24%2C%20and%20%2490%5E%5Ccirc%24.%0A%0A%20%20%20%203.%20**Projectile%20with%20terminal%20velocity.**%20A%20%2410%24%20kg%20projectile%20is%20launched%20vertically%20upward%20with%20%24v_0%3D50%5C%2C%5Cmathrm%7Bm%2Fs%7D%24.%20Compare%20the%20maximum%20height%20when%20drag%20is%20neglected%20with%20the%20height%20when%20the%20drag%20force%20is%20%24F_D%3D0.01v%7Cv%7C%24.%20Does%20the%20object%20approach%20a%20terminal%20velocity%20on%20the%20way%20down%3F%0A%0A%20%20%20%204.%20**Euler%20step-size%20audit.**%20For%20the%20football%20example%2C%20find%20a%20step%20size%20for%20Euler's%20method%20that%20makes%20the%20position%20error%20visually%20small.%20How%20small%20is%20small%20enough%20for%20your%20purpose%3F%20Explain%20your%20criterion.%0A%0A%20%20%20%205.%20**Two-dimensional%20projectile%20IVP.**%20Extend%20the%20vertical%20football%20model%20to%20two%20dimensions.%20Use%20state%20variables%20%24x%24%2C%20%24y%24%2C%20%24v_x%24%2C%20and%20%24v_y%24.%20Then%20add%20drag%20forces%20that%20oppose%20each%20velocity%20component.%0A%0A%20%20%20%206.%20**Design%20a%20mass-spring-damper%20system.**%20Use%20mass%20%24m%3D2%24%20kg%20and%20choose%20stiffness%20and%20damping%20values%20so%20the%20system%20returns%20close%20to%20equilibrium%20quickly%20but%20does%20not%20oscillate%20more%20than%20once.%20Report%20the%20values%20you%20tried%20and%20justify%20your%20final%20choice.%0A%0A%20%20%20%207.%20**One-mass%20spring%20problem.**%20A%20%241%24%20kg%20block%20is%20attached%20to%20a%20spring%20with%20stiffness%20%24k%3D100%5C%2C%5Cmathrm%7BN%2Fm%7D%24.%20At%20%24t%3D0%24%2C%20the%20block%20is%20displaced%20%240.1%24%20m%20from%20equilibrium%20and%20released%20from%20rest.%20Write%20the%20IVP%2C%20simulate%20the%20motion%20with%20%24c%3D0%24%2C%20and%20then%20find%20a%20damping%20value%20that%20removes%20visible%20oscillation.%0A%0A%20%20%20%208.%20**Compare%20solvers%20on%20an%20oscillator.**%20Compare%20semi-implicit%20Euler%20and%20%60solve_ivp%60%20on%20the%20same%20mass-spring-damper%20system.%20Which%20method%20changes%20more%20when%20you%20increase%20the%20step%20size%3F%20Which%20result%20would%20you%20use%20as%20a%20reference%20solution%3F%0A%0A%20%20%20%209.%20**Two-mass%20spring-damper%20system.**%20Two%20masses%20move%20on%20a%20horizontal%20line%20and%20are%20connected%20by%20a%20spring%20and%20damper.%20Define%20a%20state%20vector%2C%20write%20the%20first-order%20ODEs%2C%20and%20predict%20what%20happens%20to%20the%20distance%20between%20the%20masses%20when%20damping%20is%20increased.%0A%0A%20%20%20%2010.%20**Simple%20muscle-tendon%20model.**%20Treat%20the%20tendon%20as%20one%20spring%2C%20the%20muscle%20fibers%20as%20another%20spring%2C%20the%20viscous%20tissue%20as%20a%20damper%2C%20and%20the%20contractile%20element%20as%20an%20active%20force.%20Define%20a%20one-dimensional%20IVP%20for%20tendon%20length%20and%20explain%20what%20each%20parameter%20represents.%0A%0A%20%20%20%2011.%20**Force-platform%20record.**%20A%20force%20platform%20measures%20vertical%20ground%20reaction%20force%20during%20a%20jump.%20Given%20body%20mass%20and%20the%20force-time%20curve%2C%20write%20the%20ODE%20for%20the%20center-of-mass%20acceleration.%20What%20numerical%20integrations%20would%20you%20perform%20to%20estimate%20takeoff%20velocity%20and%20jump%20height%3F%0A%0A%20%20%20%2012.%20**Your%20own%20movement%20model.**%20Create%20your%20own%20IVP%20from%20a%20movement%20you%20care%20about.%20Define%20the%20state%20variables%2C%20initial%20conditions%2C%20model%20parameters%2C%20and%20one%20parameter%20you%20would%20like%20to%20estimate%20from%20data.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Final%20answers%20to%20Problems%0A%20%20%20%201.%20An%20ODE%20supplies%20the%20rate%20law%20and%20initial%20conditions%20select%20a%20trajectory.%20A%20pendulum%20released%20from%20different%20initial%20angles%20has%20different%20amplitude%20and%20energy.%0A%20%20%20%202.%20Initial%20angular%20velocity%20changes%20energy%2C%20phase%2C%20and%20possibly%20whether%20full%20rotations%20occur.%20At%2010%C2%B0%20the%20motion%20is%20approximately%20sinusoidal%3B%2045%C2%B0%20and%2090%C2%B0%20increasingly%20depart%20from%20the%20small-angle%20approximation%20and%20have%20longer%20periods.%0A%20%20%20%203.%20Heights%20above%20launch%20level%2C%20using%20%24g%3D9.8%24%20m%2Fs%C2%B2%3A%20127.55%20m%20without%20drag%3B%20113.61%20m%20with%20the%20stated%20quadratic%20drag.%20Downward%20terminal%20velocity%20is%20approximately%20%E2%88%9298.995%20m%2Fs%3B%20reaching%20it%20requires%20sufficient%20falling%20distance.%0A%20%20%20%204.%20Example%20criterion%3A%20at%20matched%20grid%20times%20through%204%20s%2C%20position%20error%20below%200.02%20m.%20For%20the%20constant-gravity%20forward%20Euler%20example%2C%20%24h%3D0.001%24%20s%20meets%20that%20criterion.%20A%20different%20accuracy%20requirement%20gives%20a%20different%20acceptable%20step.%0A%20%20%20%205.%20For%20vector%20quadratic%20drag%20with%20%24k%3D%5Crho%20C_dA%2F2%24%2C%20state%20%24(x%2Cy%2Cv_x%2Cv_y)%24%20has%20derivative%20%24(v_x%2Cv_y%2C-k%5C%7Cv%5C%7Cv_x%2Fm%2C-g-k%5C%7Cv%5C%7Cv_y%2Fm)%24.%20Specify%20four%20initial%20values.%0A%20%20%20%206.%20One%20valid%20choice%20is%20%24k%3D50%24%20N%2Fm%20and%20%24c%3D20%24%20Ns%2Fm%20for%20%24m%3D2%24%20kg%3A%20critical%20damping%20and%20no%20free%20oscillation.%20%E2%80%9CClose%E2%80%9D%20still%20needs%20a%20displacement%20tolerance%20and%20observation%20time.%0A%20%20%20%207.%20%24%5Cdot%20x%3Dv%24%2C%20%24%5Cdot%20v%3D-100x-cv%24%2C%20%24x(0)%3D0.1%24%20m%2C%20%24v(0)%3D0%24.%20With%20%24c%3D0%24%2C%20%24x%3D0.1%5Ccos(10t)%24%20m.%20Critical%20damping%20is%20%24c%3D20%24%20Ns%2Fm.%0A%20%20%20%208.%20Semi-implicit%20Euler%20changes%20more%20with%20its%20integration%20step.%20A%20%60solve_ivp%60%20result%20with%20tightened%20tolerances%20and%20demonstrated%20convergence%20is%20a%20useful%20reference%3B%20output%20spacing%20alone%20is%20not%20its%20internal%20integration%20step.%0A%20%20%20%209.%20Let%20%24q%3Dx_2-x_1-L_0%24%20and%20%24f%3Dkq%2Bc(v_2-v_1)%24.%20Then%20%24(%5Cdot%20x_1%2C%5Cdot%20x_2%2C%5Cdot%20v_1%2C%5Cdot%20v_2)%3D(v_1%2Cv_2%2Cf%2Fm_1%2C-f%2Fm_2)%24.%20Increasing%20damping%20suppresses%20relative%20oscillation%2C%20while%20excessive%20damping%20slows%20relaxation.%0A%20%20%20%2010.%20One%20massless%20series%20example%20uses%20tendon%20extension%20%24x%24%20and%20fixed%20total%20extension%20%24L%24%3A%20%24%5Cdot%20x%3D%5Bk_m(L-x)%2BF_a(t)-k_t%20x%5D%2Fc%24%2C%20%24x(0)%3Dx_0%24.%20Here%20%24k_t%2Ck_m%24%20are%20stiffnesses%2C%20%24c%24%20is%20muscle-side%20viscosity%2C%20and%20%24F_a%24%20is%20active%20force.%20This%20is%20one%20stated%20model%2C%20not%20a%20unique%20consequence%20of%20the%20verbal%20description.%0A%20%20%20%2011.%20%24%5Cdot%20v%3DF_z%2Fm-g%24%20and%20%24%5Cdot%20y%3Dv%24.%20Integrate%20from%20known%20initial%20velocity%20and%20position%3B%20after%20takeoff%2C%20ballistic%20rise%20is%20%24v_%7Bto%7D%5E2%2F(2g)%24%20when%20drag%20is%20neglected.%0A%20%20%20%2012.%20Example%3A%20damped%20knee-angle%20relaxation%20with%20state%20%24(%5Ctheta%2C%5Comega)%24%2C%20%24%5Cdot%5Ctheta%3D%5Comega%24%2C%20%24I%5Cdot%5Comega%3D-b%5Comega-k%5Ctheta%24%2C%20stated%20initial%20angle%20and%20velocity%2C%20and%20damping%20%24b%24%20estimated%20from%20the%20measured%20decay.%20Other%20well-defined%20IVPs%20are%20valid.%0A%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Go%20deeper%0A%0A%20%20%20%20Use%20these%20references%20when%20you%20want%20another%20explanation%20or%20more%20examples%3A%0A%0A%20%20%20%20-%20%5BSciPy%20%60solve_ivp%60%20documentation%5D(https%3A%2F%2Fdocs.scipy.org%2Fdoc%2Fscipy%2Freference%2Fgenerated%2Fscipy.integrate.solve_ivp.html)%20for%20solver%20options%2C%20events%2C%20tolerances%2C%20and%20examples.%0A%20%20%20%20-%20%5BSciPy%20integration%20tutorial%5D(https%3A%2F%2Fdocs.scipy.org%2Fdoc%2Fscipy%2Ftutorial%2Fintegrate.html)%20for%20a%20broader%20overview%20of%20numerical%20integration%20in%20Python.%0A%20%20%20%20-%20%5BOpenStax%20University%20Physics%3A%20Damped%20Oscillations%5D(https%3A%2F%2Fopenstax.org%2Fbooks%2Funiversity-physics-volume-1%2Fpages%2F15-5-damped-oscillations)%20for%20the%20physics%20of%20damping%2C%20driven%20oscillations%2C%20and%20energy%20loss.%0A%20%20%20%20-%20%5BMIT%20OpenCourseWare%3A%20Damped%20Harmonic%20Oscillators%5D(https%3A%2F%2Focw.mit.edu%2Fcourses%2F18-03sc-differential-equations-fall-2011%2Fpages%2Funit-ii-second-order-constant-coefficient-linear-equations%2Fdamped-harmonic-oscillators%2F)%20for%20a%20differential-equations%20treatment%20of%20spring-mass-damper%20systems.%0A%20%20%20%20-%20%5BIFAB%20Law%202%3A%20The%20Ball%5D(https%3A%2F%2Fwww.theifab.com%2Flaws%2Flatest%2Fthe-ball%2F)%20for%20current%20football%20size%20and%20mass%20ranges%20used%20in%20the%20projectile%20example.%0A%20%20%20%20-%20%5BWhat%20is%20the%20fastest%20possible%20volleyball%20serve%3F%5D(https%3A%2F%2Fuio-ccse.github.io%2Fcomputational-essay-showroom%2Fessays%2Fexampleessays%2Fvolleyball%2FVolleyball.html)%20for%20an%20example%20about%20the%20ball%20movement%20during%20a%20volleyball%20serve.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%23%20Video%20lectures%20on%20the%20Internet%0A%20%20%20%20-%20%5BEuler's%20method%5D(https%3A%2F%2Fwww.khanacademy.org%2Fmath%2Fdifferential-equations%2Ffirst-order-differential-equations%2Feulers-method-tutorial%2Fv%2Feulers-method)%20-%20Khan%20Academy.%0A%20%20%20%20-%20Numerical%20solutions%20of%20Ordinary%20Differential%20Equations%3A%20%5BIntegrating%20ODEs%5D(https%3A%2F%2Fyoutu.be%2FQBeNXHrAYns%3Fsi%3Dlf5vPgCjHs-8-Bc1)%2C%20%5BEuler's%20method%5D(https%3A%2F%2Fyoutu.be%2FMstPeOTCVzQ%3Fsi%3D_kLTgq0d1oI-fZ9O)%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20References%0A%0A%20%20%20%20-%20Bray%20K%2C%20Kerwin%20DG%20(2003)%20%5BModelling%20the%20flight%20of%20a%20soccer%20ball%20in%20a%20direct%20free%20kick%5D(https%3A%2F%2Fpeople.stfx.ca%2Fsmackenz%2Fcourses%2Fhk474%2Flabs%2Fjump%2520float%2520lab%2Fbray%25202002%2520modelling%2520the%2520flight%2520of%2520a%2520soccer%2520ball%2520in%2520a%2520direct%2520free%2520kick.pdf).%20Journal%20of%20Sports%20Sciences%2C%2021%2C%2075-85.%0A%0A%20%20%20%20-%20Downey%20AB%20(2023)%20%5BModeling%20and%20Simulation%20in%20Python%3A%20An%20Introduction%20for%20Scientists%20and%20Engineers%5D(https%3A%2F%2Fgreenteapress.com%2Fwp%2Fmodsimpy%2F).%20Green%20Tea%20Press.%0A%0A%20%20%20%20-%20IFAB%20(2026%2F27)%20%5BCurrent%20IFAB%20Law%202%3A%20The%20Ball%5D(https%3A%2F%2Fwww.theifab.com%2Flaws%2Flatest%2Fthe-ball%2F).%0A%0A%20%20%20%20-%20Kiusalaas%20J%20(2013)%20%5BNumerical%20Methods%20in%20Engineering%20with%20Python%203%5D(https%3A%2F%2Fapi.pageplace.de%2Fpreview%2FDT0400.9781139604413_A24435840%2Fpreview-9781139604413_A24435840.pdf).%203rd%20edition.%20Cambridge%20University%20Press.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20marimo%20as%20mo%0A%0A%20%20%20%20return%20(mo%2C)%0A%0A%0Aif%20__name__%20%3D%3D%20%22__main__%22%3A%0A%20%20%20%20app.run()%0A
b2d47e397a847111733de9dc147b447a