Merge branch 'master' of https://github.com/fastalgorithms/chunkie in…

…to dev-chunkerkerneval

.readthedocs.yaml

@@ -0,0 +1,32 @@
+# .readthedocs.yaml
+# Read the Docs configuration file
+# See https://docs.readthedocs.io/en/stable/config-file/v2.html for details
+
+# Required
+version: 2
+
+# Set the OS, Python version and other tools you might need
+build:
+  os: ubuntu-22.04
+  tools:
+    python: "3.12"
+    # You can also specify other tool versions:
+    # nodejs: "19"
+    # rust: "1.64"
+    # golang: "1.19"
+
+# Build documentation in the "docs/" directory with Sphinx
+sphinx:
+   configuration: docs/conf.py
+
+# Optionally build your docs in additional formats such as PDF and ePub
+# formats:
+#    - pdf
+#    - epub
+
+# Optional but recommended, declare the Python requirements required
+# to build your documentation
+# See https://docs.readthedocs.io/en/stable/guides/reproducible-builds.html
+python:
+  install:
+    - requirements: docs/requirements.txt

chunkie/+chnk/+helm2d/transmission_helper.m

@@ -151,7 +151,7 @@
                 c1 = coefs(d1);
                 c2 = coefs(d2);
                 ind1 = sum(nchs(1:i-1))*ngl*2+(1:2:2*nchs(i)*ngl);
-                ind2 = sum(nchs(1:i-1))*ngl*2+(1:2:2*nchs(i)*ngl);
+                ind2 = sum(nchs(1:i-1))*ngl*2+(2:2:2*nchs(i)*ngl);
                 targnorm = chnkrs(i).n;
                 nx = targnorm(1,:); nx = nx(:);
                 ny = targnorm(2,:); ny = ny(:);

chunkie/@chunker/chunker.m

@@ -1,8 +1,73 @@
 classdef chunker
-%CHUNKER class which describes a curve divided into chunks (panels). 
+%CHUNKER class which describes a curve divided into chunks (or "panels"). 
+%
 % On each chunk the curve is represented by the values of its position, 
-% first and second derivatives in parameter space on a Legendre grid.
+% first and second derivatives by scaled Legendre nodes.
+%
+% chunker properties:
+%   k - integer, number of Legendre nodes on each chunk
+%   nch - integer, number of chunks that make up the curve
+%   dim - integer, dimension of the ambient space in which the curve is 
+%             embedded
+%   npt - returns k*nch, the total number of points on the curve
+%   r - dim x k x nch array, r(:,i,j) gives the coordinates of the ith 
+%         node on the jth chunk of the chunker
+%   h - nch array of scaling factors for chunks. the chunk derivatives are
+%         scaled as if the coordinates in r(:,:,j) are sampled on an 
+%         interval of length 2*h(j). This affects integration formulas.
+%   d - dim x k x nch array, d(:,i,j) gives the time derivative of the 
+%         coordinate at the ith node on the jth chunk of the chunker
+%   d2 - dim x k x nch array, d(:,i,j) gives the 2nd time derivative of the 
+%         coordinate at the ith node on the jth chunk of the chunker
+%   n - dim x k x nch array of normals to the curve
+%   data - datadim x k x nch array of data attached to chunker points 
+%         this data will be refined along with the chunker
+%   adj - 2 x nch integer array. adj(1,j) is i.d. of the chunk that 
+%         precedes chunk j in parameter space. adj(2,j) is the i.d. of the
+%         chunk which follows. if adj(i,j) is 0 then that end of the chunk
+%         is a free end. if adj(i,j) is negative than that end of the chunk
+%         meets with other chunk ends in a vertex. the specific negative
+%         number acts as an i.d. of that vertex.
 %
