diff --git a/chunkie/+chnk/chunkerkerneval_smooth.m b/chunkie/+chnk/chunkerkerneval_smooth.m index 5dd487e1..ff611b19 100644 --- a/chunkie/+chnk/chunkerkerneval_smooth.m +++ b/chunkie/+chnk/chunkerkerneval_smooth.m @@ -62,7 +62,7 @@ % these are automatically ignored by the fmm and % in chunkermat corrections if ~(strcmpi(imethod,'fmm') && isempty(flag)) - flagslf = chnk.flagself(targinfo.r, chnkr.r, 1e-14*diam); + flagslf = chnk.flagself(real(targinfo.r), real(chnkr.r), 1e-14*diam); if isempty(flagslf) selfzero = sparse(opdims(1)*size(targinfo.r(:,:),2), ... opdims(2)*chnkr.npt); @@ -141,7 +141,7 @@ wts = wts(:); if strcmpi(imethod,'flam') - xflam1 = chnkr.r(:,:); + xflam1 = real(chnkr.r(:,:)); xflam1 = repmat(xflam1,opdims(2),1); xflam1 = reshape(xflam1,chnkr.dim,numel(xflam1)/chnkr.dim); diff --git a/chunkie/@chunker/flagnear.m b/chunkie/@chunker/flagnear.m index 5a4ee8e1..a8b243f8 100644 --- a/chunkie/@chunker/flagnear.m +++ b/chunkie/@chunker/flagnear.m @@ -68,8 +68,8 @@ % get all pairwise distances itarg = inbor(inbor > npt)-npt; - srcs = reshape(chnkr.rstor(:,isrc),dim,1,length(isrc)); - targs = reshape(pts(:,itarg),dim,length(itarg),1); + srcs = reshape(real(chnkr.rstor(:,isrc)),dim,1,length(isrc)); + targs = reshape(real(pts(:,itarg)),dim,length(itarg),1); dists = reshape(sqrt(sum((bsxfun(@minus,targs,srcs)).^2,1)), ... length(itarg),length(isrc)); diff --git a/chunkie/@chunkgraph/chunkgraph.m b/chunkie/@chunkgraph/chunkgraph.m index d719bd10..3ba6cc72 100644 --- a/chunkie/@chunkgraph/chunkgraph.m +++ b/chunkie/@chunkgraph/chunkgraph.m @@ -304,7 +304,7 @@ end vs = chunkends(chnkr,[1,chnkr.nch]); - vs = vs(:, [1,4]); + vs = real(vs(:, [1,4])); else msg = "CHUNKGRAPH:CONSTRUCTOR: invalid edge descriptor." + ... diff --git a/chunkie/@chunkgraph/findunbounded.m b/chunkie/@chunkgraph/findunbounded.m index 64765d0a..90f5c7d8 100644 --- a/chunkie/@chunkgraph/findunbounded.m +++ b/chunkie/@chunkgraph/findunbounded.m @@ -36,8 +36,8 @@ for jj =1:numel(edges) echnk = cgrph.echnks(abs(edges(jj))); [allrends,alltends] = chunkends(echnk); - ts1 = alltends(:,1,:); - ts2 = alltends(:,2,:); + ts1 = real(alltends(:,1,:)); + ts2 = real(alltends(:,2,:)); angs1 = atan2(ts1(2,:),ts1(1,:)); angs2 = atan2(ts2(2,:),ts2(1,:)); angdiffs = angs2-angs1; diff --git a/chunkie/@chunkgraph/vertextract.m b/chunkie/@chunkgraph/vertextract.m index b4c5e569..558d5719 100644 --- a/chunkie/@chunkgraph/vertextract.m +++ b/chunkie/@chunkgraph/vertextract.m @@ -56,7 +56,7 @@ end % compute the corresponding angles and sort them -ds = [deplus,deminus]; +ds = real([deplus,deminus]); angs = atan2(ds(2,:),ds(1,:)); [angs,isrtededges] = sort(angs); diff --git a/chunkie/chunkerkerneval.m b/chunkie/chunkerkerneval.m index 246f2794..fdcc25a8 100644 --- a/chunkie/chunkerkerneval.m +++ b/chunkie/chunkerkerneval.m @@ -162,6 +162,14 @@ error("CHUNKERKERNEVAL: input 4 is not a supported type"); end +if ~isreal(chnkobj.r) && (opts_use.accel || opts_use.forcefmm) + warning('fmm not supported for complex chunkers') +end + +if ~isreal(targinfo.r) && (opts_use.accel || opts_use.forcefmm || opts_use.flam) + warning('accel not supported for complex targets') +end + if icgrph && mk > 1 ids = chunkgraphinregion(chnkobj,targinfo.r); else @@ -282,7 +290,9 @@ if opts_use.forcepquad optsflag = []; optsflag.fac = opts_use.fac; - flag = flagnear(chnkr0,targinfo.r,optsflag); + chnkrflag = chunkerpoints(struct('r',real(chnkr0.r), ... + 'd',real(chnkr0.d),'d2',real(chnkr0.d2))); + flag = flagnear(chnkrflag,targinfo.r,optsflag); fints(irow0) = fints(irow0) + chnk.chunkerkerneval_smooth(chnkr0,kern0,opdims0,dens0,targinfo0, ... flag,opts_use); @@ -301,7 +311,9 @@ rho = 1.8; optsflag = []; optsflag.rho = rho; -flag = flagnear_rectangle(chnkr0,targinfo0.r,optsflag); +chnkrflag = chunkerpoints(struct('r',real(chnkr0.r), ... + 'd',real(chnkr0.d),'d2',real(chnkr0.d2))); +flag = flagnear_rectangle(chnkrflag,targinfo0.r,optsflag); npoly = chnkr0.k*2; nlegnew = chnk.ellipse_oversample(rho,npoly,opts_use.eps); diff --git a/chunkie/chunkerkernevalmat.m b/chunkie/chunkerkernevalmat.m index 2ec293f3..1b442d9b 100644 --- a/chunkie/chunkerkernevalmat.m +++ b/chunkie/chunkerkernevalmat.m @@ -212,7 +212,7 @@ optsadap = []; optsadap.eps = eps; -if corrections +if corrections || nonsmoothonly mat = sparse(nout,icollocs(end)-1); else mat = zeros(nout,icollocs(end)-1); @@ -267,7 +267,9 @@ % TODO: switch to oversampling version with bounding ellipse like chunkerkerneval optsflag = []; optsflag.fac = fac; -flag = flagnear(chnkr0,targinfo0.r,optsflag); +chnkrflag = chunkerpoints(struct('r',real(chnkr0.r), ... + 'd',real(chnkr0.d),'d2',real(chnkr0.d2))); +flag = flagnear(chnkrflag,real(targinfo0.r),optsflag); if forcepquad spmat = chunkerkernevalmat_pquad(chnkr0,kern0,opdims0, ... diff --git a/devtools/test/complexificationTest.m b/devtools/test/complexificationTest.m index 276fb0d3..c020a844 100644 --- a/devtools/test/complexificationTest.m +++ b/devtools/test/complexificationTest.m @@ -16,8 +16,6 @@ function complexificationTest0() zk1 = zks(1); zk2 = zks(2); -close all - a = 3; b = 1/a/2; t0 =-5; @@ -34,6 +32,11 @@ function complexificationTest0() cparams.ifclosed = 0; chnkr = chunkerfuncuni(f, nch, cparams); +rends = chunkends(chnkr,[1,chnkr.nch]); +verts = rends(:,[1,4]); +fchnk{1} = chnkr; +cgrph = chunkgraph(real(verts),[1;2],fchnk); + ddiff = kernel('helmdiff', 'd', zks); sdiff = (-1)*kernel('helmdiff', 's', zks); @@ -43,13 +46,13 @@ function complexificationTest0() K = kernel([ddiff, sdiff; ... dpdiff, spdiff]); -n = 2*chnkr.npt; -sysmat = chunkermat(chnkr, K); +n = 2*cgrph.npt; +sysmat = chunkermat(cgrph, K); sysmat = sysmat - eye(n); %% Analytic solution test -r = chnkr.r; +r = cgrph.r; r0 = [0.1; -3]; s = []; s.r = r0; @@ -58,27 +61,23 @@ function complexificationTest0() dk1 = kernel('helm', 'd', zk1); skp1 = kernel('helm', 'sp', zk1); rhs = complex(zeros(n,1)); -rhs(1:2:end) = sk1.eval(s, chnkr); -rhs(2:2:end) = skp1.eval(s, chnkr); +rhs(1:2:end) = sk1.eval(s, cgrph); +rhs(2:2:end) = skp1.eval(s, cgrph); dens = sysmat \ rhs; -rt = [0.3; 2]; +rt = [[0.3; 2] [-0.5;0.001]]; targ = []; targ.r = rt; Keval = kernel([dk1, (-1)*sk1]); opts = []; -opts.forcesmooth = true; +% opts.forcesmooth = true; opts.accel = false; -pot = chunkerkerneval(chnkr, Keval, dens, targ, opts); +pot = chunkerkerneval(cgrph, Keval, dens, targ, opts); uex = sk1.eval(s, targ); -assert(norm(uex - pot) < 1e-10); - - - - +assert(norm(uex - pot) < 1e-12); end