AdvCalc1

.docx

School

University of Waterloo *

*We aren’t endorsed by this school

Course

216

Subject

Electrical Engineering

Date

Dec 6, 2023

Type

docx

Pages

Uploaded by DrTitaniumSeal63 on coursehero.com

University of Waterloo Faculty of Engineering Department of Electrical and Computer Engineering Hello. NE 216 Laboratory 4 FALL 2011 Prepared by Kapur , Vamika UW Student ID Number: 20995888 UW User ID: v4kapur @uwaterloo.ca Prepared by Gupta , Shreeja UW Student ID Number: 20996288 UW User ID: s25gupta @uwaterloo.ca 2A Nanotechnology Engineering

14/11/23 Having watched the presentation, we will go through the following exercise: 4.1-7 Your Solution In Laboratory 2, you were required to document all aspects of your algorithm. You were also required to check all parameters for their types, etc . You must do the same for this laboratory; however, you will only need to provide the final solution for the Dormand-Prince algorithm; however, you will still be required to implement and test an implementation of the 4 th -order Runge Kutta algorithm with the signature function [t_out, y_out] = rk4( f, t_rng, y0, n ) Copy and paste your entire Matlab function for dp45 with comments and error checking here: function [t_out, y_out] = dp45( f, t_rng, y0, h, eps_abs ) if ~isa(f, 'function_handle' ) throw(MException( 'MATLAB:invalid_argument' , 'f is not a function handle' )); end if ~all(size(t_rng) == [1, 2]) throw(MException( 'MATLAB:invalid_argument' , 't_rng is not a 2-D row vector' )); end if ~isscalar(y0) throw( MException( 'MATLAB:invalid_argument' , 'y0 is not scalar' )); end if ~isscalar(h) || (h <= 0) throw(MException( 'MATLAB:invalid_argument' , 'h is not positive, scalar, or both' )); end if ~isscalar(eps_abs) || (eps_abs <= 0) throw( MException( 'MATLAB:invalid_argument' , 'eps_abs is not positive, scalar, or both' )); end A = [ 0 0 0 0 0 0 0; 1 0 0 0 0 0 0; 1/4 3/4 0 0 0 0 0; 11/9 -14/3 40/9 0 0 0 0; 4843/1458 -3170/243 8056/729 -53/162 0 0 0; 9017/3168 -355/33 46732/5247 49/176 -5103/18656 0 0; 35/384 0 500/1113 125/192 -2187/6784 11/84 0]'; by = [5179/57600 0 7571/16695 393/640 -92097/339200 187/2100 1/40]'; bz = [ 35/384 0 500/1113 125/192 -2187/6784 11/84 0]'; c = [0 1/5 3/10 4/5 8/9 1 1]'; K = zeros(1, 7); % Initialize t_out and y_out % Initialize our location to k = 1 2

t_out = t_rng(1); y_out = y0; k = 1; while t_out(k) < t_rng(2) for m = 1:7 K(m) = f(t_out(k) + h*c(m), y_out(k) + h*c(m)*K*A(:,m)); end %end of for loop y_tmp = y_out(k) + h*K*by; z_tmp = y_out(k) + h*K*bz; %error factor s = ((h.*eps_abs)./(2 * (t_rng(2) - t_rng(1)) * abs(y_tmp - z_tmp))).^(1/4); if s >= 2 y_out(k + 1) = z_tmp; t_out(k + 1) = t_out(k) + h; k = k + 1; h = 2.*h; else if s >= 1 y_out(k + 1) = z_tmp; t_out(k + 1) = t_out(k) + h; k = k + 1; else h = h./2; end %end of while loop %final check if t_out(k) + h > t_rng(2) h = t_rng(2) - t_out(k); end end end 4.8 Testing your Implementation Warning: All the questions you will do in parts a , b , and c will be repeated with the Dormand-Prince method in parts d , e , and f , respectively. Use scripts to reduce your workload! Consider using a for loop for parts c and f . We are now ready to test your code. For all questions, use format long 3

4.8 a With examples from the presentations, we had: function [dy] = f4a(t, y) dy = (y - 1).^2 .* (t - 1).^2; end function [y] = y4a_soln( t ) y = (t.^3 - 3*t.^2 + 3*t)./(t.^3 - 3*t.^2 + 3*t + 3); end and we approximated the solution with: [t2a, y2a] = rk4( @f4a, [0, 1], 0, 11 ) plot( t2a, y2a, 'or' ) hold on [t2a, y2a] = ode45( @f4a, [0, 1], 0 ); plot( t2a, y2a, 'b' ) It’s up to you to make sure your solution matches the output shown in the presentation. Save the right-hand side of the differential equation           (1) cos y t t y t y t t y t      as the function f 4 b and execute the following: [t2b, y2b] = rk4( @f4b, [0, 1], 0, 11 ) plot( t4b, y4b, 'or' ) [t2b, y4b] = ode45( @f4b, [0, 1], 0 ); hold on plot( t4b, y4b, 'b' ) title( 'uwuserid and uwuserid' ) % or title( 'uwuserid' ) 4

Your preview ends here

Eager to read complete document? Join bartleby learn and gain access to the full version