1function [y,yp,ypp] = generalwindow(x,data)
2%
3% [y,yp,ypp] = generalwindow(tol,x,data)
4%
5% evaluates blending window on [0,1] which decays to O(tol) at 0^+,
6% and 1-O(tol) at 1^-. y is vals, yp are derivs, ypp are 2nd derivs.
7% yp has support in [0,1],
8% x can be an array, and y, yp, ypp will be the same size as x.
9%
10% Note: does not return function handles!
11%
12% To test, call without arguments.
14if nargin==0, test_generalwindow; return; end
16tol = data.tol;
17gam = data.gam;
19% handle values x<0 and x>1...
20y = 0*x; yp= 0*x; ypp = 0*x; % also makes outputs correct sizes
21y(x>=1.0) = 1.0;
22% remaining x's actually need to evaluate something (not 0 or 1)...
23jj = x >0.0 & x < 1.0; % indices to eval
25if isempty(data.info) % slow evaluation
26 theta = log(1/tol); % set up window
27 W = ceil(2*theta/(pi*gam));
28 t = W*x(jj); % rescaled ordinates
29 [phit1,phitt1] = window(theta,W); % get func handles for [0,W]
30 yp(jj) = W*phit1(t); % make domain of x be [0,1]. w/ Jacobian
31 ypp(jj) = W*W*phitt1(t);
32 y(jj) = blending(phit1,t,tol); % must be scalar input
34else % use data fast eval
35 t = x(jj); % rescaled ordinates
36 y(jj) = chebEval(t,data.wei,data.info);
37 yp(jj) = chebEval(t,data.weip,data.info);
38 ypp(jj) = chebEval(t,data.weipp,data.info);
39end
41end
43 %%%%%%%%%%%%%
44function test_generalwindow
45verb = 1;
46tol = 1e-12;
47% setup the data struct
48data = struct('tol',tol,'gam',0.5,'info',[]);
50x = -1:1e-2:2.0;
51[y,yp,ypp] = generalwindow(x,data);
52if verb,
54 %figure; plot(x,[y;yp;ypp],'o-'); legend('y','yp','ypp');
56 %%% testing here
57 delta = 0.5;
58 t = x*delta;
59 [y_new,~,~] = generalwindow((t./delta),data);
60 figure;
61 plot(x,y,'o-'); xline(1); hold on;
62 plot(t,y_new,'r*'); xline(delta)
63end
65n=1e3; z = rand(n,1);
66tic;
67[Y,Yp,Ypp] = generalwindow(z,data);
68fprintf("throughput of generalwindow slow = %.3g pts/sec\n",n/toc)
70% test the fast cheb eval...
71data = setup_generalwindow(tol);
72%data
73%data.info
74[yf,ypf,yppf] = generalwindow(x,data);
76disp("test fast vs slow eval...")
77norm(y-yf,inf)
78norm(yp-ypf,inf)/norm(yp,inf)
79norm(ypp-yppf,inf)/norm(ypp,inf)
81n=1e5; z = rand(n,1);
82tic;
83[Y,Yp,Ypp] = generalwindow(z,data);
84fprintf("throughput of generalwindow cheb = %.3g pts/sec\n",n/toc)
85end