diff --git a/chunkie/+chnk/+flex2d/kern.m b/chunkie/+chnk/+flex2d/kern.m
index 7eefe85..aa1f32c 100644
--- a/chunkie/+chnk/+flex2d/kern.m
+++ b/chunkie/+chnk/+flex2d/kern.m
@@ -228,6 +228,37 @@
 
 end
 
+% clamped plate kernels for far field evaluation
+if strcmpi(type, 'clamped_plate_eval_ff')
+
+    submat = zeros(nt,2*ns);
+
+    srcnorm = srcinfo.n;
+    srctang = srcinfo.d;
+    nx = repmat(srcnorm(1,:),nt,1);
+    ny = repmat(srcnorm(2,:),nt,1);
+    dx = repmat(srctang(1,:),nt,1);
+    dy = repmat(srctang(2,:),nt,1);
+    ds = sqrt(dx.*dx+dy.*dy);
+
+    taux = dx./ds;
+    tauy = dy./ds;
+
+    [~, ~, hess, third] = chnk.flex2d.hkdiffgreen_ff(zk, src, targ);           % Hankel part
+
+    K1 = -(1/(2*zk^2).*(third(:, :, 1).*(nx.*nx.*nx) + third(:, :, 2).*(3*nx.*nx.*ny) +...
+       third(:, :, 3).*(3*nx.*ny.*ny) + third(:, :, 4).*(ny.*ny.*ny)) ) - ...
+       (3/(2*zk^2).*(third(:, :, 1).*(nx.*taux.*taux) + third(:, :, 2).*(2*nx.*taux.*tauy + ny.*taux.*taux) +...
+       third(:, :, 3).*(nx.*tauy.*tauy + 2*ny.*taux.*tauy) + third(:, :, 4).*(ny.*tauy.*tauy)));  % G_{ny ny ny} + 3G_{ny tauy tauy}
+
+    K2 =  -(1/(2*zk^2).*(hess(:, :, 1).*(nx.*nx) + hess(:, :, 2).*(2*nx.*ny) + hess(:, :, 3).*(ny.*ny)))+...
+          (1/(2*zk^2).*(hess(:, :, 1).*(taux.*taux) + hess(:, :, 2).*(2*taux.*tauy) + hess(:, :, 3).*(tauy.*tauy))); % -G_{ny ny}  + G_{tauy tauy}
+
+    submat(:,1:2:end) = K1;
+    submat(:,2:2:end) = K2;
+
+end
+
 
 %%% FREE PLATE KERNELS
 
@@ -425,6 +456,35 @@
 
 end
 
+% free plate kernels for far field evaluation
+if strcmpi(type, 'free_plate_eval_ff')           
+   srcnorm = srcinfo.n;
+   srctang = srcinfo.d;
+   nu = varargin{1};
+
+   [val,grad] = chnk.flex2d.hkdiffgreen_ff(zk,src,targ); 
+   nx = repmat(srcnorm(1,:),nt,1);
+   ny = repmat(srcnorm(2,:),nt,1);
+
+   dx = repmat(srctang(1,:),nt,1);
+   dy = repmat(srctang(2,:),nt,1);
+
+   ds = sqrt(dx.*dx+dy.*dy);
+
+   taux = dx./ds; 
+   tauy = dy./ds;
+  
+   K1 = (-1/(2*zk^2).*(grad(:, :, 1).*(nx) + grad(:, :, 2).*ny)); 
+   K1H = ((1 + nu)/2).*(-1/(2*zk^2).*(grad(:, :, 1).*(taux) + grad(:, :, 2).*tauy));                    % G_{tauy}
+   K2 = 1/(2*zk^2).*val;
+
+   submat = zeros(nt,3*ns);
+   submat(:,1:3:end) = K1;
+   submat(:,2:3:end) = K1H;
+   submat(:,3:3:end) = K2;
+
+end
+
 %%% SUPPORTED PLATE KERNELS
 
 % boundary conditions applied to a point source
@@ -727,6 +787,58 @@
 
 end
 
