Fluent Essays

Notes. When uploading your work, be sure to choose correct pages, without missing any, for each problem. Besides, make sure all pages are in order, no page is rotated, and none of your code or hand-written work is cut out. According to our homework policies, you are allowed to collaborate for homework assignments in a group of up to 5 students. If you did, please list all your collaborators at the top of your homework, below your name. Homework must be prepared using a live script (.mlx file). Submission generated using any other means (e.g., Words with MATLAB screenshots) will not be accepted. Starting from this assignment, you will be asked to write MATLAB functions. Unlike scripts, functions can be included inside your mlx file at the very end of the document. So you no longer need to write an external m-file and print its content using the type function. For more detailed instruction, please read  [Guide] How to Embed Functions in Live Scripts:   Method 1: Using Code Example Environment Suppose that you have written a function as a separate function m-file, myfoo.m and now you want to include its contents within your homework live script, HW2_Kim.mlx. At the top of MATLAB window, click on Insert > Code Example. It creates a blank box. Any text typed within this box will be formatted in monospaced fonts with syntax highlighting. Drawbacks: The box is not an executable code block.  If you make changes to myfoo.m, you need update the code example block accordingly (by copying and pasting). Method 2: Using TYPE Function Still suppose that you have a function written as a separate m-file titled myfoo.m. Then you can simply type type myfoo within a code block. The contents of myfoo.m will be printed out in your live script. Note: We have been using this to print out the contents of external script m-files. This works in the exact same way for function m-files. Drawbacks: If you forget to run that block after your final edit, your live script will show an outdated version of your function. As you complete more and more assignments, your folders will be filled with more m-files and may be difficult to manage unless you are well-organized. Method 3 (Recommended): Writing Functions at the End of Live Script Instead of writing a separate m-file, you may just write the function at the very bottom of the document. The functions are local to the live script, and so you are able to call the function anywhere inside the live script. This way, your mlx file is fully self-contained and you no longer need to go back and forth between Editor and Live Editor. See the live script accompanying Lecture 9; all functions defined for the live script are gathered in the last section titled "Functions Used". This is my recommended method! The only caveat I can think of: Suppose Problem 1(a) asks you to write a function, and so you write it at the end of your live script. You then exports it to a pdf file and find out that Problem 1 solutions are on pp. 1--2, with the function written for 1(a) appearing on p. 7. In this case, when you upload the pdf file to Gradescope, you have to select pages 7, 1, and 2 for Problem 1, in that order. Guide video link: https://www.youtube.com/watch?v=lIw-cGEDSO0 https://www.youtube.com/watch?v=JjJADydaQBk https://www.youtube.com/watch?v=iG-Pye1EcbU https://www.youtube.com/watch?v=S1tpRlp1lHI https://www.youtube.com/watch?v=Dfz8BPRykS4 https://www.youtube.com/watch?v=4BXhz_kRB2c Updated on September 26, 2021 Math 3607: Homework 4 Due: 10:00PM, Tuesday, September 28, 2021 TOTAL: 30 points • Problems marked with v are to be done by hand; those marked with � are to be solved using a computer. • Important note. Do not use Symbolic Math Toolbox. Any work done using sym or syms will receive NO credit. • Another important note. Starting from this assignment, you will be asked to write MATLAB functions. Instead of writing an external function m-file, include all your functions at the end of your live script. 1. (Sliders moving along grooves; adapted from LM 2.1–12 and Sample HW01) The mechanical device shown in Figure 1 consists of two grooves in which sliders slide. These sliders are connected to a straight rod. ` L A B Cpx,yq W O Figure 1: The bronze square is a piece of metal with two grooves cut out of it as shown. There are sliders at the points A and B which slide in these grooves. The slider at A can only slide vertically, and the one at B can only slide horizontally. There is a straight rod attached at A and B, which ex- tends to C. As the point C moves around the block, it traces out a closed curve. (a) v Analytically, determine the curve which is traced out by C in one rotation. Suggestion. Let px,yq be the coordinates of the point C. Express the variables x,y in terms of L, `, and θ, where θ P r0, 2πq is the angle from the part of the horizontal groove which is to the right of B to the rod BC. (b) � Using the previous result, plot the trajectory of C in one rotation for ` “ 2 and L “ 7. 2. (Spiral triangles to spiral polygons; adapted from LM 5.9–7, 6.8–34) � The following script1 generates spirals using equilateral triangles as shown in the figure below. 1It is slightly modified from the code included in Lecture 9 slides. Note the introduction of a new variable d_rot, which is accountable for the rotation of the innermost triangle. 1 m = 21; d_angle = 4.5; d_rot = 90; th = linspace(0, 360, 4) + d_rot; V = [cosd(th); sind(th)]; C = colormap(hsv(m)); s = sind(150 - abs(d_angle))/sind(30); R = [cosd(d_angle) -sind(d_angle); sind(d_angle) cosd(d_angle)]; hold off for i = 1:m if i > 1 V = s*R*V; end plot(V(1,:), V(2,:), ’Color’, C(i,:)) hold on end set(gcf, ’Color’, ’w’) axis equal, axis off Figure 2: Spiral triangles with m “ 21 and θ “ 4.5˝. (a) Write a function named spiralgon by modifying the script so that it generates spirals using m regular n-gons for any n ě 3. Your function must be written at the end of your homework live script (.mlx) file. Begin the function with the following header and comments. function V = spiralgon(n, m, d_angle, d_rot) % SPIRALGON plots spiraling regular n-gons % input: n = the number of vertices % m = the number of regular n-gons % d_angle = the degree angle between successive n-gons % (can be positive or negative) % d_rot = the degree angle by which the innermost n-gon % is rotated % output: V = the vertices of the outermost n-gon .... (b) Run the statements below to generate some aesthetic shapes. clf subplot(2, 2, 1), spiralgon(3, 41, 4.5, -90); subplot(2, 2, 2), spiralgon(4, 37, -2.5, 45); subplot(2, 2, 3), spiralgon(5, 61, 3, -90); subplot(2, 2, 4), spiralgon(8, 91, -4, 22.5); Note. Copy the five lines, paste them inside a single code block, and run it. This code block must precede your function(s). 3. (Machine epsilon; adapted from LM 9.3–3(a)) � Recall that the number in the computer which follows 1 is 1 ` eps , which can be verified in MATLAB by >> format long 2 >> (1 + 0.5*eps) - 1 ans = 0 >> (1 + 0.51*eps) - 1 ans = 2.220446049250313e-16 In the same manner: (a) Verify that the number in the computer which follows 8 is 8 ` 8 eps by numerically calculating 8 ` 4 eps and 8 ` 4.01 eps . (b) Verify that the number in the computer which precedes 16 is 16 ´ 8 eps by numerically calculating 16 ´ 4.01 eps and 16 ´ 4 eps . (c) What are the numbers in the computer that precedes and follows 210 “ 1024, respec- tively? Verify your claims in MATLAB by carrying out appropriate calculations. Note. Begin with format long as shown in the example above. This is needed only once before the beginning of part (a). Note. Answer each part of the problem in a single code block. No external script needs to be written. 4. (Catastrophic cancellation; LM 9.3–10) We revisit the function from Problem 3 of Homework 3. Consider the function fpxq “ $ ’ & ’ % ex ´ 1 x if x ‰ 0 1 if x “ 0, where we are interested in exploring the catastrophic cancellation which occurs as x Ñ 0 since ex Ñ 1 as x Ñ 0. (a) v Use the Taylor series expansion of ex to prove that f is continuous at 0. (b) � Now calculate fpxq numerically for x “ 10´k where k P Nr1, 20s in three slightly different ways: i. Calculate fpxq as written. ii. Calculate it as f1pxq “ ex ´ 1 log ex , for x ‰ 0. (You and I know that analytically f1pxq ” fpxq for all nonzero x – but MATLAB doesn’t.) iii. MATLAB has a function which analytically subtracts 1 from the exponential to avoid catastrophic cancellation before the result is calculated numerically. So define the function f2pxq to be the same as fpxq except that ex ´ 1 is replaced by expm1(x). Tabulate the results using disp or fprintf. The table should have four columns with the first being x, the second using fpxq, the third using f1pxq, and the fourth using f2pxq, with all shown to full accuracy. Do it as efficiently as you can, without using a loop. Note. Write your code for this part in a single code block. No external script needs to be written. 3 (c) v Comment on the results obtained in the previous part. Explain why certain methods work well while others do not. 5. (Inverting hyperbolic cosine; FNC 1.3.6) The function x “ coshptq “ et ` e´t 2 can be inverted to yield a formula for acoshpxq: t “ log ´ x ´ a x2 ´ 1 ¯ . (‹) In MATLAB, let t=-4:-4:-16 and x=cosh(t). (a) v� Find the condition number of the problem fpxq “ acoshpxq by hand. (You may use Equation (‹), or look up a formula for f 1 in a calculus book.) Then evaluate κf at the elements of x in MATLAB. (b) � Evaluate the right-hand side of Equation (‹) using x to approximate t. Record the accuracy of the answers (by displaying absolute and/or relative errors), and explain. (Warning: Use format long to get enough digits or use fprintf with a suitable format.) (c) � An alternate formula for acoshpxq is t “ ´2 log ˜ c x ` 1 2 ` c x ´ 1 2 ¸ . (:) Apply Equation (:) to x and record the accuracy as in part (b). Comment on your observation. (d) v Based on your experiments, which of the formulas (‹) and (:) is unstable? What is the problem with that formula? Note. Write your code for each of parts (a), (b), and (c) in a single code block. No external script needs to be written. 4 function x = ieee2fp(ieee) % calculates the floating-point number from its ieee representation % input: ieee = 8-character hexidecimal string or % 64-bit binary string % = (structure) when representation is in scientific notation % +- 1.xxxxx... 2^exponent or % +- 0.xxxxx... 2^(-1023) % .sign = +1 or -1 as a fp number, not text % .exponent = the exponent shown above, must be in [-1023, 1024] % .mantissa = 13-character hexidecimal string or % 52-bit binary string bin2hex = [(0:9)+'0', (0:5)+'a']; if ischar(ieee) if length(ieee) == 16 if any( (ieee < '0' | ieee > 'f') | (ieee > '9' & ieee < 'a') ) error(['***** hexidecimal string not all 0-9 or a-f: ', ieee]) end x = hex2num(ieee); elseif length(ieee) == 64 if any( ieee < '0' | ieee > '1' ) error(['***** binary string not all 0-1: ', ieee]) end Ieee = reshape(ieee, 4, 16)'; Ieee_dec = bin2dec(Ieee); ieee_hex = bin2hex(Ieee_dec+1); x = hex2num(char(ieee_hex)); else error(['***** string not correct length: ', ieee]) end else if ~isstruct(ieee) error(['***** input is not a structure variable: ', num2str(ieee)]) end if ieee.sign == +1 exponent_bin(1) = '0'; elseif ieee.sign == -1 exponent_bin(1) = '1'; else error(['***** incorrrect sign: ', num2str(ieee.sign)]) end if mod(ieee.exponent, 1) ~= 0 || ieee.exponent < -(2^10 - 1) || ieee.exponent > 2^10 error(['***** incorrect exponent: ', ... num2str(ieee.exponent)]) end exponent_bin(2:12) = dec2bin(ieee.exponent+(2^10-1), 11); if length(ieee.mantissa) == 52 if any( ieee.mantissa<'0' | ieee.mantissa>'1' ) error(['***** binary string not all 0-1: ', ieee]) end x = ieee2fp([exponent_bin, ieee.mantissa]); elseif length(ieee.mantissa) == 13 if any( (ieee.mantissa<'0'|ieee.mantissa>'f') | ... (ieee.mantissa>'9'&ieee.mantissa<'a') ) error(['***** hexidecimal mantissa not all 0-9 or a-f: ', ... ieee.mantissa]) end Ieee = reshape(exponenet_bin, 4, 3)'; Ieee_dec = bin2dec(Ieee); ieee_hex = bin2hex(Ieee_dec+1); x = hex2num([char(ieee_hex), ieee.mantissa]); else error(['***** incorrect mantissa length: ', num2str(ieee.mantissa)]) end end end function ieee = fp2ieee(x) % calculates the IEEE representation of a floating-point number in binary % input: x = floating-point number % output: ieee.bits = (structure) IEEE representation, i.e., 64 bits of 0's and 1's % .split = IEEE representation separated for easier reading % .sci = scientific notation: % +- 1.xxxxx... 2^exponent or % +- 0.xxxxx... 2^(-1023) % .sign = +1 if x >= 0 % = -1 if x < 0 if ~isfinite(x) error('***** x is not a finite number') end x_hex = num2hex(x); x_exp_dec = hex2dec(x_hex(1:3)); ieee.sign = +1; x_sign_char = '+'; x_sign_bin = '0'; if x_exp_dec >= 2^11 x_exp_dec = x_exp_dec - 2^11; ieee.sign = -1; x_sign_char = '-'; x_sign_bin = '1'; end x_real_exp = x_exp_dec - (2^10 - 1); x_mantissa = dec2bin(hex2dec(x_hex(4:16)), 52); ieee.bits = [dec2bin(hex2dec(x_hex(1:3)), 12), x_mantissa]; if x_exp_dec > 0 ieee.sci = [x_sign_char, '1.', x_mantissa]; else ieee.sci = [x_sign_char, '0.', x_mantissa]; end ieee.sci = [ieee.sci, ' x 2^', num2str(x_real_exp)]; ieee.split = [x_sign_bin, '|', split(dec2bin(x_exp_dec, 11)), '|', ... split(x_mantissa)]; end function split_string = split(string) if length(string) == 0 return end if length(string) <= 4 split_string = string; else split_string = [split(string(1:end-4)), ' ', string(end-3: ... end)]; end end function floatgui(callbackarg) %FLOATGUI Show structure of floating point numbers. % The set of positive model floating point numbers is determined % by three parameters: t, emin, and emax. It is the set of rational % numbers of the form x = (1+f)*2^e where f = (integer)/2^t, % 0 <= f < 1, e = integer, and emin <= e <= emax. % % IEEE 754 double precision has t = 52, emin = -1022, emax = 1023. % Copyright 2014 Cleve Moler % Copyright 2014 The MathWorks, Inc. % Initialize parameters if nargin == 0 t = 3; emin = -4; emax = 2; logscale = 0; Fpos=[50 300 900 250]; else t = round(get(findobj('tag','t'),'value')); emin = round(get(findobj('tag','emin'),'value')); emax = round(get(findobj('tag','emax'),'value')); logscale = get(findobj('style','check'),'value'); Fpos=get(gcf,'pos'); end % Position figure window shg clf reset set(gcf,'pos',Fpos,'name','floatgui', ... 'resize','on','defaultuicontrolunits','normalized',... 'numbertitle','off','menubar','none') % Generate and plot floating point numbers f = (0:2^t-1)/2^t; F = []; for e = emin:emax F = [F (1+f)*2^e]; end for x = F text(x,0,'|','fontunits','normalized','fontsize',0.3) end % Set axes set(gca,'pos',[.05 .6 .9 .2],'fontunits','normalized','fontsize',0.22) if logscale set(gca,'xscale','log') xmin = 1/2^(-emin+.5); xmax = 2^(emax+1.5); else set(gca,'xscale','linear') xmin = 0; xmax = 2^(emax+1); end axis([xmin xmax -1 1]) % Set tick marks fmin = min(F); fmax = max(F); xtick = 1; xticklab = {'1'}; if fmin < 1 xtick = [1/2 xtick]; xticklab = ['1/2' xticklab]; end if logscale & (fmin < 1/4) xtick = [1/4 xtick]; xticklab = ['1/4' xticklab]; end if fmin < 1/2 xtick = [fmin xtick]; xticklab = [['1/' int2str(1/fmin)] xticklab]; end if 2 < fmax xtick = [xtick 2]; xticklab = [xticklab '2']; end if 4 < fmax xtick = [xtick 4]; xticklab = [xticklab '4']; end if max(xtick) < fmax xtick = [xtick fmax]; if fmax == round(fmax) fmaxlab = int2str(fmax); else over = 2^(emax+1); fmaxlab = [int2str(over) '-1/' int2str(1/(over-fmax))]; end xticklab = [xticklab fmaxlab]; end set(gca,'xtick',xtick,'xticklabel',xticklab,'xminortick','off','ytick',[]) % Create uicontrols uicontrol('style','slider','tag','emin','value',emin, ... 'min',-8,'max',0,... 'pos',[0.15 0.26 0.13 0.07],'sliderstep',[1/8 1/8], ... 'callback','floatgui(1)'); uicontrol('style','slider','tag','t','value',t, ... 'min',0,'max',8,... 'pos',[0.435 0.26 0.13 0.07],'sliderstep',[1/8 1/8], ... 'callback','floatgui(1)'); uicontrol('style','slider','tag','emax','value',emax, ... 'min',0,'max',8,... 'pos',[0.72 0.26 0.13 0.07],'sliderstep',[1/8 1/8], ... 'callback','floatgui(1)'); uicontrol('style','text','string',['emin = ' int2str(emin)], ... 'pos',[0.15 0.35 0.13 0.07],'fontunits','normalized','fontsize',0.7) uicontrol('style','text','string',['t = ' int2str(t)], ... 'pos',[0.435 0.35 0.13 0.07],'fontunits','normalized','fontsize',0.7) uicontrol('style','text','string',['emax = ' int2str(emax)], ... 'pos',[0.72 0.35 0.13 0.07],'fontunits','normalized','fontsize',0.7) uicontrol('style','check','string','log scale','value',logscale, ... 'pos',[0.435 0.15 0.13 0.07],'fontunits','normalized','fontsize',0.7, ... 'callback','floatgui(1)'); uicontrol('style','push','pos',[0.88 0.1 0.07 0.07], ... 'fontunits','normalized','fontsize',0.7, ... 'string','close','callback','close(gcf)') % eps if fmax > 1 eps = 2^(-t); text(1,0,'|','color','r','fontunits','normalized','fontsize',0.3) text(1+eps,0,'|','color','r','fontunits','normalized','fontsize',0.3) if eps < 1 text(1.0,1.5,['eps = 1/' int2str(1/eps)], ... 'fontunits','normalized','fontsize',0.3,'fontweight','bold') else text(1.0,1.5,'eps = 1', 'fontunits','normalized','fontsize',0.3,'fontweight','bold') end end % Number of numbers % Exercise: % How many "floating point" numbers are in the set? % Complete this statement. % text(.9*xmax,2,num2str(???) Catastrophic Cancellation in log(1+x) Tae Eun Kim September 23, 2021 This note explains the catastrophic cancellation observed in Problem 4 of Homework 4. The evaluation of f(x) is severely affected by catastrophic cancellation for small x because of what is written at the beginning of the problem. Though identical to f(x) mathematically, the function f1(x) does a better job because of the following reason. Let x̂ = ̂(1 + x) − 1 = fl((1 +x) − 1), the floating-point representation of the expression (1 +x) − 1. Note that the subtraction undergoes catastrophic cancellation for small x. Also note that when log(1 + x) is evaluated in the computer, the input (1 + x) is formed first and then 1 is subtracted off from it before it is used in the algorithm based on the Taylor series log ξ = (ξ − 1) − 1 2 (ξ − 1)2 + 1 3 (ξ − 1)3 − · · · . Therefore the numerical evaluation of f1(x) can be approximated by f̂1(x) ≈ x̂ − 12x̂ 2 + 13x̂ 3 − · · · x̂ = 1 − 1 2 x̂ + 1 3 x̂2 − · · · , which resembles the series expansion used in part (a). It is clear from the right-hand side that this implementation does not involve subtraction of two nearby numbers for reasonably small x. However, when x is sufficiently small, (1 + x) is indistinguishable from 1 on the floating-point number system, in which case x̂ = ̂(1 + x) − 1 = 0. We know that it happens when x is smaller than the machine epsilon eps, which is about 2 × 10−16; in our experiment, it happens for k ≥ 16, in which case both the numerator and the denominator of f1(x) are evaluated as zeros, resulting in NaN. One way to avoid the issue with NaN is to use the series expansion obtained in part (a), namely, f(x) = 1 − 1 2 x + 1 3 x2 − · · · . Another way is to use the function log1p as suggested in the problem. This function was specifi- cally designed to avoid catastrophic cancellation occurring in the evaluation of log(1 + x) for small x by encoding log(1 + x) = x − 1 2 x2 + 1 3 x3 − · · · instead of using the series expansion for log ξ written above. 1 [Content_Types].xml _rels/.rels matlab/document.xml Homework 1 Math 3607, Spring 2021 [Your Name] Table of Contents Problem 1. (LM 2.1--6) Problem 2. (LM 2.1--8) Part (a) (i) Part (b) Part (c) Problem 3. (LM 2.1--14) Part (a) Problem 4. (Temperature Conversion) Problem 5. (Oblate Spheroid) Problem 1. (LM 2.1--6) Your answer here. % Code goes here Problem 2. (LM 2.1--8) Part (a) (i) (ii) (vi) Part (b) Part (c) Problem 3. (LM 2.1--14) Part (a) Problem 4. (Temperature Conversion) Problem 5. (Oblate Spheroid) matlab/output.xml manual code ready 0.4 metadata/coreProperties.xml 2021-01-13T21:38:03Z 2021-01-13T21:38:03Z metadata/mwcoreProperties.xml application/vnd.mathworks.matlab.code MATLAB Code R2020b metadata/mwcorePropertiesExtension.xml 9.9.0.1444674 93b378b9-5b5d-4623-b00d-0e154b3a183c metadata/mwcorePropertiesReleaseInfo.xml 9.9.0.1524771 R2020b Update 2 Nov 03 2020 2207788044