+% chunker methods: 
+%   chunker(p,t,w) - construct an empty chunker with given preferences and
+%       precomputed Legendre nodes/weights (optional)
+%   obj = addchunk(obj,nchadd) - add nchadd chunks to the structure 
+%       (initialized with zeros)
+%   obj = makedatarows(obj,nrows) - add nrows rows to the data storage.
+%   [obj,info] = sort(obj) - sort the chunks so that adjacent chunks are
+%        stored sequentially
+%   [rn,dn,d2n,dist,tn,ichn] = nearest(obj,ref,ich,opts,u,xover,aover) -
+%        find nearest point on chunker to ref
+%   obj = reverse(obj) - reverse chunk orientation
+%   rmin = min(obj) - minimum of coordinate values
+%   rmax = max(obj) - maximum of coordinate values
+%   wts = weights(obj) - scaled integration weights on curve
+%   obj.n = normals(obj) - recompute normal vectors
+%   onesmat = onesmat(obj) - matrix that corresponds to integration of a
+%             density on the curve
+%   rnormonesmat = normonesmat(obj) - matrix that corresponds to
+%             integration of dot product of normal with vector density on
+%             the curve
+%   plot(obj,varargin) - plot the chunker curve
+%   plot3(obj,idata,varargin) - 3D plot of the curve and one row of the
+%             data storage
+%   quiver(obj,varargin) - quiver plot of the chnkr points and normals
+%   scatter(obj,varargin) - scatter plot of the chnkr nodes
+%   tau = taus(obj) - unit tangents to curve
+%   obj = refine(obj,varargin) - refine the curve
+%   a = area(obj) - for a closed curve, area inside 
+%   s = arclength(obj) - get values of arclength along curve
+%   kappa = signed_curvature(obj) - get signed curvature along curve
+%   rflag = datares(obj,opts) - check if data in data rows is resolved
+%   [rc,dc,d2c] = exps(obj) - get expansion coefficients for r,d,d2
+%   ier = checkadjinfo(obj) - checks that the adjacency info of the curve
+%              is consistent
+%   [inds,adjs,info] = sortinfo(obj) - attempts to sort the curve and finds
+%              number of connected components, etc
+%   [re,taue] = chunkends(obj,ich) - get the endpoints of chunks
+%   flag = flagnear(obj,pts,opts) - flag points near the boundary
+
 
 % author: Travis Askham ([email protected])
 
@@ -49,6 +114,16 @@
 
     methods
         function obj = chunker(p,t,w)
+        %CHUNKER construct an empty chunker with some default settings
+        %
+        % syntax: chnkr = chunker(p,t,w);
+        %
+        % optional input: 
+        %   p - struct of preferences 
+        %        p.k - integer, order to be used on chunks (16)
+        %        p.dim - integer, dimension of ambient space (2)
+        %   t, w - arrays of Legendre nodes and weights of order p.k
+        % 
             if nargin < 1 || isempty(p)
                 p = chunkerpref();
             else

chunkie/@chunkgraph/chunkgraph.m

@@ -131,6 +131,72 @@
                 regions{i} = region_comp;
             end
 
+            gmat = zeros(numel(regions));
+
+            for ii=1:numel(regions)
+               rgna = regions{ii};
+               ilist = [];
+               for jj=1:numel(regions)
+                   if (ii ~=jj)
+                        rgnb = regions{jj};
+                        [inc] = regioninside(obj,rgnb,rgna);
+                        if (inc)
+                            ilist = [ilist,jj];
+                        end
+                   end
+                   gmat(ii,ilist) = 1;
+                   gmat(ilist,ii) = 1;
+               end
+               imin = min(ilist);   
+            end    
+
+            ccomp_reg = conncomp(graph(gmat));
+            [s,inds] = sort(ccomp_reg);
+            regions = regions(inds);
+
+            for ii = 1:numel(s)
+                si = s(ii);
+                for jj=1:(numel(s)-1)
+                    sj = s(jj);
+                    if (si == sj)
+                        rgna = regions{jj};
+                        rgnb = regions{jj+1};
+                        [inc] = regioninside(obj,rgna,rgnb);
+                        if (inc)
+                            regions([jj,jj+1])= regions([jj+1,jj]);
+                        end
+                    end
+                end    
+            end
+
+
+            rgns = regions;
+            rgnso= {};
+
+            for ii=1:max(s)
+                inds = find(s==ii);
+                rgnout = rgns{inds(1)};
+                for jj=2:numel(inds)
+                    indj = inds(jj);
+                    [rgnout] = mergeregions(obj,rgnout,rgns{indj});
+                end
+                rgnso{ii} = rgnout;
+            end    
+
+            regions = rgnso;
+            rgns = regions;
+            rgnout = rgns{1};
+            if (numel(rgns)>1)
+                rgn2 = rgns{2};
+                [rgnout] = mergeregions(obj,rgnout,rgn2);
+                for ii=3:numel(rgns)
+                    rgn2 = rgns{ii};
+                    [rgnout] = mergeregions(obj,rgnout,rgn2);
+                end
+            end
+
+            regions = rgnout;
+
             obj.regions = regions;
         end
         function obj = set.verts(obj,val)