+% supported plate kernels for far field evaluation
+if strcmpi(type, 'supported_plate_eval_ff')
+    srcnorm = srcinfo.n;
+    srctang = srcinfo.d;
+    srcd2 = srcinfo.d2;
+    coefs = varargin{1};
+    nu = coefs(1);
+    
+    nx = repmat(srcnorm(1,:),nt,1);
+    ny = repmat(srcnorm(2,:),nt,1);
+    
+    dx = repmat(srctang(1,:),nt,1);
+    dy = repmat(srctang(2,:),nt,1);
+    
+    ds = sqrt(dx.*dx+dy.*dy);
+    
+    taux = dx./ds;
+    tauy = dy./ds;
+    
+    dx1 = repmat(srctang(1,:),nt,1);
+    dy1 = repmat(srctang(2,:),nt,1);
+    
+    d2x1 = repmat(srcd2(1,:),nt,1);
+    d2y1 = repmat(srcd2(2,:),nt,1);
+    
+    denom = sqrt(dx1.^2+dy1.^2).^3;
+    numer = dx1.*d2y1-d2x1.*dy1;
+    
+    kappa = numer./denom; 
+    
+    kp = repmat(srcinfo.data(1,:),nt,1);
+    
+    a1 = 2-nu;
+    a2 = (-1+nu)*(7+nu)/(3 - nu);
+    a3 = (1-nu)*(3+nu)/(1+nu);
+    
+    [~, grad, hess, third] = chnk.flex2d.hkdiffgreen_ff(zk, src, targ, false); 
+    
+    K1 = -1/(2*zk^2)*(third(:,:,1).*nx.^3 + 3*third(:,:,2).*nx.^2.*ny + 3*third(:,:,3).*nx.*ny.^2 + third(:,:,4).*ny.^3) + ...
+         -a1/(2*zk^2)*(third(:,:,1).*nx.*taux.^2 + third(:,:,2).*(ny.*taux.^2 + 2*nx.*taux.*tauy) + third(:,:,3).*(nx.*tauy.^2 + 2*ny.*taux.*tauy) + third(:,:,4).*ny.*tauy.^2) + ...
+         a2*kappa./(2*zk^2).*(hess(:,:,1).*nx.^2 + 2*hess(:,:,2).*nx.*ny + hess(:,:,3).*ny.^2) + ...
+         -a3*kp./(2*zk^2).*(grad(:,:,1).*taux + grad(:,:,2).*tauy);
+
+    K2 = -1/(2*zk^2).*(grad(:,:,1).*nx + grad(:,:,2).*ny);
+
+    submat = zeros(nt,2*ns);
+
+    submat(:,1:2:end) = K1;
+    submat(:,2:2:end) = K2;
+
+end
+
 end
 
 
