Boris Haase's Homepage

Results of Research with Music

External

# Scientific Computing

The following section presupposes Set Theory and Nonstandard Analysis. Let $$\varepsilon \in [\tilde{\nu}, 1 – \tilde{\nu}]$$.

Definition: Let $$f_n^*(z) = f(\eta_nz)$$ sisters of the TS $$f(z) \in \mathcal{O}(D)$$ centred on 0 on the domain $$D \subseteq {}^{\omega}\mathbb{C}$$ where $$m, n \in {}^{\omega}\mathbb{N}^{*}$$ and $$\eta_n^m := \underline{1}^{2^{\lceil m/n \rceil}}$$. Then let $$\delta_n^*f = \check{f} – \check{f}_n^*$$ the halved sister distances of $$f$$. Let $$\mu_n^m := m!n!/(m + n)!$$ and $$\widetilde{(-n)!} := 0$$. For $$j, k \in {}^{\omega}\mathbb{N}, f \in \mathcal{C}_{\pi}^{j+2}$$ and Fourier coefficients $$c_k := \tilde{\hat{\pi}}{\uparrow}_{-\pi}^{\pi}{f(t)\tilde{\epsilon}^{\underline{k}t}{{\downarrow}t}}$$1cf. Walter, Wolfgang: Analysis 2; 5., erw. Aufl.; 2002; Springer; Berlin, p. 358 ff., the series $${\LARGE{\textbf{+}}}_{k=-\omega}^{\omega}{c_k{\epsilon}^{\underline{k}t}}$$ is called Fourier series. On the level of TS, $$\mu$$ and $$\eta$$ form a $$\mu$$-$$\eta$$-calculus.$$\triangle$$

Remark: The latter allows representing integrals and derivatives easily and finitely closed2cf. Remmert, Reinhold: Funktionentheorie 1; 3., verb. Aufl.; 1992; Springer; Berlin, p. 165 f.. A conversion to other period lengths than $$\hat{\pi }$$ is simply possible3see loc. cit., p. 364.

