acca-imperial / skprime Goto Github PK

As in the title. If the derivative is evaluated at a point z which is in the fundamental region but not in the circular domain, the result is incorrect.
differror.zip

For many reasons (design, debugging, testing e.t.c.) having simply connected support readily available would be very useful due to the ease of comparison with analytical solutions.

Matlab is much better at searching for class information using the @classdir format.

Consider

w = skprime(alpha, dv, qv);
wi = invParam(w);

Using the product formula, one can check that wi.hat(z) does not produce the correct values.

The skprime class needs a warning comment in the class definition around where the properties are defined to remind us that there is a protected copyProperties function that blindly copies all properties based on the meta-class list.

The check for an evaluation point being in the unit disk is given in skprime.evalLogXHat by abs(z) < 1 + eps(2). Should this be a <=?

Function skpDomain/ginput complains about circle intersection when called.

This is due to the new restriction on G0. See #51.

Title pretty much says it all.

A lot of duplicate code in the Green's function classes.

This is still incorrect for points abs(z) > 1.

This either needs the Gj derivatives, or a note to transform the unit domain so the parameter is not at the origin any longer.

Random test failure for unit parameter. (Seen so far with offset annulus and zero on boundary cases.) Randomness related to unit test placing unit parameter at unpredictable location every test. Need to point out what angle causes this.

Currently skprime.diff (and soon bvpFun.diff) is a function handle. Should this be an independent object?

As in title.

Classes should provide other than default behaviour for disp.

We can just use simple FD to test these.

Method skpDomain.isclose() fails on unknown variable when parameter is unit.

Current test is for circles on the order of 1e-4. How small can we get them?

Check against product formula with truncation at L=8. Use this setup:

dv = [-0.2517 + 0.3129i; 0.2307 - 0.4667i];
qv = [0.2377; 0.1557]
da = 0.01;
aarg = 1.67816099467824;
alpha = dv(1) + (qv(1) + da)*exp(1i*aarg);

Accuracy in this case for G_0 is worse, at least near some boundaries, than if, say, aarg = 0.

Most applied problems will not want the solution outside the unit disk. The extra solver step for the domain outside the unit disk adds more time in constructor-intensive applications. The solver for the outer domain, if needed by the geometry, should only be called when there is a request for values at points outside the unit disk.

Currently the software checks for FMM the first time it needs to and saves that result in a persistent variable. If FMM is removed after the first continued point evaluations, failure occurs on subsequent continued point evaluations. Should this remain persistent?

Inaccurate solutions (c.f. product formula) when parameter is near outer boundary.

Since we solve 2 BVP's when alpha is in the fundamental domain, this issue affects all cases of alpha near a boundary. If alpha is near an inner boundary, then Xhat for points outside the unit disk are inaccurate (since 1/conj(alpha) is near an inner boundary). If alpha is near an outer boundary, then Xhat for points in the unit disk are inaccurate, while Xhat for the outer points is ok.

Only give reflected boundary points in skpDomain.boundaryPts if asked. Currently one must delete points if only the bounded domain is wanted, which is the case most of the time.

Do all files have helpful help text?

In some domains at some boundary points greensC0 returns incorrect values. In the below code there is an error of 1 at two boundary points. In all cases in which I've seen errors, the error has always been 1 (but I'm not sure if that's true in general). May also apply to greensCJ.

dv = [
    -0.4-0.5i
    0.5+0.1i];
qv = [
    0.1
    0.15];

D = skpDomain(dv, qv);

zp = boundaryPts(D, 10);
alpha=0.2+0.0i;

g0 = greensC0(alpha, D);

gf=g0(zp);

w1 = skprime(alpha, D);
w2 = invParam(w1);

g0a = @(z) 1./(2i*pi).*log(w1(z)./(abs(alpha).*w2(z)));

ga=g0a(zp);

disp(ga-gf)

It appears greensC0 (and greensCj) are having the normalisation constant applied only when evaluating the full function, but not the "hat" function. The normalisation should be applied to the hat function always, since this is the BVP solution.

TL;DR: The "hat" function gives the wrong values.

Incorporate Wager's phase plot stuff under the +SKP package and add plot method to the skprime class. Will need to add some plot functions for circle domains.

Add vertcat and friends to skpParameter class. Then we can do things like [alpha; alpha - dv], etc.

Physically it makes no sense for the variable or the parameter to be outside the unit domain.

Restricting G0 means we will need to rewrite the prime function BVP for the parameter outside the unit domain. See #52.

When abs(z)>1 and abs(alpha)>1 the following can be used to obtain correct results:

dv = [
  0.052448+0.36539i
  -0.27972-0.12762i
   0.48252-0.28147i];
qv = [
  0.15197
  0.17955
  0.20956];
D = skpDomain(dv, qv);

alpha = 0.50842+0.54632i;
alpha = 1./conj(alpha);