chunkie/chunkerfunc.m

@@ -13,7 +13,7 @@
 %               [r,d] = fcurve(t);  or [r,d,d2] = fcurve(t);
 %            where d is the first derivative of r with respect to t and 
 %            d2 is the second derivative. in some situations, this will
-%            improve the convergence order. 
+%            improve the convergence order and final precision.
 %
 % Optional input:
 %	cparams - curve parameters structure (defaults)
@@ -25,9 +25,8 @@
 %           ifclosed == 0 (Inf)
 %       cparams.nover = oversample resolved curve nover
 %           times (0)
-%       cparams.eps = resolve coordinates, arclength,
-%          and first and second derivs of coordinates
-%          to this tolerance (1.0e-6)
+%       cparams.eps = tolerance to resolve coordinates and arclength 
+%           density (1.0e-6)
 %       cparams.lvlr = string, determines type of level
 %          restriction to be enforced
 %               lvlr = 'a' -> no chunk should have double the arc length 
@@ -42,6 +41,9 @@
 %       pref.nchmax - maximum number of chunks (10000)
 %       pref.k - number of Legendre nodes on chunks (16)
 %
+% Output:
+%   chnkr - a chunker object containing the discretization of the domain
+%
 % Examples:
 %   chnkr = chunkerfunc(@(t) starfish(t)); % chunk up starfish w/ standard
 %                                        % options

chunkie/demo/addpaths_loc.m

@@ -2,4 +2,4 @@ function addpaths_loc()
 
 addpath('../');
 addpath(genpath('../FLAM'));
-addpath('../fmm2d/matlab/')
+addpath('../fmm2d/matlab/');

chunkie/demo/demo_scatter.m

@@ -11,7 +11,7 @@
 
 % planewave vec
 
-kvec = 20*[1;-1.5];
+kvec = 10*[1;-1.5];
 
 %
 
@@ -101,42 +101,46 @@
 
 figure(2)
 clf
-subplot(1,3,1)
+
+t = tiledlayout(1,3,'TileSpacing','compact');
+
+nexttile
 zztarg = nan(size(xxtarg));
 zztarg(out) = uin;
 h=pcolor(xxtarg,yytarg,imag(zztarg));
 set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(brewermap([],'RdBu'));
 hold on
-plot(chnkr,'LineWidth',2)
-axis equal
-axis tight
-colormap(redblue)
-caxis([-maxu,maxu])
-title('$u_{in}$','Interpreter','latex','FontSize',24)
-
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+title('$u^{\textrm{inc}}$','Interpreter','latex','FontSize',12)
 
-subplot(1,3,2)
+nexttile
 zztarg = nan(size(xxtarg));
 zztarg(out) = uscat;
 h=pcolor(xxtarg,yytarg,imag(zztarg));
 set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(brewermap([],'RdBu'));
 hold on
-plot(chnkr,'LineWidth',2)
-axis equal
-axis tight
-colormap(redblue)
-caxis([-maxu,maxu])
-title('$u_{scat}$','Interpreter','latex','FontSize',24)
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+
+title('$u^{\textrm{scat}}$','Interpreter','latex','FontSize',12)
 
-subplot(1,3,3)
+nexttile
 zztarg = nan(size(xxtarg));
 zztarg(out) = utot;
 h=pcolor(xxtarg,yytarg,imag(zztarg));
 set(h,'EdgeColor','none')
+clim([-maxu,maxu])
+colormap(brewermap([],'RdBu'));
 hold on
-plot(chnkr,'LineWidth',2)
-axis equal
-axis tight
-colormap(redblue)
-caxis([-maxu,maxu])
-title('$u_{tot}$','Interpreter','latex','FontSize',24)
+plot(chnkr,'k','LineWidth',2)
+axis equal tight
+set(gca, "box","off","Xtick",[],"Ytick",[]);
+
+title('$u^{\textrm{tot}}$','Interpreter','latex','FontSize',12)

chunkie/guide/addpaths_loc.m

@@ -0,0 +1,4 @@
+function addpaths_loc()
+
+addpath('../');
+addpath(genpath('../FLAM'));