Skip to content

Instantly share code, notes, and snippets.

# StefanRickli/polytest.m Created Jul 17, 2017

Matlab script to test NAClab's root finder 'uvFactor'
 % Matlab script by Stefan Rickli % % This Matlab script calls uvFactor of the NAClab toolbox by Zhonggang Zeng % Get the Toolbox at http://homepages.neiu.edu/~zzeng/naclab.html %% user variables polydeg_min = 16; polydeg_max = 23; %% init f1 = figure; s1 = subplot(1,1,1); hold on; f2 = figure; s2 = subplot(1,1,1); hold on; f3 = figure; s3 = subplot(1,1,1); labels1 = {}; labels2 = {}; vertscat1 = []; vertscat2 = []; ratios = []; %% root finding for ii = polydeg_min:polydeg_max coeffs = fliplr(poly(1:ii)); ratios(end+1) = coeffs(end)/coeffs(1); fprintf('\n\n\nRatio between constant coeff and highest order coeff: %d:%d = %d\n', coeffs(end), coeffs(1), coeffs(end)/coeffs(1)); coeffs_perturbed = [coeffs(1:end-1),0]; try f1 = uvFactor(coeffs,1e-9); f2 = uvFactor(coeffs_perturbed,1e-9); roots = -f1(:,2)./f1(:,1); roots_perturbed = -f2(:,2)./f2(:,1); scatter(s1, real(roots),imag(roots)); scatter(s2, real(roots_perturbed),imag(roots_perturbed)); labels1{end+1} = ['deg(p_2) = ', num2str(ii)]; labels2{end+1} = ['deg(p_3) = ', num2str(ii)]; vertscat1(end+1) = max(abs(imag(roots))); vertscat2(end+1) = max(abs(imag(roots_perturbed))); catch err disp(getReport(err,'extended')); end end %% plots ylim(s1, [-0.15,0.15]); title(s1, 'p_2(x) for different max poly degree'); legend(s1, labels1); ylim(s2, [-0.15,0.15]); title(s2, 'p_3(x) (perturbed) for different max poly degree'); legend(s2, labels2); n = polydeg_min:polydeg_max; plot(s3, n, -0.01*log10(abs(ratios)), n, vertscat1, n, vertscat2); legend(s3, '-0.01*log10(abs(coefficient ratios))', 'vertical scattering p_2(x)', 'vertical scattering p_3(x)');
to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.