diff --git a/chunkie/demo/demo_clamped_plate_scatter_far_field.m b/chunkie/demo/demo_clamped_plate_scatter_far_field.m
new file mode 100644
index 0000000..4d3b611
--- /dev/null
+++ b/chunkie/demo/demo_clamped_plate_scatter_far_field.m
@@ -0,0 +1,183 @@
+%DEMO_CLAMPED_PLATE_SCATTER_FAR_FIELD
+%
+% Define an exterior scattering problem on a starfish-shaped domain and 
+% solve
+%
+
+% planewave vec
+
+kvec = [8;0];
+zk = norm(kvec);
+nu = 1/3;
+
+
+% discretize domain
+
+cparams = [];
+cparams.eps = 1.0e-6;
+cparams.nover = 0;
+cparams.maxchunklen = 4.0/zk; % setting a chunk length helps when the
+                              % frequency is known
+                              
+pref = []; 
+pref.k = 16;
+narms = 3;
+amp = 0.25;
+start = tic; chnkr = chunkerfunc(@(t) starfish(t,narms,amp),cparams,pref); 
+t1 = toc(start);
+
+fprintf('%5.2e s : time to build geo\n',t1)
+
+% plot geometry and data
+
+figure(1)
+clf
+plot(chnkr,'-x')
+hold on
+quiver(chnkr)
+axis equal
+
+% assembling system matrix
+
+fkern =  @(s,t) chnk.flex2d.kern(zk, s, t, 'clamped_plate');
+
+kappa = signed_curvature(chnkr);
+kappa = kappa(:);
+
+opts = [];
+opts.sing = 'log';
+
+start = tic;
+sys = chunkermat(chnkr,fkern, opts);
+sys = sys - 0.5*eye(2*chnkr.npt);
+sys(2:2:end,1:2:end) = sys(2:2:end,1:2:end) + kappa.*eye(chnkr.npt);
+
+t1 = toc(start);
+fprintf('%5.2e s : time to assemble matrix\n',t1)
+
+% building RHS
+
+[r1, grad] = planewave(kvec, chnkr.r);
+
+nx = chnkr.n(1,:); 
+ny = chnkr.n(2,:);
+
+normalderiv = grad(:, 1).*(nx.')+ grad(:, 2).*(ny.');                                % Dirichlet and Neumann BC(Clamped BC)                         
+
+firstbc = -r1;
+secondbc = -normalderiv;
+
+rhs = zeros(2*chnkr.npt, 1); rhs(1:2:end) = firstbc ; rhs(2:2:end) = secondbc;
+
+% Solving linear system
+
+start = tic; sol = gmres(sys,rhs,[],1e-12,100); t1 = toc(start);
+fprintf('%5.2e s : time for dense gmres\n',t1)    
+
+% evaluate at targets and plot
+
+rmin = min(chnkr); rmax = max(chnkr);
+xl = rmax(1)-rmin(1);
+yl = rmax(2)-rmin(2);
+nplot = 200;
+xtarg = linspace(rmin(1)-xl/2,rmax(1)+xl/2,nplot); 
+ytarg = linspace(rmin(2)-yl/2,rmax(2)+yl/2,nplot);
+[xxtarg,yytarg] = meshgrid(xtarg,ytarg);
+targets = zeros(2,length(xxtarg(:)));
+targets(1,:) = xxtarg(:); targets(2,:) = yytarg(:);
+
+start = tic; in = chunkerinterior(chnkr,{xtarg,ytarg}); t1 = toc(start);
+out = ~in;
+
+fprintf('%5.2e s : time to find points in domain\n',t1)
+
+ikern = @(s,t) chnk.flex2d.kern(zk, s, t, 'clamped_plate_eval');                              % build the kernel of evaluation          
+
+start1 = tic;
+uscat = chunkerkerneval(chnkr, ikern,sol, targets(:,out));
+t2 = toc(start1);
+fprintf('%5.2e s : time for kernel eval (for plotting)\n',t2)
+
+uin = planewave(kvec,targets(:,out));
+utot = uscat(:)+uin(:);
+
+maxu = max(abs(uin(:)));
+
+figure(2)
+clf
+
+t = tiledlayout(1,3,'TileSpacing','compact');
+
+nexttile
+zztarg = nan(size(xxtarg));
+zztarg(out) = uin;
+h=pcolor(xxtarg,yytarg,imag(zztarg),"FaceColor","interp");
+set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(redblue);
+hold on
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+title('$u^{\textrm{inc}}$','Interpreter','latex','FontSize',12)
+
+nexttile
+zztarg = nan(size(xxtarg));
+zztarg(out) = uscat;
+h=pcolor(xxtarg,yytarg,imag(zztarg),"FaceColor","interp");
+set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(redblue);
+hold on
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+title('$u^{\textrm{scat}}$','Interpreter','latex','FontSize',12)
+
+nexttile
+zztarg = nan(size(xxtarg));
+zztarg(out) = utot;
+h=pcolor(xxtarg,yytarg,imag(zztarg),"FaceColor","interp");
+set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(redblue);
+hold on
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+colorbar
+title('$u^{\textrm{tot}}$','Interpreter','latex','FontSize',12)
+title(t,"Clamped Plate BCs")
+
+
+% evaluate the far field representation
+
+ts = linspace(-pi,pi,300);
+ts = ts(1:end-1);
+
+targets_ff = [cos(ts); sin(ts)];
+
+ikern_ff = @(s,t) chnk.flex2d.kern(zk, s, t, 'clamped_plate_eval_ff', nu); 
+
+start1 = tic;
+uscat = chunkerkerneval(chnkr, ikern_ff,sol,targets_ff);
+t2 = toc(start1);
+fprintf('%5.2e s : time for kernel eval (for plotting)\n',t2)
+
+figure(3)
+clf 
+
+t = tiledlayout(1,3,'TileSpacing','compact');
+title(t,'Clamped Plate Far Field Pattern')
+
+nexttile
+plot(ts,real(uscat));
+title('Real')
+
+nexttile
+plot(ts,imag(uscat));
+title('Imag')
+
+nexttile
+plot(ts,abs(uscat));
+title('Abs')
\ No newline at end of file
diff --git a/chunkie/demo/demo_free_plate_scatter_far_field.m b/chunkie/demo/demo_free_plate_scatter_far_field.m
new file mode 100644
index 0000000..5d71d51
--- /dev/null
+++ b/chunkie/demo/demo_free_plate_scatter_far_field.m
@@ -0,0 +1,218 @@
+%DEMO_FREE_PLATE_SCATTER_FAR_FIELD
+%
+% Define an exterior scattering problem on a starfish-shaped domain and 
+% and evaluate the far field pattern
+%
+
+% planewave vec
+
+kvec = [8;0];
+zk = norm(kvec);
+nu = 1/3;
+
+% discretize domain
+
+cparams = [];
+cparams.eps = 1.0e-6;
+cparams.nover = 0;
+cparams.maxchunklen = 4.0/zk; % setting a chunk length helps when the
+                              % frequency is known
+                              
+pref = []; 
+pref.k = 16;
+narms = 3;
+amp = 0.25;
+start = tic; chnkr = chunkerfunc(@(t) starfish(t,narms,amp),cparams,pref); 
+t1 = toc(start);
+
+fprintf('%5.2e s : time to build geo\n',t1)
+
+% plot geometry and data
+
+figure(1)
+clf
+plot(chnkr,'-x')
+hold on
+quiver(chnkr)
+axis equal
+
+% defining free plate kernels
+
+fkern1 =  @(s,t) chnk.flex2d.kern(zk, s, t, 'free_plate', nu);        % build the desired kernel
+double = @(s,t) chnk.lap2d.kern(s,t,'d');
+hilbert = @(s,t) chnk.lap2d.kern(s,t,'hilb');
+
+opts = [];
+opts.sing = 'log';
+
+opts2 = [];
+opts2.sing = 'pv';
+
+% building system matrix
+
+start = tic;
+sysmat1 = chunkermat(chnkr,fkern1, opts);
+D = chunkermat(chnkr, double, opts);
+H = chunkermat(chnkr, hilbert, opts2);     
+
+sysmat = zeros(2*chnkr.npt);
+sysmat(1:2:end,1:2:end) = sysmat1(1:4:end,1:2:end) + sysmat1(3:4:end,1:2:end)*H  - 2*((1+nu)/2)^2*D*D;
+sysmat(2:2:end,1:2:end) = sysmat1(2:4:end,1:2:end) + sysmat1(4:4:end,1:2:end)*H;
+sysmat(1:2:end,2:2:end) = sysmat1(1:4:end,2:2:end) + sysmat1(3:4:end,2:2:end);
+sysmat(2:2:end,2:2:end) = sysmat1(2:4:end,2:2:end) + sysmat1(4:4:end,2:2:end);
+
+D = [-1/2 + (1/8)*(1+nu).^2, 0; 0, 1/2];  % jump matrix 
+D = kron(eye(chnkr.npt), D);
+
+sys =  D + sysmat;
+t1 = toc(start);
+fprintf('%5.2e s : time to assemble matrix\n',t1)
+
+% building RHS
+
+[~, ~, hess, third] = planewave(kvec, chnkr.r);
+
+nx = chnkr.n(1,:).'; 
+ny = chnkr.n(2,:).';
+
+dx = chnkr.d(1,:).';
+dy = chnkr.d(2,:).';
+
+ds = sqrt(dx.*dx+dy.*dy);
+taux = (dx./ds);                                                                       % normalization
+tauy = (dy./ds);
+
+kappa = signed_curvature(chnkr);
+kappa = kappa(:)';
+
+firstbc = -((hess(:,1).*(nx.*nx) + hess(:,2).*(2*nx.*ny) + hess(:,3).*(ny.*ny)) + ...
+           nu.*(hess(:,1).*(taux.*taux) + hess(:,2).*(2*taux.*tauy) + hess(:,3).*(tauy.*tauy)));
+
+secondbc = -((third(:,1).*(nx.*nx.*nx) + third(:,2).*(3*nx.*nx.*ny) +...
+       third(:,3).*(3*nx.*ny.*ny) + third(:,4).*(ny.*ny.*ny))  + ...
+        (2-nu).*(third(:,1).*(taux.*taux.*nx) + third(:,2).*(taux.*taux.*ny + 2*taux.*tauy.*nx) +...
+        third(:,3).*(2*taux.*tauy.*ny+ tauy.*tauy.*nx) +...
+        + third(:,4).*(tauy.*tauy.*ny)) + ...
+        (1-nu).*kappa'.*(hess(:,1).*taux.*taux + hess(:,2).*(2*taux.*tauy) + hess(:,3).*tauy.*tauy+...
+        -(hess(:,1).*nx.*nx + hess(:,2).*(2*nx.*ny) + hess(:,3).*ny.*ny)));
+
+rhs = zeros(2*chnkr.npt, 1); rhs(1:2:end) = firstbc ; rhs(2:2:end) = secondbc;
+
+% Solving linear system
+
+start = tic; sol = gmres(sys,rhs,[],1e-12,500); t1 = toc(start);
+fprintf('%5.2e s : time for dense gmres\n',t1)    
+
+% evaluate at targets and plot
+
+rmin = min(chnkr); rmax = max(chnkr);
+xl = rmax(1)-rmin(1);
+yl = rmax(2)-rmin(2);
+nplot = 200;
+xtarg = linspace(rmin(1)-xl/2,rmax(1)+xl/2,nplot); 
+ytarg = linspace(rmin(2)-yl/2,rmax(2)+yl/2,nplot);
+[xxtarg,yytarg] = meshgrid(xtarg,ytarg);
+targets = zeros(2,length(xxtarg(:)));
+targets(1,:) = xxtarg(:); targets(2,:) = yytarg(:);
+
+start = tic; in = chunkerinterior(chnkr,{xtarg,ytarg}); t1 = toc(start);
+out = ~in;
+
+fprintf('%5.2e s : time to find points in domain\n',t1)
+
+ikern = @(s,t) chnk.flex2d.kern(zk, s, t, 'free_plate_eval', nu); 
+
+dens_comb = zeros(3*chnkr.npt,1);
+dens_comb(1:3:end) = sol(1:2:end);
+dens_comb(2:3:end) = H*sol(1:2:end);
+dens_comb(3:3:end) = sol(2:2:end);
+
+start1 = tic;
+uscat = chunkerkerneval(chnkr, ikern,dens_comb,targets(:,out));
+t2 = toc(start1);
+fprintf('%5.2e s : time for kernel eval (for plotting)\n',t2)
+
+uin = planewave(kvec,targets(:,out));
+utot = uscat(:)+uin(:);
+
+maxu = max(abs(uin(:)));
+
+figure(2)
+clf
+
+t = tiledlayout(1,3,'TileSpacing','compact');
+
+nexttile
+zztarg = nan(size(xxtarg));
+zztarg(out) = uin;
+h=pcolor(xxtarg,yytarg,imag(zztarg),"FaceColor","interp");
+set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(redblue);
+hold on
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+title('$u^{\textrm{inc}}$','Interpreter','latex','FontSize',12)
+
+nexttile
+zztarg = nan(size(xxtarg));
+zztarg(out) = uscat;
+h=pcolor(xxtarg,yytarg,imag(zztarg),"FaceColor","interp");
+set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(redblue);
+hold on
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+title('$u^{\textrm{scat}}$','Interpreter','latex','FontSize',12)
+
+nexttile
+zztarg = nan(size(xxtarg));
+zztarg(out) = utot;
+h=pcolor(xxtarg,yytarg,imag(zztarg),"FaceColor","interp");
+set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(redblue);
+hold on
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+colorbar
+title('$u^{\textrm{tot}}$','Interpreter','latex','FontSize',12)
+title(t,"Free Plate BCs")
+
+
+
+% evaluate the far field representation
+
+ts = linspace(-pi,pi,300);
+ts = ts(1:end-1);
+
+targets_ff = [cos(ts); sin(ts)];
+
+ikern_ff = @(s,t) chnk.flex2d.kern(zk, s, t, 'free_plate_eval_ff', nu); 
+
+start1 = tic;
+uscat = chunkerkerneval(chnkr, ikern_ff,dens_comb,targets_ff);
+t2 = toc(start1);
+fprintf('%5.2e s : time for kernel eval (for plotting)\n',t2)
+
+figure(3)
+clf 
+
+t = tiledlayout(1,3,'TileSpacing','compact');
+title(t,'Free Plate Far Field Pattern')
+
+nexttile
+plot(ts,real(uscat));
+title('Real')
+
+nexttile
+plot(ts,imag(uscat));
+title('Imag')
+
+nexttile
+plot(ts,abs(uscat));
+title('Abs')
\ No newline at end of file
diff --git a/chunkie/demo/demo_supported_plate_scatter_far_field.m b/chunkie/demo/demo_supported_plate_scatter_far_field.m
new file mode 100644
index 0000000..eefbb4d
--- /dev/null
+++ b/chunkie/demo/demo_supported_plate_scatter_far_field.m
@@ -0,0 +1,215 @@
+%DEMO_SUPPORTED_PLATE_SCATTER_FAR_FIELD
+%
+% Define an exterior scattering problem on a starfish-shaped domain and 
+% evaluate the far field pattern
+%
+% this demonstration requires advanced usage in which additional geometric
+% info is passed to the kernels using chunker data. this demo also avoids
+% precision loss in the kernel evaluators by using the "native" quadrature
+% on certain components of the matrix
+
+% planewave vec
+
+kvec = [8;0];
+zk = norm(kvec);
+nu = 1/3;
+
+% discretize domain
+
+cparams = [];
+cparams.eps = 1.0e-6;
+cparams.nover = 0;
+cparams.maxchunklen = 4.0/zk; % setting a chunk length helps when the
+                              % frequency is known
+                              
+pref = []; 
+pref.k = 16;
+narms = 5;
+amp = 0.25;
+start = tic; chnkr = chunkerfunc(@(t) starfish(t,narms,amp),cparams,pref); 
+chnkr = makedatarows(chnkr,2);
+t1 = toc(start);
+
+fprintf('%5.2e s : time to build geo\n',t1)
+
+% plot geometry and data
+
+figure(1)
+clf
+plot(chnkr,'-x')
+hold on
+quiver(chnkr)
+axis equal
+
+% calculating curvature info
+
+kappa = signed_curvature(chnkr);
+kp = arclengthder(chnkr,kappa);
+kpp = arclengthder(chnkr,kp);
+
+% supported plate kernels expect (d/ds) kappa in the first data row
+% and (d^2/ds^2) kappa in the second data row
+
+chnkr.data(1,:,:) = kp;
+chnkr.data(2,:,:) = kpp;
+
+% defining supported plate kernels
+
+fkern1 =  @(s,t) chnk.flex2d.kern(zk, s, t, 'supported_plate_log',nu);           % build the desired kernel
+fkern2 =  @(s,t) chnk.flex2d.kern(zk, s, t, 'supported_plate_smooth',nu);           % build the desired kernel
+
+opts = [];
+opts.sing = 'log';
+
+opts2 = [];
+opts2.quad = 'native';
+opts2.sing = 'smooth';
+
+% building system matrix
+
+start = tic;
+M = chunkermat(chnkr,fkern1, opts);
+M2 = chunkermat(chnkr,fkern2, opts2);
+
+c0 = (nu - 1)*(nu + 3)*(2*nu - 1)/(2*(3 - nu));
+
+M(2:2:end,1:2:end) = M(2:2:end,1:2:end) + M2 + c0.*kappa(:).^2.*eye(chnkr.npt);
+M = M - 0.5*eye(2*chnkr.npt);
+
+sys =  M;
+t1 = toc(start);
+fprintf('%5.2e s : time to assemble matrix\n',t1)
+
+% building RHS
+
+[val, ~, hess] = planewave(kvec, chnkr.r);
+
+nx = chnkr.n(1,:).'; 
+ny = chnkr.n(2,:).';
+
+dx = chnkr.d(1,:).';
+dy = chnkr.d(2,:).';
+
+ds = sqrt(dx.*dx+dy.*dy);
+taux = (dx./ds);                                                                       % normalization
+tauy = (dy./ds);
+
+firstbc = -val ;
+
+secondbc = -((hess(:, 1).*(nx.*nx) + hess(:, 2).*(2*nx.*ny) + hess(:, 3).*(ny.*ny))+...
+           nu*(hess(:, 1).*(taux.*taux) + hess(:, 2).*(2*taux.*tauy) + hess(:, 3).*(tauy.*tauy)));
+    
+rhs = zeros(2*chnkr.npt, 1); rhs(1:2:end) = firstbc ; rhs(2:2:end) = secondbc;
+
+% Solving linear system
+
+start = tic; sol = gmres(sys,rhs,[],1e-12,500); t1 = toc(start);
+fprintf('%5.2e s : time for dense gmres\n',t1)    
+
+% evaluate at targets and plot
+
+rmin = min(chnkr); rmax = max(chnkr);
+xl = rmax(1)-rmin(1);
+yl = rmax(2)-rmin(2);
+nplot = 200;
+xtarg = linspace(rmin(1)-xl/2,rmax(1)+xl/2,nplot); 
+ytarg = linspace(rmin(2)-yl/2,rmax(2)+yl/2,nplot);
+[xxtarg,yytarg] = meshgrid(xtarg,ytarg);
+targets = zeros(2,length(xxtarg(:)));
+targets(1,:) = xxtarg(:); targets(2,:) = yytarg(:);
+
+start = tic; in = chunkerinterior(chnkr,{xtarg,ytarg}); t1 = toc(start);
+out = ~in;
+
+fprintf('%5.2e s : time to find points in domain\n',t1)
+
+ikern = @(s,t) chnk.flex2d.kern(zk, s, t, 'supported_plate_eval', nu); 
+
+start1 = tic;
+uscat = chunkerkerneval(chnkr,ikern,sol,targets(:,out));
+t2 = toc(start1);
+fprintf('%5.2e s : time for kernel eval (for plotting)\n',t2)
+
+uin = planewave(kvec,targets(:,out));
+utot = uscat(:)+uin(:);
+
+maxu = max(abs(uin(:)));
+
+figure(2)
+clf
+
+t = tiledlayout(1,3,'TileSpacing','compact');
+
+nexttile
+zztarg = nan(size(xxtarg));
+zztarg(out) = uin;
+h=pcolor(xxtarg,yytarg,imag(zztarg),"FaceColor","interp");
+set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(redblue);
+hold on
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+title('$u^{\textrm{inc}}$','Interpreter','latex','FontSize',12)
+
+nexttile
+zztarg = nan(size(xxtarg));
+zztarg(out) = uscat;
+h=pcolor(xxtarg,yytarg,imag(zztarg),"FaceColor","interp");
+set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(redblue);
+hold on
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+title('$u^{\textrm{scat}}$','Interpreter','latex','FontSize',12)
+
+nexttile
+zztarg = nan(size(xxtarg));
+zztarg(out) = utot;
+h=pcolor(xxtarg,yytarg,imag(zztarg),"FaceColor","interp");
+set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(redblue);
+hold on
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+colorbar
+title('$u^{\textrm{tot}}$','Interpreter','latex','FontSize',12)
+title(t,"Supported Plate BCs")
+
+
+% evaluate the far field representation
+
+ts = linspace(-pi,pi,300);
+ts = ts(1:end-1);
+
+targets_ff = [cos(ts); sin(ts)];
+
+ikern_ff = @(s,t) chnk.flex2d.kern(zk, s, t, 'supported_plate_eval_ff', nu); 
+
+start1 = tic;
+uscat = chunkerkerneval(chnkr, ikern_ff,sol,targets_ff);
+t2 = toc(start1);
+fprintf('%5.2e s : time for kernel eval (for plotting)\n',t2)
+
+figure(3)
+clf 
+
+t = tiledlayout(1,3,'TileSpacing','compact');
+title(t,'Supported Plate Far Field Pattern')
+
+nexttile
+plot(ts,real(uscat));
+title('Real')
+
+nexttile
+plot(ts,imag(uscat));
+title('Imag')
+
+nexttile
+plot(ts,abs(uscat));
+title('Abs')
\ No newline at end of file