Example: $$\text{Re} \; c \in [\tilde{\nu}, 1 + \tilde{\nu}[, c \in {}^{\omega}\mathbb{C}$$ and $$\varepsilon := \tilde{2}^j, j \in {}^{\omega}\mathbb{N}^{*}$$ imply $$\zeta (n+c)=a_n{\LARGE{\textbf{+}}}_{k=1}^{n}{\delta_n^* v_c(\varepsilon u^k)}$$ for $$z \in {}^{1-\tilde{\nu}}\dot{\mathbb{C}} \subset \mathbb{D}$$4cf. Remmert, Reinhold: Funktionentheorie 2; 1. unveränd. Nachdruck der 1. Aufl.; 1992; Springer; Berlin, p. 42 and $$v_c(z):={\LARGE{\textbf{+}}}_{m=1}^{\omega }{\zeta (m+c){{z}^{m}}}=z{\LARGE{\textbf{+}}}_{m=1}^{\omega }{{{{\widetilde{m}}}^{c}}\widetilde{z-m}}.\square$$

Example: A Richardson extrapolation determines the digamma function $$\psi$$ and$\zeta(\grave{n}) =\widetilde{\varepsilon^n\hat{n}}{\LARGE{\textbf{+}}}_{k=1}^n{(\psi(\varepsilon u^k \underline{1}^{2\tilde{n}}) – \psi(\varepsilon u^k))} + \mathcal{O}(\varepsilon^{\hat{n}})$ for $$n = 2$$ and $$\varepsilon = \tilde{2}^{13}$$ as $$\zeta(3) = 2^{24}{\LARGE{\textbf{+}}}_{k=1}^2{(\psi(\underline{\varepsilon} u^k) – \psi(\varepsilon u^k))} + \mathcal{O}(10^{-16})$$.

Representation theorem for integrals: The TS (see below) $$f(z) \in \mathcal{O}(D)$$ centred on 0 for $$D \subseteq {}^{\omega}\mathbb{C}$$ gives for $$\grave{m}, n \in {}^{\omega}\mathbb{N}^*$$${\uparrow}_0^z…{\uparrow}_0^{\zeta_2}{f(\zeta_1){\downarrow}\zeta_1\;…\;{\downarrow}\zeta_n} = \widetilde{n!} f(z\mu_n) z^n.\square$Example: For the TS $$f(x), g(x) \in {}^{\omega}\mathbb{R}$$, it holds that $${\uparrow}_0^x{f(v){\downarrow}v}{\uparrow}_0^x{\uparrow}_0^{y}{g(v){\downarrow}v{\downarrow}y} = f(x\mu_1)\check{g}(x\mu_2)x^3.$$

Representation theorem for derivatives: For $${}^{\widetilde{\nu}}\dot{\mathbb{C}} \subset \mathbb{D} \subseteq {}^{\omega}\mathbb{C}, n$$-th unit roots, TS$f(z):=f(0) + {\LARGE{\textbf{+}}}_{m=1}^{\omega }{\widetilde{m!}\,{{f}^{(m)}}(0){z^m}},$$$\varepsilon := \tilde{2}^j\tilde{r}, j \in {}^{\omega}\mathbb{Z}, n = {\epsilon}^{\sigma} \in {}^{\omega}\mathbb{N}^{*}, u :={\epsilon}^{\tilde{n} \hat{\underline{\pi}}}$$ and $$f$$’s radius of convergence $$r \in {}^{\nu}{\mathbb{R}}_{>0}$$ imply${{f}^{(n)}}(0)=2^{jn}\acute{n}!{\LARGE{\textbf{+}}}_{k=1}^{n}{\delta_n^* f(\tilde{2}^j u^k)}.\square$Universal multistep theorem: For $$n \in {}^{\nu}\mathbb{N}_{\le p}, k, m, p \in {}^{\nu}\mathbb{N}^{*}, {\downarrow}\overset{\rightharpoonup}{x} \in\, ]0, 1[, x \in [a, b] \subseteq {}^{\omega}\mathbb{R}, y : [a, b] \rightarrow {}^{\omega}\mathbb{R}^q,$$ $$f : [a, b]\times{}^{\omega}\mathbb{R}^{q \times n} \rightarrow {}^{\omega}\mathbb{R}^q, g_k(\overset{\rightharpoonup}{x}) := g_{\acute{k}}(x)$$, and $$g_0(a) = f(\overset{\leftharpoonup}{a}, y_0, … , y_{\acute{n}})$$, the TS of the initial value problem $$y^\prime(x) = f(x, y((\rightharpoonup)^0 x), … , y((\rightharpoonup)^{\acute{n}} x))$$ of order $$n$$ implies$y(\overset{\rightharpoonup}{x}) = y(x) + {\downarrow}\overset{\rightharpoonup}{x}{\LARGE{\textbf{\pm}}}_{k=1}^{p}{\left (g_{p-k}(\overset{\rightharpoonup}{x}){\LARGE{\textbf{+}}}_{m=k}^{p}{\widetilde{m!}\tbinom{\acute{m}}{\acute{k}}}\right )} + \mathcal{O}(({\downarrow}\overset{\rightharpoonup}{x})^{\grave{p}}).\square$Theorem for (anti-) derivatives of TS: For $$m \in {}^{\omega}\mathbb{Z}$$, $$q = \tilde{\varepsilon}(z-a)$$, $$a \in \mathbb{D}$$ and $$j, k \in \mathbb{N}_{<n}$$, modular arithmetic5cf. Knuth, Donald Ervin: The Art of Computer Programming Volume 2; 3rd Ed.; 1997; Addison Wesley; Reading, p. 302 – 311 and $$n$$-th unit roots result in the corresponding DFT form:${\updownarrow}^mf_n(z) := \tilde{n}(q^j)^T(\delta_{jk}\widetilde{\varepsilon q}^jj!/(j-m)!)({\tilde{u}}^{jk})(f({\varepsilon u}^k+a))+\mathcal{O}(\varepsilon^n).\square$Reversion theorem of TS: For $$y \in f(\mathbb{D}), y(a) = b$$ and $$y^{\prime}(a) \ne 0$$, Bürmann’s theorem yields6see Whittaker, Edmund T.; Watson, George N.: A Course of Modern Analysis; 5th Ed.; 2021; Cambridge University Press; Cambridge, p. 129 ff. :$f^{-1}(y) = a + \tilde{n} {\LARGE{\textbf{+}}}_{m=1}^n{\widetilde{m}{\tilde{\varepsilon}}^{\acute{m}}(y – b)^m({\tilde{u}}^{\acute{m}k})^T(f(\varepsilon u^k + a)^{-m})}+\mathcal{O}(\varepsilon^n).\square$Remark: Put $$y = 0$$ for a sufficiently small $$|a – z_0|$$ such that $$z_0 \in {}^{\omega}\mathbb{C}$$ is a zero of $$f(z) \in {}^{\omega}\mathbb{C}$$.

Theorem for derivatives of Fourier series: For $$f \in \mathcal{C}_{\pi}^{m+2}$$7cf. Walter, loc. cit., $$m \in {}^{\omega}\mathbb{N}, j \in [-\acute{n}, n] \cap \mathbb{Z}, t \in [-\pi, \pi]$$ and $$k \in \mathbb{N}_{<\hat{n}}$$, the following DFT form is implied allowing integration of individual terms:${\downarrow}^mf_n(t):=\tilde{n}(u^{\tilde{\pi}jnt})^T(\delta_{(j+\acute{n})k}\underline{j}^m)(\tilde{u}^{jk})(\check{f}(\tilde{n}\pi k – \pi))+\mathcal{O}(\tilde{n}).\square$Conclusion: Supporting points $$kr := \tilde{n}\pi k$$ of the smooth $$f(kr)$$ yield for $$j$$ like $$k$$ interpolating in $$\mathcal{O}(\sigma n)$$:${\downarrow}^mf_n(t):=\tilde{n}(u^{\tilde{r}jt})^T(\delta_{jk}\underline{j}^m)(\tilde{u}^{jk})(\check{f}(kr)).\square$Theorem (DFT fixed-point method): Every $$z \in {}^{\omega}\mathbb{C}$$ of an arbitrary $$m$$-polynomial $$p(z) = 0$$ with $$m \in [2, n] \cap \mathbb{N}$$ for $$n := 2^r, r \in {}^{\nu}\mathbb{N}^*$$ and coefficients from $${}^{\nu}\mathbb{C}$$ can be (initiating) determined also in $$\mathcal{O}(\sigma n).$$

Proof and algorithm: Let $$U = (\tilde{u}^{jk})$$ for $$j, k \in \mathbb{N}_{<n}, u :={\epsilon}^{\tilde{n} \hat{\underline{\pi}}}, q := \hat{z}$$ and $$s_k := p(\tilde{2}u^k)$$. A simple transform achieves $$|q| < \tilde{2}$$ for all zeros $$\zeta$$ of $$p(z)$$ and $$p(0) = 1$$. It follows from $$p(z) = \tilde{n}(q^j)^TUs = \tilde{n}\mu^Ts = 0$$ the simplified iteration $$\mu^* = U_{1}^{-T}\mu U((\delta_{jk}\tilde{u}^j)U^{-1}\mu-(U_{\acute{n}}^{-T}\mu+\beta s^T\mu, \beta s^T\mu, …, \beta s^T\mu)^T)$$. Here let $$\delta_{jk}$$ the Kronecker delta, starting point $$q := \tilde{2}$$ and $$\beta \in {}^{\nu}\mathbb{C}^*$$ such that each time $$||\mu^*-\mu||$$ is roughly halved and $$\mu^Ts = 0$$ holds. Finish by polynomial division where $$m > 2.\square$$

Remark: The identity instead of $$\delta_n^*$$ provides arbitrarily precise approximations for the $$f^{(n)}$$. The last theorems are equally valid for multidimensional TS and Laurent series. The error is $$\mathcal{O}(\varepsilon^n)$$ instead of $$\mathcal{O}(\varepsilon)$$ for analogously defined $$m$$-dimensional DFT forms of the TS with $$\binom{m+n}{n}$$ derivatives where the effort is comparable. A good choice for the crown points $$\varepsilon u^m + a$$ around the crown centre $$a$$ is $$\varepsilon = \tilde{2}$$ and $$n = 64$$.

Approximation theorem: The $$f^{(s)}(x) \in {}^{\omega}\mathbb{R}$$ for $$x \in {}^{\omega}\mathbb{R}$$ as above allow computing the interpolant

$$g(x) := {\LARGE{\textbf{+}}}_{r=0}^{\acute{m}}{\chi_{]x_r, x_{\grave{r}}[}(x)((x_{\grave{r}}-x)p_r(x)+(x-x_r)p_{\grave{r}}(x))/(x_{\grave{r}}-x_r)}+$$$${\LARGE{\textbf{+}}}_{r=0}^m{\chi_{\{x_r\}}(x)p_r(x)}$$

for $$m, n \in {}^{\nu}\mathbb{N}$$ and $$p_r(x) := {\LARGE{\textbf{+}}}_{s=0}^n{f^{(s)}(x_r){(x-x_r)}^s/s!}$$ in $$\mathcal{O}(\sigma mn)$$ where $$f^{(s)}(x_r) = g^{(s)}(x_r)$$ holds for every $$x_r \in {}^{\omega}\mathbb{R}$$. Replace in the complex case $${}^{\omega}\mathbb{R}$$ by $${}^{\omega}\mathbb{C}$$ and put $$x = \gamma(t) \in {}^{\omega}\mathbb{C}$$ for the path $$\gamma(t)$$ where $$t \in {}^{\omega}\mathbb{R}.\square$$

Theorem (differential equations): Let $$d(z) \in {}^{\omega}\mathbb{C}^{n+2}$$ vector of DFT forms for $$n \in {}^{\omega}\mathbb{N}$$. If the system $$y^{(\grave{n})} = \left(y^{(n)},y^{(\acute{n})}, …, y, 1\right)d(z) \in {}^{\omega}\mathbb{C}$$ where (not) $$\left(y^{(n)}(a),\ y^{(\acute{n})}(a), …, y(a)\right)^T\in {}^{\omega}\mathbb{C}^{\grave{n}}$$ and $$a, z \in {}^{\omega}\mathbb{C}$$ has the solution $$y_0$$, this can be determined in $$\mathcal{O}(n^3)$$ giving the error mentioned above.$$\square$$

Remark: If function values of crown points are evaluated and terminating TS are used, non-linear or partial (systems of) differential equations may be solved via DFT forms.

Example: For $$b = (f(\varepsilon u^0 + a), …, f(\varepsilon u^{\acute{n}} + a), 1)^T \in {}^{\omega}\mathbb{C}^{\grave{n}}$$ and an upper triangle matrix $$T(z) \in {}^{\omega}\mathbb{C}^{\grave{n}\times\grave{n}}$$, ordinary differential equations have only solutions if $$T(z)b = 0$$ holds (Cholesky decomposition!).

Definition: To increase precision, let according to the trapezoidal rule $${}^{T}{\uparrow}_{z\in A}{f(z){\downarrow}z:={\LARGE{\textbf{+}}}_{z\in A}{\tilde{2}(f(z)+f(\overset{\rightharpoonup}{z}))(\overset{\rightharpoonup}{z}-z)}}$$ and according to the midpoint rule $${}^{M}{\uparrow}_{z\in A}{f(z){\downarrow}z:={\LARGE{\textbf{+}}}_{z\in A}{f(\tilde{2}(z\,+\overset{\rightharpoonup}{z}}))(\overset{\rightharpoonup}{z}-z)}.\triangle$$

Series theorem for integrals: For $$c \in [0, a], f \in {}^{\nu}\mathcal{C}^{n}$$ and $$n, \hat{p}, r\,(= 2^p), t\,(=2^v), \hat{v} \in {}^{\nu}2\mathbb{N}^*$$, DFT forms of TS show

$${\uparrow}_0^a f(w){\downarrow}w = {\uparrow}_0^{\tilde{t}a}{\LARGE{\textbf{+}}}_{s=1}^t{f(x+\acute{s}\tilde{t}a){\downarrow}x} = {\uparrow}_0^1 g(y){\downarrow}y$$$$= {\LARGE{\textbf{+}}}_{\check{q}=1}^{\check{r}}{{\LARGE{\textbf{+}}}_{\check{m}=0}^{\check{n}-1}{\widetilde{\grave{m}!}\hat{g}^{(m)}\left (\tilde{r}\acute{q} \right){\tilde{r}}^{\grave{m}}}}+\mathcal{O}\left(\widetilde{\grave{n}!}g^{(n)}(\tilde{a}c){\tilde{r}}^n\right).\square$$

Remarks: DFT equivalents may replace any $$g^{(m)}$$ for $$g(y) := af(ay)$$ and $$a > 0$$. Transforming $$z := {\epsilon}^{\pm x}$$ makes here infinite bounds of integration finite. A finite integral decreases the remainder’s modulus sufficiently. The midpoint rule is here more advantageous than the trapezoidal rule.

Example: For $$g_{\mp}(\vartheta) := \check{\pi}(1 – {\varepsilon}^2 {\sin}^2(\check{\pi}\vartheta))^{\mp\tilde{2}}$$, complete elliptic integrals of the first and second kind are ${\uparrow}_0^{1} g_{\mp}(\vartheta) {\downarrow}\vartheta = g_{\mp}(\check{1}) + \widetilde{24} g_{\mp}^{(2)}(\check{1}) + \widetilde{1920} g_{\mp}^{(4)}(\check{1}) + \mathcal{O}\left(\widetilde{9!}g_{\mp}^{(6)}(\check{1})\right).$Second Euler-Maclaurin formula: For $$f(\check{q}) = g(\tilde{r}\acute{q})$$ and $$k = \grave{a} = \grave{r}$$, the preceding theorem yields in $$\mathcal{O}(\sigma n)$$${\LARGE{\textbf{+}}}_{\check{q}=1}^{\check{r}} f(\check{q}) = {\uparrow}_{\check{1}}^{\check{k}} f(x){\downarrow}x + {\LARGE{\textbf{+}}}_{\check{m}=1}^{\check{n}-1} H_m \left (f^{(\acute{m})}(\check{k}) – f^{(\acute{m})}(\check{1}) \right ) + \mathcal{O}\left (H_n \left (f^{(\acute{n})}(\check{k}) – f^{(\acute{n})}(\check{1}) \right )\right ).$Proof: For $$h(x) = x/$$sin $$x$$ and $$H_m := \underline{1}^m\widetilde{m!}h^{(m)}(0)=\widetilde{m!}B_m(2-2^m) \rightarrow 2\tilde{\underline{\pi}}\text{-}^m$$,${\LARGE{\textbf{+}}}_{\check{q}=1}^{\check{r}} g(\tilde{r}\acute{q}) = \check{r}{\uparrow}_0^1 g(y){\downarrow}y – {\LARGE{\textbf{+}}}_{\check{q}=1}^{\check{r}}{\LARGE{\textbf{+}}}_{\check{m}=1}^{\check{n}-1} {\widetilde{\grave{m}!}g^{(m)}\left (\tilde{r}\acute{q} \right)} + \mathcal{O}\left(\widetilde{\grave{n}!}g^{(n)}(\tilde{a}c){\tilde{r}}^n\right)$results in the claim by inserting and combining the factorials that are counted depending on partitions where $$B_m$$ are Bernoulli numbers and $$B_{\grave{m}} = 0, B_0=1$$ as well as $$B_1=-\tilde{2}.\square$$

First Euler-Maclaurin formula: Mathematical induction for $$n$$ shows also that8cf. Mollin, Richard A.: Advanced Number Theory with Applications; 1st Ed.; 2017; Routledge; Milton Park, p. 193 f.

$${\LARGE{\textbf{+}}}_{\check{q}=0}^{\check{r}}f(\check{q})=\uparrow_0^{\check{r}}{f\left(x\right){\downarrow}x}+\check{f}(\check{r})+\check{f}(0)+{\LARGE{\textbf{+}}}_{\check{m}=1}^{\check{n}-1}{\widetilde{m!}}B_m\left(f^{(\acute{m})}\left(\check{r}\right)-f^{(\acute{m})}\left(0\right)\right)+$$$$\mathcal{O}\left({\widetilde{n!}}B_n\left(f^{(\acute{n})}\left(\check{r}\right)-f^{(\acute{n})}\left(0\right)\right)\right).\square$$

code of the FFT form

© 2010-2024 by Boris Haase

top

References

↑1 cf. Walter, Wolfgang: Analysis 2; 5., erw. Aufl.; 2002; Springer; Berlin, p. 358 ff. cf. Remmert, Reinhold: Funktionentheorie 1; 3., verb. Aufl.; 1992; Springer; Berlin, p. 165 f. see loc. cit., p. 364 cf. Remmert, Reinhold: Funktionentheorie 2; 1. unveränd. Nachdruck der 1. Aufl.; 1992; Springer; Berlin, p. 42 cf. Knuth, Donald Ervin: The Art of Computer Programming Volume 2; 3rd Ed.; 1997; Addison Wesley; Reading, p. 302 – 311 see Whittaker, Edmund T.; Watson, George N.: A Course of Modern Analysis; 5th Ed.; 2021; Cambridge University Press; Cambridge, p. 129 ff. cf. Walter, loc. cit. cf. Mollin, Richard A.: Advanced Number Theory with Applications; 1st Ed.; 2017; Routledge; Milton Park, p. 193 f.