summaryrefslogtreecommitdiff
path: root/Build/source/utils/asymptote/base/math.asy
diff options
context:
space:
mode:
Diffstat (limited to 'Build/source/utils/asymptote/base/math.asy')
-rw-r--r--Build/source/utils/asymptote/base/math.asy50
1 files changed, 50 insertions, 0 deletions
diff --git a/Build/source/utils/asymptote/base/math.asy b/Build/source/utils/asymptote/base/math.asy
index 85d41d86d07..eef703cb7ee 100644
--- a/Build/source/utils/asymptote/base/math.asy
+++ b/Build/source/utils/asymptote/base/math.asy
@@ -38,6 +38,15 @@ bool polygon(path p)
return cyclic(p) && piecewisestraight(p);
}
+// Return the intersection time of the point on the line through p and q
+// that is closest to z.
+real intersect(pair p, pair q, pair z)
+{
+ pair u=q-p;
+ real denom=dot(u,u);
+ return denom == 0 ? infinity : dot(z-p,u)/denom;
+}
+
// Return the intersection time of the extension of the line segment PQ
// with the plane perpendicular to n and passing through Z.
real intersect(triple P, triple Q, triple n, triple Z)
@@ -185,6 +194,30 @@ real[][] operator /(real[][] a, real b)
return a*(1/b);
}
+private string incommensurate=
+ "Multiplication of incommensurate matrices is undefined";
+
+pair[][] operator * (pair[][] a, pair[][] b)
+{
+ int n=a.length;
+ int nb=b.length;
+ int nb0=b[0].length;
+ pair[][] m=new pair[n][nb0];
+ for(int i=0; i < n; ++i) {
+ pair[] ai=a[i];
+ pair[] mi=m[i];
+ if(ai.length != nb)
+ abort(incommensurate);
+ for(int j=0; j < nb0; ++j) {
+ pair sum;
+ for(int k=0; k < nb; ++k)
+ sum += ai[k]*b[k][j];
+ mi[j]=sum;
+ }
+ }
+ return m;
+}
+
bool square(real[][] m)
{
int n=m.length;
@@ -371,3 +404,20 @@ pair[] quarticroots(real a, real b, real c, real d, real e)
return roots;
}
+
+pair[][] fft(pair[][] a, int sign=1)
+{
+ pair[][] A=new pair[a.length][];
+ int k=0;
+ for(pair[] v : a) {
+ A[k]=fft(v,sign);
+ ++k;
+ }
+ a=transpose(A);
+ k=0;
+ for(pair[] v : a) {
+ A[k]=fft(v,sign);
+ ++k;
+ }
+ return transpose(A);
+}