hw2p2 - % minus the trivial root k = 0, which we didn't...

Info iconThis preview shows page 1. Sign up to view the full content.

View Full Document Right Arrow Icon
function k = hw2p2(omega, h) % Assuming that kh << 1 so that we can make the approximation that % tanh = kh + epsilon for epsilon <<< 1 % Find the root to (T_p)*h*k^4 + g*h*k^2 - w^2 = 0 as the starting point of % newton's method T_p = 7.2E-5; g = 9.8; tolerance = 1.E-12; % root of the above approximation k0 = sqrt((-g*h + sqrt((g*h)^2+4*omega^2*T_p*h))/(2*T_p*h)); drawline(k0); k = k0; y = abs(k * (g+T_p * k ^ 2) * tanh(h * k) - omega^2); % since let g(x) be the original function k(g+T_p*k^2)tanh(k*h)-w^2, % then f(x) = g + t_p*k^2 - w^2/(k*tanh(k*h)) will have the same root
Background image of page 1
This is the end of the preview. Sign up to access the rest of the document.

Unformatted text preview: % minus the trivial root k = 0, which we didn't want from the start. while abs(y) &gt;= tolerance k = k - gap(k, omega, h, T_p, g); drawline(k); y = k * (g+T_p * k ^ 2) * tanh(h * k) - omega^2; end hold on ks = k0:0.002:k+0.01; plot(ks, ks .* (g+T_p .* ks .^ 2) .* tanh(h .* ks)-omega^2), grid on; hold off end function dk = gap(k, w, h, T_p, g) % gap = |k[n+1]-k[n]| within the newton iterations dk = (g + T_p*k^2 - w^2/(k*tanh(k*h))) / (2*T_p*k^2+w^2/k^2*(coth(h*k) +h*k*csch(h*k)^2)); end function drawline(k) hold on line([k,k], [-1, 1]) hold off end...
View Full Document

Ask a homework question - tutors are online