zin = -0.36035+0.32187i;
zin = 1./conj(zin);

%%
% Analytic version.

w = skprime(zin, D);
dw = diff(w);

dag0a = 1./(2i*pi)*(dw(alpha)/w(alpha)-0.5/alpha); 

disp(dag0a)

%%
% BVP version of derivative.

dag0 = greensC0Dp(alpha, D);

disp(dag0(zin));

%%
% Modified BVP version.

dag0 = greensC0Dp(1./conj(alpha), D);

F = conj(dag0(1./conj(zin)));
F = (2i*pi).*F./alpha.^2+0.5./alpha;
dag0=1./(2i*pi).*(F-0.5./alpha);

disp(dag0);

Thus, when the case abs(z)>1 and abs(alpha)<1 is fixed, all cases will be catered for.

The function diff(g0) seems to now be giving erroneous results? (This wasn't previously the case).

In the below example, diff(g0) is the 'odd one out'.

dv = [
  0.052448+0.36539i
  -0.27972-0.12762i
   0.48252-0.28147i];
qv = [
  0.15197
  0.17955
  0.20956];
D = skpDomain(dv, qv);

alpha = 0.50842+0.54632i;
% alpha = 1./conj(alpha);

zin = -0.36035+0.32187i;

%%
% v1

w = skprime(alpha, D);
dw = diff(w);

wi = skprime(1./conj(alpha), D);
dwi = diff(wi);

dg0a = 1./(2i.*pi).*(dw(zin)./w(zin)-dwi(zin)./wi(zin)); 

disp(dg0a)

%%
% v2

g0 = greensC0(alpha, D);

dg0 = diff(g0);

disp(dg0(zin));

%% 
% v3

h = 1e-6;

g0a = @(z) 1./(2i*pi)*log(w(z)./(abs(alpha).*wi(z)));

disp((g0a(zin+h)-g0a(zin))/h)

When evaluating the function greensC0 at specific locations the 'undetermined real constant' from the Schwarz problem is not evaluated correctly. Solutions using this function differ from those in which G0 is formed using the Schottky Klein Prime function by a real constant (which for practical applications should not be the case).

Probably also applied to GJ (but have not checked)
g0constant.zip

.

Formulation needed for these.

If the third argument is a bvpFun object, the constructor fails.

Class constructor signature has changed since help text was last updated.

How about this: Let the derivative function signature be diff(w, n) where w is the prime object and n is the order of the derivative. In other words, can we get away with a recursive implementation of this?

Derivatives w.r.t. to the variable and the parameter are needed for many practical applications.

Allow skprime to take a vector of centres and radii to construct skpDomain automatically. Not everyone wants to set one of those up every time.

See title.

Promote derivative function skprime.diff to parent class bvpFun. Should allow spectral derivatives for the first-kind integrals and the Green's function.

There is evidence the skprimeconj class does not give the same quality output as skprime. I'll attach the script that shows this as soon as I find it ...

probably related to #39 .
When abs(alpha) > 1 the real part of greensC0(zeta) is incorrect.
If abs(alpha) <= 1 the result is correct for zeta anywhere in the fundamental domain.

dv = [
  0.052448+0.36539i
  -0.27972-0.12762i
   0.48252-0.28147i];
qv = [
  0.15197
  0.17955
  0.20956];
D = skpDomain(dv, qv);

alpha = 0.50842+0.54632i;
alpha = 1./conj(alpha);

zin = -0.36035+0.32187i;
% zin = 1./conj(zin);

%%

g0 = greensC0(alpha, D);

disp(g0(zin));

%% 

w = skprime(alpha, D);
wi = skprime(1./conj(alpha), D);

g0a = @(z) 1./(2i*pi)*log(w(z)./(abs(alpha).*wi(z)));

disp(g0a(zin))

For parameter a, if 1/conj(a) is "near enough" an inner circle, then the solution to the Schwarz problem for G_0 and the prime function are unstable.

The fix for this, using G_j for z on C_j, will not work until #9 is resolved.

Current unit tests for Greens functions aren't complete or helpful. Should check points and parameter in fundamental domain. Should check derivatives.

Is it currently possible to evaluate 2nd and 3rd order derivatives of the prime function (w.r.t. zeta) with spectral accuracy?

If not, we should get this in asap.

For points z on a boundary C_j, the equality

G_j(z,a) = G_0(z,a) - v_j(z) + v_j(a) + angle(a/(a - d_j))/(2*pi)

does not work numerically. If one compares the values of Green's function from the Schwarz problem with what comes out of the product formula, it will be noticed the values differ by a rotational constant (there is a real value c such that the rotational constant is exp(1i*c)).

The class vjFirstKind needs dependent parameter taujj to compute tau_jj = vj(theta_j(z)) - vj(z).

The skpDomain.isin function only checks that points are in the open, punctured unit domain.

Should there be multiple functions? One for all of the fundamental domain as well? One to include the closure?

See title.